A Parallel-Kinetic-Perpendicular-Moment Model for Magnetized Plasmas

TL;DR

This work formulates the parallel-kinetic-perpendicular-moment (PKPM) model, a reduced kinetic framework for weakly collisional, magnetized plasmas that preserves parallel dynamics while representing perpendicular velocity structure with a spectral Laguerre-Fourier closure. By transforming velocity coordinates to a frame moving with the local flow and then to a field-aligned (CGL) frame, the authors optimize the spectral basis and derive a minimal yet physically rich set of equations: a gyrotropic distribution expanded in Laguerre in $v_ perp$ and a couple of Laguerre moments (through $F_0$ and $F_1$) coupled to a parallel-kinetic coordinate $v_ parallel$, plus a momentum equation closed by the gyrotropic pressure tensor. The lowest-order PKPM system reduces the full Vlasov-Maxwell problem to two coupled 4D kinetic equations plus a momentum equation, enabling efficient simulation of phenomena such as parallel electrostatic shocks and moderate guide-field magnetic reconnection, while capturing essential finite-Larmor-radius effects through a limited number of Fourier harmonics. The paper situates PKPM within the historical KMHD/Ramos framework, contrasts it with other spectral/hybrid approaches, and demonstrates its physics fidelity and numerical robustness through two nonlinear benchmarks, highlighting its potential as a cost-effective, scalable tool for magnetized plasma dynamics and its role in a broader multi-part PKPM research program.

Abstract

We describe a new model for the study of weakly-collisional, magnetized plasmas derived from exploiting the separation of the dynamics parallel and perpendicular to the magnetic field. This unique system of equations retains the particle dynamics parallel to the magnetic field while approximating the perpendicular dynamics through a spectral expansion in the perpendicular degrees of freedom, analogous to moment-based fluid approaches. In so doing, a hybrid approach is obtained which is computationally efficient enough to allow for larger-scale modeling of plasma systems while eliminating a source of difficulty in deriving fluid equations applicable to magnetized plasmas. We connect this system of equations to historical asymptotic models and discuss advantages and disadvantages of this approach, including the extension of this parallel-kinetic-perpendicular-moment beyond the typical region of validity of these more traditional asymptotic models. This paper forms the first of a multi-part series on this new model, covering the theory and derivation, alongside demonstration benchmarks of this approach including shocks and magnetic reconnection.

Figures (11)

• Figure 1: Electron (left column) and proton (right column) mass density (top row), momentum density (middle top row), total energy density (middle bottom row), and parallel pressure (bottom row) at $t=1500\omega_{pe}^{-1}$ for upstream Mach number $M_s = 3.0$ (black) and $M_s = 5.0$ (red). All quantities are normalized to their upstream values for ease of comparison between the $M_s = 3.0$ and $M_s = 5.0$ cases since their upstream flows and energies are different. The characteristics of a shock wave are clearly identifiable in the $M_s = 3.0$ simulation: a sharp pile-up of the density, a rapid stagnation of the flow, significant electron heating over the same length scale, and a decrease in the ion energy from the rapid conversion of ion energy into both electron heating and electromagnetic energy. On the other hand, the $M_s = 5.0$ case exhibits no such sharp transitions, with a smooth gradient up to a total mass density $\rho \sim 2$ and total momentum $\rho u_x \sim 0.0$ for both the electrons and protons, corresponding to two interpenetrating beams of plasma.
• Figure 2: The $F_0$ distribution function in the local fluid flow frame for the electrons (left column) and protons (left middle column), and the $F_0$ distribution function in the lab frame for the electrons (right middle column) and protons (right column) for the $M_s = 3.0$ (top row) and $M_s = 5.0$ (bottom row) simulations. In the lab frame, the incoming proton beam is centered at the upstream velocity, $u_{x_0} = 6.0 v_{th_i}$ and $u_{x_0} = 10.0$, as we expect, and the characteristics of the shock with the trapped electron and ion populations are identifiable in the $M_s = 3.0$ simulation, while the $M_s = 5.0$ simulation shows only two distinct ion beams propagating through each other. We also draw attention to the form of the proton distribution function in the local fluid flow frame and emphasize that these are the distribution functions which are are directly solved for by the numerical method. As expected, the distribution function adjusts in the local fluid flow frame to preserve the identity that the first moment is zero. Note that these distribution function plots are normalized to their respective maximum values on the grid, e.g., $F_{0_s} = F_{0_s}/\max(F_{0_s})$.
• Figure 3: Evolution of the out-of-plane current density, $J_z$ with contours of the in-plane magnetic field super-imposed by computing $A_z$, the out-of-plane vector potential, from the in-plane $B_x$ and $B_y$ for the $B_g = B_0$ (left) and $B_g = 0.5 B_0$ (right) lower resolution simulations. We observe morphologies of the current layer consistent with Le:2013, which found at lower electron $\beta_e$ a transition from a regime at lower guide field in which an extended current layer forms from the magnetized electrons developing strong anisotropy and driving a perpendicular current across field lines, to a regime in which the magnetic tension in the guide field causes the current and density to peak near the diagonally opposed separator field lines and negate the impact of the electron anisotropy on the magnetic field's overall tension---see Figure \ref{['fig:pkpm_comp_temp_aniso']}. This contrast is especially clear at $t=20 \Omega_{ci}^{-1}$ and $30 \Omega_{ci}^{-1}$ as the reconnection rate reaches its peak values---see Figure \ref{['fig:recon_rate']}---and we can see a more concentrated current layer in the $B_g = B_0$ simulation compared to the extended current layer in the $B_g = 0.5 B_0$ simulation.
• Figure 4: Zoom in of the $B_g = 0.5 B_0$ simulation with lower resolution and larger hyper-diffusion (left) and higher resolution and smaller hyper-diffusion (right). While the mode is identifiable in the lower resolution simulation, the secondary instability is especially prominent at increased resolution.
• Figure 5: Comparison of the electron parallel temperature normalized to the initial electron temperature (top), electron perpendicular temperature normalized to the initial electron temperature (middle top), electron temperature anisotropy (middle bottom), and electron firehose criteria (bottom) at $t = 20 \Omega_{ci}^{-1}$ for the $B_g = B_0$ simulation (left) and $B_g = 0.5 B_0$ simulation (right). In both cases, a significant electron anisotropy from an excess of parallel pressure develops in the current layer, but a depletion of electron perpendicular pressure in the $B_g = 0.5 B_0$ simulation further increases the electron anisotropy in the layer. Combined with the lower guide field and thus a weaker magnetic field at the X-point, the electron firehose criteria is much closer to marginal stability $\Lambda_{firehose} \sim 0$ for the $B_g = 0.5 B_0$ simulation. Thus, the electrons more significantly modify the tension in the magnetic field at the reconnecting X-point compared to the higher guide field simulation, driving a perpendicular current that spreads the current layer into the exhaust.