r'''A collection of plasma physics related formulae.
.. todo::
Extended implementation should be done in
:mod:`irfpy.util.physics`.
Collection of formulae which are really frequently used.
- plasma frequency (:func:`plasma_freq_mksa`)
- gyro frequency
- gyro radius
- gyro period
- Debye length
- Ion acoustic velocity
- Electron-volt to kelvin conversion (:func:`eV2K`)
The unit system is complecated.
MKSA is very useful for theory, but not for practical.
In this module, MKSA is used only when _mksa is specified in the method.
Practically, we uses ``km``, ``amu``, ``sec``, ``nT``, ``eV``, and so on.
Thus, most of functions uses the :abbr:`PPUS ('plasma physisist unit system)`
(plasma physisist unit system).
For example, there are two methods to calculate Debye length.
:func:`debye_length_mksa` is for MKSA, and :func:`debye_length` is in PPUS.
This is for one of the typical solar wind conditon.
Typical debye length is 10m.
>>> print('%.2f' % debye_length_mksa(1.6e-19, 3e6, 50000))
8.92
This is for similar condition in plasma physist unit,
i.e. q=1, n=3 cm\ :superscript:`-3`, and T=5eV. Returned unit is in km.
>>> print('%.4f' % debye_length(1, 3, 5))
0.0096
List of formulation
--------------------
.. math::
\omega_\mathrm{pe} \mathrm{[/s]} = 5.641\times 10^4
\times \sqrt{n \mathrm{[/cc]}}
f_\mathrm{pe} \mathrm{[/s]} = 8.979\times 10^3
\times \sqrt{n \mathrm{[/cc]}}
\omega_\mathrm{pi} \mathrm{[/s]} = 1.317\times 10^3
\times \sqrt{n \mathrm{[/cc]}}
f_\mathrm{pi} \mathrm{[/s]} = 2.095\times 10^2
\times \sqrt{n \mathrm{[/cc]}}
\omega_\mathrm{ge} \mathrm{[/s]} = 1.759\times 10^2
\times B \mathrm{[nT]}
f_\mathrm{ge} \mathrm{[/s]} = 2.799\times 10^1
\times B \mathrm{[nT]}
\omega_\mathrm{gi} \mathrm{[/s]} = 9.579\times 10^{-2}
\times B \mathrm{[nT]}
f_\mathrm{gi} \mathrm{[/s]} = 1.525\times 10^{-2}
\times B \mathrm{[nT]}
T_\mathrm{ge} \mathrm{[s]} = 3.572\times 10^{-2}
\times \frac{1}{B\mathrm{[nT]}}
T_\mathrm{gi} \mathrm{[s]} = 6.559\times 10^1
\times \frac{1}{B \mathrm{[nT]}}
R_\mathrm{ge} \mathrm{[km]} = 5.686\times 10^{-3} \times
\frac{v \mathrm{[km/s]}}{B \mathrm{[nT]}}
R_\mathrm{gi} \mathrm{[km]} = 1.044\times 10^1 \times
\frac{v \mathrm{[km/s]}}{B \mathrm{[nT]}}
\lambda_\mathrm{De} \mathrm{[km]} = 7.434\times 10^{-3}
\times \sqrt{\frac{T_\mathrm{e} \mathrm{[eV]}}{n \mathrm{[/cc]}}}
'''
from irfpy.util import constant as k
import numpy as _np
import warnings
warnings.warn('''The irfpy.util.formula module is deprecated. Use irfpy.util.physics module instead''', DeprecationWarning)
_twopi = _np.pi * 2
[docs]def plasma_afreq_mksa(mass, charge, density):
''' Returns the plasma angular frequency, wp=sqrt(n e**2)/(m eps0)
'''
return _np.sqrt(density / mass / k.eps0) * charge
[docs]def plasma_afreq(mass_amu, charge_qe, density__cm3):
''' Returns the plasma angular frequency (wp) in plasma unit.
'''
mass = mass_amu * k.amu
charge = charge_qe * k.qe
density = density__cm3 * 1e6
return plasma_afreq_mksa(mass, charge, density)
[docs]def plasma_afreq_el(density__cm3):
''' Returns the electron's plasma angular frequency (wpe).
'''
mass = k.me
charge = k.qe
density = density__cm3 * 1e6
return plasma_afreq_mksa(mass, charge, density)
[docs]def plasma_afreq_pr(density__cm3):
''' Returns the proton's plasma angular frequency (wppi).
'''
mass = k.mp
charge = k.qe
density = density__cm3 * 1e6
return plasma_afreq_mksa(mass, charge, density)
[docs]def plasma_freq_mksa(mass, charge, density):
''' Returns the plasma frequency, fp = wp / 2 pi.
'''
return plasma_afreq_mksa(mass, charge, density) / _twopi
[docs]def plasma_freq(mass_amu, charge_qe, density__cm3):
return plasma_afreq(mass_amu, charge_qe, density__cm3) / _twopi
[docs]def plasma_freq_el(density__cm3):
return plasma_afreq_el(density__cm3) / _twopi
[docs]def plasma_freq_pr(density__cm3):
return plasma_afreq_pr(density__cm3) / _twopi
[docs]def gyro_afreq_mksa(mass, charge, magfield):
''' Return the gyro angular frequency in Hz, wg = qB / m.
'''
return charge * magfield / mass
[docs]def gyro_freq_mksa(mass, charge, magfield):
''' Return the gyro angular frequency in Hz, fg = qB / 2 pi m.
'''
return gyro_afreq_mksa(mass, charge, magfield) / _twopi
[docs]def gyro_afreq(mass_amu, charge_qe, magfield_nT):
m = mass_amu * k.amu
q = charge_qe * k.qe
B = magfield_nT * (1e-9)
return gyro_afreq_mksa(m, q, B)
[docs]def gyro_freq(mass_amu, charge_eq, magfield_nT):
return gyro_afreq_mksa(mass_amu, charge_eq, magfield_nT) / _twopi
[docs]def gyro_afreq_el(magfield_nT):
return gyro_afreq_mksa(k.me, k.qe, magfield_nT * 1e-9)
[docs]def gyro_freq_el(magfield_nT):
return gyro_freq_mksa(k.me, k.qe, magfield_nT * 1e-9)
[docs]def gyro_afreq_pr(magfield_nT):
return gyro_afreq_mksa(k.mp, k.qe, magfield_nT * 1e-9)
[docs]def gyro_freq_pr(magfield_nT):
return gyro_freq_mksa(k.mp, k.qe, magfield_nT * 1e-9)
[docs]def gyro_period_mksa(mass, charge, magfield):
return 1. / gyro_freq_mksa(mass, charge, magfield)
[docs]def gyro_period(mass_amu, charge_qe, magfield_nT):
return 1. / gyro_freq(mass_amu, charge_qe, magfield_nT)
[docs]def gyro_period_el(magfield_nT):
return 1. / gyro_freq_mksa(k.me, k.qe, magfield_nT * 1e-9)
[docs]def gyro_period_pr(magfield_nT):
return 1. / gyro_freq_mksa(k.mp, k.qe, magfield_nT * 1e-9)
[docs]def gyro_radius_mksa(mass, charge, magfield, velocity):
''' Returns a gyro_radius, rg = v / wg = mv / qB
'''
return mass * velocity / (charge * magfield)
[docs]def gyro_radius(mass_amu, charge_qe, magfield_nT, velocity_km_s):
''' Returns a gyro_radius in the unit system of plasma scientist.
'''
rg_m = gyro_radius_mksa(mass_amu * k.amu, charge_qe * k.qe,
magfield_nT * 1e-9, velocity_km_s * 1e3)
return rg_m * 1e-3 # unit of km/s
[docs]def gyro_radius_el(magfield_nT, velocity_km_s):
rg_m = gyro_radius_mksa(k.me, k.qe, magfield_nT * 1e-9,
velocity_km_s * 1e3)
return rg_m * 1e-3
[docs]def gyro_radius_pr(magfield_nT, velocity_km_s):
rg_m = gyro_radius_mksa(k.mp, k.qe, magfield_nT * 1e-9,
velocity_km_s * 1e3)
return rg_m * 1e-3
[docs]def debye_length_mksa(charge, density, temperature):
''' Returns the Debye length, Ld = sqrt( eps0 kB T / n e**2)
'''
return _np.sqrt(k.eps0 * k.kB * temperature / density) / charge
[docs]def debye_length(charge_qe, density__cc, temperature_eV):
''' Debye length in PPUS.
'''
q = charge_qe * k.qe
n = density__cc * 1e6
T = eV2K(temperature_eV)
dm = debye_length_mksa(q, n, T)
return dm * 1e-3 # in km
[docs]def debye_length_el(density__cc, temperature_eV):
return debye_length(1, density__cc, temperature_eV)
[docs]def debye_length_pr(density__cc, temperature_eV):
return debye_length(1, density__cc, temperature_eV)
[docs]def ion_acoustic_speed_mksa(temperature_el, mass_ion=k.mp, gamma=3):
''' Returns the ion acoustic speed in m/s.
The ion acoustic speed, S0, is calculated by
S0 = sqrt( gamma * kB * Te / mi )
@param temperature_el electron temperature in K
@param mass_ion Ion mass in kg
@param gamma The ratio of specific heat for electrons.
'''
S0 = _np.sqrt(gamma * k.kB * temperature_el / mass_ion)
return S0
[docs]def ion_acoustic_speed(temperature_el_eV, mass_ion_amu, gamma=3):
''' Returns the ion acoustic speed in km/s.
The ion acoustic speed, S0, is calculated by
S0 = sqrt( gamma * kB * Te / mi )
@param temperature_el electron temperature in K
@param mass_ion Ion mass in kg
@param gamma The ratio of specific heat for electrons.
'''
t = eV2K(temperature_el_eV)
m = mass_ion_amu * k.amu
return ion_acoustic_speed_mksa(t, m, gamma) * 1e-3 # in km/s
[docs]def eV2K(temp_eV):
r''' Convert from eV to Kelvin
Electron volt is used to specify the temperature.
Here the function provides the conversion from eV to K.
The factor is :math:`q_e/k_B\sim 1.16\times 10^4`.
:Parameters:
temp_eV : float
The temperature in electron volts.
:Returns:
temp_K : float
Corresponding temperature in Kelvins.
'''
return k.qe / k.kB * temp_eV # More precise values are needed.