Saha¶
-
plasmapy.formulary.ionization.
Saha
(g_j, g_k, n_e: Unit("1 / m3"), E_jk: Unit("J"), T_e: Unit("K")) -> Unit(dimensionless)¶ The Saha equation, derived in statistical mechanics, gives an approximation of the ratio of population of ions in two different ionization states in a plasma. This approximation applies to plasmas in thermodynamic equilibrium where ionization and recombination of ions with electrons are balanced.
\[\frac{N_j}{N_k} = \frac{1}{n_e} \frac{g_j}{4 g_k a_0^{3}} \left( \frac{k_B T_e}{\pi E_H} \right)^{\frac{3}{2}} \exp\left( \frac{-E_{jk}}{k_B T_e} \right)\]Where \(k_B\) is the Boltzmann constant, \(a_0\) is the Bohr radius, \(E_H\) is the ionization energy of Hydrogen, \(N_j\) and \(N_k\) are the population of ions in the j and k states respectively. This function is equivalent to Eq. 3.47 in Drake.
Parameters: Warns: UnitsWarning
– If units are not provided, SI units are assumed.Raises: TypeError
– TheT_e
,E_jk
, orn_e
is not aQuantity
and cannot be converted into a ~astropy.units.Quantity.UnitConversionError
– If theT_e
,E_jk
, orn
not in appropriate units.
Returns: ratio – The ratio of population of ions in ionization state j to state k.
Return type: Examples
>>> import astropy.units as u >>> T_e = 5000 * u.K >>> n = 1e19 * u.m ** -3 >>> g_j = 2 >>> g_k = 2 >>> E_jk = 1 * u.Ry >>> Saha(g_j, g_k, n, E_jk, T_e) <Quantity 3.299...e-06> >>> T_e = 1 * u.Ry >>> n = 1e23 * u.m ** -3 >>> Saha(g_j, g_k, n, E_jk, T_e) <Quantity 1114595.586...>
Notes
For reference to this function and for more information regarding the Saha equation, see chapter 3 of R Paul Drake’s book, “High-Energy-Density Physics: Foundation of Inertial Fusion and Experimental Astrophysics” (DOI: 10.1007/978-3-319-67711-8_3).