thermal_speed
- plasmapy.formulary.speeds.thermal_speed(
- T: Quantity,
- particle: str | int | integer | Particle | CustomParticle | Quantity,
- method: Literal['most_probable', 'rms', 'mean_magnitude', 'nrl'] = 'most_probable',
- mass: Quantity = None,
- ndim: int = 3,
Calculate the speed of thermal motion for particles with a Maxwellian distribution. (See the Notes section for details.).
Aliases:
vth_
Lite Version:
thermal_speed_lite
- Parameters:
T (
Quantity
) – The temperature of the particle distribution, in units of kelvin or energy.particle (
Particle
) – Representation of the particle species (e.g.,"p"
for protons,"D+"
for deuterium, or"He-4 +1"
for singly ionized helium-4).method (
str
, optional) – (Default"most_probable"
) Method to be used for calculating the thermal speed. Valid values are"most_probable"
,"rms"
,"mean_magnitude"
, and"nrl"
.mass (
Quantity
) – Mass override in units convertible to kg. If given, thenmass
will be used instead of the mass value associated withparticle
.ndim (
int
) – (Default3
) Dimensionality (1D, 2D, 3D) of space in which to calculate thermal speed. Valid values are1
,2
, or3
.
- Returns:
vth – Thermal speed of the Maxwellian distribution.
- Return type:
- Raises:
UnitConversionError – If the particle temperature is not in units of temperature or energy per particle.
ValueError – The particle temperature is invalid or particle cannot be used to identify an isotope or particle.
- Warns:
RelativityWarning
– If the ion sound speed exceeds 5% of the speed of light.UnitsWarning
– If units are not provided, SI units are assumed.
Notes
There are multiple methods (or definitions) for calculating the thermal speed, all of which give the expression
\[v_{th} = C_o \sqrt{\frac{k_B T}{m}}\]where \(T\) is the temperature associated with the distribution, \(m\) is the particle’s mass, and \(C_o\) is a constant of proportionality determined by the method in which \(v_{th}\) is calculated and the dimensionality of the system (1D, 2D, 3D). The \(C_o\) used for the
thermal_speed
calculation is determined from the input argumentsmethod
andndim
, and the values can be seen in the table below: ↓ method
ndim →
1
2
3
"most_probable"
\[0\]\[1\]\[\sqrt{2}\]"rms"
\[1\]\[\sqrt{2}\]\[\sqrt{3}\]"mean_magnitude"
\[\sqrt{2/π}\]\[\sqrt{π/2}\]\[\sqrt{8/π}\]"nrl"
\[1\]The coefficients can be directly retrieved using
thermal_speed_coefficients
.The Methods
In the following discussion the Maxwellian distribution \(f(\mathbf{v})\) is assumed to be 3D, but similar expressions can be given for 1D and 2D.
Most Probable
method = "most_probable"
This method expresses the thermal speed of the distribution by expressing it as the most probable speed a particle in the distribution may have. To do this we first define another function \(g(v)\) given by
\[\int_{0}^{∞} g(v) dv = \int_{-∞}^{∞} f(\mathbf{v}) d^3\mathbf{v} \quad \rightarrow \quad g(v) = 4 π v^2 f(v)\]then
\[\begin{split}g^{\prime}(v_{th}) = \left.\frac{dg}{dv}\right|_{v_{th}} = 0\\ \implies v_{th} = \sqrt{\frac{2 k_B T}{m}}\end{split}\]Root Mean Square
method = "rms"
This method uses the root mean square to calculate an expression for the thermal speed of the particle distribution, which is given by
\[v_{th} = \left[\int v^2 f(\mathbf{v}) d^3 \mathbf{v}\right]^{1/2} = \sqrt{\frac{3 k_B T}{m}}\]Mean Magnitude
method = "mean_magnitude"
This method uses the mean speed of the particle distribution to calculate an expression for the thermal speed, which is given by
\[v_{th} = \int |\mathbf{v}| f(\mathbf{v}) d^3 \mathbf{v} = \sqrt{\frac{8 k_B T}{π m}}\]NRL Formulary
method = "nrl"
The NRL Plasma Formulary [Richardson, 2019] uses the square root of the Normal distribution’s variance as the expression for thermal speed.
\[v_{th} = σ = \sqrt{\frac{k_B T}{m}} \quad \text{where} \quad f(v) \sim e^{v^2 / 2 σ^2}\]
Examples
>>> import astropy.units as u >>> thermal_speed(5*u.eV, 'p+') <Quantity 30949.6... m / s> >>> thermal_speed(1e6*u.K, particle='p+') <Quantity 128486... m / s> >>> thermal_speed(5*u.eV, particle='e-') <Quantity 132620... m / s> >>> thermal_speed(1e6*u.K, particle='e-') <Quantity 550569... m / s> >>> thermal_speed(1e6*u.K, "e-", method="rms") <Quantity 674307... m / s> >>> thermal_speed(1e6*u.K, "e-", method="mean_magnitude") <Quantity 621251... m / s>
For user convenience
thermal_speed_coefficients
andthermal_speed_lite
are bound to this function and can be used as follows.>>> from plasmapy.particles import Particle >>> mass = Particle("p").mass.value >>> coeff = thermal_speed.coefficients(method="most_probable", ndim=3) >>> thermal_speed.lite(T=1e6, mass=mass, coeff=coeff) 128486...