Maxwellian_1D
- plasmapy.formulary.distribution.Maxwellian_1D(v, T, particle: str | Integral | Particle | CustomParticle | Quantity = 'e', v_drift=0, vTh=nan, units='units', *, mass_numb=None, Z=None)
Probability distribution function of velocity for a Maxwellian distribution in 1D.
Returns the probability density function at the velocity
v
in m/s to find a particleparticle
in a plasma of temperatureT
following the Maxwellian distribution function.- Parameters:
v (
Quantity
) – The velocity in units convertible to m/s.T (
Quantity
) – The temperature in kelvin.particle (
str
, optional) – Representation of the particle species(e.g.,'p'
for protons,'D+'
for deuterium, or'He-4 +1'
for singly ionized helium-4), which defaults to electrons.v_drift (
Quantity
, optional) – The drift velocity in units convertible to m/s.vTh (
Quantity
, optional) – Thermal velocity (most probable velocity) in m/s. This is used for optimization purposes to avoid re-calculatingvTh
, for example when integrating over velocity-space.units (
str
, optional) – Selects whether to run function with units and unit checks (when equal to “units”) or to run as unitless (when equal to “unitless”). The unitless version is substantially faster for intensive computations.mass_numb (integer, keyword-only, optional) – The mass number associated with
particle
.Z (real number, keyword-only, optional) – The charge number associated with
particle
.
- Returns:
f – Probability density in units of velocity-1, normalized so that \(\int_{-∞}^{+∞} f(v) dv = 1\).
- Return type:
- Raises:
TypeError – The parameter arguments are not Quantities and cannot be converted into Quantities.
UnitConversionError – If the parameters are not in appropriate units.
ValueError – If the temperature is negative, or the particle mass or charge state cannot be found.
Notes
In one dimension, the Maxwellian distribution function for a particle of mass m, velocity v, a drift velocity V and with temperature T is:
\[f = \sqrt{\frac{m}{2 \pi k_B T}} e^{-\frac{m}{2 k_B T} (v-V)^2} \equiv \frac{1}{\sqrt{\pi v_{Th}^2}} e^{-(v - v_{drift})^2 / v_{Th}^2}\]where \(v_{Th} = \sqrt{2 k_B T / m}\) is the thermal speed
Examples
>>> from astropy import units as u >>> v=1*u.m/u.s >>> Maxwellian_1D(v=v, T=30000 * u.K, particle='e', v_drift=0 * u.m / u.s) <Quantity 5.9163...e-07 s / m>