## CHAPTER B - B3 Development of Concept Improvements and Advances in Fundamental Understanding of Fusion Plasmas

## 3.2 Understanding of plasma characteristics for improved concepts

### 3.2.1.Analysis of turbulent structures in fusion plasmas

Background and Objectives: Coherent turbulence structures play an important in the dynamics of turbulence and modify the transport of (passive) scalars and contribute to the preferential concentration of transported particles. In the case of edge turbulence in fusion plasmas, it is important to be able to model the effect of coherent structures on density and electric potential fields. The turbulence structure dimensionality tensor is a one-point tensor that carries non- local information providing a statistical description of coherent structures. This activity is aimed at further developing the dispersed phase structure dimensionality tensor diagnostic so that it can become integrated in the electrostatic edge-SOL (ESEL) turbulence code.

Work performed in 2011:

(i)In the first part of 2011, initial simulations for structure tensors in edge plasmas were computed using the edge-SOL code ESEL. The structure tensors were computed a posteriori using the results generated from the code. Larger simulations carried out in the second half of 2011 are currently being analysed, while further simulations and DPSD computations will be done in 2012 in an effort to understand the correlation between the computed edge plasma parameters and the structure tensors.

### 3.2.2.Magnetic field modelling in plasma multiphysics code

Background and Objectives: A two dimensional axisymmetric and a three-dimensional code, employing the transient solution of the continuity equations for charged particles coupled with the solution of the electric field, have been previously developed. The three dimensional code has the potential to simulate much more complex and realistic Tokamak configurations. However, there is a further need to incorporate the influence of the magnetic field, an effect that had not been included in the three-dimensional code. As a result our effort has been first of all to develop further the fully three-dimensional code in order to incorporate magnetic field effects.

Work performed in 2011:

(i) Two-dimensional axisymmetric plasma simulations were initially performed during the first part of 2011. The fully three-dimensional code was then developed further in order to include magnetic field effects. The augmented three-dimensional code was initially validated against the two-dimensional axisymmetric code for a simple parallel plane configuration, and subsequently, the code was utilised in the solution of plasma related problems in progressively complex geometries. The code development will be completed in early 2012.

## 3.4 Theory and modelling

### 3.4.1.Stationary MHD modes in magnetically confined plasmas

Background and Objectives: This activity aims at constructing equilibria and relaxed states of fusion plasmas with flow or/and finite electrical conductivity and investigating their linear and nonlinear stability in connection with the ITER project.

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

(i) Construction of equilibrium configurations with incompressible flow parallel to the magnetic field and ITER pertinent shaping (lower X-point) as well as confinement figures of merit on the basis of generic linear analytic solutions has been carried out in the framework of the EFDA Physics MHD Topical Group. Application of a sufficient condition for linear stability of equilibria with parallel flow [Throumoulopoulos, Tasso, Phys. Plasmas 14, 122104 (2007)] for the equilibria constructed showed that this condition is not satisfied. However, since the condition is sufficient it does not imply that the equilibrium is unstable. This result confirms a conjecture that equilibrium nonlinearity may activate flow stabilising effects (see **Annex 25**) (in co-operation with MEdC and IPP).

(ii) A class of exact axisymmetric tokamak equilibria with sheared flows parallel to the magnetic field has been constructed, generalising previous work on the subject [Ap Kuiroukidis, Plasma Phys. Control. Fusion 52, 015002 (2010)]. The additional free parameters associated with new terms in the solution make it possible to construct up-down asymmetric configurations with a divertor X-point and desirable values of confinement figures of merit as the safety factor on the magnetic axis and plasma betas. In particular, we construct a number of ITER-pertinent equilibria. Their stability with respect to linear MHD perturbations is also examined by applying the sufficient condition mentioned in item (i). Details are given in **Annex 26**.

(iii) In the framework of the EFDA Tokamak Integrated Modelling Task Force we have extended the HELENA code to incompressible plasma rotation parallel to the magnetic field. In this case the generalised Grad Shafranov equation by means of an integral transformation involving the Alfvén Mach function becomes identical in form with the usual static equation. Thus, the extension has been accomplished by using the existing static version of the code in conjunction with the Alfvén Mach function (see **Annex 27**) (in co-operation with IPP).

(iv) A general stability condition for plasma-vacuum systems with resistive walls has been derived by using the Frieman Rotenberg Lagrangian stability formulation [Rev. Mod. Phys. 32, 898 (1960)]. It is shown that the Lyapunov stability limit for external modes does not depend upon the gyroscopic term but upon the sign of the perturbed potential energy only. In the absence of dissipation in the plasma such as viscosity, it is expected that the flow cannot stabilise the system. Details about this study can be found in [H. Tasso and G.N. Throumoulopoulos, Phys. Plasmas 18, 070702 (2011)] (in co-operation with IPP).

### 3.4.2.Turbulence and anomalous transport phenomena

Background and Objectives: This activity addresses several tasks:

In a first part of the activity, we study particle and heat transport in turbulent environments, with the use of realistic models of turbulence in toroidal topologies, and with particular focus on the edge region, impurities, fast particles, and anomalous transport behaviour. Several tools are used in parallel in this study: (1) We perform test-particle simulations in turbulent electromagnetic field environments. (2) A code for the solution of the Fokker-Planck equation is developed, as an alternative tool to the test-particle simulations and the Langevin equation, and moreover for comparison to the ETS code. (3) Random walk models and fractional diffusion equations are developed as tools to describe anomalous transport, with particular aim to describe in parallel the combined diffusion in position and velocity space. (4) Eulerian- Lagrangian CFD models for the study of near-wall motion of particulate matter and its interaction with stochastic turbulent fluid flow in the presence of deterministic/stochastic magnetic fields.

In a second part of the activity, we explore particular turbulence models that are based on the concept of Self-Organised Criticality (SOC), where the turbulent plasma is considered as a complex system, in which localised instabilities are relaxed in small-scale diffusive events. SOC is usually modelled in the form of Cellular Automata and the sand-pile analogy is made use of, yet it is a particular aim of our approach to use the natural physical variables and that the SOC models have a consistent physical interpretation. In particular, we develop the Extended Cellular Automaton (X-CA) model for the magnetic field, which is fully MHD compatible, as well as a SOC model for turbulence driven by micro instabilities. In parallel, we build the MHD code MYDAS2, as an alternative to the X-CA model and in order to validate results of the latter. Also, gyrokinetic simulations are performed, with the gyrokinetic code GENE, in order to study turbulence driven by micro-instabilities (ITG, TEM, and ETG mode driven instabilities).

The last part of this activity utilises Hamiltonian methods for obtaining evolution equations of the particle distribution functions under the presence of turbulent or other perturbing modes. These equations describe anomalous transport phenomena more accurately than standard theories, such as the quasilinear theory, since they incorporate the actual underlying particle dynamics and are not based on commonly used statistical assumptions with questionable validity in realistic cases.

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

(i) The module of Lagrangian particle tracking (LPT) has been modified by including the electromagnetic force. The new LPT module has been successfully verified in certain simple test cases. A parametric study has been initiated to investigate the interaction of charged particles with turbulent MHD flows. The transport and deposition of charged particles in wall-bounded turbulent flows with and without heat transfer in the presence of uniform magnetic fields has been investigated. It is found that the electromagnetic force substantially affects the trajectories of the charged particles. A wall-normal magnetic field does not affect directly the motion towards the walls, while the streamwise and spanwise magnetic fields produce forces that enable them to reach the walls. For a magnetic field normal to the channel walls, particle deposition is reduced, while it is increased for streamwise and spanwise magnetic fields. In the latter case, particle deposition is also larger than that found in the corresponding cases with neutral particles and/or ‘bubbles’. Particle dispersion was also investigated in anisotropic wall-bounded turbulent flows with electromagnetic fields produced by the numerical solution of the 3D full MHD equations (see **Annex 28**).

(ii) Data obtained from EMEDGE3D code have been found time consuming. The code did not involved beyond its first version since the developer (now in Marseille) is now fully concentrated in a different scientific area. Thus, we decided to build a simpler, though more focused to the mechanism behind the blob generation, physical model aiming towards its numerical implementation. Therefore, we succeeded in formulating an analytical model for evaluating the statistics of the drift wave turbulence. This effort was a collaboration with IST-Portugal through the Portuguese Mobility Contract. A journal publication was the outcome of this effort: J. T. Mendonca, K. Hizanidis, “Improved model of quasi-particle turbulence (with applications to Alfven and drift wave turbulence)”, Phys. Plasmas 18, 112306 (2011). The model will be numerically implemented in order to provide spectral results in the framework of our Work Plan 2012-2013.

(iii) A Matlab code based on the FP model is operational. The code, for the time being, involves two actions and the inhomogeneous, time dependent diffusion tensor can handle turbulent fields as well. Extension to a 3D model (involving 3 actions), in FORTRAN is under way. We also plan to couple this code, as far as the fields are concerned, with the one which will be based on the model already developed in task (ii).

(iv) No work was done in the period in subject, since priority was given to the further development of the MHD code MYDAS, which will yield useful information that will allow to further develop the X-CA SOC model. (AUTH).

(v) No work was done in the period in subject, since priority was given to the modeling of micro-instabilities with a standard tool, the gyrokinetic code GENE, which will yield transport coefficients that can be compared with those from the SOC model (in collaboration with Univ. of Bayreuth). (AUTH)

(vi) Transport coefficients of ions, as well as distributions of final positions in radial direction, have been determined in a magnetic island topology (NTM), for particles near the X- and near the O-point, and they have been compared to the transport coefficients in the unperturbed topology. Since the Lorentz force was integrated in the simulations, the computation times were large for obtaining good enough statistics, and the test-particle code has been parallelised over the particles, which move independently of each other, leading to a very large reduction in computational cost. First runs have been made for electrons, they turned though out to be still too demanding on computational cost when using the Lorentz force, so an appropriate guiding center approximation will have to be made. Also, with a variant of the test-particle code, magnetic field lines can be traced, and Poincare surfaces of sections in the poloidal plane have been determined, in order to visualise the magnetic topology and to locate the X- and the O-point. For more details, see **Annex 29** (in collaboration with, MEdC, Univ. of Craiova).

(vii) In order to study the current drive in magnetic islands (ECRH/CD), the Fokker-Planck code CHET has been extended from 1 to 2 dimensions in velocity space, and has successfully been tested. In a next step, the coordinate system has been changed from Cartesian to spherical, taking into account the radial (v) and azimuthal θ (angle with the background magnetic field) components. In order to treat the coordinate singularities, the radial grid was extended over the pole to negative radii, with even number of grid points so that the pole is not part of the grid. The azimuthal grid was extended from half-circle to full circle, such that it becomes periodic, and the equi-spaced grid points have a small constant phase such that the θ = 0, θ = p axis is excluded. The radial derivative is treated with Chebyshev expansion (and applying parity corrections at the negative radii), whereas the azimuthal derivatives by Fourier expansion. The grid and pseudospectral differentiations were successfully tested with simple diffusion test cases. The code CHET was also parallelised, and the computation times could indeed be reduced. Last, a collision operator is being implemented. For more details, see **Annex 30**. Concerning the ETS, an alternative basic integrator had been written in the previous year, based on the pseudo-spectral method. In this year, the time-integration has been modified from 4^{th} order Runge-Kutta to an implicit Euler method, based on a matrix formulation of the pseudospectral, Chebyshev-based derivatives, and exploiting the linearity of the equations solved by the ETS. The implementation has been completed and been benchmarked against a different integrator (the Progonka solver, number 3), with our own driving routine, though. It remains to include our solver within the ETS itself and to benchmark it there again. Finally, no work has been done on determining the transport coefficients in the SOC model for ITG mode driven turbulence (with aim their implementation in CHET) (in collaboration with Univ. of Bayreuth; ITM IMP3).

(viii) The numerical solution of the FTE (with the Grunwald Letnikov definition of fractional derivatives in combination with the method of lines) had already been established. In the period in subject, it has been clarified that one cause for anomalous transport in the FTE is an assumed asymmetry between the left and the right fractional derivatives. So-far, the effect is though too small to explain e.g. profile stiffness under conditions of off-axis heating, the mechanism has still to be further investigated (in collaboration with MEdC, Univ. of Craiova).

(ix) Several modifications and improvements of the MHD code MYDAS were made: (1) Time-stepping has been modified and is done for the variables in Fourier space (transformed in the poloidal direction), and (2) several time-stepping schemes have been evaluated for stability and speed, and finally fixed step-size 4^{th} order Rung-Kutta was considered most adequate and is being used. (3) The equations have been re-implemented in flux- conservative form, and (4) the resistive term has been added to the induction equation. (5) The pole problem (coordinate singularity) has been treated in a new way, no pole equations are used anymore, but the grid has been modified in the radial direction to include also negative radii, with an even number of grid points such that the pole is excluded (and taking care of the correct parity at the negative radii). (6) Moreover, the possibility to solve the MHD equations in toroidal coordinates has been added, with a simple switch that allows changing between cylindrical and toroidal coordinates. (7) The code has been parallelised with OpenMP. (8) A Perfectly Matched Layer (PML) had been implemented as absorbing boundary condition, and it worked in the tests made, it is though not very practical in that it uses auxiliary variables and equations that have to be tailored to the original equations. For this reason, the PML was replaced by a simpler Absorbing Layer, in which a spatial (dissipative) filter gradually imposes the boundary conditions towards the edge (9) The de-aliasing was improved by replacing the simpler 2/3-rule by an exponential low pass filter in Fourier-Chebyshev-space. Finally, a new general family of analytical solutions (equilibria, and propagating and standing waves) in polar coordinates has been derived as test-cases, and MYDAS reproduced them successfully, in ring-shaped domains as well as in the entire disk, and, with the new de-aliasing, the non-linear test case of the Orszag-Tang vortex could be reproduced (colliding shock-waves in Cartesian coordinates). For more details, see **Annex 31**. Because of the various technical and physics improvements of MYDAS, work on the study of edge turbulence (ELMs) has not been done yet.

(x) Towards the end of the year in subject, the study of turbulence driven by micro-instabilities (ITG, TEM, ETG modes) has been initiated, with focus on electron heat transport under conditions of multi-scale turbulence. The 3^{rd} party gyrokinetic code GENE (IPP Garching) has been installed on a local machine, and is running in parallel execution (OpenMP and MPI), as needed. A number of small test runs have been performed, in order to familiarise with the code, its initialisation parameters, its potential use, as well as the diagnostic tools for the analysis of the results of the code. (In collaboration with IPP Garching.)

### 3.4.3.Discrete kinetic and stochastic models for transport (UoThly)

Background and Objectives: Discrete kinetic simulations may provide an alternative mesoscale type of description of the transport properties of complex physical systems. The main unknown is a particle distribution function, which obeys a suitably chosen kinetic equation, while the bulk quantities of practical interest are obtained by taking moments of the distribution function. The benefit of this description, compared to the more traditional micro- and macro-scopic approaches, is due to the fact that physics, at particle level, may be investigated with moderate computational effort. This activity aims to develop kinetic computer codes capable to simulate in a computationally efficient manner non-linear resistive MHD flows and vacuum systems in fusion devices.

Work performed in year 2011 (in cooperation with KIT):

(i) In several occasions the vacuum and pumping systems of DT fusion reactors may operate under transient conditions, e.g. due to starting or ending pumping operation or due to unsteadiness in plasma scenarios. The starting gas flow in a channel has been investigated in the whole range of the Knudsen number by numerically solving the governing time dependent BGK kinetic equations in a fully deterministic manner. In particular, a time- dependent kinetic code has been developed and implemented to the transient solution for starting fully developed gas flows in circular and rectangular channels. The gas is initially at rest and then due to a suddenly imposed uniform pressure gradient, is starting to flow. The motion is time-dependent up to the point where the steady-state flow conditions are recovered. The solution provides a detailed description of the evolution of the flow field with regard to time from the starting point, where the gas is at rest up to a certain time where almost steady-state conditions are recovered. Based on the results some insight of how rapidly a vacuum flow will respond to a sudden change, related to an externally imposed pressure gradient coming from a vacuum pump or a valve, is obtained (see **Annex 32**). The developed algorithm may also be applied to more complex time-dependent systems including leak detection applications

(ii) The recently developed deterministic nonlinear kinetic algorithm has been implemented to solve vacuum flows through circular and rectangular channels of very small up to moderate lengths for various pressure ratios. As the ratio of the length over the hydraulic diameter of the channel is gradually increased the nonlinear solution tends to the corresponding linear fully developed solution. In addition, as the difference between the inlet and outlet pressure is decreased the nonlinear results recover the corresponding linear results for any channel length. A comparison with DSMC results is also performed in order to further validate and benchmark the deterministic code. Very good agreement in all cases has been observed (see **Annex 33**). In addition, good agreement with experimental results available in the literature has been observed. This work is useful in order to simulate all types of vacuum flows through pipe elements as well to develop a complete data base to be integrated in the pipe network code for modelling gas distribution systems of arbitrary complexity under any vacuum conditions. In addition to the above an investigation has been performed in order to quantify the influence of the channel ends by applying the effective length concept. The main idea is to extend the geometry of the channel by additional length increments for the purposes of calculating the correct pressure profile, justifying the name of the approach. The inlet and outlet parts of the channel are considered separately and the deviation from the fully developed profile is expressed via the effective length concept. The values extracted in this work can be used to efficiently obtain practically important quantities, such as the mass flow rate through the channel, by a simple procedure and within only a few minutes. Therefore, the well-known linear fully developed analysis accordingly modified by the input of the present work is applicable not only for long channels where L/R>80 as in the traditional solution of the problem but it is extended also for the case of moderately long channels where L/R>10. The results clearly indicate that for each channel cross section the effective length depends only on the rarefaction parameter (see **Annex 34**).

(iii) Vacuum gas flows through single pipe elements of various cross sections have been investigated over the years both numerically and experimentally. In fusion related vacuum applications, including the vacuum systems of DT reactors, these single pipe elements are combined together in order to form a pipe network. Computational algorithms for solving gas pipe networks in the viscous regime are well developed but corresponding tools for solving pipe networks under any degree of gas rarefaction have not been developed so far. Very recently such an algorithm has been developed by the UoThly team in the case of pipe networks consisting of long circular pipes. This algorithm has been generalised and extended to pipe networks consisting of channels of arbitrary cross-section as is the case at the vacuum pumping network of ITER. The analysis is valid and the results are accurate in the whole range of the Knudsen number, while the involved computational effort is very small. This is achieved by successfully integrating the kinetic results for the corresponding single channels into the general solver for designing the gas pipe network. To demonstrate the feasibility of the approach two typical systems consisting of long rectangular and trapezoidal micro-channels are solved (see **Annex 35**). The extension of the algorithm to networks consisting of pipes of arbitrary length although initially planned for 2011, it will be performed during the next period. The completion of this latter step will lead to the development of an advanced computational tool for the design of the vacuum systems of DT fusion reactors.

### 3.4.4.2-D MHD code for plasma trapping in fusion devices including specific modules on magnetic field topology, nuclear reaction products, neutronics and neutron-pellet interaction and “Development of an alternative source and neutralizer technique aiming at increased efficiency of neutral beam heating systems”. [Task Agreement HCD-02-02: Neutral Beam Advanced Technologies]

Background and Objectives:

The main aim of our investigations and numerical simulations the last few years concern the development of a Compact Magnetic Fusion (CMF) device in open magnetic field topology for high neutron flux production and study fusion process in high density and high temperature plasmas. Details on the physics and the 2-D MHD resistive numerical code used to study the CMF device, the plasma production in the open magnetic field topology and the plasma trapping in the externally applied high magnetic field were presented in previous annual reports and the corresponding annexes. At present we investigate (a) the initial plasma formation in the CMF device using non-linear laser propagation and filamentation effects during the laser- plasma interaction in order to increase the interaction volume and improve the energy transfer from the laser beam to plasma, (b) the refuelling of the compact magnetic fusion device using pellet techniques. One of the future promising schemes for Tokamak heating is the use of neutral high energy and high current hydrogen (or D) beam interaction with the torus plasma. The main objective is the development of a negative hydrogen beam (or D beam) and the neutralisation using laser based techniques.

Work performed in year 2011:

(i) Numerical investigations, performance and development of a high neutron flux source: Significant progress on numerical studies on fusion process in high density and high temperature plasma produced by a high intensity laser pulse with clusters in a high external magnetic field, enable us to propose the development of a compact CMF device for high flux neutron production. The non-linear propagation of a laser beam in plasma allows increasing the interaction volume by an important factor. A 2-D MHD code enables to calculate the spatio-temporal evolution of plasma parameters for different initial spatial plasma distribution and calculate the plasma trapping and neutron production. The results from numerical simulations and the use of advanced techniques concerning the development of the compact CMF device confirm (i) the efficiency of the proposed configuration which combines laser self-guiding propagation in plasmas and pulsed power techniques for high magnetic generation in mirror-like topology and (ii) the high neutron flux production. The work is oriented (see **Annex 36**) to propose the assembling of the proposed investigations, technologies and calculations for the development of an experimental device for neutron material testing for Tokamak applications. The development of the neutron source, with a rep.rate of 10-100 Hz, requests an efficient refueling scheme (or system). The injection of deuterated gas (or clusters) is very expensive because only a small amount of the volume of the injected gas interacts with the focused beam. In order to improve the operation of the CMF device we initiate calculations on the potential injection of a small pellet with approximately zero initial velocity (due to small radial dimension of the CMF device) in the device. The volume in the mirror-like magnetic configuration is up to 0,5-1 cm^{3}. Our investigation shows that the plasma density can be maintained up to 10^{17}- 10^{18} cm^{-3} in the device for the injection of a pellet of a diameter of 10^{-2} mm and a solid density of 10^{28} m^{-3}. For these calculations we do not take into account the plasma losses due to the open magnetic field topology of the device. The fabrication of such small pellets injected with zero initial velocity does not request to resolve serious technical problems for potential future applications within the proposed CMF device. Near future laser systems allow operation of the high intensity laser beam at 100 Hz with MW average power which enables to couple the proposed CMF device with a laser beam and the transfer of complementary amount of energy to plasma, produced by bigger pellets, in order to cover the losses.

(ii) Investigations on Negative Ions: The proposed work for the negative ion production, acceleration and neutralisation is based on the study of magnetically isolated diodes coupled with a pulsed power generator and the neutralisation of negative ions (H or D) by laser photo-detachment techniques. We initiate studies on (a) numerical calculation of a pulsed power generator working from 20 ns up to µs in order to produce 500 kV and 200 A and (b) the development of a 1-D, two fluid model in cylindrical geometry describing the spatio-temporal evolution of negative ions (H or D) and electrons in a magnetically isolated diode which is between two concentric cylinders. The model equations include conservation of particles, momentum and energy of electrons and negative ions, coupled with the Poisson’s equations for the space charge separation. The code will allow describing the final negative ion extraction current for different externally applied electric and magnetic fields.

Last Updated (Friday, 11 January 2013 12:47)