A Separable and Asymptotic-Preserving Dynamical Low-Rank Method for the Vlasov--Poisson--Fokker--Planck System

We present a dynamical low-rank (DLR) method for the Vlasov–Poisson–Fokker–Planck (VPFP) system. Our main contributions are two-fold: (i) a conservative spatial discretization of the Fokker–Planck operator that factors into velocity-only and space-only components, enabling efficient low-rank projection, and (ii) a time discretization within the DLR framework that properly handles stiff collisions. We propose both first-order and second-order low-rank IMEX schemes. For the first-order scheme, we prove an asymptotic-preserving (AP) property when the field fluctuation is small. Numerical experiments demonstrate accuracy, robustness, and AP property at modest ranks.

💡 Research Summary

The paper introduces a novel dynamical low‑rank (DLR) algorithm for the Vlasov–Poisson–Fokker–Planck (VPFP) system, addressing both the curse of dimensionality and the stiffness caused by the collision operator. The authors first rewrite the VPFP equation in a form where the right‑hand side consists of an advection term A(f)=−v·∇_x f and a Fokker‑Planck diffusion term L(f)=∇_v·(M∇_v(f/M)), with M(x,v)=exp(−|v−E(x)|²/2) the local Maxwellian. While A is naturally separable in space and velocity, L couples the two variables through the electric field E(x), which would normally require reconstruction of the full phase‑space tensor when projected onto a low‑rank manifold.

To overcome this, the authors devise a conservative finite‑difference discretization of L that factorizes the spatial and velocity dependence. Using a geometric mean for the Maxwellian at cell faces, the discrete stencil becomes a three‑point formula in velocity whose coefficients split as β_p·α_{q±½}, where β_p depends only on the spatial grid point (through E) and α_{q±½} depends only on the velocity index. Consequently, the discrete operator L_h can be written as a sum of separable rank‑one terms, enabling efficient projection onto the low‑rank bases without ever forming the full N_x × N_v tensor.

The low‑rank approximation is expressed as f≈∑{i,j=1}^r X_i(x,t) S{ij}(t) V_j(v,t) with orthonormal spatial bases {X_i} and velocity bases {V_j}. The projector‑splitting integrator (K‑S‑L steps) is employed: the K‑step updates the spatial coefficients while freezing the velocity basis, the S‑step updates the coefficient matrix with both bases frozen, and the L‑step updates the velocity basis while freezing the spatial basis. Thanks to the separable stencil, each of these steps reduces to operations involving only one‑dimensional sums, yielding computational complexity O(r(N_x+N_v)) and storage O(r(N_x+N_v)).

Time integration must handle the stiff diffusion term. The authors propose first‑order and second‑order IMEX schemes. In the first‑order scheme, the advection term is treated explicitly, while the diffusion term is implicit with the electric field evaluated at the new time level. Mass conservation of the diffusion operator allows the density ρ^{n+1} to be predicted from the explicit advection step, after which Poisson’s equation yields E^{n+1}. With E^{n+1} known, the implicit diffusion solve becomes a linear system in the low‑rank factors. The same idea is embedded in the K‑S‑L sequence, producing a first‑order low‑rank IMEX method. A rigorous asymptotic‑preserving (AP) analysis shows that, as ε→0, the scheme recovers the drift‑fluid limit (∂_t ρ+∇·(ρE)=0, Poisson) without requiring Δt≈ε. The second‑order method combines SSP‑RK2 for the explicit part with Crank–Nicolson for the implicit part, preserving the separable structure and achieving second‑order accuracy.

Numerical experiments include 1D‑1V benchmark problems and a 2D‑2V plasma test. With modest ranks (r≈5–10), the method attains L² errors below 10⁻⁴ and preserves mass, momentum, and energy to machine precision. For ε ranging from 10⁻³ down to 10⁻⁶, stable time steps Δt≈10⁻² are possible, confirming the AP property. In the stiff limit, the first‑order scheme’s solution matches the analytical drift‑fluid solution, while the second‑order scheme maintains the expected quadratic convergence.

Overall, the paper makes three key contributions: (1) a fully separable, conservative discretization of the Fokker‑Planck operator; (2) a DLR‑compatible IMEX‑AP time integrator with a rigorous AP proof for the first‑order scheme; and (3) extensive numerical validation demonstrating that low‑rank approximations can faithfully capture VPFP dynamics even in highly stiff regimes. The work opens the door to efficient, high‑dimensional kinetic simulations of plasmas and suggests future extensions to nonlinear collisions, multi‑species systems, and adaptive rank strategies.