utils¶

Utility modules: quick-start presets, batched parameter sweeps, and Matplotlib styles for figures.

Presets¶

Quick-start factory functions for common bubble configurations.

from jbubble.utils.presets import free_bubble, lipid_bubble, thick_shell_bubble

BubblePreset ¶

Bases: tuple

Assembled (eom, pulse) pair that the preset factory functions return.

You can unpack it directly:

eom, pulse = free_bubble()

Or you can access the fields by name:

preset = free_bubble()
result = run_simulation(preset.eom, preset.pulse)

Attributes:

Name	Type	Description
`eom`	`EquationOfMotion`	Assembled equation of motion.
`pulse`	`Pulse`	Driving pulse.

free_bubble ¶

free_bubble(
    R0: float = 2e-06,
    freq: float = 1000000.0,
    pressure: float = 100000.0,
    cycle_num: int = 5,
    *,
    gamma: float = 1.4,
    sigma: float = 0.072,
    mu: float = 0.001,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset

Uncoated gas bubble in a Newtonian liquid (Keller-Miksis).

The simplest physically meaningful preset: a polytropic gas bubble with Laplace surface tension only and no shell coating. Suitable for free cavitation bubbles, and as a baseline before you add a shell.

Physics: KellerMiksis + PolytropicGas + NoShell + NewtonianMedium, driven by a sine ToneBurst.

Parameters:

Name	Type	Description	Default
`R0`	`float`	Equilibrium radius [m]. Default: `2e-6` (2 µm).	`2e-06`
`freq`	`float`	Driving frequency [Hz]. Default: `1e6` (1 MHz).	`1000000.0`
`pressure`	`float`	Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).	`100000.0`
`cycle_num`	`int`	Number of tone-burst cycles. Default: `5`.	`5`
`gamma`	`float`	Polytropic exponent. Default: `1.4`, the adiabatic exponent of a diatomic gas such as air.	`1.4`
`sigma`	`float`	Surface tension [N/m]. Default: `0.072` (surface tension of water at 20 °C: 0.0727 N/m, NIST).	`0.072`
`mu`	`float`	Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C: 1.0016e-3 Pa s, NIST).	`0.001`
`P_amb`	`float`	Ambient pressure [Pa]. Default: `101325.0` (one standard atmosphere).	`101325.0`
`rho_L`	`float`	Liquid density [kg/m³]. Default: `998.0` (water at 20 °C: 998.2 kg/m³, NIST).	`998.0`
`c_L`	`float`	Speed of sound in the liquid [m/s]. Default: `1500.0`, the round value that ultrasound modelling commonly uses for water (NIST gives 1482 m/s at 20 °C and 1497 m/s at 25 °C).	`1500.0`

Returns:

Type	Description
`BubblePreset`	`(eom, pulse)` pair, ready for `run_simulation`.

References

NIST Chemistry WebBook, NIST Standard Reference Database 69, Thermophysical Properties of Fluid Systems: water at 0.101325 MPa. doi:10.18434/T4D303

Source code in jbubble/utils/presets.py

def free_bubble(
    R0: float = 2e-6,
    freq: float = 1e6,
    pressure: float = 100e3,
    cycle_num: int = 5,
    *,
    gamma: float = 1.4,
    sigma: float = 0.072,
    mu: float = 1e-3,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset:
    """Uncoated gas bubble in a Newtonian liquid (Keller-Miksis).

    The simplest physically meaningful preset: a polytropic gas bubble with
    Laplace surface tension only and no shell coating. Suitable for free
    cavitation bubbles, and as a baseline before you add a shell.

    Physics: [`KellerMiksis`][jbubble.bubble.eom.KellerMiksis] +
    [`PolytropicGas`][jbubble.bubble.gas.PolytropicGas] +
    [`NoShell`][jbubble.bubble.shell.NoShell] +
    [`NewtonianMedium`][jbubble.bubble.medium.NewtonianMedium], driven by
    a sine [`ToneBurst`][jbubble.pulse.tone_burst.ToneBurst].

    Parameters
    ----------
    R0 : float
        Equilibrium radius [m]. Default: `2e-6` (2 µm).
    freq : float
        Driving frequency [Hz]. Default: `1e6` (1 MHz).
    pressure : float
        Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).
    cycle_num : int
        Number of tone-burst cycles. Default: `5`.
    gamma : float
        Polytropic exponent. Default: `1.4`, the adiabatic exponent of a
        diatomic gas such as air.
    sigma : float
        Surface tension [N/m]. Default: `0.072` (surface tension of water
        at 20 °C: 0.0727 N/m, NIST).
    mu : float
        Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C:
        1.0016e-3 Pa s, NIST).
    P_amb : float
        Ambient pressure [Pa]. Default: `101325.0` (one standard
        atmosphere).
    rho_L : float
        Liquid density [kg/m³]. Default: `998.0` (water at 20 °C:
        998.2 kg/m³, NIST).
    c_L : float
        Speed of sound in the liquid [m/s]. Default: `1500.0`, the round
        value that ultrasound modelling commonly uses for water (NIST gives
        1482 m/s at 20 °C and 1497 m/s at 25 °C).

    Returns
    -------
    BubblePreset
        `(eom, pulse)` pair, ready for
        [`run_simulation`][jbubble.simulation.run_simulation].

    References
    ----------
    NIST Chemistry WebBook, NIST Standard Reference Database 69,
    Thermophysical Properties of Fluid Systems: water at 0.101325 MPa.
    [doi:10.18434/T4D303](https://doi.org/10.18434/T4D303)
    """
    eom = KellerMiksis(
        gas=PolytropicGas(gamma=gamma),
        shell=NoShell(sigma=sigma),
        medium=NewtonianMedium(mu=mu, rho_L=rho_L, c_L=c_L),
        R0=R0,
        P_amb=P_amb,
    )
    pulse = ToneBurst(freq=freq, pressure=pressure, shape=Sine(), cycle_num=cycle_num)
    return BubblePreset(eom=eom, pulse=pulse)

lipid_bubble ¶

lipid_bubble(
    R0: float = 2e-06,
    freq: float = 1000000.0,
    pressure: float = 100000.0,
    cycle_num: int = 5,
    *,
    kappa_s: float = 7.5e-09,
    chi: float = 0.5,
    sigma_rupture: float = 0.072,
    R_buckle_ratio: float = 0.98058,
    smoothing: float = 0.01,
    gamma: float = 1.095,
    mu: float = 0.001,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset

Lipid-shelled ultrasound contrast agent: the Marmottant (2005) model.

Models an SF6 bubble coated with a thin lipid monolayer that buckles under compression, stretches elastically, and ruptures under expansion, with surface-dilatational viscosity. The surface tension follows SmoothMarmottantSurfaceTension, which matches the piecewise Marmottant law to within smoothing * ln(2) * sigma_rupture and keeps gradients smooth for fitting.

The shell and gas defaults are the SonoVue parameters of Gümmer et al. (2021), who take them from earlier characterisations of SonoVue. The preset solves the Keller-Miksis equation with a polytropic gas, not the equations of motion and the hard-core gas that Gümmer et al. use.

Physics: KellerMiksis + PolytropicGas + LipidShell with SmoothMarmottantSurfaceTension + NewtonianMedium, driven by a sine ToneBurst.

Parameters:

Name	Type	Description	Default
`R0`	`float`	Equilibrium radius [m]. Default: `2e-6` (2 µm, the top of the 1-2 µm range that Gümmer et al. study).	`2e-06`
`freq`	`float`	Driving frequency [Hz]. Default: `1e6` (1 MHz).	`1000000.0`
`pressure`	`float`	Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).	`100000.0`
`cycle_num`	`int`	Number of tone-burst cycles. Default: `5`.	`5`
`kappa_s`	`float`	Shell surface-dilatational viscosity [N s/m]. Default: `7.5e-9` (Gümmer et al., who consider 5e-9 to 1e-8).	`7.5e-09`
`chi`	`float`	Shell elasticity [N/m]. Default: `0.5` (Gümmer et al.).	`0.5`
`sigma_rupture`	`float`	Surface tension after rupture [N/m]. Default: `0.072`, the surface tension of clean water (Gümmer et al.).	`0.072`
`R_buckle_ratio`	`float`	Buckling radius as a fraction of `R0`. Default: `0.98058`, which gives the surface tension at `R0` of Gümmer et al., \(\sigma_0 = 0.020\) N/m, through their Eq. 11: \((1 + \sigma_0/\chi)^{-1/2} = 0.98058\) at \(\chi = 0.5\) N/m. If you change `chi`, \(\sigma_0\) changes with it.	`0.98058`
`smoothing`	`float`	Corner width of the smoothed Marmottant law as a fraction of `sigma_rupture`. Default: `0.01`; see `SmoothMarmottantSurfaceTension`.	`0.01`
`gamma`	`float`	Polytropic exponent. Default: `1.095`, the value Gümmer et al. use for the SF6 in SonoVue.	`1.095`
`mu`	`float`	Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C: 1.0016e-3 Pa s, NIST).	`0.001`
`P_amb`	`float`	Ambient pressure [Pa]. Default: `101325.0` (one standard atmosphere).	`101325.0`
`rho_L`	`float`	Liquid density [kg/m³]. Default: `998.0` (water at 20 °C: 998.2 kg/m³, NIST).	`998.0`
`c_L`	`float`	Speed of sound in the liquid [m/s]. Default: `1500.0`, the round value that ultrasound modelling commonly uses for water (NIST gives 1482 m/s at 20 °C and 1497 m/s at 25 °C).	`1500.0`

Returns:

Type	Description
`BubblePreset`	`(eom, pulse)` pair, ready for `run_simulation`.

References

Gümmer, J., Schenke, S., & Denner, F. (2021). Modelling lipid-coated microbubbles in focused ultrasound applications at subresonance frequencies. Ultrasound Med. Biol. 47(10), 2958-2979. doi:10.1016/j.ultrasmedbio.2021.06.012

Marmottant, P., van der Meer, S., Emmer, M., Versluis, M., de Jong, N., Hilgenfeldt, S., & Lohse, D. (2005). A model for large amplitude oscillations of coated bubbles accounting for buckling and rupture. J. Acoust. Soc. Am. 118(6), 3499-3505. doi:10.1121/1.2109427

NIST Chemistry WebBook, NIST Standard Reference Database 69, Thermophysical Properties of Fluid Systems: water at 0.101325 MPa. doi:10.18434/T4D303

Source code in jbubble/utils/presets.py

def lipid_bubble(
    R0: float = 2e-6,
    freq: float = 1e6,
    pressure: float = 100e3,
    cycle_num: int = 5,
    *,
    kappa_s: float = 7.5e-9,
    chi: float = 0.5,
    sigma_rupture: float = 0.072,
    R_buckle_ratio: float = 0.98058,
    smoothing: float = 0.01,
    gamma: float = 1.095,
    mu: float = 1e-3,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset:
    r"""Lipid-shelled ultrasound contrast agent: the Marmottant (2005) model.

    Models an SF6 bubble coated with a thin lipid monolayer that buckles
    under compression, stretches elastically, and ruptures under
    expansion, with surface-dilatational viscosity. The surface tension
    follows
    [`SmoothMarmottantSurfaceTension`][jbubble.bubble.shell.SmoothMarmottantSurfaceTension],
    which matches the piecewise Marmottant law to within
    `smoothing * ln(2) * sigma_rupture` and keeps gradients smooth for
    fitting.

    The shell and gas defaults are the SonoVue parameters of Gümmer et al.
    (2021), who take them from earlier characterisations of SonoVue. The
    preset solves the Keller-Miksis equation with a polytropic gas, not
    the equations of motion and the hard-core gas that Gümmer et al. use.

    Physics: [`KellerMiksis`][jbubble.bubble.eom.KellerMiksis] +
    [`PolytropicGas`][jbubble.bubble.gas.PolytropicGas] +
    [`LipidShell`][jbubble.bubble.shell.LipidShell] with
    [`SmoothMarmottantSurfaceTension`][jbubble.bubble.shell.SmoothMarmottantSurfaceTension] +
    [`NewtonianMedium`][jbubble.bubble.medium.NewtonianMedium], driven by
    a sine [`ToneBurst`][jbubble.pulse.tone_burst.ToneBurst].

    Parameters
    ----------
    R0 : float
        Equilibrium radius [m]. Default: `2e-6` (2 µm, the top of the
        1-2 µm range that Gümmer et al. study).
    freq : float
        Driving frequency [Hz]. Default: `1e6` (1 MHz).
    pressure : float
        Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).
    cycle_num : int
        Number of tone-burst cycles. Default: `5`.
    kappa_s : float
        Shell surface-dilatational viscosity [N s/m]. Default: `7.5e-9`
        (Gümmer et al., who consider 5e-9 to 1e-8).
    chi : float
        Shell elasticity [N/m]. Default: `0.5` (Gümmer et al.).
    sigma_rupture : float
        Surface tension after rupture [N/m]. Default: `0.072`, the surface
        tension of clean water (Gümmer et al.).
    R_buckle_ratio : float
        Buckling radius as a fraction of `R0`. Default: `0.98058`, which
        gives the surface tension at `R0` of Gümmer et al.,
        $\sigma_0 = 0.020$ N/m, through their Eq. 11:
        $(1 + \sigma_0/\chi)^{-1/2} = 0.98058$ at $\chi = 0.5$ N/m. If you
        change `chi`, $\sigma_0$ changes with it.
    smoothing : float
        Corner width of the smoothed Marmottant law as a fraction of
        `sigma_rupture`. Default: `0.01`; see
        [`SmoothMarmottantSurfaceTension`][jbubble.bubble.shell.SmoothMarmottantSurfaceTension].
    gamma : float
        Polytropic exponent. Default: `1.095`, the value Gümmer et al. use
        for the SF6 in SonoVue.
    mu : float
        Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C:
        1.0016e-3 Pa s, NIST).
    P_amb : float
        Ambient pressure [Pa]. Default: `101325.0` (one standard
        atmosphere).
    rho_L : float
        Liquid density [kg/m³]. Default: `998.0` (water at 20 °C:
        998.2 kg/m³, NIST).
    c_L : float
        Speed of sound in the liquid [m/s]. Default: `1500.0`, the round
        value that ultrasound modelling commonly uses for water (NIST gives
        1482 m/s at 20 °C and 1497 m/s at 25 °C).

    Returns
    -------
    BubblePreset
        `(eom, pulse)` pair, ready for
        [`run_simulation`][jbubble.simulation.run_simulation].

    References
    ----------
    Gümmer, J., Schenke, S., & Denner, F. (2021). Modelling lipid-coated
    microbubbles in focused ultrasound applications at subresonance
    frequencies. *Ultrasound Med. Biol.* 47(10), 2958-2979.
    [doi:10.1016/j.ultrasmedbio.2021.06.012](https://doi.org/10.1016/j.ultrasmedbio.2021.06.012)

    Marmottant, P., van der Meer, S., Emmer, M., Versluis, M., de Jong, N.,
    Hilgenfeldt, S., & Lohse, D. (2005). A model for large amplitude
    oscillations of coated bubbles accounting for buckling and rupture.
    *J. Acoust. Soc. Am.* 118(6), 3499-3505.
    [doi:10.1121/1.2109427](https://doi.org/10.1121/1.2109427)

    NIST Chemistry WebBook, NIST Standard Reference Database 69,
    Thermophysical Properties of Fluid Systems: water at 0.101325 MPa.
    [doi:10.18434/T4D303](https://doi.org/10.18434/T4D303)
    """
    shell = LipidShell(
        sigma=SmoothMarmottantSurfaceTension(
            R_buckle_ratio=R_buckle_ratio,
            chi=chi,
            sigma_rupture=sigma_rupture,
            smoothing=smoothing,
        ),
        kappa_s=kappa_s,
    )
    eom = KellerMiksis(
        gas=PolytropicGas(gamma=gamma),
        shell=shell,
        medium=NewtonianMedium(mu=mu, rho_L=rho_L, c_L=c_L),
        R0=R0,
        P_amb=P_amb,
    )
    pulse = ToneBurst(freq=freq, pressure=pressure, shape=Sine(), cycle_num=cycle_num)
    return BubblePreset(eom=eom, pulse=pulse)

thick_shell_bubble ¶

thick_shell_bubble(
    R0: float = 2e-06,
    freq: float = 1000000.0,
    pressure: float = 100000.0,
    cycle_num: int = 5,
    *,
    d_s: float | None = None,
    G_s: float = 11700000.0,
    mu_s: float = 0.45,
    sigma: float = 0.04,
    gamma: float = 1.4,
    mu: float = 0.001,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset

Polymer-shelled bubble: the Church (1995) thick-shell model.

Models a gas bubble inside a viscoelastic solid shell of finite thickness, with ThickShell. The defaults describe the air-filled, polymer-shelled agent that Hoff et al. (2000) characterised, whose stiff shell dominates the bubble's response and damps it heavily. With the defaults, the undamped natural frequency is about 6.9 MHz, against about 2.0 MHz for free_bubble, and the quality factor is about 0.6. Damping this heavy leaves the linear response without a resonance peak. The bubble therefore oscillates far less than lipid_bubble at the same drive: its peak expansion at 100 kPa and 1 MHz is about 1%.

Protein-shelled agents, such as Albunex and Optison, have an albumin shell. For those, Qin and Ferrara (2010, Tables 1 and 2) use G_s=88.8e6, mu_s=1.77, and d_s=15e-9 at R0=1.5e-6.

Physics: KellerMiksis + PolytropicGas + ThickShell + NewtonianMedium, driven by a sine ToneBurst.

Parameters:

Name	Type	Description	Default
`R0`	`float`	Equilibrium (outer) radius [m]. Default: `2e-6` (2 µm).	`2e-06`
`freq`	`float`	Driving frequency [Hz]. Default: `1e6` (1 MHz).	`1000000.0`
`pressure`	`float`	Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).	`100000.0`
`cycle_num`	`int`	Number of tone-burst cycles. Default: `5`.	`5`
`d_s`	`float or None`	Shell thickness at `R0` [m]. Default: `None`, which uses `0.05 * R0`, the thickness-to-radius ratio of Hoff et al.	`None`
`G_s`	`float`	Shell shear modulus [Pa]. Default: `11.7e6` (11.7 MPa), within the 10.6-12.9 MPa that Hoff et al. estimate. Qin and Ferrara (2010, Table 1) use the same value for a polymer shell.	`11700000.0`
`mu_s`	`float`	Shell shear viscosity [Pa s]. Default: `0.45`, within the 0.39-0.49 Pa s that Hoff et al. estimate. Qin and Ferrara (2010, Table 1) use the same value for a polymer shell.	`0.45`
`sigma`	`float`	Surface tension [N/m]. Default: `0.04`, the gas-shell interfacial tension in Qin and Ferrara (2010, Table 1).	`0.04`
`gamma`	`float`	Polytropic exponent. Default: `1.4`, the adiabatic exponent of air, the gas in the Hoff et al. agent.	`1.4`
`mu`	`float`	Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C: 1.0016e-3 Pa s, NIST).	`0.001`
`P_amb`	`float`	Ambient pressure [Pa]. Default: `101325.0` (one standard atmosphere).	`101325.0`
`rho_L`	`float`	Liquid density [kg/m³]. Default: `998.0` (water at 20 °C: 998.2 kg/m³, NIST).	`998.0`
`c_L`	`float`	Speed of sound in the liquid [m/s]. Default: `1500.0`, the round value that ultrasound modelling commonly uses for water (NIST gives 1482 m/s at 20 °C and 1497 m/s at 25 °C).	`1500.0`

Returns:

Type	Description
`BubblePreset`	`(eom, pulse)` pair, ready for `run_simulation`.

References

Church, C. C. (1995). The effects of an elastic solid surface layer on the radial pulsations of gas bubbles. J. Acoust. Soc. Am. 97(3), 1510-1521. doi:10.1121/1.412091

Hoff, L., Sontum, P. C., & Hovem, J. M. (2000). Oscillations of polymeric microbubbles: Effect of the encapsulating shell. J. Acoust. Soc. Am. 107(4), 2272-2280. doi:10.1121/1.428557

Qin, S., & Ferrara, K. W. (2010). A model for the dynamics of ultrasound contrast agents in vivo. J. Acoust. Soc. Am. 128(3), 1511-1521. doi:10.1121/1.3409476

NIST Chemistry WebBook, NIST Standard Reference Database 69, Thermophysical Properties of Fluid Systems: water at 0.101325 MPa. doi:10.18434/T4D303

Source code in jbubble/utils/presets.py

def thick_shell_bubble(
    R0: float = 2e-6,
    freq: float = 1e6,
    pressure: float = 100e3,
    cycle_num: int = 5,
    *,
    d_s: float | None = None,
    G_s: float = 11.7e6,
    mu_s: float = 0.45,
    sigma: float = 0.04,
    gamma: float = 1.4,
    mu: float = 1e-3,
    P_amb: float = 101325.0,
    rho_L: float = 998.0,
    c_L: float = 1500.0,
) -> BubblePreset:
    """Polymer-shelled bubble: the Church (1995) thick-shell model.

    Models a gas bubble inside a viscoelastic solid shell of finite
    thickness, with [`ThickShell`][jbubble.bubble.shell.ThickShell]. The
    defaults describe the air-filled, polymer-shelled agent that Hoff et
    al. (2000) characterised, whose stiff shell dominates the bubble's
    response and damps it heavily. With the defaults, the undamped natural
    frequency is about 6.9 MHz, against about 2.0 MHz for
    [`free_bubble`][jbubble.utils.presets.free_bubble], and the quality
    factor is about 0.6. Damping this heavy leaves the linear response
    without a resonance peak. The bubble therefore oscillates far less than
    [`lipid_bubble`][jbubble.utils.presets.lipid_bubble] at the same drive:
    its peak expansion at 100 kPa and 1 MHz is about 1%.

    Protein-shelled agents, such as Albunex and Optison, have an albumin
    shell. For those, Qin and Ferrara (2010, Tables 1 and 2) use
    `G_s=88.8e6`, `mu_s=1.77`, and `d_s=15e-9` at `R0=1.5e-6`.

    Physics: [`KellerMiksis`][jbubble.bubble.eom.KellerMiksis] +
    [`PolytropicGas`][jbubble.bubble.gas.PolytropicGas] +
    [`ThickShell`][jbubble.bubble.shell.ThickShell] +
    [`NewtonianMedium`][jbubble.bubble.medium.NewtonianMedium], driven by
    a sine [`ToneBurst`][jbubble.pulse.tone_burst.ToneBurst].

    Parameters
    ----------
    R0 : float
        Equilibrium (outer) radius [m]. Default: `2e-6` (2 µm).
    freq : float
        Driving frequency [Hz]. Default: `1e6` (1 MHz).
    pressure : float
        Peak acoustic pressure amplitude [Pa]. Default: `100e3` (100 kPa).
    cycle_num : int
        Number of tone-burst cycles. Default: `5`.
    d_s : float or None
        Shell thickness at `R0` [m]. Default: `None`, which uses
        `0.05 * R0`, the thickness-to-radius ratio of Hoff et al.
    G_s : float
        Shell shear modulus [Pa]. Default: `11.7e6` (11.7 MPa), within the
        10.6-12.9 MPa that Hoff et al. estimate. Qin and Ferrara (2010,
        Table 1) use the same value for a polymer shell.
    mu_s : float
        Shell shear viscosity [Pa s]. Default: `0.45`, within the
        0.39-0.49 Pa s that Hoff et al. estimate. Qin and Ferrara (2010,
        Table 1) use the same value for a polymer shell.
    sigma : float
        Surface tension [N/m]. Default: `0.04`, the gas-shell interfacial
        tension in Qin and Ferrara (2010, Table 1).
    gamma : float
        Polytropic exponent. Default: `1.4`, the adiabatic exponent of air,
        the gas in the Hoff et al. agent.
    mu : float
        Liquid dynamic viscosity [Pa s]. Default: `1e-3` (water at 20 °C:
        1.0016e-3 Pa s, NIST).
    P_amb : float
        Ambient pressure [Pa]. Default: `101325.0` (one standard
        atmosphere).
    rho_L : float
        Liquid density [kg/m³]. Default: `998.0` (water at 20 °C:
        998.2 kg/m³, NIST).
    c_L : float
        Speed of sound in the liquid [m/s]. Default: `1500.0`, the round
        value that ultrasound modelling commonly uses for water (NIST gives
        1482 m/s at 20 °C and 1497 m/s at 25 °C).

    Returns
    -------
    BubblePreset
        `(eom, pulse)` pair, ready for
        [`run_simulation`][jbubble.simulation.run_simulation].

    References
    ----------
    Church, C. C. (1995). The effects of an elastic solid surface layer on
    the radial pulsations of gas bubbles. *J. Acoust. Soc. Am.* 97(3),
    1510-1521. [doi:10.1121/1.412091](https://doi.org/10.1121/1.412091)

    Hoff, L., Sontum, P. C., & Hovem, J. M. (2000). Oscillations of
    polymeric microbubbles: Effect of the encapsulating shell. *J. Acoust.
    Soc. Am.* 107(4), 2272-2280.
    [doi:10.1121/1.428557](https://doi.org/10.1121/1.428557)

    Qin, S., & Ferrara, K. W. (2010). A model for the dynamics of
    ultrasound contrast agents in vivo. *J. Acoust. Soc. Am.* 128(3),
    1511-1521. [doi:10.1121/1.3409476](https://doi.org/10.1121/1.3409476)

    NIST Chemistry WebBook, NIST Standard Reference Database 69,
    Thermophysical Properties of Fluid Systems: water at 0.101325 MPa.
    [doi:10.18434/T4D303](https://doi.org/10.18434/T4D303)
    """
    if d_s is None:
        d_s = 0.05 * R0
    eom = KellerMiksis(
        gas=PolytropicGas(gamma=gamma),
        shell=ThickShell(sigma=sigma, d_s=d_s, G_s=G_s, mu_s=mu_s),
        medium=NewtonianMedium(mu=mu, rho_L=rho_L, c_L=c_L),
        R0=R0,
        P_amb=P_amb,
    )
    pulse = ToneBurst(freq=freq, pressure=pressure, shape=Sine(), cycle_num=cycle_num)
    return BubblePreset(eom=eom, pulse=pulse)

GridSweep¶

Batched Cartesian-product parameter sweeps. GridSweep runs chunks of the grid in parallel on worker threads, one per CPU core by default (up to 31 per CPU device), and returns NumPy arrays in grid order. The guide Parameter sweeps explains how to choose workers and a batch size, and how to keep sweeps reproducible.

from jbubble.utils.gridsweep import GridSweep

GridSweep ¶

GridSweep(
    fn: Callable[..., Any],
    search_space: Mapping[str, Any],
    batch_size: int = 512,
    progress: bool = True,
    devices: int | Sequence[Device] | None = None,
    workers: int | None = None,
)

Batched, parallel sweep over a Cartesian product of named parameter axes.

GridSweep splits the grid into chunks of per_worker points and evaluates jax.vmap(fn) on each chunk. It traces fn once, compiles one executable per device ahead of time, and runs the chunks concurrently from workers threads. Results come back to the host as NumPy arrays, in grid order.

Parameters:

Name	Type	Description	Default
`fn`	`callable`	`fn(**params) -> PyTree`. `GridSweep` calls it with scalar keyword arguments, one per axis name, so it must be JAX-vmappable and jit-able.	required
`search_space`	`dict[str, Array]`	`{param_name: 1-D array of values}`. The sweep covers the Cartesian product of all axes, in `sorted()` order of the names (Unicode code-point order, so uppercase names sort before lowercase ones).	required
`batch_size`	`int`	Target number of grid points in each batch that `batches` yields. A batch holds one chunk per worker, and each worker evaluates `ceil(min(batch_size, total_points) / workers)` points per call. The batch size that `GridSweep` uses, in the `batch_size` attribute, is therefore at least `min(batch_size, total_points)` and less than that plus `workers`, and a grid with at most `batch_size` points comes in one batch. Default: `512`.	`512`
`progress`	`bool`	Whether to show a `tqdm` progress bar during iteration. Default: `True`.	`True`
`devices`	`int or sequence of jax.Device`	Devices to run on. `None` uses every local device on GPU or TPU. On CPU, it uses one local CPU device per 31 workers: the default CPU device for up to 31 workers, and more only if you set `JAX_NUM_CPU_DEVICES`. An `int` `n` uses the first `n` local devices, and a sequence uses exactly those devices. Default: `None`.	`None`
`workers`	`int`	Number of chunks evaluated at the same time, each from its own thread. `None` uses one worker per CPU core available to the process (respecting its CPU affinity), up to 31 per CPU device, on CPU, and one worker per device on GPU or TPU. `1` evaluates the chunks one after another on the calling thread. `workers` must be at least the number of devices. On CPU, `GridSweep` lowers a value above 31 per device to that limit; see the parallelism note. The default follows the core count, so it also sets `per_worker`; see the reproducibility note. Default: `None`.	`None`

Attributes:

Name	Type	Description
`devices`	`list[Device]`	Devices that the sweep runs on. A grid with fewer chunks than devices uses only as many devices as it has chunks.
`workers`	`int`	Number of chunks evaluated at the same time. It's at most the number of chunks in the grid.
`per_worker`	`int`	Number of grid points in each chunk, that is, in each call to the compiled `jax.vmap(fn)`.
`batch_size`	`int`	Number of grid points in each batch: `per_worker * workers`.

Raises:

Type	Description
`ValueError`	If `search_space` is empty or one of its axes isn't a non-empty 1-D array, if `batch_size` or `workers` is less than 1, if `devices` asks for more devices than are available or repeats a device, or if `workers` is less than the number of devices.
`TypeError`	If `fn` isn't callable, a name in `search_space` isn't a string, `batch_size` or `workers` isn't an integer, or `devices` has the wrong type.

Warns:

Type	Description
`UserWarning`	If `workers` is more than 31 per CPU device. `GridSweep` then runs 31 workers per CPU device.

Notes

Parallelism. On CPU, one JAX device runs computations from several threads in parallel, so up to 31 workers need no JAX_NUM_CPU_DEVICES or XLA_FLAGS setup. XLA runs at most 32 computations at a time on each CPU device, and a host callback that runs a JAX operation needs a free slot, so GridSweep runs at most 31 chunks at a time on each CPU device. To use more than 31 cores, set JAX_NUM_CPU_DEVICES to ceil(cores / 31) before you import JAX; devices=None then uses as many CPU devices as the workers need. Below 32 workers, extra virtual CPU devices only add compiles. To leave cores free on a shared machine, or in a container whose CPU quota Python can't see, set workers.
Compilation. Every call has the same shape, per_worker points, so fn is traced once and compiled once per device, and a ragged last chunk doesn't trigger a recompile. GridSweep pads that chunk with copies of its last grid point, never with synthetic values such as zeros, and drops the padding from the output.
Load balancing. Up to 2 * workers chunks are in flight at a time, across batch boundaries, so a slow chunk holds up only its own worker. Within a chunk, jax.vmap runs the solver loop until its slowest point finishes, so a point that hits max_steps still slows its whole chunk. Check converged in the results.
Reproducibility. For a given per_worker, results don't depend on workers or devices: a sweep is bitwise identical to a serial sweep (workers=1) with the same per_worker. per_worker is ceil(min(batch_size, total_points) / workers), so by default it depends on the core count and the grid size. A different chunk size compiles a different program, which can shift saved trajectories near violent collapses by more than the solver tolerances. To get the same chunks on another machine, pass workers and batch_size explicitly. A different CPU or JAX version can still change the results slightly.
Stopping early. If you stop iterating over batches early, for example with break, an exception, or Ctrl+C, GridSweep cancels the chunks that haven't started and waits for the running ones to finish.
Memory. Device memory holds at most workers chunks at a time, and host memory holds at most 2 * workers finished chunks plus the current batch.
Thread safety. Host callbacks inside fn, such as jax.debug.callback or jax.debug.print, can run at the same time from several threads, so they must be thread-safe. A host callback must not run JAX operations, such as comparing or reducing the jax.Array arguments of a jax.debug.callback; convert them with np.asarray first. Such an operation needs a free computation slot on the device, so it's slow, and it can deadlock if other JAX computations in the process fill the slots. Worker threads don't inherit JAX settings that a with block sets for the calling thread only, such as jax.enable_x64(...) or jax.debug_nans(...), so with more than one worker, the chunks run without them. Set them with jax.config.update instead.
Known issue. From JAX 0.11.1, XLA on CPU can deadlock when at least as many large FFTs inside a loop run at the same time as XLA has threads in its pool (jax-ml/jax#41265), and chunks that run at the same time add to that count. If fn runs FFTs inside the ODE solve and the sweep stops making progress, set XLA_FLAGS=--xla_cpu_multi_thread_eigen=false before you import JAX, as the issue suggests; each FFT then runs on one thread. To keep FFTs multithreaded, use the issue's other workaround: before you import JAX, set PJRT_NPROC to more than the number of FFTs in flight, which gives the pool more threads. Fewer workers also means fewer FFTs in flight, but that doesn't help when a single chunk runs enough FFTs at the same time.

Examples:

Sweep the peak expansion ratio of a free bubble over radius and pressure on every available core:

import jax.numpy as jnp

from jbubble import run_simulation
from jbubble.utils.gridsweep import GridSweep
from jbubble.utils.presets import free_bubble


def peak_ratio(R0, pressure):
    eom, pulse = free_bubble(R0=R0, pressure=pressure)
    result = run_simulation(eom, pulse)
    return {"ratio": result.radius.max() / R0, "ok": result.converged}


search_space = {
    "R0": jnp.linspace(1e-6, 5e-6, 20),
    "pressure": jnp.array([50e3, 200e3]),
}
grid = GridSweep(peak_ratio, search_space).run()
print(grid["ratio"].shape)  # (20, 2), a NumPy array
print(grid["ok"].all())  # True if every simulation converged

To leave cores free on a shared machine, limit the worker count:

grid = GridSweep(peak_ratio, search_space, workers=4).run()

Source code in jbubble/utils/gridsweep.py

def __init__(
    self,
    fn: Callable[..., PyTree],
    search_space: Mapping[str, Any],
    batch_size: int = 512,
    progress: bool = True,
    devices: int | Sequence[jax.Device] | None = None,
    workers: int | None = None,
) -> None:
    if not callable(fn):
        raise TypeError(f"fn must be callable, got {type(fn).__name__}.")
    if not search_space:
        raise ValueError("search_space must have at least one axis.")
    if not all(isinstance(k, str) for k in search_space):
        raise TypeError("search_space names must be strings.")
    batch_size = operator.index(batch_size)
    if batch_size < 1:
        raise ValueError(f"batch_size must be at least 1, got {batch_size}.")

    self.fn = fn
    self.search_space = search_space
    self.progress = progress

    # Stable ordering so multi-index resolution is deterministic
    self._keys = sorted(search_space)
    self._axes = [jnp.asarray(search_space[k]) for k in self._keys]
    for k, a in zip(self._keys, self._axes, strict=True):
        if a.ndim != 1 or a.shape[0] == 0:
            raise ValueError(
                f"search_space[{k!r}] must be a non-empty 1-D array, "
                f"got shape {a.shape}."
            )
    # Host copies: GridSweep gathers each chunk's parameters with NumPy.
    self._host_axes = [np.asarray(a) for a in self._axes]
    self._sizes = [a.shape[0] for a in self._host_axes]
    self._N = math.prod(self._sizes)

    resolved, workers = _resolve_devices_and_workers(devices, workers)
    self.per_worker = math.ceil(min(batch_size, self._N) / workers)
    self.workers = min(workers, math.ceil(self._N / self.per_worker))
    self.batch_size = self.per_worker * self.workers
    # Don't compile for devices that a small grid leaves idle.
    self.devices = resolved[: self.workers]

    self._vmapped = jax.jit(jax.vmap(self._call))
    self._executables: dict[jax.Device, Any] = {}

axes `property` ¶

axes

Parameter axes in sweep order, that is, sorted() order of the names.

grid_shape `property` ¶

grid_shape

Shape of the full parameter grid, one int per axis.

num_batches `property` ¶

num_batches

Number of batches that batches yields.

total_points `property` ¶

total_points

Total number of grid points.

batches ¶

batches() -> Generator[tuple[dict[str, ndarray], Any]]

Iterate lazily over the grid, one batch at a time.

The first iteration compiles fn. Batches arrive in grid order; each has batch_size points except possibly the last.

Yields:

Name	Type	Description
`params`	`dict[str, ndarray]`	Batch of parameter values: one array per axis name, with one value per grid point in the batch.
`outputs`	`Any`	Corresponding outputs of `fn`, as NumPy arrays with a leading batch axis.

Source code in jbubble/utils/gridsweep.py

def batches(
    self,
) -> Generator[tuple[dict[str, np.ndarray], PyTree]]:
    """Iterate lazily over the grid, one batch at a time.

    The first iteration compiles `fn`. Batches arrive in grid order;
    each has `batch_size` points except possibly the last.

    Yields
    ------
    params : dict[str, numpy.ndarray]
        Batch of parameter values: one array per axis name, with one
        value per grid point in the batch.
    outputs : PyTree
        Corresponding outputs of `fn`, as NumPy arrays with a leading
        batch axis.
    """
    buffer: list[PyTree] = []
    start = 0
    chunks = self._progress_chunks()
    try:
        for out in chunks:
            buffer.append(out)
            end = min(start + self.per_worker * len(buffer), self._N)
            if len(buffer) == self.workers or end == self._N:
                params = self._gather(np.arange(start, end))
                outputs = jax.tree.map(lambda *xs: np.concatenate(xs), *buffer)
                buffer, start = [], end
                yield params, outputs
    finally:
        chunks.close()

run ¶

run() -> Any

Run the full sweep and return a grid-shaped PyTree.

Returns:

Type	Description
`Any`	The structure that `fn` returns, with NumPy-array leaves of shape `(grid_shape, leaf_shape)`. `grid_shape` follows the axes in `search_space` in `sorted()` order of the names (Unicode code-point order, so uppercase names sort before lowercase ones). The grid is row-major: the last axis varies fastest.

Source code in jbubble/utils/gridsweep.py

def run(self) -> PyTree:
    """Run the full sweep and return a grid-shaped PyTree.

    Returns
    -------
    PyTree
        The structure that `fn` returns, with NumPy-array leaves of
        shape `(*grid_shape, *leaf_shape)`. `grid_shape` follows the axes
        in `search_space` in `sorted()` order of the names (Unicode
        code-point order, so uppercase names sort before lowercase
        ones). The grid is row-major: the last axis varies fastest.
    """
    # Fill preallocated arrays chunk by chunk, so the host never holds the
    # result twice.
    grid: list[np.ndarray] = []
    treedef = None
    start = 0
    chunks = self._progress_chunks()
    try:
        for c, out in enumerate(chunks):
            leaves, treedef = jax.tree.flatten(out)
            if c == 0:
                grid = [np.empty((self._N, *x.shape[1:]), x.dtype) for x in leaves]
            end = start + self._chunk_size(c)
            for dst, src in zip(grid, leaves, strict=True):
                dst[start:end] = src
            start = end
    finally:
        chunks.close()
    shape = self.grid_shape
    return jax.tree.unflatten(
        treedef, [x.reshape(*shape, *x.shape[1:]) for x in grid]
    )

Example¶

import jax.numpy as jnp
from jbubble.utils.gridsweep import GridSweep
from jbubble import run_simulation, SaveSpec
from jbubble.utils.presets import free_bubble


def peak_ratio(R0, pressure):
    preset = free_bubble(R0=R0, pressure=pressure)
    result = run_simulation(
        preset.eom,
        preset.pulse,
        save_spec=SaveSpec(500),
        t_max=10e-6,
    )
    return result.radius.max() / R0


sweep = GridSweep(
    fn=peak_ratio,
    search_space={
        "R0": jnp.linspace(1e-6, 5e-6, 20),
        "pressure": jnp.array([50e3, 100e3, 200e3, 400e3]),
    },
    batch_size=256,
)

grid = sweep.run()  # NumPy array of shape (20, 4)
print(sweep.grid_shape)  # (20, 4)
print(sweep.total_points)  # 80

Figure styles¶

Matplotlib style sheets in a light and a dark theme, with the colours of the jbubble docs and figures.

The jbubble figure palette in the light theme: eight colour swatches, four damped radius curves, a grey drive trace, and the Blues colour map

style ¶

Matplotlib styles for jbubble figures, in a light and a dark theme.

The package ships two Matplotlib style sheets as package data. To use one, pass its dotted name to matplotlib.pyplot.style.use:

import matplotlib.pyplot as plt

plt.style.use("jbubble.style.light")  # or "jbubble.style.dark"

Matplotlib 3.7 and later find a style sheet inside an importable package, so you need no file paths. Neither theme changes anything outside Matplotlib's rcParams, and this package doesn't import Matplotlib, so jbubble keeps Matplotlib an optional dependency (the examples extra).

Both themes draw on a solid background (#ffffff for light, #0d1117 for dark), hide the top and right spines, and add a hairline grid. The colour cycle has eight categorical slots, chosen for contrast against each background and for colour-vision deficiencies:

Slot	Name	Light	Dark	Typical role
`C0`	azure	`#0d77ca`	`#3e96ea`	the bubble radius
`C1`	coral	`#e35a2d`	`#e6653c`	the first comparison
`C2`	teal	`#129483`	`#1ea28f`	the second comparison
`C3`	plum	`#9f1169`	`#b51c79`	the third comparison
`C4`	green	`#438f31`	`#4f9b3e`
`C5`	violet	`#7447c8`	`#8f6ce0`
`C6`	amber	`#da9516`	`#c0851f`	below 3:1 contrast on white
`C7`	crimson	`#af2938`	`#c83846`

Every pair among the first four slots stays distinguishable in both themes. With more than four series on one set of axes, add direct labels, line styles, or small multiples, or use a sequential colour map for ordered data (Blues in the light theme and Blues_r in the dark theme, which are the image defaults). Draw the acoustic drive in the neutral grey #8c959f (light) or #6e7681 (dark), ideally on its own axes.

The style sheets are cosmetic. Their colours and sizes can change in any release and aren't covered by jbubble's versioning policy.

utils¶

Presets¶

BubblePreset ¶

free_bubble ¶

lipid_bubble ¶

thick_shell_bubble ¶

GridSweep¶

GridSweep ¶

axes property ¶

grid_shape property ¶

num_batches property ¶

total_points property ¶

batches ¶

run ¶

Example¶

Figure styles¶

style ¶

axes `property` ¶

grid_shape `property` ¶

num_batches `property` ¶

total_points `property` ¶