def banana_width(rho: Numeric, q: Numeric, epsilon: Numeric) -> Numeric:
    r"""Radial width of a thermal banana orbit, the conventional estimate.

    $$\Delta_b = \frac{q\rho}{\sqrt\epsilon}$$

    Parameters
    ----------
    rho : float or np.ndarray
        Larmor radius in the toroidal field, positive [m].
    q : float or np.ndarray
        Safety factor, positive [-].
    epsilon : float or np.ndarray
        Inverse aspect ratio $r/R_0$ of the surface, in (0, 1) [-].

    Returns
    -------
    float or np.ndarray
        Banana width [m].

    Raises
    ------
    ValueError
        A non-positive input or ``epsilon`` outside (0, 1).

    Convention
    ----------
    Conservation of $P_\phi$ gives the excursion
    $\Delta r = q\rho\,(v_\parallel/v)/\epsilon$ for the midplane parallel
    speed of the orbit; a thermal trapped particle has
    $v_\parallel/v \sim \sqrt\epsilon$, hence this estimate, with the
    pitch-dependent order-unity factor set to one. The Larmor radius is
    $mv/(|Z|eB_\phi)$ at whatever speed the caller chooses.

    Physical interpretation
    -----------------------
    The step of the banana-regime random walk: larger than the Larmor radius
    by $q/\sqrt\epsilon$, which is why neoclassical transport exceeds
    classical. It is also the finite-orbit width that makes a local
    transport description break down near steep gradients.

    Assumptions
    -----------
    Large aspect ratio, circular surface, banana width small compared with
    the minor radius (not a potato orbit near the axis).

    References
    ----------
    .. [1] J. Wesson, *Tokamaks*, 4th ed., Oxford University Press (2011),
           Sec. 3.12.
    .. [2] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 7.
    """
    rho = _validate_positive("rho", rho)
    q = _validate_positive("q", q)
    epsilon = _validate_epsilon(epsilon)
    return _maybe_scalar(q * rho / np.sqrt(epsilon))

def collisions_per_transit(nu: Numeric, v: Numeric, q: Numeric, R0: Numeric) -> Numeric:
    r"""Collision frequency in units of the transit frequency, $\hat\nu = qR_0\nu/v$.

    $$\hat\nu = \frac{\nu}{\omega_t} = \frac{qR_0\,\nu}{v}$$

    Parameters
    ----------
    nu : float or np.ndarray
        Collision frequency, non-negative; the caller states which one (a
        deflection frequency for the trapped-particle orderings) [1/s].
    v : float or np.ndarray
        Particle speed, positive [m/s].
    q : float or np.ndarray
        Safety factor, positive [-].
    R0 : float or np.ndarray
        Major radius, positive [m].

    Returns
    -------
    float or np.ndarray
        $\hat\nu$: connection length over mean free path, per radian [-].

    Raises
    ------
    ValueError
        A negative ``nu`` or a non-positive ``v``, ``q`` or ``R0``.

    Convention
    ----------
    Not a $\nu_*$. The trapped-particle collisionality is
    $\nu_\mathrm{eff}/\omega_b \sim \hat\nu/\epsilon^{3/2}$, and VAFT's
    existing $\nu_*$ functions (Sauter's among them, #353) each fix their own
    prefactors and thermal-speed conventions; they share the symbol but not
    the value. $\hat\nu$ is kept explicit so its regime boundaries,
    ``neoclassical_regime_boundaries``, carry the $\epsilon$ dependence in
    the open.

    Physical interpretation
    -----------------------
    Below $\epsilon^{3/2}$ a trapped particle completes its banana before it
    is scattered out; above one, a particle collides before it has seen the
    poloidal variation of $B$ at all.

    Assumptions
    -----------
    Large aspect ratio, circular surfaces.

    References
    ----------
    .. [1] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 7.
    """
    nu = _validate_non_negative("nu", nu)
    v = _validate_positive("v", v)
    q = _validate_positive("q", q)
    R0 = _validate_positive("R0", R0)
    return _maybe_scalar(q * R0 * nu / v)

def coulomb_logarithm_electron_sauter(n_e: Numeric, T_e: Numeric) -> Numeric:
    r"""Electron Coulomb logarithm $\ln\Lambda_e$ in the Sauter convention.

    $$\ln\Lambda_e = 31.3 - \ln\!\left(\frac{\sqrt{n_e}}{T_e}\right)$$

    Parameters
    ----------
    n_e : float or np.ndarray
        Electron density, strictly positive [m^-3].
    T_e : float or np.ndarray
        Electron temperature, strictly positive [eV].

    Returns
    -------
    float or np.ndarray
        Electron Coulomb logarithm [-].

    Raises
    ------
    ValueError
        For non-finite or non-positive input.

    Convention
    ----------
    This is Sauter Eq. (18d), and it is **not** the NRL expression already in
    this package: :func:`vaft.formula.equilibrium.coulomb_logarithm_from_n_T`
    uses 30.9, this uses 31.3.  The difference is about 1.3 percent of a typical
    value, but the collisionality coefficients here were fitted with this one,
    so mixing them biases $\nu_e^*$ consistently.  Use this function wherever a
    Sauter or Redl coefficient is downstream.

    Validity
    --------
    Thermal electrons well above the ionisation stage; the same
    $T_e \gtrsim 10$ eV floor as the NRL form.

    References
    ----------
    .. [1] O. Sauter, C. Angioni and Y. R. Lin-Liu, Phys. Plasmas 6 (1999)
           2834, Eq. (18d).

    See Also
    --------
    vaft.formula.equilibrium.coulomb_logarithm_from_n_T : the NRL convention.
    """
    n_array = _validate_positive("n_e", n_e)
    t_array = _validate_positive("T_e", T_e)
    return _maybe_scalar(31.3 - np.log(np.sqrt(n_array) / t_array))

def coulomb_logarithm_ion_sauter(n_i: Numeric, T_i: Numeric, Z: Numeric) -> Numeric:
    r"""Ion-ion Coulomb logarithm $\ln\Lambda_{ii}$ in the Sauter convention.

    $$\ln\Lambda_{ii} = 30.0 - \ln\!\left(\frac{Z^3\sqrt{n_i}}{T_i^{3/2}}\right)$$

    Parameters
    ----------
    n_i : float or np.ndarray
        Ion density, strictly positive [m^-3].
    T_i : float or np.ndarray
        Ion temperature, strictly positive [eV].
    Z : float or np.ndarray
        Ion charge number, strictly positive [-].

    Returns
    -------
    float or np.ndarray
        Ion-ion Coulomb logarithm [-].

    Raises
    ------
    ValueError
        For non-finite or non-positive input.

    Convention
    ----------
    Sauter Eq. (18e).  The charge enters cubed, so an impurity species changes
    this substantially more than it changes $\ln\Lambda_e$.

    Validity
    --------
    Thermal ions of a single charge state.

    References
    ----------
    .. [1] O. Sauter, C. Angioni and Y. R. Lin-Liu, Phys. Plasmas 6 (1999)
           2834, Eq. (18e).
    """
    n_array = _validate_positive("n_i", n_i)
    t_array = _validate_positive("T_i", T_i)
    z_array = _validate_positive("Z", Z)
    return _maybe_scalar(30.0 - np.log(z_array**3 * np.sqrt(n_array) / t_array**1.5))

def deeply_trapped_bounce_frequency(v_perp: Numeric, q: Numeric, R0: Numeric, epsilon: Numeric) -> Numeric:
    r"""Bounce frequency of a deeply trapped particle, the small-amplitude limit of the banana orbit.

    $$\omega_b = \sqrt{\frac{\epsilon}{2}}\,\frac{v_\perp}{qR_0}$$

    Parameters
    ----------
    v_perp : float or np.ndarray
        Perpendicular speed at the outboard midplane, positive [m/s].
    q : float or np.ndarray
        Safety factor, positive [-].
    R0 : float or np.ndarray
        Major radius, positive [m].
    epsilon : float or np.ndarray
        Inverse aspect ratio $r/R_0$ of the surface, in (0, 1) [-].

    Returns
    -------
    float or np.ndarray
        Angular frequency of the harmonic bounce motion [rad/s].

    Raises
    ------
    ValueError
        A non-positive input or ``epsilon`` outside (0, 1).

    Convention
    ----------
    With $B = B_0(1 - \epsilon\cos\theta)$ and $s = qR_0\theta$ along the
    field, the mirror force $-\mu\,\partial B/\partial s$ is harmonic about
    the outboard midplane for small excursions, with this frequency exactly.
    Barely trapped particles bounce more slowly (the period diverges
    logarithmically at the trapped-passing boundary). The textbook ordering
    $\omega_b \sim \sqrt\epsilon\,v/(qR_0)$ drops the $1/\sqrt 2$; this is
    the exact limit it estimates.

    Physical interpretation
    -----------------------
    $\omega_b/\omega_t = \sqrt{\epsilon/2}$ at $v_\perp = v$: trapped particles
    are slower than passing ones by $\sqrt\epsilon$, which is why the banana
    regime needs a collision frequency below $\epsilon^{3/2}\omega_t$, not
    below $\omega_t$.

    Assumptions
    -----------
    Large aspect ratio, circular surface, deeply trapped
    ($v_\parallel \ll \sqrt\epsilon\,v_\perp$ at the midplane), drifts
    neglected in the parallel motion.

    References
    ----------
    .. [1] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 7.
    .. [2] J. Wesson, *Tokamaks*, 4th ed., Oxford University Press (2011),
           Sec. 3.12.
    """
    v_perp = _validate_positive("v_perp", v_perp)
    q = _validate_positive("q", q)
    R0 = _validate_positive("R0", R0)
    epsilon = _validate_epsilon(epsilon)
    return _maybe_scalar(np.sqrt(epsilon / 2.0) * v_perp / (q * R0))

def neoclassical_regime_boundaries(epsilon: Numeric):
    r"""Banana--plateau and plateau--Pfirsch--Schl\"uter orderings in $\hat\nu$.

    $$\hat\nu_\mathrm{b|p} = \epsilon^{3/2},\qquad \hat\nu_\mathrm{p|PS} = 1$$

    Parameters
    ----------
    epsilon : float or np.ndarray
        Inverse aspect ratio $r/R_0$ of the surface, in (0, 1) [-].

    Returns
    -------
    tuple of float or np.ndarray
        ``(banana_plateau, plateau_pfirsch_schlueter)`` in units of
        ``collisions_per_transit`` [-].

    Raises
    ------
    ValueError
        ``epsilon`` outside (0, 1).

    Convention
    ----------
    Asymptotic orderings, not transitions: $\hat\nu \ll \epsilon^{3/2}$ is
    the banana regime ($\nu_\mathrm{eff} \ll \omega_b$),
    $\epsilon^{3/2} \ll \hat\nu \ll 1$ the plateau, $\hat\nu \gg 1$
    Pfirsch--Schl\"uter. The coefficients are set to one; real transport
    coefficients cross over smoothly across a decade around each value
    (Sauter's fits interpolate across them), so there is deliberately no
    regime classifier.

    Physical interpretation
    -----------------------
    The plateau exists only because trapped particles are slow: its width,
    $\epsilon^{-3/2}$ in $\hat\nu$, shrinks as the aspect ratio falls, and
    at spherical-tokamak $\epsilon$ the plateau nearly disappears.

    Assumptions
    -----------
    Large aspect ratio, circular surfaces.

    References
    ----------
    .. [1] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 8.
    .. [2] F. L. Hinton and R. D. Hazeltine, Rev. Mod. Phys. 48 (1976) 239.
    """
    epsilon = _validate_epsilon(epsilon)
    return _maybe_scalar(epsilon**1.5), _maybe_scalar(np.ones_like(epsilon))

def redl_neoclassical_conductivity(
    sigma_spitzer: Numeric, f_trap: Numeric, nu_e_star: Numeric, Z_eff: Numeric
) -> Numeric:
    r"""Neoclassical parallel conductivity $\sigma_{\mathrm{neo}}$, Redl 2021.

    $$\frac{\sigma_{\mathrm{neo}}}{\sigma_{\mathrm{Sptz}}}
      = 1 - \left(1 + \frac{0.21}{Z}\right) X_{33}
        + \frac{0.54}{Z} X_{33}^2 - \frac{0.33}{Z} X_{33}^3$$

    with $X_{33} = f_t / \left[1 + 0.25(1 - 0.7 f_t)\sqrt{\nu_e^*}
    (1 + 0.45\sqrt{Z-1}) + 0.61(1 - 0.41 f_t)\nu_e^*/\sqrt{Z}\right]$.

    Parameters
    ----------
    sigma_spitzer : float or np.ndarray
        Reference Spitzer conductivity from
        :func:`sauter_spitzer_conductivity` [S m^-1].
    f_trap : float or np.ndarray
        Trapped-particle fraction of the surface, in [0, 1] [-].
    nu_e_star : float or np.ndarray
        Electron collisionality from
        :func:`electron_collisionality_sauter`, non-negative [-].
    Z_eff : float or np.ndarray
        Effective ion charge, at least one [-].

    Returns
    -------
    float or np.ndarray
        Neoclassical parallel conductivity [S m^-1].

    Raises
    ------
    ValueError
        For non-finite input, *f_trap* outside [0, 1], negative *nu_e_star*, or
        *Z_eff* below one.

    Convention
    ----------
    The reference conductivity is still Sauter Eq. (12), so
    :func:`sauter_spitzer_conductivity` is the right input here despite the
    name; Redl refits the correction, not the normalisation.  As in the Sauter
    case the collisionality must be the Eq. (18b) one (issue #353).

    Validity
    --------
    Empirical fit.  Fitted to NEO solutions spanning tight aspect ratio, so it
    stays usable at the trapped fractions where the 1999 fit does not.
    $Z_{\mathrm{eff}} \ge 1$ is required because $\sqrt{Z-1}$ appears.

    References
    ----------
    .. [1] A. Redl, C. Angioni, E. Belli and O. Sauter, Phys. Plasmas 28 (2021)
           022502.
    """
    sigma_array = _validate_positive("sigma_spitzer", sigma_spitzer)
    trapped = _validate_fraction("f_trap", f_trap)
    collisionality = _validate_non_negative("nu_e_star", nu_e_star)
    z_array = np.asarray(Z_eff, dtype=float)
    if np.any(~np.isfinite(z_array)):
        raise ValueError(f"Z_eff must be finite. Got {Z_eff!r}")
    if np.any(z_array < 1.0):
        raise ValueError(f"Z_eff must be >= 1. Got {Z_eff!r}")
    x33 = trapped / (
        1.0
        + 0.25
        * (1.0 - 0.7 * trapped)
        * np.sqrt(collisionality)
        * (1.0 + 0.45 * np.sqrt(z_array - 1.0))
        + 0.61 * (1.0 - 0.41 * trapped) * collisionality / np.sqrt(z_array)
    )
    ratio = (
        1.0
        - (1.0 + 0.21 / z_array) * x33
        + (0.54 / z_array) * x33**2
        - (0.33 / z_array) * x33**3
    )
    return _maybe_scalar(sigma_array * ratio)

def sauter_spitzer_conductivity(
    T_e: Numeric, Z_eff: Numeric, ln_Lambda_e: Numeric
) -> Numeric:
    r"""Spitzer parallel conductivity $\sigma_{\mathrm{Sptz}}$ in the Sauter normalisation.

    $$\sigma_{\mathrm{Sptz}} = 1.9012\times10^{4}\,
      \frac{T_e^{3/2}}{Z_{\mathrm{eff}} N_Z(Z_{\mathrm{eff}}) \ln\Lambda_e},
      \qquad N_Z = 0.58 + \frac{0.74}{0.76 + Z_{\mathrm{eff}}}$$

    Parameters
    ----------
    T_e : float or np.ndarray
        Electron temperature, strictly positive [eV].
    Z_eff : float or np.ndarray
        Effective ion charge, strictly positive [-].
    ln_Lambda_e : float or np.ndarray
        Electron Coulomb logarithm, strictly positive; use
        :func:`coulomb_logarithm_electron_sauter` [-].

    Returns
    -------
    float or np.ndarray
        Classical parallel conductivity, without trapped-particle correction
        [S m^-1].

    Raises
    ------
    ValueError
        For non-finite or non-positive input.

    Convention
    ----------
    This is the reference conductivity the neoclassical correction multiplies,
    Sauter Eq. (12), and it is not the same object as
    :func:`vaft.formula.equilibrium.spitzer_resistivity_from_T_e_Z_eff_ln_Lambda`:
    that one applies the charge dependence linearly, this one through the
    fitted $N_Z$, so their reciprocals differ by tens of percent at
    $Z_{\mathrm{eff}} > 1$.  Pair this one with the neoclassical corrections in
    this module.

    Validity
    --------
    Empirical fit.  $N_Z$ is a fit to the Spitzer-Harm charge dependence, valid
    for $1 \le Z_{\mathrm{eff}} \lesssim 5$.

    References
    ----------
    .. [1] O. Sauter, C. Angioni and Y. R. Lin-Liu, Phys. Plasmas 6 (1999)
           2834, Eq. (12).
    .. [2] L. Spitzer and R. Harm, Phys. Rev. 89 (1953) 977.

    See Also
    --------
    vaft.formula.equilibrium.spitzer_resistivity_from_T_e_Z_eff_ln_Lambda
    """
    t_array = _validate_positive("T_e", T_e)
    z_array = _validate_positive("Z_eff", Z_eff)
    log_array = _validate_positive("ln_Lambda_e", ln_Lambda_e)
    charge_factor = 0.58 + 0.74 / (0.76 + z_array)
    value = _SPITZER_CONDUCTIVITY_COEF * t_array**1.5 / (
        z_array * charge_factor * log_array
    )
    return _maybe_scalar(value)

def transit_frequency(v: Numeric, q: Numeric, R0: Numeric) -> Numeric:
    r"""Poloidal transit frequency of a passing particle on a large-aspect-ratio circular surface.

    $$\omega_t = \frac{v}{qR_0}$$

    Parameters
    ----------
    v : float or np.ndarray
        Speed along the field line of the particle considered, positive; the
        caller chooses which speed (a thermal speed is a convention of its
        own) [m/s].
    q : float or np.ndarray
        Safety factor, positive [-].
    R0 : float or np.ndarray
        Major radius, positive [m].

    Returns
    -------
    float or np.ndarray
        Rate at which the particle advances in poloidal angle [rad/s].

    Raises
    ------
    ValueError
        A non-positive or non-finite input.

    Convention
    ----------
    One poloidal circuit is a parallel path of length $2\pi qR_0$, so
    $\omega_t$ is per radian of poloidal angle; the transit *time* is
    $2\pi/\omega_t$. The speed is the parallel speed of the orbit, not a
    thermal-speed convention (#353 tracks those).

    Physical interpretation
    -----------------------
    The reference rate of the neoclassical orderings: a collision frequency
    compared with $\omega_t$ says how many collisions a particle suffers in
    one connection length (``collisions_per_transit``).

    Assumptions
    -----------
    Large aspect ratio, circular concentric surfaces, $B$ nearly uniform along
    the path.

    References
    ----------
    .. [1] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 7.
    """
    v = _validate_positive("v", v)
    q = _validate_positive("q", q)
    R0 = _validate_positive("R0", R0)
    return _maybe_scalar(v / (q * R0))

def trapped_particle_effective_collision_frequency(nu_deflection: Numeric, epsilon: Numeric) -> Numeric:
    r"""Rate at which pitch-angle scattering detraps a trapped particle.

    $$\nu_\mathrm{eff} = \frac{\nu_D}{\epsilon}$$

    Parameters
    ----------
    nu_deflection : float or np.ndarray
        Pitch-angle-scattering (deflection) frequency $\nu_D$, non-negative
        [1/s].
    epsilon : float or np.ndarray
        Inverse aspect ratio $r/R_0$ of the surface, in (0, 1) [-].

    Returns
    -------
    float or np.ndarray
        Effective collision frequency of the trapped population [1/s].

    Raises
    ------
    ValueError
        A negative ``nu_deflection`` or ``epsilon`` outside (0, 1).

    Convention
    ----------
    The input is the *deflection* frequency, not a momentum-loss or
    energy-exchange frequency; the trapped region is $\sim\sqrt\epsilon$ wide
    in pitch, and diffusion across it takes $\sim\epsilon/\nu_D$. This is the
    conventional definition, with the order-unity factor set to one.

    Physical interpretation
    -----------------------
    Trapped particles are detrapped $1/\epsilon$ times faster than a
    particle is deflected through a radian. The banana regime needs
    $\nu_\mathrm{eff} < \omega_b$, and the ratio of the two is the
    trapped-particle collisionality, the quantity every $\nu_*$ convention
    estimates.

    Assumptions
    -----------
    Large aspect ratio, a diffusive (small-angle) scattering operator.

    References
    ----------
    .. [1] P. Helander and D. J. Sigmar, *Collisional Transport in Magnetized
           Plasmas*, Cambridge University Press (2002), Ch. 7.
    """
    nu_deflection = _validate_non_negative("nu_deflection", nu_deflection)
    epsilon = _validate_epsilon(epsilon)
    return _maybe_scalar(nu_deflection / epsilon)

def trapped_particle_fraction(epsilon: Numeric) -> Numeric:
    r"""Trapped-particle fraction $f_t$ of a circular surface of inverse aspect ratio.

    $$f_t = 1 - \frac{(1-\epsilon)^2}{\sqrt{1-\epsilon^2}\,(1 + 1.46\sqrt{\epsilon})}$$

    Parameters
    ----------
    epsilon : float or np.ndarray
        Inverse aspect ratio $r/R_0$ of the surface, in [0, 1) [-].

    Returns
    -------
    float or np.ndarray
        Fraction of particles trapped on the surface [-].

    Raises
    ------
    ValueError
        For non-finite input, or input outside [0, 1).

    Convention
    ----------
    $f_t$ is the flux-surface quantity that enters every coefficient in this
    module.  This function is the *circular* approximation to it, written in
    terms of the inverse aspect ratio alone.  A shaped equilibrium's trapped
    fraction is an integral over the field-strength distribution on the surface
    and differs from this at the tens-of-percent level at strong elongation, so
    a caller holding a real equilibrium should compute $f_t$ from it and pass
    that instead of calling this.  NEO writes its own value to
    ``out.neo.diagnostic_geo`` as ``f_trap``.

    Physical interpretation
    -----------------------
    The share of the local Maxwellian whose parallel energy is too small to
    cross the magnetic well on the outboard side.  It rises steeply with
    $\epsilon$, which is why bootstrap current matters far more in a spherical
    tokamak than in a conventional one.

    Validity
    --------
    Empirical fit.  Concentric circular surfaces; the coefficient 1.46 is a fit
    to the exact integral rather than a derived value.

    Limitations
    -----------
    At $\epsilon \to 1$ the expression tends to 1 but the underlying expansion
    has long since stopped being controlled.  VEST reaches
    $\epsilon \approx 0.6$, where this returns about 0.91: nearly every
    particle counted as trapped, which is precisely where treating a fitted
    neoclassical coefficient as reliable stops being safe.

    References
    ----------
    .. [1] Y. R. Lin-Liu and R. L. Miller, Phys. Plasmas 2 (1995) 1666.
    .. [2] J. Wesson, *Tokamaks*, 4th ed., Oxford University Press (2011),
           Sec. 4.9 (trapped particles).
    """
    array = np.asarray(epsilon, dtype=float)
    if np.any(~np.isfinite(array)):
        raise ValueError(f"epsilon must be finite. Got {epsilon!r}")
    if np.any(array < 0.0) or np.any(array >= 1.0):
        raise ValueError(f"epsilon must lie in [0, 1). Got {epsilon!r}")
    numerator = (1.0 - array) ** 2
    denominator = np.sqrt(1.0 - array**2) * (1.0 + 1.46 * np.sqrt(array))
    return _maybe_scalar(1.0 - numerator / denominator)