Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Use more generalized keyword argument for gyroradius #1210

Merged
merged 7 commits into from Jul 25, 2021
Merged
Show file tree
Hide file tree
Changes from 5 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
2 changes: 2 additions & 0 deletions changelog/1210.removal.rst
@@ -0,0 +1,2 @@
Use more generalized keyword argument T instead of T_i in gyroradius functionality
Modified gyroradius tests i.e. replaced all T_i in parameters to T
namurphy marked this conversation as resolved.
Show resolved Hide resolved
82 changes: 49 additions & 33 deletions plasmapy/formulary/parameters.py
Expand Up @@ -60,7 +60,7 @@
check_relativistic,
validate_quantities,
)
from plasmapy.utils.exceptions import PhysicsWarning
from plasmapy.utils.exceptions import PhysicsWarning, PlasmaPyFutureWarning

__all__ += __aliases__

Expand Down Expand Up @@ -1094,6 +1094,7 @@ def gyrofrequency(B: u.T, particle: Particle, signed=False, Z=None) -> u.rad / u
@validate_quantities(
Vperp={"can_be_nan": True},
T_i={"can_be_nan": True, "equivalencies": u.temperature_energy()},
T={"can_be_nan": True, "equivalencies": u.temperature_energy()},
validations_on_return={"equivalencies": u.dimensionless_angles()},
)
def gyroradius(
Expand All @@ -1102,6 +1103,7 @@ def gyroradius(
*,
Vperp: u.m / u.s = np.nan * u.m / u.s,
T_i: u.K = np.nan * u.K,
T: u.K = np.nan * u.K,
) -> u.m:
r"""Return the particle gyroradius.

Expand All @@ -1122,8 +1124,12 @@ def gyroradius(
The component of particle velocity that is perpendicular to the
magnetic field in units convertible to meters per second.

T : `~astropy.units.Quantity`, optional, keyword-only
The particle temperature in units convertible to kelvin.

T_i : `~astropy.units.Quantity`, optional, keyword-only
The particle temperature in units convertible to kelvin.
Note: Deprecated. Use T instead.

Returns
-------
Expand Down Expand Up @@ -1152,9 +1158,9 @@ def gyroradius(

Notes
-----
One but not both of ``Vperp`` and ``T_i`` must be inputted.
One but not both of ``Vperp`` and ``T`` must be inputted.

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

Expand All @@ -1173,69 +1179,79 @@ def gyroradius(
Examples
--------
>>> from astropy import units as u
>>> gyroradius(0.2*u.T, particle='p+', T_i=1e5*u.K)
>>> gyroradius(0.2*u.T, particle='p+', T=1e5*u.K)
<Quantity 0.002120... m>
>>> gyroradius(0.2*u.T, particle='p+', T_i=1e5*u.K)
>>> gyroradius(0.2*u.T, particle='p+', T=1e5*u.K)
<Quantity 0.002120... m>
>>> gyroradius(5*u.uG, particle='alpha', T_i=1*u.eV)
>>> gyroradius(5*u.uG, particle='alpha', T=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)
>>> gyroradius(B=0.01*u.T, particle='e-', T=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)
>>> gyroradius(0.2*u.T, 'e-', T=1e5*u.K)
<Quantity 4.94949...e-05 m>
>>> gyroradius(5*u.uG, 'e-', T_i=1*u.eV)
>>> gyroradius(5*u.uG, 'e-', T=1*u.eV)
<Quantity 6744.25... m>
>>> gyroradius(400*u.G, 'e-', Vperp=1e7*u.m/u.s)
<Quantity 0.001421... m>
"""

isfinite_Ti = np.isfinite(T_i)
# Backwards Compatibility and Deprecation check for keyword T_i
if not np.isnan(T_i):
warnings.warn(
"keyword T_i will is deprecated, use T instead",
namurphy marked this conversation as resolved.
Show resolved Hide resolved
PlasmaPyFutureWarning,
)
# Assign T_i to T, if T is not assigned already else use T value instead
namurphy marked this conversation as resolved.
Show resolved Hide resolved
if np.isnan(T):
T = T_i

isfinite_T = np.isfinite(T)
isfinite_Vperp = np.isfinite(Vperp)

# check 1: ensure either Vperp or T_i invalid, keeping in mind that
# check 1: ensure either Vperp or T invalid, keeping in mind that
# the underlying values of the astropy quantity may be numpy arrays
if np.any(np.logical_and(isfinite_Vperp, isfinite_Ti)):
if np.any(np.logical_and(isfinite_Vperp, isfinite_T)):
raise ValueError(
"Must give Vperp or T_i, but not both, as arguments to gyroradius"
"Must give Vperp or T, but not both, as arguments to gyroradius"
)

# check 2: get Vperp as the thermal speed if is not already a valid input
if np.isscalar(Vperp.value) and np.isscalar(
T_i.value
): # both T_i and Vperp are scalars
T.value
): # both T and Vperp are scalars
# we know exactly one of them is nan from check 1
if isfinite_Ti:
# T_i is valid, so use it to determine Vperp
Vperp = thermal_speed(T_i, particle=particle)
if isfinite_T:
# T is valid, so use it to determine Vperp
Vperp = thermal_speed(T, particle=particle)
# else: Vperp is already valid, do nothing
elif np.isscalar(Vperp.value): # only T_i is an array
# this means either Vperp must be nan, or T_i must be array of all nan,
elif np.isscalar(Vperp.value): # only T is an array
# this means either Vperp must be nan, or T must be array of all nan,
# or else we couldn't have gotten through check 1
if isfinite_Vperp:
# Vperp is valid, T_i is a vector that is all nan
# Vperp is valid, T is a vector that is all nan
# uh...
Vperp = np.repeat(Vperp, len(T_i))
Vperp = np.repeat(Vperp, len(T))
else:
# normal case where Vperp is scalar nan and T_i is valid array
Vperp = thermal_speed(T_i, particle=particle)
elif np.isscalar(T_i.value): # only Vperp is an array
# this means either T_i must be nan, or V_perp must be array of all nan,
# normal case where Vperp is scalar nan and T is valid array
Vperp = thermal_speed(T, particle=particle)
elif np.isscalar(T.value): # only Vperp is an array
# this means either T must be nan, or V_perp must be array of all nan,
# or else we couldn't have gotten through check 1
if isfinite_Ti:
# T_i is valid, V_perp is an array of all nan
if isfinite_T:
# T is valid, V_perp is an array of all nan
# uh...
Vperp = thermal_speed(np.repeat(T_i, len(Vperp)), particle=particle)
# else: normal case where T_i is scalar nan and Vperp is already a valid array
Vperp = thermal_speed(np.repeat(T, len(Vperp)), particle=particle)
# else: normal case where T is scalar nan and Vperp is already a valid array
# so, do nothing
else: # both T_i and Vperp are arrays
else: # both T and Vperp are arrays
# we know all the elementwise combinations have one nan and one finite, due to check 1
# use the valid Vperps, and replace the others with those calculated from T_i
# use the valid Vperps, and replace the others with those calculated from T
Vperp = Vperp.copy() # avoid changing Vperp's value outside function
Vperp[isfinite_Ti] = thermal_speed(T_i[isfinite_Ti], particle=particle)
Vperp[isfinite_T] = thermal_speed(T[isfinite_T], particle=particle)

omega_ci = gyrofrequency(B, particle)

Expand Down
77 changes: 41 additions & 36 deletions plasmapy/formulary/tests/test_parameters.py
Expand Up @@ -52,6 +52,7 @@
from plasmapy.utils.exceptions import (
PhysicsError,
PhysicsWarning,
PlasmaPyFutureWarning,
RelativityError,
RelativityWarning,
)
Expand Down Expand Up @@ -798,13 +799,13 @@ def test_gyrofrequency():
def test_gyroradius():
r"""Test the gyroradius function in parameters.py."""

assert gyroradius(B, "e-", T_i=T_e).unit.is_equivalent(u.m)
assert gyroradius(B, "e-", T=T_e).unit.is_equivalent(u.m)

assert gyroradius(B, "e-", Vperp=25 * u.m / u.s).unit.is_equivalent(u.m)

# test for possiblity to allow nan for input values
assert np.isnan(gyroradius(np.nan * u.T, particle="e-", T_i=1 * u.K))
assert np.isnan(gyroradius(1 * u.T, particle="e-", T_i=np.nan * u.K))
assert np.isnan(gyroradius(np.nan * u.T, particle="e-", T=1 * u.K))
assert np.isnan(gyroradius(1 * u.T, particle="e-", T=np.nan * u.K))
assert np.isnan(gyroradius(1 * u.T, particle="e-", Vperp=np.nan * u.m / u.s))

Vperp = 1e6 * u.m / u.s
Expand All @@ -830,32 +831,36 @@ def test_gyroradius():
assert np.isnan(gyroradius(np.nan * u.T, "e-", Vperp=1 * u.m / u.s))

with pytest.raises(ValueError):
gyroradius(3.14159 * u.T, "e-", T_i=-1 * u.K)
gyroradius(3.14159 * u.T, "e-", T=-1 * u.K)

with pytest.warns(u.UnitsWarning):
assert gyroradius(1.0, "e-", Vperp=1.0) == gyroradius(
1.0 * u.T, "e-", Vperp=1.0 * u.m / u.s
)

with pytest.warns(u.UnitsWarning):
assert gyroradius(1.1, "e-", T_i=1.2) == gyroradius(
1.1 * u.T, "e-", T_i=1.2 * u.K
assert gyroradius(1.1, "e-", T=1.2) == gyroradius(
1.1 * u.T, "e-", T=1.2 * u.K
)

with pytest.raises(ValueError):
gyroradius(1.1 * u.T, "e-", Vperp=1 * u.m / u.s, T_i=1.2 * u.K)
gyroradius(1.1 * u.T, "e-", Vperp=1 * u.m / u.s, T=1.2 * u.K)

with pytest.raises(u.UnitTypeError):
gyroradius(1.1 * u.T, "e-", Vperp=1.1 * u.m, T_i=1.2 * u.K)
gyroradius(1.1 * u.T, "e-", Vperp=1.1 * u.m, T=1.2 * u.K)

assert gyroradius(B, particle="p", T_i=T_i).unit.is_equivalent(u.m)
# Check for Deprecation warning when using T_i instead of T
with pytest.warns(PlasmaPyFutureWarning):
gyroradius(1.1 * u.T, "e-", T_i=1.2 * u.K)

assert gyroradius(B, particle="p", T=T_i).unit.is_equivalent(u.m)

assert gyroradius(B, particle="p", Vperp=25 * u.m / u.s).unit.is_equivalent(u.m)

# Case when Z=1 is assumed
assert np.isclose(
gyroradius(B, particle="p", T_i=T_i),
gyroradius(B, particle="H+", T_i=T_i),
gyroradius(B, particle="p", T=T_i),
gyroradius(B, particle="H+", T=T_i),
atol=1e-6 * u.m,
)

Expand All @@ -876,82 +881,82 @@ def test_gyroradius():
particle2 = "alpha"
Vperp2 = thermal_speed(T2, particle=particle2)
gyro_by_vperp = gyroradius(B2, particle="alpha", Vperp=Vperp2)
assert gyro_by_vperp == gyroradius(B2, particle="alpha", T_i=T2)
assert gyro_by_vperp == gyroradius(B2, particle="alpha", T=T2)

explicit_positron_gyro = gyroradius(1 * u.T, particle="positron", T_i=1 * u.MK)
assert explicit_positron_gyro == gyroradius(1 * u.T, "e-", T_i=1 * u.MK)
explicit_positron_gyro = gyroradius(1 * u.T, particle="positron", T=1 * u.MK)
assert explicit_positron_gyro == gyroradius(1 * u.T, "e-", T=1 * u.MK)

with pytest.raises(TypeError):
gyroradius(u.T, particle="p", Vperp=8 * u.m / u.s)

with pytest.raises(ValueError):
gyroradius(B, particle="p", T_i=-1 * u.K)
gyroradius(B, particle="p", T=-1 * u.K)

with pytest.warns(u.UnitsWarning):
gyro_without_units = gyroradius(1.0, particle="p", Vperp=1.0)
gyro_with_units = gyroradius(1.0 * u.T, particle="p", Vperp=1.0 * u.m / u.s)
assert gyro_without_units == gyro_with_units

with pytest.warns(u.UnitsWarning):
gyro_t_without_units = gyroradius(1.1, particle="p", T_i=1.2)
gyro_t_with_units = gyroradius(1.1 * u.T, particle="p", T_i=1.2 * u.K)
gyro_t_without_units = gyroradius(1.1, particle="p", T=1.2)
gyro_t_with_units = gyroradius(1.1 * u.T, particle="p", T=1.2 * u.K)
assert gyro_t_with_units == gyro_t_without_units

with pytest.raises(ValueError):
gyroradius(1.1 * u.T, particle="p", Vperp=1 * u.m / u.s, T_i=1.2 * u.K)
gyroradius(1.1 * u.T, particle="p", Vperp=1 * u.m / u.s, T=1.2 * u.K)

with pytest.raises(u.UnitTypeError):
gyroradius(1.1 * u.T, particle="p", Vperp=1.1 * u.m, T_i=1.2 * u.K)
gyroradius(1.1 * u.T, particle="p", Vperp=1.1 * u.m, T=1.2 * u.K)

with pytest.raises(u.UnitTypeError):
gyroradius(1.1 * u.T, particle="p", Vperp=1.2 * u.m, T_i=1.1 * u.K)
gyroradius(1.1 * u.T, particle="p", Vperp=1.2 * u.m, T=1.1 * u.K)


class Test_gyroradius:

# some custom numpy array tests here, because of the T_i / Vperp situation
# some custom numpy array tests here, because of the T / Vperp situation
def test_handle_numpy_array(self):
# Tests to verify that can handle Quantities with numpy array as the value:
assert gyroradius(B_arr, "e-", Vperp=V_arr)[0] == gyroradius(
B_arr[0], "e-", Vperp=V_arr[0]
)
assert gyroradius(B_arr, "e-", T_i=T_arr)[0] == gyroradius(
B_arr[0], "e-", T_i=T_arr[0]
assert gyroradius(B_arr, "e-", T=T_arr)[0] == gyroradius(
B_arr[0], "e-", T=T_arr[0]
)

def test_handle_mixed_Qarrays(self):
# If both Vperp or Ti are input as Qarrays, but only one of the two is valid
# If both Vperp or T are input as Qarrays, but only one of the two is valid
# at each element, then that's fine, the function should work:
assert gyroradius(B_arr, "e-", Vperp=V_nanarr, T_i=T_nanarr2)[0] == gyroradius(
B_arr[0], "e-", Vperp=V_nanarr[0], T_i=T_nanarr2[0]
assert gyroradius(B_arr, "e-", Vperp=V_nanarr, T=T_nanarr2)[0] == gyroradius(
B_arr[0], "e-", Vperp=V_nanarr[0], T=T_nanarr2[0]
)

def test_raise_two_valid_inputs(self):
# If both Vperp or Ti are nan-less, Qarrays or not, should raise ValueError:
# If both Vperp or T are nan-less, Qarrays or not, should raise ValueError:
with pytest.raises(ValueError):
gyroradius(B_arr, "e-", Vperp=V, T_i=T_arr)
gyroradius(B_arr, "e-", Vperp=V, T=T_arr)
with pytest.raises(ValueError):
gyroradius(B_arr, "e-", Vperp=V_arr, T_i=T_i)
gyroradius(B_arr, "e-", Vperp=V_arr, T=T_i)

def test_all_valid_and_one_valid(self):
# If one of (Vperp, Ti) is a valid and one is Qarray with at least one valid, ValueError:
# If one of (Vperp, T) is a valid and one is Qarray with at least one valid, ValueError:
with pytest.raises(ValueError):
gyroradius(B_arr, "e-", Vperp=V, T_i=T_nanarr)
gyroradius(B_arr, "e-", Vperp=V, T=T_nanarr)
with pytest.raises(ValueError):
gyroradius(B_arr, "e-", Vperp=V_nanarr, T_i=T_i)
gyroradius(B_arr, "e-", Vperp=V_nanarr, T=T_i)

def test_scalar_and_nan_qarray(self):
# If either Vperp or Ti is a valid scalar and the other is a Qarray of all nans,
# If either Vperp or T is a valid scalar and the other is a Qarray of all nans,
# should do something valid and not raise a ValueError
assert np.all(np.isfinite(gyroradius(B_arr, "e-", Vperp=V, T_i=T_allnanarr)))
assert np.all(np.isfinite(gyroradius(B_arr, "e-", Vperp=V_allnanarr, T_i=T_i)))
assert np.all(np.isfinite(gyroradius(B_arr, "e-", Vperp=V, T=T_allnanarr)))
assert np.all(np.isfinite(gyroradius(B_arr, "e-", Vperp=V_allnanarr, T=T_i)))

def test_keeps_arguments_unchanged(self):
Vperp1 = u.Quantity([np.nan, 1], unit=u.m / u.s)
Vperp2 = u.Quantity([np.nan, 1], unit=u.m / u.s) # an exact copy
T_i = u.Quantity([1, np.nan], unit=u.K)

gyroradius(B_arr, "e-", Vperp=Vperp1, T_i=T_i)
gyroradius(B_arr, "e-", Vperp=Vperp1, T=T_i)
assert_quantity_allclose(Vperp1, Vperp2)


Expand Down