Quasilinear Diffusion Theory

• Quasilinear diffusion theory is a framework that models effective transport via resonant wave–particle interactions and state-dependent coefficients.
• It employs kinetic formulations, Hamiltonian reformulations, and gauge-invariant methods to derive diffusion equations in velocity, momentum, or action space.
• The theory extends to nonlinear PDEs and includes numerical, variational, and stochastic generalizations, ensuring conservation and positivity in complex systems.

Quasilinear diffusion theory denotes a family of diffusion frameworks in which the effective transport law is not fixed a priori but is generated by the state of the system or by an averaged interaction mechanism. In plasma kinetics, it describes the slow evolution of a distribution function under weak, resonant wave–particle interactions, yielding diffusion equations in velocity, momentum, or action space. In nonlinear PDE theory, it designates equations whose highest-order part depends on the unknown or its derivatives, such as $$-\operatorname{div}\!\big(\mu(|\nabla u|)\nabla u\big)=f$$, together with nonlocal, stochastic, kinetic, and weighted reaction–diffusion variants (Brizard et al., 2022, Heid et al., 2021, Lo et al., 2022).

1. Foundational kinetic formulation

In the standard kinetic setting, quasilinear diffusion arises from a weakly perturbed Vlasov system after separating a slowly evolving background distribution from fast oscillatory fluctuations. For a uniform magnetized plasma, the quasilinear evolution of a gyrotropic distribution $f_0(v_\parallel,v_\perp)$ can be written in divergence form,

$$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$

with diffusion tensor

$$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$

The resonance condition

$\omega-k_\parallel v_\parallel-\ell\Omega=0$

is therefore the organizing principle of the theory: only resonant particles contribute to secular diffusion (Brizard et al., 2022).

The classical assumptions are weak perturbations, random phases, separation of time scales, and dominance of resonant interactions. In this regime, the cumulative effect of many small wave-induced kicks is diffusive rather than ballistic. In invariant velocity variables for a uniform field, the natural coordinates are $p_\parallel=Mv_\parallel$ and $\mu=Mv_\perp^2/(2B_0)$, and the standard Kennel–Engelmann equation becomes a symmetric $2\times2$ diffusion system in $(p_\parallel,\mu)$ or, equivalently, in $(\mathcal E,\mu)$, where the energy kick depends only on the electric field (Brizard et al., 2022).

An analogous quasilinear closure appears for long-range Vlasov systems. For the Hamiltonian Mean Field model, one writes

$f_0(v_\parallel,v_\perp)$0

linearizes the Vlasov equation around the slowly varying homogeneous part $f_0(v_\parallel,v_\perp)$1, and derives a diffusion equation

$f_0(v_\parallel,v_\perp)$2

For even single-humped distributions in the HMF model, only the $f_0(v_\parallel,v_\perp)$3 modes matter, the unstable eigenfrequency is purely imaginary, and the diffusion coefficient reduces to

$f_0(v_\parallel,v_\perp)$4

This quasilinear theory works reasonably well for weakly unstable initial conditions and predicts the energy marking the out-of-equilibrium phase transition between unmagnetized and magnetized quasi-stationary states (Campa et al., 2016).

2. Hamiltonian and invariant-space reformulations

A central development in modern plasma quasilinear theory is its Hamiltonian reformulation. Instead of expressing the wave forcing directly through $f_0(v_\parallel,v_\perp)$5, one introduces the gauge-invariant effective potential

$f_0(v_\parallel,v_\perp)$6

and rewrites the Vlasov equation with a noncanonical Poisson bracket. In this representation, the second-order quasilinear evolution is generated by

$f_0(v_\parallel,v_\perp)$7

which makes gauge invariance, entropy production, and the dyadic structure of the diffusion tensor explicit (Brizard et al., 2022).

In a uniform magnetic field, the natural invariant space is two-dimensional. In a nonuniform axisymmetric magnetic field, the appropriate invariants become

$f_0(v_\parallel,v_\perp)$8

so quasilinear transport is intrinsically three-dimensional in invariant space. The corresponding evolution equation is

$f_0(v_\parallel,v_\perp)$9

and the diffusion tensor has the explicit $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$0 form

$$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$1

This tensor couples gyro/pitch-angle transport, energy diffusion, and radial diffusion in a single object (Brizard et al., 2022).

The same structure is obtained abstractly in canonical action–angle variables: $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$2 This canonical form shows that diffusion is along resonance directions $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$3 in action space and that the tensor is a sum of dyads weighted by resonant quasilinear potentials (Brizard et al., 2022).

A closely related conservation-law refinement appears in the extension of the Kennel–Engelmann tensor from two to four dimensions. The standard magnetized tensor respects the energy–parallel-momentum relation but neglects perpendicular momentum absorption. Enforcing four-momentum conservation leads to

$$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$4

which, after transformation to constants-of-motion variables $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$5, yields the diffusion path

$$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$6

This matches the form required by action-angle Hamiltonian theory after bounce averaging (Ochs, 12 Nov 2025).

3. Positivity, conservation, and numerical structure

For toroidal plasmas, bounce averaging of the classical Kennel–Engelmann coefficients can destroy positive definiteness because the parallel inhomogeneity makes the resonance kernel asymmetric under the bounce average. A positive-definite alternative is obtained by evaluating the phase integral along the trajectory before averaging. The resulting bounce-averaged diffusion tensor can be written as

$$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$7

which is manifestly positive semidefinite. In this form, resonant contributions are expressed through Airy-function factors near resonance points, while nonresonant contributions arise near outer-midplane, inner-midplane, and trapping-tip locations where the phase curvature is small. The construction includes both resonant and non-resonant contributions, and the correlations between the consecutive resonances and in many poloidal periods (Lee et al., 2017).

A related finite-Larmor-radius reduction has been developed for ion cyclotron heating. Starting from the kinetic energy change $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$8 and matching it to the quasilinear Fokker–Planck operator, one derives reduced coefficients that preserve the diffusion directions, wave polarizations, and H-theorem of the full Kennel–Engelmann model. For the fundamental damping, the lowest-order reduction corresponds to $$\frac{1}{\Omega}\frac{\partial f_0}{\partial \tau} = \frac{\partial}{\partial \mathbf{v}}\cdot \left(\mathbf{D}\cdot \frac{\partial f_0}{\partial \mathbf{v}}\right),$$9 and neglect of higher Bessel terms; for the second harmonic damping, the lowest nonzero contribution is $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$0 and corresponds to $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$1 (Lee et al., 2017).

Conservation can also be enforced at the discretization level. In a three-dimensional momentum / three-dimensional spectral quasilinear model with cylindrical symmetry, the particle pdf and wave sed satisfy coupled bilinear equations that admit an unconditionally conservative weak form. A conservative Galerkin discretization with continuous basis functions for the particle pdf, discontinuous basis functions for the wave sed, and a consistent quadrature rule preserves particle number, momentum, and energy rigorously, independently of the singular transition probability. The resonance manifold is integrated by a marching simplex algorithm, which converts the singular kernel into quadrature over simplices on the manifold (Huang et al., 2022).

These developments show that positivity and conservation are not merely numerical conveniences. They encode the structural content of quasilinear diffusion: entropy production, compatibility with wave action, and preservation of the geometric diffusion directions determined by resonance.

4. Variational and elliptic-parabolic PDE formulations

In nonlinear PDE theory, quasilinear diffusion refers to equations whose principal part depends on the unknown or on its gradient. A prototypical elliptic model is

$$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$2

with mixed boundary conditions on $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$3. In the variational setting,

$$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$4

and the weak problem is equivalent to minimizing $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$5 over the admissible affine set $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$6 (Heid et al., 2021).

The classical Kačanov scheme freezes the diffusion law at the previous iterate: $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$7 or, equivalently,

$$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$8

Its traditional convergence theory assumes that $$\mathbf{D} = \sum_{\ell}\mathrm{Re}(-i\Delta_\ell)\, \tilde{\mathbf{v}}_\ell^\ast\tilde{\mathbf{v}}_\ell, \qquad \mathrm{Re}(-i\Delta_\ell) = \pi\Omega\,\delta(\omega-k_\parallel v_\parallel-\ell\Omega).$$9 is continuously differentiable, decreasing, uniformly bounded above and below, and associated with a uniformly convex energy. The key limitation is the monotonicity requirement $\omega-k_\parallel v_\parallel-\ell\Omega=0$0, which excludes shear-thickening and non-monotone laws even though the iteration often converges numerically (Heid et al., 2021).

A modified Kačanov iteration introduces damping,

$\omega-k_\parallel v_\parallel-\ell\Omega=0$1

and proves convergence from an energy decrease estimate

$\omega-k_\parallel v_\parallel-\ell\Omega=0$2

For quasilinear diffusion, this estimate follows from the two-sided structural condition

$\omega-k_\parallel v_\parallel-\ell\Omega=0$3

which yields strong monotonicity and Lipschitz continuity of $\omega-k_\parallel v_\parallel-\ell\Omega=0$4 without requiring $\omega-k_\parallel v_\parallel-\ell\Omega=0$5 to be decreasing or differentiable. The resulting solver converges for decreasing, non-monotone, and strongly non-monotone viscosity laws, and adaptive damping can both accelerate convergence and recover convergence where the classical Kačanov iteration fails (Heid et al., 2021).

A time-dependent quasilinear reaction–diffusion model displays a different aspect of the theory: $\omega-k_\parallel v_\parallel-\ell\Omega=0$6 For initial data satisfying

$\omega-k_\parallel v_\parallel-\ell\Omega=0$7

the solution remains uniformly bounded above and below,

$\omega-k_\parallel v_\parallel-\ell\Omega=0$8

For compactly supported data, the positivity set expands with finite speed and obeys

$\omega-k_\parallel v_\parallel-\ell\Omega=0$9

while in outer sets one has

$p_\parallel=Mv_\parallel$0

These results identify $p_\parallel=Mv_\parallel$1 as a threshold between grow-up regimes and globally bounded regimes (Iagar et al., 7 Jun 2026).

5. Nonlocal, stochastic, and kinetic generalizations

A broad abstract framework for quasilinear diffusion systems replaces the local elliptic operator by a symmetric, coercive operator $p_\parallel=Mv_\parallel$2 on $p_\parallel=Mv_\parallel$3 and allows the coefficients to depend on both the state and nonlocal derivatives: $p_\parallel=Mv_\parallel$4 Here $p_\parallel=Mv_\parallel$5 may represent classical gradients, Riesz fractional gradients, nonlocal gradients, or higher-order fractional derivatives, with order $p_\parallel=Mv_\parallel$6. Under boundedness and coercivity of $p_\parallel=Mv_\parallel$7 and of the matrix $p_\parallel=Mv_\parallel$8, together with linear or sublinear growth conditions on $p_\parallel=Mv_\parallel$9, one obtains global existence in the maximal regularity space

$\mu=Mv_\perp^2/(2B_0)$0

for local, fractional, and anisotropic nonlocal diffusion operators alike (Lo et al., 2022).

In the stochastic setting, the abstract quasilinear equation

$\mu=Mv_\perp^2/(2B_0)$1

is treated by combining deterministic quasilinear parabolic theory with evolution semigroups. Because the evolution family $\mu=Mv_\perp^2/(2B_0)$2 is not adapted in the Itô sense, the mild formulation is replaced by a pathwise mild representation in which the stochastic convolution is rewritten by integration by parts. Under sectoriality, Lipschitz continuity of $\mu=Mv_\perp^2/(2B_0)$3 in suitable fractional-domain scales, and local Lipschitz conditions on $\mu=Mv_\perp^2/(2B_0)$4 and $\mu=Mv_\perp^2/(2B_0)$5, the equation has a unique local pathwise mild solution and a maximal local solution characterized by the blow-up alternative

$\mu=Mv_\perp^2/(2B_0)$6

(Kuehn et al., 2018).

A further hypoelliptic generalization arises in kinetic theory. For the Kolmogorov operator

$\mu=Mv_\perp^2/(2B_0)$7

kinetic maximal $\mu=Mv_\perp^2/(2B_0)$8-regularity is established in anisotropic spaces $\mu=Mv_\perp^2/(2B_0)$9 adapted to the kinetic scaling, and the trace space is identified with anisotropic Besov space

$2\times2$0

This linear theory supports a local well-posedness result for quasilinear kinetic diffusion,

$2\times2$1

where $2\times2$2 is uniformly elliptic. The solution theory is formulated in weighted kinetic maximal-regularity spaces and yields local existence, uniqueness, continuous dependence, and instantaneous smoothing (Niebel et al., 2020).

6. Regimes, thresholds, and broader significance

A recurring theme in quasilinear diffusion theory is that the formal diffusion equation remains meaningful only within a specific dynamical regime. For a single charged particle interacting with a discrete spectrum of electrostatic waves, the perturbative quasilinear regime is characterized by free-streaming-based diffusion in velocity space with diffusion coefficient

$2\times2$3

In that regime, diffusion occurs only when wave–particle interaction is local in phase velocity; conversely, numerical results indicate that chaotic diffusion can occur even when wave–particle interaction is not local. KAM tori bound the accessible velocity interval, so diffusion is intrinsically finite in extent, and a renormalized Gaussian kernel can model the time evolution of the velocity distribution while accounting for phase-space boundaries (Bénisti, 20 Aug 2025).

The long-range Vlasov setting shows a different threshold phenomenon. Quasilinear diffusion around a weakly unstable homogeneous state drives the angle-averaged distribution toward a marginally stable quasi-stationary state. For Gaussian and semi-elliptical initial data in the HMF model, the quasilinear theory works reasonably well for weakly unstable initial conditions and predicts the energy marking the out-of-equilibrium phase transition between unmagnetized and magnetized quasi-stationary states; at lower energies, the disagreement grows, and the quasi-stationary states are remarkably well fitted by polytropic distributions with index $2\times2$4 in the Gaussian case or $2\times2$5 in the semi-elliptical case (Campa et al., 2016).

Taken together, these results show that quasilinear diffusion theory is not a single theorem but a class of asymptotic closures. In plasma physics it links resonant wave spectra, invariant-space transport, and conservation laws; in nonlinear PDEs it organizes existence theory, variational structure, iterative linearization, and asymptotic rates for equations with state-dependent principal part. The common content is the replacement of a complicated nonlinear or oscillatory microscopic dynamics by an effective diffusion process whose coefficients are themselves determined by resonance, invariants, or nonlinear constitutive structure (Brizard et al., 2022, Heid et al., 2021).