two_fluid_dispersion_solution

plasmapy.dispersion.two_fluid_dispersion.two_fluid_dispersion_solution(*, B: Unit(‘T’), ion: Union[str, plasmapy.particles.particle_class.Particle], k: Unit(‘rad / m’), n_i: Unit(‘1 / m3’), T_e: Unit(‘K’), T_i: Unit(‘K’), theta: Unit(‘deg’), gamma_e: Union[float, int] = 1, gamma_i: Union[float, int] = 3, z_mean: Union[float, int] = None)

Using the solution provided by Bellan 2012, calculate the analytical solution to the two fluid, low-frequency (\(\omega/kc \ll 1\)) dispersion relation presented by Stringer 1963. This dispersion relation also assummes a uniform magnetic field \(\mathbf{B_o}\), no D.C. electric field \(\mathbf{E_o}=0\), and quasi-neutrality. For more information see the Notes section below.

Aliases: tfds_

Parameters
  • B (Quantity) – The magnetic field magnitude in units convertible to \(T\).

  • ion (str or Particle) – Representation of the ion species (e.g., 'p' for protons, 'D+' for deuterium, 'He-4 +1' for singly ionized helium-4, etc.). If no charge state information is provided, then the ions are assumed to be singly ionized.

  • k (Quantity, single valued or 1-D array) – Wavenumber in units convertible to \(rad / m\). Either single valued or 1-D array of length \(N\).

  • n_i (Quantity) – Ion number density in units convertible to \(m^{-3}\).

  • T_e (Quantity) – The electron temperature in units of \(K\) or \(eV\).

  • T_i (Quantity) – The ion temperature in units of \(K\) or \(eV\).

  • theta (Quantity, single valued or 1-D array) – The angle of propagation of the wave with respect to the magnetic field, \(\cos^{-1}(k_z / k)\), in units must be convertible to \(deg\). Either single valued or 1-D array of size \(M\).

  • gamma_e (float or int, optional) – The adiabatic index for electrons, which defaults to 1. This value assumes that the electrons are able to equalize their temperature rapidly enough that the electrons are effectively isothermal.

  • gamma_i (float or int, optional) – The adiabatic index for ions, which defaults to 3. This value assumes that ion motion has only one degree of freedom, namely along magnetic field lines.

  • z_mean (float or int, optional) – The average ionization state (arithmetic mean) of the ion composing the plasma. Will override any charge state defined by argument ion.

Returns

omega – A dictionary of computed wave frequencies in units \(rad/s\). The dictionary contains three keys: 'fast_mode' for the fast mode, 'alfven_mode' for the Alfvén mode, and 'acoustic_mode' for the ion-acoustic mode. The value for each key will be a \(N x M\) array.

Return type

Dict[str, Quantity]

Raises
Warns

PhysicsWarning – When the computed wave frequencies violate the low-frequency (\(\omega/kc \ll 1\)) assumption of the dispersion relation.

Notes

The complete dispersion equation presented by Springer 1963 2 (equation 1 of Bellan 2012 1) is:

\[\begin{split}\left( \cos^2 \theta - Q \frac{\omega^2}{k^2 {v_A}^2} \right) & \left[ \left( \cos^2 \theta - \frac{\omega^2}{k^2 {c_s}^2} \right) - Q \frac{\omega^2}{k^2 {v_A}^2} \left( 1 - \frac{\omega^2}{k^2 {c_s}^2} \right) \right] \\ &= \left(1 - \frac{\omega^2}{k^2 {c_s}^2} \right) \frac{\omega^2}{{\omega_{ci}}^2} \cos^2 \theta\end{split}\]

where

\[\begin{split}Q &= 1 + k^2 c^2/{\omega_{pe}}^2 \\ \cos \theta &= \frac{k_z}{k} \\ \mathbf{B_o} &= B_{o} \mathbf{\hat{z}}\end{split}\]

\(\omega\) is the wave frequency, \(k\) is the wavenumber, \(v_A\) is the Alfvén velocity, \(c_s\) is the sound speed, \(\omega_{ci}\) is the ion gyrofrequency, and \(\omega_{pe}\) is the electron plasma frequency. This relation does additionally assume low-frequency waves \(\omega/kc \ll 1\), no D.C. electric field \(\mathbf{E_o}=0\) and quasi-neutrality.

Following section 5 of Bellan 2012 1 the exact roots of the above dispersion equation can be derived and expressed as one analytical solution (equation 38 of Bellan 2012 1):

\[\frac{\omega}{\omega_{ci}} = \sqrt{ 2 \Lambda \sqrt{-\frac{P}{3}} \cos\left( \frac{1}{3} \cos^{-1}\left( \frac{3q}{2p} \sqrt{-\frac{3}{p}} \right) - \frac{2 \pi}{3}j \right) + \frac{\Lambda A}{3} }\]

where \(j = 0\) represents the fast mode, \(j = 1\) represents the Alfvén mode, and \(j = 2\) represents the acoustic mode. Additionally,

\[\begin{split}p &= \frac{3B-A^2}{3} \; , \; q = \frac{9AB-2A^3-27C}{27} \\ A &= \frac{Q + Q^2 \beta + Q \alpha + \alpha \Lambda}{Q^2} \; , \; B = \alpha \frac{1 + 2 Q \beta + \Lambda \beta}{Q^2} \; , \; C = \frac{\alpha^2 \beta}{Q^2} \\ \alpha &= \cos^2 \theta \; , \; \beta = \left( \frac{c_s}{v_A}\right)^2 \; , \; \Lambda = \left( \frac{k v_{A}}{\omega_{ci}}\right)^2\end{split}\]

References

1(1,2,3)

PM Bellan, Improved basis set for low frequency plasma waves, 2012, JGR, 117, A12219, doi: 10.1029/2012JA017856.

2

TE Stringer, Low-frequency waves in an unbounded plasma, 1963, JNE, Part C, doi: 10.1088/0368-3281/5/2/304

Examples

>>> from astropy import units as u
>>> from plasmapy.dispersion import two_fluid_dispersion
>>> inputs = {
...     "k": 0.01 * u.rad / u.m,
...     "theta": 30 * u.deg,
...     "B": 8.3e-9 * u.T,
...     "n_i": 5e6 * u.m ** -3,
...     "T_e": 1.6e6 * u.K,
...     "T_i": 4.0e5 * u.K,
...     "ion": "p+",
... }
>>> omegas = two_fluid_dispersion_solution(**inputs)
>>> omegas
{'fast_mode': <Quantity 1520.57... rad / s>,
 'alfven_mode': <Quantity 1261.75... rad / s>,
 'acoustic_mode': <Quantity 0.688152... rad / s>}
>>> inputs = {
...     "k": [1e-7, 2e-7] * u.rad / u.m,
...     "theta": [10, 20] * u.deg,
...     "B": 8.3e-9 * u.T,
...     "n_i": 5e6 * u.m ** -3,
...     "T_e": 1.6e6 * u.K,
...     "T_i": 4.0e5 * u.K,
...     "ion": "He+",
... }
>>> omegas = two_fluid_dispersion_solution(**inputs)
>>> omegas['fast_mode']
<Quantity [[0.00767..., 0.00779... ],
           [0.01534..., 0.01558...]] rad / s>