Stellarator ideal magnetohydrodynamic (MHD) codes that assume nested flux surfaces such as \texttt{VMEC} and \texttt{DESC} solve the equilibrium problem by prescribing the total toroidal magnetic flux, plasma profiles and the last closed flux surface (LCFS), which analytically determines the unique field in vacuum, whereas at finite $β$ bifurcations and distinct equilibria sharing the same boundary have been reported in the literature. In this paper, we propose prescribing the Poincaré cross-section of the field at a single toroidal plane, and implement it in \texttt{DESC}, where the new condition enters only through the linear constraints and therefore costs no more than a fixed-LCFS solve. Fixing the geometry on one plane rather than on a full toroidal surface generally leaves more of the spectral coefficients free, and the ones it frees are those carrying the toroidal variation of the boundary flux surface, which a prescribed LCFS holds fixed at every toroidal angle; together these allow better-converged numerical solutions. We solve equilibria with the cross-section held fixed, starting either from the axisymmetric shape obtained by revolving that cross-section toroidally, or from an existing fixed-LCFS solution. In the latter case, the volume-averaged normalized force error falls by an order of magnitude while the configuration stays close to the original one, and re-solving the resulting boundary with the conventional fixed-LCFS solver recovers the same equilibrium. We further show that the Poincaré coefficients can be used directly as design variables by optimizing a quasi-helical configuration that maintains high-fidelity force balance throughout the process.
Analytical studies of zonal-flow (ZF) excitation by fluctuating fields, such as drift-wave turbulence and Alfvén eigenmodes (AEs) driven by energetic particles (EPs) in tokamaks, are usually carried out via intricate gyrokinetic calculations and simulations whose results are notoriously difficult to interpret. Here, we report a transparent result that is also not limited to any particular modes or frequency range. Using oscillation-center theory that captures both ponderomotive forces and quasilinear diffusion, we derive a compact formula for the ZF drive produced by any fluctuating field in terms of the field's canonical momentum and the dissipation power density. For ZFs generated by drift waves, our model subsumes the local relation between the zonal velocity and the drift-wave energy density that was previously derived ad hoc. The simplicity and generality of our result also opens a path toward optimization of ZF excitation with external waves and bridges the physics of AE-EP interactions with that of poloidal flows driven by radiofrequency waves in fusion plasmas.
Quasi-poloidal (QP) magnetic fields have desirable properties for confining plasma: no radial drift of guiding centres (with positive implications for neoclassical transport), zero Pfirsch–Schlüter current and a lower level of damping for poloidal flows. Despite their attractive properties, QP fields are not amenable to the near-axis expansion, a major theoretical tool for understanding toroidal fields. In this paper we provide a novel framework for defining and understanding QP flux surfaces. This framework relies on a simplification that transforms the task of finding a QP flux surface from a three-dimensional problem to a two-dimensional (2D) problem. This simplification also applies to asymmetric magnetic mirrors with desirable properties. We sketch how this 2D problem can form the basis of an efficient optimisation problem for finding QP flux surfaces. We leverage this 2D problem for theoretical understanding: for instance, we identify a route to finding QP flux surfaces that are naturally flat mirrors (Velasco et al. 2023, Nucl. Fusion , 63, 126038). The reduced model is qualitatively checked against numerically optimised QP equilibria. These numerical solutions only satisfy QP approximately, but we predictably find that local discrepancies with the reduced model correspond to significant local QP errors, anomalous parallel currents and field lines deviating from geodesics.
Theories of ion-scale microinstabilities in tokamaks and stellarators typically assume that the passing electrons respond adiabatically due to their fast propagation speed. However, when the magnetic shear becomes sufficiently small, ion-scale modes can extend far along the magnetic field and the non-adiabatic response of passing electrons becomes important. We derive a theory of extended modes at low magnetic shear through a multiscale expansion of the gyrokinetic equation. The theory elucidates the physics of the geodesic extended mode, a new type of microinstability. The new mode couples the non-adiabatic physics of both electrons and ions, unlike extended modes at magnetic shear of order unity. The theory is validated against gyrokinetic simulations and the parameter dependences of the new mode are studied.
Integrated codes simulating interactions between Alfvén waves and fast ions in tokamak plasmas use perturbative models for the relatively slow processes of instability growth, saturation, chirping and bursting, and transport. Faster processes by which an Alfvén mode's spatiotemporal structure forms are assumed to have been completed within the mode's oscillation period, $τ_0 \equiv 2π/ω_0$. This separation of time scales underlies the computational efficiency of perturbative models, where the Alfvén mode's time-dependence is reduced to that of a scalar signal $s(t) = A(t)\sin(-ω_0 t - φ(t))$ with variable amplitude $A(t)$ and phase $φ(t)$. For this, accurate input data in the form of a mode's spatial structure $δΦ({\mathbf x})$, damping rate $γ_{\rm d}$, and initial frequency $ω_0$ are required. For modes residing in dense or continuous spectra, $δΦ$ and $γ_{\rm d}$ could be estimated from the form of the continua and fast ion orbits, but it is difficult to guess the seed frequency $ω_0$. Here, we report results of numerical experiments showing that it is possible to find $ω_0$ using a prompt frequency shift that occurs during the first few $100$ time steps of a simulation. Restarts with the shifted frequency iteratively converge to a value of $ω_0$ that seems to maximize the resonant drive, suggesting an auto-optimization process. The need for iteration is attributed to the fact that the terms required for rapid frequency adjustments were truncated when deriving the perturbative model. Meanwhile, the fact that partial auto-optimization is possible at all is attributed to the fact that series truncation alone (without filter) does not strictly enforce slowness. Remnants of and cross-talk with faster dynamics still occur in numerical implementations. This entails potential for both uncertainty and utility.
R. Michael Churchill, Matt Landreman, Jong Youl Choi, Byoungchan Jang, Rory Conlin, Noah Mandell, Anima Anandkumar, Valentin Duruisseaux, Jeffrey Larson, Dario Panici, et al.
Previous work built AI-based surrogates for a nonlinear gyrokinetic simulation code with the goal of using them for fast, direct calculations of turbulent ion heat flux in stellarator design optimizations and scenario planning for experiments. These AI surrogates were trained on data from >200k nonlinear, adiabatic electron gyrokinetic simulations with the gyrokinetic flux-tube code GX, using a wide range of stellarator magnetic configurations ($\sim$23k), positions in the plasma, and gradient scale lengths. In this paper, we demonstrate the use of the AI-based turbulence surrogate in the optimization of stellarator magnetic equilibrium and to speed up stellarator transport solvers. Due to its speed ($\sim$ms), the AI-based surrogate enables previously unattainable optimization objectives, such as full radial profiles of ion turbulent heat flux, or directly optimizing to maximize the turbulent critical gradient at multiple locations across the plasma. These direct calculations provide a potentially more accurate optimization target and reduce reliance on ad-hoc heuristics that may not accurately capture the variation of turbulent transport with magnetic configuration. By including the AI-based surrogate for turbulent heat flux in a transport solver, we can quickly postprocess and confirm the improved ion temperature resulting from the optimized equilibrium. Finally, we demonstrate the use of AI agents with strong reasoning AI models to automate the outer loop, exploring many objective and hyperparameter configurations with this AI-based turbulence surrogate to discover improved turbulence optimized magnetic equilibria.
A rigorous implementation of energy-dependent ion- and electron-induced electron emission in a continuum-kinetic framework is used to reveal their effects on the scaling of plasma properties in the sheath with an applied bias potential to the walls. The approach comes with a novel methodology for modeling particle-induced electron emission (PIEE) that includes 1) improved fitting functions for the yield and spectra, 2) the use of SRIM and a summation of the Lindhard formula and a modified Bethe formula for obtaining accurate stopping powers, and 3) a binding energy correction to the PIEE spectra and ion-induced yield based on density functional theory (DFT) calculations. The emission models are implemented as a fully energy dependent and dynamic boundary condition. For this investigation, tungsten and graphite walls are studied for their relevance in magnetic fusion experiments. Equations are derived from fluid theory that predict the relative importance of ion- and electron-induced emission on the structure of the sheath. The simulations provide evidence for the theoretical predictions, showing that a transition from a classical to space-charge limited (SCL) sheath depends primarily on electron-induced emission. Furthermore, claims in previous literature of increased heat and particle loads to the walls due to electron emission are supported, however differing mechanisms for the increase are observed. The increased thermal and particle fluxes due to PIEE are primarily driven by collisional transfer of energy from emitted electrons in the presheath. Finally, quantitative predictions for the device modeled in this study coincide with previous modeling efforts and experimental measurements.
It is shown how evaporative cooling in hot plasma is essentially different from evaporative cooling in other media, such as neutral gases or Bose-Einstein condensates. The fundamental difference in plasmas arises both from the large mass ratio between electrons and ions in fully ionized plasma and the unusually sensitive dependence of plasma collisionality on speed. Thus, a hot ion mode plasma ($T_e < T_i$) is shown to support a distinctive evaporative cooling regime, where a strong reduction in the ion evaporation rate appears. This new regime may have application to approaches to economical nuclear fusion, where the ion tail plays an outsized role in the fusion reaction rate.
Fish-scale divertor plates are becoming popular in the design of next-generation magnetic fusion reactors to make the edge of one divertor tile sits safely behind the shadowed profile of the preceding tile, hence to protect the leading edge from extreme heat deposition. Recently, ST40 tokamak has observed an extremely narrow and peaked divertor heat-load footprint near the separatrix strike point on the fish-scale plates, on top of the usual ion-drift-width scale footprint. This observation raised concern over the severely localized burn even when the well-known ion-drift-width scale burn issue can be resolved [X. Zhang et al., Nucl. Mater. Energy 41, 101772 (2024)]. In this report we demonstrate that the lost ion gyro-orbits to the tilted tile edges can lead to a significant positive-potential sheath at a practical magnetic field incidence angle that is much greater than the usual incidence angle for the ion gyro-sheath formation on flat strike-surface. This could make the kinetic electron heat-flow spilled over from the confined region down the divertor legs to gain a significant kinetic energy amplification and allow a highly localized heat-load peak under a low edge turbulence condition. Possible intrinsic mitigation mechanism is discussed. This issue should be considered carefully in the design of next generation magnetic fusion reactors.
Fusion reactions can release a substantial fraction, and in aneutronic reactions nearly all, of their energy as the kinetic energy of fusion-product ions. In magnetized plasmas, these energetic ions can form ring distributions in velocity space. Such distributions can drive Dory-Guest-Harris (DGH) electrostatic instabilities, but these self-instabilities require a finite energetic-ion density and can be stabilized by a thermal background. We show that the same background can instead enable a distinct instability mechanism: a cyclotron pole of the minority energetic-ion susceptibility destabilizes a stable ion Bernstein wave (IBW) eigenmode. When the energetic-ion harmonic is distinct from the thermal-ion harmonics, exact root-pole resonance is thresholdless in the ideal collisionless limit. A species-resolved power balance shows that the energetic ions supply the free energy while the thermal plasma receives it. In the LAPD proton--alpha example, $99.96\%$ of the alpha-particle power loss enters the coherent proton response. This self-excited, ion-directed transfer provides a possible linear building block for alpha-particle energy channeling in proton-Boron11 fusion.
J. Griff-McMahon, V. Valenzuela-Villaseca, C. A. Walsh, S. Malko, B. McCluskey, K. Lezhnin, H. Landsberger, L. Berzak Hopkins, G. Fiksel, M. J. Rosenberg, et al.
Self-generated magnetic fields are commonly produced in high-power laser–plasma interactions. These fields can inhibit plasma heat-flow, which makes them important in inertial fusion and controlled laboratory astrophysics experiments. In this work, we characterize the time evolution of self-generated magnetic fields using multi-view proton tomography at two timings. Tomographic reconstructions of the magnetic field show a clear transition from fields located close to the target at early time to more extended coronal fields at later time. The tomographic inversion and mesh radiography also enable a direct measurement of the magnetic flux evolution. Comparisons with extended-magnetohydrodynamic simulations show only moderate agreement in field structure but good agreement in magnetic flux. This suggests that the field generation model is largely correct under these conditions, while the magnetic transport model requires additional development to reproduce the observed field structure.
Why do nonaxisymmetric stellarators avoid ballooning crashes that afflict tokamaks? Three-dimensional geometry induces Anderson localization of ballooning modes, converting a global instability into a Ginzburg-Landau network of isolated wave packets. Global stability reduces to a percolation problem: Below a critical threshold, instability is arrested; above it, a crash occurs. This explains benign stellarator saturation, predicts vulnerability in quasisymmetric designs, and introduces the critical threshold as a nonlinear stability metric for reactor optimization, pending experimental validation.
For over eighteen years, the electron current sheet measured in the Magnetic Reconnection Experiment (MRX) has stood a factor of 2--5 wider than predicted by MRX-like kinetic simulations, and the measured electron force balance has not closed with the classical terms alone. As a consequence, the dominant nonideal terms responsible for breaking the frozen-in condition have remained unexplained. Here, two-dimensional kinetic simulations with binary Coulomb collisions are performed in a cylindrical geometry representative of MRX, at the realistic hydrogen mass ratio $m_i/m_e = 1836$ and at MRX-relevant collisionality. For the first time, these simulations reproduce the measured electron current sheet half-width. The simulated value, $δ_{BT} = 0.744 \pm 0.054$~cm or $6.26 \pm 0.45$ electron skin depths ($d_e$), lies within the experimental range of 5.5--7.5~$d_e$. The electron force balance closes through the classical channels alone: the pressure-tensor divergence supports 76\% of the nonideal electric field and collisional friction the remainder. The historical force-balance deficit reappears only when the simulated layer is sampled at the experimental 3~cm outflow resolution, suggesting that the deficit reflects probe resolution rather than anomalous dissipation. Beyond this reproduction, an analytic model of the layer width is developed that orders the meandering electrons by the coherence of their orbits against collisions. In this model, the Dreicer ratio $E_D/|E_y|$ selects the electrons whose current-carrying motion survives, and the resulting width prediction brackets the measured values across a wide collisionality scan. A discrepancy remains in the width normalized to the local electron gyroradius ($ρ_e$), whose measured value lies a factor of 3--6 above both the model and the simulations.
Edge localized mode (ELM) pacing via vertical plasma oscillations—or "jogging"—has been successfully demonstrated on the DIII-D tokamak. Rapid vertical movement of the plasma toward the X-point effectively triggers ELMs. By vertically oscillating the plasma at 10 Hz, the ELM frequency increased from its natural rate of ∼5 Hz in similar DIII-D discharges to 10 Hz. Notably, downward jogs have been observed to trigger multiple ELMs in a single cycle. ELMs triggered at these higher-than-natural frequencies lead to smaller decreases in stored energy, dropping from 10% to below 1%. Consequently, the peak heat flux to the divertor is reduced by a factor of ∼2, alongside an observed reduction in carbon impurity concentration. During downward jogs in the lower single null (LSN) configuration, the X-point movement is slower and smaller than the top of the plasma, resulting in a reduced plasma cross-section and volume. To understand the mechanism of ELM triggering by jogging, a model of the edge toroidal current was developed and tested against DIII-D experimental data. Both the data and the model suggest that moving the plasma toward the X-point locally induces a net positive toroidal current in the edge region. ELITE stability analysis indicates that this induced current pushes the plasma state across the peeling side of the peeling-ballooning stability boundary into the unstable region, triggering the ELMs.
The CENTAUR Collaboration, Samuel W. Freiberger, Evan Bursch, Javier Chiriboga, Hiro J. Farre-Kaga, Eliot Felske, Sophia Guizzo, John Labbate, Shreyas Seethalla, Frederick Sheehan, et al.
This work presents the compact experimental negative triangularity reactor (CENTAUR), a low overnight cost, high-field tokamak, breakeven reactor design, achieving a predicted total fusion power of 40MW and scientific energy gain of 1.3. Ballooning stability calculations confirm that the device's pedestal is within the first stability regime, which is consistent with the expected ELM-free operation associated with negative triangularity (NT) plasmas. The geometry of the NT divertor allows for high fraction of radiated power (13.5$\%$) between the separatrix and plasma facing components. Heat transport modeling based on simulations of the edge region show heat loads into plasma facing components well below material limits. The magnet system employs rare-earth barium copper oxide (REBCO) high-temperature superconductors in 18 toroidal field coils, an hourglass-shaped central solenoid, and six poloidal field coils to support high-field ($B_0=10.9$ T) plasma confinement, shaping, and current drive. Neutronics analysis shows that a 12 cm $B_4C$ shield keeps superconducting magnet heating below the 33~K quench limit during 10 s, 40 MW DT pulses. With this shielding, the modeled fluence indicates HTS components can survive more than ten times the 3000-pulse design lifetime. Iteration of economic analysis in tandem with the technical design process allows CENTAUR to achieve its overnight cost goal of $\$$2B determined using a custom costing model that predicts a total overnight cost of $1.6$B$\pm0.2$B.
Plasma jets are formed in various astrophysical systems as plasma is rapidly ejected from a source, with a subset of these jets being magnetized and having a helical structure. Here, we demonstrate that helical jets may be formed using a ring of laser pulses that arrive on planar foils sequentially with increasing energy. The formation of the jets and their properties, including kinetic helicity, are studied through a set of three-dimensional magneto-hydrodynamics simulations with conditions informed by the parameters of the OMEGA laser facility. We find that jets with a higher degree of helicity may be generated under realistic experimental conditions when compared to a uniform jet. Synthetic x-ray and Thomson scattering diagnostics computed from simulated data demonstrate that the helical jet provides a unique fingerprint in both its morphology and plasma parameters. This laboratory helical jet platform may allow for controlled experimental study of the dynamics of helical plasma structures and, through interaction with other jets or targets, can allow for studies of shear-driven turbulence and mixing relevant to interactions between astrophysical jets and ambient clouds or crosswind.
Quasisymmetry, omnigenity and piecewise omnigenity confine trapped particles by making the bounce action independent of the field-line label. Recent optimizations produce mixed-symmetry stellarators that confine alpha particles well without them. We propose a general theory for them. From Whitham modulation theory we define iso-action, which requires only that the drift surface close and allows misalignment with flux surfaces. A solvable model supplies an exact relation between trapped segments while branch actions vary. We develop a proxy $Γ_W$ for the reach that misalignment costs.
Minseok Kim, Young-Ho Lee, SangKyeun Kim, Minwoo Kim, Sang-hee Hahn, Hiro J. Farre-Kaga, Ricardo Shousha, Juhyeok Jang, SooHyun Son, Yoon Seong Han, et al.
The edge-localized mode (ELM) frequency ($f_{\mathrm{ELM}}$) was successfully controlled in real time on KSTAR using a proportional-integral (PI) feedback controller, employing a $\mathrm{D}_2$ divertor gas puff as the actuator under tungsten lower-divertor conditions. The controller accurately tracked a two-step target---a 30 Hz increase in $f_{\mathrm{ELM}}$ for 4 s, followed by a 30 Hz decrease for 3 s---yielding mean and median absolute percentage errors of approximately 13% and 12%, respectively. Compared to a reference discharge, the actively controlled shot did not exhibit a significant drop in volume-integrated core radiation, confirming that excessive gas use merely degrades overall plasma performance. However, when contrasted with the exponential increase in core radiation observed in the absence of divertor gas puffing, these results underscore the critical need for real-time optimization. Specifically, divertor gas commands must be actively managed to maintain an $f_{\mathrm{ELM}}$ sufficient for flushing tungsten from the core while maximizing global plasma performance.
At the core of molecular dynamics (MD) simulations of plasma–surface interactions is the interatomic potential that predicts the energy and forces of atomic configurations. Recently, machine-learned interatomic potentials (MLIPs) have become popular in related fields. These MLIPs, developed for near-equilibrium calculations, are challenged when used for the relatively high-energy, chaotic conditions of plasma–surface interactions. In this paper, active learning is used to produce a large dataset of density functional theory calculations featuring C, H, O, and Ar in configurations relevant to simulations of plasma–surface interactions. These data are then used to train both an MLIP and a classical interatomic potential (reactive force field, ReaxFF) for direct comparison. Both potentials are trained using typical machine learning methods, namely, optimization of a loss function via automatic differentiation with respect to the interatomic potential parameters. Both models performed well on a test dataset, producing comparable errors. However, MD simulations using the MLIP were not consistent with published experiments. In contrast, the trained ReaxFF potential appears to perform well on these tasks. Active learning accompanied by machine-learning-style parameter fitting appears promising as a method for producing transferable interatomic potentials for simulations of plasma–surface interactions.
Laser wakefield acceleration promises compact electron accelerators for applications in medicine, industry, and fundamental science. Yet, despite rapid progress, accurately predicting the electron energy attainable in a given experimental configuration and the acceleration length required to reach it remains an open challenge. Here we use Bayesian optimization combined with advanced particle-in-cell simulation techniques to determine the maximum electron energy that a self-guided laser wakefield accelerator driven by a laser of a given energy and wavelength can produce. By systematically optimizing the accelerator performance across a range of laser energies and wavelengths, we derive energy-optimized scaling laws. These scaling laws yield the highest electron energy over the shortest acceleration length possible, are expressed solely in terms of laser energy and wavelength, and are accompanied by the complete set of laser and plasma parameters required to enable the scaling. The resulting scaling laws provide practical guidance for designing state-of-the-art laser wakefield acceleration experiments operating at their fundamental performance limits.
