Session BO05: Fundamental Plasmas: Computational Techniques

Live

Computer simulations are one of the main ways in studying plasma systems. Employing simulation enable researcher to perform indirect measurement of systems from observation as well as to explore parameter space. However, accurate simulations are often slow to execute, which limits their applicability to extensive parameter exploration, large-scale data analysis, and uncertainty quantification. A promising route to accelerate simulations by building fast emulators with machine learning requires large training datasets, which can be prohibitively expensive to obtain with slow simulations. Here we present a method based on neural architecture search to build accurate emulators even with a limited number of training data. The method successfully accelerates simulations by up to 2 billion times in 10 scientific cases including 4 cases in plasma physics. Our approach also inherently provides emulator uncertainty estimation, adding further confidence in their use. We anticipate this work will accelerate research involving expensive simulations, allow more extensive parameters exploration, and enable new, previously unfeasible computational discovery.   

Live

The onset of major disruptions is an important issue for advanced tokamak plasmas such as the ITER experiment. While advanced statistical methods have been used to address the problem of tokamak disruption prediction and control, recent approach based on deep learning have proven particularly compelling. In this presentation, we will introduce new improvements to the fusion recurrent neural network (FRNN) software suite, which had recently delivered disruption predictions with unprecedented accuracy for measured signals from the EUROFUSION/JET and DIII-D tokamaks. Up to now, FRNN was based on the long short-term memory (LSTM) variant of recurrent neural networks. Here, we introduce and implement the “temporal convolutional neural network (TCN)” architecture for the representation of the time-dependent input signals, thus rendering the FRNN architecture fully convolutional. Our results demonstrate that FRNN models based on TCN achieve improved computational performance and prediction results when compared with the LSTM architecture for a representative fusion database. This represents a step forward in establishing a plasma device performance prediction platform with flexible architecture, capable of being tuned and adapted for different prediction needs for experimental datasets.   

Live

The nonlinear Fokker-Planck (FP) collision operator is an important component of high-fidelity plasma codes that contain non-Maxwellian plasmas (caused by neutral beam injection, alpha particles, or edge effect), yet when simulating multiple impurity charge states, numerically solving for the FP operator becomes prohibitively expensive. We present here a machine learning method, based on an encoder-decoder neural network, which accelerates the numerical solution of a multi-species, fully nonlinear FP operator. This neural network is trained on particle distribution function data generated from the XGC code, which is a massively parallel, gyrokinetic turbulence code focused on simulating the edge of tokamak plasmas. The network is trained such that its objective includes not only matching the output distribution function from the collision operator, but also a penalty to enforce physical conservation properties (mass, momentum, energy) of the FP operator. Relative conservation error on the order 10^-4 has been so-far achieved in the model training, and still improving. Implementation details including comparisons of the Picard iteration solver versus the neural network approach in terms of computational speed, memory usage, and accuracy in the XGC code will be shown.   

Live

The kinetic Landau-fluid (LF) closure with collision and periodic boundary condition is used in the development of deep learning (DL) surrogate model. Classical neural network, namely feedforward neural network (FNN), is constructed and trained to learn the kinetic LF closure with static limit and arbitrary mean-free-path in configuration space. The preliminary relation between best hyperparameters and critical parameters for data generation is found. Comparing with the numerical approach (non-Fourier method) of the LF closure, the deep learning surrogate model shows an order of magnitude of improvement in terms of accuracy. The surrogate model closure has been integrated for the first time with fluid simulations. Our DL-enabled fluid simulations give the correct Landau damping rate for a wide range of wavevectors, while Hammett-Perkins closure cannot produce correct damping rate. We correctly connect collisionless Hammett-Perkins closure and collisional Braginskii closure to reproduce the intrinsic nonlocal feature of heat flux with DL techniques. The simulations with the deep learning surrogate model are as good as simulations with the analytic closure in terms of long-term numerical stability in the linear Landau damping test.   

Live

Townsend discharge has been simulated using the speed-limited particle-in-cell (SLPIC) algorithm [1]. This algorithm limits the speed of the fastest electrons in the simulation to enable larger time steps and, therefore, faster computing times. The SLPIC algorithm facilitates a straightforward, fully-kinetic treatment of effects such as secondary emission and collisions. In this application, we simulated the amplification of a current accelerating across an argon-filled voltage gap, including electron-impact ionization, electron-neutral elastic collisions, and ion-induced secondary-electron emission. SLPIC provided more than an order of magnitude speed-up. SLPIC correctly determined the breakdown voltage (where the electron avalanche begins). [1]Werner, G. R., Jenkins, T. G., Chap, A. M., & Cary, J. R. (2018). Speeding up simulations by slowing down particles: Speed-limited particle-in-cell simulation. Physics of Plasmas, 25(12), 123512.   

Live

The guiding center approximation is commonly adopted in fusion particle codes to resolve the motion of charged particles in the core and at the plasma edge of a magnetic fusion device. However, the guiding center approximation does not hold at regions where the field gradients are large, such as the plasma sheath and presheath formed in front of a material surface. Here we present a hybrid particle method, which converts the guiding center particles into finite-orbit particles once they reach a region close to the boundary, where the field gradients become large. The domain is divided in a guiding center region and a full-orbit region, and the conversion of the particle is performed both ways. The conversion requires the adoption of a ring-charge method to avoid discontinuities in the potentials. We present numerical results from an implementation of this new scheme within the hPIC Particle-in-Cell developed at Illinois. We show that the approach can accurately resolve the full-orbit physics of the plasma sheath, still allowing to treat the vast majority of the plasma bulk with the usual guiding-center approximation.   

Live

The necessity for modeling collisional processes in plasmas is becoming ever more important as experimental efforts using higher density plasmas and solid targets come to fruition. We present progress towards an efficient module for simulating collisional ionization events within a particle-in-cell (PIC) framework. Our model has been tested rigorously for physical accuracy and does not suffer from statistical noise, thus decreasing the number of particles needed for a given simulation. This is done by calculating the rate of ionization deterministically and then adjusting the species densities within the simulation accordingly, which acts as a ``smoothing'' process reducing noise generated. Our model also includes proper momentum transfer due to the collisional process. This module has been integrated into the PIC code OSIRIS and has been benchmarked against other PIC codes, such as EPOCH and Smilei. We also use our model to simulate a variety of physical situations including electron beam propagation through air, electron stopping power through collisional ionization, and fast electron propagation through solids.   

Live

Symplectic integration methods are being investigated to determine if sufficiently smooth force interpolations can suppress the finite-grid instability and diffusive energy growth. The finite-grid instability arises in particle-in-cell (PIC) simulations when the Debye length is not resolved by the grid, while diffusive energy growth can occur even once the Debye length is resolved. Both numerical effects cause unphysical heating of the plasma over time. Sufficiently smooth symplectic methods preserve phase space structure, and explicit symplectic methods with smooth force interpolations have been found to have superior energy conservation properties than other explicit symplectic integrators [1,2]. Consequently, smooth symplectic algorithms may mitigate the growth of the grid instability, and enable explicit time integration that allows the Debye length to be under-resolved without incurring the computational cost of an implicit method. The effects of grid resolution on the unphysical energy growth rate in electrostatic PIC simulations will be presented for both standard and smooth symplectic time-integrators. [1] JR Cary and I Doxas, J. Comput. Phys., 107 (1), 98-104. [2] I Doxas and JR Cary, Physics of Plasmas, 4 (7), 2508-2518.   

Live

Experiments have shown that stellarators and heliotrons are subject to soft limits and loss of equilibrium due to changes in magnetic topology as plasma-beta is increased. We expect that non-ideal, time-dependent simulations, such as those in Refs. [1-2], will contribute to the understanding of this behavior. An efficient numerical representation is important for achieving resolution at realistic values of physical parameters, so we are developing a variant of the NIMROD code that features 3D geometry with the logical-to-physical coordinate mapping using a 2D spectral-element/1D Fourier representation. The logical periodic coordinate is a generalized toroidal angle. The only restrictions are that the mapping is continuous and periodic. Improved convergence is demonstrated through anisotropic heat-transport computations. We also compare the numerical ideal-MHD spectra and divergence-constraint properties of different formulations for the magnetic advance, using magnetic field or vector potential with either H1 or edge elements. [1] K. Ichiguchi, et al., PPCF 55 014009 (2013); [2] T. A. Bechtel, C. C. Hegna, and C. R. Sovinec, this meeting.   

Live

A key component in particle simulation models of plasma is solving for the trajectories of charged particles in electromagnetic fields. For strongly magnetized plasmas, this so-called particle pushing problem presents numerical difficulties due to the wide disparity in time scales between the fast timescale gyromotion and the long timescale macroscopic behavior, such as drift waves, which are the actual focus of interest in computational simulations. In other words, the problem is numerically stiff. Exponential integration is a lesser known time integration technique with excellent accuracy and numerical stability properties, and hence, offers a promising alternative to standard numerical particle pushers. For the specific case of E x B drift problems, we found that exponential integrators outperform the Boris and Buneman algorithms in terms of accuracy and computation time, even in the presence of an electric field gradient. We also discuss our current work on exploring the Hamiltonian formulation of the problem and in the development of symplectic and approximately symplectic exponential integrators for these problems.   

Live

We perform non-adiabatic simulations of warm dense aluminum based on the electron-force field (EFF) variant of wave-packet molecular dynamics (WPMD) [1]. Comparison of the ion-ion structure factor and dispersion relation with density functional theory (DFT) is used to validate the technique across a range of temperatures and densities spanning the warm dense matter regime. At 3.5eV and 5.2 g/cm$^{\mathrm{3}}$ we find a dispersion relation in close agreement with the more robust and adiabatic Kohn-Sham DFT. However, the approximations within EFF begin to breakdown at higher densities. To extend the region of applicability we suggest two improvements to the EFF model: (1) Improved approximation of the exchange and correlation energy and (2) a more physical basis set [2, 3]. \begin{enumerate} \item A. Jaramillo-Botero et al. J.Comput. Chem. 32, 497 (2011). \item R. A. Davis, W. A. Angermeier, R. Hermsmeier, and T. G. White. 2020. doi: ArXivID:2003.05566. \item H. Xiao et al., Mechanics of Materials 90, 243--252 (2015). \end{enumerate}   

Live

High-fidelity plasma simulations involving atomic kinetics requires state-space completeness often accompanied by many types of atomic transitions. The computational expense associated with this requirement drove the development of a collisional-radiative (CR) modeling package simulating reduced-order atomic kinetics through an automatic level-grouping process. Starting from the novel Boltzmann grouping methodology devised by Le et al.\footnote{Le et al. \textit{Phys Plasmas} 20, 1-19 (2013)}, spectral clustering techniques taken from machine learning automated the construction of appropriate level groups for time-dependent simulations. While the results from these simulations captured global plasma evolution well for various plasma conditions\footnote{Abrantes et al. \textit{J Comput Phys} 407, (2020)}, obtaining accurate radiative spectra remained problematic throughout the investigation. Further refinement of the clustering methodology to better account for spectral accuracy is therefore needed. In this work, a preliminary investigation adapting the clustering method to better capture these radiative effects will be introduced to facilitate the rapid development of automatically-reduced CR models that better approximate radiative effects for a wider range of plasma regimes.   

Arc Simulation with a Tightly Coupled Near-Wall Non-Equilibrium Sheath-Layer Model

Representing all the length and time scales of discharge in a numerical simulation can be prohibitively expensive for large-scale discharges such as in arc heaters. Simplifications, especially a local thermal equilibrium approximation, can reduce computational cost and is appropriate over much of the domain in such an application. However, equilibrium fails near electrodes. We develop and test a near-wall closure model to capture these non-equilibrium and integrate it with a bulk local equilibrium approximation to limit total computational cost. A boundary layer approximation is made such that the near-wall region is locally one dimensional, and multi-fluid non-equilibrium plasma simulations are solved using appropriate methods for coping with their numerical stiffness (due to electron transport) and are then tabulated for inclusion as boundary conditions in the full domain. The approach is demonstrated in an arc discharge with tungsten electrodes at atmospheric pressures. Agreement for total arc voltage depends upon inclusion of the non-equilibrium sheath representation and its coupling with the overall equilibrium simulation.