3D MHD modelling of plasmoid drift following massive material injection in a tokamak

Mechanisms of plasmoid drift following massive material injection are studied via 3D non-linear MHD modelling with the JOREK code, using a transient neutral source deposited at the low field side midplane of a JET H-mode plasma to clarify basic processes and compare with existing theories. The simulations confirm the important role of the propagation of shear Alfvén wave (SAW) packets from both ends of the plasmoid (SAW braking'') and the development of external resistive currents along magnetic field lines (Pégourié braking’’) in limiting charge separation and thus the $\mathbf{E}\times \mathbf{B}$ plasmoid drift, where $\mathbf{E}$ and $\mathbf{B}$ are the electric and magnetic fields, respectively. The drift velocity is found to be limited by the SAW braking on the few microseconds timescale for cases with relatively small source amplitude while the Pégourié braking acting on a longer timescale is shown to set in earlier with larger toroidal extent of the source, both in good agreement with existing theories. The simulations also identify the key role of the size of the $\mathbf{E}\times \mathbf{B}$ flow region on plasmoid drift and show that the saturated velocity caused by dominant SAW braking agrees well with theory when considering an effective pressure within the $\mathbf{E}\times \mathbf{B}$ flow region. The existence of SAWs in the simulations is demonstrated and the 3D picture of plasmoid drift is discussed.

💡 Research Summary

This paper presents a comprehensive three‑dimensional non‑linear magnetohydrodynamic (MHD) study of plasmoid drift that follows massive material injection (MMI) in a tokamak, using the JOREK code to model a transient neutral source deposited on the low‑field‑side (LFS) mid‑plane of a JET H‑mode discharge. The authors aim to clarify the basic physical processes governing plasmoid motion and to benchmark the simulation results against existing analytical theories that describe two principal “braking” mechanisms: shear Alfvén wave (SAW) braking and Pégourié braking.

The physical picture underlying the theory is that a dense, high‑pressure plasmoid created by pellet ablation experiences a vertical polarization due to the ∇B drift, which generates a charge separation. This separation drives an E×B drift toward the LFS. The drift is initially accelerated by the pressure imbalance (the first term in Eq. 1) but is rapidly limited by currents that act to neutralize the charge separation. SAW braking arises from Alfvénic wave packets launched from both ends of the plasmoid; these packets propagate at the Alfvén speed C_A and carry current j_A = –∇⊥²Φ/(μ₀ C_A), thereby providing an external path for the charge to flow and reducing the net polarization. Pégourié braking, on the other hand, involves the development of external, purely parallel resistive currents that connect the top and bottom of the plasmoid after the SAW packets have reflected back. This mechanism operates on a longer timescale that depends on the toroidal extent of the plasmoid (through the external connection length L_con).

The JOREK simulations reproduce both mechanisms. A series of runs were performed with neutral deuterium injection rates ranging from 2 × 10²⁶ to 1 × 10²⁸ atoms s⁻¹ (total injected atoms 10¹⁹–5 × 10²⁰) and with toroidal source widths Δϕ from 0.5 rad to 2 rad. The baseline equilibrium corresponds to JET pulse #96874 (R ≈ 2.96 m, I_p ≈ 3 MA, B_ϕ₀ ≈ 2.8 T, central T_e ≈ 7 keV, n_e ≈ 0.85 × 10²⁰ m⁻³). The neutral cloud has a Gaussian shape with half‑widths Δr = L_θ = 4 cm in the radial and poloidal directions; the toroidal width is varied as described.

Key findings are:

1. SAW Braking Dominates at Low Source Amplitude – For injection rates ≤ 1.4 × 10²⁷ s⁻¹, the saturated drift velocity V_D,lim scales linearly with the plasmoid pressure, exactly as predicted by Eq. (2) (V_D,lim ≈ μ₀ C_A p₀ Z₀ R / B_ϕ²). The simulations show that the drift accelerates within a few microseconds and then plateaus at a value that matches the theoretical estimate when the effective pressure inside the E×B flow region is used.

2. Pégourié Braking Becomes Important for Large Amplitude or Large Toroidal Extent – When the injection rate is increased to 1 × 10²⁸ s⁻¹, the observed V_D,lim is roughly half of the value expected from simple SAW scaling. The authors attribute this reduction to an earlier onset of Pégourié braking, which is facilitated by stronger magnetic field line stochasticity and a shorter external connection length L_con. The toroidal width study confirms this: increasing Δϕ from 0.5 rad to 2 rad leads to an earlier and more rapid decline of V_D, indicating that the longer toroidal extent in the simulation artificially delays Pégourié braking, a known limitation of present 3D MHD codes.

3. Size of the E×B Flow Region Controls the Saturated Velocity – By analyzing the spatial distribution of the electrostatic potential U (which determines the E×B flow), the authors demonstrate that a larger flow region encloses more plasma pressure, raising the effective pressure p_eff that enters Eq. (2). Consequently, simulations with broader plasmoids (larger Δϕ) exhibit higher V_D,lim even when the injected neutral amount is the same, provided SAW braking remains dominant.

4. Direct Evidence of SAW Packets – The paper presents spatio‑temporal maps of the parallel current density and the electrostatic potential, showing wave‑like structures propagating from the plasmoid ends at a speed consistent with the local Alfvén velocity. These structures match the theoretical description of SAW packets and confirm that the current closure path involves both parallel and perpendicular components as depicted in the schematic Fig. 1(b).

5. Implications for ITER and Future Experiments – The authors discuss that realistic SPI scenarios in ITER will involve plasmoids with toroidal extents comparable to the machine circumference, whereas current 3D MHD simulations must often use enlarged Δϕ for numerical stability. This artificial enlargement can underestimate the impact of Pégourié braking, leading to over‑prediction of drift velocities and, consequently, of the mitigation efficiency of LFS‑deuterium SPI.

Overall, the study validates the two‑braking theory within a full 3D non‑linear MHD framework, quantifies the parameter regimes where each mechanism dominates, and highlights the importance of correctly representing toroidal geometry and E×B flow region size for predictive modeling of disruption mitigation. Future work suggested includes higher toroidal resolution, incorporation of kinetic neutral models, and coupling to impurity radiation physics to capture the full dynamics of massive material injection in next‑generation devices.