Session ZO04: Magnetic Confinement: Simulation and Whole Device Modeling

Live

Global gyrokinetic simulation of micro-turbulence without a source does not reach a steady state due to density/temperature profile relaxation. Numerical sources can be used to prevent profile relaxation. Several common forms of energy/particle source have been implemented in the gyrokinetitc $\delta \!f$-PIC code GEM. The rate of the source is often chosen somewhat arbitrarily, e.g a fraction of the maximum linear growth rate, or a fraction of the inverse eddy turnover time. We will compare simulations using different forms of source, and simulations using the same form of source but varying source rates, and answer the question whether there exists a regime where the steady state heat flux is insensitive to the exact form and value of source. This insensitivity is needed to DEFINE a steady state turbulence and anomalous transport. We will also investigate the closely related problem of the evolution of the marker density in phase-space. It is shown that, in partially-linearized PIC simulations, the marker density necessarily evolves away from the initially loaded marker distribution. Such evolution invalidates the typically used weight evolution equation. Techniques to mitigate this problem will be discussed.   

Live

To exchange particle distribution function between coupled codes in the ECP-WDM project and to evaluate the dissipative operations, such as plasma collisions, the total-$f$ gyrokinetic particle-in-cell code {\sc xgc}~[S.~Ku et al. {\em Phys. Plasmas} {\bf 25} 056107 (2018)] implements a velocity grid and a bilinear mapping between marker particles (continuous in velocity space) and velocity grid [Yoon, Chang {\em Phys. Plasmas} {\bf 21} 032503 (2014)]. The bilinear operation ensures conservation of particle density and momentum, but fails to conserve particle energy with enough accuracy. In the present work we have updated {\sc xgc}~to instead use a novel mapping technique, recently included in the {{\sc pets}c}~library [S.~Balay et al. {{{\sc pets}c}~Users Manual} {ANL-95/11~-~Revision~3.13} (2020)], which employs a pseudo-inverse to preserve moments up to the order of the discretization space [Hirvijoki et al. arXiv:1802.05263]. We demonstrate the functionality and that $2^{\mathrm{nd}}$-order elements, in addition to particle and momentum conservation, also conserve energy.   

Live

This project develops fast, parallel numerical methods to resolve Kelvin-Helmholtz and drift-wave instabilities with high-order (HO) accuracy in space and time. A drift-reduced extended magnetohydrodynamic (XMHD) model is used to describe the effects of macroscopic transport phenomena in plasma at the edge of a confinement device, which we then solve using high-order (HO) approximations, including HO curvilinear meshes discretized with HO finite elements and HO time integration schemes. HO methods offer unmatched resolution of the stiff nonlinear behavior of edge plasma and drift instabilities, but introduce numerous difficulties in solving the resulting equations. Here, we present a broad framework for the parallel numerical solution of HO XMHD models. The spatial problem is discretized using HO finite elements in the MFEM library, yielding a semidiscrete set of differential algebraic equations in time. For each time step, we apply a new framework for the fast parallel solution of fully implicit Runge-Kutta methods, coupled with an Anderson-accelerated nonlinear iteration. Each inner linear iteration is then solved implicitly using block preconditioning techniques, and a new nonsymmetric algebraic multigrid method called AIR is applied to the highly advective variables.   

Live

The high-order MFEM finite element framework, which provides state-of-the-art numerical techniques for discretization and meshing, is used to efficiently solve the drift-reduced extended magnetohydrodynamics physics equations. By eliminating the fast compressional Alfven and sound waves, these models efficiently describe turbulence driven by drift wave instabilities and magnetic reconnection mediated by the shear Alfven wave. The discretization is based on finite element exterior calculus, which mimics the fundamental theorems of vector calculus. By utilizing the extensive set of high-order finite element spaces, such as H1, Nedelec, Raviart-Thomas, and Discontinuous Galerkin, one can achieve desirable numerical properties such as the ability to prescribe divergence-free magnetic fields, incompressible flows, and conservation of energy to machine accuracy. In addition to providing greater numerical accuracy, the high-order elements also provide curved meshes that can be aligned with magnetic flux surfaces. We illustrate these capabilities through applications to edge plasma physics examples and explore the use of advanced meshing techniques for improving accuracy near X-points caused by divertors and magnetic islands. This work was performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344 and was supported by LLNL Laboratory Directed Research and Development project 20-ERD-038.   

Live

Modeling and control of plasmas is a notoriously challenging, yet vital topic in modern physics. This work focuses on the development of several novel reduced-order modeling frameworks for compressible plasmas, leveraging decades of progress in first-principles and data-driven modeling of fluids. These theoretical frameworks enable the development of sparse and interpretable nonlinear reduced-order models from data that are intrinsically connected to the underlying physics. We demonstrate the effectiveness of these approaches on data from high-fidelity numerical simulations. These techniques prove promising for the prediction, estimation, and control in industrial and laboratory plasmas.   

Live

A new, fully data-driven transport model has been developed that uses neural networks to predict plasma profiles on a scale of the energy confinement time into the future given actuators and the present plasma state. The model was trained and tested on DIII-D data from the 2010-2018 experimental campaigns. The model is accurate on average and is shown to scale properly with actuators, with rotational transform predictions the worst and pressure predictions the best. The model can run in milliseconds and is very simple to use. This makes it a potentially useful tool for operators and physicists when planning plasma scenarios. It also is a candidate for doing phase-space exploration without going through the DIII-D database or complicated and expensive simulation codes. A reduced model using only realtime diagnostics has also been developed and formed the basis for a model-predictive control algorithm implemented and successfully tested on DIII-D.   

Live

Accurately calculating the heat flux in fusion plasmas is computationally prohibitive using first-principles gyrokinetic codes. The trapped gyro-Landau-fluid (TGLF) code addresses this problem by solving a reduced set of gyrokinetic equations. To accurately model the nonlinear saturation of turbulence, the TGLF code employs so-called saturation rules SAT0 or SAT1. To validate the TGLF model, we built a database containing $2500$ plasma discharges in the DIII-D tokamak, for which we have generated a corresponding database of $1.8\times10^5$ time and space slices. The data was filtered to eliminate unphysical cases with negative energy fluxes and MHD unstable cases. Moreover, we have eliminated cases close to the thresholds of kinetic ballooning modes and drift wave turbulence. The two saturation models SAT0 and SAT1 were subsequently validated with the filtered dataset of $10^5$ cases. Lastly, to help in the validation efforts we applied machine learning tools to the filtered dataset. As a consistency check for our neural network, we find that we are able to accurately reproduce the free parameters of the saturation rules SAT0 and SAT1, which have previously been calibrated by GYRO. These tools will help identify any promising areas for improvement of these saturation rules.   

Live

The physics of fusion plasmas includes multiscale phenomena that can be prohibitively expensive to simulate. Projective integration offers a computational framework for bridging divided timescales. Toward this end, we demonstrate the data-driven construction of linear time advance operators, which are devised to be separable into the fast and slow components of simulated dynamics. Our work is based on a least-squares approximation obtained via singular value and dynamic mode decomposition. We filter the components of these procedures to specifically enable large timesteps, resulting in a spatially reduced rank and temporally partitioned operator for use in standard integration schemes. We verify the capability of these operators to reconstruct dynamics extracted from simulated data. In addition, we test their stability characteristics for first order methods. Our results are detailed for linear cases of 1d diffusion and advection. We then investigate modeling scenarios relevant to fusion plasma physics by applying this approach to 2D SOLPS steady-state simulations. From these benchmarks, we conclude by extending the developed data-driven algorithm to a projective integration framework for multiscale simulation.   

Live

Anew integrated modeling approach has been developed allowing the prediction of the kinetic profiles of tokamak plasmas from magnetic axis to separatrix only using global parameters as inputs. In particular, anew pedestal transport model, based on empirical observations from multiple devices, is included in the ASTRA transport code and applied in combination with the TGLF and NCLASS modules for core turbulent and neoclassical transport. Asimple but realisticscrape-off layer model computesthe separatrix boundary conditions asfunction of the main engineering parameters. In this way, no information fromkinetic profile measurementsis required as input of the integrated modeling workflow, and the only inputs of the model are the magnetic field, the plasma current, the heating power, the fueling rate, and the plasma geometry. The pedestal top pressure is determined using the MISHKA MHD stability code. This model is applied to 50 stationary ASDEX Upgrade H-mode plasmas. Changes inpedestal structure and coregradients, produced by variations in many operational parameters, are well captured by the model. The predicted stored energies are in better agreement with the experimental observations than those obtained bythe IPB98(y,2) scaling law.   

Live

A new workflow in the OMFIT integrated-modelling framework has been developed to predict profiles and energy confinement, based on zero-dimensional (0D) tokamak quantities. This workflow addresses one of the present limitations of systems studies, which currently relies on the experimental energy confinement scaling law $\tau_{98,y2}$. We seek to obtain a fully theoretical prediction by progressively dropping assumptions and replacing simple scaling laws with state-of-the-art theory-based physics models. In OMFIT the PRO-create (profiles creator) module generates physically plausible plasma profiles and a consistent equilibrium using the same $\tau_{98,y2}$ 0D parameters. This result forms the starting point for the STEP (Stability, Transport, Equilibrium, and Pedestal) module which iterates between equilibrium, sources, core transport, and pedestal calculations to obtain a self-consistent solution. We will report on the validation of this workflow, as has been carried out on the ITER H-98P(y,2) database and a series of DIII-D plasmas, yielding a theory-based energy confinement scaling, and its applications towards the evaluation of potential DIII-D upgrades and the design optimization of next generation fusion devices.   

Live

Extensive experiments were conducted in TFTR with deuterium and tritium plasmas 1993 - 1997. We report results from recent reanalysis updating previous results using the TRANSP integrated analysis code. The purpose is to better understand neutral beam injection (NBI) effects in DT plasmas. These results are relevant to upcoming JET DT experiments, and to improve understanding of Toroidal Alfven mode effects seen in TFTR. Improvements to the TRANSP code since the early analysis include more comprehensive and accurate atomic physics data, improved Monte Carlo simulation of fast ion parameters, and improved equilibrium solutions. The D and T species mix in the NBI correlates with the core ion and electron temperatures, the thermal hydrogenic isotopic mix, and the energy confinement times. The shape of the carbon density was well diagnosed, but shapes of other impurities were not measured. We discuss evidence that in the core early in the NBI phase Z$_{eff}$ is lower than that implied by only the carbon density (consistent with dilution from Li injection) and higher late in the NBI (consistent with higher Z impurities accumulating). Comparisons of NBI and alpha heating are discussed.   

Live

To resolve large gradients in the plasma sheath region, the mesh of the hPIC Particle-In-Cell code has been modified from a classical uniform mesh to a non-uniform block-structured implicit mesh using the Parallel Unstructured Mesh Infrastructure (PUMI) library. The implicit nature allows to define a mesh with a minimal number of parameters and to generate all mesh quantities on-the-fly. The algorithm allows to split the entire domain into a number of submeshes either of uniform or boundary layer type. A boundary layer submesh employs a geometric gradation in element size starting from one side. The performance measurement of the multi-block PUMI mesh incorporated in hPIC has been done for different domain sizes and mesh configurations. For a small domain size of 500 Debye lengths, a speed-up of up to 16 times with respect to a uniform grid is achieved, maintaining the global error on the solution at about 1{\%}. For a large plasma domain (1.5 m), a speed-up of more than 100 times is achieved with respect to the uniform grid. This demonstrates the ability, thanks to the non-uniform mesh, of simulating a large plasma domain (meters long) at fusion relevant conditions with a finite-orbit PIC in a reasonable computational time.   

Live

The Exascale High-Fidelity Whole-Device-Modeling project aims at delivering an application composed of many physics components coupled together at the first-principles level. We study the spatial coupling of two gyrokinetic codes, one for the core (such as GENE or GEM) and one for the edge (XGC). The new generalized coupling scheme, which is introduced, combines the coupling of Poisson equation [Dominski et al, Phys. Plasmas 072308 (2018)] with the new kinetic coupling of particle distribution functions. The coupling of particle distribution functions is performed only once in multiple time-steps, by interfacing the core and edge simulations with a flexible 5D grid. This generalized interfacing enables the coupling of different models, such as delta-f and total-f, as well as the coupling of particle-in-cell and continuum codes. The transfer of the particle distribution function between marker particles and the 5D grid is based on a new resampling technique. A first demonstration of first-principles spatial coupling in a DIII-D like plasma will be shown, by using the XGC suite of codes for core and edge sides   

Live

Recent fusion planning activities recommend that the U.S. should pursue innovative science and technology to enable construction of a pilot plant that produces net electricity from fusion at reduced capital cost. Such a mission requires discovery, development, prototyping, and integration of multiple physics and technology innovations. In this work, tokamak configurations are explored to determine the potential synergies and conflicts between proposed innovations. For example, lower aspect ratio (A=1.8-2.6) tokamaks are potentially advantageous for maximizing self-driven current fraction and plasma performance per unit magnet cost but would have reduced space for a central solenoid for long-pulse partial-inductive operation. Performance characteristics are studied as a function of device size, magnet capability, aspect ratio, and stability and confinement assumptions. Realistic shaping poloidal field coils using free-boundary equilibrium calculations and superconducting magnet layout and stress analysis are also investigated while incorporating long-leg and liquid metal divertors. Initial device configurations, physics scenarios, and engineering studies for a tokamak-based next-step facility dedicated to the exploration and integration of these innovations are described.   

Perspectives for the High Field, Compact Machine Approach and Advances Within the Ignitor Program

The confirmed confinement and purity properties of the high density plasmas produced by the high field, compact line of experiments and the fact that the reactivity of D-T plasmas increases as $B_{P}^{4} $, continue to indicate that this line is the most promising in the effort to approach ignition conditions ($B_{P} $ is the poloidal field). The Ignitor Program, started from the Alcator and the Frascati Torus programs, has been the first to develop a complete machine design with the objective of approaching ignition, and has incorporated, systematically, advances made in relevant structural engineering, materials (including high field superconducting magnets), etc. in order to maintain it at the forefront of fusion research. Ongoing activities include involving outstanding departments of the Sapienza University to advance the existing design and, in the near term, the fabrication of the machine central post, of the central solenoid and of the large vertical field coils, the collaboration with Rosatom on the cryogenic system, etc.