Choose timezone
Your profile timezone:
→
Europe/London
HOW room (Culham Campus)
HOW room
Culham Campus
Abingdon, OX14 3DB, UK
Description
Welcome to the ECW2026 event page!
The UKAEA Exhaust Code Workshop (ECW) is a space where both developers and users of exhaust-relevant codes can share progress that might not normally make it to a conference talk or a paper: implementation details, challenges, half-started (but interesting) ideas, what didn't work and what worked well!
Building on the huge success of the Hermes-3 workshop in 2024 and its 2025 edition under the ECW name, we are now making the workshop fully open and inviting broader contributions. We will still give updates on our projects, including Hermes-3 and the GPU kinetic neutral model VANTAGE.
ECW2026 is a free event. Lunch will be provided.
There will be a conference dinner - the pre-payment page will open nearer to the event.
Who can attend?
ECW is open to anyone involved in or interested in fusion device exhaust code development. This is not exclusive to tokamaks, and we welcome submissions relevant to stellarators, dipoles or any other configurations. We want to create an open, informal and welcoming atmosphere where there are no stupid questions or stupid answers - non-experts and early-stage students are very welcome!
All external visitors must bring valid photo ID, such as a driving licence or passport.
Format
We have four focus topics this year: mean-field, turbulence, kinetic neutrals and research software engineering.
The program will include:
-
Keynote talks from world-leading experts (30mins + 10mins questions)
-
Moderated panel discussions
-
Contributed talks (15 mins + 5 mins questions)
-
Lightning talks with posters (2 slides in 4 mins + A0 poster)
The event will be hybrid apart from the poster session, but we encourage you to come in person - the most useful discussions often happen over a coffee or over a beer in the evening.
Key dates, registration, submission and other details
Please see the relevant tabs on the left hand side of this webpage.
Visas / Electronic Travel Authorisation (ETA)
Please check the UK Government website to see if you need a visa / ETA to attend this conference / training. We recommend you apply for a visa at least 16 weeks ahead of the event.
If you require a conference invitation letter to support your application, please email the support team with your full name, organisation, and nationality as stated in your passport, along with your request.
Contact
-
-
13:00
Welcome and Introduction
-
1
European effort towards edge fluid modelling tools for self-consistent reactor-relevant simulations: main outcomes and remaining gaps
Predictive modelling of plasma edge and divertor physics is a key challenge for ITER and future fusion power plants, where narrow operational margins require self-consistent treatment of turbulence, neutrals, magnetic geometry, and plasma-wall interactions. The TSVV3 project was established to advance European edge fluid turbulence modelling tools towards reactor-relevant predictive capability. We report here on its main outcomes and discuss identified remaining gaps.
The TSVV3 project coordinated developments across the main European edge turbulence codes (GBS, GRILLIX, FELTOR and SOLEDGE3X), combining advances in numerical methods, physical models and high-performance computing. The project addressed four major challenges: 1- realistic magnetic and wall geometries, 2- self-consistent neutral dynamics in detached divertor conditions, 3- turbulence modelling in high-performance plasma regimes, and 4- scalability towards reactor-sized simulations. Code developments were accompanied by validation against experiments on several European tokamaks and by applications to both tokamak and stellarator configurations.
The project enabled arbitrary axisymmetric magnetic geometries in all contributing European edge turbulence codes and extended simulations to three-dimensional configurations, including stellarators. Self-consistent neutral models based on fluid, kinetic and coupled Monte Carlo approaches were implemented, allowing first turbulence simulations in detached divertor regimes and demonstrating significant turbulence-induced scrape-off layer broadening in line with experimental results. Electromagnetic turbulence and improved collisional closures extended the applicability of fluid models towards high-β and low-collisionality plasmas, eventually enabling first self-consistent edge turbulence simulations in H-mode conditions. In parallel, substantial advances in numerical algorithms and GPU acceleration improved computational performance. Validation activities across multiple experimental devices have demonstrated the growing flexibility and predictive capability of these models.
While TSVV3 has significantly expanded the physics fidelity and applicability of European edge fluid turbulence codes, significant challenges remain on the path to self-consistent reactor-scale simulations. These include improved wall geometry treatment, scalable kinetic neutral modelling, impurity physics in three-dimensional turbulence simulations, and numerical strategies capable of addressing the spatial and temporal scales of fusion power plants. We will browse through some of these issues and highlight possible solutions currently explored in the TSVV-B project which took over TSVV3.Speaker: Dr Patrick Tamain (CEA) -
2
Numerical challenges in SOLPS-ITER modelling for STEP divertor design
We have been using the SOLPS-ITER code to design and optimise the divertor for STEP (Spherical Tokamak for Energy Production). While this modelling has yielded important physical insights — such as the impacts of fuelling locations and dome shape — it has also highlighted significant numerical challenges. These issues introduce additional uncertainties on top of existing variances driven by assumed input parameters and reduced physics models. This presentation focuses on the two most significant numerical challenges encountered: unexpected up-down asymmetries under nominally symmetric conditions, and particle balance errors. We will detail the possible causes of these numerical issues and discuss actionable strategies to achieve more robust, reliable SOL/Divertor modelling within the timescales demanded by STEP.
Speaker: Ryoko Osawa (UKAEA) -
3
Global fluid turbulence simulations in a general frame
Toroidal magnetic confinement fusion experiments require a rotational transform
to achieve adequate confinement. This rotational transform can be achieved
via any combination of internal plasma currents, torsion of the magnetic axis,
and the toroidal rotation of the plasma cross section. Internal currents can lead
to current-driven instabilities, and therefore the first stellarators employed magnetic
axis torsion in a figure-8 form to provide their rotational transform [1].
These early configurations were unstable, however, and the figure-8 stellarator
was abandoned in favor of the toroidally-shaped tokamaks and classical stellarators.
As a result, nearly every simulation and analysis framework for magnetic
fusion experiments has been written assuming a toroidal geometry. Recently,
however, advancements in near-axis theory [2] and equilibrium reconstruction
codes [3] have renewed the interest in configurations unrealizable by traditional
simulation and analysis frameworks [4].
Here we demonstrate the capability of the BSTING project [5, 6] to simulate
global fluid transport and turbulence a generalized frame. We utilize a hot-ion
turbulence model within the Hermes-3 family of multifluid models [7], which provides
an inherent flexibility of model fidelity. As a proof-of-principle, transport
characteristics of an optimized quasi-isodynamic figure-8 stellarator [4] are presented.
It is determined that the figure-8 indicates strong ballooning transport
in high-curvature regions, with radial heat and particle flux bands extending
into zero-curvature regions. A comparison of the turbulent flux with known
turbulence scalings, as well as how the flux compares to that in conventional
stellarator geometries is provided.Speaker: Dr Brendan Shanahan (MIT Plasma Science and Fusion Center) -
4
Dynamic system identification simulations of MAST-U using SOLPS-ITER
The safe and controlled exhaust of heat and particles from magnetically confined fusion plasmas requires having dynamic information about the system. On present-day reactors, this information is obtained using system identification experiments, which involve perturbing the gas injection rate and observing the Scrape-Off-Layer (SOL) response in the frequency domain. However, this strategy is cumbersome for future reactors, as accidental disruptions or reattachment of the SOL plasma could lead to significant damage to the device. Therefore, numerical dynamic models of the exhaust plasma, validated against experiments, are needed. To this end, we present a multi-sine gas injection perturbation approach that uses time-dependent SOLPS-ITER simulations, with a grid extending to the vessel wall, to quantify the plasma response. Comparison to an existing set of experimental data on MAST-U is done to validate this approach. The simulation grid was generated using the Grid Optimization and Adaptation Toolbox (GOAT), and both the Advanced Fluid Neutral (AFN) model and the kinetic neutral model Eirene available in SOLPS-ITER are evaluated. In steady state, it is possible to simulate divertor states ranging from attached to deeply detached using both fluid and kinetic neutrals. The comparison in the frequency domain shows that the AFN SOLPS-ITER simulations predict faster response times to divertor fuelling perturbations than the experiments at low frequencies, likely in part due to the fluid approximation used for the neutrals. Simulations employing time-dependent Eirene for kinetic treatment of neutrals show good agreement with the experimental data throughout the measured frequency range of around 9 to 30 Hz. Simulations with Eirene using the Quasi-Steady-State (QSS) approximation result in different dynamic behaviour due to the non-linear coupling between the plasma solver B2.5 and Eirene. Proper time-dependent kinetic treatments of the neutrals is therefore essential for accurately simulating plasma evolution in response to changes in gas injection.
Speaker: Stijn Kobussen (DIFFER, Dutch Institute for Fundamental Energy Research) -
15:00
Break
-
5
Hermes-3 simulation of plasmas confined by a levitated dipole magnet
Laboratory dipole confinement was proposed by Hasegawa in 1987 as a possible route to fusion-relevant plasmas and has been explored in levitated-dipole experiments such as LDX, RT-1, and the recent Junior device. A dipole magnetic field is characterized by its strong magnetic-flux expansion, i.e., the magnetic field strength varies strongly along a field line, which poses challenges for numerical simulations. In this work, we use Hermes-3 to study three-dimensional plasma dynamics near the edge of a dipole-confined plasma.
Fluctuations in the closed-field-line region in dipole-confined plasmas are often believed to be flute-like with $k_\parallel\approx 0$. We find that this picture can break down near a marginal stationary profile. In slightly supercritical profiles, three-dimensional simulations reveal that non-flute fluctuations with finite parallel structure naturally emerge due to the strong variation of the local magnetic geometry along the field line. Furthermore, these fluctuations produce a differential radial transport: along the same field line, particles are pinched inward near the outer midplane while being driven outward near the inner midplane. This poloidal-angle-dependent transport is not captured by purely flute-like, flux-tube-averaged descriptions. This result is further supported by eigenmode solutions from a simplified drift-reduced fluid model, whose quasilinear flux reproduces the differential radial transport. Future work will extend the simulation domain to include the X-point and open-field-line regions in a dipole field.
Speaker: Yichen Fu (Columbia University) -
6
Implementation and Analysis of Transport Simulations with External 3D Fields in Hermes-3
Despite tokamaks being primarily treated as axisymmetric beasts, many aspects of actually running such a machine are inherently nonaxisymmetric. This includes error fields, magnetic ripple, and even Edge Localized Mode-suppressing (ELM-suppressing) Resonant Magnetic Perturbations (RMPs). Each of these effects are extremely important to consider while the fusion field is moving towards next generation machines. Specifically, recent SPARC studies have addressed error field mitigation [1], and magnetic ripple-driven alpha particle losses [2], additionally many ITER studies have explored the implementation of RMPs to suppress machine-killing ELMs [3].
Upon adding 3D magnetic fields, transport is strongly affected. Loss of axisymmetry breaks integrability of field line trajectories, resulting in a chaotic equilibrium magnetic field. This chaotic field changes the structure of the differential operators within the governing transport equations, and provides additional transport in and of itself [4]. In this work, these effects are studied through a transport simulation of an RMP-applied DIII-D shot, although the methodology could nominally be applied to any other machine, and any other mode of externally-applied asymmetry. The implementation methodology is walked through, highlighting difficulties, and possible pitfalls as well as interesting physics effects which become necessary when considering a chaotic equilibrium field. Finally, the results of the study, such as the effect of 3D fields on anomalous transport coefficients, and divertor heat flux patterns are presented.
[1] S. Munaretto et.al., Nuclear Fusion, (2025).
[2] S. D. Scott et.al., Journal of Plasma Physics, (2020).
[3] L. Zhou et.al., Plasma Phys. Control., (2016).
[4] A. B. Rechester and M. N. Rosenbluth, Physical Review Letters, (1978).Speaker: Sidney Williams (University of California San Diego) -
7
Progress of the drift and turbulence capabilities of Hermes-3 for non-axisymmetric geometries
Particle drifts play a crucial role in understanding of the edge physics of current and future fusion
devices. While tokamak simulations routinely include drift effects, they have been largely neglected
in stellarator simulations until now. Advancements in BOUT++ enable the efficient simulation
of stellarator geometries. Hermes-3 is a plasma fluid model for the edge of fusion devices and
covers a wide range of levels of fidelity. In this work, we report on the progress of the Hermes3 implementation for non-axisymmetric geometries using the flux-coordinates-independent (FCI)
approach. We elucidate the implementation for steady-state transport simulations with drifts and
show a first example simulation of W7-X. Additionally, we demonstrate the turbulence capabilities
by simulating a W7-AS-like stellarator geometry. In addition to recent progress, we show some of
the remaining challenges such as neutral physics or the creation of closed divertor grids.Speaker: Tobias Tork (IPP Greifswald) -
8
Immersed Boundary Extension to BOUT++/Hermes-3
Realizing commercially viable fusion power plants (FPPs) demands whole-device modeling efforts which integrate high‐fidelity plasma simulations with engineering systems and design capabilities. Among the most pressing challenges are the accurate prediction of heat and particle exhaust in the divertor region and the self‐consistent treatment of plasma–neutral and plasma–wall interactions in the edge and scrape‐off layer. In this work, the BOUT++/Hermes-3 [1] fluid framework is extended with a Cartesian flux-coordinate-independent (FCI) discretization [2], allowing for full-device transport and turbulence studies from the core to the first wall. To fully characterize the plasma, the device wall is represented as an immersed boundary [3], from which boundary conditions are imposed on the perpendicular anomalous diffusion and ExB advection terms. For terms employing a finite volume method, plasma cells intersected by the boundary are further treated using a cut-cell discretization [4]. The parallel discretization and associated boundary conditions are presently inherited from the quasi-FCI implementation BSTING [5].
This work was performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344. LLNL-ABS-2020022
[1] B. Dudson, M. Kryjak, H. Muhammed, P. Hill, and J. Omotani Comp. Phys. Comm. 296 108991 (2024)
[2] F. Hariri and M. Ottaviani Comp. Phys. Comm. 184 2419-2429 (2013)
[3] R. Ghias, R. Mittal, and H. Dong Journal of Computational Physics 225 528–553 (2007)
[4] H. Johansen and P. Colella Journal of Computational Physics 147, 60-85 (1998)
[5] B. Shanahan, B. Dudson, and P. Hill Plasma Phys. Control. Fusion 61 025007 (2019)Speaker: Stefan Tirkas (Lawrence Livermore National Laboratory) -
17:00
Day 1 Close
-
13:00
-
-
08:55
Day 2 Open
-
9
On the challenge of kinetic neutrals modelling in reactor relevant conditions: characterization of the issues and possible solutions
The particle, momentum and energy exchanges between the plasma and the neutral gas from recycling are key aspects of divertor physics. Atoms and molecules are often not in a collisional regime and a kinetic description is required to get transport right. In detached regimes, where the temperature in the divertor drops below a few eVs, the complexity of the reaction channels for neutrals increases. The Monte Carlo approach implemented in EIRENE allows to solve the kinetic problem on a 3D grid, and can accommodate both the geometrical complexity and the relevant species and reactions channels. However, in large machines such as ITER, in which the mean free path of neutrals can be very short compared to the size of the divertor, regions where neutrals become collisional can appear. Such regions make the Monte Carlo approach computationally much more costly, since atoms or molecules may undergo tens of thousands of collisions before being ionized or pumped. It also makes fluid approaches accurate at least in these regions. This presentation will focus on EIRENE simulations in ITER conditions, showing that the computational cost increase is related to physics and that sacrificing the physics is not an option if the goal is to assess e.g. peak heat fluxes in the divertor. Advanced fluid models, fully consistent with the underlying kinetic model used in EIRENE, have been developed for atoms and shown promising results. However, they are generally not valid everywhere in the simulation domain and hybrid kinetic-fluid models are shown to provide a way to combine accuracy and efficiency.
Speaker: Dr Yannick Marandet (Aix-Marseille University / CNRS) -
10
The KINetic Deterministic NEutral Solver for turbulence simulations in the boundary of magnetic confinement devices
We present KINDNES (KINetic Deterministic NEutral Solver), a tool for solving the kinetic Boltzmann equation applied to neutral particle dynamics in the boundary plasma of tokamaks and stellarators. Originally developed as the neutral model within the GBS turbulence code, KINDNES is now developed as a standalone library.
Neutral particles play a central role in edge plasma physics, governing key processes such as recycling, momentum and energy exchange, and the onset of plasma detachment. While existing neutral models typically rely either on fluid approximations or on Monte Carlo kinetic methods, KINDNES provides a fully kinetic, deterministic framework, free of statistical noise. The model discretizes the Boltzmann equation for each neutral species and integrates it along its characteristics, reducing the problem to the inversion of a discrete linear system.
The fundamental physical mechanisms and numerical implementation of KINDNES are first presented, followed by its first benchmark against the Monte Carlo code EIRENE. Using a deuterium plasma background derived from a SOLPS simulation, the two codes are compared in an attached divertor configuration, with EIRENE configured to match the physical assumptions currently available in KINDNES. Within this common physics framework, the two codes show good agreement in atomic deuterium density, temperature and velocities, representing an encouraging validation of KINDNES as a deterministic alternative for neutral particle modeling.
Building on this validation, the most recent developments of KINDNES are then presented: the extension to molecular deuterium dynamics, enabling the simulation of detached divertor conditions; the introduction of flexible wall geometry support, broadening the range of device configurations that can be modeled; the development of a full 3D solver, that includes the proper tratment of toroidal velocities; the relaxation of the adiabatic approximation, rendering the solver aware of the plasma evolution in time and the adoption of a Hierarchical Matrices method to substantially accelerate the solution of the linear system in terms of both memory and computational cost.Speaker: Davide Mancini (EPFL) -
11
Using SPARTA as a standalone DSMC test bed for divertor neutral transport and exhaust boundary conditions
Kinetic neutral transport remains one of the least reducible components of fusion exhaust modelling. While coupled plasma-edge tools such as SOLPS-ITER/EIRENE remain the standard route for integrated divertor simulations, standalone kinetic-neutral calculations provide a complementary way to isolate rarefied-gas effects, test boundary-condition assumptions, and interpret neutral transport patterns without the full complexity of plasma feedback.
In this contribution, we discuss the use of SPARTA, a direct simulation Monte Carlo code for rarefied-gas dynamics, as a standalone tool for divertor-relevant neutral transport studies. The emphasis is not on replacing established plasma-edge workflows, but on identifying where a DSMC treatment can provide useful diagnostic and modelling insight. We first outline the practical advantages and limitations of this approach relative to EIRENE-like neutral models, including the use of binary neutral--neutral collisions rather than BGK-type closures, the treatment of neutral--wall interaction models such as CLL and TRIM-informed reflection, and the construction of plasma-mimicking boundary conditions from prescribed neutral sources and fluxes.
We then present a set of simplified private-flux-region and toroidal-plane test cases motivated by divertor exhaust applications. These cases show that neutrals need not behave as an isotropic diffusive background: ordered sources, wall interaction, geometry, and neutral--neutral collisions can generate structured flows, preferential capture, backflow, and non-trivial momentum transport. Such behaviour has direct implications for reduced exhaust models, including the common treatment of pump surfaces through fixed albedo or absorption coefficients. In particular, DSMC calculations allow pump capture, return flux, and effective exhaust conductance to be estimated from the local rarefied flow rather than prescribed as constants.
The talk will focus on lessons for kinetic-neutral and exhaust-code developers: which assumptions are exposed by standalone DSMC, which quantities can be passed back to integrated edge codes, and where systematic comparison against EIRENE-like models is still needed.
Speaker: Matteo Moscheni (Gauss Fusion GmbH, Parkring 29, 85748 Garching bei München, Germany) -
10:20
Break
-
12
A fluid model for the deuterium molecules in the plasma edge based on the AMJUEL database and benchmark with EIRENE
In plasma boundary modelling, fluid neutrals offer several practical disadvantages compared to their more accurate kinetic Monte Carlo counterparts. They are typically much faster, and their deterministic nature is a major benefit when doing optimization or when coupling to plasma turbulence codes.
Besides missing kinetic effects, fluid neutral models are typically also lagging behind with regard to the physics processes they include. The vast majority of present-day fluid neutral models consider only atoms.
The first commonly referenced fluid molecular models in the context of plasma edge models were developed in UEDGE. The UEDGE work focused primarily on the details of the molecular physics (Collisional Radiative Models) and their impact on detachment in DIII-D, while the benchmark with EIRENE was more restricted. In this contribution, we aim for a complementary approach, considering only the default molecular physics model from EIRENE, but performing for the first time a thorough benchmark between the fluid model and EIRENE for all macroscopic parameters (density, velocity and temperature) and all exchange terms with the plasma and atom populations (ion/atom particle sources, ion/atom parallel momentum sources, electron energy sources and ion/atom energy sources), on a plasma background representative of the detached regime.
We show that the particle velocity distribution of the kinetic molecules approaches a Maxwellian in a large part of the relevant domain, due to the elastic collisions with the detached plasma background. We also show that the bulk of the molecular physics can be captured reasonably well with the presented fluid molecular model, showcasing it to be an excellent starting point for further research in this direction.
Speaker: Wim Van Uytven (KU Leuven (EU)) -
13
Propagator-based Multi-level Monte-Carlo for Kinetic Neutral Species
In recent work, we introduced a novel multi-level Monte-Carlo scheme for kinetic simulations of neutral particles in plasmas near the isotropic cross-section limit. Of particular interest is the case of edge plasmas in detatched or near-detached regimes. The multi-level scheme in this case is based on a propagator formalism enumerating discrete charge-exchange collision events. This scheme reproduces the results of known algorithms in the large particle number and fine mesh limit. More importantly, the new formalism shows promise for differentiable and deterministic particle-based models. At the end, we will discuss preliminary work implementing this method for fully-implicit coupled plasma-kinetic neutral models. This talk/poster is based on joint work with Maxim Umansky and Ben Dudson.
Speaker: Gregory Parker (UC Berkeley) -
14
NESO-Particles: A Performance Portable Particle Framework for Fusion Plasma Simulation
Simulation of fusion plasma often requires representations of non-Maxwellian distributions. A common method is to represent the distribution function as a collection of individual particles where the state of each particle evolves over time. In this talk we describe NESO-Particles (NP) which is a performance portable framework for describing particle based algorithms.
NP provides users with abstractions for describing particle data and looping operations for the fusion plasma use case. These abstractions form a separation of concerns between the description of particle based algorithms and the execution of these algorithms on modern HPC hardware. Our implementation is realised as a MPI+SYCL library that executes on both CPU and GPU hardware.
Speaker: Will Saunders (UKAEA) -
12:00
Lunch
-
12:45
📸 Group Photo
-
15
The kinetic neutral model in MHD code JOREK
A few years ago, the nonlinear extended MHD code JOREK [2] has been extended with a dedicated kinetic neutral model [3], opening up the possibility to transiently model the interplay between edge-plasma and plasma-wall interactions and the plasma processes for which JOREK is well known such as ELMs, disruptions and RE beams, where previously these worlds could not be self-consistently combined in any code. Since then, the kinetic neutral model has been further developed and used as exhaust code by a growing part of the JOREK community and has attracted increasingly more attention from the exhaust physics field.
The MHD model in JOREK uses a flux surface aligned poloidal finite element grid up to the true first wall with the option to model variation in the toroidal direction using Fourier expansion. The kinetic neutral model is built on JOREK’s hybrid OMP-MPI particle-in-cell framework [4] and follows an array of variable weight particles through time, accumulating sources and sinks that are fed into the plasma equations. Both the MHD side and the kinetic side are inherently suited and optimized for time dependent modelling, making the model interesting for investigating dynamic physics and for control. Included in the neutral model physics are ionisation, charge exchange, wall – and volume recombination, line radiation and puffing, while molecular interactions are currently still missing.
As part of recent developments, various kinetic applications in JOREK (neutrals, impurities, runaway electrons and fast particles) have been combined into a single flexible kinetic framework, with object oriented control of the simulation setup, such that the included particle species, the particle pushers used for each species and the kinetic-MHD coupling scheme can be selected. Additionally, the neutral model is now compatible with a two temperature plasma description; direct binary neutral self-collisions have been added; and capability to resolve simulation actions (such as neutral collisions) on timescales in between the MHD timestep and
the kinetic particle timestep have been implemented. Furthermore, there is work underway for GPU acceleration and making the kinetic neutral model compatible with stellarator geometry.
This talk will discuss JOREK’s kinetic neutral model, how it is implemented and what is and is not currently included in the physics, along with details on the recent upgrades, encountered issues, and a selection of application and validation simulations of small ELM burn-through, XPR formation and transition into a MARFE, QCE access, and dynamic detachment studies.
[1] Hoelzl, Nucl. Fusion 66 (2026) 116006.
[2] Hoelzl, Nucl. Fusion 61 (2021) 065001.
[3] Korving, Phys. Plasmas 30 (2023) 042509.
[4] van Vugt, PhD Thesis, TU Eindhoven, 978-90-386-4811-8 (2019).
Speaker: Daniël Maris (DIFFER) -
16
VANTAGE-Reactions – performance-portable neutral reactive physics for exhaust modelling
UKAEA (United Kingdom Atomic Energy Authority), Culham Campus, Abingdon, Oxfordshire, OX14 3DB, UK.
The accurate simulation of neutral particles in the tokamak exhaust has been a topic of research across many decades, with the primary numerical approach being Markov Chain Monte Carlo methods, widely used in linear transport problems[1,2]. Historically, neutral-neutral collisions have been treated using the BGK method in order to fit in with the linear transport constraints[3], with reduced models recently gaining traction[4].
Recent work at UKAEA has been focused on performance-portable particle libraries built on top of NESO-Particles[5], with VANTAGE-Reactions[6] providing abstractions and general implementations for plasma-neutral, particle-surface, and neutral-neutral interactions. It is actively being integrated into finite volume (Hermes-3[7]) and finite element (PENKNIFE[8]) codes. Following a less common path to kinetic neutral modelling for tokamaks[9,10,11], we build on methods from rarefied gas dynamics, specifically the Direct Simulation Monte Carlo method[12], with a focus on performance-portability and flexibility in both steady state and time-dependent scenarios. In this talk, we will present the general weighted particle approach used by VANTAGE-Reactions, benchmarking results so far, as well as early progress in development of binary collision support through an event-splitting stochastic weighted particle method[13,14].
This work has been part-funded by the EPSRC Fusion Grant 2022/27 [grant number EP/W006839/1].
References:
[1] D. Reiter, Journal of Nuclear Materials, 196–198 80–89 (1992)
[2] D.P. Stotler et al. “DEGAS 2 neutral transport modeling of high density, low temperature plasmas.” (1997).
[3] V. Kotov, et al., Plasma Physics and Controlled Fusion, 50 10 (2008)
[4] D. V. Borodin et al., Nuclear Fusion, 62 8 (2022)
[5] https://github.com/ExCALIBUR-NEPTUNE/NESO-Particles
[6] https://github.com/UKAEA-Edge-Code/VANTAGE-Reactions
[7] B. Dudson et al. Computer Physics Communications 296 2024
[8] https://github.com/ExCALIBUR-NEPTUNE/PENKNIFE
[9] S. Varoutis et al., Fusion Engineering and Design, 121 13–21 (2017)
[10] S. Q. Korving et al, Physics of Plasmas, 30 4 (2023)
[11] K. Kvist et al., Physics of Plasmas, 31 3 (2024)
[12] G. Bird, Molecular Gas Dynamics and the Direct Simulation of Gas Flows, Oxford University Press (1994)
[13] S. Rjasanow et al., Journal of Computational Physics, 124. 2 243–253 (1996)
[14] G. Oblapenko et al. Journal of Computational Physics 466 (2022)
Speaker: Stefan Mijin (UKAEA) -
⚡ Lightning Talks / Poster Sessions: Set A
-
17
Kinetic Electron Effects in Expanding Scrape-Off-Layer Flux Tubes
Magnetic flux expansion in the scrape-off layer (SOL) is a critical mechanism for spreading exhaust heat loads across divertor target surfaces, and its accurate modelling is essential for predicting performance in present and next-generation devices. This phenomenon is typically modelled using fluid codes, which incorporate flux expansion geometrically through area factors in the governing equations [1,2]. However, fluid models often rely on assumptions such as those in the Braginskii closure [3], which break down in the plasma edge. To accurately capture the non-local and non-equilibrium behaviour displayed by the divertor plasma during detached conditions and transient events like edgelocalised modes (ELMs), kinetic treatment of the SOL is required [4,5,6].
A finite volume discretisation of the kinetic electron equation that includes magnetic flux expansion has been developed and implemented within the ReMKiT1D framework, a reduced 1D2V Vlasov-Fokker-Planck (VFP) code [7]. Unlike higher-dimensional models, where these effects emerge naturally from the coordinate geometry, the 1D reduction requires their explicit reconstruction at the level of the governing terms and their spatial discretisation. This approach shares the mathematical structure of the cosmic ray transport formulation of Bell et al. [8]. However, while that work embeds the method within a full multi-physics model, here the contribution of the kinetic flux expansion terms is isolated and validated numerically via velocity-space moment analysis of the electron distribution function.
To clearly expose the role of kinetic electrons, the model intentionally employs simplified treatments of ions and neutrals so that any observed effects can be attributed to non-Maxwellian electron dynamics. Within this controlled framework, systematic parameter scans over flux expansion factor and input power for both steady-state and transient scenarios reveal how kinetic electron effects modify SOL profiles and target heat flux under conditions relevant to detachment and ELM-driven transients. Beyond characterising these effects, this implementation enables future work on the development of reduced kinetic–fluid coupling strategies for transient SOL modelling.
References
[1] Dudson, B. D. et al. PPCF 61(6) (2019)
[2] Derks, G. L. et al. PPCF 64 (2022)
[3] Braginskii, S. I. Rev. Plasma Phys. 1, 205 (1965)
[4] Tskhakaya, D. et al. Contrib. Plasma Phys., 48(1–3), 89–93 (2008)
[5] Chankin, A. V. et al. PPCF 60 (2018)
[6] Mijin, S. et al. PPCF 62(9) (2020)
[7] Mijin, S. et al. Comput. Phys. Commun., 300 (2024)
[8] Bell, A. R. et al. MNRAS 539, 1236–1247 (2025)This work has been part-funded by the EPSRC Energy Programme [grant number EP/W006839/1] and used the ARCHER2 UK National Supercomputing Service as well as resources provided by the Cambridge Service for Data Driven Discovery (CSD3) operated by the University of Cambridge Research Computing Service.
Speaker: Vera Oberhauser (UKAEA) -
18
Electromagnetic Fluid Turbulence Simulations of Edge Tokamak Plasmas with Hermes-3
Hermes-3 is a multi-component plasma fluid code based on the BOUT++ framework [1] which simulates turbulence in the edge region and scrape-off-layer of tokamak plasmas with full 3D diverted tokamak geometry. The effect of turbulent magnetic fluctuations due to local turbulent currents has been implemented in Hermes-3, including ‘magnetic flutter’, how local magnetic fluctuations impact the evaluation of parallel gradients.
Flux-driven turbulence simulations were performed with evolving electron and ion density, momentum, and pressure. The scenario is based on a DIII-D plasma with magnetic fluctuations observed in the dynamics of L-H and H-L transitions [2]. The plasma evolves self-consistently from a given grid and input power to the core. The model has upper single null diverted geometry, unfavourable ion grad-B drift direction, deuterium fuel, toroidal magnetic field BT = 2.0T, plasma current IP =1.0MA, major radius R = 1.72m, and minor radius a = 0.6m. The radial extent of the domain is between normalised radius 0.9<𝜌<1.07, where 𝜌 is the normalised toroidal flux radial coordinate.
Three simulations were performed with varying levels of electromagnetism active in the model: an electrostatic simulation (ES), an electromagnetic simulation without flutter enabled (EM), and an electromagnetic simulation with flutter enabled (Flut). A comparison of the saturated equilibria shows EM has increased turbulent transport than ES, and the Flut simulation shows that inclusion of magnetic flutter reduces this turbulent transport back to levels comparable to the ES simulation.
Magnetic flutter has been seen to be crucial for the GRILLIX team to model their H-mode-like edge profiles [2]. The implementation of flutter in Hermes-3 paves the way for H-mode plasmas and L-H transitions to be simulated.
Speaker: Mr Tom Ashton-Key (Imperial College London) -
19
The effect of fluid neutrals on 1D detachment burn-through for STEP
\noindent The STEP program run by UK Fusion Energy Ltd (UKFE, UKAEA group) aims to deliver a tritium self sufficient prototype fusion power plant generating 100 MW net electric power based on the spherical tokamak concept. The compact spherical geometry of STEP's SPP-2 design raises significant exhaust challenges, with powers crossing the separatrix reaching $P_\text{sep} \approx110 - 140\,\text{MW}$, corresponding to $P_\text{sep}/R_0 \approx 25-32\,\text{MW/m}$. To overcome these challenges advanced solutions, such as a double null with an extended outer leg divertor, will be employed to achieve sufficient detachment and protect plasma-facing components [1]. While these measures are predicted to achieve acceptable steady-state power to the divertor targets, transient power loads from fluctuations in core fusion/radiation power or vertical displacement-induced disconnected double null effects may lead to burn-through of the detachment front, and unacceptable power loading to material surfaces.
High fidelity (SOLPS-ITER) simulations of transient burn-through remain prohibitively expensive to explore the wide parameter space of reactor scale devices. The main cause of this expense is the kinetic treatment of neutral species. To provide initial screening for STEP and comparison with analytical models [2], it is attractive to adopt the more computationally tractable method offered by reducing dimensionality to 1D and adopting a fluid treatment of the neutral species.
We use the multi-fidelity Hermes-3 fluid code to study power transient burn-through in STEP-relevant 1D exhaust scenarios, including a netural reservoir model based on the form implemented in DIV1D [3] to capture neutral cross field transport. Results show that while the inclusion of neutral reservoirs significantly improves steady state agreement to 2D SOLPS-ITER simulations, strong reservoir action can lead to non-physical to detachment front transient response. In addition to this, at the point of reattachment, strong recycling leads to extreme buffering of the 5eV temperature front, requiring strong pumping to resolve. Achieving physical results with 1D fluid neutrals therefore requires careful consideration of both reservoir implementation and target recycling. We assess both effects against available MAST-U transient SOLPS-ITER simulations.
Speaker: John Lloyd Baker (University of York, UKAEA) -
20
Benchmarking 1D plasma-neutral transport between a ReMKiT1D exhaust code and SOLPS-ITER
We present the benchmarking of one-dimensional (1D) simulations comparing a newly developed 1D exhaust code RMK_REX [1] with the modelling suite SOLPS-ITER [2]. RMK_REX is developed in the ReMKiT1D framework [3] for fluid-kinetic and collisional-radiative modelling, with the aim of simulating detached plasmas with non-local electron effects. The benchmark serves to calibrate RMK_REX to SOLPS-ITER as these features are developed, using grids and equilibria from the latter to initialise the former. Boundaries in SOLPS-ITER are configured to mimic elements of the 1D treatment including the neutral recycling and the inherent absence of side walls.
Three cases are examined in a linearly expanding flux tube geometry: a plasma-only case, and two plasma-atom cases without and with a cross-field pressure diffusion contribution to the 1D parallel neutral transport. The plasma-only case gives very close agreement, differing only slightly at the downstream boundary. In the plasma-atom case without the cross-field approximation, there are qualitative similarities in electron density, temperature and energy flux, however the target plasma is hotter and less dense, and energy content is overestimated compared to the SOLPS-ITER data. Introducing the cross-field approximation to neutral transport improves agreement in the density and temperature profiles, while underestimating ion and neutral energies compared to SOLPS-ITER. The neutral pressure downstream of the front is found to be much more sensitive to the energetic efficiency with which neutrals reflect from side walls than the target boundary condition itself. Results from this benchmark will help identify improvements and limitations of the 1D reduction of neutral transport as a step towards more complex time-dependent simulations of the detached exhaust.
Finally, future work seeks to incorporate kinetic electron transport with flux expansion. Kinetics alone will introduce deviations from the classical fluid transport coefficients [4], interactions with neutrals [5] and the boundary fluxes [6]. These deviations are likely to be significant for energetic transients in current experiments like MAST Upgrade, or during steady state in plant scale devices.
[1] Leigh, Moulton and Mijin (2026); 27th International Conference on Plasma Surface Interactions, Regensburg, Germany. Submitted to J. Nuclear Materials and Energy.
[2] Wiesen et al 2015 Journal of Nuclear Materials 463 480–484
[3] Mijin et al 2024 Comp. Phys. Comms 300 (2024) 109195
[4] Power et al 2021 Eur. Phys. J. Plus, 136
[5] Chankin et al 2018 Plasma Phys. Control. Fusion 60
[6] Mijin et al 2020 Plasma Phys. Control. Fusion 62 095004Speaker: Sid Leigh (UKAEA) -
21
Adding fluid molecules to Hermes-3.
Molecules are expected to play a critical role in the physics of detachment in tokamak divertors.
Fluid modelling of molecular species has been a longstanding and conspicuous absentee from the Hermes-3 feature list, but recent changes to the code have greatly simplified its implementation.
The handling of chemical reactions has been substantially overhauled, introducing a framework that automatically computes many of the associated source terms and allows rate data to be retrieved from arbitrary sources. These improvements have reduced the development cost of adding new reactions, paving the way for a variety of mechanisms to be added.
Here, we report on the addition of the widely-used Kotov mechanism and examine the possibility of implementing more advanced reaction networks, as well as verifying results using other exhaust codes and validating against data from devices likes MAST-U.Speaker: Owen Parry (UKAEA) -
22
First Snowflake Divertor Simulations in Hermes-3
The divertor problem is one of the biggest ongoing challenges in the fusion energy sector. The magnetic geometry and topology of the divertor region in tokamaks significantly influences plasma edge dynamics as well as heat and particle exhaust.[1] Advanced divertor configurations, like the snowflake (SF) divertor, offer advantages in mitigating heat loads and enhancing plasma performance.[2] However, simulating these geometries can be extremely challenging due to the challenging meshing constraints and the complex physics of the edge region.
BOUT++ is a widely used open-source framework for simulating plasmas at the edge regions. [3] However, its capability to handle complex divertor geometries and topologies was limited by the previous mesh, which was unsuitable for non-standard divertor geometries. To address this, we have upgraded BOUT++'s mesh, thus giving it the ability to simulate complex divertor configurations, including the snowflake divertor. The upgraded mesh allows users to define complex divertor geometries (up to two X-Points) with flexibility and precision. We have also expanded INGRID[4] (a meshing code which can already handle complex divertor topologies) to produce BOUT++ grid files, which gives us a new tool that can already create grids for most SF geometries. Hermes-3 is one of the leading plasma edge codes; it is a physics model implemented on the BOUT++ framework. It has the ability to create 1, 2, and 3D simulations, making it the ideal candidate for understanding the SF. Ideal SFs, however, are hard to maintain experimentally, and they can easily evolve into another category of SF (SF+, and SF-) depending on where the secondary X-Point is. So we have started a study on the sensitivity of heat and particle fluxes at the divertor target to the X-point separation; this distinguishes the topology from ideal, SF+, and SF- configurations by implementing the upgraded BOUT++ mesh into Hermes-3.
The expanded capability to simulate complex divertor geometries within BOUT++ and Hermes-3 opens new ways of investigating advanced divertor concepts and optimizing divertor design for future fusion reactors. This development provides a valuable resource for theoretical and computational exploration of plasma behavior in complex divertor configurations. Understanding the behaviour and advantages of the advanced divertor geometries will be crucial for any fusion powerplant in the future.Bibliography:
[1] A. Loarte, “Effects of divertor geometry on tokamak plasmas,” Plasma Phys. Control. Fusion, vol. 43, no. 6, p. R183, Jun. 2001, doi: 10.1088/0741-3335/43/6/201.
[2] D. D. Ryutov and V. A. Soukhanovskii, “The snowflake divertor,” Phys. Plasmas, vol. 22, no. 11, p. 110901, Nov. 2015, doi: 10.1063/1.4935115.
[3] B. Dudson et al., BOUT++. (Oct. 10, 2025). Zenodo. doi: 10.5281/zenodo.17313945.
[4] B. Garcia, M. Umansky, J. Watkins, J. Guterl, and O. Izacard, INGRID: an interactive grid generator for 2D edge plasma modeling. 2021. doi: 10.48550/arXiv.2102.07040.Speaker: Sebastian Ruiz Gonzalez (University of York) -
23
Modelling 2D edge transport in MAST-U with Hermes-3
Understanding the mechanism of edge transport in magnetic confinement fusion devices is critical for future reactors, as it determines how effectively high heat loads can be exhausted to achieve the sustainability of the devices.
In this work, the multi-fidelity exhaust code Hermes-3 [1] is employed to model scrape-off layer (SOL) transport using a fluid framework, solving the Braginskii equations [2] for the plasma and the advanced fluid equations for neutrals [3,4]. This approach contrasts with conventional exhaust modelling codes, such as SOLPS-ITER, which describe neutral transport kinetically at high computational cost. Replacing kinetic neutral models with fluid descriptions could accelerate simulations, whereas it is valid only in regions dominated by high charge-exchange collisions. Despite this limitation, fluid neutral models remain attractive due to their efficiency. This work compares fluid and kinetic treatments of neutral particles and assesses their impact on modelling in MAST-U.
[1] B.Dudson, M.Kryjak, H.Muhammed, P.Hill, J,Omotani Hermes-3: Multi-component plasma simulations with BOUT++ Comp. Phys. Comm. 2023 108991
[2] Braginskii S.I. 1965 Transport processes in a plasma Reviews of Plasma Physics vol 1
[3] W. Van Uytven, et. al. ‘Assessment of advanced fluid neutral models for the neutral atoms in the plasma edge and application in ITER geometry’, Nucl. Fusion, vol. 62, no. 8, p. 086023, June 2022
[4] et. al. N. Horsten, ‘Development and assessment of 2D fluid neutral models that include atomic databases and a microscopic reflection model’, Nucl. Fusion, vol. 57, no. 11, p. 116043, Aug. 2017
Speaker: Lin Shih (University of York) -
24
Drifting cold-ion flow boundary condition at the magnetic presheath entrance
When using a fluid model to simulate the plasma in the Scrape-Off Layer, boundary conditions are imposed at the entrance of the boundary layer next to the solid targets, where the model assumptions fail. This boundary layer reflects most electrons away from the target. With finite ion temperature, the largest length scale in the boundary layer is the collisional mean free path (along the magnetic field), associated with the collisional layer [1]. In the cold-ion limit, the largest length scale is instead the ion sound gyroradius, associated with the magnetic presheath [2]. Here, the ion polarisation drift grows significantly as the electric field strongly increases towards the wall, until the ion flow is bent across the magnetic field towards the target.
In the absence of gradients tangential to the target, the Bohm-Chodura condition requires that cold ions enter the magnetic presheath flowing at the sound speed parallel to the magnetic field [2,3]. However, due to the small angle between the magnetic field and the target, the small transport in the direction normal to the target caused by tangential gradients perpendicular to the magnetic field competes with the small projection of the ion parallel streaming [4]. This effect is included in some fluid codes via perturbative corrections to the sonic parallel flow [5]. Here, we explore the cold-ion flow boundary condition when cross-field transport causes large corrections, as is often the case.
[1] M. Abazorius, University of Oxford, PhD thesis (2025).
[2] R. Chodura, Phys. Fluids 25(9), 1628-1633 (1982).
[3] K.-U. Riemann, Physics of plasmas 1(3), 552-558 (1994).
[4] P. C. Stangeby, A. V. Chankin, Physics of Plasmas 2(3), 707-715 (1995).
[5] J. Loizu, P. Ricci, F. D. Halpern, S. Jolliet, Physics of plasmas 19, 12307 (2012).Speaker: Alessandro Geraldini (University of York) -
25
The collisional layer: connecting the bulk plasma to the kinetic sheath
In magnetic confinement fusion devices, the divertor is often negatively charged because electrons are more mobile and are preferentially captured by the divertor. This negative charge is shielded by a positively charged Debye sheath. In magnetised plasmas, another layer of width of the order of the ion gyroradius, known as the magnetic presheath, forms on top of the Debye sheath. In this work, we consider a simplified 1D problem in space, focusing on a distribution function that only has gradients perpendicular to the wall. We show that for a sufficiently collisional plasma, there exists another layer with a width that scales with the ion mean free path, the collisional layer. This layer connects the plasma far away from the wall, where fluid equations are used to model the plasma, to the collisionless magnetic presheath. The collisional layer is described by the steady state electrostatic ion drift kinetic equation in one spatial dimension, together with quasineutrality and adiabatic electrons. We show that the kinetic Chodura condition must be satisfied at the magnetic presheath entrance and that the potential at the magnetic presheath entrance diverges. Finally, we introduce a semi-Lagrangian finite element code CLOVER (Collisional Layer Solver) developed to solve our system of equations numerically. This code provides the distribution function at the magnetic presheath entrance and the potential drop across the collisional layer.
Speaker: Mantas Abazorius (United Kingdom Atomic Energy Authority) -
14:40
Poster Session / Break
-
17
-
26
Hermes-3 performance and neutral model update
Hermes-3 is an open-source fluid edge plasma code written using the BOUT++ framework. Building on the features of SD1D, Hermes-2 and STORM, it enables self-consistent simulation of scrape-off layer plasmas in 1D, 2D and 3D for both steady-state and unsteady/turbulent regimes. In this work, we present an update on several activities contributing towards improved neutral model capability in the code.
Hermes-3 is a fully implicit code. In 2D mean-field, it relies on the SNES solver in PETSc with backward Euler time integration and an LU-type preconditioner. The implementation of separate flux limiters for cross-field neutral advection, conduction and viscosity led to severe performance degradation. Upon investigation, the limiter implementation was found to lead to several numerical problems due to instances of non-differentiability and insufficient regularisation. Resolving these has led to a marked speedup and improvement in robustness across a range of levels of limiter saturation.
Following this, implementing a fully differentiable slope limiter has led to reduced solver failures, steadier timesteps and improved performance. A further significant speed enhancement was found by selectively lagging terms in the neutral flux limiter formulation, revealing their large impact on the stiffness of the system.
We also present an update on coupling Hermes-3 to VANTAGE (Versatile, Accurate Neutral Transport and Analysis for GPU Execution), a new performance-portable kinetic neutral model. Mass conservation is demonstrated in a simple slab simulation with an ionization-recombination balance between the fluid plasma and kinetic neutrals.
This work has been funded by the Fusion Futures Programme. As announced by the UK Government in October 2023, Fusion Futures aims to provide holistic support for the development of the fusion sector.
Performed in part under the auspices of the U.S. DOE by LLNL under contract DE-AC52-07NA27344.
Speaker: Mike Kryjak -
16:00
💬 Discussion - Kinetic neutrals and reaction physics
-
17:00
Day 2 Close
-
08:55
-
-
08:55
Day 3 Open
-
27
SOLEDGE-HDG: a novel unstructured approach for transport simulations in tokamak-relevant conditions
Fluid transport simulations based on mean-field quantities, with the averaged out turbulence fluctuations being modelled, remain a standard in the international community when considering engineering studies of tokamak relevant configurations (for a review of state-of-the-art transport codes, see Schwander et al. 2024).
In the perspective of enhancing computational efficiency and expanding codes capabilities, we started developing SOLEDGE-HDG, a high-order, finite-element code within the SOLEDGE code suite approximately ten years ago [Giorgiani et al. 2018 ; Capasso et al. 2025].
The hybridization of the DG method enables a special choice of the numerical traces (or numerical fluxes) for spatial discretization, which makes the HDG stand out thanks to its stability features, reduced number of degrees of freedom, and super-convergence properties (the error converges faster than the degree of the polynomials used to represent approximate the solution) [Giorgiani et al. 2018]. The innovation of this HDG method for tokamak simulations lies in its mesh flexibility which eliminates the need for alignment with magnetic field lines or flux surfaces. This enables accurate discretization of complex plasma-facing components and singularities (e.g., X-points) as well as the simulation of non-steady phases (e.g., start-up) without costly remeshing, addressing a gap left by current codes. Relieving the constraint of aligment with magnetic field comes at the price of additional spurious numerical diffusion in the perpendicular direction that can be bounded as long as high-order interpolations are used [Giorgiani et al. 2020]. The HDG scheme also promotes implicit time integration that relies on rapid Newton-Raphson convergence, thus allowing large time steps and fast computational access to steady states, leading to superior performance of SOLEDGE-HDG in achieving 2D transport equilibria when compared to semi-implicit codes.
In this presentation, we will detail the algorithm, highlighting both its attractive properties and its current limitations. Over recent years, the code has been successfully applied to simulate tokamak-relevant configurations [Scotto et al. 2022, Kudashev et al. 2026]. Selected results will be presented to illustrate the code’s new capabilities and its strong potential for addressing transport and heat exhaust challenges in existing devices, as well as in next-generation machines such as ITER [Scotto et al. 2024] or SPARC.Speakers: Dr Frederic Schwander (Aix-Marseille University / M2P2), Dr Eric Serre (Aix-Marseille University / M2P2) -
28
Advancing BIT1 Towards Exascale: Multi-GPU Hybrid Particle-in-Cell Monte Carlo Simulations at Scale
Particle-in-Cell (PIC) Monte Carlo (MC) simulations of the plasma edge play an important role in magnetic confinement fusion research, providing insights into plasma-surface interactions and plasma, impurity and neutral particle transport, relevant to both present and future fusion devices. As plasma simulations continue to increase in size and complexity, efficiently exploiting modern heterogeneous supercomputing systems presents significant challenges, including data movement, load imbalance, scalability, and resilience. This talk presents recent developments in advancing the Berkeley Innsbruck Tbilisi 1D3V (BIT1) PIC MC code [1,2] towards exascale computing through a portable hybrid MPI+OpenMP implementation capable of leveraging large-scale multi-GPU systems across both NVIDIA and AMD architectures. The work focuses on practical aspects of code development and optimisation, including GPU offloading strategies, overlapping communication and computation, memory management approaches, particle load balancing, and scalable checkpoint/restart capabilities. The integration of openPMD and ADIOS2 for standardised high-performance I/O is also discussed, enabling parallel I/O [3], enhanced diagnostics [4], in-memory data streaming [5], and in-situ analysis and visualisation workflows [4,5]. Performance and scalability results from leading pre-exascale and exascale systems, including Frontier, LUMI-G and MareNostrum 5 ACC are presented, together with lessons learned from profiling [6], porting [7,8], and enabling portability [9]. The talk concludes with ongoing efforts to improve resilience and workflow efficiency [10], together with future research extending hybrid BIT1 to Intel GPUs towards exascale, targeting Aurora and Europe's first exascale system, JUPITER Booster, to assess portability and performance across NVIDIA, AMD, and Intel GPUs and identify any remaining architecture-specific limitations at scale.
References
[1] D. Tskhakaya, et al., “Optimization of PIC codes by improved memory management,” Journal of Computational Physics, vol. 225, no. 1, pp. 829–839, 2007. doi:10.1016/j.jcp.2007.01.002
[2] D. Tskhakaya, et al., “PIC/MC code BIT1 for plasma simulations on hpc,” in 2010 18th Euromicro, pp. 476–481, IEEE, 2010. doi:10.1109/PDP.2010.47
[3] J.J. Williams, et al., "Enabling high-throughput parallel I/O in particle-in-cell Monte Carlo simulations with OpenPMD and darshan I/O monitoring." in 2024 IEEE International Conference on Cluster Computing Workshops (CLUSTER Workshops), IEEE, 2024. doi:10.1109/CLUSTERWorkshops61563.2024.00022
[4] J.J. Williams, et al., “Understanding the Impact of OpenPMD on BIT1, a Particle-in-Cell Monte Carlo Code, Through Instrumentation, Monitoring, and In-Situ Analysis,” in European Conference on Parallel Processing, pp. 214–226, Springer Nature Switzerland, 2024. doi:10.1007/978-3-031-90200-0_18
[5] J.J. Williams, et al., “Integrating High Performance In-Memory Data Streaming and In-Situ Visualization in Hybrid MPI+ OpenMP PIC MC Simulations Towards Exascale,” The International Journal of High Performance Computing Applications, p. 10943420251409229, 2025. doi:10.1177/10943420251409229
[6] J.J. Williams, et al., “Leveraging HPC Profiling and Tracing Tools to Understand the Performance of Particle-in-Cell Monte Carlo Simulations,” in European Conference on Parallel Processing, pp. 123–134, Springer Nature Switzerland, 2023. doi:10.1007/978-3-031-50684-0_10
[7] J.J. Williams, et al., “Optimizing BIT1, a Particle-in-Cell Monte Carlo Code, with OpenMP/OpenACC and GPU Acceleration,” in International Conference on Computational Science, pp. 316–330, Springer Nature Switzerland, 2024. doi:10.1007/978-3-031-63749-0_22
[8] J.J. Williams, et al., “Accelerating Particle-in-Cell Monte Carlo Simulations with MPI, OpenMP/OpenACC and Asynchronous Multi-GPU Programming,” Journal of Computational Science, vol. 88, p. 102590, 2025. doi:10.1016/j.jocs.2025.102590
[9] J.J. Williams, et al., “Multi-GPU Hybrid Particle-in-Cell Monte Carlo Simulations for Exascale Computing Systems,” in International Conference on Computational Science, pp. 32–47, Springer Nature Switzerland, 2026. doi:10.1007/978-3-032-29921-5_3
[10] J.J. Williams, et al., “High-Performance Resilient Multi-GPU Hybrid Particle-in-Cell Monte Carlo Simulations at Scale,” in Euro-Par 2026: Parallel Processing Workshops: Euro-Par 2026 International Workshops (BIGHPC), Pisa, Italy, August 24–28, Accepted For Publication, 2026 doi:10.48550/arXiv.2606.28534
Speaker: Mr Jeremy J. Williams (KTH Royal Institute of Technology) -
29
Drifts and Turbulence in Hermes-3
We report on progress in modeling of plasma drifts and 3D turbulence using Hermes-3, in tokamak and stellarator configurations. During the past year improvements have been made to models for the diamagnetic and polarization currents, ion viscosity and energy conservation. Tokamak applications include simulation of turbulence in NSTX and NSTX-U configurations, with comparison to experimental data.
Stellarator boundary modeling is challenging due to the complex magnetic field topology and wall geometry, in addition to the many interactions between plasma, neutrals and surfaces present in modeling of tokamak divertors. The Flux-Coordinate Independent (FCI) method has emerged as a promising approach to plasma discretization [1,2] that has sufficient flexibility to model real 3D devices. Recent advances in BOUT++ include relaxation methods to find steady-state plasma potentials[3], support-operator method (SOM) based FCI operators that conserve fluxes to machine precision on complex meshes, and nonlinear solvers and preconditioners to accelerate the solution of steady-state transport calculations.
Progress has been made to port BOUT++ and Hermes-3 to GPUs. Field operators have been rewritten to use template expressions to enable kernel fusion. The design and tradeoffs and future plans will be discussed.
Speaker: Dr Ben Dudson (LLNL) -
10:20
Break
-
30
Overview of solver, Kinetics for Transport
The dynamics of the edge of tokamak plasmas is understood to be sensitive to the electric field profile. Meanwhile modelling the long-wavelength component of the radial electric field is challenging for core physics, but essential to understand the rotation profile.
The moment-kinetics approach to evaluating transport was introduced under the Excalibur project (report TN-11). This aims to evolve a suitably normalised distribution function on the transport timescale consistently with momentum conservation, across the separatrix and in the presence of the strong changes in parameter values typical of detaching plasmas.
We present the work to date on the Kinetics for Transport (KfT) solver, which implements the moment-kinetics approach.
The solver is being developed in Julia, in a finite element framework with a Chebyshev basis, and supports various configurations for physics studies. Benchmarks and physics studies in 1D1V and 1D2V have been completed, while 2D1V simulations with adiabatic electrons demonstrated the presence of an instability relying on finite pitch angle of the magnetic field.
We have focused on developing a numerical approach that can simulate the sheath boundary conditions with kinetic ions and electrons with scalable performance, to enable 2D2V simulations. The implementation of the moment-kinetics structure and approach to electron preconditioning will be discussed.
Speaker: Sarah Newton (UKAEA) -
31
Some aspects of the finite element method applied to edge plasma modelling
The finite element method (FEM) concerns the numerical discretization of continuum partial differential equations. This talk briefly covers some applications of FEM to problems in plasma physics, including strongly orthotropic transport and plasma turbulence. We will emphasize the potential advantages of FEM in terms of numerical properties – stability, rapid convergence, and the suppression of numerical diffusivity (even without a fieldline-following computational geometry) – and geometrical fidelity – higher-order geometries, including faithful representations of curved surfaces. Finally, we will discuss how certain types of discrete representation, known as compatible - or mimetic- discretizations, preserve features of the continuum model.
Speaker: Ed Threlfall (UKAEA) -
32
Plasma Edge Numerics using Kinetic Neutrals Integrated into Finite Elements
UKAEA's edge code PENKNIFE is under development, incorporating the finite element framework Nektar++ and the reactions library VANTAGE. The code is designed to be used with modern hardware, exploiting the higher arithmetic intensity of finite element methods and the GPU parallelism of NESO-Particles for improved performance. A range of features have been implemented, encompassing a hierarchy of complexity from simple mean-field transport to 3D turbulence. Results from the latest simulations are presented, including comparisons with SOLPS.
Speaker: James Edgeley -
12:00
Lunch
-
13:00
💬 Discussion - Turbulence and mean-field codes
-
⚡ Lightning Talks / Poster Sessions: Set B
-
33
Transient SOLPS-ITER simulations detachment burn-through events on MAST-U
We present transient simulations of the MAST-U tokamak using the SOLPS-ITER model, operated in time-dependent mode and including time-dependent kinetic neutrals. Given the expected audience, we focus on the numerical requirements and challenges to run SOLPS-ITER in time-dependent mode, as well as the separate validation of the Eirene and B2.5 codes.
Speaker: David Moulton (UKAEA) -
34
Absolutely calibrated edge neutral inferences on MAST-U and first comparisons with edge models
Steven Thomas1,2,*, Jerry W. Hughes1, Alex Tookey2, Davis Easley3, Amal Reeja Biju4, Yacopo Damizia4, Jamie Dunsmore1, Bart Lomanowski3, Saskia Mordijck4, Ekin Öztürk4, Scott Silburn2, and the MAST Upgrade Team2,†
1MIT Plasma Science and Fusion Center, Cambridge, MA 02139, USA
2UKAEA, Culham Campus, Abingdon, Oxfordshire, OX14 3DB, UK
3Oak Ridge National Laboratory, Oak Ridge, TN 37831-6169, USA
4William & Mary, Williamsburg, VA 23185, USA
*email: sthoma@mit.edu; †See author list of [1]Recent work [2] has optimised the high-speed video (HSV) diagnostic on MAST Upgrade to obtain absolutely calibrated line-of-sight-integrated emission from neutral deuterium. The HSV is a wide-angle, high sampling rate (≥ 30 kHz) optical camera, fitted with an optical interference filter to isolate Dα line emission (n = 3 → 2). HSV data, interpreted using a collisional-radiative model, is used to infer radial profiles of ionization source rate, Sion, and neutral density, n0. Furthermore, Sion is used to obtain radial ion flux profiles, Γ.
Future fusion pilot plants will operate with burning plasmas to achieve high fusion power output, which requires high core pressures. As core transport is stiff, the maximum achievable pressure in the core relies on large edge pedestal gradients [3]. Therefore, predictions of plasma core performance require good predictions of the pedestal density and gradients, ne and ∇ne, respectively, which are set by a balance between Γ and Sion through the continuity equation:
∂t ne = −∇ · Γ + Sion. In the edge, sources from cold neutrals are plentiful and are the dominant contribution to the flux [2].The data from HSV provides an excellent resource for validating codes and testing models on MAST-U. For example, Γ is used in a diffusive-convective ansatz to infer diffusive and convective flux transport coefficients which are used as inputs to SOLPS-ITER and the Sion output directly compared to the experimental Sion inference. Preliminary attempts to model MAST-U discharges with KN1D [4] show disagreement between the simulated neutral penetration length and the inferred n0 profile from HSV. We use this opportunity to discuss the new edge neutral inferences in MAST Upgrade and how they can be used to constrain and test modelling of the plasma exhaust.
Acknowledgements: Work supported by DOE Awards DE-SC0023289, DE-SC0023372, and DE-AC05-00OR22725, and by the Engineering and Physical Sciences Research Council [grant number EP/W006839/1].
References[1] J. R. Harrison et al., Nucl. Fusion, vol. 64, p. 112017, 2024.
[2] S. Thomas et al., Plasma Phys. Control. Fusion, vol. 68, p. 105033 2026
[3] J. E. Kinsey et al., Nucl. Fusion, vol. 51, p. 083001, 2011.
[4] B. LaBombard, "KN1D: A 1-D Space, 2-D Velocity, Kinetic Transport Algorithm for Atomic and Molecular Hydrogen in an Ionizing Plasma," Tech. Rep. PSFC Research Report PSFC/RR-01-3, 2001.Speaker: Steven Thomas (MIT PSFC) -
35
Accelerating EDGE2D-EIRENE Studies with Deep Learning Surrogate Models
Surrogate enable the rapid evaluation of plasma quantities across different divertor operating conditions. Building on previous work on surrogate modelling for geometry, fuelling and impurity seeding studies, this work focuses on the continued development and evaluation of deep learning-based surrogates trained to predict 2D plasma variables from simulation inputs (both 1d and 2D).
We present an investigation into the impact of model design choices, including input feature selection, multi-output versus single-output training, model size, embedding dimension and dataset scaling. We will also demonstrate the surrogate's performance across conventional, Super-X and expanded divertor configurations, highlighting its ability to generalise across magnetic geometries and outlining downstream applications where rapid surrogate predictions can accelerate analysis and parameter studies.
Speaker: Naomi Carey (UKAEA) -
36
Physical preconditioning for neutrals in Hermes-3
Hermes-3 uses the pressure-diffusion model to represent the cross-field transport of fluid neutrals. This enables neutral advection, conduction and viscosity in the perpendicular direction according to diffusion coefficients which are each subject to flux limitation.
Hermes-3 has a fully implicit solver. Due to the inherent nonlinearity in the neutral equations and high neutral speeds through small cells, achieving good performance requires preconditioning.
While algebraic preconditioners such as LU are highly effective in 2D systems, they do not scale well to the large domains present in 3D turbulence or FCI applications of the code.
In this work, we consider an approach for physical preconditioning the neutral equations. Through an approximate block factorisation based on the physical blocks of neutral density, momentum and pressure, we derive a preconditioning approach that aims to preserve the underlying physics as much as possible while requiring an efficient approximation of the Schur complement that arises. Nonetheless, many questions and challenges still remain from the extreme scales and desire for overall computational efficiency in realistic simulations.Speaker: Niall Bootland (STFC) -
37
Assessing the Hermes-3 neutral model for stellarator scenarios
As part of the effort to establish Hermes-3 as a simulation tool for the stellarator edge and scrape-off layer, we examine the current implementation of the advanced fluid neutral model.
Due to the toroidal asymmetry of the strike line, stellarator scenarios could be more sensitive to artificial anisotropies in the neutral transport than tokamaks. We therefore evaluate the neutral anisotropy introduced by the pressure-diffusion model by means of simple test cases. This effort leads to an estimation under which conditions a full velocity fluid neutral model might be more opportune.
To further substantiate trust in the Hermes-3 neutral model, we attempt a qualitative comparison of steady-state results to the kinetic neutral model EIRENE, commonly used in conjunction with the fluid plasma model EMC3 as the state-of-the-art simulation suite for the stellarator edge.Speaker: Annika Stier (IPP-Greifswald) -
38
Understanding the 3D island divertor with a simplified model
The island divertor concept was proposed for power and particle exhaust in low-shear stellarators, such as Wendelstein 7-X. So far low downstream density has been experimentally measured. Stronger scaling of the downstream density is necessary to ensure the reactor-relevance of the island divertor in terms of particle exhaust. Numerical studies of the island divertor SOL require 3D codes, among which EMC3-Eirene is the state of the art. In this work, a two-point model is used to investigate the recycling regimes in simplified island divertor geometries to interpret the results of 3D simulations.
The two-point model used for tokamaks has been adapted into a stellarator two-point model (STPM). This contribution extends the STPM to include volumetric parameters: a convected power fraction, a target-localized dissipated power fraction, and a more general parametrization of the momentum loss factor. The extended STPM is validated against EMC3-Eirene simulations, using an island divertor geometry that is limited to a few flux surfaces of the powercarrying layer (PCL) in order to neglect flux-surface-perpendicular transport. The importance of the volumetric parameters is demonstrated and, when extracted from the 3D simulations, a reasonable agreement is found between both models. This work also shows that the STPM predicts a detrimental ”diffusionlimited” transport regime which suppresses high recycling. The PCL geometry retains typical 3D features of the island divertor such as the non-axisymmetry of the divertor targets and the presence of a target-shadowed region (TSR). The impact of these 3D features on parallel and perpendicular profiles is presented. Finally, the extended STPM provides the ability to do two-point-model formatting which allows to interpret 3D simulations results. Open and closed island divertor geometries are compared using such a formatting.
Speaker: Nassim Maaziz (IPP-Greifswald) -
39
Performance and architectural improvements in Hermes-3
Hermes-3 is a multifluid plasma simulation code, built on BOUT++. It is used for simulating transport and turbulence at the edge of magneticically-confined plasmas. It uses a system of modular components which make it highly flexible. I have provide research software engineering skills to support the ongoing development and optimisation of Hermes-3 and BOUT++. I will be presenting the work I have done on this over the past year, such as developing a new normalisation scheme for the BOUT++ SNES solver and improvements to the overall software architecture.
Speaker: Chris MacMackin (UK Atomic Energy Authority) -
40
Unstructured meshes for Hermes-3 using the FCI approach
In order to model the scrape-off layer of a fusion device,
traditionally field aligned codes have been used. For stellarators, this is not feasible, and the flux coordinate independent approach (FCI) is used.
Aligning the grid with the boundary is challenging, and generating
body fitted, structured grids restricts the homogeneity and
orthogonality of the grid, making numerical integration challenging.
Unstructured grids are well established and allow high quality, body fitted meshes and there should not be any road blocks to enable FCI codes to use unstructured meshes.For the FCI approach, the mesh is not field or flux surfaces aligned, thus removing the motivation why tokamak grids use the structured grid approach.
It is sufficient to interpolate the field for the parallel slices, but there is no fundamental reason why this should not be possible for unstructured meshes.
Currently a significant amount of time is spent on generating meshes. Closed divertors grids are especially challenging, restricting the grid quality significantly, as well as requiring a significant amount of work to generate body fitted, structured meshes. Generating unstructured meshes is trivially done, using state of the art meshing codes such as gmsh.Unstructured meshes allow to a much larger degree to refine the mesh, where steps in the plasma pressure are to be resolved. Especially if griding is easy and fast, as in the case for unstructured meshes, iterating with different meshes might be a viable approach for running steady state transport simulations.
It might also become more interesting to generate different meshes for neutrals and plasma quantities, as the neutral mesh could be much coarser in the core, which is not possible using structured meshes.This contribution will discuss the advantages, some alternatives and steps to a Hermes-3 version using unstructured meshes.
Speaker: David Bold (IPP Greifswald) -
41
Development and Validation of a Closed Island Divertor Concept for reactor-relevant exhaust in Wendelstein 7-X
The development of a viable plasma exhaust solution is one of the principal challenges on the path towards stellarator fusion reactors. While the island divertor of Wendelstein 7-X (W7-X) has demonstrated excellent performance in terms of detachment and steady-state operation, its present open geometry limits neutral compression and particle exhaust. A closed divertor configuration is therefore being investigated as a promising route towards enhanced density buildup, improved impurity retention, and increased radiative power dissipation.
In this work, we present the physics basis and numerical assessment of closed island divertor concepts for W7-X. The study is performed primarily with the EMC3-EIRENE code, a three-dimensional mean-field plasma-neutral transport solver, to investigate the influence of divertor closure on plasma and neutral dynamics. Guided by simplified heat transport modelling and realistic engineering constraints, a range of closed divertor geometries is explored. The simulations demonstrate promising trends, including enhanced downstream density buildup and favourable localization of plasma-neutral interactions, while maintaining acceptable heat loads on the divertor targets. The underlying physical mechanisms governing these improvements are identified and discussed.
Building on this understanding, we present simulations of a tungsten closed divertor configuration with impurity seeding to assess its exhaust performance under reactor-relevant conditions. These results provide important insight into the potential of closed island divertors and establish a physics-based foundation for future divertor upgrades and experimental implementation on W7-X.
Speaker: Amit Kharwandikar (Max-Planck-Institut für Plasmaphysik, Greifswald) -
14:40
Poster Session / Break
-
33
-
42
KN1DPy: A new Python translation of the KN1D neutral transport code
We present KN1DPy, a Python translation of KN1D, a spatially 1D kinetic neutral transport code that was originally written in IDL to model gas fuelling on the Alcator C-Mod tokamak. Given input profiles of electron density, $n_e$, and temperature, $T_e$, in the scrape-off-layer and on closed flux surfaces, KN1D solves the Boltzmann Equation to obtain the distribution functions for atomic and molecular neutral hydrogen. While its 1D geometry represents a step down in fidelity compared to 2D codes like EIRENE and DEGAS-2, this simplification also vastly reduces the computation time. Each KN1D run takes only a few seconds, and as a result KN1D can easily be included within integrated transport models to give an estimate for the edge ionisation source from gas fuelling. Furthermore, because KN1D returns the full distribution function for the neutrals, it also yields detailed information about the effective neutral temperature, the fraction of impinging neutrals reflected by the plasma, and the fraction of neutrals that have undergone charge-exchange.
We show some examples from DIII-D and Alcator C-Mod, and include comparisons with results from 2D codes on these devices. We also demonstrate how KN1DPy can be used to provide a quick estimate for the neutral source in integrated transport models.
Speaker: Jamie Dunsmore (MIT) -
43
GUERNICA: A continuum kinetic code for neutral dynamics in fusion plasmas
In this talk, GUERNICA, a new continuum kinetic neutral code, is presented for the simulation of atomic hydrogen transport in fusion edge plasmas, solving the Boltzmann equation with charge exchange, ionization and elastic neutral-neutral collisions. Neutral atoms strongly influence tokamak edge and divertor physics through fueling, momentum exchange and volumetric power loss [1]. Accurate neutral modeling is increasingly important for future devices, particularly for detachment which is central to managing power exhaust [2]. In contrast to Monte Carlo approaches, which suffer from noise and high computational cost in strongly collisional regimes, GUERNICA combines a discontinuous Galerkin discretization in configuration space with a discrete velocity formulation, enabling noise-free, high-fidelity solutions whose cost is independent of plasma-neutral collisionality. Verification and validation are demonstrated using analytic test problems and comparisons against DEGAS2 in one-dimensional configurations. Performance results show favorable scaling with velocity resolution and high throughput on single-GPU hardware. These results are extended to initial two-dimensional simulations in a MAST-U geometry driven by a static plasma background, underscoring the potential of deterministic kinetic neutral modeling for next-generation edge studies.
[1] S. Mordijck, Nucl. Fusion 60, 082006 (2020)
[2] U. Fantz et al, J. Nucl. Mater. 290–293, 367–73 (2001)Work supported by US DOE under DE-SC0007880.
Speaker: Jack Gabriel (William & Mary) -
44
Non-Local Transport of Vibrationally Excited D2 and its Impact on Divertor Detachment in MAST Upgrade Super-X Plasmas
Molecular processes are known to influence divertor detachment, but conventional effective-rate models assume vibrationally excited molecules remain in local equilibrium and neglect their transport [1, 2]. In this work, vibrationally resolved D2 simulations using the X1EXT molecular database are applied to both simplified divertor-leg and full-device MAST Upgrade Super-X SOLPS-ITER simulations to investigate the role of vibrational transport,electron cooling, and plasma-surface interactions during detachment [3, 4].
In isolated divertor-leg simulations, vibrationally excited molecules were found to survive for tens of microseconds and travel distances of several tens of centimetres before reacting. More than 55% of molecular trajectories exceeded the local electron-temperature gradient length, demonstrating strongly non-local behaviour. This transport enables excited molecules generated in warmer regions to penetrate detached plasmas, enhancing molecular charge exchange (MCX), molecular-activated recombination (MAR), and molecular-activated dissociation (MAD) far from their point of origin. Electron-impact excitation of D2 was also identified as an important power-loss channel, accounting for approximately 10% of total divertor power losses. These excitation losses (“electron quenching”) reduce electron temperature, promote low-temperature plasma-molecular interactions, and accelerate the onset of detachment.
The impact of these mechanisms was assessed in full-device simulations of MAST Upgrade Super-X plasmas and compared with spectroscopic measurements. Relative to effective-rate approaches, the vibrationally resolved model predicts detachment rollover at lower upstream density through enhanced MCX-driven MAR and stronger divertor cooling. The simulations reproduces experimentally observed Dα emission profiles, Fulcher-band emission,ionisation-front migration, and target ion-flux rollover more accurately than either AMJUEL or vibrationally unresolved X1EXT simulations.
Sensitivity studies show that plasma-surface interaction models have a strong impact on the vibrational state distributions of recycled molecules. The Eley-Rideal recycling model produces an overpopulation of high vibrational states and leads to unrealistically early detachment rollover, whereas the partial-thermalisation approach gives the closest agreement with experimental observations [5, 6]. These results highlight that predictive modelling of detached divertor plasmas requires a self-consistent treatment of vibrational state distributions, non-local molecular transport, and electron quenching, alongside plasma-surface interactions.
[1] K Verhaegh, et.al.
The role of plasma-atom and molecule interactions on power & particle balance during
detachment on the MAST Upgrade Super-X divertor. Nuclear Fusion, 63(126023), 2023.
[2] S. Kobussen, et.al.
Collisional radiative modelling with improved cross sections to investigate
plasma molecular interactions in divertor plasmas. Technical report, Masters Thesis, 2023.
[3] J Bryant, et.al.
Impact of yacora evaluated molecular effective rate coefficients
on detached solps-iter simulations. Nuclear Fusion, 2 2025.
[4] D Wunderlich, et.al.
Yacora on the Web: Online collisional
radiative models for plasmas containing H, H2 or He. Journal of Quantitative Spectroscopy
and Radiative Transfer, 240(106695), 1 2020.
[5] Anthony J.H.M. Meijer, et.al.
Isotope effects in the
formation of molecular hydrogen on a graphite surface via an eley-rideal mechanism. Journal
of Physical Chemistry A, 106, 2002.
[6] K Verhaegh, et.al.
Divertor shaping with neutral baffling as a solution to the tokamak power
exhaust challenge. Communications Physics, 8, 5 2025Speaker: Joseph Bryant (Tokamak Energy) -
45
CRFAX: Collisional-Radiative modelling Framework in JAX
The plasma behaviour in the exhaust region of tokamaks is strongly influenced by atomic and molecular processes, as experimentally demonstrated on, e.g., MAST-U [1], JET, AUG [2], Alcator C-Mod [3]. Faithful representation of these processes is therefore necessary for robust exhaust simulation capability. However, fully-resolved atomic and molecular reaction systems typically include hundreds and even thousands of excitation and rovibration levels across ionisation stages [4], rendering direct inclusion into exhaust fluid or kinetic codes computationally unfeasible. At the same time, an asymptotic treatment is not applicable, as neither local-thermodynamic nor coronal equilibriums hold in the typical divertor conditions. This, combined with high stiffness due to extreme timescale ranges, warrants for reduced (effective) models accounting for both collisional and radiative processes.
CRFAX (Collisional Radiative modelling Framework in JAX), being developed at UKAEA, is a Python package for modelling, analysis, reduction and post-processing of atomic and molecular systems in plasmas. A non-exhaustive list of present capabilities includes reduction methods for linearised reaction systems (a method by Greenland [5] and quasi-steady-state assumption), automatic reduced-model identification method [5], calculation of validity timescales and metrics of reduced models, solvers for linear and non-linear models, as well as calculation of rate contributions (emission spectrum and energy balance). Multi-dimensional scanning, automatic differentiation [6] and gradient-based fitting applied to parametrised reaction systems allow for applications such as tabulation of effective rate coefficients, sensitivity analysis and spectrum fitting. The provenance of CRFAX objects is tracked using W3C Prov data model [7]. The future development of CRFAX is aimed at reduction techniques for non-linear systems [8, 9] and implementation of uncertainty quantification.
[1] K. Verhaegh et al., Nucl. Fusion 63, 016014 (2023).
[2] M. Bernert et al., Nucl. Mater. Energy 12, 111-118 (2017).
[3] A. Mathews et al., Rev. Sci. Instrum. 93, 063504 (2022).
[4] K. Sawada, M. Goto, Atoms 4(4), 29 (2016).
[5] P. T. Greenland, Proc. R. Soc. A: Math. Phys. Eng. Sci., 457, 1821-1839 (2001).
[6] R. Frostig, M. J. Johnson, C. Leary, SysML Conf. 2018, hal-05188750, (2019).
[7] L. Moreau et al., PROV-DM: The PROV data model, http://www.w3.org/TR/prov-dm/ (2013).
[8] U. Maas, S.B. Pope, Proc. Combust. Inst. 24, 103-112 (1992).
[9] S. H. Lam, D. A. Goussis, Int. J. Chem. Kinet. 26(4), 461-486 (1994).Speaker: Andrei Ludvig-Osipov (UKAEA) -
17:00
Day 3 Close
-
08:55
-
-
08:55
Day 4 Open
-
46
Overview of the GRILLIX code
An overview of the GRILLIX [1] code is presented, a three-dimensional turbulence simulation code for the self-consistent modeling of the plasma edge and scrape-off layer. The code is based on a global drift-reduced fluid model coupled to a three-moment fluid neutral model. It employs the Flux-Coordinate Independent (FCI) approach, enabling efficient treatment of complex magnetic geometries, and has recently been extended to non-axisymmetric stellarator configurations [2]. Key scientific results obtained with GRILLIX are summarized, including simulations of advanced confinement regimes such as the X-point radiator [3] and Quasi-Continuous Exhaust regime [4], along with their experimental validation, studies of confinement transitions (L–H) [5], and first simulations of island divertor configurations in Wendelstein 7-X.
We present the modeling, numerical, and technical developments required to obtain these results. On the modeling side, in particular the self-consistent treatment of electromagnetic effects and extensions to low-collisionality regimes proved essential. While the Flux-Coordinate Independent (FCI) approach provides an efficient framework for accurately capturing strongly anisotropic dynamics, it also introduces several challenges, including the treatment of boundary conditions, the coupling to neutral models, and the achievement of good parallel scalability. Furthermore, the HPC aspects of GRILLIX are discussed. Next-generation supercomputing systems increasingly require efficient execution on GPU-accelerated hardware. We therefore present the porting of the elliptic field solver to GPUs [6], which represents both the most algorithmically complex and the most computationally demanding component of GRILLIX, while also featuring a well-defined interface that enables modular porting. Significant performance improvements are achieved with the GPU-accelerated implementation across different architectures.
The presentation concludes with an overview of ongoing work on open problems. In particular, the breakdown of the fluid assumptions may necessitate a kinetic description, for which the GENE-X code is being developed in close collaboration with GRILLIX. To address the challenges associated with boundary conditions, we present an approach that departs from the currently employed immersed boundary method and offers the potential for a more accurate and flexible implementation of boundary conditions within the FCI framework. Finally, lessons learned from the porting of the elliptic field solver have important implications for the GPU porting of the remaining GRILLIX components, which will be discussed.
[1] A. Stegmeir, D. Coster, A. Ross et al., PPCF 60:035005 (2018).
[2] A. Stegmeir, M. Finkbeiner, C.Pitzal et al., CPC 318:109874 (2026).
[3] K. Eder, W. Zholobenko, A. Stegmeir et al., NF 65:096029 (2025).
[4] K. Zhang, W. Zholobenko, A. Stegmeir et al., PRL submitted (2026).
[5] W. Zholobenko, F. Jenko, K. Zhang et al., PRL 136:075101 (2026).
[6] A. Stegmeir, C. Lalescu, M. Lin et al., PASC26 (accepted) (2026)Speaker: Andreas Stegmeir (Max Planck Institute for Plasma Physics) -
47
Porting Hermes-3 to GPUs: Current Status, Challenges, and Roadmap
Hermes-3, a multi-fluid plasma simulation framework built on BOUT++, is currently being ported to GPUs using the RAJA performance portability layer and CUDA. This work builds upon the GPU infrastructure previously developed at Lawrence Livermore National Laboratory (LLNL) for BOUT++, enabling the reuse of established programming models and optimization strategies.
The porting effort is being validated using the Annulus plasma turbulence model and three-dimensional edge plasma simulations of the TCV tokamak. We present the current status of GPU-enabled kernels and preliminary performance results, highlighting key challenges related to memory management, solver portability, and data movement. The roadmap toward full GPU capability, including further physics module acceleration and optimization for future exascale systems, is also discussed.
Speaker: Jony Castagna (UKRI-STFC Hartree Centre) -
48
Old Solver, New Tricks: modernizing and developing the UEDGE codebase in the age of AI
The 2D multi-fluid edge-plasma code UEDGE employs the fully time-implicit NKSOL Jacobian-Free Newton-Krylov (JFNK) solver (P.N. Brown, Y. Saad, SIAM J. Sci. Comput. 11, 1990). UEDGE solves a set of coupled plasma-neutral equations delivering 10–100$\times$ speedup compared explicitly-coupled fluid-plasma/kinetic-neutral code and numerical convergence to steady-state when including magnetic and $\mathbf{E}\times\mathbf{B}$ drift flows. Maintaining compatibility with the NKSOL solver while extending the physics model has shaped the development of UEDGE over three decades. The codebase, tracing its lineage to Braam’s B2 code, has grown organically into a tightly integrated Fortran architecture with nested physics loops with a stencil extending beyond local cells due to higher-order terms and the staggered grid, and order-dependent sections. This heritage of hard-won institutional expertise, encoded in its numerics, presents a challenge to code refactoring efforts: modernization is in perpetual competition with the interdependencies that provides robustness to the solver.
UEDGE’s time-tested but long-challenged fluid neutral model is one of its most notable features. The model has employed a 9-point stencil from its inception, as needed for non-orthogonal meshes: a choice since validated as critical for physics accuracy (Dekeyser et al., Nucl. Mater. Energy 18, 2019), dispelling much of the early criticism of fluid neutral models and paving the way for advanced fluid neutral (AFN) models. The focus of AFN model development has turned to capturing kinetic effects and the role of plasma chemistry involving minority species, such as hydrogenic molecules, and there are recent UEDGE development efforts addressing both topics. To capture molecule-assisted reaction pathways, a self-consistent model of diffusive fluid molecules, coupled to the plasma and atoms via effective collisional-radiative rates calculated with the bespoke CRUMPET code, was developed and implemented into UEDGE. To preserve compatibility with the JFNK solver, an expedient vacuum transport model using pre-computed DEGAS2 neutral trajectories through the halo plasma was developed and implemented, avoiding the significant refactoring associated with an unstructured, wall-conformal grid implementation. The inherent uncertainties carried by the effective dissociation rates, arising from the inevitable truncation of the reaction chains, and intricacies of the boundary conditions will be presented and discussed.
Several recent numerical and computational developments have focused on improving solver performance. OpenMP parallelization of the preconditioning Jacobian assembly and right-hand side evaluation achieves 80% reduction in wall-clock computational time on 32 CPU threads. The optimal decomposition is non-trivial and Random Forest regression models have been applied to identify the dominant performance factors. A continuation solver, leveraging the NKSOL JFNK solver to trace a continuous path through parameter space, producing a dense set of steady-state solutions, with ~22 s computational time per solution for a benchmark scan, was developed and implemented, enabling the creation of large databases needed for data-driven applications. Regression models trained on solver timings failed to identify meaningful improvements to the solver hyperparameters, and the underlying challenges will be reviewed.
The discussion will conclude with an overview of ongoing UEDGE developments. These include modernization efforts to decompose the monolithic Fortran source into modular, Python-level interfaces and potentially reimplement physics submodules in JAX; domain extensions to reach the vessel wall using semi-structured grids and API-like recursive calls; and a multi-layer agentic workflow for autonomous, accessible simulation. A demonstration of simulation capabilities will serve as an opportunity to introduce the Python-installable UEDGE and UETOOLS interface packages, both openly available on GitHub and PyPI.
This work was supported by US DOE under contract nos. DE-AC52-07NA27344 and was supported by the LLNL-LDRD Program under Projects No. 25-ERD-014 and23-ERD-015. LLNL-ABS-2020769
Speaker: Andreas Holm (Lawrence Livermore National Laboratory) -
10:20
Break
-
49
A One-Dimensional Integrated Workflow for Tungsten Sputtering, Transport, and Redeposition in the Divertor Region
Tungsten sputtering, redeposition, and transport in the divertor region are key processes governing impurity migration and contamination of the plasma core in fusion devices. Predictive modelling of these processes requires a coupled description of edge-plasma transport, sheath acceleration, plasma–surface interaction, sputtering, ionization, recombination, and impurity redeposition. In this work, we present a one-dimensional integrated workflow for modelling tungsten impurity production and transport in the divertor. The framework couples Hermes-3 for edge-plasma transport with a kinetic/analytical sheath model, a Monte Carlo impurity ion model, and a Direct Simulation Monte Carlo model for neutral tungsten atoms. Tungsten sputtering yields are calculated using a Monte Carlo ion tracer together with the Eckstein formula. Ionisation and recombination processes are included for tungsten charge states from neutral tungsten to $\text{W}^{20+}$. The model also supports multi-species incident ion populations contributing to tungsten sputtering.
The developed workflow provides self-consistent estimates of tungsten redeposition rates, sputtering yields, and charge-state-resolved tungsten impurity distributions. In addition, the one-dimensional model can supply sputtering and redeposition data for two-dimensional impurity transport simulations. Development of a fully self-consistent two-dimensional model is currently underway.
Speaker: Tengfei Tang (Agency of Science, Technology and Research) -
50
2D transport mechanisms near the X-points in snowflake divertors
Snowflake divertor (SFD) experiments show enhanced radial transport near the X-points, more than can be explained by standard 2D edge transport modelling (i.e. with poloidally uniform diffusion coefficients and without drifts). While 3D effects may be relevant here, several 2D transport mechanisms have been proposed in the literature which could explain this phenomenon: 1) a region of enhanced ExB drifts driven by large poloidal gradients in plasma profiles close to the primary X-point [1], and 2) the ‘churning mode’ [2], a toroidally-symmetric plasma vortex in the region of weak poloidal magnetic field near the X-points.
We have investigated these effects through two separate modelling approaches. Firstly, we have carried out interpretive modelling of recent SFD experiments at MAST-U using the code UEDGE with equilibrium drifts and currents included in the model. Secondly, we have developed a model of convective transport in the vicinity of the X-points which includes both the effect of equilibrium drifts and the churning mode. This latter model, implemented in BOUT++, allows us to directly compare the two mechanisms at conditions relevant to current and future SFD experiments.
From UEDGE modelling with drifts, we have found that ExB drifts can drive transport across the X-point region leading to activation of secondary strike points (SPs) and reduced peak heat fluxes at the primary SPs, consistent with the predictions in [1]. Direct heat flux measurements at the SPs were not available on the experiments modelled, but the UEDGE results show good agreement in midplane $T_e$, $n_e$ profiles and radiated power profiles in the divertor region.
From BOUT++ modelling, we find that the equilibrium drifts effect is dominant when the ratio of plasma pressure at the X-points to magnetic pressure at the midplane is low,$\beta_{pm} \lesssim 1\%$. Above this, electromagnetic effects become increasingly important in driving additional transport. Such conditions are attainable in future SFD experiments at MAST-U and NSTX-U, as well as during ELMs in current devices.
This work was carried out under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC5207NA27344.
[1] Canal et al., Nucl. Fusion 55 (2015)
[2] Ryutov et al., Phys. Scr. 89 (2014)Speaker: Dominic Power (LLNL / UKAEA) -
51
A large EDGE2D-EIRENE MAST-U database and learned surrogate models for diverse tasks
Accurate scrape-off layer (SOL) modelling is essential for tokamak heat and particle exhaust studies, supporting experiment interpretation, scenario optimisation, and reactor design. Existing methods require a compromise between speed and fidelity: analytical models are fast but simplified, while comprehensive simulations capture plasma and neutral physics but are too costly for rapid-turnaround applications.
This work presents a machine-learning framework for accelerating high-fidelity SOL modelling. A scalable automation workflow for EDGE2D-EIRENE was developed to generate tens of thousands of simulations spanning input power, fuelling and impurity seeding rates and locations, magnetic geometry, and transport coefficients. Neural networks, ensemble methods, tabular foundation models and neural operators were trained for tasks including two-dimensional plasma state reconstruction, derived quantity prediction (e.g., divertor peak heat flux), and boundary condition generation for core-edge integrated modelling.
The resulting surrogate models provide millisecond inference while retaining the fidelity of the underlying simulations, enabling high-fidelity SOL modelling for experimental support and broad design-space exploration.
Speaker: George Holt -
12:00
Lunch
-
13:00
💬 Discussion - Research Software Engineering and AI/agentic coding
-
14:00
End of Workshop
-
08:55