gyroradius¶
-
plasmapy.physics.parameters.gyroradius(B, particle='e-', *, Vperp=<Quantity nan m / s>, T_i=<Quantity nan K>)¶ Return the particle gyroradius.
Parameters: - B (Quantity) – The magnetic field magnitude in units convertible to tesla.
- 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. If no charge state information is provided, then the particles are assumed to be singly charged. - Vperp (Quantity, optional) – The component of particle velocity that is perpendicular to the magnetic field in units convertible to meters per second. Must be input as a keyword argument.
- T_i (Quantity, optional) – The particle temperature in units convertible to kelvin. Must be input as a keyword argument.
Returns: r_Li – The particle gyroradius in units of meters. This ~astropy.units.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: Raises: TypeError– The arguments are of an incorrect typeUnitConversionError– The arguments do not have appropriate unitsValueError– If any argument contains invalid values
Warns: ~astropy.units.UnitsWarning – If units are not provided, SI units are assumed
Notes
One but not both of
VperpandT_imust be inputted.If any of
B,Vperp, orT_iis a number rather than aQuantity, 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}}{omega_{ci}}\]where \(V_{\perp}\) is the component of particle velocity that is perpendicular to the magnetic field and \(\omega_{ci}\) is the particle gyrofrequency. If a temperature is provided, then \(V_\perp\) 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.00212087 m> >>> gyroradius(0.2*u.T,particle='p+',T_i=1e5*u.K) <Quantity 0.00212087 m> >>> gyroradius(5*u.uG,particle='alpha',T_i=1*u.eV) <Quantity 288002.38837768 m> >>> gyroradius(400*u.G,particle='Fe+++',Vperp=1e7*u.m/u.s) <Quantity 48.23129811 m> >>> gyroradius(B=0.01*u.T,T_i=1e6*u.K) <Quantity 0.00313033 m> >>> gyroradius(B=0.01*u.T,Vperp=1e6*u.m/u.s) <Quantity 0.00056856 m> >>> gyroradius(0.2*u.T,T_i=1e5*u.K) <Quantity 4.94949252e-05 m> >>> gyroradius(5*u.uG,T_i=1*u.eV) <Quantity 6744.2598183 m> >>> gyroradius(400*u.G,Vperp=1e7*u.m/u.s) <Quantity 0.00142141 m>