coupling_parameter
- plasmapy.formulary.collisions.dimensionless.coupling_parameter(T: Unit("K"), n_e: Unit("1 / m3"), species, z_mean: ~numbers.Real = nan, V: Unit("m / s") = <Quantity nan m / s>, method='classical') -> Unit(dimensionless)
Ratio of the Coulomb energy to the kinetic (usually thermal) energy.
Classical plasmas are weakly coupled (\(Γ ≪ 1\), where \(Γ\) is the coupling parameter). Dense plasmas tend to have significant to strong coupling (\(Γ ≥ 1\)). For more details, see the notes section below.
- Parameters:
T (
Quantity
) – Temperature in units of temperature or energy per particle, which is assumed to be equal for both the test particle and the target particle.n_e (
Quantity
) – The electron number density in units convertible to per cubic meter.species (
tuple
) – A tuple containing string representations of the test particle (listed first) and the target particle (listed second).z_mean (
Quantity
, optional) – The average ionization (arithmetic mean) of a plasma for which a macroscopic description is valid. This parameter is used to compute the average ion density (given the average ionization and electron density) for calculating the ion sphere radius for non-classical impact parameters.z_mean
is a required parameter ifmethod
is"ls_full_interp"
,"hls_max_interp"
, or"hls_full_interp"
.V (
Quantity
, optional) – The relative velocity between particles. If not provided, thermal velocity is assumed: \(μ V^2 \sim 2 k_B T\) where \(μ\) is the reduced mass.method (
str
, optional) – The method by which to compute the coupling parameter: either"classical"
or"quantum"
. The default method is"classical"
. The Notes section of this docstring has more information about these two methods.
- Returns:
coupling – The coupling parameter for a plasma.
- Return type:
- Raises:
ValueError – If the mass or charge of either particle cannot be found, or any of the inputs contain incorrect values.
UnitConversionError – If the units on any of the inputs are incorrect.
RelativityError – If the input velocity is same or greater than the speed of light.
- Warns:
UnitsWarning
– If units are not provided, SI units are assumed.RelativityWarning
– If the input velocity is greater than 5% of the speed of light.
Notes
The coupling parameter is given by
\[Γ = \frac{E_{Coulomb}}{E_{Kinetic}}\]The Coulomb energy is given by
\[E_{Coulomb} = \frac{Z_1 Z_2 q_e^2}{4 π ε_0 r}\]where \(r\) is the Wigner-Seitz radius, and 1 and 2 refer to particle species 1 and 2 between which we want to determine the coupling.
In the classical case the kinetic energy is the thermal energy:
\[E_{kinetic} = k_B T_e\]The quantum case is more complex. The kinetic energy is dominated by the Fermi energy, modulated by a correction factor based on the ideal chemical potential. This is obtained more precisely by taking the thermal kinetic energy and dividing by the degeneracy parameter, modulated by the Fermi integral [Gericke et al., 2002]:
\[E_{kinetic} = 2 k_B T_e / χ f_{3/2} (μ_{ideal} / k_B T_e)\]where \(χ\) is the degeneracy parameter, \(f_{3/2}\) is the Fermi integral, and \(μ_{ideal}\) is the ideal chemical potential.
The degeneracy parameter is given by
\[χ = n_e Λ_{de Broglie} ^ 3\]where \(n_e\) is the electron density and \(Λ_{de Broglie}\) is the thermal de Broglie wavelength.
See equations 1.2, 1.3 and footnote 5 in Bonitz [1998] for details on the ideal chemical potential.
Examples
>>> import astropy.units as u >>> n = 1e19 * u.m**-3 >>> T = 1e6 * u.K >>> species = ('e', 'p') >>> coupling_parameter(T, n, species) <Quantity 5.8033...e-05> >>> coupling_parameter(T, n, species, V=1e6 * u.m / u.s) <Quantity 5.8033...e-05>