Plasma Physics Simulations Research Articles

Verification methods for drift–diffusion reaction models for plasma simulations

Compared to other computational physics areas such as codes for general computational fluid dynamics, the documentation of verification methods for plasma fluid codes remains under developed. Current analytical solutions for plasma are often highly limited in terms of testing highly coupled physics, due to the harsh assumptions needed to derive even simple plasma equations. This work highlights these limitations, suggesting the method of manufactured solutions (MMSs) as a potential option for future verification efforts. To demonstrate the flexibility of MMS in verifying these highly coupled systems, the Multiphysics Object-Oriented Simulation Environment (MOOSE) framework was utilized. Thanks to the MOOSE framework’s robustness and modularity, as well as to its physics module capabilities and ecosystem applications (i.e. Zapdos and the chemical reaction network) developed for plasma physics modeling and simulation, this report lays the groundwork for a structured method of conducting plasma fluid code verification.

Compression of tokamak boundary plasma simulation data using a maximum volume algorithm for matrix skeleton decomposition

This report demonstrates satisfactory data compression of SOLPS-ITER simulation output ranging from 2D fields, 1D profiles, and 0D scalar variables with a novel matrix decomposition approach. The singular value decomposition (SVD) scales poorly for large matrix sizes and is unsuited to the application on high dimensional data common to fusion plasma physics simulation. We employ the columns-submatrix-rows (CUR) matrix factorization technique in order to compute a low-rank approximation up to two orders of magnitude faster than the SVD, but within a nominal L2-norm relative error of ϵ=10−2. In addition, the CUR approach maintains the original format of the data, in its extracted columns and rows, allowing for interpretable data storage at the original resolution of the simulation. We utilize an iterative algorithm to compute the CUR decomposition of simulation output by maximizing the volume, or linearly independent information content, of a low-rank submatrix contained within the data. Experiments over n×n randomized test matrices with embedded rank-deficient features show that this maximum volume implementation of CUR matrix approximation has reduced asymptotic computational complexity on the order of n compared to the SVD, which scales approximately as n3. These results show that the CUR technique can be used to effectively select time step snapshots (columns) of over 140 SOLPS-ITER output variables and the associated discretized coordinate timeseries (rows) allowing for reconstruction of the complete simulation dynamics.

Stochastic Galerkin particle methods for kinetic equations of plasmas with uncertainties

The study of uncertainty propagation is of fundamental importance in plasma physics simulations. To this end, in the present work we propose a novel stochastic Galerkin (sG) particle method for collisional kinetic models of plasmas under the effect of uncertainties. This class of methods is based on a generalized polynomial chaos (gPC) expansion of the particles' position and velocity. In details, we introduce a stochastic particle approximation for the Vlasov-Poisson system with a BGK term describing plasma collisions. A careful reformulation of such dynamics is needed to perform the sG projection and to obtain the corresponding system for the gPC coefficients. We show that the sG particle method preserves the main physical properties of the problem, such as conservations and positivity of the solution, while achieving spectral accuracy for smooth solutions in the random space. Furthermore, in the fluid limit the sG particle solver is designed to possess the asymptotic-preserving property necessary to obtain a sG particle scheme for the limiting Euler-Poisson system, thus avoiding the loss of hyperbolicity typical of conventional sG methods based on finite differences or finite volumes. We tested the schemes considering the classical Landau damping problem in the presence of both small and large initial uncertain perturbations, the two stream instability and the Sod shock tube problems under uncertainties. The results show that the proposed method is able to capture the correct behavior of the system in all test cases, even when the relaxation time scale is very small.

Data-driven linear time advance operators for the acceleration of plasma physics simulation

We demonstrate the application of data-driven linear operator construction for time advance with a goal of accelerating plasma physics simulation. We apply dynamic mode decomposition (DMD) to data produced by the nonlinear SOLPS-ITER (Scrape-off Layer Plasma Simulator - International Thermonuclear Experimental Reactor) plasma boundary code suite in order to estimate a series of linear operators and monitor their predictive accuracy via online error analysis. We find that this approach defines when these dynamics can be represented by a sequence of approximate linear operators and is essential for providing consistent projections when compared to an unconstrained application. For linear diffusion and advection–diffusion fluid test problems, we construct and apply operators within explicit and implicit time advance schemes, demonstrating that stability can be robustly guaranteed in each case. We further investigate the use of the linear time advance operators within several integration methods including forward Euler, backward Euler, and the matrix exponential. The application of this method to simulation data from SOLPS-ITER, with varying levels of Markov chain Monte Carlo numerical noise, shows that constrained DMD operators yield a capability to identify, extract, and integrate a (slow) subset of the present timescales. Example applications show that for projected speedup factors of 2×, 4×, and 8×, a mean relative error of 3%, 5%, and 8% and maximum relative error less than 20% are achievable, which appears acceptable for typical SOLPS-ITER steady-state simulations.

Combined molecular and atomic potentials for elastic cross sections of electrons scattering off diatomic molecules at intermediate energies

A simple model is proposed to compute electron-diatomic molecule elastic differential cross sections at intermediate energies within the framework of an analytical local optical potential. In a spherical harmonic expansion of the molecular potential we treat the isotropic term with the partial-wave decomposition and apply the independent atom model (IAM) on all higher orders. This model is seen to properly converge to the IAM at high energies, while bringing significant improvement at lower energies. We compare the results with a well-tested program called ELSEPA, tailored for high-energy electron-atom scattering and address its further extension to molecules with the method proposed. The simplicity of the calculations and the encouraging agreement in shape with experimental data could promote attractiveness among plasma physics simulations in need of coherent and well-resolved differential cross sections.

Modelling of plasma gamma ray sources in large tokamaks

• Methodology for creation of realistic plasma gamma ray source for Monte Carlo simulations. • Connection between plasma physics simulations and gamma ray Monte Carlo transport. • Methodologies for generation of position, energy and direction sampling distributions from plasma parameters. • Methodology testing with a simplified MCNP tokamak model. Understanding the physics of fast ions in a fusion plasma is widely considered as one of the crucial tasks for the reliable operation of fusion tokamak reactors. Measurements on tokamaks have shown that gamma rays are produced when fast ions react with either plasma fuel ions or with the plasma impurities such as beryllium, carbon and oxygen. The spectroscopy of these gamma rays can be used for measurement of fusion rates in plasma or behaviour and confinement of fusion reaction products, such as the alpha particles. To computationally support experiments a methodology has been developed for creation of a realistic plasma gamma ray source for Monte Carlo simulations. The methodology presented in this paper consists of several steps. First the generation of realistic plasma parameters with the code TRANSP on the basis of experimental measurements. Second the generation of sampling distribution functions from the calculated plasma parameters such as the calculation of the source position distribution function based on the rate densities for the gamma ray emitting reaction. Third the creation of realistic plasma gamma ray source for the code MCNP gamma ray transport simulation. The developed methodology was analysed on a MCNP tokamak model to demonstrate that the methodology generates a realistic plasma gamma ray source for Monte Carlo simulations.

A new method to dispatch split particles in Particle-In-Cell codes

Particle-In-Cell codes are widely used for plasma physics simulations. It is often the case that particles within a computational cell need to be split to improve the statistics or, in the case of non-uniform meshes, to avoid the development of fictitious self-forces. Existing particle splitting methods are largely empirical and their accuracy in preserving the distribution function has not been evaluated in a quantitative way. Here we present a new method specifically designed for codes using adaptive mesh refinement. Although we point out that an exact, distribution function preserving method does exist, it requires a large number of split particles and its practical use is limited. We derive instead a method that minimizes the cost function representing the distance between the assignment function of the original particle and that of the sum of split particles. Depending on the interpolation degree and on the dimension of the problem, we provide tabulated results for the weight and position of the split particles. This strategy represents no overhead in computing time and for a large enough number of split-particles it asymptotically tends to the exact solution.

Evolution and observational signatures of the cosmic ray electron spectrum in SN 1006

ABSTRACT Supernova remnants (SNRs) are believed to be the source of Galactic cosmic rays (CRs). SNR shocks accelerate CR protons and electrons which reveal key insights into the non-thermal physics by means of their synchrotron and γ-ray emission. The remnant SN 1006 is an ideal particle acceleration laboratory because it is observed across all electromagnetic wavelengths from radio to γ-rays. We perform 3D magnetohydrodynamics (MHD) simulations where we include CR protons and follow the CR electron spectrum. By matching the observed morphology and non-thermal spectrum of SN 1006 in radio, X-rays, and γ-rays, we gain new insight into CR electron acceleration and magnetic field amplification. (1) We show that a mixed leptonic–hadronic model is responsible for the γ-ray radiation: while leptonic inverse-Compton emission and hadronic pion-decay emission contribute equally at GeV energies observed by Fermi, TeV energies observed by imaging air Cherenkov telescopes are hadronically dominated. (2) We show that quasi-parallel acceleration (i.e. when the shock propagates at a narrow angle to the upstream magnetic field) is preferred for CR electrons and that the electron acceleration efficiency of radio-emitting GeV electrons at quasi-perpendicular shocks is suppressed at least by a factor ten. This precludes extrapolation of current 1D plasma particle-in-cell simulations of shock acceleration to realistic SNR conditions. (3) To match the radial emission profiles and the γ-ray spectrum, we require a volume-filling, turbulently amplified magnetic field and that the Bell-amplified magnetic field is damped in the immediate post-shock region. Our work connects microscale plasma physics simulations to the scale of SNRs.

Nine Outstanding Questions of Solar Wind Physics.

In situ measurements of the solar wind have been available for almost 60 years, and in that time plasma physics simulation capabilities have commenced and ground‐based solar observations have expanded into space‐based solar observations. These observations and simulations have yielded an increasingly improved knowledge of fundamental physics and have delivered a remarkable understanding of the solar wind and its complexity. Yet there are longstanding major unsolved questions. Synthesizing inputs from the solar wind research community, nine outstanding questions of solar wind physics are developed and discussed in this commentary. These involve questions about the formation of the solar wind, about the inherent properties of the solar wind (and what the properties say about its formation), and about the evolution of the solar wind. The questions focus on (1) origin locations on the Sun, (2) plasma release, (3) acceleration, (4) heavy‐ion abundances and charge states, (5) magnetic structure, (6) Alfven waves, (7) turbulence, (8) distribution‐function evolution, and (9) energetic‐particle transport. On these nine questions we offer suggestions for future progress, forward looking on what is likely to be accomplished in near future with data from Parker Solar Probe, from Solar Orbiter, from the Daniel K. Inouye Solar Telescope (DKIST), and from Polarimeter to Unify the Corona and Heliosphere (PUNCH). Calls are made for improved measurements, for higher‐resolution simulations, and for advances in plasma physics theory.

Versatile fusion source integrator AFSI for fast ion and neutron studies in fusion devices

ASCOT Fusion Source Integrator AFSI, an efficient tool for calculating fusion reaction rates and characterizing the fusion products, based on arbitrary reactant distributions, has been developed and is reported in this paper. Calculation of reactor-relevant D–D, D–T and D–3He fusion reactions has been implemented based on the Bosch–Hale fusion cross sections. The reactions can be calculated between arbitrary particle populations, including Maxwellian thermal particles and minority energetic particles. Reaction rate profiles, energy spectra and full 4D phase space distributions can be calculated for the non-isotropic reaction products. The code is especially suitable for integrated modelling in self-consistent plasma physics simulations as well as in the Serpent neutronics calculation chain. Validation of the model has been performed for neutron measurements at the JET tokamak and the code has been applied to predictive simulations in ITER.

Computer Simulation in Low‐Temperature Plasma Physics: Future Challenges

Computer simulations can be carried out with various aims. Perhaps the most challenging is prediction under conditions where experiments are difficult or inaccessible, especially when failure to predict adequately may have unhappy consequences. There is, probably, not much confidence at present in the capability of low‐temperature plasma physics simulations in such a context. Other fields have attempted to meet this challenge using a collection of techniques collectively known as “Verification and Validation,” or “V&V.” These are methods for enhancing confidence in the correctness and fidelity of computer simulations. This paper surveys these techniques and discusses their application to improvements of simulation capability in low‐temperature plasma physics.

InfiniCloud 2.0: distributing High Performance Computing across continents

InfiniCloud 2.0 is World’s first native InfiniBand High Performance Cloud distributed across four continents, spanning Asia, Australia, Europe and North America. The project provides researchers with instant access to computational, storage and network resources distributed around the globe. These resources are then used to build a geographically distributed, virtual supercomputer, complete with globally-accessible parallel file system and job scheduling.This paper describes high level design and the implementation details of InfiniCloud 2.0. A gene sequencing pipeline as well as plasma physics simulation code are used to demonstrate system’s capabilities.

A high-order Boris integrator

This work introduces the high-order Boris-SDC method for integrating the equations of motion for electrically charged particles in electric and magnetic fields. Boris-SDC relies on a combination of the Boris-integrator with spectral deferred corrections (SDC). SDC can be considered as preconditioned Picard iteration to compute the stages of a collocation method. In this interpretation, inverting the preconditioner corresponds to a sweep with a low-order method. In Boris-SDC, the Boris method, a second-order Lorentz force integrator based on velocity-Verlet, is used as a sweeper/preconditioner. The presented method provides a generic way to extend the classical Boris integrator, which is widely used in essentially all particle-based plasma physics simulations involving magnetic fields, to a high-order method. Stability, convergence order and conservation properties of the method are demonstrated for different simulation setups. Boris-SDC reproduces the expected high order of convergence for a single particle and for the center-of-mass of a particle cloud in a Penning trap and shows good long-term energy stability.

Spectral Solver for Multi-scale Plasma Physics Simulations with Dynamically Adaptive Number of Moments

A spectral method for kinetic plasma simulations based on the expansion of the velocity distribution function in a variable number of Hermite polynomials is presented. The method is based on a set of non-linear equations that is solved to determine the coefficients of the Hermite expansion satisfying the Vlasov and Poisson equations. In this paper, we first show that this technique combines the fluid and kinetic approaches into one framework. Second, we present an adaptive strategy to increase and decrease the number of Hermite functions dynamically during the simulation. The technique is applied to the Landau damping and two-stream instability test problems. Performance results show 21% and 47% saving of total simulation time in the Landau and two-stream instability test cases, respectively.

An asymptotic preserving automatic domain decomposition method for the Vlasov–Poisson–BGK system with applications to plasmas

In this work we present an efficient strategy to deal with plasma physics simulations in which localized departures from thermodynamical equilibrium are present. The method is based on the introduction of intermediate regions which allows smooth transitions between kinetic and fluid zones. In this paper we extend Domain Decomposition techniques, obtained through dynamic coupling and buffer zones, to the study of plasmas and, moreover, we combine them with Asymptotic Preserving and Asymptotically Accurate strategies for the time integration. We use a hybrid scheme in which both kinetic and fluid descriptions are considered and coupled together while the kinetic model is solved by asymptotic preserving and accurate methods, in order to guarantee high efficiency and accuracy in all regimes. The numerical scheme is validated and its performances are analyzed by numerical simulations.

MATHEMATICAL STUDY OF A HYPERBOLIC REGULARIZATION TO ENSURE GAUSS' LAW CONSERVATION IN MAXWELL–VLASOV APPLICATIONS

This paper studies a hyperbolic modification of Maxwell's equations to ensure Gauss' law. This correction was obtained by adding a parameter-dependent new unknown and is of great interest for the numerical simulation in plasma physics since the discretization of the Maxwell–Vlasov system does not grant straightforwardly the physical conservation of the charge. Such problems are encountered while using Particle-In-Cell schemes. In this paper the new proposed system has the interest of still being a Friedrichs' one. Its asymptotic behavior with respect to the parameter and the link between modified and original Maxwell's systems are thus investigated. At last, we look for some boundary conditions, granting the well-posedness of the system. Generalizations of the Silver–Müller condition, perfect electric and magnetic conductors, as well as impedance and admittance representation of materials are detailed.

Entropy Based Information Visualization of Scientific Data

Visualization of large scale time-varying scientific data has been a challenging problem due to their ever-increasing size. Identifying and presenting the most informative (or important) aspects of the data plays an important role in facilitating an efficient visualization. In this paper, an information assisted method is presented to locate temporal and spatial data containing salient physical features and accordingly accelerate the visualization process. Based on information measures, the method adaptively picks up important time steps and sub-regions with the maximum information content so that the time-varying data can be effectively visualized in limited time or using limited resources without loss of potential useful physical features. The experiments on the data of radiation diffusion dynamics and plasma physics simulation demonstrate the effectiveness of the proposed method. The method can remarkably improve the way in which scientists analyze and understand large scale time-varying scientific data

Fast parallel Particle-To-Grid interpolation for plasma PIC simulations on the GPU

Particle-In-Cell (PIC) methods have been widely used for plasma physics simulations in the past three decades. To ensure an acceptable level of statistical accuracy relatively large numbers of particles are needed. State-of-the-art Graphics Processing Units (GPUs), with their high memory bandwidth, hundreds of SPMD processors, and half-a-teraflop performance potential, offer a viable alternative to distributed memory parallel computers for running medium-scale PIC plasma simulations on inexpensive commodity hardware. In this paper, we present an overview of a typical plasma PIC code and discuss its GPU implementation. In particular we focus on fast algorithms for the performance bottleneck operation of Particle-To-Grid interpolation.

Unstructured finite element solvers in gyrokinetic turbulence simulations of burning plasmas

Many plasma physics simulation applications utilize an unstructured finite element method (FEM) discretization of an elliptic operator, of the form −∇2u + αu = f; when α is equal to zero a "pure" Laplacian or Poisson equation results and when α is greater than zero a Helmholtz equation is produced. Discretized equations of this form, often 2D, occur in many tokamak fusion plasma, or burning plasma, applications — from MHD to Gyrokinetic codes. This report investigates the performance characteristics of basic classes of linear solvers (ie, direct, one-level iterative, and multilevel iterative methods) on 2D unstructured FEM problems of the form −∇2u + αu = f, with both α = 0 and α ≠ 0. The purpose of this work is to provide computational physicists guidelines as to appropriate linear solvers for their problems via detailed performance analysis, in terms of both the scalability and the constants in the solution times. We show, as expected, almost perfect scalability of multilevel methods and quantify the solution costs on a common computational platform — the IBM SP Power3.

Boris push with spatial stepping

The Boris push is commonly used in plasma physics simulations because of its speed and stability. It is second-order accurate, requires only one field evaluation per time step, and has good conservation properties. However, for accelerator simulations it is convenient to propagate particles in z down a changing beamline. A ‘spatial Boris push’ algorithm has been developed which is similar to the Boris push but uses a spatial coordinate as the independent variable, instead of time. This scheme is compared to the fourth-order Runge–Kutta algorithm, for two simplified muon beam lattices: a uniform solenoid field, and a ‘FOFO’ lattice where the solenoid field varies sinusoidally along the axis. Examination of the canonical angular momentum, which should be conserved in axisymmetric systems, shows that the spatial Boris push improves accuracy over long distances.