gyroradius

plasmapy.formulary.parameters.gyroradius(B: Unit("T"), particle: plasmapy.particles.particle_class.Particle, *, Vperp: Unit("m / s") = <Quantity nan m / s>, T_i: Unit("K") = <Quantity nan K>) -> Unit("m")

Return the particle gyroradius.

Aliases: rc_, rhoc_

Parameters
  • B (Quantity) – The magnetic field magnitude in units convertible to tesla.

  • particle (Particle) – Representation of the particle species (e.g., 'p' for protons, 'D+' for deuterium, or 'He-4 +1' for singly ionized helium-4). If no charge state information is provided, then the particles are assumed to be singly charged.

  • Vperp (Quantity, optional, keyword-only) – The component of particle velocity that is perpendicular to the magnetic field in units convertible to meters per second.

  • T_i (Quantity, optional, keyword-only) – The particle temperature in units convertible to kelvin.

Returns

r_Li – The particle gyroradius in units of meters. This Quantity will be based on either the perpendicular component of particle velocity as inputted, or the most probable speed for an particle within a Maxwellian distribution for the particle temperature.

Return type

Quantity

Raises
Warns

UnitsWarning – If units are not provided, SI units are assumed.

Notes

One but not both of Vperp and T_i must be inputted.

If any of B, Vperp, or T_i is a number rather than a Quantity, then SI units will be assumed and a warning will be raised.

The particle gyroradius is also known as the particle Larmor radius and is given by

\[r_{Li} = \frac{V_{\perp}}{ω_{ci}}\]

where \(V_⟂\) is the component of particle velocity that is perpendicular to the magnetic field and \(ω_{ci}\) is the particle gyrofrequency. If a temperature is provided, then \(V_⟂\) will be the most probable thermal velocity of an particle at that temperature.

Examples

>>> from astropy import units as u
>>> gyroradius(0.2*u.T, particle='p+', T_i=1e5*u.K)
<Quantity 0.002120... m>
>>> gyroradius(0.2*u.T, particle='p+', T_i=1e5*u.K)
<Quantity 0.002120... m>
>>> gyroradius(5*u.uG, particle='alpha', T_i=1*u.eV)
<Quantity 288002.38... m>
>>> gyroradius(400*u.G, particle='Fe+++', Vperp=1e7*u.m/u.s)
<Quantity 48.23129... m>
>>> gyroradius(B=0.01*u.T, particle='e-', T_i=1e6*u.K)
<Quantity 0.003130... m>
>>> gyroradius(0.01*u.T, 'e-', Vperp=1e6*u.m/u.s)
<Quantity 0.000568... m>
>>> gyroradius(0.2*u.T, 'e-', T_i=1e5*u.K)
<Quantity 4.94949...e-05 m>
>>> gyroradius(5*u.uG, 'e-', T_i=1*u.eV)
<Quantity 6744.25... m>
>>> gyroradius(400*u.G, 'e-', Vperp=1e7*u.m/u.s)
<Quantity 0.001421... m>