cold_plasma_permittivity_SDP¶
-
plasmapy.formulary.dielectric.
cold_plasma_permittivity_SDP
(B: Unit("T"), species, n, omega: Unit("rad / s"))¶ Magnetized Cold Plasma Dielectric Permittivity Tensor Elements.
Elements (S, D, P) are given in the “Stix” frame, ie. with B // z.
The \(\exp(-i \omega t)\) time-harmonic convention is assumed.
Parameters: - B (Quantity) – Magnetic field magnitude in units convertible to tesla.
- species (list of str) – List of the plasma particle species e.g.: [‘e’, ‘D+’] or [‘e’, ‘D+’, ‘He+’].
- n (list of ~astropy.units.Quantity) –
list
of species density in units convertible to per cubic meter The order of the species densities should follow species. - omega (Quantity) – Electromagnetic wave frequency in rad/s.
Returns: - sum (~astropy.units.Quantity) – S (“Sum”) dielectric tensor element.
- difference (~astropy.units.Quantity) – D (“Difference”) dielectric tensor element.
- plasma (~astropy.units.Quantity) – P (“Plasma”) dielectric tensor element.
Notes
The dielectric permittivity tensor is expressed in the Stix frame with the \(\exp(-i \omega t)\) time-harmonic convention as \(\varepsilon = \varepsilon_0 A\), with \(A\) being
\[\begin{split}\varepsilon = \varepsilon_0 \left(\begin{matrix} S & -i D & 0 \\ +i D & S & 0 \\ 0 & 0 & P \end{matrix}\right)\end{split}\]where:
\[ \begin{align}\begin{aligned}S = 1 - \sum_s \frac{\omega_{p,s}^2}{\omega^2 - \Omega_{c,s}^2}\\D = \sum_s \frac{\Omega_{c,s}}{\omega} \frac{\omega_{p,s}^2}{\omega^2 - \Omega_{c,s}^2}\\P = 1 - \sum_s \frac{\omega_{p,s}^2}{\omega^2}\end{aligned}\end{align} \]where \(\omega_{p,s}\) is the plasma frequency and \(\Omega_{c,s}\) is the signed version of the cyclotron frequency for the species \(s\).
References
- T.H. Stix, Waves in Plasma, 1992.
Examples
>>> from astropy import units as u >>> from numpy import pi >>> B = 2*u.T >>> species = ['e', 'D+'] >>> n = [1e18*u.m**-3, 1e18*u.m**-3] >>> omega = 3.7e9*(2*pi)*(u.rad/u.s) >>> permittivity = S, D, P = cold_plasma_permittivity_SDP(B, species, n, omega) >>> S <Quantity 1.02422...> >>> permittivity.sum # namedtuple-style access <Quantity 1.02422...> >>> D <Quantity 0.39089...> >>> P <Quantity -4.8903...>