def test_upper_hybrid_frequency(): r"""Test the upper_hybrid_frequency function in parameters.py.""" omega_uh = upper_hybrid_frequency(B, n_e=n_e) omega_ce = gyrofrequency(B) omega_pe = plasma_frequency(n=n_e) assert omega_ce.unit.is_equivalent(u.rad / u.s) assert omega_pe.unit.is_equivalent(u.rad / u.s) assert omega_uh.unit.is_equivalent(u.rad / u.s) left_hand_side = omega_uh**2 right_hand_side = omega_ce**2 + omega_pe**2 assert np.isclose(left_hand_side.value, right_hand_side.value) with pytest.raises(ValueError): upper_hybrid_frequency(5 * u.T, n_e=-1 * u.m**-3) with pytest.warns(u.UnitsWarning): assert upper_hybrid_frequency(1.2, 1.3) == upper_hybrid_frequency( 1.2 * u.T, 1.3 * u.m**-3) with pytest.warns(u.UnitsWarning): assert upper_hybrid_frequency(1.4 * u.T, 1.3) == upper_hybrid_frequency( 1.4, 1.3 * u.m**-3) assert_can_handle_nparray(upper_hybrid_frequency)
def test_lower_hybrid_frequency(): r"""Test the lower_hybrid_frequency function in parameters.py.""" ion = 'He-4 1+' omega_ci = gyrofrequency(B, particle=ion) omega_pi = plasma_frequency(n=n_i, particle=ion) omega_ce = gyrofrequency(B) omega_lh = lower_hybrid_frequency(B, n_i=n_i, ion=ion) assert omega_ci.unit.is_equivalent(u.rad / u.s) assert omega_pi.unit.is_equivalent(u.rad / u.s) assert omega_ce.unit.is_equivalent(u.rad / u.s) assert omega_lh.unit.is_equivalent(u.rad / u.s) left_hand_side = omega_lh**-2 right_hand_side = 1 / (omega_ci**2 + omega_pi**2) + omega_ci**-1 * omega_ce**-1 assert np.isclose(left_hand_side.value, right_hand_side.value) with pytest.raises(ValueError): lower_hybrid_frequency(0.2 * u.T, n_i=5e19 * u.m**-3, ion='asdfasd') with pytest.raises(ValueError): lower_hybrid_frequency(0.2 * u.T, n_i=-5e19 * u.m**-3, ion='asdfasd') with pytest.raises(ValueError): lower_hybrid_frequency(np.nan * u.T, n_i=-5e19 * u.m**-3, ion='asdfasd') with pytest.warns(u.UnitsWarning): assert lower_hybrid_frequency(1.3, 1e19) == lower_hybrid_frequency( 1.3 * u.T, 1e19 * u.m**-3) assert_can_handle_nparray(lower_hybrid_frequency)
def cold_plasma_LRP(B, species, n, omega ,flag= True ): """ @param B: Magnetic field in z direction @param species: ['e','p'] means electrons and protons @param n: The density of different species @param omega: wave_frequency @param flag: if flag == 1 then omega means to be a number, if flag is not True the omega means a symbol @return: (l, R, P) """ if flag: wave_frequency = omega w = wave_frequency.value else: wave_frequency = symbols('w') w = wave_frequency L, R, P = 1, 1, 1 for s, n_s in zip(species, n): omega_c = parameters.gyrofrequency(B=B, particle=s, signed=True) omega_p = parameters.plasma_frequency(n=n_s, particle=s) L += - omega_p.value ** 2 / (w * (w - omega_c.value)) R += - omega_p.value ** 2 / (w * (w + omega_c.value)) P += - omega_p.value ** 2 / w ** 2 return L, R, P
def __init__(self, name, density, magnetic_field, Tx, Tz): self.name = name self.density = density self.gyro_f = parameters.gyrofrequency(B=magnetic_field, particle=name, signed=True) self.plasma_f = parameters.plasma_frequency(density, particle=name) self.Tx = Tx self.Tz = Tz self.vth_x = parameters.thermal_speed(Tx * u.K, name).value self.vth_z = parameters.thermal_speed(Tz * u.K, name).value
def test_plasma_frequency(): r"""Test the plasma_frequency function in parameters.py.""" assert plasma_frequency(n_e).unit.is_equivalent(u.rad / u.s) assert np.isclose(plasma_frequency(1 * u.cm**-3).value, 5.64e4, rtol=1e-2) with pytest.raises(TypeError): plasma_frequency(u.m**-3) with pytest.raises(u.UnitConversionError): plasma_frequency(5 * u.m**-2) with pytest.raises(ValueError): plasma_frequency(np.nan * u.m**-3) with pytest.warns(u.UnitsWarning): assert plasma_frequency(1e19) == plasma_frequency(1e19 * u.m**-3) assert plasma_frequency(n_i, particle='p').unit.is_equivalent(u.rad / u.s) # Case where Z=1 is assumed assert plasma_frequency(n_i, particle='H-1+') == plasma_frequency(n_i, particle='p') assert np.isclose(plasma_frequency(mu * u.cm**-3, particle='p').value, 1.32e3, rtol=1e-2) with pytest.raises(ValueError): plasma_frequency(n=5 * u.m**-3, particle='sdfas') with pytest.warns(u.UnitsWarning): plasma_freq_no_units = plasma_frequency(1e19, particle='p') assert plasma_freq_no_units == plasma_frequency(1e19 * u.m**-3, particle='p') plasma_frequency(1e17 * u.cm**-3, particle='p') # testing for user input z_mean testMeth1 = plasma_frequency(1e17 * u.cm**-3, particle='p', z_mean=0.8).si.value testTrue1 = 333063562455.4028 errStr = f"plasma_frequency() gave {testMeth1}, should be {testTrue1}." assert np.isclose(testMeth1, testTrue1, atol=0.0, rtol=1e-15), errStr assert_can_handle_nparray(plasma_frequency)
def time_plasma_frequency(self): plasma_frequency(1e19*u.m**-3, particle='p', to_hz=True)
def cold_plasma_permittivity_SDP(B, species, n, omega): r""" Magnetized Cold Plasma Dielectric Permittivity Tensor Elements. Elements (S, D, P) are given in the "Stix" frame, ie. with B // z. The :math:`\exp(-i \omega t)` time-harmonic convention is assumed. Parameters ---------- B : ~astropy.units.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 : ~astropy.units.Quantity Electromagnetic wave frequency in rad/s. Returns ------- S : ~astropy.units.Quantity S ("Sum") dielectric tensor element. D : ~astropy.units.Quantity D ("Difference") dielectric tensor element. P : ~astropy.units.Quantity P ("Plasma") dielectric tensor element. Notes ----- The dielectric permittivity tensor is expressed in the Stix frame with the :math:`\exp(-i \omega t)` time-harmonic convention as :math:`\varepsilon = \varepsilon_0 A`, with :math:`A` being .. math:: \varepsilon = \varepsilon_0 \left(\begin{matrix} S & -i D & 0 \\ +i D & S & 0 \\ 0 & 0 & P \end{matrix}\right) where: .. math:: 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} where :math:`\omega_{p,s}` is the plasma frequency and :math:`\Omega_{c,s}` is the signed version of the cyclotron frequency for the species :math:`s`. References ---------- - T.H. Stix, Waves in Plasma, 1992. Examples -------- >>> from astropy import units as u >>> from plasmapy.constants 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) >>> S, D, P = cold_plasma_permittivity_SDP(B, species, n, omega) >>> S <Quantity 1.02422902> >>> D <Quantity 0.39089352> >>> P <Quantity -4.8903104> """ S, D, P = 1, 0, 1 for s, n_s in zip(species, n): omega_c = parameters.gyrofrequency(B=B, particle=s, signed=True) omega_p = parameters.plasma_frequency(n=n_s, particle=s) S += -omega_p**2 / (omega**2 - omega_c**2) D += omega_c / omega * omega_p**2 / (omega**2 - omega_c**2) P += -omega_p**2 / omega**2 return S, D, P
def permittivity_1D_Maxwellian(omega, kWave, T, n, particle, z_mean=None): r""" The classical dielectric permittivity for a 1D Maxwellian plasma. This function can calculate both the ion and electron permittivities. No additional effects are considered (e.g. magnetic fields, relativistic effects, strongly coupled regime, etc.) Parameters ---------- omega : ~astropy.units.Quantity The frequency in rad/s of the electromagnetic wave propagating through the plasma. kWave : ~astropy.units.Quantity The corresponding wavenumber, in rad/m, of the electromagnetic wave propagating through the plasma. This is often modulated by the dispersion of the plasma or by relativistic effects. See em_wave.py for ways to calculate this. T : ~astropy.units.Quantity The plasma temperature - this can be either the electron or the ion temperature, but should be consistent with density and particle. n : ~astropy.units.Quantity The plasma density - this can be either the electron or the ion density, but should be consistent with temperature and particle. particle : str The plasma particle species. z_mean : str The average ionization of the plasma. This is only required for calculating the ion permittivity. Returns ------- chi : ~astropy.units.Quantity The ion or the electron dielectric permittivity of the plasma. This is a dimensionless quantity. Notes ----- The dielectric permittivities for a Maxwellian plasma are described by the following equations [1]_ .. math:: \chi_e(k, \omega) = - \frac{\alpha_e^2}{2} Z'(x_e) \chi_i(k, \omega) = - \frac{\alpha_i^2}{2}\frac{Z}{} Z'(x_i) \alpha = \frac{\omega_p}{k v_{Th}} x = \frac{\omega}{k v_{Th}} :math:`chi_e` and :math:`chi_i` are the electron and ion permittivities respectively. :math:`Z'` is the derivative of the plasma dispersion function. :math:`\alpha` is the scattering parameter which delineates the difference between the collective and non-collective Thomson scattering regimes. :math:`x` is the dimensionless phase velocity of the EM wave propagating through the plasma. References ---------- .. [1] J. Sheffield, D. Froula, S. H. Glenzer, and N. C. Luhmann Jr, Plasma scattering of electromagnetic radiation: theory and measurement techniques. Chapter 5 Pg 106 (Academic press, 2010). Example ------- >>> from astropy import units as u >>> from plasmapy.constants import pi, c >>> T = 30 * 11600 * u.K >>> n = 1e18 * u.cm**-3 >>> particle = 'Ne' >>> z_mean = 8 * u.dimensionless_unscaled >>> vTh = parameters.thermal_speed(T, particle, method="most_probable") >>> omega = 5.635e14 * 2 * pi * u.rad / u.s >>> kWave = omega / vTh >>> permittivity_1D_Maxwellian(omega, kWave, T, n, particle, z_mean) <Quantity -6.72809257e-08+5.76037956e-07j> """ # thermal velocity vTh = parameters.thermal_speed(T=T, particle=particle, method="most_probable") # plasma frequency wp = parameters.plasma_frequency(n=n, particle=particle, z_mean=z_mean) # scattering parameter alpha. # explicitly removing factor of sqrt(2) to be consistent with Froula alpha = np.sqrt(2) * (wp / (kWave * vTh)).to(u.dimensionless_unscaled) # The dimensionless phase velocity of the propagating EM wave. zeta = (omega / (kWave * vTh)).to(u.dimensionless_unscaled) chi = alpha**2 * (-1 / 2) * plasma_dispersion_func_deriv(zeta.value) return chi.to(u.dimensionless_unscaled)
def cold_plasma_permittivity_LRP(B: u.T, species, n, omega: u.rad / u.s): r""" Magnetized Cold Plasma Dielectric Permittivity Tensor Elements. Elements (L, R, P) are given in the "rotating" basis, ie. in the basis :math:`(\mathbf{u}_{+}, \mathbf{u}_{-}, \mathbf{u}_z)`, where the tensor is diagonal and with B // z. The :math:`\exp(-i \omega t)` time-harmonic convention is assumed. Parameters ---------- B : ~astropy.units.Quantity Magnetic field magnitude in units convertible to tesla. species : list of str 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 : ~astropy.units.Quantity Electromagnetic wave frequency in rad/s. Returns ------- L : ~astropy.units.Quantity L ("Left") Left-handed circularly polarization tensor element. R : ~astropy.units.Quantity R ("Right") Right-handed circularly polarization tensor element. P : ~astropy.units.Quantity P ("Plasma") dielectric tensor element. Notes ----- In the rotating frame defined by :math:`(\mathbf{u}_{+}, \mathbf{u}_{-}, \mathbf{u}_z)` with :math:`\mathbf{u}_{\pm}=(\mathbf{u}_x \pm \mathbf{u}_y)/\sqrt{2}`, the dielectric tensor takes a diagonal form with elements L, R, P with: .. math:: L = 1 - \sum_s \frac{\omega_{p,s}^2}{\omega\left(\omega - \Omega_{c,s}\right)} R = 1 - \sum_s \frac{\omega_{p,s}^2}{\omega\left(\omega + \Omega_{c,s}\right)} P = 1 - \sum_s \frac{\omega_{p,s}^2}{\omega^2} where :math:`\omega_{p,s}` is the plasma frequency and :math:`\Omega_{c,s}` is the signed version of the cyclotron frequency for the species :math:`s`. References ---------- - T.H. Stix, Waves in Plasma, 1992. Examples -------- >>> from astropy import units as u >>> from plasmapy.constants 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) >>> L, R, P = cold_plasma_permittivity_LRP(B, species, n, omega) >>> L <Quantity 0.63333549> >>> R <Quantity 1.41512254> >>> P <Quantity -4.8903104> """ L, R, P = 1, 1, 1 for s, n_s in zip(species, n): omega_c = parameters.gyrofrequency(B=B, particle=s, signed=True) omega_p = parameters.plasma_frequency(n=n_s, particle=s) L += -omega_p**2 / (omega * (omega - omega_c)) R += -omega_p**2 / (omega * (omega + omega_c)) P += -omega_p**2 / omega**2 return L, R, P