{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Frontiers: Copper Particle in Viscous Oil\n\nThis example corresponds to section 3.3 in\n`our publication <CitingOsaft>`.\nIn this example we study the acoustic radiation force (ARF) on a copper\nparticle suspended in a viscous oil. We compare the theory by Doinikov\n(rigid, 1994) and the theory by Settnes & Bruus (2012).\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "As always we start off by importing the necessary Python modules. For our\nexample we are going to need the osaft library, and the packages NumPy and\nMatplotlib.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy as np\nfrom matplotlib import pyplot as plt\nfrom matplotlib.patches import Circle\n\nimport osaft"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The next step is to define the properties for our example,\nthese include the material properties, the properties of the acoustic field\nand the radius. We always assume SI-units.\n\nThe wave type is set using the ``osaft.WaveType`` enum. Currently, there are\ntwo options:\n``osaft.WaveType.STANDING`` and ``osaft.WaveType.TRAVELLING`` for a plane\nstanding wave and a plane travelling wave, respectively.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# --------\n# Geometry\n# --------\n# Radius\nR_0 = 5e-6  # [m]\n\n# --------------------\n# Properties of Copper\n# --------------------\n# Speed of sound\nrho_cu = 8_930  # [m/s]\n# Density\nc_cu = 5_100  # [kg/m^3]\n\n# -------------------\n# Properties of Oil\n# -------------------\n# Speed of sound\nc_oil = 1_445  # [m/s]\n# Density\nrho_oil = 922.6  # [kg/m^3]\n# Viscosity\neta_oil = 0.03  # [Pa s]\nzeta_oil = 0  # [Pa s]\n\n# --------------------------------\n# Properties of the Acoustic Field\n# --------------------------------\n# Frequency\nf = 5e5  # [Hz]\n# Pressure\np_0 = 1e5  # [Pa]\n# Wave type\nwave_type = osaft.WaveType.STANDING\n# Position of the particle in the field\nposition = np.pi / 4  # [rad]\n\n# Once all properties are defined we can initialize the solution classes.\n# In this example, we use the classes ``osaft.doinikov1994rigid.ARF()``, and\n# ``osaft.settnes2012.ARF()``.\n\ndoinikov = osaft.doinikov1994rigid.ARF(\n    f=f,\n    R_0=R_0,\n    rho_s=rho_cu,\n    rho_f=rho_oil,\n    c_f=c_oil,\n    eta_f=eta_oil,\n    zeta_f=zeta_oil,\n    p_0=p_0,\n    wave_type=wave_type,\n    position=position,\n    long_wavelength=True,\n)\n\nsettnes = osaft.settnes2012.ARF(\n    f=f,\n    R_0=R_0,\n    rho_s=rho_cu,\n    c_s=c_cu,\n    rho_f=rho_oil,\n    c_f=c_oil,\n    eta_f=eta_oil,\n    p_0=p_0,\n    wave_type=wave_type,\n    position=position,\n)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Next, we want to compare the boundary layer thickness and the particle\nradius for the given parameter. Both quantities are properties of our\nsolution classes and can easily be evaluated, and we can compute the ratio.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(f\"{settnes.delta / settnes.R_0 = :.2f}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "With the model from Doinikov it is also possible to compute the scattering\nfield. To plot the scattering field we need to initialize a\n``osaft.FluidScatteringPlot()`` instance.\nThe class takes the model as an argument and the radius range ``r_max`` that\nis plotted. For more options see the documentation.\nThe method ``plot()`` will then generate the plot.\nAs always with Matplotlib, we need to call ``plt.show()`` to display the\nplot.\n``plot()`` returns two Matplotlib objects, a ``Figure`` and an ``Axes``\nobject.\nThese can be used to further manipulate the plot and to save it.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Scattering plot small viscosity\nscattering_plot = osaft.FluidScatteringPlot(doinikov, r_max=5 * doinikov.R_0)\nfig, ax = scattering_plot.plot_velocity(inst=False, incident=False)\n\n# Adding a circle to illustrate the boundary layer thickness\ncircle = Circle(\n    (0, 0),\n    doinikov.R_0 + doinikov.delta,\n    fill=False,\n    edgecolor=\"white\",\n    linestyle=\"--\",\n)\nax.add_artist(circle)\nplt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Finally, we want to compare the ARF in the different models. To plot the\nARF we need to initialize an ``osaft.ARFPlot()`` instance. With the method\n``add_solutions()`` we can add our models to the plotter.\nWith ``set_abscissa()`` we define the variable that we want to\nplot the ARF against. Here we select the fluid viscosity `eta_f`.\nFinally, ``osaft.ARFPlot.plot_solutions()`` will generate the plot. Again,\nwe call ``plt.show()`` to display the plot.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# sphinx_gallery_thumbnail_number = -1\n\n# Initializing plotting class\narf_plot = osaft.ARFPlot()\n\n# Add solutions to be plotted\narf_plot.add_solutions(settnes, doinikov)\n\n# Define independent plotting variable (in this case the fluid viscosity)\narf_plot.set_abscissa(np.linspace(1e-5, 0.1, 300), \"eta_f\")\nfig, ax = arf_plot.plot_solutions()\n\n# Setting the axis labels\nax.set_xlabel(r\"$\\eta \\, \\mathrm{[Pa s]}$\")\n\nplt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "<div class=\"alert alert-info\"><h4>Note</h4><p>It is only possible to plot the ARF  against\n      properties that wrap an underlying ``PassiveVariable``, i.e. an input\n      parameter of the model.\n      You can get a list of all input variables using the method\n      ``input_variables()``.</p></div>\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(f\"{doinikov.input_variables() = }\")\nprint(f\"{settnes.input_variables() = }\")"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.9.15"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}