Enhanced Preconditioner for JOREK MHD Solver

The JOREK extended magneto-hydrodynamic (MHD) code is a widely used simulation code for studying the non-linear dynamics of large-scale instabilities in divertor tokamak plasmas. Due to the large scale-separation intrinsic to these phenomena both in space and time, the computational costs for simulations in realistic geometry and with realistic parameters can be very high, motivating the investment of considerable effort for optimization. In this article, a set of developments regarding the JOREK solver and preconditioner is described, which lead to overall significant benefits for large production simulations. This comprises in particular enhanced convergence in highly non-linear scenarios and a general reduction of memory consumption and computational costs. The developments include faster construction of preconditioner matrices, a domain decomposition of preconditioning matrices for solver libraries that can handle distributed matrices, interfaces for additional solver libraries, an option to use matrix compression methods, and the implementation of a complex solver interface for the preconditioner. The most significant development presented consists in a generalization of the physics based preconditioner to "mode groups", which allows to account for the dominant interactions between toroidal Fourier modes in highly non-linear simulations. At the cost of a moderate increase of memory consumption, the technique can strongly enhance convergence in suitable cases allowing to use significantly larger time steps. For all developments, benchmarks based on typical simulation cases demonstrate the resulting improvements

Probing non-linear MHD stability of the EDA H-mode in ASDEX Upgrade

Regimes of operation in tokamaks that are devoid of large ELMs have to be better understood to extrapolate their applicability to reactor-relevant devices. This paper describes non-linear extended MHD simulations that use an experimental equilibrium from an EDA H-mode in ASDEX Upgrade. Linear ideal MHD analysis indicates that the operational point lies slightly inside of the stable region. The non-linear simulations with the visco-resistive extended MHD code, JOREK, sustain non-axisymmetric perturbations that are linearly most unstable with toroidal mode numbers of n = \{6 \dots 9\}, but non-linearly higher and lower n become driven and the low-n become dominant. The poloidal mode velocity during the linear phase is found to correspond to the expected velocity for resistive ballooning modes. The perturbations that exist in the simulations have somewhat smaller poloidal wavenumbers (k_{\theta} \sim 0.1 to 0.5 cm^{-1} ) than the experimental expectations for the quasi-coherent mode in EDA, and cause non-negligible transport in both the heat and particle channels. In the transition from linear to non-linear phase, the mode frequency chirps down from approximately 35 kHz to 13 kHz, which corresponds approximately to the lower end of frequencies that are typically observed in EDA H-modes in ASDEX Upgrade

MHD stability in X-point geometry:simulation of ELMs

A non-linear MHD code, named JOREK, is under development with the aim of studying the non-linear evolution of the MHD instabilities thought to be responsible for edge localized modes (ELMs): external kink (peeling) and medium-n ballooning modes. The full toroidal X-point geometry is taken into account including the separatrix, open and closed field lines. Analysis of the influence of the separatrix shows a strong stabilization of the ideal and resistive MHD external kink/peeling modes. One instability remains unstable in the presence of the X-point, characterized by a combination of a tearing and a peeling mode. The so-called peeling-tearing mode shows a much weaker dependence on the edge q. Non-linearly the n = 1 peeling-tearing mode saturates at a constant amplitude yielding a mostly kink-like perturbation of the boundary with an island-like structure close to the X-point. The non-linear evolution of a medium-n ballooning mode shows the formation of density filaments. The density filaments are sheared off from the main plasma by an n = 0 flow non-linearly induced by the Maxwell stress. The amplitude of the ballooning mode is limited by this n = 0 flow and multiple (in time) density filaments can develop to bring the plasma below the stability boundary

Equilibrium flows in non-linear MHD simulations of x-point plasmas

In non-linear MHD simulations of ELMs [1], a radially localised, toroidally symmetric, poloidal flow layer exists in the H-mode pedestal region. This sheared flow layer could have a significant influence on the linear stability properties of MHD instabilities and their non-linear evolution. Using the non-linear MHD simulation code JOREK [1] with reduced resistive MHD equations, we study the edge-localised poloidal flow in both circular and X-point tokamak plasmas at equilibrium (toroidal symmetry). For the circular case, an analytical interpretation is derived. In the simulations of X-point plasmas, the flow can have both m = 0 and m = 1 components. In fact, abrupt transitions take place between the two equilibrium states, accompanied by a strong increase in the kinetic energy. Similar transitions between equilibrium flow states have been predicted by Strauss [2] for m = 0 poloidal flow patterns. Scalings are obtained for both the m = 1 and m = 0 flows

Progress in understanding disruptions triggered by massive gas injection via 3D non-linear MHD modelling with JOREK

3D non-linear MHD simulations of a D 2 massive gas injection (MGI) triggered disruption in JET with the JOREK code provide results which are qualitatively consistent with experimental observations and shed light on the physics at play. In particular, it is observed that the gas destabilizes a large m/n  =  2/1 tearing mode, with the island O-point coinciding with the gas deposition region, by enhancing the plasma resistivity via cooling. When the 2/1 island gets so large that its inner side reaches the q  =  3/2 surface, a 3/2 tearing mode grows. Simulations suggest that this is due to a steepening of the current profile right inside q  =  3/2. Magnetic field stochastization over a large fraction of the minor radius as well as the growth of higher n modes ensue rapidly, leading to the thermal quench (TQ). The role of the 1/1 internal kink mode is discussed. An I p spike at the TQ is obtained in the simulations but with a smaller amplitude than in the experiment. Possible reasons are discussed

Linear MHD stability analysis of post-disruption plasmas in ITER

Most of the plasma current can be replaced by a runaway electron (RE) current during plasma disruptions in ITER. In this case the post-disruption plasma current profile is likely to be more peaked than the pre-disruption profile. The MHD activity of such plasma will affect the runaway electron generation and confinement and the dynamics of the plasma position evolution (Vertical Displacement Event), limiting the timeframe for runaway electrons and disruption mitigation. In the present paper, we evaluate the influence of the possible RE seed current parameters on the onset of the MHD instabilities. By varying the RE seed current profile, we search for subsequent plasma evolutions with the highest and the lowest MHD activity. This information can be applied to a development of desirable ITER disruption scenario

Magnetohydrodynamics modelling of H-mode plasma response to external resonant magnetic perturbations

The response of an H-mode plasma to Resonant Magnetic Perturbations (RMPs) generated by so-called I-coils in DIII-D experiments on type I edge localized modes suppression is modelled using the nonlinear reduced magnetohydrodynamics (with zero-β, i.e. zero plasma temperature, in the version used here) code JOREK in X-point geometry. JOREK self-consistently advances in time the magnetic flux, vorticity, and plasma density in the presence of the RMPs. Without any toroidal rotation, the magnetic response from the plasma does not significantly modify the islands widths. A radial convective E⃗×B⃗ plasma transport is observed to occur in the presence of the RMPs. The possibility that this mechanism could explain the enhanced density transport observed experimentally in DIII-D is discussed. Simulations with a rigid-body-like rotation at a fixed velocity shows evidence of a screening of the RMPs. The extension of our results to realistic parameters is discussed

Performance analysis and optimization of the JOREK code for many-core CPUs

This report investigates the performance of the JOREK code on the Intel Knights Landing and Skylake processor architectures. The OpenMP scaling of the matrix construction part of the code was analyzed and improved synchronization methods were implemented. A new switch was implemented to control the number of threads used for the linear equation solver independently from other parts of the code. The matrix construction subroutine was vectorized, and the data locality was also improved. These steps led to a factor of two speedup for the matrix construction

Effect of divertor plasma conditions & drifts on ELM power fluxes at the ITER divertor targets

\u3cp\u3eThe effect of divertor recycling, in/out recycling asymmetries and ion ∇B drift direction on in/out divertor power and particle flux asymmetries for stationary plasma conditions and during ELMs have been modelled with the 2-D PARASOL PIC kinetic code. The direction of the ion ∇B drift has a strong effect on the steady-state in/out heat/particle flux divertor asymmetries and this effect is even larger for ELMs. The modelled changes of the in/out divertor asymmetries with ∇B for steady state conditions are contrary to experimental findings (the inner divertor heat flux does not become similar to outer one when ∇B direction is reversed). Simulations of ELMs find that the energy load to the inner divertor is largest for normal ∇B and smallest for reversed ∇B for an ELM energy loss of δW\u3csub\u3eELM\u3c/sub\u3e/W-12%. This finding is robust to modelling assumptions (recycling value, in/out recycling ratio, ELM energy loss magnitude and plasma collisionality). This is good qualitative agreement with experiment, although magnitude of the predicted changes is much larger than in experiment. The magnitude of δW\u3csub\u3eELM\u3c/sub\u3e/W itself is also found to affect the in/out ELM energy deposition asymmetry and E\u3csub\u3ein\u3c/sub\u3e/E\u3csub\u3eout\u3c/sub\u3e increases with δW\u3csub\u3eELM\u3c/sub\u3e/W for normal ∇B while it decreases with reversed ∇B for low collisionalities. Further 2-D PARASOL simulations to study the role of drifts, recycling and thermoelectric currents are in progress to refine these findings.\u3c/p\u3