The aim of the conference is to provide a discussion forum covering all areas of fusion-oriented theoretical activities.
Its main topics are:
To reach ITER’s mission goals, a thorough understanding of the underlying physical processes, plant systems and operational scenarios is mandatory. Preparation for ITER operations is supported by state-of-the-art modelling implemented in the ITER Integrated Modelling & Analysis Suite (IMAS). This modular framework has been developed to standardize the communication between the various codes describing different ITER subsystems, thus enabling modelling studies to optimize the performance of ITER. An important example of the advanced physics modelling capability recently implemented in IMAS is the Heating and Current Drive (H&CD) workflow, capable to describe all of ITER’s heating systems and synergies between them. The standardization imposed by the IMAS platform also allows a flexible coupling of the H&CD workflow with various transport solvers, particularly important to study burning plasmas in ITER with dominant alpha particle heating. This presentation starts by introducing briefly the ITER Research Plan, together with ongoing modelling efforts. The talk will illustrate various capabilities of the physics workflows in IMAS towards the development of a high-fidelity plasma simulator for ITER. Finally, examples of physics studies for scenario optimization in the non-active and D-T operational phases of ITER will be shown.
Good confinement of the bulk plasma and fusion-generated alpha particles are two basic design properties of a fusion reactor. First, small radial energy fluxes are necessary for the plasma to achieve fusion-relevant conditions. In turn, fusion-born alpha particles are expected to contribute to heat the plasma, which requires their confinement time to be sufficiently long. In stellarators, neoclassical processes are a major concern with respect to bulk and energetic ion confinement. In a generic three-dimensional device, trapped particles move back and forth along the field lines with a non-zero average radial drift. This, in combination with collisions, gives place to a variety of stellarator-specific neoclassical transport regimes of which the most relevant ones for low-collisionality bulk plasma transport are the 1/ν, sqrt(ν) and superbananaplateau regimes. The same processes and regimes can be observed in tokamaks with broken symmetry. In most fusion-revelant cases, bulk ions are typically in the sqrt(v) regime, and the superbanana orbits are the mechanism behind prompt energetic ion losses. Stellarator optimization is the numerical procedure by which the magnetic configuration is tailored to meet several design criteria, and neoclassical transport of bulk and energetic ions are two of them. The effective ripple, a figure of merit for the level of transport in the 1/ν regime, is usually addressed; the loss fraction of energetic ions is targeted by means of simplified proxies. In this talk, we present the recently developed code KNOSOS (KiNetic Orbit-averaging SOlver for Stellarators), a freely-available open-source code that provides a fast computation of low collisionality neoclassical transport in three-dimensional magnetic confinement devices by rigorously solving the radially local bounce-averaged drift kinetic equation coupled to the quasineutrality equation. In the first part of the talk, we show that KNOSOS reproduces the results of the standard neoclassical code DKES and can be orders of magnitude faster. For this reason, KNOSOS can provide new figures of merit for stellarator optimization for bulk neoclassical transport that could not be considered before, such as the level of transport in the sqrt(ν) regime and the superbana-plateau regime. The latter can be described because the component of the magnetic drift tangent to flux surfaces is rigorously included in the equations.In the second part of the talk, we go one step further and show that, by also keeping the radial component of the magnetic drift in the particle orbits, the resulting radially global bounceaverageddrift-kinetic equation can be employed to accurately describe the neoclassical transport of energetic ions in stellarators. This is demonstrated by comparing the solutions of this equation against guiding-center and full-orbit simulations using the code ASCOT. The bounce average reduces the dimensionality of the equation and could lead to potentially faster calculations, that would be specially suited for stellarator optimization and parameter scans.
Simulating global plasma instabilities and turbulence in stellarators with gyrokinetic codes is significantly more complicated and computationally expensive than in tokamaks and, partly because of this, the field is significantly less developed in stellarators than in the tokamak counterpart. The main reason is that the three-dimensional geometry of stellarators makes the flux tube model unsuitable for these devices. Consequently, simulating instabilities and turbulence in stellarators requires covering at least a full surface or a (full-surface and) radially global computational domain. The latter is the most complete computational domain and the most expensive in terms of the computational cost of the simulations. In order to reach a predictive capability, simulation codes require their verification through comparisons between different codes and against theoretical models and their validation against experimental measurements. In this contribution, we address the verification of two global gyrokinetic delta-f codes specifically designed for three-dimensional geometries, EUTERPE and GENE-3D. EUTERPE was developed as a linear code and afterward was extended incorporating non-linear saturation mechanisms, different models for kinetic electrons and impurities, and electromagnetic effects. GENE-3D has recently been extended to include electromagnetic effects. For this verification activity, we use two magnetic configurations, the standard configuration of the helical device LHD and a configuration of the optimized stellarator W7-X. We tackle several problems with different levels of complexity, from linear simulations with adiabatic electrons to non-linear or multi-species simulations. In addition, we also study two problems that are of particular relevance for stellarators and for which the global codes are well suited: the spatial localization of instabilities and turbulence and the influence of a long-wavelength radial electric field on them. A thorough comparison of results from both codes is performed, finding a very good agreement between them.
Three dimensional magnetic equilibria are in general composed of nested flux surfaces, magnetic islands and chaotic field lines, although it is possible to design stellarator coil configurations that produce vacuum fields with nested flux surfaces (Pedersen, S. T. et al. 2016, Nature comm.). At finite $\beta$ however, currents self-generated by the plasma, such as diamagnetic, Pfirsch-Schlüter or bootstrap, perturb the magnetic field, thus breaking nested flux surfaces and ultimately impairing confinement. To date, there is no theory, nor extensive numerical study that characterizes the maximum achievable $\beta$ above which magnetic surfaces are destroyed, nor a theory on the dependency of this critical $\beta$ on other relevant operational parameters. We propose using the Stepped Pressure Equilibrium Code (SPEC) (Hudson, S. R. et al. 2012, Phys. of Plasmas), which can compute 3-dimensional stepped-pressure equilibria with magnetic islands and chaos, to study the effect of finite $\beta$ on the magnetic topology of stellarators. Recent numerical work significantly improved the speed and robustness of SPEC (Qu, Z. et al, 2020, Plasma Phys. Cont. Fusion), which allows large parameter scans in a reasonable amount of time (Loizu, J. et al. 2017, J. Plasma Phys.). In addition, SPEC has recently been extended to allow free-boundary calculations (Hudson, S. R. et al. 2020, Plasma Phys. Cont. Fusion) with prescribed net toroidal current profiles (Baillod, A. et al. 2021, J. Plasma Phys.). Leveraging these new capabilities, we present here the first extensive and comprehensive study of the equilibrium $\beta$-limit with bootstrap current. We consider a number of representative configurations, such as classical, quasi-axisymmetric and quasi-helically symmetric stellarators.
Local gyrokinetic simulations use a field-aligned domain that twists due to the magnetic shear of the background magnetic equilibrium. However, if the magnetic shear is strong and/or the domain is long, the twist can become so extreme that it fails to properly resolve the turbulence. In this work, we derive and implement the non-twisting flux tube, a local simulation domain that remains rectangular at all parallel locations. Convergence and run time tests indicate that it can calculate the heat flux more efficiently than the conventional flux tube. For one test case, it was 30 times less computationally expensive and we found no case for which it was more expensive. It is most advantageous when the magnetic shear is high and the domain includes at least two regions of turbulent drive, which makes it potentially useful for pedestal simulations, stellarator simulations, and tokamak simulations with several poloidal turns. It also more accurately models the inboard midplane when the magnetic shear is large. Additionally, we show how the non-twisting flux tube can be generalized to allow further optimization and control of the simulation domain. Lastly, we will explore the possibility of using the non-twisting flux tube to elegantly and efficiently include arbitrary radial variation in the safety factor profile, thereby blurring the distinction between local and global simulations.
Given the large anisotropy of transport processes in magnetized plasmas, the magnetic field structure can strongly impact the heat diffusion: magnetic surfaces and cantori form barriers to transport while chaotic layers and island structures can degrade confinement. When a small but finite amount of perpendicular diffusion is included, the structure of the magnetic field becomes less important, allowing finite pressure gradients to be supported across chaotic regions and island chains. We introduce a metric for the effective volume of non-integrability based on the solution to the anisotropic heat diffusion equation. To validate this metric, we consider model fields with a single island chain and a strongly chaotic layer for which analytic predictions of the relative parallel and perpendicular transport can be made. We also analyze critically chaotic fields produced from different sets of perturbations, highlighting the impact of the mode number spectrum on the heat transport. We propose that this metric be used to assess the impact of non-integrability on the heat transport in stellarator equilibria.
Minimising the presence and size of magnetic islands in the core of stellarators is necessary to retain the advantageous confinement properties of nested flux surfaces. Therefore, efficient optimisation schemes to reduce island size in stellarator configurations are desirable. However, a configuration with small islands may not be good enough: it may be sensitive to the exact positions of the coils, which also tend to have complex shapes that make them more difficult to manufacture precisely. Therefore, it is also crucial to reliably and efficiently quantify the sensitivity of island size with respect to coils. For example, in NCSX tight coil tolerances — also related to island sensitivity — inflated the construction cost, which was one of the reasons that led to the experiment being cancelled.An adjoint approach to calculating the gradient of two quantities related to magnetic island size is presented. The first quantity is the full radial width of an island calculated using a method developed by Cary and Hanson, which relies on the island size being small compared to the system size. The island width calculation is verified using an analytical magnetic field configuration. The second quantity, known as Greene’s residue, is simpler to calculate and can be calculated for any periodic field line; for an island centre (O point) the residue is strongly correlated with the island width. Convergence tests for the gradients of island width and residue are presented. The result of a gradient-based optimisation of residues in a magnetic configuration produced by a pair of helical coils are presented. To our knowledge this is the first optimisation for good flux surfaces using analytic derivatives. Furthermore, the shape gradient— a measure of sensitivity — of the island width with respect to coils is calculated and numerically verified for a magnetic island in the NCSX configuration. The applications and extensions of this work in the area of stellarator optimisation and design are discussed.
Direct numerical solution of the nonlinear Boltzmann transport equation remains elusive nearly 150 years after its discovery. Appropriate approximations continue to serve as the foundation of aerodynamics and plasma physics. However, some important problems in fusion don’t lend themselves to such approximations, such as runaway electron avalanche and regimes of dense neutral populations in scrape-off layers. Routine and robust solution of the nonlinear Boltzmann equation will be needed to close these gaps in computational capability. Several approaches have been successfully developed throughout the 20th century to numerically solve the Boltzmann equation in various limits, particularly in the case of linear collision operators. These include the Monte Carlo technique used in EIRENE and DEGAS2, and the weakly anisotropic expansion used in BOLSIG+. Advances in applied mathematics now make routine numerical solution of the nonlinear equation feasible. For example, a Fourier representation has been recently applied to find new fluid closures using machine learning. This talk focuses on a new computational framework (LightningBoltz), which is based on a generalization of the conservative spectral method of Gamba, Rjasanow, and Keßler. This framework takes ample advantage of pre-computing the discrete collision operator in a distributed computing system. Depending on the demands of the specific problem, this capability allows the full bilinear Boltzmann operator to be calculated in fractions of a millisecond on a laptop, while time-dependent and spatially-inhomogeneous kinetic problems can be readily solved on a workstation or small cluster. Other features include: force-field acceleration, implicit time-stepping, and arbitrary cross sections read from atomic physics databases. LightningBoltz has been verified against constructed analytic solutions to the Boltzmann equation and against the Chapman-Cowling expansion at low Knudsen number. Cross- verification with the DEGAS2, Gkeyll, and BOLSIG+ codes is presented. A proof-of-principle for a one-dimensional model of a tokamak scrape-off layer is shown, which includes rigorous treatment of elastic scattering: a unique capability. The spectral technique is also being applied to the coupling between DEGAS2 and the gyrokinetic edge code XGC to rigorously account for elastic scattering and enforce conservation in charge-exchange interactions. Runaway electron distributions pose unique challenges due to the disparate velocity scales, and potential paths to tackling the nonlinear kinetic avalanche problem will be discussed.
Magnetic fusion is now approaching next-step D-T burning plasma experiments in ITER, which will be mostly self-heated by fusion-born alpha-particles. Burning plasmas will be a highly nonlinear medium, and predicting with confidence the alpha-particle heating and alpha-particle losses in such plasmas is a very challenging task. To advance it, dedicated experimental studies are being performed on present-day machines worldwide with fast ions matching the key ITER dimensionless parameters βfast/βthermal, Vfast/VA, R0·d(βfast)/dr. Progress with experimental studies is closely coupled with the development and validation of theoretical analyses and numerical simulations, together with their extrapolation to future burning plasmas. Among all alpha-particle effects, instabilities excited by wave-particle resonances with the super-Alfvénic alpha-particles constitute the largest uncertainty in the predictions for ITER. Thus, studying of fast-ion driven instabilities is one of the highest priority experimental and theoretical research tasks. In Europe, fast particle experimental studies are mostly performed on the largest tokamak JET, and on several medium-size tokamaks: ASDEX-Upgrade (AUG), TCV, MAST, and on stellarators Wendelstein 7-X and TJ-II. The JET machine leads the research on fast ions in the MeV energy range. In view of the upcoming D-T campaign (2021), JET has performed dedicated experiments on scenarios for alpha-particle-driven Alfvén eigenmodes (AEs), and on investigating alpha-particles born in aneutronic D-3He fusion. A rich variety of AEs with frequencies ranging from ~50 kHz to ~500 kHz were excited in recent D-3He experiments on JET, part of which were attributed to alpha-particles with a bump-on-tail distribution self-sustained by periodic modulation of the fusion source with sawteeth. The medium-size tokamaks AUG and TCV have contributed to the development of novel methods for controlling AEs with ECRH/ECCD, which could be used eventually for burning plasma control. Finally, both the spherical tokamak MAST and the TJ-II stellarator have observed beam-driven Alfvén instabilities with fast frequency sweeping (the “chirping modes”), successfully interpreted in the framework of the non-linear Berk-Breizman theory as discussed in this talk. First observations of fast-ion driven instabilities have been also reported from Wendelstein 7-X. The wealth of fast-ion experimental data described above was validated with fast-ion modelling tools. High-confidence exists now for the predictive capabilities of linear spectral MHD codes used in EUROfusion such as CASTOR-K, MISHKA, and HALO to compute the spatial structure of AEs and their frequencies. This very accurate computational analysis has resulted in a new technique ‘Alfvén spectroscopy’, providing valuable information on the evolution of the q-profile and plasma density from the observed AE spectrum. The MHD codes are routinely used to identify the wave-particle resonances and compute energetic-ion drive mechanisms. Concurrently, kinetic codes such as LIGKA have been advanced to provide fairly accurate estimates of the kinetic damping effects, some of which are exponentially sensitive to plasma parameters. In this talk, we illustrate that strongly nonlinear codes have become an accurate tool to describe the dynamics of energetic-particle modes well above the threshold of the mode excitation. We also discuss further developments and needs in fast-ion theory, modelling and experiment towards a better understanding of AEs in burning plasmas. These include studies of multiple fast-ion driven modes with high toroidal mode numbers as expected for ITER. The resonance overlap and global stochasticity of fast-ion orbits could play a dominant role in the case of multiple modes, but such conditions have not been sufficiently studied yet. This advancement will require a dedicated effort both in developing appropriate diagnostics and modelling tools.
Exploiting mechanisms to actively control and reduce turbulence is essential for maximizing the performance of fusion devices. Already in 2010, it became clear that low levels of turbulent transport observed in several JET experiments could not be explained by the usual turbulence reduction mechanisms such as strong plasma rotation and magnetic shear. This talk summarizes main analysis and modeling steps that have been undertaken in the period 2010-2020 to reveal new physical mechanisms responsible for the observed turbulence modification and transport reduction on JET. Initial results clearly indicated that strong pressure gradients, notably from fast ions, can have a beneficial impact on the reduction of heat transport in high-beta JET plasmas. As a next step the stabilization of turbulent transport in low-rotation JET plasmas heated mostly by moderately energetic (~100 keV) ICRF-generated fast ions was demonstrated. This progress has recently led to the important insight that under certain conditions alpha particles can also significantly contribute to stabilize ITG turbulence and reduce heat transport in D-T plasmas of ITER. This result is particularly important for achieving and maintaining high ion temperatures in a fusion reactor, characterized by dominant core electron heating from fusion-born alphas. Furthermore, highly energetic alpha particles can also excite Alfvénic instabilities, commonly considered as detrimental for plasma confinement. In this talk, we will discuss a novel type of turbulence suppression and improved energy confinement in JET plasmas with MeV-range fast ions and strong fast-ion driven Alfvénic eigenmodes (AEs). The detailed gyrokinetic analysis of recent dedicated fast-ion experiments at JET shows that MeV ions trigger nonlinearly large-scale zonal flows both in the electrostatic and electromagnetic potentials, leading to the suppression of the electrostatic transport and leaving electromagnetic fluctuations as the main source of weak electron transport. This mechanism is rather efficient for raising ion temperature in fusion plasmas with a large population of fast ions and collisional electron heating, and hold promises for D-T plasmas with alpha particle heating. The talk is concluded with the discussion of the implication of these results for the upcoming D-T campaign on JET and future ITER operations.
The hybrid tokamak scenario is characterized by low magnetic shear in the plasma core and a central value of the safety factor close to unity. It represents a hybrid between standard scenarios and advanced scenarios and is a candidate scenario for ITER and DEMO. The hybrid scenario allows for high-performance, sawtooth-free operation with extended discharge lengths and has the advantage that its characteristic safety factor profile is automatically maintained. The latter property of hybrid discharges is due to a self-regulating current redistribution mechanism called magnetic flux pumping which had not been understood yet. Current diffusion calculations that have been performed for hybrid discharges falsely predict values of the central safety factor below unity and hence sawtoothing. Understanding this effect is crucial in order to extrapolate the accessibility and properties of the scenario to future tokamaks. Based on 3D nonlinear MHD simulations of tokamak plasmas performed with the M3D-C1 code, we propose an explanation for magnetic flux pumping. In these simulations, a saturated quasi-interchange instability creates helical (m=1,n=1) perturbations of the magnetic and velocity fields in the central region of the plasma. The perturbations combine via an MHD dynamo to give an effective loop voltage flattening the background current density profile in the plasma core. This mechanism is self-regulating and prevents sawtoothing by keeping the central safety factor profile flat and close to unity. Since the quasi-interchange instability is pressure-driven, the maximal amount of flux pumping that can be provided by the dynamo loop voltage effect scales with the pressure. The beta threshold for the avoidance of sawteeth depends on how much the current density is being peaked centrally, e.g. by central current drive. In ASDEX Upgrade tokamak discharges which have been set up to test this model, measurement results qualitatively agree with these theoretical predictions. In these discharges, positive ECCD has been applied in several steps to drive the safety factor on axis to lower values, while at the same time an NBI power scan has been performed to increase beta, resulting in an alternation between sawtoothing and sawtooth-free phases. During the sawtooth-free phases, experimental evidence for anomalous current redistribution by a continuous (m=1,n=1) mode, leading to q≈1 in the core, is found in accordance with the theoretical model. A quantitative comparison between theory and experiment is ongoing work. To this end, nonlinear MHD simulations using the new full MHD model of the JOREK code have been set up based on selected time points during the different phases of one of the described ASDEX Upgrade discharges.
The trajectories of collisionless charged particles in static magnetic fields are Hamiltonian flows, hence may be understood using tools of Hamiltonian dynamics. In the presence of three independent constants of the motion, it is known that a trajectory is regular. A low-energy particle in an axis-symmetric device possesses two exact constants of motion: the total energy and the momentum along the direction of symmetry; and a third one, adiabatic invariant in reality: the magnetic moment µ. Thus, its trajectory is regular. By increasing the particle energy up to values such that its Larmor radius becomes no longer negligible with respect to ambient scales, one may reach a scenario where µ is no longer conserved. Accordingly, the trajectory may become chaotic. In this work we present a study of the dynamics of fast particles in an axis-symmetric magnetic geometry patterned after that of the NSTX Spherical Tokamak. Trajectories are computed from full Hamilton’s equations of motion, not of the guiding center ones. Thus, no a-priori hypothesis is done about the conservation of the magnetic moment µ, which is instead computed and monitored at the post-processing stage. We find—consistently with earlier studies — that, within a sizable domain of the machine, µ encompasses three regimes: an adiabatic one, where µ is roughly constant with superposed small periodic oscillation at the Larmor frequency, and which is characteristic of the lower end of the energy range; a superadiabatic one, where µ still conserves a constant mean value, but large-scale oscillations at the bounce frequency are superposed; finally, a genuinely chaotic regime, where µ experiences sudden irregular jumps between widely differing mean values. These two latter regimes are characteristic of the energy range of beam ions for NSTX scenario, or alpha particle in fusion-scale devices, with appropriately rescaled fields. The boundary between these two regimes is not sharp: we find instances of trajectories which appear superadiabatic even over times very long with respect to bounce frequency, yet ultimately turn into non-adiabatic. An analysis of the details of the trajectories shows that the breakup of adiabaticity occurs when an orbit with a large enough velocity crosses a region of low magnetic field amplitude and with a large variation of the orientation of the field. We highlight also an interesting crosstalk between the different degrees of freedom: even when µ is not conserved and chaotic orbits are present, the conservation of the toroidal momentum sets strong constraints about the volume available to the particles, and their radial diffusion stays bounded. Our findings have both practical and physical consequences. First, when µ has large oscillations (not necessarily chaotic jumps), the use of guiding-center or gyrokinetic calculations for such orbits is questionable. Second, the capability of these particle to excite Alfvènic instabilities decreases, since the velocity of the particle fluctuates with respect to the phase-velocity of the Alfvèn wave, which imposes a fluctuating sign to the energy exchanges with this wave. Signatures of this effect might have already been evidenced in experiments.
The Guiding Center (GC) theory has been widely used for more than four decades as the basis for the study of single and collective particle dynamics in toroidal magnetic fields utilized in fusion devices. The Hamiltonian formulation of the theory has been originally given in terms of non-canonical variables and then extended to canonical ones. The canonical Hamiltonian description, apart from its elegant structure, is accompanied by an arsenal of powerful mathematical methods of practical importance, which are readily applicable when a final transformation to Action-Angle variables is performed. This transformation is possible for the case of GC motion in axisymmetric equilibria, where the corresponding equations are integrable. All symmetry-breaking perturbations related to non-axisymmetric and/or time-dependent modes and waves can be expressed in terms of the Action-Angle variables of the unperturbed axisymmetric system, facilitating the analytical and numerical study of the collective particle behaviour. More specifically, the transformation provides the Orbital Spectrum of the GC orbits in terms of the bounce/transit frequencies as well as the bounce/transit-averaged toroidal precession and gyration frequencies. The respective Orbital Spectrum Analysis gives an a priori knowledge of the exact location and the strength of the resonant particle interactions with multi-scale spatio-temporal perturbations in the GC phase space. This is of particular importance for understanding and predicting the modifications of the respective distribution functions of bulk and energetic particles and the synergetic effects between different types of perturbations. In addition, the unique advantage of the Action-Angle formulation, related to the clear separation of the different time scales of the particle motion in different degrees of freedom, enables the systematic dynamical reduction of the kinetic description to a hierarchy of evolution equations for the distribution functions which extends standard bounce-averaged and quasilinear kinetic theories. The practical use of the aforementioned advantages has been hindered mostly due to complications on the calculation of the transformation to Action-Angle variables. In fact, analytical results have been reported only for the case of a Large Aspect Ratio (LAR) equilibrium under the Zero Drift Width (ZDW) approximation and only for deeply trapped particles, which is the case described in many textbooks. Our work focuses on raising these restrictions in two ways: We obtain analytical results beyond the deeply trapped approximation, that can be further utilized in analytical perturbation theory. Moreover, we study Finite Drift Width (FDW) effects, which are important for energetic particles, and we consider the transformation to Action-Angle variables and the calculation of the Orbital Spectrum for realistic numerical and experimental equilibria.
The passing electron response to ion-gyroradius-scale instabilities is often considered to be adiabatic: on irrational flux surfaces, passing electrons are assumed to be able to transit the entire flux surface on time scales faster than the mode e-folding growth time as a result of the small electron-to-ion mass ratio. This argument fails on mode-rational flux surfaces where magnetic field lines only cover a subset of the surface. In fact, gyrokinetic simulations have revealed modes driven by the passing electron response in narrow radial layers near rational flux surfaces in both electrostatic and electromagnetic scenarios, with possible impacts on transport. In the ballooning representation, these modes appear with giant tails in extended poloidal angle. The small electron-to-ion mass ratio limit of linear electrostatic gyrokinetics is presented, for ion-gyroradius-scale modes, in the ballooning representation, including the nonadiabatic response of passing electrons. This theory reveals novel modes driven solely by the nonadiabatic passing electron response, and recovers the usual ion and trapped-electron driven instabilities. The collisionless and collisional limits of the theory are considered, and simple scaling predictions for basic properties of the eigenmodes are obtained. The scaling predictions are shown to be in agreement with simulation results from the gyrokinetic code GS2. The extension of the theory to the electromagnetic case is discussed: the extended theory may be relevant to the study of micro-tearing modes. Finally, it is observed that the novel electrostatic passing-electron-response-driven modes may be insensitive to equilibrium flow shear, since these modes form narrow radial layers. Hence, it may be important to consider the stability of passing-electron-response driven modes in projected scenarios where turbulence is assumed to be suppressed by equilibrium flow shear.
In tokamak the parallel divergence of the parallel current is non-zero due to the time variation of the vorticity. This has an interesting connection with the axial anomaly usually invoked in baryogenesis. In two-dimensions (good approximation for the plasma in strong magnetic field) the field of vorticity can be seen as a discrete set of positive and negative elementary vortices, of fixed magnitude and two orientations. Their density is the “matter”, their interaction is mediated by the field of velocity and the formalism is similar to a classical field theory. Consider the presence of an incipient, large scale vortex. The time variation of the topological “winding” of velocity field in the region around it is equivalent to generation of new elementary vortices. They are attracted by the initial vortex and this will grow. The current of vorticity along the main vortex has non-zero divergence which equals the time variation of the rotational of the velocity around (i.e. the time variation of the topology of the velocity field is equivalent to creation of new elementary vortices). One notes the possibility of a mapping to the formalism of the quantum axial anomaly: the main vortex is a scalar string, elementary vortices are fermions and the velocity field is the gauge field with nontrivial topology. In the quantum axial anomaly the topological winding of the gauge field is converted into fermions running along the scalar string. The model described above is parallel and purely classical. It may be useful to describe the processes of concentration of vorticity, in vortices or in the sheared rotation layers of barriers, including the H-mode layer. The two-dimensional model (mentioned above) consisting of discrete elementary vortices has two formulations: (1) ideal fluid and (2) plasma/planetary atmosphere. In (1) the governing equation (Euler) does not exhibit any intrinsic length, it is conformally invariant. In (2) the equation is Charney Hasegawa Mima (CHM) and there is an intrinsic length (Larmor radius/ Rossby radius). The difference is essential. In the CHM case the Ertel’s theorem states that (ω+Ωc)/n is a Lagrangian invariant. In the ideal Euler fluid the connection between the vorticity ω and and density n is lost. The distinction between ideal Euler (conformal) fluid and CHM plasma (fixed Larmor radius) becomes very important for the poloidally rotating plasma. The reason is a factor which occurs systematically in the perturbative treatment of the ion rotation , (ρeff)^(-2) = (ρs)^(-2) (1-vdia/u). We call it “effective Larmor radius” because it replaces in some formula the basic Larmor radius, as if the system would base its dynamics on this modified parameter. It exhibits the fundamental change that occurs when the CHM dynamics reaches asymptotically the Euler dynamics. When vdia is close to u the plasma behaves as Euler (rho-effective is very large) and the connection between the vorticity (shear of the rotation velocity) and the density is reduced. This places a particular limit on the pedestal evolution. In addition, since the ion “effective Larmor radius” is very large, a class of electron vortices (dependent on this condition) is excited, possibly as precursors of ELMs.
Relativistic effects in astrophysical objects and fusion plasmas do not necessarily require extremely high temperatures and energies. They appear to be non-negligible even for electron temperatures $T_e$ of the order of tens of keV, i.e. when $T_e\ll m_ec^2$. Relativistic effects in transport physics appear due to macroscopic features of the relativistic thermodynamic equilibrium given by the Maxwell-J{\"u}ttner distribution function. In fusion devices such as ITER, DEMO and especially in future aneutronic fusion schemes, relativistic effects are surely not important for ions, but electrons require more careful attention. In particular, it was shown recently that relativistic effects in electron transport for temperatures characteristic of fusion reactor leads to noticeably different results for both tokamaks and stellarators compared with the corresponding non-relativistic calculations. Nevertheless, practically all transport codes applied for simulation of the fusion reactor scenarios are based on the non-relativistic approach. In the literature on relativistic plasmas one can find two different approaches. The most rigorous one is based on the covariant formulations of the kinetics and hydrodynamics of plasmas. This formalism is required for the correct treatment of astrophysical objects. Another approach, aimed at the description of the laboratory fusion plasma, does not satisfy Lorentz invariance but takes features of relativistic thermodynamic equilibrium into account. For different reasons, both these approaches are not convenient enough for implementation into reactor transport codes. The present work is focused on transport processes in hot plasmas with relativistic electrons with characteristic velocities of fluxes $V\ll v_{te}$ and $v_{te}< c$. The main goal is to derive Braginskii equations in a weakly relativistic approach for the fluxes, i.e. neglecting terms of the higher order than $V^2/c^2$ and $v_{te}V/c^2$, while the bulk electrons are described as fully relativistic. The final equations are mathematically very similar to the non-relativistic ones and have a transparent physical interpretation. For closure of the transport equations, it is proposed to use generalized Laguerre polynomials of order $3/2+{\mathcal{R}}(T_e)$, where ${\mathcal{R}}(T_e)=(15/8)(T_e/m_e c^2)+{\mathcal{O}}(T_e^2/m_e^2c^4)$. Just like Sonine polynomials in non-relativistic approach, these polynomials give a finite number of terms for the right-hand-side of the relativistic linearized kinetic equation. In contrast to earlier attempts to use Sonine polynomials for relativistic transport, this makes the proposed presentation more convenient and practical.
While the curse of high dimensionality associated with the Fokker-Planck equation can be resolved by Monte Carlo solution to the underlying stochastic process, the distribution of discretization points for the resulting random walk-based method follows the evolution of the target distribution function, leading to the weight-spreading phenomenon in the $\delta f$ method, and therefore increase in statistical sampling error. In order to overcome this issue, we propose a Monte Carlo solution algorithm to the Fokker-Planck equation with fixed particles in phase space using the method of fundamental solutions. The numerical challenges introduced by the integral solution in resolving small time step sizes are addressed and dealt with by approximating Dirac delta functions with Gaussian kernels of finite variance around each particle. The numerical scheme is tested against the benchmark by studying the relaxation of non-equilibrium distribution functions following Brownian motion and the linearized Landau equation. For the relevant and finite range of time step sizes, the proposed convolution solution method provides a reasonable accuracy, yet the limit of its capability is investigated.
In order to achieve sustainable confinement in fusion plasmas, it is crucial to understand and mitigate all transport mechanisms. In recent years observations show that there is a strong evidence that the overall transport of heat and particles is to a large part caused by intermittency (or bursty events) related to coherent structures. In this work a novel approach where a global heat flux model based on a fractional derivative of plasma pressure is proposed for the heat transport in fusion plasmas. A simplified transport coefficient is assumed to capture all physical properties of the plasma, including neoclassical and anomalous transport. This assumption yields an exponent of the main derivative to be fractional instead of integer. The fractional degree of the heat transport is defined through the power balance analysis. In the proposed fractional model, a single constant fractionality index, α, is used as the dominant global scale dependence of the transport which is modified as compared to a diffusive model where α=2. Our aim with this work is to find a reduced (i.e., with the least number of parameters involved) transport model that can predict most plasmas, therefore ignoring the detail nature and the classifications of the transport processes involved, and bundle their average (time/radial) effect into one constant parameter, α. The method was used to study the heat transport in a selected set of JET plasmas, including C-Wall and ITER like Wall, L-mode, H-mode, with many different heating and fueling schemes from a wide range of experimental programs and plasmas with and without ELMs and various MHD modes active. The average fractional degree of the heat flux over the dataset was found as α ~ 0.8. Note that the model will instantly yield necessary profile evolutions and could thus be used as a feedback control of the plasma stability and control in real time by predicting profiles and providing a tool to detect and perhaps prevent or mitigate destructive transport events. It should be noted, that in some cases there is a wider range of α parameters over the database because these plasmas are in essence very different with one another on many factors such as the NBI or ICRH input power, fueling scheme, ELM control, etc. What we are observing however, is that a significant number of these plasmas fall into a similar range for α parameter specially as the input power is increased yielding a transport model with predictive power in a wide parameter regime. Finally, we would like to make a note that this study is the first of its kind and its findings are expected to encourage further discussion on the validity and the mathematical limitations of our current models to address global properties of transport in fusion plasmas.
The pedestal, and transport barriers in general, play an important role in tokamak performance and thus it is desirable to find a comprehensive model for these regions. In transport barriers, the applicability of standard neoclassical theory is limited because of sharp gradients of temperature, density, and radial electric field. We have developed a new neoclassical approach that sets the scale length in transport barriers to be the poloidal gyroradius. This ordering implies that the poloidal component of the E×B-drift becomes of the order of the thermal velocity. Using large aspect ratio and low collisionality expansions, we define a new set of variables based on conserved quantities, which simplifies the drift kinetic equation whilst keeping finite orbit width effects. Previous work, which accounted for strong gradients in density and electric field, is extended by allowing the temperature gradient to have the same scale length as the density gradient, and by including a poloidally varying part in the electric potential, which modifies the trapping condition for particles in the transport barrier. Studying contributions from both passing particles and trapped particles in the banana regime, we find that the resulting transport equations of particles, parallel momentum and energy are dominated by trapped particles. The poloidally varying part of the electric potential does not only contribute to transport but plays an important role in trapping particles as it grows across the transport barrier. We manage to reproduce profiles of density and temperature that are similar to those measured in real devices and find that it is necessary to choose a finite particle flux as well as a non-zero parallel momentum input to do so.
We provide a gyrofluid model of a collisionless and magnetized plasma, valid for finite βe, finite parallel magnetic perturbations and electron finite Larmor radius effects. This model is used to study the linear and non-linear evolution of magnetic reconection and magnetic islands. Gyrofluid models provide an effective tool, complementary to kinetic models, for studying such effects.
To understand plasma behaviour in the scrape-off layer (SOL), we need to know the boundary conditions for the plasma and electromagnetic fields near a divertor. At the boundary, in the direction perpendicular to the wall, there are four length scales of interest. These are the Debye length $\lambda_D$, the ion gyroradius $\rho_i$, the projection of the collisional mean free path in the direction normal to the wall $\lambda_{\perp}$ and the device size L. Assuming that the plasma near the divertor satisfies the scale separation $\lambda_D \ll \rho_i \ll \lambda_\perp \ll L$, we can split the plasma-wall boundary into three separate layers, the layer closest to the wall is the Debye sheath of width $\lambda_D$, then follows magnetic presheath of width $\rho_i$ and then the collisional layer of width $\lambda_{\perp}$. Plasma dynamics in the first two layers are well understood. In the SOL at distances much greater than $\lambda_{\perp}$ from the wall collisionality is high and Braginskii fluid equations are used to model the plasma behaviour, the ion and electron distribution functions are assumed to be approximately Maxwellian. The collisional layer connects this region of high collisionality with the collisionless magnetic presheath, where the ion distribution function is far from Maxwellian. The distribution function must satisfy the Chodura condition at the entrance to the magnetic presheath. We have also found that at the entrance of the collisional layer the flow of ions has to be supersonic. The numerical analysis should recover these results. To analyse the collisional layer we use the Galerkin method to solve the drift kinetic equation in one spatial dimension with the full Fokker-Planck collision operator, together with the quasineutrality equation and the assumption of adiabatic electrons. For our boundaries we assume that all ions that reach the wall are absorbed and we set the distribution function far away from the wall to be approximately Maxwellian.
The transport of particles and energy from the core plasma in a divertor tokamak to the reactor walls, via the scrape-off layer (SOL), occurs largely parallel to the magnetic field lines. Experimental evidence and theoretical considerations suggest that a fluid approach to modelling this transport may miss some important behaviour. In particular, temperatures at the target may be modified by nonlocal heat flow, and reaction rates with atomic and molecular species may be influenced by strongly non-Maxwellian electron distributions. Both of these effects are closely related to the onset of detachment. Here, an approach to kinetic modelling of parallel transport using the code SOL-KiT will be presented, which features fully kinetic electrons as well as a self-consistent fluid model for comparison. Extensions to the SOL-KiT model will also be discussed, including a quasi-2D fluid model for hydrogenic neutrals and allowing for independent ion and electron temperatures. Results will be presented of comparisons between fluid and kinetic simulations, showing modifications to temperature gradients, plasma-neutral reaction rates and electron-ion energy transfer.
Particle-in-cell (PIC) codes are one of the workhorses for numerically exploring plasma dynamics across a vast parameter space. While explicit discretizations of PIC systems allow for a straightforward time integration, stability requirements set strict limitations on both the maximal time-step and maximal grid spacing. Implicit PIC methods on the other hand put laxer restrictions on the time-step and the grid size, but require the solution to demanding linear systems for each time step. In particular, recent work shows that such linear systems are effectively solved using Jacobian-Free Newton-Krylov methods. In this contribution we demonstrate that machine learning can accelerate the GMRES solver, used in time-stepping of the implicit PIC method. To accelerate numerical simulations with machine learning one must determine where predictions from machine learning models can interface with numerical solvers. Also, machine learning models are known to interpolate well and extrapolate poorly. Therefore one must also explore how well a machine learning model can accelerate a simulation, depending on the training data set. To address these fundamental questions, we consider the simplest possible case, 1d electron plasma oscillations. We train neural networks to propose vectors that augment the Krylov space on which GMRES operates. We present an extensive hyperparameter scan over model architectures and explore the performance of multiple loss functions and regularizations for this task. Furthermore, we define a matrix of simulation initial conditions and evaluate model performance for disjoint training and validation sets. As a result, we find that even moderate augmentation GMRES’ Krylov space with proposal vectors from machine learning models can reduce the residuals by a factor of 10.
Drift-reduced plasma fluid models are commonly used to model the edge-SOL region of L-mode discharges in tokamaks. It is often observed that electrostatic simulations of plasma turbulence are restricted to a very small explicit timestep, or the implicit timestep is very poorly conditioned. The origins of this restriction can be traced to the properties of the linear waves supported by the model. We present the results of our recently submitted paper, where we review the impact on the system dispersion relation of different approximations that can be made in the model Ohm's law. The analysis uses a simplified set of equations which nevertheless contain the key linear features of a wide class of drift-reduced models used in the plasma community to model SOL turbulence. The electrostatic limit features dispersion relations that diverge at small perpendicular wave number $k_{\perp}$, giving rise to fast dynamics that may limit timestep size in simulations, or result in a poorly conditioned implicit timestep, as has been noted in gyrokinetic simulations. An electromagnetic Ohm's law removes the difficulty at small $k_{\perp}$, but in the zero electron mass case introduces a divergence at large $k_{\perp}$ due to the kinetic Alfven wave. Including in addition finite electron mass results in a well behaved dispersion relation, with parallel wave speeds between the Alfven speed at low $k_{\perp}$ and the electron thermal speed at high $k_{\perp}$. The practical significance of our conclusions from the dispersion relation is demonstrated in non-linear seeded filament simulations using STORM. Even in the very low \beta conditions of the SOL, including both electromagnetic effects and finite electron mass results in the fewest iterations per unit simulation time being taken by the implicit time-solver.
Numerical codes for electromagnetic wave propagation in fusion plasmas are mainly based on frequency-domain asymptotic methods, which provide a fast solution and are thus valuable for experiment design and control applications. However, in several cases of practical interest (like O-X-B mode conversion, mm-diagnostics) these tools run close to their limits of validity and should be compared to analytic or numerical full-wave solutions. The Fortran code RFFW (successor of FWTOR) solves Maxwell’s equations for the evolution of electromagnetic waves in plasmas. The propagation is computed with the FDTD method, whereas the medium response is described in terms of the generated plasma current. The plasma equilibrium can be analytic or experimental, including non-axisymmetric perturbations and density fluctuations, and the wave electric field spectrum can be arbitrary. In fusion-relevant applications, the code may conduct investigations of RF propagation and absorption in plasma geometry relevant to H&CD, reflectometry and MHD control. The results may provide a more detailed picture of the physics involved, and also benchmark the validity of the currently established tools in this field, which are based on simpler models. The code has been parallelized under EUROfusion HLST support, which has allowed the exploitation of the much larger CPU power and memory of supercomputers. The first step was to obtain a hybrid parallel (MPI + OpenMP) version while retaining the original code structure, with an overall speedup factor of 2500 as compared to the serial code version. Afterwards, the work was focused on MPI communication improvements like e.g., piggy-backing from multiple array exchanges in one and finer grain exchanges, further optimization of the memory usage from loop structures, and introduction of parallel I/O and restart functionality. All these improvements, apart from the simulation of cases where a small grid scale is required (e.g., EC absorption in high-density plasma), allow to upgrade the physics model of the code with the inclusion of processes that require themselves a very large computational burden, like the convolutional computation of the plasma current response in time domain based on the Fourier transform of the frequency-domain dielectric tensor. The latter may open the way for the exploration of EC physics problems for which asymptotic methods cannot be applied. In this work, first we present the main aspects of the code physics and technical setup, and also refer to details of the parallelization scheme. Then, we show results that exhibit the strong scaling of the code, as well as examine cases relevant to ECRH/ECCD application in medium-sized tokamaks (AUG, TCV).
Wave-particle interactions are ubiquitous in space and laboratory plasma systems and have been the subject of intense research interest for many decades. From a theoretical point of view, the nonlinear motion of a charged particle with an electrostatic wave has been one of the basic paradigms of complex and chaotic Hamiltonian dynamics. However, although single particle dynamics have been extensively studied, the collective dynamics of an ensemble of particles with a distribution of initial momenta, has not been systematically investigated. Analytical results are well-known for two limiting cases: Landau damping for particle interactions with low amplitude wavepackets under the quasilinear theory and particle interactions with a monochromatic wave having finite amplitude. In this work we systematically investigate the collective particle dynamics, in terms of momentum and energy transport, under interaction with finite amplitude electrostatic Localized Wavepackets of various spatial widths. Both resonant and ponderomotive effects are shown along with their crucial dependence on the parameters of the wavepackets and the initial momentum distribution functions. The collective particle dynamics are described in terms of the time evolution of the momentum distribution function as well as its first moments. Our numerical findings are systematically compared to analytical results for the aforementioned limiting cases and extend to more general configurations. Moreover, we consider an analytical kinetic description of reduced kinetic models [5] and discuss its paramount importance for higher-dimensional systems such as those describing particle motion in toroidal plasmas where numerical particle tracing is computationally expensive.
The assumption of a small relative deviation from background f_0 under the delta-f scheme often used to simulate the plasma core will not be valid when simulating the plasma edge, characterized by low density and temperature and strong gradients. In order to retain the noise reduction benefit of the delta-f scheme as compared to the full-f approach, a study of a transition scheme by means of a time-dependent temperature profile in the background Maxwellian is done. This profile is adapted by locally relaxing the kinetic energy accumulating in the deviation to the background component of the distribution function. The background distribution f_0 becomes time-dependent, with an adaptation rate that is specified as a free parameter. To this end, simulations of a simplified system mimicking the plasma edge are run using GK-engine, which is a delta-f PIC code solving the nonlinear electrostatic gyrokinetic equations in sheared slab geometry. Initial radial density and temperature profiles exhibiting high logarithmic gradients are used. Radially dependent sources in the form of a Krook operator are introduced to achieve quasi-steady state with sustained heat flux. All simulations are subsequently switched to flux-driven with specified heat sources. Convergence studies and ensemble runs are done against increasing number of markers and under different background temperature adaptation rates. Signal-to-noise diagnostics of initial results have shown an improvement as high as seven times the SNR value, when comparing the adaptive scheme against cases with a stationary background Maxwellian.
The development of the EU-DEMO reactor is at the pre-conceptual design phase. At this stage, close attention is paid to the heating mix necessary to fulfill all the plasma requirements: breakdown, ramp-up, L-H transition, burn control, NTM stabilization, sawteeth pacing, radiative instability control and ramp-down. Integrated modeling is an effective tool to compare the impact of dominant electron vs. ion heating on turbulence and plasma kinetic profile evolution. Thus, the ability of a given heating mix to meet the aforementioned requirements can be systematically studied. We have utilized the ASTRA transport code, coupled to the TGLF turbulent transport model for the core region of DEMO in order to compare the plasma response to dominant electron and ion heating mixes representative of ECRH, NBI and ICRF in L- and H-mode plasmas. Suitable boundary conditions are applied at the separatrix via a 0D 2-point model. Initial estimations point to the feasibility of the L-H transition and significant fusion power production even in purely electron-heated plasmas. The L-H transition work has been revisited and extended through the inclusion of Xe and W impurities, an expansion of the scanned parameter space in terms of density and ECRH power, a scan of the ECRH power deposition location and width as well as TGLF-predicted electron density profiles. The inclusion of impurities sets a strict concentration threshold for the L-H transition at high densities, while at low to intermediate Greenwald fraction the L-H transition becomes accessible at much higher impurity concentrations. Simulation in H-mode are performed by modeling the pedestal with scaling laws coupled to the core parameters. The fraction of electron heating and Xe concentration is scanned and we determine its impact on the obained fusion power and separatrix heat fluxes.
The MHD dynamo effect is an intrinsic and fundamental feature of reversed-field pinch (RFP) plasmas. It plays an important role in the tokamak as well (commonly as referred to as the "flux pumping" mechanism) in particular for the hybrid scenario with central safety factor close to one. In this contribution, we review results based on the above-mentioned nonlinear 3D MHD theory and simulations, and related experiments. Such results allow to identify the underlying physics of the MHD dynamo effect common to tokamak and RFP configurations: a helical core displacement modulates parallel current density along flux tubes, which requires a helical electrostatic potential to build up, giving rise to a helical MHD dynamo flow. Similarities between the MHD dynamo at play in the reversed-field pinch and tokamak configuration will be discussed, with the aim of providing a common theoretical framework for the two configurations. Both the quasi-periodic sawtoothing regime and in the stationary helical regime (obtained either with application of magnetic perturbations or at high plasma pressure will be considered as result of the nonlinear 3D MHD codes SpeCyl and PIXIE3D.
JT-60SA is a fully superconducting new tokamak device, designed, built and exploited jointly by Europe and Japan. It is the largest tokamak ever built before ITER and it is now in its commissioning phase. JT-60SA will exploit and extend the legacy both of JET and of the superconducting tokamaks presently in operation (WEST, EAST, KSTAR). It is expected to be at the forefront of the international fusion programme for many years, both before and during the D-T phase of ITER operation. The main missions of JT-60SA are: i) to support the ITER experimental programme as a satellite machine; ii) to pave the way to the next step of the international fusion programme, i.e., the demonstration fusion reactor (DEMO); iii) to promote the growth of a new generation of physicists and engineers who will exploit the next step fusion devices. In this tutorial, the main characteristics of JT-60SA and the highlights of its scientific programme will be presented, with particular focus on the opportunities for physics studies and validation of theory and models that the machine will offer.
We present results from a variety of dedicated L-H transition studies at JET-ILW, emphasizing the discrepancies between experimental data and accepted models of the transition. From earlier experiments in JET-ILW it is known that as plasma isotopic composition changes from deuterium, through varying deuterium/hydrogen concentrations, to pure hydrogen, the value of the density at which the threshold is minimum, ne,min, increases, leading us to expect ne,min(T) < ne,min(D). Preliminary analysis of the first JET-ILW Tritium L-H experiments, shows transient ohmic L-H transitions for ne< ne,min(D), as expected. At higher densities, with NBI heating, we see hints of PLH(T) < PLH(D). An analysis of Doppler reflectometer measurements of the radial electric field in D and He plasmas has been carried out. We do not find a critical radial electric field value or vExB rotation before the transition. Instead, it appears that the diamagnetic velocity, proportional to ∇p, may be a better indicator of the required conditions for an L-H transition. In H vs D it has been shown that the reason for the increased PLH in H is that lower confinement in H implies higher fuelling and power are required to match the edge pressure profiles before the transition. This also tells us that ∇p before the transition is important, and is a reminder that PLH is in fact determined by plasma transport characteristics in L-mode. Planned confinement studies in T may help elucidate this connection.
This presentation will focus on recent progress towards the validation of 3D non-linear MHD disruption simulations with the JOREK code. Simulations of a disruption triggered by an argon massive gas injection in JET pulse 85943 have been compared in detail to experimental data. Synthetic diagnostics have been used for the purpose, including interferometry, bolometry and saddle loops. A good global match has been found for these diagnostics as well as for the evolution of the plasma current, including the characteristic spike, giving confidence in the model. Analysis of the simulations reveals that the drive of the m=2/n=1 tearing mode via radiative cooling inside the island is a key mechanism in the disruption process. This drive leads the 2/1 island to grow to a very large size, causing a relaxation all the way to the plasma center. Other cases, including with shattered pellet injection, will also be discussed.
Shattered pellet injection (SPI) is the current concept for the ITER disruption mitigation system (DMS) to prevent disruption-related damage from thermal loads, electromagnetic forces or runaway electron (RE) beams. Compared with impurity or mixed deuterium-impurity pellets that contain large quantity of impurities like neon, pure deuterium (D2) SPI is expected to strongly dilute the plasma without immediately triggering a thermal quench and could be instrumental for RE avoidance in ITER. Detailed simulations of D2-SPI-induced disruptions, especially the (pre-)thermal quench phase in a JET discharge (#96874) with the JOREK code will be presented. In this discharge, a pure D2 pellet with a diameter of 12.5mm (barrel A) was shattered before entering the 3MA/2.9T H-mode plasma, exhibiting a much longer cooling time (~9ms between the arrival of the first shards at the plasma edge and the onset of thermal quench) than that with impurity pellets (typically 1~2ms). 3D non-linear MHD simulations show that the plasma is not fully cooled during the penetration and ablation of the D2 shards, allowing the current density (j) to evolve on a much slower time scale than the penetration of the shards and no evident global contraction of j has been observed. Detailed comparison between simulations and experimental measurements will be presented to validate the model as well as clarify the underlying physics. In particular, the role of background impurities, pointed out as a concern in the predictions for ITER, will be examined in detail. Possible effects of plasma rotation on (pre-)thermal quench dynamics will also be discussed.
The resistive tearing mode instability is a well-known phenomenon that is often linked to magnetic reconnection. When magnetic field lines reconnect, this leads to a conversion of magnetic energy into kinetic or thermal energy. In turn, this can result in interesting events such as solar flares in the solar corona or the disruption of plasma confinement in tokamaks. Therefore, understanding how this tearing instability is influenced by various physical effects can lead to new insights on how such events are triggered. Since the literature regularly claims that flow has a stabilizing effect on the resistive tearing mode, we use the new linear 1D magnetohydrodynamic (MHD) spectral code Legolas (https://legolas.science) to combine the effects of equilibrium flow and resistivity to examine the influence of the flow on the tearing mode parametrically. This allows us to investigate the influence of the flow speed, the wavelength, the angle between flow and wavevector, etc. in different magnetic field configurations. This parametric study serves as a first step in identifying the combined effects of equilibrium flow, magnetic shear and finite resistivity on the linear MHD spectrum.
The ideal-MHD theory of axisymmetric modes with toroidal mode number n = 0 in tokamak plasmas is developed. These modes are resonant at the magnetic X-points of the tokamak divertor separatrix. Consequently, current sheets form along the separatrix, which profoundly affect the stability of vertical plasma displacements. In particular, current sheets at the magnetic separatrix lead to stabilization of n = 0 modes, at least on the ideal-MHD time scale, adding an essential ingredient to the mechanism of passive feedback stabilization. The theory discussed here presents analogies with the physics of current sheet formation from the evolution of internal kink modes and the magnetic island coalescence problem. In these works, the ideal-MHD constraint causes magnetic flux to pile up near the X-points, leading to perturbed localized currents and a stabilizing effect in the ideal-MHD limit. For the island coalescence problem, it was found that a chain of magnetic islands becomes ideal-MHD unstable when the island width exceeds a critical threshold. In any case, flux pile-up prevents any further nonlinear evolution for unstable internal kinks and island coalescence unless the ideal-MHD constraint is relaxed, e.g., by resistivity. In our problem, vertical displacements are linearly stable in the ideal-MHD limit when the mode resonance at the equilibrium X-points of the divertor separatrix is appropriately taken into account. We remark that the resonant behavior of n=0 modes at divertor X-points is an ideal-MHD phenomenon. Therefore, it is reasonable that it must be treated first according to ideal-MHD. Future work will consider extended-MHD effects. We have found that, when the plasma density extends to the magnetic separatrix and n = 0 perturbations resonate at the magnetic X-points, vertical displacements are stable, at least on ideal-MHD time scales, without any need for passive stabilization elements. The stabilization mechanism is a direct consequence of the ideal-MHD flux-freezing constraint on the X-points, generating current sheets localized along the magnetic separatrix, exerting a force capable of pushing back the plasma in its vertical motion. This also suggests that plasma electrical resistivity in a thin boundary layer along the magnetic separatrix, in addition to wall resistivity, may have a profound impact on the stability of n = 0 vertical displacements.
There is considerable interest in studying plasma physics problems that will be solvable to error-correcting quantum computers. An interesting class of such problems is the propagation of electromagnetic waves based on Maxwell equations, with the plasma physics determining the dielectric properties of the medium. In order to develop the fundamental concepts for casting a classical wave propagation equation into a form suitable for quantum computers, we consider a medium that can be described by a scalar dielectric. We formulate an interleaved sequence of collide-stream operators acting on a lattice of qubits, in which the collision operators entangle the local on-site qubits, while the streaming operators spread this entanglement throughout the lattice. We study the interaction of a plane electromagnetic pulse (propagating in the x-direction) with a 2D cylindrical dielectric in the x-z plane. Initially, Ey and Bz are functions only of x, t. Following the interaction with the dielectric, the fields become 2D with a Bx(x,z,t) generated to preserve div B = 0. Multiple internal reflections create a series of circular wavefronts that emanate from the dielectric cylinder. Our mesoscopic QLA codes are ideally parallelized on classical supercomputers and can be readily encoded onto a quantum computer. For inhomogeneous media one finds both Hermitian and anti-Hermitian evolution operators: the anti-Hermitian operators require a doubling of the number of qubits/lattice site in order to yield a unitary representation. The QLA is easily extended to 3D by taking tensor products of the individual 1D QLAs.
The solar wind carries a broadband of fluctuations in density, velocity and magnetic fields that, at the large scales, have been interpreted in terms of an ongoing magnetohydrodynamic turbulent cascade. Alfvénic fluctuations have indeed been commonly observed in the solar wind since the first in-situ measurements, and they are thought to provide a possible mechanism to heat the solar corona at temperatures in excess of one million degrees and to accelerate the solar wind. Parker Solar Probe (PSP) was launched in August, 2018. It will be the first spacecraft to fly into the sun’s corona, to within about 10 solar radii from the sun’s surface, with the goal to understand what heats the corona and accelerates the solar wind. Early measurements from PSP have already provided a glimpse of the “young” solar wind in regions never explored before. Closer to the sun, the wind appears to be permeated by magnetic field lines which are strongly perturbed, to the point that they produce local inversions of the radial magnetic field, known as switchbacks. The corresponding signature of switchbacks in the velocity field is that of local enhancements in the radial speed (or jets) that display the typical velocity/magnetic field correlation that characterizes Alfvén waves propagating away from the sun. Switchbacks are thus an extreme case of the Alfvénic fluctuations that dominate the solar wind energy spectrum further away, and may be the remnant of coronal processes leading to solar wind formation –although their origin is still open to debate. After reviewing the main properties of Alfvénic fluctuations and switchbacks in the solar wind, we will address how their stability and evolution is affected by nonlinearities, kinetic effects and solar wind expansion. We will discuss what are the implications for models of switchback generation, and we will conclude by outlining remaining open issues.
The performance of magnetic confinement devices is determined mainly by turbulent transport inducing particle and energy losses, and limiting plasma confinement. Among the different experimental actuators of turbulence, supra-thermal particles – generated via external heating – are typically considered one of the most efficient in suppressing ion-temperature-gradient (ITG) driven turbulence. In this context, valuable steps into understanding the underlying physical mechanisms responsible for this turbulence regulation have been taken from first principle gyrokinetic simulations. In particular, recent results have shown that fast particles can resonate with the bulk ion-driven ITG micro-instabilities when their magnetic drift frequency is close to the ITG frequencies. When this resonant condition is fulfilled, supra-thermal particles can deplete the free energy content of the ITG instabilities if the fast particle temperature gradient exceeds the respective density gradient. Signatures of this potentially beneficial resonant interaction have been observed at ASDEX Upgrade, JET and are also expected for W7-X and the ramp-up phase of an ITER standard scenario. In this contribution, we present results of radially global GENE simulations showing that this wave-particle resonant mechanism might trigger the formation of a new type of transport barrier, called F-ATB (fast ion-induced anomalous transport barrier). These numerical findings guided the design of a recent ASDEX Upgrade experiment, maximizing the beneficial role of the wave-particle resonance interaction. Features of transport reduction and the formation of a central region of improved confinement were observed together with no degradation of the energy confinement during the external heating power ramp-up.
Solar flares are the most powerful explosions in the solar system, and are widely accepted to result from release of stored magnetic energy through magnetic reconnection. Many aspects of flare physics remain poorly understood: in particular, the large numbers of high-energy electrons and ions, forming a non-thermal element of the energy distribution, which can be the predominant energy release channel. Understanding the acceleration and transport of these non-thermal particles is a challenge for plasma physics theory, and a comprehensive model must account for both the global magnetic field structure on length-scales of around 10 – 100 Mm, and the plasma kinetic processes associated with energy dissipation and particle acceleration, on scales as small as 1 cm (electron gyro-radius). I will outline our current knowledge of energetic particles in flares based on Hard X-ray and radio observations, and the mechanisms which are proposed to explain their origin. Then I will focus on a scenario for confined solar flares, in which reconnection and energy release are triggered by the kink instability in twisted coronal loops, with distributed particle acceleration in fragmented current sheets. A forward-modelling approach predicts observable emissions associated with both thermal and non-thermal plasma in this scenario, and provides a means potentially to detect twisted magnetic fields in the corona. Building on this approach, it has recently been shown that an individual flaring event can be effectively modelled using a combination of 3D magnetohydrodynamic and kinetic (test particle) simulations, giving good agreement with observations. This allows modelling of the escape of energetic particles into the heliosphere, an important issue for space weather prediction. Many other stars also exhibit flares, and flares in young stars, known as T-Tauri stars, can be especially vigorous. I will describe a recent model for the radio and X-ray emission in T-Tauri flares, building on current understanding of solar flares, in which non-thermal particles generated by star-disk interaction populate large loops. Finally, I will discuss the outstanding unsolved problems in this area, and suggest some directions for future work, including some synergies with fusion plasmas.
Metallic PFC melting events in contemporary machines as well as future reactors fall under a rather unique regime; (i). Due to the limited wetted area, the liquid pools are surrounded by progressively colder solid surfaces so that once the melt is accelerated out of the pool, under the action of plasma-induced forces, it promptly solidifies. This necessitates modelling of the full coupling between the fluid dynamics and heat transfer, in order to capture the constantly evolving melt depth and the strength of the viscous damping correctly. (ii). The liquid pools are characterized by a vast scale separation, since they feature extents up to an order of a meter versus 100’s micron depths. The 3D description of the severely deformed free surface and moving solid-liquid boundary at extended spatiotemporal scales is computationally prohibitive. The use of the shallow water approximation, justified by (ii), offers a suitable compromise, reducing the dimensionality of fluid equations by one owing to depth-integration of the Navier-Stokes equations. Such averaging makes computations feasible, while retaining sufficient details for an adequate description of the macroscopic melt motion. The MEMOS-U model is based on applying the shallow water approximation to metallic melts induced by hot magnetized plasmas where phase transitions & electromagnetic responses are pivotal. The code solves the incompressible resistive thermoelectric MHD equations within the magnetostatic limit together with the convection-diffusion equation for the temperature. In order to formulate boundary conditions on the free surface, not only the liquid but also the plasma-plus-vapour ambient medium are implicitly treated as immiscible TEMHD fluids. Ultimately, the plasma-vapour effects appear exclusively through the pressure and drag force in the stress balance conditions as well as the incident plasma heat flux and the vapor cooling flux in the boundary condition for the temperature equation. Those plasma properties are quantities extracted from external models or experimental observations; these inputs and their underlying assumptions are critically analyzed. In addition, a rigorous description for the bulk current density (dependent on the plasma scenario and material properties) entering the free surface needs to be formulated. In the case of tungsten, owing to its high melting point and moderate work function, the only non-trivial boundary condition for the current continuity equation refers to the thermionic current density that escapes to the pre-sheath. Its estimation constitutes the most sensitive aspect of ‘plasma-effects’ modelling. This is due to the fact that the escaping thermionic current density simultaneously affects the heat balance (as an effective cooling channel) and the dynamics (through the dominant JxB acceleration) leading to a rather involved overall effect on the results. The MEMOS-U model employs results of dedicated PIC studies which revealed that the transition to the space-charge limited regime and the subsequent strict limitation of the escaping thermionic current density are global characteristics of strongly emitting sheaths also in the presence of inclined magnetic fields. The MEMOS-U model has been tested against multiple dedicated experiments in ASDEX-Upgrade and JET featuring different PFC materials, exposure geometries and plasma scenarios. The achieved quantitative agreement confirmed that MEMOS-U accurately describes the main physical mechanisms responsible for shallow melting and the final surface deformation profiles, the latter being the principal quantity of interest in predictive simulations of wall damage in future fusion reactors.
The suppression of turbulence in fusion plasmas, crucial to the success of next-generation tokamaks such as ITER, depends on a variety of physical mechanisms including the shearing of turbulent eddies via zonal flow and possibly the generation of intrinsic rotation. The turbulence exhibits interesting features such as avalanche structures and self-organisation, and its absence is associated with the formation of internal transport barriers. In order to successfully capture all of these effects in gyrokinetic simulation, one may need to allow for the inclusion of global effects (such as radial profile variation), as well as other often-neglected effects that are small in rho. A careful numerical treatment is necessary to ensure that both the global and local physics are calculated accurately at reasonable expense. To that end, we develop a novel approach to gyrokinetics where multiple flux-tube simulations are coupled together in a way that consistently incorporates global profile variation while allowing the use of Fourier basis functions, thus retaining spectral accuracy. By doing so, the need for Dirichlet boundary conditions typically employed in global gyrokinetic simulation, where fluctuations are nullified at the simulation boundaries, is obviated. This results in a smooth convergence to the local periodic limit as rho -> 0. In addition, our scale-separated approach allows the use of transport-averaged sources and sinks, offering a more physically motivated alternative to the standard sources based on Krook-type operators. Having implemented this approach in the flux-tube code stella, we study the role of transport barriers and avalanche formation in the transition region between the quiescent core and the turbulent pedestal, as well as the efficacy of intrinsic momentum generation by radial profile variation. Finally, we show that near-marginal plasmas can exhibit a radially localized Dimits shift, where strong coherent zonal flows give way to flows which are more turbulent and smaller scale.
A simplified local model of a tokamak plasma is derived in the low-beta limit of gyrokinetics in a slab of constant magnetic field curvature and gradient. The ordering adopted was chosen in order to retain Alfvénic perturbations to the magnetic field, while ordering out compressive perturbations, in a similar manner to previous work. In the electromagnetic regime, we demonstrate the existence of the novel “Thermo-Aflvénic instability” that arises due to a deviation from isothermality of the total temperature along the (exact) perturbed field line. This instability both destabilises kinetic Alfvén waves and enhances the conventional curvature-mediated ETG instability, driving turbulence on scales above the electron skin depth. Assuming critical balance, it is shown that the resultant turbulent heat flux is larger than that due to the electrostatic ETG modes, presenting a significant departure from the expected picture of the electron turbulent heat transport.
This communication introduces the concept of plasma-wall self-organization (PWSO) in magnetic fusion. The basic idea is the existence of a time delay in the feedback loop relating radiation and impurity production on divertor plates. Both a zero and a one-dimensional description of PWSO are provided. They lead to an iterative equation whose equilibrium fixed point is unstable above some threshold. This threshold corresponds to a radiative density limit, which can be reached for a ratio of total radiated power to total input power as low as 1/2. When detachment develops and physical sputtering dominates, this limit is progressively pushed to very high values if the radiation of non-plate impurities stays low. Therefore, PWSO comes with two basins for this organization: the usual one with a density limit, and a new one with density freedom, in particular for machines using high-Z materials. Two basins of attraction of PWSO are shown to exist for the tokamak during start-up, with a high density one leading to this freedom. This basin might be reached by a proper tailoring of ECRH assisted ohmic start-up in present middle-size tokamaks, mimicking present stellarator start-up. In view of the impressive tokamak DEMO wall load challenge, it is worth considering and checking this possibility, which comes with that of more margins for ITER and of smaller reactors. Important facts and results:The L-mode density limit increases with rising heating power P like P^{0.4} It is shown with a data base of 5 tokamaks that a scaling like (IP/a^4)^{4/9}, with a the small radius and I the total current, organizes much better the data than Greenwald's one in I/a^2; especially device per device. This scaling is a part of those derived in earlier works, which are in much better agreement with the tokamak and RFP databases. They apply to the stellarator too. They describe the density limit as a radiative one, and therefore include naturally a clear explicit dependence on P. A large part of the radiative density limit comes from the radiation of impurities. Their amount is governed by plasma-wall interaction. We provide a self-consistent description of this interaction and introduce the concept of PWSO. Both zero and one-dimensional descriptions lead to a delay equation whose simplest expression is R+ = \alfa (P-R), where P is the total input power in the plasma, R is the total radiated power, and R+ is its delayed value. This makes the plasma-wall system unstable for \alfa >1. Since \alfa is proportional to the density below detachment, this threshold defines a density limit. It can be reached for a ratio of total radiated power to total input power as low as 1/2. When detachment develops, the plasma temperature at the plates decreases, which makes \alfa to vanish. This pushes the radiative density limit to very high values when physical sputtering dominates, in particular for tungsten. Hence density freedom. The 0D and 1D models of PWSO apply to the stellarator and to the reversed field pinch as well.
We present the results of 3D, flux-driven, global, two-fluid electrostatic turbulence simulations in a 5-field period stellarator with an island divertor. The numerical simulations are carried out with the GBS code, which solves the two-fluid drift-reduced Braginskii equations and has been extended recently to simulate plasma turbulence in non-axisymmetric magnetic equilibria. The vacuum magnetic field used in the simulations is carefully constructed using Dommaschk potentials in order to describe a configuration with a central region of nested flux surfaces, surrounded by a chain of magnetic islands. In a similar way to the diverted configurations of W7-AS and W7-X, particles and heat, transported radially outwards from the core region, reach the island region, which effectively acts as a scrape-off-layer with the open field lines striking the walls at specific toroidal locations of the device wall. We find that the radial particle and heat transport is mainly driven by a field-aligned mode with low poloidal wavenumber, whose origin is investigated theoretically. The equilibrium radial electric field in the core is found to be in the ion-root regime, Er<0. Transport is observed to be larger on the high-field-side of the device, where the amplitude of fluctuations is larger, despite the cross-phase between density and potential being smaller. A very good agreement is obtained when comparing the radial ExB flux of the simulation with the one predicted theoretically due to a single (dominant) coherent mode. In contrast to tokamak simulations and experiments, we do not observe radial propagation of coherent filamentary structures (blobs) that usually contribute to intermittent transport events in axisymmetric configurations, thus shedding light on the surprising differences between transport mechanisms in stellarator and tokamak configurations.
Understanding plasma-wall interactions produced by the huge particle and heat fluxes reaching the vessel walls in nuclear fusion devices is of uttermost importance for the next generation of reactors. For example, the divertor of ITER is expected to withstand heat loads of around 10MW m-2 in steady state operation. To recreate these conditions, the linear plasma device Magnum-PSI is currently operated at DIFFER. It allows exposing a target material to both steady state or pulsed plasma relevant for the ones expected in ITER. However, simulations are indispensable to understand the Atomic and Molecular (A&M) processes that take place in Magnum-PSI, and to extrapolate these observations to tokamak divertors. In this regard, SOLPS-ITER and B2.5-Eunomia are currently being employed to model Magnum-PSI in a wide range of parameters such as plasma current, magnetic field and plasma density. Both packages include the same CFD plasma code (B2.5) but different Mote-Carlo modules for A&M neutrals: Eirene and Eunomia, respectively. Eunomia was designed specifically for linear plasma devices like Magnum-PSI. Therefore, an effort is currently undertaken to port some of its functionality into Eirene, primarily designed for tokamaks. Although both neutral modules have similar capabilities, the implementation of some collision process lead to relevant differences. Due to the distinct multi-step implementation of Molecular Assisted Recombination (MAR) or the different collision rates used for electron impact ionization/excitation for H atoms, disparate sources of energy and particles for B2.5 are calculated by the Monte-Carlo codes. Elastic collisions of ions with molecules, although reading the same collision rate, are incorporated into the neutral codes with a distinct formulation for the an-isotropic post-collision angle, which result in different neutral distributions. B2.5-Eunomia can simulate Magnum-PSI for a wide variety of plasma scenarios. However, the lack of measured data for the potential profile at the source demands additional assumptions leading to uncertainties in the comparison. Moreover, the appropriate transport coefficients for plasma required to model Magnum-PSI are still under discussion, due to the differing from tokamaks range of plasma temperature and density in which the machine can operate. The present work is focused on code-code comparison aimed to reproduce Magnum-PSI Thomson Scattering measurements for radial electron density and temperature profiles few centimeters in front of the target. A special focus is put on how the particular implementation of different A&M processes in the two neutral modules lead to different coupled results.
Experimental results in the first campaigns of Wendelstein 7-X (W7-X) have shown that, due to the optimization of the magnetic configuration with respect to neoclassical transport, turbulence is essential to understand and predict the total particle and energy fluxes. This has motivated much work on gyrokinetic modelling in order to interpret the already available experimental results and to prepare the next experimental campaigns. Thus, it is desirable to have a sufficiently complete, documented and well-verified set of linear and nonlinear gyrokinetic simulations in W7-X geometry against which new codes or upgrades of existing codes can be tested and benchmarked. This work is an attempt to provide such a set of simulations through a comprehensive benchmark between the recently developed code stella and the well-established code GENE in W7-X geometry. It consists of linear and nonlinear collisionless electrostatic flux-tube simulations, organized into five different ‘tests’. They include stability analyses of linear ITGs and density-gradient-driven TEMs, computation of the collisionless relaxation of zonal potential perturbations and calculation of ITG-driven heat fluxes. As different magnetic field lines are not equivalent in stellarator geometry, simulations in two different flux tubes are provided, clarifying the similarities and differences between the stability features in both of them.
The excitation of toroidicity induced Alfven eigenmodes (TAEs) using an electromagnetic antenna acting on a confined toroidal plasma is studied. The antenna is described by an electrostatic potential resembling the target TAE mode structure along with its corresponding parallel electromagnetic potential computed from Ohm's law. Stable long-time linear simulations are achieved by integrating the antenna within the framework of a mixed representation and pullback scheme. By decomposing the plasma electromagnetic potential into symplectic and Hamiltonian parts and using Ohm's law, the destabilizing contribution of the potential gradient parallel to the magnetic field is canceled in the equations of motion. Besides evaluating frequency as well as growth/damping rates of excited modes compared to referenced TAEs, we study the interaction of antenna-driven modes with fast particles and indicate their margins of instability. Furthermore, we show preliminary results for nonlinear simulations in the presence of a TAE mode excited with the antenna.
Axisymmetric modes (i.e. with toroidal mode number n=0) destabilized by fast ions have been observed in recent JET experiments. Motivated by these experimental results, we have reconsidered the dispersion relation of macroscopic n=0 vertical displacements in shaped tokamak plasmas. Vertical displacements are normally stable thanks to a combination of passive and active feedback stabilization. Passive feedback stabilization, consisting of currents flowing along the wall and/or in plasma facing components, is required in order to stabilize the ideal instability. Considering the presence of a nearby resistive wall, the n=0 mode dispersion relation is cubic. Under relevant conditions two roots are oscillatory and weakly damped with $\omega \sim e_0^{1/2}\omega_A$, where $\omega_A$ is the poloidal Alfvén frequency and the ellipticity $e_0 = (b^2-a^2)/(b^2+a^2)$, with $b$ and $a$ major and minor semi-axis of the ellipse describing the plasma boundary. Since $e_0<1$, the mode frequency lies below the minimum of the Alfvén continuum, therefore the mode is not affected by continuum damping. The third root is unstable with growth rate related to the wall resistivity. Much work has been done regarding the third unstable root, which can be stabilized with active feedback. However, not enough attention has been given to the two oscillatory roots. In fact, even though these modes are damped by the wall resistivity, due to their oscillatory character they can be driven unstable by resonant interaction with energetic particles. In this work a fully analytic derivation of the cubic dispersion relation for the n=0 modes forms within the so-called reduced ideal MHD model is presented. Thereafter, effects of energetic particles are considered on the basis of the hybrid kinetic-MHD model, where thermal plasma is treated using ideal-MHD, while fast particles are described in terms of the collisionless drift-kinetic equation. The mode-particle resonant condition for these modes is $\omega =p\omega_{b/t}$, where $p$ is an integer number labelling harmonics over particle orbit periodicity and $\omega_{b/t}$ is the bounce (or transit) frequency of magnetically confined fast particles. In order to drive the mode unstable a fast particles distribution function with $\partial F_h/\partial E > 0$ is necessary. It is possible to show that such a distribution function can be obtained considering fast particles losses or a modulation of their source. Both trapped and passing particles can resonate. In particular for typical values of for JET and fast particles energy of the order of the MeV, the largest contributions to the mode resonance involve the first harmonic ($p=\pm 1$) for transit particles and the second harmonic ($p=\pm 2$) for trapped particles. The growth rate introduced by the resonant interaction can overcome the damping introduced by wall resistivity. This theory presents a possible explanation to the observed n=0 modes in presence of fast particles, alternative to the one proposed in previous work, which involves global Alfvén eigenmodes (GAEs).
While the achievement of high confinement regime is ensured through the formation and sustainment of edge transport barriers associated with sheared flows, the effect of non-axisymmetric perturbations of the magnetic field, like non-resonant magnetic perturbations or ripple, on the transition remains an open issue. The underlying loss of axisymmetry is responsible for a toroidal torque, sometimes called Neoclassical Toroidal Viscosity (NTV), which, as predicted by neoclassical theory, leads to magnetic braking. This constraint lifts the degeneracy between the mean toroidal velocity and the radial electric field which results from toroidal symmetry in the standard neoclassical theory. As a result, predictions on meaningful quantities regarding transport barriers, i.e. mean flows and radial electric field, are available in the low collisionality limit. The implementation of non-axisymmetric perturbations in the gyrokinetic code GYSELA shows quantitative agreement between neoclassical theory and simulations for experimentally relevant magnetic perturbation magnitude and collisionalities. In realistic plasmas, mean flows result from a competition/synergy between neoclassical effects and turbulence. In ion temperature gradient driven turbulence, this competition reveals that even a small perturbation amplitude, typically a few percents of the unperturbed magnetic field magnitude, tends to change the plasma mean toroidal rotation toward the counter-current direction. This effect has already been observed experimentally in JET, JT-60U and Tore Supra. The radial electric field increases by up to +80%, thus showing that a non-axisymmetric perturbation may impact the transition toward transport barriers. To understand and quantify this competition, the turbulent Reynold’s stress is compared with the theoretical magnetic drag. Simulations with realistic plasma parameters show that turbulent transport prevails for a ripple amplitude below a critical value. Above this threshold, sheared flows are controlled by the NTV.
An improvement of the boundary conditions scheme of the 3D nonlinear MHD numerical code SpeCyl is presented. Boundary conditions have been shown to play a key role in the helical self-organization both in Reversed Field Pinch and tokamak plasmas. Two different sets of boundary conditions have been extensively tested against ubiquitous relaxation phenomena induced by plasma current in toroidal devices: ideal kinks and tearing modes. The role of wall position and resistivity on linear perturbations profiles and their exponential growth rates was tested, motivating the need for a reformulation of fluid boundary conditions as well. Preliminary results of such new boundary conditions are also presented, along with a summary of the relevant theoretical framework underlying linear MHD instabilities.
In this work, an implicit scheme for electromagnetic particle-in-cell/Fourier simulations is developed using the v_∥ formula and applied to studies of Alfvén waves in one dimension and in tokamak plasmas on structured meshes. While the “particle enslavement” scheme has been introduced for reducing the degree of freedom of particles in the field-particle system, in this work, we focus on the theoretical analyses of the convergence of the system. An analytical treatment is introduced to achieve efficient convergence of the iterative solution of the implicit field-particle system. The essence of this scheme is to represent the particle moments such as the density and the parallel current as functions of field variables. Then the correction matrix of the implicit field-particle system can be obtained analytically. As a result, this treatment is termed with “moment enslavement”. Its application to the one-dimensional uniform plasma demonstrates the applicability in a broad range of β/m_e values. The simulation results and the theoretical results agree well in our studied regime where β_e m_i/m_e∈[1/16,32], k_⊥ ρ_ti=0.2, where β_e is the electron beta, k_⊥ is the perpendicular (to B) wave vector, ρ_ti is the Larmor radius of thermal ions, m_i and m_e are ion and electron masses. The toroidicity induced Alfvén eigenmode (TAE) is simulated using the widely studied case defined by the ITPA Energetic particle (EP) Topical Group. The real frequency and the growth (or damping) rate of the TAE with (or without) EPs agree with previous results reasonably well. The full f electromagnetic particle scheme established in this work provides a natural choice for EP transport studies where large profile variation and arbitrary distribution need to be captured in kinetic simulations. By combining recent developments with the applications to mode structure symmetry breaking studies, the ongoing work focuses on the simulations of kinetic ballooning mode simulations and EP driven Alfvén modes in realistic tokamak geometry, with the capability of including the plasma separatrix. The full f and delta f simulations are compared, demonstrating their features in performance and the capabilities of the full f scheme in treating specific problems.
Filaments in the boundary of magnetically confined fusion plasmas lead to enhanced erosion of the main chamber walls. These high-density, coherent structures can be thought of as intermittent fluctuations described by the Filtered Poisson Process (FPP) as a superposition of pulses with a fixed shape and a constant duration. Additionally, these fluctuations have large amplitudes compared to the background level which can arise in far-from equilibrium and turbulent systems as well as due to sudden, large-amplitude forcing. Conditional averaging has been a much-used tool for finding amplitudes of intermittent events, the waiting time between them as well as the average waveform. However, conditional averaging requires both significant amplitude thresholding and a large minimal distance between events, significantly limiting the number of found events. This reduces the accuracy of reconstructing intermittent data time series and the statistical analyses that can be derived in order to monitor and predict filaments. In this contribution, we study a variant of the Richardson-Lucy deconvolution algorithm. This is an iterative method converging to a least squares solution which can be used to recover event amplitudes and arrival times from an intermittent time series for a known typical event shape. The method was applied to synthetically generated data time series consisting of a superposition of one-sided exponential pulses. Signal reconstruction and recovery of event amplitudes and arrivals is excellent for low to moderate overlap between events. As event overlap increases, signal recovery remains excellent, but an empirical threshold for reconstruction of amplitudes and arrival times is found. The sampling time must be 10 times lower than the average waiting time or less. In the presence of noise with the same correlation function as the signal, spurious events are observed, and some thresholding or filtering must be used to separate the noise from the data time series. The deconvolution method requires event thresholding related to the noise to signal ratio and the degree of event overlap. Events separated by as little as two sampling times may be distinguished. Lastly, the deconvolution algorithm will be applied to MAST data from SOL probes in order to recover filament amplitudes and arrival times. We will compare and discuss these results to those extracted from the conditional averaging technique.
In a tokamak, radio frequency (RF) electromagnetic waves that propagate through low density plasma (n_l) enter the strongly turbulent edge region (n_e) before passing into the fusion plasma (n_p). Whether used for diagnostics or for heating and current drive, it is important to quantify the spectral properties of these waves. The magnetized n_l, n_p and n_e (via homogenization) regions are defined through the cold plasma, anisotropic permittivity tensor. Experimental evidence suggest that drift waves and rippling modes are present in the n_e region. Thus it is assumed that the n_e region is separated from the n_l and n_p regions by periodic density interfaces (plasma gratings) formed as a superposition of spatial modes with varying periodicity and random amplitudes. The ScaRF full-wave, 3D electromagnetic code has been developed for analyzing scattering scenarios of this form. ScaRF can be used for scattering analysis of any cold plasma RF wave and consequently for the scattering of electron cyclotron waves in ITER-type and medium-sized tokamaks. Since the density interfaces are random, the power reflection coefficient (R), obtained by ScaRF, is a random variable and is calculated for different realizations of the density interface. In this work, the uncertainty of R is rigorously quantified by use of the Polynomial Chaos Expansion method using Sparse Basis (SB-PCE) and Hyperbolic truncation schemes. The SB-PCE method is proven accurate, faster that other methods, and much more efficient, requiring only 11 R samples compared to thousands used by reference methods such as the Monte Carlo (MC) approach.
Three-dimensional force-free states describing toroidal plasmas with D-shaped cross section, are constructed. The construction is carried out by perturbing two-dimensional axisymmetric single-Beltrami states with translationally symmetric ones. The perturbation and the unperturbed magnetic field have a common Beltrami parameter λ, thus their superposition still satisfies the Beltrami equation. The boundary was imposed on the axisymmetric state upon using proper conditions for specific boundary points according to the shaping method. The addition of the translationally symmetric component as a small perturbation has a noticeable impact on the equilibrium state, i.e. one can observe helical magnetic islands in Poincaré maps, which appear in certain rational magnetic surfaces. Furthermore, the conjecture of previous work is confirmed, according to which the surfaces of the resulting 3D configuration remain closed (toroidal) in the vicinity of the magnetic axis if the axisymmetric field has sufficiently high weight in the superposition.
H-mode plasma scenarios, with enhanced energy and particle confinement at plasma edge, are a viable option for fusion energy production and represent the core operational regimes for present day and next generation fusion experiments such as ITER. Good confinement however comes at the price of potentially large pressure and current density gradients in the edge region, leading to the so-called Edge Localized Modes. Type-I ELMs in particular are large bursts that can damage the plasma facing components causing large heat and particle fluxes. Applying 3D resonant magnetic perturbations (RMPs) with non-axisymmetric saddle coils is a promising method to mitigate or suppress type-I ELMs. The DTT experiment, presently under realization, will be built in Frascati (Italy) with the main mission of developing reactor-relevant power exhaust technologies. ELM control is therefore of particular importance for DTT operational scenarios. Non-axisymmetric in-vessel coils are being considered for the purpose of ELM mitigation and suppression. Linear plasma response modelling is exploited in this work to assess the effect of different coil geometries on ELM stability in a full power DTT scenario. Peeling-like plasma response in particular is found to be correlated with ELM control and is computed by solving single fluid MHD equations in toroidal geometry with the MARS-F code, which includes the effect of plasma flow. The effects on the magnetic field spatial spectrum due to the three-dimensional geometry of the non-axisymmetric in-vessel coils are evaluated with the CARIDDI code. In order to evaluate the RMP effect on ELM stability, local plasma displacement in the x-point region is used as the main metric. This criterion has been correlated with ELM mitigation thresholds in MAST and ASDEX-Upgrade. By taking as reference a displacement value of 3 mm at the x-point, the optimal operational space for non-axisymmetric coils in DTT is sketched in terms of coil current and phasing between independent toroidal arrays. Configurations with either two or three sets of 9 coils each are used in the model. While n=1 RMPs induce an important core response in the reference scenario, n=2,3 response shows a main peeling-like component. The effect of edge q-profile variations and comparison with other metrics is considered, such as with the Chirikov parameter criterion correlated with the RMP-induced edge stochasticity.
It has been experimentally observed in both tokamaks and stellarators that peaked density profiles lead to enhanced confinement regimes. The reduction in transport is believed to be related to the stabilization of ion-scale turbulence. In this conference contribution, we perform gyrokinetic simulations with the gyrokinetic code stella focusing on the effect of the density gradient on nonlinear heat fluxes. The influence of the magnetic geometry is investigated by means of an inter-machine study that includes the W7-X, LHD, TJ-II and NCSX stellarators, as well as the Cyclone Base Case tokamak. The
simulations are collisionless and a vanishing electron temperature gradient is assumed. For the devices listed above, we have computed the ion heat flux as a function of the normalized density gradient, a/L n , for a fixed value of the normalized ion temperature gradient, a/L Ti =3, by means of nonlinear stella simulations with kinetic electrons. We show that, in a broad range of the scanned a/L n values, W7-X and NCSX exhibit a strong reduction of the ion heat flux with increasing a/L n . In TJ-II the reduction is more modest and in LHD the ion heat flux has a weak dependence on a/L n . In contrast to the stellarators, we have found that the ion heat flux of the tokamak increases strongly with the density gradient. By comparing the results, it is clear that, in stellarator geometry, the behavior of the linear growth rates as a function of a/L n does not correlate with the behavior of the ion heat flux. Finally, in this conference contribution we will also discuss the effect on the ion heat flux of treating the electrons adiabatically or kinetically, as well as the effect of taking zero or finite a/L Ti.
Verification and validation of the internal kink instability in tokamak have been performed for both gyrokinetic (GTC) and kinetic-MHD codes (GAM-solver, M3D-C1-K, NOVA, XTOR-K). Using realistic magnetic geometry and plasma profiles from the same equilibrium reconstruction of the DIII-D shot #141216, these codes exhibit excellent agreements for the growth rate and mode structure of the n=1 internal kink mode in ideal MHD simulations by suppressing all kinetic effects. The simulated radial mode structure agrees quantitatively with the electron cyclotron emission measurement after adjusting, within the experimental uncertainty, the q=1 flux-surface location in the equilibrium reconstruction. Equilibrium plasma pressure gradient and compressible magnetic perturbation strongly destabilize the kink, while poloidal variations of the equilibrium current density stabilize the kink. Furthermore, kinetic effects of thermal ions are found to decrease the kink growth rate in kinetic-MHD simulations, but increase the kink growth rate in gyrokinetic simulations, possibly due to the additional drive of the ion temperature gradient and parallel electric field. Kinetic thermal electrons are found to have negligible effects on the internal kink instability.
In the Wendelstein 7-X stellarator, the vacuum rotational transform, iota, has a flat radial profile and does not cross any major rational resonance. Nevertheless, during plasma operation the iota profile can be strongly modified by electron cyclotron current drive (ECCD) in such a way that the resulting iota profile passes through low-order rational values, and this can trigger magnetohydrodynamic (MHD) events. Indeed, W7-X plasmas are sometimes subject to repetitive collapses of core confinement, which can be observed regardless of the direction in which the ECCD current is driven. These phenomena are periodic, rapid (on the order of Alfvénic times) and large “crashes” of electron temperature, which are reminiscent of tokamak “sawtooth” instabilities. There are also examples of discharges where related events lead to complete termination of the entire plasma. Even though the origin of these MHD instabilities is not yet clear, the fast crashes are likely to be connected to the formation of magnetic islands and magnetic reconnection. In the present work, we discuss the MHD events happening during the course of sawtooth cycles in W7-X using the combination of slow current diffusion with fast relaxation that conserves the corresponding helical flux. Utilising Taylor relaxation theory, we predict the nonlinear redistribution of plasma current caused by the largest of the observed events and obtain a 3D post-crash state. To study different types of crashes — not only fast and large ones, where we demonstrate a good agreement with Taylor relaxed state, but also smaller ones — we combine the Taylor theory approach and resistive evolution aspects of the current profile development.
Radio frequency waves (RF) are scattered by filamentary structures which exist in the edge region of a tokamak plasma. The waves are reflected, refracted, and diffracted leading to a change in their spectral properties. The spatial profile of the launched power gets fragmented and part of the launched power can be coupled to an unwanted cold plasma wave. Consequently, the efficiency of heating and current drive by RF waves in the core plasma can be negatively impacted. The scattering depends on the frequency and the wavelength of the launched RF wave, the size of the turbulent structure, the density of the ambient plasma, and the density inside the filament. We will present results on studies of wave scattering by RF waves in three different frequency regimes – lower hybrid waves, helicon waves, and ion cyclotron waves. While the filamentary structures scatter RF waves, they can also be affected by the radiation pressure of the waves. The filamentary structures can be pulled towards the RF source or pushed away depending on a variety of plasma and wave parameters. We will present results on RF scattering and on the RF induced radiation pressure in the three different frequency regimes.
Understanding the field excitation and power losses in the plasma edge from an ICRH antenna is of paramount importance to avoid impurity release from the edge in reactor conditions. The semi-analytical code ANTITER IV provides a complete description in plane geometry (z along the total B0 field in front of the antenna, x the radial component) in the cold plasma approximation and with appropriate boundary conditions of the excitation of Fast (F) and Slow (S) waves in the inhomogeneous plasma in front of the antenna. For usual non-inverted heating scenarios (i.e. ω>ωci of majority ions) the Lower Hybrid (LH) resonance is present in the SOL region of the plasma density profiles considered in this study, representative for large machines like ITER or DEMO. The system of 4 first order differential equations describing the waves is singular at the position of the LH and integration through the LH layer is performed by analytical continuation. The S waves are excited by their confluence with the F waves or directly by the Ez component radiated by the antenna. The case of the 24 straps antenna array of ITER is modeled using reference profiles for low and high coupling conditions. An antenna based on this design is presently considered for DEMO. The analysis is made for an antenna and Faraday screen that are non field-aligned with respect to B0 and compared to the field-aligned cases. The results show that, for the considered density profile in front of the antenna: (i) A surface is wave excited along the LH layer together with TEM/z coaxial modes propagating between the wall and the LH region. These are decaying waves at both sides of the antenna but are of significant amplitude far from the antenna. (ii) The radial electric field component Ex is dominant in front of the antenna and along the wall for all antenna phasing cases. The Ez component is much smaller even if the Faraday screen is not aligned with B0. (iii) Ex is resonant at the LH resonance and is the only RF electric field component entering into the power loss expression at the LH. (iv) The minimization of the low |kz| components in the antenna excitation spectrum, mainly in the region |kz|=k0, leads to a reduction in the Ex amplitude in front and on the edges of the antenna, the amplitude of the coaxial and surface waves and the edge power loss. It leads also to a reduction of the z component of the surface current density Jz=-Hy at both sides of the antenna array. (v) The large Ex component extending outside the antenna region and diverging at the lower hybrid layer could be a main cause for the impurity production by an ICRH antenna. (vi) When the density profile is amputated to avoid the LH density the edge wave behavior is strongly modified, the edge power losses cancelled and the coupling results to the plasma center can be perturbed. Such a modeling simplification can lead to wrong conclusions.
In the pedestal, large bootstrap and Pfirsch-Schluter currents, arising from the steep pressure gradient, drive kink instabilities and enhance the drive for the peeling-ballooning modes. The latter triggers periodic plasma eruptions at the edge of a tokamak that significantly influence the pedestal properties, limiting the pedestal pressure gradient and degrading core confinement. The conventional gyrokinetic theory usually employed to describe the physics of the pedestal is not sufficient to recover the full spectrum of ideal MHD results, e.g. it cannot capture the physics of kink modes. Indeed, it adopts a Maxwellian equilibrium and therefore does not retain the current density gradient along the magnetic field lines. These current density gradient effects can typically be ignored in a core plasma of a conventional tokamak, but must be retained in the core of a spherical tokamak plasma or in the vicinity of the pedestal of a conventional tokamak where the pressure gradient is the strongest. We have developed a new electromagnetic gyrokinetic theory that captures the strong gradients in the pedestal. The particle distribution function is written in the form: $f=F + \delta f$, where $F$ describes the equilibrium piece (in the gyrokinetic ordering), i.e. it varies on transport time scales, and $\delta f$ is the fluctuating piece (the particle species indices are dropped for simplicity). To recover ideal MHD, we must retain neoclassical effects in the equilibrium distribution function and thus solve perturbatively the drift kinetic equation that describes $F$, expanding it in powers of $\delta = \rho_{\vartheta}/L \ll 1$, where $\rho_{\vartheta}$ is the poloidal Larmor radius and $L$ is the characteristic size of the system. Here the leading order corresponds to the Maxwellian, employed in conventional gyrokinetics. The following order introduces the neoclassical flows (e.g. bootstrap current flows), required to recover kink physics. The gyro-angle dependent piece of the second order correction in a $\delta$ expansion must be retained as well to ensure a consistent ordering. This $F$ is then employed to obtain the extended gyrokinetic equation for fluctuations. To solve the gyrokinetic equation for $\delta f$, we employ an expansion in $\Delta = \rho/L$, where $\rho$ is the particle Larmor radius. We provide two sets of equations: (a) global, integro-differential gyrokinetic equations that allow one to capture the global properties of modes and (b) local gyrokinetic equations which employ an eikonal ansatz for the fluctuating electromagnetic field. The latter provides the theoretical formalism in a form ready for implementing in existing local gyrokinetic codes.
The Debye sheath is shown to vanish completely in magnetised plasmas for a sufficiently small electron gyroradius and small angle between the magnetic field and the wall. This angle depends on the current onto the wall. When the Debye sheath vanishes, there is still a potential drop between the wall and the plasma across the magnetic presheath. The magnetic field angle corresponding to sheath collapse is shown to be much smaller than previous estimates, scaling with the electron-ion mass ratio and not with the square root of the mass ratio. This is shown to be a consequence of the finite ion orbit width effects, which are not captured by fluid models. The wall potential with respect to the bulk plasma at which the Debye sheath vanishes is calculated. Above this wall potential, it is possible that the Debye sheath will invert.
What if you are interested in performing a particle simulation with kinetic electrons and gyrokinetic ions? In principle, the kinetic motion of electrons is described in terms of electric and magnetic fields, while the standard gyrokinetic motion of ions is described in terms of electric and magnetic potentials. The dependence of standard gyrokinetic theory on perturbed potentials, instead of perturbed fields, introduces the requirement of specifying a choice of gauge. Of course, you might decide to use the potential representations of the electromagnetic fields in your kinetic description, but that seems wrong to be moving away from physical fields. In this tutorial talk, I will show how standard gyrokinetic theory can be transformed into a gauge-free gyrokinetic theory, which is entirely expressed in terms of the perturbed electric and magnetic fields. A gauge-free gyrokinetic model can be used for hybrid-kinetic simulations of magnetized plasmas in which particle species can be represented in terms of either a Vlasov kinetic description or a gauge-free Vlasov gyrokinetic description. Explicit conservation laws of energy and toroidal angular momentum associated with an axisymmetric background magnetic field can also easily be constructed.
We report our numerical observation of seeded blob dynamics at the SOL region of tokamak. We have used Gkeyll computational plasma framework to perform our 5D gyrokinetic simulation for helical open magnetic field lines. We simulate for plasma parameters similar to ASDEX upgrade experimental shots. Toroidally elongated coherent density structures (known as plasma blobs) are seeded just outside the LCFS. The blobs propagate radially outward and further break into smaller structures depending on the different plasma parameters. The dynamics of the blobs are found to vary as plasma density and temperature is varied. The blob dynamics is also found to be sensitive to the perpendicular and parallel extent of the blobs. A detailed numerical study is performed for various plasma parameters to measure the radial and poloidal blob-velocity as the blobs propagate within the SOL region.