Impact of the improved parallel kinetic coefficients on the helium and neon transport in SOLPS-ITER for ITER

New Grad's-Zhdanov module is implemented in the SOLPS-ITER code and applied to ITER impurity transport simulations. Significant difference appears in the helium transport due to improved parallel kinetic coefficients. As a result 30\% decrease of the separatrix-averaged helium relative concentration is observed for the constant helium source and pumping speed. Change of the impurity behaviour is discussed. For the neon changes are less pronounced. For the first time the ion distribution functions are studied in the ITER Scrape-off layer conditions to reveal the origin of the kinetic coefficient improvements and theory limitations.


Introduction
Impurity transport is a crucial problem in fusion research. Impurity accumulation in the confinement region of fusion devices should be avoided. Edge and Scrape-off layer (SOL) transport affects significantly global impurity transport. The SOLPS-ITER code has been developed to study transport in the SOL and Edge [1,2]. In the standard SOLPS-ITER model Zhdanov-Yushmanov expressions [3] are used for the thermal and friction forces coefficients calculation [4]. However, for the ions with close masses the complete Grad's-Zhdanov method should be used [3]. It is important to mention that multi-temperature generalisation of this method is in progress for the Soledge3x-EIRENE code [5,6].
In the present paper a new Grad's-Zhdanov module based on the improved approach [7] and complete Grad's-Zhdanov method [3,7] has been implemented in the SOLPS-ITER code (Sec. 2). Impact on the helium and neon transport is studied (Sec. 3). The new approach leads to changes in the parallel impurity momentum balance, which causes the poloidal impurity stagnation point to shift in the divertor region. The solution is very sensitive to the stagnation point variations. A simple model is proposed to explain how small changes of the stagnation point lead to the significant difference in the impurity balance. For the first time for the ITER SOL, the ion distribution functions are studied using this method to examine the origin of the kinetic coefficients changes (Sec. 4). The applicability of the method is also discussed.

Grad's-Zhdanov module
In the paper [7] the reduced form of the improved analytical expressions was implemented into the SOLPS-ITER code and tested for the neon transport in ITER. Now a complete Grad's-Zhdanov module is implemented into the SOLPS-ITER code for the kinetic coefficients calculation using the 21N-moment method [3]. This module operates in two regimes: 1) Improved analytical method (IAM) [7]; 2) Explicit matrix inversion method (EMIM) [3,7]. The new module calculates thermal and friction forces coefficients between all ions, whereas the standard approach considered the thermal force only due to main ion interaction with impurities. Moreover, the heat flux according to the Grad's-Zhdanov closure is applied instead of the Braginskii ion heat conductivity [8]. Also, the ion velocity difference dependent heat flux was taken into account. Even though no significant impact from the new heat flux was observed on the solution, now all vector parallel moments are considered according to the 21N-moment Grad's-Zhdanov method. Therefore, the thermal and friction forces are treated self-consistently with the heat flux. The Grad's-Zhdanov stress-viscosity implementation is underway.
3 Impurity transport simulation results in ITER and sensitivity tests  In this section, helium and neon transport for 500MW baseline ITER scenario (Energy flux through the separatrix P SOL ≈ 100M W ) using the new Grad's-Zhdanov module is considered. The simulation has been done for the deuterium, helium and neon mixture with full drifts and currents and coupled with EIRENE. Investigation of the case, which has been presented in our previous paper [7], was continued. It is catalogued in the ITER integrated modeling analysis suite (IMAS) database as 123081.
First of all, we consider helium transport. After the implementation of the new module, the separatrix-averaged helium concentration has dropped from c He = 0.86% for the standard approach (close to results in [9]) to c He = 0.61% for the analytical (IAM) and c He = 0.57% for the numerical (EMIM) mode correspondingly. Also, one can see the difference of the helium distribution in the divertor and SOL ( figure 1). For the new approach helium is better compressed in the divertor. In figure 1 it can be seen that the helium density above the X-point is lower for this method than for the standard one. Moreover, the high helium density region is pushed closer to the target in the inner divertor (figure 1). Indeed, for the helium the new approach provides larger kinetic coefficient for the friction force and smaller kinetic coefficient for the thermal force than the standard SOLPS-ITER approach. In other words, the helium thermal diffusion (balance between friction and thermal forces), which drags helium upstream, is weaker for given parallel temperature gradient for the improved method. It can be illustrated in terms of impurity-deuterium velocity difference, which in the Scrape-off layer is mostly determined by the interplay between thermal and friction forces of impurity with deuterium [10]. Since the difference between analytical and numerical methods is not significant for helium, we will further focus on the comparison between Grad's-Zhdanov module application using analytical method and the standard SOLPS-ITER model. In figure 2 one can see the velocity difference between helium and deuterium normalised by local deuterium thermal velocity. In the regions where the deuterium flow velocity is directed towards the target and the ion temperature gradient is sufficiently high, the helium velocity is directed from the divertor towards upstream. This shifts the helium stagnation point closer to the target. Its relative position to the impurity ionisation source is one of the major effects that defines impurity leakage from the divertor [10]. The neutrals, which are ionized above the stagnation point, leak from the divertor and contribute significantly in the ion leakage flow. Indeed, for the standard approach the helium ionisation source further upstream than the stagnation point in each flux tube in the divertor region is 9.4·10 20 particles/s (total for both inner and outer targets), whereas the helium ion leakage upstream from the divertor region is 10.5 · 10 20 particles/s (the difference is due to a small radial ion flow inside the stagnation points region). After the Grad's-Zhdanov module application the helium stagnation point is shifted further upstream (figure 3), and the ionisation source further upstream than the stagnation point is 5.9 · 10 20 particles/s, whereas the helium ion leakage upstream from the divertor region is 7.2 · 10 20 particles/s. However, the significant change of the stagnation point appears only in the inner divertor, whereas a much larger contribution in the helium leakage comes from the outer divertor, where the stagnation point does not change significantly (figure 3). To give the reader an understanding why small changes in the stagnation point significantly affect particle flows, a simple model is proposed: where Γ in is the entering impurity ion poloidal flow into the divertor region through the Far-SOL; Γ out is the impurity ion poloidal leakage flow from the divertor region closer to the separatrix, and α is the fraction of entering ion poloidal flow which leaves the divertor region as ions; Γ neut is the impurity neutral leakage flow from the divertor region, which penetrates into the SOL above X-point and into the confinement region, and β is the fraction of entering ion poloidal flow, which leaves the divertor region as neutrals; Γ source is the flow due to the imprity source ( figure 4). Note, the exact choice of the divertor region boundary (red dash in the figure 3 and 4) is arbitrary, however it only slightly affects distribution between α and β. Eq. (2) is a global particle balance. As a result: Note, fractions of the entering impurity flow α and β are measured directly according to the modeling results. The ion leakage fraction is α = 0.76 and β = 0.18 for the standard approach and α = 0.74 and β = 0.16 for our improved approach. The sum α + β is equal to 0.94 (or 0.90 for the new module). It makes the denominator in (3) small and the result (3) sensitive to the denominator variation. Therefore, for chosen helium source (close to the fusion ash production) Γ source = 1.0 · 10 20 particles/s, the ion helium leakage is Γ out = 12.7 · 10 20 particles/s and Γ out = 7.4 · 10 20 particles/s for the standard and improved (IAM) approach correspondingly. It is close to the results of the simulation, which is discussed earlier. Thus, a 4% change of the distribution of the entering flow leads to a >30% decrease of the ion leakage flow. This is the result of the feedback response of the entrance ion flow, which is described by Eq (1)- (2). Indeed, the less flow leaks from the divertor, the less flow returns back into the divertor and this affects the leakage flow, as a result. Similar feedback occur in the recycling phenomena on the target. This significant reduction of the impurity leakage flow then leads to the separatrixaveraged relative helium concentration drop. Finally it is worth to mention that this change of α is a result of the stagnation point shifting. Indeed, the fraction of the entering impurity flow, which penetrates above the stagnation point as neutrals and then got ionized there, is 69%. If the stagnation point boundary is simply shifted from the standard case (red solid line in figure 3) to the case with Grad's-Zhdanov module applied (red dotted line in figure 3) without changing the ionisation source profile, the ionisation source above the new boundary becomes 1.2 · 10 20 particles/s less or 8% in terms of the entering flow fraction.
The same effect, however less pronounced, takes place for neon, as well (figure 5) (here the source is a seeding valve above the mid-plane). Before application of the Grad's-Zhdanov module the separatrix-averaged relative neon concentration is c N e = 1.0% 1 . After application of the Grad's-Zhdanov module it drops to c N e = 0.8%. Indeed, a smaller (than for the standard method) fraction of the impurity flow returns upstream and a larger fraction of the impurity flow goes into the pumping system, as it was mentioned in our earlier work [7].
It is also worth to mention that in such ITER baseline scenario simulations the velocity difference for the helium can be up to 20% of the deuterium thermal velocity (figure 2), and for the neon up to 25% of the deuterium thermal velocity ( figure 5) . This means, we are close to the limit of applicability of the linear collisional theory [3]. This will be discussed in more detail in section 4.

Ion distribution functions
In this section, ion distribution functions for deuterium, helium and neon are presented for the ITER simulation with the aim of showing why the helium kinetic coefficients are affected more than the neon ones and discussing the method applicability. The ion distribution function, which is used to obtain kinetic coefficients, can be reconstructed based on its moments (they are calculated in the code). Indeed, let us consider the Hermite polynomial approximation of the distribution function that is used in this method [3]: where moments, which we output from the code, are: w a is the diffusive velocity; h T a /h w a and r T a /r w a are the temperature gradient/velocity difference dependent parts of the heat flux and additional vector moment correspondingly (for details see [3]). And: a n a u a a m a n a , w a = u a − u Two spatial points are chosen: one in the inner divertor and one in the outer divertor (blue circles in figures 2b, 5b). In figure 6 predicted ion distribution functions at these points are shown. The parallel momentum transfer is studied, thus distribution functions are integrated over the perpendicular velocities.
One of the main improvements for the thermal and friction forces calculation is the inclusion of the correction for the impurity distribution function due to the impact from the impurity heat flux and impurity additional vector moment r, whereas for the standard SOLPS-ITER model only light main ion heat flux and additional vector moment r affect the momentum collision term. The shaping of the impurity distribution function becomes important, when the m D /m impurity ratio is not small, since the widths of distribution functions become comparable (figures 6) and impurity velocity contributes significantly to the relative velocity in the deuterium-impurity collisions. Therefore additional corrections for the helium distribution function provide a larger effect than for the neon distribution function.
It is interesting to mention that, for the deuterium, the main contribution into the shaping of the distribution function comes from the ion temperature gradient (figure 6), which leads to the so-called "tail": overpopulated (with respect to the non-modified Maxwell function) high velocity phase space and underpopulated low velocity phase space down the temperature gradient and opposite up the temperature gradient. In figure 6b deuterium parallel flow velocity is positive (small), however the tail leads to a more effective momentum transfer to the impurity in the negative direction, due to a velocity dependence of the Coulomb collision cross-section. This phenomenon is the thermal force. It leads to the shift of the impurity distribution function. This shift leads to the friction force, which balances the thermal force. As a result, impurity parallel flow velocity is negative. This is the nature of the impurity thermal diffusion, such diffusion nature for electrons is described in [8].
In figure 6b one can see a negative part of the neon distribution. Such negative parts of distribution functions are present in this method or in [8], however inside the area of theory applicability they are small and do not contribute to the result. In some places, where the ion temperature gradient is sufficiently strong (figure 6b), the neon diffusive velocity w N e is ≈60% of the neon thermal velocity T i /m N e . In such a case we are on the boundary of the applicability of linear Grad's-Zhdanov theory. Moreover, in the ITER simulations such distribution function approximation can be not sufficiently accurate, when it is applied for the ions heavier than neon, since their thermal velocity is even smaller. The applicability of this method for heavy impurities requires further investigations.

Conclusions
Improved analytical parallel kinetic coefficients are used for the helium and neon transport simulations for the ITER baseline scenario. New thermal and friction forces kinetic coefficients lead to weaker thermal diffusion, which drags impurity upstream. Significant impact on the helium transport is observed compared to the standard SOLPS-ITER approach. Relative helium concentration (separatrix-averaged) decreases by 30%. For neon, changes are less pronounced. A simple leakage model is presented for illustration of the parallel impurity transport impact on the global impurity balance. Using linear Grad's-Zhdanov method ion distribution functions are analyzed for ITER Scrape-off layer predictive simulation for the first time. For helium this method is well applicable, whereas for neon ITER cases are on the boundary of applicability.