B1 Provision of Support to the Advancement of the ITER Physics Basis

1.2 Energy and particle confinement/transport

1.2-1. Turbulence and transport phenomena

Background and Objectives: In this activity we study drift-type turbulence, large-scale flows and transport phenomena in tokamaks. The analysis of novel experimental results of high performance tokamak plasma discharges has lead to the conclusion that plasma turbulence driven from spatial gradients is the dominant mechanism for the anomalous transport in tokamaks. Furthermore, macroscopic flows, coherent structures or large-scale events (generated by plasma turbulence) may lead plasma to a state of intermittency regulating the transport. Therefore, the description of the dynamics of plasma instabilities and the associated turbulence becomes crucial for the understanding of transport in tokamaks. The main goal is to develop, extend and investigate (both analytically and numerically) advanced models of plasma turbulence, which will be able to describe the complex dynamics of plasma turbulence and justify the features of the observed transport in tokamaks.

Work performed in year 2007 (in co-operation with the Institutes indicated):

1. A series of numerical investigations on the evolution of the Resistive Ballooning Mode (RBM) instability has been performed, using the global fluid code EMEDGE3D which includes parallel dynamics. Extensive analysis of the results was performed using standard and advanced tools for the study of turbulence phenomenology. In this task, we investigated the influence of the electromagnetic effects in the intermittency phenomena of RBM turbulence which are attributed to the existence of coherent structures and large scale flows. We showed that increase of plasma β may lead to multi-fractal scaling of plasma turbulence. This work is almost completed and details can be found in Annex I (in co-operation with CEA).
2. The nonlinear tokamak microturbulence code GYRO was installed in the 64bit cluster of NTUA. The code was developed at General Atomics and uses a fixed (Eulerian) grid to solve the 5-D gyrokinetic-Maxwell equations. A series of numerical calculations using the GYRO code have been initiated in order to extend our previous investigations regarding the stability and the dispersion properties of the coupled Ion Temperature Gradient (ITG) and Trapped Electron Mode (TEM) instabilities, which are mainly responsible for the observed drift-type turbulence and transport. This work is in progress (in co-operation with UofUppsala).
3. Theoretical and numerical tools have been developed for the further and more complete understanding of tokamak transport. To this extent we have performed further numerical gyrokinetic calculations and constructed a transparent gyrofluid turbulence model that predicts the toroidal momentum transport including a new pinch due to the Coriolis force. This model was also used to clarify discrepancies between the gyrokinetic and gyrofluid approaches with respect to the Prandtl number calculations. Comparison of gyrokinetic theoretical calculations (using the codes GS2, GYRO, GKW and GENE) with JET experiments has also been made. Having performed a large number of linear as well as nonlinear gyrokinetic calculations, we have proved that the observed ion temperature profiles are close to the threshold of the linear ITG mode instability. Details on this work can be found in Annex II and Annex III (in co-operation with IPP and UKAEA).

1.3 MHD stability and plasma control

1.3-1. ECRH for MHD control in TEXTOR; NTM avoidance in ASDEX-Upgrade by early application of ECCD

Background and Objectives: The localised application of Electron Cyclotron Resonant Heating (ECRH) and Current Drive (ECCD) has been the most effective technique for the suppression of (classical and neoclassical) tearing modes, which degrade the energy confinement in fusion plasmas. Preliminary experiments on tearing mode control by ECRH in TEXTOR (conducted during 2004-2005) provided substantial evidence that the early application of ECRH (i.e. before the growth of the magnetic island) could be more effective than the late application (i.e. inside a pre-existing saturated island) for the suppression of the tearing mode. The objective of this activity is the theoretical and experimental investigation of the advantage of the early application of ECRH for the suppression of tearing modes, which are generated in TEXTOR by the facility of the Dynamic Ergodic Divertor.

Work performed in year 2007 (in co-operation with IIP-FOM and IIP-FZJ):

1. A theoretical model for the suppression of tearing modes by ECRH was developed, in which the advantage of the early application was attributed to the multiple solutions of the Rutherford equation which are present when localised heating is combined with non-inductive co-current drive. An experimental investigation was conducted (for the assessment of the above theoretical model) in TEXTOR, using dedicated experimental time on the TEXTOR tokamak and the (recently repaired) FOM gyrotron. More specifically we investigated the behaviour of the tearing modes, which are destabilised by perturbation fields (generated by the Dynamic Ergodic Divertor), during early and late application of ECRH, at various powers (200-400 kW), injection angles (co-counter CD), and DED currents (2-3 kA). The major result of the project (with important implications for the ECRH program of ITER) was the confirmation that early application of ECRH has a significant advantage for the prevention of classical tearing modes when it is combined with non-inductive co-current drive, as predicted by the recently developed theory. (For more details, see Annex IV.)
2. Since the defective power generator of ASDEX Upgrade could not be repaired (contrary to what was originally planned) and ASDEX Upgrade will shut down until March 2008, it was agreed with IPP that the activity on “NTM avoidance in ASDEX-Upgrade by early application of ECCD” will be moved to year 2008.

1.3-2. MHD flows and turbulence

Background and Objectives:

The main goal of this subtask is to develop the necessary mathematical models and numerical codes and to carry out parametric studies in order understand the effects of magnetic fields on the flow and heat transfer of conducting fluids as related to the ITER and DEMO fusion machines. The work for 2007 focused on the study of MHD laminar and turbulent flow and heat transfer in cavities and ducts, using finite volume computational fluid dynamics (CFD) codes and direct numerical simulations (DNS) for viscous MHD flows. Stability analysis of MHD flows via finite elements and bifurcation theory were also carried out. Part of these computations was performed on clustered PC systems developed for this purpose.

Work performed in year 2007 (in co-operation with the Institutes indicated):

1. We performed direct numerical simulations (DNS) of natural convection flow of a liquid metal placed between two coaxial vertical cylinders with electrically insulating walls, under the combined effect of a homogeneous magnetic field and internal heat sources. In particular, DNS simulations were performed for the following cases: (a) Transverse (θ = 0o) magnetic field, top and bottom walls adiabatic, and the cylindrical sidewalls kept at equal (constant) temperature or at constant temperature difference, (b) Adiabatic cylindrical sidewalls, constant temperature at the bottom wall, sinusoidal temperature at the top wall, and magnetic field in the transverse or axial direction. Depending on the external Rayleigh number representing the strength of natural convection, the internal Rayleigh number expressing the effect of the internal heat sources, and on the Hartmann number representing the strength of the magnetic field, the flow remains laminar, becomes periodic or transitional, or develops turbulent regions. The increase of the internal and external Rayleigh numbers promotes turbulence and enhances heat transfer. The increase of the magnetic field suppresses turbulence in all cases (see Annex V). (In co-operation with ULB.)
2. All suitable numerical facilities for further development and run of the code SFELES of ULB at UoThly were integrated. The creation of a new computer cluster at the Fluids Laboratory of UoThly helped in this task. More specifically, the toroidal, non-linear and full MHD version of code SFELES was developed. The model is validated against known MHD instabilities (mainly the kink and tearing instability), which are common in fusion plasmas. The simulations are conducted in a cylindrical domain, where well-known analytical solutions exist. Starting from this solution, and increasing the number of modes in the axial direction, a fully resolved transitional and turbulent simulation is performed. Preliminary results of the transitional regime as well as of the turbulent evolution and important statistics were obtained, see Annex VI. The simulations were extended to include transient and turbulent natural convection of liquid metals in cylindrical containers in the presence of mean shear and particles. (This activity will be finalised in early 2008 and reported then. II is performed in collaboration with ULB and CEA Cadarache.)
3. Interpolation methods were incorporated into the computational fluid dynamics (CFD) MHD models in order to calculate the electromagnetic and viscous torque on conducting and insulating rotating walls. A recently developed 3D CFD collocated code for viscous MHD flows, includes new total variation diminishing (TVD) discretisation schemes based on 2nd order interpolation, which were corrected by 3rd order estimations of convective fluxes. The SIMPLE-family of algorithms was tested on fine grids. Natural and MHD natural convection in a cavity were taken as examples. It was found that the SIMPLEX always requires the longest CPU time, the SIMPLER came next, while SIMPLE and SIMPLEC require the least. As far as robustness is concerned, the SIMPLE algorithm is worst, SIMPLER comes next, while SIMPLEX and SIMPLEC are better. The SIMPLEC algorithm is more appropriate, especially for the computation at a fine grid (see Annex VII). A numerical study was also performed to assess 3D effects in the cavity flow, as reported in the PhD thesis of Fidaros (U of Thessaly, July 2007). (In collaboration with FZK)
4. Stability analysis of MHD flows: Development of the 3D stability code has been initiated. More specifically, the governing equations and boundary conditions have been recast in finite element format. The Arnoldi algorithm that calculates the eigenvalues of complex matrices has been implemented. The developed code has been optimised in order to achieve minimum matrix storage requirements. The three-dimensional stability code was tested in a number of simplified cases, and subsequently a parametric study was conducted over a wide parameter range. In the absence of magnetic field and at low Pr = 0.032 the base flow is seen to be more unstable to 3D perturbations, having a smaller threshold value of Gr in comparison to 2D disturbances, leading to a stationary arrangement that is periodic in the longitudinal direction and dominated by small wave-numbers (see Annex VIII). (In collaboration with FZK and ULB.)
5. CFD modelling of MHD rotating flows in shells: The CFD modelling of tangential layers in steady state magneto-hydrodynamic rotating flows was concluded in 2007. The flow in a cylinder with electrically perfectly conducting walls driven by a small rotating disk and the flow in a spherical shell with an inner rotating sphere were studied. The fluid was subjected to an external axial uniform magnetic field. Two types of cores, separated from each other by a tangential layer parallel to the axis of rotation are observed, an inner following a solid-body rotation and an outer quasi stagnant. A counter-rotating jet forms in the tangential layer between the cores. The characteristics of the tangential layer and the properties of the meridional motion are determined for a wide range of Hartmann numbers (see Annex IX). Further numerical tools were developed such as TVD schemes for the convective terms and Lorentz force discretisation, and fast solution algorithms for linear systems of equations. A preliminary study of the corrosive behaviour of conductive fluids on metallic walls was also attempted, as reported in the PhD thesis of Fidaros (U of Thessaly, July 2007). (In co-operation with FZK.)
6. CFD modelling of LiPb flow in the DCLL blanket. A CFD model was developed for the simulation of laminar flow of LiPb in a straight duct of the Dual Coolant Lithium Lead blanket. An extensive literature review was also carried out for both laminar and turbulent MHD flow in ducts. Finally, the study was extended to include turbulent heat transfer under the effect of various magnetic fields in channels. (See Annex X; the work is also presented in the Final Report to ENEA February 2008 for the work performed under a separate agreement and in collaboration with ENEA). The work performed for ENEA showed that the development of very thin Hartmann and side layers in duct flows under the presence of very strong magnetic fields does not allow for CFD simulations to be performed at Ha>1000. However, attempts will be made to perform parametric studies for channel optimisation in 2008 using data available in the open literature.

1.3-3. Use and development of MHD code in relation to liquid metal blanket

Background and Objectives: The liquid metal circulating through the reactor blanket is exposed to strong magnetic fields. While magnetic fields tend to suppress turbulence, curvature and system rotation can, under certain conditions, destabilise the flow and lead to turbulence augmentation. It is therefore important to understand the combined effects of mean shear, strong magnetic fields and system rotation on heat transfer effects. The goal of this activity is to examine the effects of thermal gradients (modelled as a passive scalar) in the presence of strong magnetic fields, mean shear, and system rotation using a fully-3D spectral code. (This is a multi-annual activity performed in co-operation with ULB).

Work performed in year 2007 (in co-operation with ULB):

1. A fully-3D spectral code was modified to include passive scalar transport and the effects of frame rotation. The code has the capability to solve both the quasi-static and full MHD equations on a deforming grid. The new code has been benchmarked and verified against existing literature with the following specific results:

• The code was validated for the case of passive scalar transport in hydrodynamic shear flow under the influence of frame rotation, for which existing data were available, and was found to be in excellent agreement.
• Parametric studies (see Annex XI) had shown that a much higher resolution (10243 instead of 2563) was required for certain combinations of magnetic interaction number (N) and Rosby number (Ro) values. These grid sizes were not achievable with the amount of distributed memory available on the computational cluster.
• The issues identified in (b) were successfully resolved. An efficient memory management module was implemented in the MHD code. At the same time, the computational cluster memory was upgraded from 128 Gb to 256 Gb total distributed memory. The successful completion of these two tasks allowed the execution of the first simulations at a 10243 resolution for testing purposes. Testing has been completed, and production simulation runs on large size computational grids are now underway. These will allow us to establish stability limits for all combinations of magnetic interaction numbers (N) and Rosby numbers (Ro) and to clarify the role of the eddy structure on passive scalar transport. (Workto be continued into next period.)

2. As originally planned and in collaboration with ULB, we initiated a task on Large-Eddy Simulations (LES) of turbulent and pseudo-laminar Hartmann flow near the transitional regime. Previous measurements and works on stability analysis have shown that a critical modified Reynolds number, based on the Hartmann layer thickness R = Re/Ha, exists at approximately R = 380. The LES confirmed that R is the relevant parameter that describes the transition and a critical Reynolds number for re-laminarisation was observed at R  500. However, for turbulent flows at moderate values of the Ha and Re as considered in this study, the similarity with respect to N appears to be better than for R, especially for the first order statistics.

1.4 Power and particle exhaust, plasma-wall interaction

1.4-2. Plasma multiphysics fluid analysis

Background and Objectives: For the accurate analysis of Finite-Element simulations, it is imperative that one uses good quality mesh elements, ideally equilateral; otherwise one violates the assumption embedded in the Finite Element method of linear variation of variables between mesh nodes, which leads to numerical instabilities in the simulations. At the start of year 2007, the algorithms for the adaptive mesh generator were already available emanating from the work of the authors. The objective was to implement the adaptive mesh generator algorithm using software language C++ in an attempt to reduce computational times and improve the accuracy of the results. This is critical for future multiphysics fluid analysis of plasma (electric field, charged and neutral gas particle dynamics analysis) within the fusion reactor. Furthermore, at the start of year 2007, a Finite-Element Flux Corrected Transport algorithm had already been developed for solving drift-diffusion equations. This algorithm was implemented in C++ and numerically solves the conservation of mass equation for neutral gas particles. The objective is to include also the conservation equations for momentum and energy of the neutral gas particles in order to develop the capability to study the full neutral gas particle dynamics, including heating effects in the fusion reactors.

1. To reduce the computational power necessary for multiphysics plasma analysis without compromising the quality of the results, an adaptive mesh generator was developed with the following specific attributes:
   • Normalised error indicators, varying from 0 to 1, for all carriers involved in the analysis, to guarantee optimal mesh resolution.
   • An h- and r- refinement algorithm that uses mesh jiggling, edge swapping and node addition/removal techniques has been developed and successfully tested in two dimensional Cartesian and axisymmetric coordinates, resulting in significant element quality improvements.
   • For an adaptive mesh generator to be efficient, results need to be interpolated from one mesh to the other. Thereby, a fast interpolation algorithm between meshes has been developed. It uses a rectangular grid, equally divided in both directions as a reference grid, to efficiently link the geometries of the two meshes to be interpolated.

The adaptive mesh generator has been tested for highly transient phenomena, such as the avalanche and streamer propagation in regular geometries to start with, such as parallel plates in atmospheric air in order to assess its accuracy and the computational improvement it offers when compared to non-adaptive mesh techniques (see Annex XII for details). Driven by the satisfactory results of the adaptive mesh generator, we have initiated an extension of the successful operation of the adaptive mesh algorithm to irregular geometries and problems such as point planes, where the real advantages of adaptive meshing techniques will be utilised. (Work to be continued into next period.)

2. A compressible Navier-Stokes numerical model that includes the conservation of mass, momentum and energy for the neutral gas particles has been developed in the two-dimensional Cartesian and two-dimensional cylindrical axisymmetric coordinates and coupled to the existing plasma model. The following tasks have been addressed:
   • Formulation of the equations and the various terms that comprise the conservation of momentum and energy equations.
   • Finite Element formulation using the Finite-Element Flux Corrected Transport method (FE-FCT) of the two equations treating each term separately.
   • Implementation of the equations in software programming language C++.
   • Coupling of the Navier-Stokes equations with the charge particle continuity and Poisson’s equations through production and loss processes, such as Joule heating, ionisation, attachment, recombination and momentum transfer.

Testing of the compressible Navier-Stokes solver has been performed through source terms emulating Joule heating effects both in Cartesian and cylindrical axisymmetric coordinates along with validation against published results (see Annex XIII for details).

1.5 Physics of plasma heating and current drive

1.5-1. Transport and chaos in fusion plasmas during ECRH

Background and Objectives: In this activity, we study the wave propagation, absorption and current drive, and the wave-particle interaction in the EC frequency regime. The activity has three parts:

In the first part, we focus on the importance of nonlinear effects, hot plasma effects and (in collaboration with IPP) effects due to the beam shape, in the propagation and absorption of EC waves in tokamak plasmas. Our target is to upgrade the physical schemes used to study the wave propagation, absorption and current drive. In the field of wave propagation in plasmas, the mainstream in theory and applications is oriented to frequency-domain asymptotic methods, where the determination of the plasma response presents less difficulty. However, in many cases of interest (like e.g. mode conversion) this approach breaks down, the solution becomes questionable and therefore a full-wave treatment is called for. We implement, as a wave solver, the Finite Difference-Time Domain method (FDTD), which is nowadays accepted as a reliable tool in numerical electromagnetics. Furthermore, in collaboration with IPP, we consider the propagation and absorption of non-Gaussian EC beams in tokamak plasmas. In experiments involving EC waves, beams with a non-Gaussian amplitude profile can be generated by the launching system or during the propagation in the plasma. The study of non-Gaussian beams is also crucial in forming a theory for the description of the change in the beam shape due to localised absorption, in terms of the generation of higher-order non-Gaussian modes.

In the second part, which has been initiated in this period, we analyse the propagation of EC waves (using ray tracing techniques), their absorption, the current drive and the transport of particles when the NTM instability is excited inside a tokamak. The general idea for the stabilisation of the NTM is to raise the current in the vicinity of the island’s O-point, by utilising the current driven by the EC wave. So far, the analysis of the EC propagation and absorption was done for a stable magnetic topology, ignoring the presence of islands. Our analysis starts with a stable magnetic configuration, which is perturbed and gradually forms a magnetic island. In this unstable geometry, we study the propagation near the reconnection layer, focusing on the effect of the island topology on the efficiency of the EC absorption and ECCD. The plasma response is studied with the use of detailed particle dynamics in the perturbed magnetic topology. The presence of the island affects both the transport of the EC waves in the vicinity of the island and the transport of particles.

Finally, in the third part nonlinear dynamics of wave-particle interactions in fusion plasmas are considered from the point of view of applications to Electron Cyclotron Resonant Heating (ECRH) and its role in the stabilization of the Neoclassical Tearing Modes (NTM), Electron Cyclotron Current Drive (ECCD) and plasma diagnostics through interactions with rf pulses. The Hamiltonian formalism along with a set of accompanying tools such as the phase space analysis, the symplectic integration, the Canonical Perturbation Method and Lie transforms are used in order to extend the theory beyond the quasilinear approximation for realistic fusion plasma configurations and general form of the wave spectra.

Work performed in year 2007 (in co-operation with Institutes indicated):

1. We continued the development of a full-wave code for the description of the propagation and absorption of EC beams in tokamak geometry, based on the Finite Difference-Time Domain method (FDTD). As a benchmark, we studied the perpendicular EC propagation in simplified geometry under different physics models for the dielectric response of the plasma. In general, since FDTD is in the time-domain, conversion from the frequency domain is needed in order to be able to exploit the existing know-how. However, in the case of a constant-frequency wave propagating in stationary plasma, the FDTD method can be applied directly using the frequency-domain dielectric tensor.
2. We initiated the analysis of the propagation, absorption and current drive of EC waves when the NTM instability is excited in a tokamak. In a first attempt, we studied the possibility to use as background magnetic field the perturbed magnetic topology yielded by mapping models like the tokamap. We found though that this method needs to be systematically modified in order to yield the desired information on the magnetic field, and we investigate it in a separate activity [see 3.4.2(iii) and pertinent Annex cited therein]. We then proceeded to start from a standard equilibrium topology in analytical form (simple s-α geometry), which later will be perturbed on the base of the NTM instability. We almost completed the part of the code that calculates the stable background magnetic field, both in Cartesian and curvilinear coordinates.
3. Studies on the self-consistent Vlasov–Maxwell system have been performed. For the case of weak wave-particle interactions, a reduced Vlasov equation for the particle phase-averaged distribution function (having the form of an quasilinear action diffusion equation), obtained in our previous studies, has been considered for its appropriateness for replacing the original Vlasov equation, resulting in a system with lower number of dimensions. These studies have been extended for the cases of stronger wave-particle interactions, with the utilisation of Lie transform techniques. It has been shown that a hierarchy of higher-order action diffusion equations can replace the original Vlasov equation and describe electron dynamics under stronger interactions. These equations have terms with higher-order derivatives of the phase-averaged distribution functions with respect to the action. These terms are shown to be appropriate for describing not only interactions with stronger wave fields, but also effects related to nonlinear resonances, where the particles resonate with the beat of more that one waves. The selection of the appropriate member of the hierarchy of action diffusion equations that can replace the Vlasov equation can be based on the physical parameters determining the strength of the wave-particle interactions and the importance of higher-order (beyond quasilinear) terms. (In co-operation with MIT.)
4. An explicit near-symplectic Hamiltonian mapping, describing the electron motion under interaction with a wave, has been obtained with the utilisation of higher-order Lie perturbation theory. This mapping can be used for arbitrary interaction strength between waves and particles and has several advantages in comparison with previously obtained implicit mappings. The performance of the mapping has been tested against other schemes for a specific Hamiltonian system, describing electron dynamics in an ECRH configuration (see Annex XIV). (In co-operation with U. Craiova and TEKES-Uof Latvia.)
5. Investigation of the application of the Hamiltonian methods and results to other types of wave-particle interactions in fusion plasmas has been performed. We have concluded that both higher-order action diffusion equations [see task (iii)] describing collective particle behaviour, as well as explicit near-symplectic mappings [see task (iv) (PLEASE correct ref.)] describing single particle motion, are applicable for describing wave-particle interactions in ICRH configurations. In comparison to ECRH, qualitatively different interaction features are expected, due to the nonrelativistic character of the ions. (In co-operation with MIT.)
6. Studies of wave-particle interactions in axisymmetric toroidal magnetic field configurations have been initiated. In order to apply Hamiltonian methods, a canonical set of variables have been used for describing the particle motion in the absence of wave fields. The resulting Hamiltonian system is shown to be integrable for the case of an axisymmetric equilibrium. The wave fields, as well as any deviation of the static magnetic field from axisymmetry, can be considered as a perturbation and standard methods can be applied in order to derive an action diffusion equation describing particle and momentum transport, under the presence of ECRH and magnetic island perturbations. (This task will be continued in year 2008, in co-operation with MIT.)

1.6 Energetic particle physics

1.6-1. Alfvén wave-particle interactions at sub-Alfvénic velocities

Background and Objectives: The aim of this project is to develop a coherent theory capable of describing the experimental results on sub Alfvénic resonances observed in JET. In particular, the interest lies in the understanding of the behaviour of α-particles originating from nuclear reactions in deuteron–deuteron or deuteron–triton thermonuclear reactions, or due to injection of neutral beams.

Work performed in 2007 (Partly in co-operation with JET):

1. Due to unforeseen circumstances related to the co-operation with JET, the work outlined in the “Background and Objectives” section had to be redirected to the related question of the role of the geodesic acoustic mode in the broadening of the profile of the otherwise sharp Alfvén modes It is being shown that besides the usual geodesic mode there exist another pressure (or, rather, ‘gradient of pressure’) mode perpendicular to the direction of the geodesic curvature. Furthermore it has been proved that the ‘geodesic’ mode is the usual ExB ambipolar drift. Details are presented in Annex XV (Work on the subject is still in progress.)
2. Work initiated during the previous visits at JET on the equilibrium of toroidal plasma configurations like JET/ITER, when the ‘Shafranov’ shift is not constant, culminated in the derivation of the equivalent of Grad-Shafranov equation, as presented in Annex XVI. (The work has been presented to EPS-34, Warsaw 2007.)
3. Work has been performed on the energy transfer through Landau resonances in non-thermalised plasmas satisfying distributions functions of the Weibull or kappa- type. A presentation of those results, initially scheduled for December 2007 (eventually realized in January 2008) was made at Culham at the invitation of TF Team. Some of the ideas were considered worth pursuing further in the experimental JET Campaign of February 2008.