Repository: PlasmaPy/PlasmaPy
Pull request #396 (2018-04-23 20:54:01)

Use helper method to fix up drift units
Avoids quite a bit of code repetition in distributions.py


Patch:
diff --git a/plasmapy/classes/plasma.py b/plasmapy/classes/plasma.py
index 614e9d9d..431bbbac 100644
--- a/plasmapy/classes/plasma.py
+++ b/plasmapy/classes/plasma.py
@@ -1,7 +1,6 @@
 """
 Defines the core Plasma class used by PlasmaPy to represent plasma properties.
 """
-import warnings
 
 import numpy as np
 import astropy.units as u
@@ -24,7 +23,7 @@
                              integer_charge,
                              )
 
-from plasmapy.utils import call_string, CouplingWarning
+from plasmapy.utils import call_string
 
 
 class Plasma3D:
@@ -226,15 +225,11 @@ def coupling(self):
         are important. This compares Coulomb potential energy to thermal
         kinetic energy.
         """
-        couple = coupling_parameter(self.T_e,
-                                    self.n_e,
-                                    (self.particle, self.particle),
-                                    self.Z)
-        if couple < 0.01:
-            warnings.warn(f"Coupling parameter is {couple}, you might have strong coupling effects",
-                          CouplingWarning)
-
-        return couple
+        return coupling_parameter(self.T_e,
+                                  self.n_e,
+                                  (self.particle, self.particle),
+                                  self.Z)
+
     def quantum_theta(self):
         """
         Quantum theta parameter, which compares Fermi kinetic energy to
diff --git a/plasmapy/examples/plot_distribution.py b/plasmapy/examples/plot_distribution.py
index 33bd6639..815994c2 100644
--- a/plasmapy/examples/plot_distribution.py
+++ b/plasmapy/examples/plot_distribution.py
@@ -2,8 +2,7 @@
 1D Maxwellian distribution function
 ===================================
 
-We import the usual modules, and the hero of this notebook,
-the Maxwellian 1D distribution:
+We import the usual and the hero of this notebook, the Maxwellian 1D distribution:
 """
 
 
@@ -12,41 +11,43 @@
 import plasmapy
 import matplotlib.pyplot as plt
 from plasmapy.constants import (m_p, m_e, c, mu0, k_B, e, eps0, pi, e)
-
 from plasmapy.physics.distribution import Maxwellian_1D
 
 
 ############################################################
-# Given we'll be plotting, import astropy's quantity support:
+# Given we'll be plotting:
 
 
 from astropy.visualization import quantity_support
 quantity_support()
 
 ############################################################
-# As a first example, let's get the probability density of
-# finding an electron with a speed of 1 m/s if we have a
-# plasma at a temperature of 30 000 K:
+# Let's get the probability density of finding an electron at 1 m/s if we
+# have a plasma at 30 000 K:
 
 
-p_dens = Maxwellian_1D(v=1 * u.m / u.s, T=30000 * u.K,
-                       particle='e', V_drift=0 * u.m / u.s)
-print(p_dens)
+Maxwellian_1D(
+    v=1 * u.m / u.s,
+    T=30000 * u.K,
+    particle='e',
+    V_drift=0 * u.m / u.s)
 
 ############################################################
-# Note the units! Integrated over speed, this will give us a
+# Note the units! Integrated over velocities, this will give us a
 # probability. Let's test that for a bunch of particles:
 
-T = 3e4 * u.K
-dv = 10 * u.m / u.s
-v = np.arange(-5e6, 5e6, 10) * u.m / u.s
 
 ############################################################
-# Check that the integral over all speeds is 1
-# (the particle has to be somewhere):
+# With a vector (from Numpy)
+start = -500000
+stop = -start
+v = np.arange(start, stop) * 10 * u.m / u.s
+dv = v[1] - v[0]
 
+############################################################
+# Test the normalization to 1 (the particle has to be somewhere)
 for particle in ['p', 'e']:
-    pdf = Maxwellian_1D(v, T=T, particle=particle)
+    pdf = Maxwellian_1D(v, T=30000 * u.K, particle=particle)
     integral = (pdf).sum() * dv
     print(f"Integral value for {particle}: {integral}")
     plt.plot(v, pdf, label=particle)
@@ -54,12 +55,29 @@
 
 
 ############################################################
-# The standard deviation of this distribution should give us back the
+# The standard deviation of this distribution should give us back our
 # temperature:
 
 
+T = 30000 * u.K
 std = np.sqrt((Maxwellian_1D(v, T=T, particle='e') * v**2 * dv).sum())
 T_theo = (std**2 / k_B * m_e).to(u.K)
 
-print('T from standard deviation:', T_theo)
-print('Initial T:', T)
+print(T_theo / T)
+
+
+############################################################
+# And the center of the distribution is, as can be seen below:
+
+
+T_e = 30000 * u.K
+V_drift = 10 * u.km / u.s
+
+start = -5000
+stop = - start
+dv = 10000 * u.m / u.s
+
+v_vect = np.arange(start, stop, dtype='float64') * dv
+
+print(v_vect[Maxwellian_1D(v_vect, T=T_e,
+                           particle='e', V_drift=V_drift).argmax()])
diff --git a/plasmapy/physics/distribution.py b/plasmapy/physics/distribution.py
index 96a2a153..87d88e96 100644
--- a/plasmapy/physics/distribution.py
+++ b/plasmapy/physics/distribution.py
@@ -6,6 +6,15 @@
 from scipy.special import gamma
 
 
+def _v_drift_units(v_drift):
+    # Helper method to assign units to  v_drift if it takes a default value
+    if (v_drift == 0 and
+        not isinstance(v_drift, astropy.units.quantity.Quantity)):
+        v_drift = v_drift * u.m / u.s
+    else:
+        v_drift = v_drift.to(u.m / u.s)
+    return v_drift
+
 def Maxwellian_1D(v,
                   T,
                   particle="e",
@@ -13,8 +22,6 @@ def Maxwellian_1D(v,
                   vTh=np.nan,
                   units="units"):
     r"""
-    pdf of velocity for a Maxwellian distribution.
-
     Returns the probability at the velocity `v` in m/s
     to find a particle `particle` in a plasma of temperature `T`
     following the Maxwellian distribution function.
@@ -28,15 +35,16 @@ def Maxwellian_1D(v,
         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 :math:`He_4^{+1}`
-        (singly ionized helium-4), which defaults to electrons.
+        Representation of the particle species(e.g., `'p'` for protons, `'D+'`
+        for deuterium, or `'He-4 +1'` for :math:`He_4^{+1}`
+        (singly ionized helium-4),
+        which defaults to electrons.
 
     V_drift: ~astropy.units.Quantity, optional
         The drift velocity in units convertible to m/s.
 
     vTh: ~astropy.units.Quantity, optional
-        Thermal velocity (most probable velocity) in m/s. This is used for
+        Thermal velocity (most probable) in m/s. This is used for
         optimization purposes to avoid re-calculating vTh, for example
         when integrating over velocity-space.
 
@@ -49,7 +57,7 @@ def Maxwellian_1D(v,
     Returns
     -------
     f : ~astropy.units.Quantity
-        Probability in units of Velocity^-1, normalized so that
+        Probability in Velocity^-1, normalized so that
         :math:`\int_{-\infty}^{+\infty} f(v) dv = 1`.
 
     Raises
@@ -73,9 +81,9 @@ def Maxwellian_1D(v,
     .. math::
 
         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}
+        f = (\pi * v_Th^2)^{-1/2} e^{-(v - v_{drift})^2 / v_Th^2}
 
-    where :math:`v_{Th} = \sqrt{2 k_B T / m}` is the thermal speed
+    where :math:`v_Th = \sqrt(2 k_B T / m)` is the thermal speed
 
     Examples
     --------
@@ -85,18 +93,13 @@ def Maxwellian_1D(v,
     >>> Maxwellian_1D(v=v, T= 30000*u.K, particle='e',V_drift=0*u.m/u.s)
     <Quantity 5.91632969e-07 s / m>
     """
-
     if units == "units":
         # unit checks and conversions
         # checking velocity units
         v = v.to(u.m / u.s)
-        # catching case where drift velocities have default values, they
+        # Catching case where drift velocities have default values, they
         # need to be assigned units
-        if V_drift == 0:
-            if not isinstance(V_drift, astropy.units.quantity.Quantity):
-                V_drift = V_drift * u.m / u.s
-        # checking units of drift velocities
-        V_drift = V_drift.to(u.m / u.s)
+        V_drift = _v_drift_units(V_drift)
         # convert temperature to Kelvins
         T = T.to(u.K, equivalencies=u.temperature_energy())
         if np.isnan(vTh):
@@ -137,10 +140,8 @@ def Maxwellian_velocity_3D(vx,
                            vTh=np.nan,
                            units="units"):
     r"""
-    pdf of velocity for a 3D Maxwellian distribution.
-
     Return the probability of finding a particle with velocity components
-    `vx`, `vy`, and `vz` in m/s in an equilibrium plasma of temperature
+    `vx`, `vy`, and `vz`in m/s in an equilibrium plasma of temperature
     `T` which follows the 3D Maxwellian distribution function. This
     function assumes Cartesian coordinates.
 
@@ -159,8 +160,8 @@ def Maxwellian_velocity_3D(vx,
         The temperature, preferably in Kelvin.
 
     particle: str, optional
-        Representation of the particle species (e.g., ``'p'`` for protons,
-        ``'D+'`` for deuterium, or ``'He-4 +1'`` for :math:`He_4^{+1}`
+        Representation of the particle species (e.g., `'p'` for protons, `'D+'`
+        for deuterium, or `'He-4 +1'` for :math:`He_4^{+1}`
         (singly ionized helium-4), which defaults to electrons.
 
     Vx_drift: ~astropy.units.Quantity, optional
@@ -210,10 +211,9 @@ def Maxwellian_velocity_3D(vx,
 
     .. math::
 
-        f = (\pi v_{Th}^2)^{-3/2} \exp \left [-(\vec{v} -
-        \vec{V}_{drift})^2 / v_{Th}^2 \right ]
+        f = (\pi * v_Th^2)^{-3/2} \exp(-(\vec{v} - \vec{V_{drift}})^2 / v_Th^2)
 
-    where :math:`v_{Th} = \sqrt{2 k_B T / m}` is the thermal speed.
+    where :math:`v_Th = \sqrt(2 k_B T / m)` is the thermal speed.
 
     See also
     --------
@@ -244,19 +244,9 @@ def Maxwellian_velocity_3D(vx,
         vz = vz.to(u.m / u.s)
         # catching case where drift velocities have default values, they
         # need to be assigned units
-        if Vx_drift == 0:
-            if not isinstance(Vx_drift, astropy.units.quantity.Quantity):
-                Vx_drift = Vx_drift * u.m / u.s
-        if Vy_drift == 0:
-            if not isinstance(Vy_drift, astropy.units.quantity.Quantity):
-                Vy_drift = Vy_drift * u.m / u.s
-        if Vz_drift == 0:
-            if not isinstance(Vz_drift, astropy.units.quantity.Quantity):
-                Vz_drift = Vz_drift * u.m / u.s
-        # checking units of drift velocities
-        Vx_drift = Vx_drift.to(u.m / u.s)
-        Vy_drift = Vy_drift.to(u.m / u.s)
-        Vz_drift = Vz_drift.to(u.m / u.s)
+        Vx_drift = _v_drift_units(Vx_drift)
+        Vy_drift = _v_drift_units(Vy_drift)
+        Vz_drift = _v_drift_units(Vz_drift)
         # convert temperature to Kelvins
         T = T.to(u.K, equivalencies=u.temperature_energy())
         if np.isnan(vTh):
@@ -293,8 +283,6 @@ def Maxwellian_speed_1D(v,
                         vTh=np.nan,
                         units="units"):
     r"""
-    pdf of speed for a Maxwellian distribution.
-
     Return the probability of finding a particle with speed `v` in m/s
     in an equilibrium plasma of temperature `T` which follows the
     Maxwellian distribution function.
@@ -353,9 +341,9 @@ def Maxwellian_speed_1D(v,
 
     .. math::
 
-       f(v) = 4\pi v^2 \frac{1}{(\pi v_{Th}^2)^{3/2}} \exp(-(v - V_{drift})^2 / v_{Th}^2 )
+       f(v) = 4 \pi v^2 (\pi * v_Th^2)^{-3/2} \exp(-(v - V_{drift})^2 / v_Th^2)
 
-    where :math:`v_{Th} = \sqrt{2 k_B T / m}` is the thermal speed.
+    where :math:`v_Th = \sqrt(2 k_B T / m)` is the thermal speed.
 
     Example
     -------
@@ -370,13 +358,9 @@ def Maxwellian_speed_1D(v,
         # unit checks and conversions
         # checking velocity units
         v = v.to(u.m / u.s)
-        # catching case where drift velocity has default value, and
+        # Catching case where drift velocity has default value, and
         # needs to be assigned units
-        if V_drift == 0:
-            if not isinstance(V_drift, astropy.units.quantity.Quantity):
-                V_drift = V_drift * u.m / u.s
-        # checking drift velocity units
-        V_drift = V_drift.to(u.m / u.s)
+        V_drift = _v_drift_units(V_drift)
         # convert temperature to Kelvins
         T = T.to(u.K, equivalencies=u.temperature_energy())
         if np.isnan(vTh):
@@ -418,10 +402,8 @@ def Maxwellian_speed_3D(vx,
                         vTh=np.nan,
                         units="units"):
     r"""
-    pdf of speed for a 3D Maxwellian distribution.
-
     Return the probability of finding a particle with speed components
-    `vx`, `vy`, and `vz` in m/s in an equilibrium plasma of temperature
+    `vx`, `vy`, and `vz`in m/s in an equilibrium plasma of temperature
     `T` which follows the 3D Maxwellian distribution function. This
     function assumes Cartesian coordinates.
 
@@ -491,10 +473,10 @@ def Maxwellian_speed_3D(vx,
 
     .. math::
 
-       f = 4 \pi \vec{v}^2 (\pi v_{Th}^2)^{-3/2} \exp(-(\vec{v} -
-       \vec{V}_{drift})^2 / v_{Th}^2)
+       f = 4 \pi \vec{v}^2 (\pi * v_Th^2)^{-3/2} \exp(-(\vec{v} -
+       \vec{V_{drift}})^2 / v_Th^2)
 
-    where :math:`v_{Th} = \sqrt{2 k_B T / m}` is the thermal speed.
+    where :math:`v_Th = \sqrt(2 k_B T / m)` is the thermal speed.
 
     See also
     --------
@@ -522,20 +504,11 @@ def Maxwellian_speed_3D(vx,
         vx = vx.to(u.m / u.s)
         vy = vy.to(u.m / u.s)
         vz = vz.to(u.m / u.s)
-        # catching case where drift velocities have default values, they
+        # Catching case where drift velocities have default values, they
         # need to be assigned units
-        if Vx_drift == 0:
-            if not isinstance(Vx_drift, astropy.units.quantity.Quantity):
-                Vx_drift = Vx_drift * u.m / u.s
-        if Vy_drift == 0:
-            if not isinstance(Vy_drift, astropy.units.quantity.Quantity):
-                Vy_drift = Vy_drift * u.m / u.s
-        if Vz_drift == 0:
-            if not isinstance(Vz_drift, astropy.units.quantity.Quantity):
-                Vz_drift = Vz_drift * u.m / u.s
-        Vx_drift = Vx_drift.to(u.m / u.s)
-        Vy_drift = Vy_drift.to(u.m / u.s)
-        Vz_drift = Vz_drift.to(u.m / u.s)
+        Vx_drift = _v_drift_units(Vx_drift)
+        Vy_drift = _v_drift_units(Vy_drift)
+        Vz_drift = _v_drift_units(Vz_drift)
         # convert temperature to Kelvins
         T = T.to(u.K, equivalencies=u.temperature_energy())
         if np.isnan(vTh):
@@ -654,10 +627,10 @@ def kappa_velocity_1D(v,
     .. math::
 
        f = A_\kappa \left(1 + \frac{(\vec{v} -
-       \vec{V_{drift}})^2}{\kappa v_{Th},\kappa^2}\right)^{-\kappa}
+       \vec{V_{drift}})^2}{\kappa v_Th,\kappa^2}\right)^{-\kappa}
 
-    where :math:`v_{Th},\kappa` is the kappa thermal speed
-    and :math:`A_\kappa = \frac{1}{\sqrt{\pi} \kappa^{3/2} v_{Th},\kappa^2
+    where :math:`v_Th,\kappa` is the kappa thermal speed
+    and :math:`A_\kappa = \frac{1}{\sqrt{\pi} \kappa^{3/2} v_Th,\kappa^2
     \frac{\Gamma(\kappa + 1)}{\Gamma(\kappa - 1/2)}}`
     is the normalization constant.
 
@@ -811,10 +784,10 @@ def kappa_velocity_3D(vx,
     .. math::
 
        f = A_\kappa \left(1 + \frac{(\vec{v} -
-       \vec{V_{drift}})^2}{\kappa v_{Th},\kappa^2}\right)^{-(\kappa + 1)}
+       \vec{V_{drift}})^2}{\kappa v_Th,\kappa^2}\right)^{-(\kappa + 1)}
 
-    where :math:`v_{Th},\kappa` is the kappa thermal speed
-    and :math:`A_\kappa = \frac{1}{2 \pi (\kappa v_{Th},\kappa^2)^{3/2}}
+    where :math:`v_Th,\kappa` is the kappa thermal speed
+    and :math:`A_\kappa = \frac{1}{2 \pi (\kappa v_Th,\kappa^2)^{3/2}}
     \frac{\Gamma(\kappa + 1)}{\Gamma(\kappa - 1/2) \Gamma(3/2)}` is the
     normalization constant.
 
@@ -851,21 +824,11 @@ def kappa_velocity_3D(vx,
         vx = vx.to(u.m / u.s)
         vy = vy.to(u.m / u.s)
         vz = vz.to(u.m / u.s)
-        # catching case where drift velocities have default values, they
+        # Catching case where drift velocities have default values, they
         # need to be assigned units
-        if Vx_drift == 0:
-            if not isinstance(Vx_drift, astropy.units.quantity.Quantity):
-                Vx_drift = Vx_drift * u.m / u.s
-        if Vy_drift == 0:
-            if not isinstance(Vy_drift, astropy.units.quantity.Quantity):
-                Vy_drift = Vy_drift * u.m / u.s
-        if Vz_drift == 0:
-            if not isinstance(Vz_drift, astropy.units.quantity.Quantity):
-                Vz_drift = Vz_drift * u.m / u.s
-        # checking units of drift velocities
-        Vx_drift = Vx_drift.to(u.m / u.s)
-        Vy_drift = Vy_drift.to(u.m / u.s)
-        Vz_drift = Vz_drift.to(u.m / u.s)
+        Vx_drift = _v_drift_units(Vx_drift)
+        Vy_drift = _v_drift_units(Vy_drift)
+        Vz_drift = _v_drift_units(Vz_drift)
         # convert temperature to Kelvins
         T = T.to(u.K, equivalencies=u.temperature_energy())
         if np.isnan(vTh):
diff --git a/plasmapy/physics/transport/__init__.py b/plasmapy/physics/transport/__init__.py
index 5cc4550b..de8b6c5e 100644
--- a/plasmapy/physics/transport/__init__.py
+++ b/plasmapy/physics/transport/__init__.py
@@ -2,7 +2,7 @@
                          collision_frequency,
                          collision_rate_electron_ion,
                          collision_rate_ion_ion,
-                         impact_parameter_perp,
+                         b_perp,
                          impact_parameter,
                          mean_free_path,
                          mobility,
diff --git a/plasmapy/physics/transport/collisions.py b/plasmapy/physics/transport/collisions.py
index b986e6ec..82d82b60 100644
--- a/plasmapy/physics/transport/collisions.py
+++ b/plasmapy/physics/transport/collisions.py
@@ -138,45 +138,6 @@ def Coulomb_logarithm(T,
     This means the Coulomb logarithm will not break down for Lambda < 0,
     which occurs for dense, cold plasmas.
 
-    Methods
-    ---
-    Classical
-        classical Landau-Spitzer approach. Fails for large coupling
-        parameter where Lambda can become less than zero.
-    GMS-1
-        1st method listed in Table 1 of reference [3]
-        Landau-Spitzer, but with interpolated bmin instead of bmin
-        selected between deBroglie wavelength and distance of closest
-        approach. Fails for large coupling
-        parameter where Lambda can become less than zero.
-    GMS-2
-        2nd method listed in Table 1 of reference [3]
-        Another Landau-Spitzer like approach, but now bmax is also
-        being interpolated. The interpolation is between the Debye
-        length and the ion sphere radius, allowing for descriptions
-        of dilute plasmas. Fails for large coupling
-        parameter where Lambda can become less than zero.
-        3rd method listed in Table 1 of reference [3]
-        classical Landau-Spitzer fails for argument of Coulomb logarithm
-        Lambda < 0, therefore a clamp is placed at Lambda_min = 2
-    GMS-4
-        4th method listed in Table 1 of reference [3]
-        Spitzer-like extension to Coulomb logarithm by noting that
-        Coulomb collisions take hyperbolic trajectories. Removes
-        divergence for small bmin issue in classical Landau-Spitzer
-        approach, so bmin can be zero. Also doesn't break down as
-        Lambda < 0 is now impossible, even when coupling parameter is large.
-    GMS-5
-        5th method listed in Table 1 of reference [3]
-        Similar to GMS-4, but setting bmin as distance of closest approach
-        and bmax interpolated between Debye length and ion sphere radius.
-        Lambda < 0 impossible.
-    GMS-6
-        6th method listed in Table 1 of reference [3]
-        Similar to GMS-4 and GMS-5, but using interpolation methods
-        for both bmin and bmax.
-    Lambda < 0 impossible.
-
     Examples
     --------
     >>> from astropy import units as u
@@ -211,26 +172,61 @@ def Coulomb_logarithm(T,
                                   z_mean=z_mean,
                                   V=V,
                                   method=method)
-    if method in ["classical", "GMS-1", "GMS-2"]:
+    if method == "classical":
+        # classical Landau-Spitzer approach. Fails for large coupling
+        # parameter where Lambda can become less than zero.
+        ln_Lambda = np.log(bmax / bmin)
+    elif method == "GMS-1":
+        # 1st method listed in Table 1 of reference [3]
+        # Landau-Spitzer, but with interpolated bmin instead of bmin
+        # selected between deBroglie wavelength and distance of closest
+        # approach. Fails for large coupling
+        # parameter where Lambda can become less than zero.
+        ln_Lambda = np.log(bmax / bmin)
+    elif method == "GMS-2":
+        # 2nd method listed in Table 1 of reference [3]
+        # Another Landau-Spitzer like approach, but now bmax is also
+        # being interpolated. The interpolation is between the Debye
+        # length and the ion sphere radius, allowing for descriptions
+        # of dilute plasmas. Fails for large coupling
+        # parameter where Lambda can become less than zero.
         ln_Lambda = np.log(bmax / bmin)
     elif method == "GMS-3":
+        # 3rd method listed in Table 1 of reference [3]
+        # classical Landau-Spitzer fails for argument of Coulomb logarithm
+        # Lambda < 0, therefore a clamp is placed at Lambda_min = 2
         ln_Lambda = np.log(bmax / bmin)
         if ln_Lambda < 2:
             ln_Lambda = 2 * u.dimensionless_unscaled
-    elif method in ["GMS-4", "GMS-5", "GMS-6"]:
+    elif method == "GMS-4":
+        # 4th method listed in Table 1 of reference [3]
+        # Spitzer-like extension to Coulomb logarithm by noting that
+        # Coulomb collisions take hyperbolic trajectories. Removes
+        # divergence for small bmin issue in classical Landau-Spitzer
+        # approach, so bmin can be zero. Also doesn't break down as
+        # Lambda < 0 is now impossible, even when coupling parameter is large.
+        ln_Lambda = 0.5 * np.log(1 + bmax ** 2 / bmin ** 2)
+    elif method == "GMS-5":
+        # 5th method listed in Table 1 of reference [3]
+        # Similar to GMS-4, but setting bmin as distance of closest approach
+        # and bmax interpolated between Debye length and ion sphere radius.
+        # Lambda < 0 impossible.
+        ln_Lambda = 0.5 * np.log(1 + bmax ** 2 / bmin ** 2)
+    elif method == "GMS-6":
+        # 6th method listed in Table 1 of reference [3]
+        # Similar to GMS-4 and GMS-5, but using interpolation methods
+        # for both bmin and bmax.
+        # Lambda < 0 impossible.
         ln_Lambda = 0.5 * np.log(1 + bmax ** 2 / bmin ** 2)
     else:
         raise ValueError("Unknown method! Choose from 'classical' and 'GMS-N', N from 1 to 6.")
     # applying dimensionless units
     ln_Lambda = ln_Lambda.to(u.dimensionless_unscaled).value
-    if ln_Lambda < 0:
-        raise utils.PhysicsError(f"A Coulomb logarithm = {ln_Lambda} < 0 is nonphysical.")
-    elif ln_Lambda < 2 and method in ["classical", "GMS-1", "GMS-2"]:
-        warnings.warn(f"Coulomb logarithm is {ln_Lambda} and {method} relies on weak coupling.",
-                      utils.CouplingWarning)
-    elif ln_Lambda < 4:
+    if ln_Lambda < 4:
         warnings.warn(f"Coulomb logarithm is {ln_Lambda}, you might have strong coupling effects",
-                      utils.CouplingWarning)
+                      utils.PhysicsWarning)
+    # if ln_Lambda < 1:
+    #     raise utils.PhysicsError(f"Coulomb logarithm is {ln_Lambda}, less than 1")
     return ln_Lambda
 
 
@@ -255,8 +251,6 @@ def _boilerPlate(T, particles, V):
     # obtaining reduced mass of 2 particle collision system
     reduced_mass = atomic.reduced_mass(*particles)
 
-    if V == 0:
-        raise utils.exceptions.PhysicsError("You cannot have a collision for zero velocity!")
     # getting thermal velocity of system if no velocity is given
     if np.isnan(V):
         V = parameters.thermal_speed(T, mass=reduced_mass)
@@ -266,9 +260,9 @@ def _boilerPlate(T, particles, V):
 
 @check_quantity({"T": {"units": u.K, "can_be_negative": False}
                  })
-def impact_parameter_perp(T,
-                          particles,
-                          V=np.nan * u.m / u.s):
+def b_perp(T,
+           particles,
+           V=np.nan * u.m / u.s):
     """
     Distance of closest approach for a 90 degree Coulomb collision.
 
@@ -291,7 +285,7 @@ def impact_parameter_perp(T,
 
     Returns
     -------
-    impact_parameter_perp : float or numpy.ndarray
+    b_perp : float or numpy.ndarray
         The distance of closest approach for a 90 degree Coulomb collision.
 
     Raises
@@ -321,7 +315,7 @@ def impact_parameter_perp(T,
 
     Notes
     -----
-    The distance of closest approach, impact_parameter_perp, is given by [1]_
+    The distance of closest approach, b_perp, is given by [1]_
 
     .. math::
         b_{\perp} = \frac{Z_1 Z_2}{4 \pi \epsilon_0 m v^2}
@@ -332,7 +326,7 @@ def impact_parameter_perp(T,
     >>> from astropy import units as u
     >>> T = 1e6*u.K
     >>> particles = ('e', 'p')
-    >>> impact_parameter_perp(T, particles)
+    >>> b_perp(T, particles)
     <Quantity 8.35505011e-12 m>
 
 
@@ -480,9 +474,9 @@ def impact_parameter(T,
     # deBroglie wavelength
     lambdaBroglie = hbar / (2 * reduced_mass * V)
     # distance of closest approach in 90 degree Coulomb collision
-    bPerp = impact_parameter_perp(T=T,
-                                  particles=particles,
-                                  V=V)
+    bPerp = b_perp(T=T,
+                   particles=particles,
+                   V=V)
     # obtaining minimum and maximum impact parameters depending on which
     # method is requested
     if method == "classical":
@@ -633,6 +627,14 @@ def collision_frequency(T,
     taken as the thermal velocity), and :math:`\ln{\Lambda}` is the Coulomb
     logarithm accounting for small angle collisions.
 
+    The collisional cross-section is obtained by
+
+    .. math::
+        \sigma = \pi \rho_{\perp}^2
+
+    where :math:`\rho_{\perp}` is the distance of closest approach for
+    a 90 degree Coulomb collision.
+
     See eq (2.14) in [2]_.
 
     Examples
@@ -664,9 +666,9 @@ def collision_frequency(T,
             V = V_reduced
         # electron-electron collision
         # impact parameter for 90 degree collision
-        bPerp = impact_parameter_perp(T=T,
-                                      particles=particles,
-                                      V=V_reduced)
+        bPerp = b_perp(T=T,
+                       particles=particles,
+                       V=V_reduced)
         # Coulomb logarithm
         cou_log = Coulomb_logarithm(T,
                                     n,
@@ -684,9 +686,9 @@ def collision_frequency(T,
             V = np.sqrt(2 * k_B * T / m_e)
         # need to also correct mass in collision radius from reduced
         # mass to electron mass
-        bPerp = impact_parameter_perp(T=T,
-                                      particles=particles,
-                                      V=V) * reduced_mass / m_e
+        bPerp = b_perp(T=T,
+                       particles=particles,
+                       V=V) * reduced_mass / m_e
         # Coulomb logarithm
         # !!! may also need to correct Coulomb logarithm to be
         # electron-electron version !!!
@@ -702,9 +704,9 @@ def collision_frequency(T,
         if np.isnan(V):
             V = V_reduced
         # ion-ion collision
-        bPerp = impact_parameter_perp(T=T,
-                                      particles=particles,
-                                      V=V)
+        bPerp = b_perp(T=T,
+                       particles=particles,
+                       V=V)
         # Coulomb logarithm
         cou_log = Coulomb_logarithm(T,
                                     n,
@@ -713,7 +715,7 @@ def collision_frequency(T,
                                     V=np.nan * u.m / u.s,
                                     method=method)
     # collisional cross section
-    sigma = Coulomb_cross_section(bPerp)
+    sigma = np.pi * (2 * bPerp) ** 2
     # collision frequency where Coulomb logarithm accounts for
     # small angle collisions, which are more frequent than large
     # angle collisions.
@@ -721,51 +723,6 @@ def collision_frequency(T,
     return freq.to(u.Hz)
 
 
-@check_quantity({
-    'impact_param': {'units': u.m,
-                    'can_be_negative': False}
-    })
-def Coulomb_cross_section(impact_param: u.m):
-    """
-    Cross section for a large angle Coulomb collision
-
-    Parameters
-    ----------
-    impact_param : ~astropy.units.Quantity
-        Impact parameter for the collision.
-
-    Examples
-    --------
-    >>> Coulomb_cross_section(7e-10*u.m)
-    <Quantity 6.1575216e-18 m2>
-    >>> Coulomb_cross_section(0.5*u.m)
-    <Quantity 3.14159265 m2>
-
-    Notes
-    -----
-    The collisional cross-section (see [1]_ for a graphical demonstration)
-    for a 90 degree Coulomb collision is obtained by
-
-    .. math::
-        \sigma = \pi (2 * \rho_{\perp})^2
-
-    where :math:`\rho_{\perp}` is the distance of closest approach for
-    a 90 degree Coulomb collision. This function is a generalization of that
-    calculation. Please note that it is not guaranteed to return the correct
-    results for small angle collisions.
-
-    Returns
-    -------
-    ~astropy.units.Quantity
-        The Coulomb collision cross section area.
-
-    References
-    ----------
-    .. [1] https://en.wikipedia.org/wiki/Cross_section_(physics)#Collision_among_gas_particles
-    """
-    sigma = np.pi * (2 * impact_param) ** 2
-    return sigma
-
 @utils.check_quantity({
     'T_e': {'units': u.K, 'can_be_negative': False},
     'n_e': {'units': u.m ** -3, 'can_be_negative': False}
diff --git a/plasmapy/utils/__init__.py b/plasmapy/utils/__init__.py
index bb47e8b2..d4115d62 100644
--- a/plasmapy/utils/__init__.py
+++ b/plasmapy/utils/__init__.py
@@ -15,7 +15,6 @@
                          InvalidParticleError,
                          PlasmaPyWarning,
                          PhysicsWarning,
-                         CouplingWarning,
                          RelativityWarning,
                          AtomicWarning,
                          MissingAtomicDataWarning)
diff --git a/plasmapy/utils/exceptions.py b/plasmapy/utils/exceptions.py
index cefb4b26..4c08ac80 100644
--- a/plasmapy/utils/exceptions.py
+++ b/plasmapy/utils/exceptions.py
@@ -109,12 +109,6 @@ class RelativityWarning(PhysicsWarning):
     pass
 
 
-class CouplingWarning(PhysicsWarning):
-    """
-    A warning for functions that rely on a particular coupling regime to be valid.
-    """
-
-
 class AtomicWarning(PlasmaPyWarning):
     """The base warning for the `~plasmapy.atomic` subpackage."""
     pass


Test patch:
diff --git a/plasmapy/physics/transport/tests/test_collisions.py b/plasmapy/physics/transport/tests/test_collisions.py
index 53c0721e..5a2ceaa6 100644
--- a/plasmapy/physics/transport/tests/test_collisions.py
+++ b/plasmapy/physics/transport/tests/test_collisions.py
@@ -3,7 +3,7 @@
 from astropy import units as u
 
 from plasmapy.physics.transport import (Coulomb_logarithm,
-                                        impact_parameter_perp,
+                                        b_perp,
                                         impact_parameter,
                                         collision_frequency,
                                         mean_free_path,
@@ -38,12 +38,6 @@ def setup_class(self):
         self.gms6 = 3.635342040477818
         self.gms6_negative = 0.030720859361047514
 
-    def test_symmetry(self):
-        lnLambda = Coulomb_logarithm(self.temperature1, self.density2, self.particles)
-        lnLambdaRev = Coulomb_logarithm(self.temperature1, self.density2, self.particles[::-1])
-        assert lnLambda == lnLambdaRev
-
-
     def test_Chen_Q_machine(self):
         """
         Tests whether Coulomb logarithm gives value consistent with
@@ -66,7 +60,6 @@ def test_Chen_Q_machine(self):
                   f"{lnLambdaChen} and not {lnLambda}.")
         assert testTrue, errStr
 
-
     def test_Chen_lab(self):
         """
         Tests whether Coulomb logarithm gives value consistent with
@@ -183,13 +176,20 @@ def test_GMS1_negative(self):
         Murillo, and Schlanges PRE (2002). This checks for when
         a negative (invalid) Coulomb logarithm is returned.
         """
-        with pytest.raises(exceptions.PhysicsError, match="< 0"):
+        with pytest.warns(exceptions.PhysicsWarning, match="strong coupling effects"):
             methodVal = Coulomb_logarithm(self.temperature2,
                                           self.density2,
                                           self.particles,
                                           z_mean=np.nan * u.dimensionless_unscaled,
                                           V=np.nan * u.m / u.s,
                                           method="GMS-1")
+        testTrue = np.isclose(methodVal,
+                              self.gms1_negative,
+                              rtol=1e-15,
+                              atol=0.0)
+        errStr = (f"Coulomb logarithm for GMS-1 should be "
+                  f"{self.gms1_negative} and not {methodVal}.")
+        assert testTrue, errStr
 
     def test_GMS2(self):
         """
@@ -217,13 +217,20 @@ def test_GMS2_negative(self):
         Murillo, and Schlanges PRE (2002). This checks for when
         a negative (invalid) Coulomb logarithm is returned.
         """
-        with pytest.raises(exceptions.PhysicsError, match="< 0"):
+        with pytest.warns(exceptions.PhysicsWarning, match="strong coupling effects"):
             methodVal = Coulomb_logarithm(self.temperature2,
                                           self.density2,
                                           self.particles,
                                           z_mean=self.z_mean,
                                           V=np.nan * u.m / u.s,
                                           method="GMS-2")
+        testTrue = np.isclose(methodVal,
+                              self.gms2_negative,
+                              rtol=1e-15,
+                              atol=0.0)
+        errStr = (f"Coulomb logarithm for GMS-2 should be "
+                  f"{self.gms2_negative} and not {methodVal}.")
+        assert testTrue, errStr
 
     def test_GMS3(self):
         """
@@ -466,7 +473,7 @@ def test_invalid_particle_error(self):
     particles = ('e', 'p')
 
 
-class Test_impact_parameter_perp:
+class Test_b_perp:
     @classmethod
     def setup_class(self):
         """initializing parameters for tests """
@@ -475,24 +482,19 @@ def setup_class(self):
         self.V = 1e4 * u.km / u.s
         self.True1 = 7.200146594293746e-10
 
-    def test_symmetry(self):
-        result = impact_parameter_perp(self.T, self.particles)
-        resultRev = impact_parameter_perp(self.T, self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
         """
-        methodVal = impact_parameter_perp(self.T,
-                                          self.particles,
-                                          V=np.nan * u.m / u.s)
+        methodVal = b_perp(self.T,
+                           self.particles,
+                           V=np.nan * u.m / u.s)
         testTrue = np.isclose(self.True1,
                               methodVal.si.value,
                               rtol=1e-1,
                               atol=0.0)
         errStr = ("Distance of closest approach for 90 degree Coulomb "
-                  f"collision, impact_parameter_perp, should be {self.True1} and "
+                  f"collision, b_perp, should be {self.True1} and "
                   f"not {methodVal}.")
         assert testTrue, errStr
 
@@ -502,14 +504,14 @@ def test_fail1(self):
         value comparison by some quantity close to numerical error.
         """
         fail1 = self.True1 + 1e-15
-        methodVal = impact_parameter_perp(self.T,
-                                          self.particles,
-                                          V=np.nan * u.m / u.s)
+        methodVal = b_perp(self.T,
+                           self.particles,
+                           V=np.nan * u.m / u.s)
         testTrue = not np.isclose(methodVal.si.value,
                                   fail1,
                                   rtol=1e-16,
                                   atol=0.0)
-        errStr = (f"impact_parameter_perp value test gives {methodVal} and "
+        errStr = (f"b_perp value test gives {methodVal} and "
                   f"should not be equal to {fail1}.")
         assert testTrue, errStr
 
@@ -530,11 +532,6 @@ def setup_class(self):
         self.V = 1e4 * u.km / u.s
         self.True1 = np.array([7.200146594293746e-10, 2.3507660003984624e-08])
 
-    def test_symmetry(self):
-        result = impact_parameter(self.T, self.n_e, self.particles)
-        resultRev = impact_parameter(self.T, self.n_e,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -604,11 +601,6 @@ def setup_class(self):
         self.True_protons = 44450104815.91857
         self.True_zmean = 1346828153985.4646
 
-    def test_symmetry(self):
-        result = collision_frequency(self.T, self.n, self.particles)
-        resultRev = collision_frequency(self.T, self.n, self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -718,11 +710,6 @@ def setup_class(self):
         self.V = 1e4 * u.km / u.s
         self.True1 = 4.4047571877932046e-07
 
-    def test_symmetry(self):
-        result = mean_free_path(self.T, self.n_e, self.particles)
-        resultRev = mean_free_path(self.T, self.n_e,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -776,11 +763,6 @@ def setup_class(self):
         self.True1 = 1.2665402649805445e-3
         self.True_zmean = 0.00020264644239688712
 
-    def test_symmetry(self):
-        result = Spitzer_resistivity(self.T, self.n, self.particles)
-        resultRev = Spitzer_resistivity(self.T, self.n,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -848,11 +830,6 @@ def setup_class(self):
         self.True1 = 0.13066090887074902
         self.True_zmean = 0.32665227217687254
 
-    def test_symmetry(self):
-        result = mobility(self.T, self.n_e, self.particles)
-        resultRev = mobility(self.T, self.n_e,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -923,11 +900,6 @@ def setup_class(self):
         self.V = 1e4 * u.km / u.s
         self.True1 = 440.4757187793204
 
-    def test_symmetry(self):
-        result = Knudsen_number(self.length, self.T, self.n_e, self.particles)
-        resultRev = Knudsen_number(self.length, self.T, self.n_e,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.
@@ -984,11 +956,6 @@ def setup_class(self):
         self.True_zmean = 10.689750083758698
         self.True_quantum = 0.3334662805238162
 
-    def test_symmetry(self):
-        result = coupling_parameter(self.T, self.n_e, self.particles)
-        resultRev = coupling_parameter(self.T, self.n_e,  self.particles[::-1])
-        assert result == resultRev
-
     def test_known1(self):
         """
         Test for known value.