Firstly, a brief overview will be given on different models that are able to describe the behaviour of wave propagation as a function of specific frequency ranges. Each range corresponds to different heating systems, namely, 20–100 MHz for the ion cyclotron resonant heating, 2–20 GHz for lower-hybrid heating or current drive, and 100–250 GHz for electron cyclotron resonant heating or current drive systems. The specification of every system will be explained in detail, including the typical set of equations and the assumptions needed to describe the properties of these heating or current drive systems, as well as their specific domains of validity. In these descriptions, special attention will be paid to the boundary conditions. A review of specific physical problems associated with the wave heating systems will also be provided. The review will detail the role of simulation in answering questions that arise from experiments on magnetized plasma devices devoted to fusion. A few examples that will be covered are the impact of edge turbulence on wave propagation and its consequences on heating system performance, the effects of fast particles and ponderomotive effects, among others. A study that is more focused on radio-frequency sheath effects will also be discussed. It shows that such simulations require very sophisticated tools to gain a partial understanding of the observations undertaken in dedicated experiments. To conclude this review, an overview will be given about the requirements and progress necessary to obtain relevant predictive simulation tools able to describe the wave heating systems used in fusion devices.

These notes are based on lectures given by one of us (J.B.T.) at the University of Texas in Austin in 1991. Part I concerns some basic features of plasma confinement by magnetic fields as an introduction to an account of plasma relaxation in Part II. Part III discusses confinement by magnetic mirrors, especially minimum-$B$ systems. It also includes a general discussion of adiabatic invariants and of the principle of maximal ordering in perturbation theory. Part IV is devoted to the analysis of perturbations in toroidal plasmas and the stability of ballooning modes.

The transport of energetic electrons is sensitive to magnetic perturbations. By using three-dimensional numerical simulation of test particle drift orbits we show that the transport of untrapped electrons through an open region with magnetic perturbations cannot be described by a diffusive process. Based on our test particle simulations, we propose a model that leads to an exponential loss of particles.

In this work we numerically demonstrate both significant transient (i.e. non-modal) linear amplification and sustained nonlinear turbulence in a kinetic plasma system with no unstable eigenmodes. The particular system considered is an electrostatic slab with magnetic shear, kinetic electrons and ions, weak collisions and a density gradient, but with no temperature gradient. In contrast to hydrodynamic examples of non-modal growth and subcritical turbulence, here there is no sheared flow in the equilibrium. Significant transient linear amplification is found when the magnetic shear and collisionality are weak. It is also demonstrated that nonlinear turbulence can be sustained if initialized at sufficient amplitude. We prove that these two phenomena are related: when sustained turbulence occurs without unstable eigenmodes, states that are typical of the turbulence must yield transient linear amplification of the gyrokinetic free energy.

In this paper, we present the guiding-centre transformation of the radiation–reaction force of a classical point charge travelling in a non-uniform magnetic field. The transformation is valid as long as the gyroradius of the charged particles is much smaller than the magnetic field non-uniformity length scale, so that the guiding-centre Lie-transform method is applicable. Elimination of the gyromotion time scale from the radiation–reaction force is obtained with the Poisson-bracket formalism originally introduced by Brizard (Phys. Plasmas, vol. 11, 2004, 4429–4438), where it was used to eliminate the fast gyromotion from the Fokker–Planck collision operator. The formalism presented here is applicable to the motion of charged particles in planetary magnetic fields as well as in magnetic confinement fusion plasmas, where the corresponding so-called synchrotron radiation can be detected. Applications of the guiding-centre radiation–reaction force include tracing of charged particle orbits in complex magnetic fields as well as the kinetic description of plasma when the loss of energy and momentum due to radiation plays an important role, e.g. for runaway-electron dynamics in tokamaks.

The particle-in-cell (PIC) algorithm is the most popular method for the discretisation of the general 6D Vlasov–Maxwell problem and it is widely used also for the simulation of the 5D gyrokinetic equations. The method consists of coupling a particle-based algorithm for the Vlasov equation with a grid-based method for the computation of the self-consistent electromagnetic fields. In this review we derive a Monte Carlo PIC finite-element model starting from a gyrokinetic discrete Lagrangian. The variations of the Lagrangian are used to obtain the time-continuous equations of motion for the particles and the finite-element approximation of the field equations. The Noether theorem for the semi-discretised system implies a certain number of conservation properties for the final set of equations. Moreover, the PIC method can be interpreted as a probabilistic Monte Carlo like method, consisting of calculating integrals of the continuous distribution function using a finite set of discrete markers. The nonlinear interactions along with numerical errors introduce random effects after some time. Therefore, the same tools for error analysis and error reduction used in Monte Carlo numerical methods can be applied to PIC simulations.

We construct nonlinear toroidal equilibria of fixed diverted boundary shaping with reversed magnetic shear and flows parallel to the magnetic field. The equilibria have hole-like current density and the reversed magnetic shear increases as the equilibrium nonlinearity becomes stronger. Also, application of a sufficient condition for linear stability implies that the stability is improved as the equilibrium nonlinearity correlated to the reversed magnetic shear gets stronger with a weaker stabilizing contribution from the flow. These results indicate synergetic stabilizing effects of reversed magnetic shear, equilibrium nonlinearity and flow in the establishment of Internal Transport Barriers (ITBs).

In this work we outline an approach to the problem of local equilibrium in non-axisymmetric configurations that adheres closely to Miller's original method for axisymmetric plasmas (Miller et al. 1998 Phys. Plasmas5, 973). Importantly, this method is novel in that it allows not only specification of 3D shape, but also explicit specification of the shear in the 3D shape. A spectrally-accurate method for solution of the resulting nonlinear partial differential equations is also developed. We verify the correctness of the spectral method, in the axisymmetric limit, through comparisons with an independent numerical solution. Some analytic results for the two-dimensional case are given, and the connection to Boozer coordinates is clarified.

Gyrokinetic simulations of electromagnetic turbulence in magnetically confined torus plasmas including tokamak and heliotron/stellarator are reviewed. Numerical simulation of turbulence in finite beta plasmas is an important task for predicting the performance of fusion reactors and a great challenge in computational science due to multiple spatio-temporal scales related to electromagnetic ion and electron dynamics. The simulation becomes further challenging in non-axisymmetric plasmas. In finite beta plasmas, magnetic perturbation appears and influences some key mechanisms of turbulent transport, which include linear instability and zonal flow production. Linear analysis shows that the ion-temperature gradient (ITG) instability, which is essentially an electrostatic instability, is unstable at low beta and its growth rate is reduced by magnetic field line bending at finite beta. On the other hand, the kinetic ballooning mode (KBM), which is an electromagnetic instability, is destabilized at high beta. In addition, trapped electron modes (TEMs), electron temperature gradient (ETG) modes, and micro-tearing modes (MTMs) can be destabilized. These instabilities are classified into two categories: ballooning parity and tearing parity modes. These parities are mixed by nonlinear interactions, so that, for instance, the ITG mode excites tearing parity modes. In the nonlinear evolution, the zonal flow shear acts to regulate the ITG driven turbulence at low beta. On the other hand, at finite beta, interplay between the turbulence and zonal flows becomes complicated because the production of zonal flow is influenced by the finite beta effects. When the zonal flows are too weak, turbulence continues to grow beyond a physically relevant level of saturation in finite-beta tokamaks. Nonlinear mode coupling to stable modes can play a role in the saturation of finite beta ITG mode and KBM. Since there is a quadratic conserved quantity, evaluating nonlinear transfer of the conserved quantity from unstable modes to stable modes is useful for understanding the saturation mechanism of turbulence.

Understanding the key processes occurring in the tokamak scrape-off layer (SOL) is becoming of the outmost importance while we enter the ITER era and we move towards the conception of future fusion reactors. By controlling the heat exhaust, by playing an important role in determining the overall plasma confinement, and by regulating the impurity level in tokamak core, the dynamics of the fusion fuel in the SOL is, in fact, related to some of the most crucial issues that the fusion program is facing today. Because of the limited diagnostic access and in view of predicting the SOL dynamics in future devices, simulations are becoming crucial to address the physics of this region. The present paper, which summarizes the lecture on SOL simulations that was given during the 7th ITER international school (August 25–29, 2014, Aix-en-Provence, France), provides a brief overview of the simulation approaches to the SOL dynamics. First, disentangling the complexity of the system, the key physics processes occurring in the SOL are described. Then, the different simulation approaches to the SOL dynamics are presented, from first-principles kinetic and fluid models, to the phenomenological analysis.

The dynamics of burning plasma is very complicated physics, which is dominated by multi-scale and multi-physics phenomena. To understand such phenomena, numerical simulations are indispensable. Fundamentals of numerical methods used in fusion science numerical modeling are briefly discussed in this paper. In addition, the parallelization technique such as open multi processing (OpenMP) and message passing interface (MPI) parallel programing are introduced and the loop-level parallelization is shown as an example.

To rigorously model fast ions in fusion plasmas, a non-Maxwellian equilibrium distribution must be used. In this work, the response of high-energy alpha particles to electrostatic turbulence has been analyzed for several different tokamak parameters. Our results are consistent with known scalings and experimental evidence that alpha particles are generally well confined: on the order of several seconds. It is also confirmed that the effect of alphas on the turbulence is negligible at realistically low concentrations, consistent with linear theory. It is demonstrated that the usual practice of using a high-temperature Maxwellian, while previously shown to give an adequate order-of-magnitude estimate of the diffusion coefficient, gives incorrect estimates for the radial alpha particle flux, and a method of correcting it in general is provided. Furthermore, we see that the timescales associated with collisions and transport compete at moderate energies, calling into question the assumption that alpha particles remain confined to a flux surface that is used in the derivation of the slowing-down distribution.

We present a comprehensive survey of various computational methods in CEDRES++ (Couplage Equilibre Diffusion Résistive pour l'Etude des Scénarios) for finding equilibria of toroidal plasma. Our focus is on free-boundary plasma equilibria, where either poloidal field coil currents or the temporal evolution of voltages in poloidal field circuit systems are given data. Centered around a piecewise linear finite element representation of the poloidal flux map, our approach allows in large parts the use of established numerical schemes. The coupling of a finite element method and a boundary element method gives consistent numerical solutions for equilibrium problems in unbounded domains. We formulate a new Newton method for the discretized nonlinear problem to tackle the various nonlinearities, including the free plasma boundary. The Newton method guarantees fast convergence and is the main building block for the inverse equilibrium problems that we can handle in CEDRES++ as well. The inverse problems aim at finding either poloidal field coil currents that ensure a desired shape and position of the plasma or at finding the evolution of the voltages in the poloidal field circuit systems that ensure a prescribed evolution of the plasma shape and position. We provide equilibrium simulations for the tokamaks ITER and WEST to illustrate the performance of CEDRES++ and its application areas.

Superbanana and superbanana plateau transport processes are critical to plasma confinement in tokamaks with broken symmetry. The transport is caused by the superbanana resonance, which occurs at a pitch angle that makes the toroidal drift speed vanish, i.e. the tips of the superbananas. The physics consequences of the resonance on the symmetry breaking induced toroidal momentum damping and on the energetic alpha particle transport have been demonstrated using large aspect ratio expansion. Here, the existing theory for the superbanana and superbanana plateau transport is extended for finite aspect ratio tokamaks with broken symmetry. The effects of finite plasma β, and magnetic field shear are naturally included. Here, β is the ratio of the thermal plasma pressure to the magnetic field pressure. The explicit expressions for the transport fluxes in these regimes in terms of the equilibrium quantities are presented. It is shown that the main effects are to modify the resonance function G(k) and the expression for the pitch angle parameter k in the existing theory.

BOUT++ is a 3D nonlinear finite-difference plasma simulation code, capable of solving quite general systems of Partial Differential Equations (PDEs), but targeted particularly on studies of the edge region of tokamak plasmas. BOUT++ is publicly available, and has been adopted by a growing number of researchers worldwide. Here we present improvements which have been made to the code since its original release, both in terms of structure and its capabilities. Some recent applications of these methods are reviewed, and areas of active development are discussed. We also present algorithms and tools which have been developed to enable creation of inputs from analytic expressions and experimental data, and for processing and visualisation of output results. This includes a new tool Hypnotoad for the creation of meshes from experimental equilibria. Algorithms have been implemented in BOUT++ to solve a range of linear algebraic problems encountered in the simulation of reduced Magnetohydrodynamics (MHD) and gyro-fluid models: A preconditioning scheme is presented which enables the plasma potential to be calculated efficiently using iterative methods supplied by the PETSc library (the Portable, Extensible Toolkit for Scientific Computation) (Balay et al. 2014), without invoking the Boussinesq approximation. Scaling studies are also performed of a linear solver used as part of physics-based preconditioning to accelerate the convergence of implicit time-integration schemes.

Ripples in magnetic or electrostatic confinement fields give rise to trapping separatrices, and conventional neoclassical transport theory describes the collisional trapping/detrapping of particles with fractured distribution function. Our experiments and novel theory have now characterized a new kind of neoclassical transport processes arising from chaotic (nominally collisionless) separatrix crossings, which occur due to E × B plasma rotation along θ−ruffled or wave-perturbed separatrices. This chaotic neoclassical transport becomes dominant at low collisionality when the collisional spreading of particle energy during the dynamical period is less than the separatrix energy ruffle.

In previous works, we derived stationary density distribution functions (DDF) where the local equilibrium is determined by imposing the maximum entropy (MaxEnt) principle, under the scale invariance restrictions, and the minimum entropy production theorem. In this paper we demonstrate that it is possible to reobtain these DDF solely from the MaxEnt principle subject to suitable scale invariant restrictions in all the variables. For the sake of concreteness, we analyse the example of ohmic, fully ionized, tokamak-plasmas, in the weak-collisional transport regime. In this case we show that it is possible to reinterpret the stationary distribution function in terms of the Prigogine distribution function where the logarithm of the DDF is directly linked to the entropy production of the plasma. This leads to the suggestive idea that also the stationary neoclassical distribution functions, for magnetically confined plasmas in the collisional transport regimes, may be derived solely by the MaxEnt principle.

A reduced four-dimensional (integrated over perpendicular velocity) gyrokinetic model of slab ion temperature gradient-driven turbulence is used to study the phase-space scales of free energy dissipation in a turbulent kinetic system over a broad range of background gradients and collision frequencies. Parallel velocity is expressed in terms of Hermite polynomials, allowing for a detailed study of the scales of free energy dynamics over the four-dimensional phase space. A fully spectral code – the DNA code – that solves this system is described. Hermite free energy spectra are significantly steeper than would be expected linearly, causing collisional dissipation to peak at large scales in velocity space even for arbitrarily small collisionality. A key cause of the steep Hermite spectra is a critical balance – an equilibration of the parallel streaming time and the nonlinear correlation time – that extends to high Hermite number n. Although dissipation always peaks at large scales in all phase space dimensions, small-scale dissipation becomes important in an integrated sense when collisionality is low enough and/or nonlinear energy transfer is strong enough. Toroidal full-gyrokinetic simulations using the Gene code are used to verify results from the reduced model. Collision frequencies typically found in present-day experiments correspond to turbulence regimes slightly favoring large-scale dissipation, while turbulence in low-collisionality systems like ITER and space and astrophysical plasmas is expected to rely increasingly on small-scale dissipation mechanisms. This work is expected to inform gyrokinetic reduced modeling efforts like Large Eddy Simulation and gyrofluid techniques.