Self-consistent, global, neoclassical radial-electric-field calculations of electron-ion-root transitions in the W7-X stellarator

The radial electric field in the Wendelstein 7-X stellarator is computed by means of self-consistent, global, neoclassical simulations using the gyrokinetic particle-in-cell code EUTERPE. The simulation results are compared with local predictions obtained from a transport code using locally computed neoclassical transport coefficients. The analysis focuses on ion-electron-root transitions and investigates their dependence on collisionality, normalised ion gyroradius, and the electron-ion temperature ratio. Several of the results cannot be reproduced using conventional, local neoclassical transport theory. An approximate criterion for root transitions is derived, which results in an analytical scaling law that is useful for understanding how the position of the transition layer varies with plasma parameters.


Introduction
Neoclassical transport processes play an important role in stellarators for several different reasons.At high enough temperatures, neoclassical transport fluxes become very large and can exceed the turbulent ones.For instance, in the 1/ν-regime, which usually pertains to the electrons, the particle and energy transport coefficients scale with the temperature T as T 7/2 , which is stronger than the gyro-Bohm scaling (T 3/2 ) typical Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. of turbulent transport.The ion transport is reduced by the appearance of a radial electric field, which limits the radial excursions of trapped particles through poloidal E × B rotation.This electric field is determined by the requirement that the particle transport be ambipolar.It is important to note that, even if the turbulent transport fluxes exceed the neoclassical ones, the radial electric field is still expected to be determined by neoclassical transport.The reason is that in stellarator configurations, to leading order in the normalised gyroradius ρ * = ρ/L ≪ 1 (ρ is the gyroradius and L is the scale length of the plasma profiles), the latter depends on the radial electric field whereas turbulent transport does not.Neoclassical transport is intrinsically ambipolar only in axisymmetric or quasisymmetric magnetic configurations, for which other transport processes determine the radial electric field [1].Therefore, the radial electric field and, hence, poloidal plasma rotation is expected to be determined by neoclassical rather than turbulent transport in most stellarators, at least on length scales significantly exceeding the ion gyroradius [2].
Standard neoclassical theory is local in minor radius, i.e. it makes predictions on the basis of the local density, temperature, and their gradients.Many of these predictions have been confirmed experimentally, including the bootstrap current and the radial electric field [3,4].Sometimes terms that are customarily neglected need to be retained in the local theory.For instance, the observation of an impurity hole in the Large Helical Device (LHD) has led to the recognition of the importance of electric potential variation on flux surfaces in impurity transport studies [5].Moreover, in certain limits, for omnigenous or large aspect ratio stellarators, a local theory can be applied that includes the drifts tangential to the flux surface [6].These become important for very low collision frequencies (below the 1/ν-regime) [7].
However, there are also phenomena where a radially local treatment of neoclassical transport breaks down.A particularly important example, which is the topic of this article, occurs when the electric field varies rapidly with the plasma radius.The equation for neoclassical-transport ambipolarity is nonlinear and can have several roots [8][9][10].The so-called ionroot is usually realised, corresponding to a negative (inwardpointing) radial electric field, E r < 0, but if the electrons are much hotter than the ions, the electron-root, E r > 0, can arise.This typically occurs in the centre of the plasma, and, as a result, the electric field then undergoes a sharp transition at some intermediate radius.This transition cannot be described by local neoclassical theory but can be simulated by selfconsistent, global, particle-in-cell-code simulations.Such simulations are the subject of the present paper.
Hitherto, global radial electric field calculations were undertaken in [11] to study impurity hole plasmas, but therein only the ions were treated kinetically.This methodology was adopted even though the utilised code, FORTEC3D [12], is capable of treating electrons kinetically [13], primarily to save on computation time.
An earlier attempt to simulate the electron-ion-root transition self-consistently was undertaken in [14] using the code GTC [15], wherein buffer zones were introduced in the field equation.Although the simulations are global within the scope of their computational domain, the presence of these buffer zones at the boundaries could potentially influence the solution.In addition to calculating the radial electric field, turbulence simulations were carried out in the presence of this field in order to study whether the sheared plasma rotation suppresses turbulent transport.It has been hypothesised that large values of E ′ r in the transition layer strongly shear the plasma flow and thus reduce the turbulent transport [16].
In the work presented here, we treat both the electrons and the ions kinetically, and carefully calculate the radial electric field self-consistently in the entire plasma volume.We compare the global electric field, neoclassical particle-and energy fluxes with their counterparts computed using local theory, and we explore how the results depend on the normalised gyroradius, collisionality, and electron-to-ion temperature ratio.
This article is structured as follows: in section 2 we discuss the details of the numerical implementation in EUTERPE, show a comparison between the local and global neoclassical orderings, and provide an analytical insight into the problem based on the local theory.In section 3 we give a brief description of the types of magnetic configurations of W7-X employed in this study.In the following section 4, we provide a discussion focusing on the benchmark of the results of global EUTERPE simulations against results from Neotransp, a transport code based on local neoclassical theory, which performs the integration of the mono-energetic transport coefficients in the same way as the NTSS code [17].Finally, in section 5 we analyse the dependence of the transition layer on the ion collisionality ν * i , the electron-ion temperature ratio τ , and the normalised ion gyroradius ρ * i .We focus on phenomena that cannot be explained using the conventional, local theory which in turn highlights the importance of self-consistent global electric field simulations.

Methodology and numerical techniques used in EUTERPE
EUTERPE [18] is a gyrokinetic Monte Carlo particle-incell δf code capable of treating any ideal magnetohydrodynamic (MHD) equilibrium calculated with the VMEC code [19].Since intra-and inter-species collisions are implemented in EUTERPE, it is capable of simulating neoclassical transport processes.In the past, the electrostatic potential and its first-order perturbation (small variation within the flux surface) have been computed with EUTERPE in the local approximation [20].In the present work, we extend the code to solve the self-consistent, neoclassical equations globally.We consider the following equations of motion for each species, where we display the species index only if necessary: where R is the gyrocentre position, v ∥ = b • Ṙ the velocity parallel to the field lines, and µ = v 2 ⊥ 2B the specific magnetic moment.Furthermore, m is the species mass, q the charge, and The code solves the drift kinetic equation for each species, where is the Lorentz operator for pitch-angle scattering with ζ = v ∥ /v and ν the species collision frequency, which is a function of the radial coordinate r = a √ s (where a is the minor radius and s ∈ (0, 1) the normalised toroidal flux) and the velocity.The collision frequency of species α colliding against other species, which are indexed by β (α = β for self-collisions) is denoted by: where: Here, we defined a reference collision frequency with K α/β = (m β /m α )(T α /T β )K, K = mv 2 /2T the normalised kinetic energy, n = ´f d 3 v the species' density, log(Λ) the Coulomb logarithm, ϵ 0 the permittivity of free space, and the error function.For electrons, ν e is the sum of electronelectron and electron-ion collision frequencies.On the other hand, for ions, ν i includes only ion-ion collisions, since ionelectron collisions effectively act as a relatively small friction force, which can be neglected.A detailed description of the implementations of collisions in EUTERPE is given in [21].
The omission of the energy scattering term in equation ( 5) can be motivated by the following observations.First, from the viewpoint of the local theory, in the long-mean-free-path regime, neoclassical transport in stellarators is dominated by contributions from 'locally trapped' particles that have their trapping state changed (or modified) by pitch-angle scattering.Since energy scattering plays no role in determining whether a particle is locally trapped or 'locally passing', it can therefore be neglected.In the global picture though, the simulation results could be influenced if the redistribution of markers during the simulation leads to an energy distribution that is considerably different from what is expected from the initial density and temperature profiles.However, we observed that the relative changes in density and temperature profiles which occur during the simulation are small.This shows that confinement for W7-X is sufficiently good so that energy scattering would have negligible influence on the results.
The electric potential Φ is determined by the drift-kinetic limit of the gyrokinetic quasineutrality equation: where n 0,i is the unperturbed ion density profile.On length scales L exceeding the ion gyroradius ρ i , this equation implies quasineutrality, since the variations in the electrostatic potential are expected to be of order ϕ ∼ T/q.Equation ( 10) has already been implemented in EUTERPE and is more natural than using a timedependent field equation presented in [14], since the boundary conditions are unambiguous, removing the need for introducing buffer zones.In addition, it allows simulating poloidal and toroidal Fourier modes of the electric field.For simplicity, we only resolve Φ 0 ≡ Φ (0,0) -the (m, n) = (0, 0) Fourier component of the electric potential Φ, i.e. the higher order variation of the electric potential Φ 1 is neglected.This is done by applying a corresponding Fourier filter on Φ after each time step.The boundary conditions are Φ = 0 at the plasma boundary and the natural regularity condition at the magnetic axis, that is Φ ′ (0,n) = 0, and Φ (m,n) = 0, for m ̸ = 0 [18].EUTERPE uses a δf-method to reduce particle noise, where equations ( 1)- (10) are solved for each species by introducing δf with where v T = T/m is the species' thermal velocity.Starting with δf = 0 the thermodynamic forces drive the system which, after a quick initial transient phase, converges to a steady-state solution after a few collision times.Once the steady-state solution is reached, the species' total particle flux F Γ through the surface S is calculated using: where r is the radial coordinate, ⟨•⟩ is the flux surface average, and γ = ´δf v d 3 v is the particle flux density.For the total energy flux, F Q = Q ´S √ g dθ dϕ, where Q = ⟨q • ∇r⟩, and q = ´m 2 v 2 δf v d 3 v is the particle energy flux density.Similarly, the bootstrap current density is calculated using: In addition to the electric field profile, these three quantities will be used to benchmark EUTERPE results with local theory using the codes DKES [22,23] and Neotransp.

Predictions of local neoclassical transport theory
According to the ordering scheme used in the local neoclassical theory [24][25][26], the deviation from a local Maxwellian is small where Moreover, δf is assumed to vary on the same length scale as f 0 , so that as opposed to the global theory, where smaller-scale variations are allowed Such localised, narrow variations are indeed expected in root transitions, necessitating a global treatment.With these ordering differences in mind, we can summarise both approaches in a single equation: with the corresponding definitions in table 1.The consequence of neglecting v B in Ṙ0 is that local theory effectively decouples the flux surfaces.Note that local simulations are linear whereas global simulations must be nonlinear because the backreaction of the trajectories on the self-consistent field is essential for establishing ambipolarity.
In either case, only the radial part of the electric field, E r = −dΦ 0 /dr, which causes poloidal rotation, is included in v Φ .In local simulations each flux surface can thus be simulated separately since no radial derivatives act on δf in equation (18).The radial electric field enters as a parameter in this equation, and regulates the particle flux.By performing multiple simulations with different E r , the field E r that satisfies the ambipolarity constraint can be determined at each radius: where Γ i and Γ e now refer to the ion and electron particle fluxes in a plasma with a single ion species with charge number Z.This local approximation allows for an analytical treatment providing valuable insight into the structure of the electric field.Here, we describe a simplified model that gives rise to the electric field root transition phenomenon.
For given electric field, equation ( 18) is linear in δf and the right hand side is a linear combination of the so-called thermodynamic forces, i.e. the radial density and temperature gradients, and the electric field.As a result, all moments of δf, including the particle and energy fluxes and the bootstrap current, are linearly dependent on the thermodynamic  (18).The methods of calculating the radial electric field are also compared between the local and global theories.While in general, Φ = Φ 0 + Φ 1 , where Φ 1 are higher order corrections to the flux-surface-averaged potential Φ 0 , in this work we consider Φ = Φ 0 .

Transport Er calculation
forces driving them.When calculating the particle flux, particularly in stellarators, the parallel electric field E ∥ can usually be neglected, i.e. the Ware pinch is negligible.Furthermore, since the collision operator is approximated by pure pitchangle scattering in equation ( 18), it follows that the particle flux is fully determined by the so-called D 11 mono-energetic transport coefficient [25].In particular, the following relations hold: where primes represent differentiation with respect to r, and Traditionally, different regimes of transport are distinguished based on their collision frequency dependence and analytical estimates for D 11 are available in each case.For our consideration, the most important regimes are: 1.The 1/ν-regime: 2. The √ ν-regime: where ω E = E r /rB is the characteristic E × B precession frequency, ϵ eff is the effective ripple, b 10 characterises the principal variation of B in the poloidal direction (refer to equation ( 29)), ϵ t = r R with R the major radius, and v d = mv 2 /2qBR is the characteristic drift velocity.For these 'pure' regimes, D For instance, if at a certain radius in the plasma, the ions are in the √ ν-regime and the electrons are in the 1/νregime, it is possible to obtain multiple solutions for E r that satisfy the ambipolarity constraint (equation ( 19)).This result is discussed in [8].This situation arises in a transition layer, where depending on the sign of the electric field (i.e. the root of the ambipolarity condition), we distinguish between the ion-root solution, satisfying E r < 0 and the electron-root solution E r > 0, along with a third intermediate solution, which is always unstable, as follows from the analysis of the linear stability of these solutions [10].
We now derive an approximate criterion for the existence of the electron-root based on the local Ansatz in terms of three important dimensionless parameters: the collisionality, the normalised ion gyroradius, and the electron-ion temperature ratio: where the subscript 'th' indicates that the corresponding quantities are evaluated at K = 1, and Ω i = ZeB/m i is the ion gyrofrequency.
Typically, the plasma profiles in W7-X satisfy n ′ /n < 0 and T ′ /T < 0. It then follows from equation (19), that L 11,e > L 11,i is a necessary condition for the electron root to exist, and conversely L 11,e < L 11,i is a necessary condition for the ion root existence.It is thus instructive to analyse the ratio L 11,e /L 11,i .In the 1/ν-regime, it becomes L i/e τ 7/2 , where µ i/e is the ion to electron mass ratio, which shows that for an equithermal plasma the ions diffuse faster than electrons.This would contradict ambipolarity, and a radial electric field therefore arises that causes the ions to enter the √ ν-regime whilst the electrons stay in the 1/ν-regime.It follows from equations ( 21)-( 23) that the ratio L 11,i is given by: On the other hand, we expect the radial electric field to be of the order of E r ∼ T i /(ea), where a denotes the minor radius of the plasma, so that Therefore, An electron-root is expected if the electrons are so hot that this ratio exceeds unity, which happens for τ greater than a critical value τ * : Since τ is usually a monotonically decreasing function of r, we expect that any electron-root first appears close to the magnetic axis, and that its extent increases with τ , ρ * i and ϵ eff , and decreases with ν * i and b 10 .This expected behaviour will be compared with the simulation results below.For this comparison we will solve the equation τ = γτ * with a proportionality constant γ, for a radial position r * .We emphasise that ν * i , ρ * i , τ , b 10 and ϵ eff all depend on r.
Finally, we note that in global simulations, one cannot decouple the transport properties between neighbouring flux surfaces.It is then necessary to calculate the self-consistent radial electric field with the gyrokinetic Poisson's equation (equation ( 19)) or by further combining it with the guiding centre continuity equation [14], since the ambipolarity constraint relies on the local Ansatz.While EUTERPE is in principle able to run local simulations, we have decided to instead benchmark the global electric field against Neotransp.This was done in order to exclude the possibility of any persistent errors within the local and global simulations in EUTERPE.

Description of the W7-X stellarator and its confinement properties
W7-X is often described as an 'optimised stellarator'.This term pertains to its design era, reflecting its adherence to multiple optimisation criteria for effective plasma confinement.Crucially for our simulations, it is characterised by favourable neoclassical confinement properties [27].In our simulations, we have observed that the self-consistent, global radial electric field converges on the timescale of approximately the collision time of the slowest colliding ions τ νi .If during this time a significant portion of the particles escapes from the confinement region, then naturally, the choice of particle boundary conditions at the plasma edge (particle reinsertion) will affect the solution.However, in our analysis, we have found that the choice of particle boundary conditions affects the solution only in a small layer close to the plasma edge for simulations in W7-X.This is expected, since even in the absence of the electric field, the particle orbits are sufficiently well confined on the timescale τ νi .On the other hand, in the case of classical or poorly optimised stellarators, particle reinsertion schemes may affect the steady-state solution, and hence require careful treatment.

Overview of equilibrium configurations and plasma parameters
The non-planar, modular coil system of W7-X allows for great flexibility, providing the possibility of comparing different magnetic configurations [28].In this work, we consider the following reference configurations of W7-X [25]: • The low-mirror configuration (b 01 = 0 on axis), which tries to minimise the variation of the field strength along the magnetic axis.
• The standard configuration (b 01 = 0.046 on axis), characterised by equal currents in all non-planar coils.• The high-mirror configuration (b 01 = 0.1 on axis), characterised by improved confinement of the deeply trapped particles at volume averaged β > 2%.This is the configuration closest to the target configuration of the optimisation which led to W7-X.
Here b 01 characterises the principal variation of B in the toroidal direction, as obtained from the decomposition of a VMEC [19,29] equilibrium into Fourier modes in Boozer coordinates: where b m,n = B m,n /B 0,0 , and W7-X has N = 5 field periods.
For the sake of simplicity, we are exclusively examining the vacuum fields, i.e. the finite-β effects are neglected.We aim to demonstrate that EUTERPE is able to calculate the selfconsistent global electric field for different MHD equilibria, and that this code can adequately resolve the radial structure of the electron-to-ion-root transition layer.

Benchmarking EUTERPE and Neotransp simulations
In our investigation, we use density and temperature profiles, n exp and T exp , of a pure hydrogen plasma Z = 1 obtained by interpolating experimental results from [3], shown in figure 1(a).The simulation time is equal to one ion collision time of the ions in the core of the plasma volume, which is approximately 0.005 s.The time-step is set to a fixed value of ∆t = 0.5/Ω i , and the particle trajectories are calculated using the fourth order-Runge Kutta integrator.For collisions, being modelled as a stochastic process [21], the Euler-Maruyama method is utilised.The total number of macro-particles is N p,e = 10 9 for electrons and N p,i = 10 8 for ions.The requirement for such a large quantity of macro-particles may seem surprising.However, it is worth noting that there have been previous instances of global neoclassical electron transport simulations employing a similar magnitude of macro-particles [13].Electrons, being lighter, exhibit a quicker response to the electric field, which in turn induces fluctuations in the electric field.About 10 6 CPU hours were required to complete the simulation on a high-performance computer.The spatial grid was divided into N r = 32 radial, N θ = 16 poloidal, and N ϕ = 16 toroidal points.EUTERPE simulates only one field period of W7-X, employing periodic boundary conditions in the toroidal direction.
Even for these high macro-particle numbers, the electrostatic potential shows strong noise-induced fluctuations.It was thus necessary to average the simulation results over a time window of length 0.5τ νi .Near the magnetic axis, the Jacobian of the coordinate systems that is used for charge assignment in EUTERPE goes to zero, inducing large fluctuations there.The work-around used in our simulations was to significantly increase the macro-particle number near r = 0.As a result, the noise-induced fluctuations in the particle-and energy fluxes could be reduced.Consequently, convergence of the timeaveraged radial electric field to a stationary state was employed as a criterion for determining the reliability of simulations.
In figure 1(b) we show the radial electric field obtained with EUTERPE for the standard configuration of W7-X.Also shown in the figure are the ion-and electron-root values of the radial electric field corresponding to locally ambipolar neoclassical particle fluxes computed by the local code Neotransp.At small radii, r < 0.12 m, there is only an electron-root (E r > 0), and at large radii, r > 0.2 m, only an ion-root (E r < 0).At intermediate radii both roots are present, and the EUTERPE solution undergoes a transition from one to the other.Outside the transition layer, the electric field from EUTERPE is in good agreement with the prediction of Neotransp.This is in accordance with theoretical expectations for a well-optimised stellarator with small orbit widths and few unconfined trajectories.We observe that the inflection point and the zero of the global electric field lie within the region of the overlap of the local roots.Finally, we notice that around the inflection point, the derivative of the global field has a symmetric, parabolic profile.Throughout this article, the term 'inflection point', denoted as r i , refers to the global minimum of E ′ r .Figure 2 shows the bootstrap currents, neoclassical energy flux and neoclassical particle flux computed by EUTERPE.These are compared with the predictions of the local theory.The electron-related quantities generally exhibit high fluctuations but their averages align well with the local predictions.We note the following results: • The total particle fluxes are ambipolar.
• The ion-root solution of the local neoclassical theory significantly overestimates the particle and energy fluxes in the transition layer.This is especially visible for the ions.We observe that the global solution is closer to the value predicted by the local electron-root in the transition layer.• The electron bootstrap current agrees with the local solution, but the particle noise makes it difficult to assess its behaviour in the transition layer.On the other hand, the ion bootstrap current has not reached convergence within the timescale of this simulation.To accurately determine its profile, a longer simulation time would be necessary.We note that both the local and global results do not include the momentum correction term in the collision operator.
In order to save computing time, for the remaining simulations presented in this work, we decided to reduce the number of markers to N p,e = 10 8 and N p,i = 10 7 respectively.Figure 3(a) shows the resulting electric field in the standard configuration.It is in qualitative agreement with the field shown in figure 1(b).The most significant difference is seen in the ion-root regime, however this region is not of interest, since the local theory accurately describes the electric field there.We have decided that this level of detail is sufficient for our analysis.A further, quantitative analysis of the effect of reducing the number of markers in the standard configuration is discussed in section 4.2 and summarised in table 3.In order to demonstrate that a variety of magnetic configurations can be simulated, we also present results from simulations in the low-and high-mirror configurations with the same density and temperature profiles as before.Figures 3(b) and (c) show the results of these simulations.In both cases, the electric field is in good agreement with the local theory in the electron-root region.Similarly to the standard configuration results, a small deviation from the local results is seen in the ion-root regime.However, this discrepancy is expected to disappear with increased macro-particle number.The inflection point and zero of the global electric field reside within the intersection region of the local roots.In the vicinity of the inflection point, the derivative of the global field exhibits a parabolic behaviour.For the high-mirror case, the transition region is significantly shifted outwards, and the electric field possesses a much stronger electron-root.This is unlike the low-mirror case, where the electric field profile is comparable with the standard configuration.Table 2 summarises this information at certain radii of interest.

Methodology for the analysis of electron-ion-root transitions
We now explore the transition region in greater detail and introduce several quantities in order to characterise its properties.This is necessary in order to study the variations of the transition layer width and its position with the dimensionless plasma parameters in section 5.
We begin by defining three radii, r 0 , r max and r i , by the relations E r (r 0 ) = 0, E ′ r (r max ) = 0 and E ′ ′ r (r i ) = 0, respectively.If Comparison of the location of the maximum of the electric field rmax, its inflection point r i , its zero r 0 and the scaling law prediction r * for different magnetic configurations (all in cm).Note that the proportionality constant γ of the scaling law was chosen so that r 0 = r * for the experimental profiles in the W7-X standard configuration.The vertical lines from left to right mark the locations of rmax, r i , and r 0 , accordingly.The orange, solid line denotes the best parabola fit at the minimum of E ′ r , i.e. at the point r i .For this example, the zeros of the parabola appear in the vicinity of rmax and r 0 .This observation lead to the definition of the submaximal transitions.the electron-root is realised in the plasma centre, E(r) is positive for small r and necessarily reaches a maximum before changing sign, so we must have r max < r 0 .Empirically, we also observe that usually r max < r i < r 0 .Figure 1(b) illustrates the three radial locations r max , r i and r 0 .In order to characterise the width of the transition region, we can use either of the quantities 2(r i − r max ), 2(r 0 − r i ) or r 0 − r max .A fourth possibility, shown in figure 4 is to define the width w as the distance between the two zero-crossings of the best parabola fit of E ′ r around r i .In the case of figure 4 all these measures of the transition width are comparable to each other.As can be seen in table 3, they range between 4.3 cm and 4.8 cm in the best resolved simulations.This table summarises the variation of the described quantities with macro-particle number for the simulations in the standard configuration with N p,e = 10 9 , N p,i = 10 8 and N p,e = 10 8 , N p,i = 10 7 , correspondingly.As we can see, r max r i and r 0 vary only by about 1 cm, and since these quantities are of primary interest in our research, in the next section we perform a parameter scan with a reduced number of markers.However, in general, we expect the quantities related to the width of the transition zone to have a larger variation, since they rely on first order differences of the radial positions.
Another quantity of interest given in table 3 is , where v th = 2T 0,i /m i is the ion thermal speed.The quantity E ′ r * measures the shear of the E × B flow and will be considered in the next section.
We also find a special type of transitions in which r max < r max,e , where r max,e is the location of maximum E r,e , the local electron-root solution.We refer to such cases as submaximal transitions.In this case, a clear indication of the beginning of the transition is visible, since at a certain radius sign(E ′ r ) ̸ = sign(E ′ r,e ) causing a significant deviation between the two curves, E r and E r,e .We observe that this approximately happens at r max .However, if the transition is not submaximal, that is r max = r max,e , then the deviation from the local electron-root begins at r > r max and can develop smoothly, making it difficult to define the exact location of the beginning of the transition.A similar difficulty always arises at the outer end of the transition layer, where the solution gradually merges with the local ion-root, irrespective of whether the transition is submaximal.To identify the end of the transition, one can instead consider the point where the electric field changes sign, i.e. the point r 0 .Finally, the aforementioned parabolic profile of E ′ r around r i seems to be a universal property of the root-transitions, which has been observed in every simulation we performed.Therefore, r i lies approximately in the middle of the transition.This leads us to the following criterion for identifying submaximal transitions: r 0 − r max ∼ 2(r i − r max ) ∼ 2(r 0 − r i ).For non-submaximal transitions, the transition begins at r > r max , and only r i and r 0 characterise the transition region.In such instances 2(r 0 − r i ) and w are the most reasonable definitions for the width of the transition.

Dependence on dimensionless plasma parameters
Our final discussion relates to the dependence of the electronion-root transition layer on the dimensionless plasma parameters introduced in equation ( 24).The findings in this section pertain to the W7-X standard configuration.Validation of the claims regarding different magnetic geometries is left for future research.First, the ion collision frequency can be approximated as: where for simplicity we neglect any variation of the Coulomb logarithm (i.e.log(Λ) ≈ const.).Thus, the considered dimensionless plasma parameters can be written explicitly in terms of the T i , T e and n profiles: We explore the solution space by introducing a parameter α and altering the experimental profiles, n exp , T exp,e , and T exp,i , in such a way that two dimensionless parameters are held fixed and the remaining one varies linearly with α.We choose the following dependencies: .75, 0.875, 1, 1.125, 1.25, 1.5} and: • For ν * i , α = α ν * ∈ {0.5, 1, 1.5, 2, 2.5, 3} and: • For τ , α = α τ ∈ {0, 0.5, 1, 1.5, 2, 2.5} and: Note that this parameterisation ensures T i = T e at the edge independently of α τ .
As mentioned previously, to compare the simulation results with the scaling law provided in equation ( 28), we solve the equation τ = γτ * for the radial position r * .The constant of proportionality γ, which is a free parameter, is chosen so that r * = r 0 for the profiles with α = 1, i.e. for T i = T exp,i , T e = T exp,e , and n = n exp .Figure 5 shows the plots of ρ * i , ν * i , τ , and τ * for this case.We proceed with investigating the dependence of the transition location, its width, and the maximum shear of the E × B flow on α ρ * , α ν * , and α τ .
Figure 6 shows how the location of the transition layer depends on the considered plasma parameters in the W7-X standard configuration.Figure 6(b) indicates that the transition layer is shifted inwards as ν * i is increased, and figures 6(a) and (c) show that the transition layer moves towards the plasma edge as ρ * i and/or τ increases.The analytical scaling law is in excellent agreement with the simulation results for the ν * i and ρ * i parameterisations.For the τ parameterisation, the two points in figure 6(c) located at α τ = 0 and α τ = 0.5 correspond to ion-root plasmas, even though r * indicates a transition.This illustrates that, when varying τ , the prediction of r * equation ( 28) becomes unreliable in weak electron-root plasmas.Finally, in all cases explored, the scaling law (28) agrees well with the r 0 dependence on the magnetic geometry, with ϵ eff being the quantity that varies most significantly between the reference configurations.This result is summarised in table 2.
It is important to mention that, in the limit of ρ * i → 0, one expects r 0 = r i from the local theory, but it is computationally expensive to simulate such low ρ * i plasmas in EUTERPE and this verification is therefore left for future research.
The level of plasma turbulence and the transport it causes is known to be reduced by sheared E × B flows.It is thus of interest to find magnetic configurations that have a strongly sheared radial electric field over a significant portion of the plasma volume.Figure 7 shows how the transition layer width depends on the considered plasma parameters.For a given value of α, a large discrepancy between 2(r 0 − r i ) and either of the r 0 − r max and 2(r i − r max ) suggest that a transition is not submaximal.This is the case for α ν * ⩾ 2 in the ν * i parameterisation, for α τ ⩾ 1.5 in the τ parameterisation, and for α ρ * ⩾ 1.125 in the ρ * i parameterisation.On the other hand, for submaximal transitions, the four employed width measures agree within a difference of order 1cm.Overall, the transition zone width increases with ν * i , τ and ρ * i with τ having the most significant impact.However, it appears that there exists a certain τ beyond which the transition zone width remains constant.Regarding the shear of the E × B flow, in figure 8    electric field was sheared enough to cause a significant turbulence suppression in these scenarios remains to be investigated in the future.
Finally, we have found that there exist scenarios when the global simulations are essential for predicting the structure of the electric field.Figure 9(a) shows that in the case of a weak electron-root, a global simulation is essential to determine whether a transition layer exists at all.It arises in the τ scan with α τ = 0.5.The ion-root exists in almost the entire plasma and, according to our results, the electric field remains in the ion-root everywhere.In other words, although an electron root is possible according to the local theory, it does not materialise in the global simulation.Furthermore, figure 9(b) shows that it is possible that the zero crossing of the electric field in the global simulation can appear outside of the local electron-ion root overlap region.This appears for α ρ * = 1.25 and α ρ * = 1.5 in the ρ * i scan.As mentioned earlier, the above results were obtained from simulations in the W7-X standard configuration.While it is important to approach these findings with care when applying them to different magnetic configurations, the local theory suggests that the location of the electric field transition should still follow the scaling law presented in equation (28).Additionally, for W7-X, table 4 indicates a substantial variation in E ′ r * across different configurations, suggesting a positive correlation with ϵ eff , in contrast to the transition layer width, which varies relatively insignificantly.Our final comment in this section pertains to the quantitative analysis of the ν * i and ρ * i exponents appearing in the critical temperature ratio τ * , as seen in equation (28).From our fits, we measured that the ρ * i exponent for global simulations is −0.26, differing from the predicted −3/7.Similarly, the ν * i exponent is measured to be 0.18, as opposed to the expected  3/7.This suggests that while equation (28) may not be quantitatively precise, it still provides useful qualitative insights into the behavior of the global, radial electric field.

Conclusions
We have presented self-consistent, global neoclassical radialelectric-field simulations with EUTERPE for the optimised stellarator W7-X.The simulations were performed in the standard, high-mirror and low-mirror configurations of this device.Outside the transition region between the electron and ion roots, the radial electric field and the neoclassical transport fluxes are in excellent agreement with output from the local code Neotransp.However, global simulations are necessary to resolve the transition region and to determine its radial location.An example was found where the transition extends beyond the overlap region of the local electronand ion-root solutions, which cannot be understood from conventional local theory.In some cases, global simulations are necessary even to establish the very existence of a transition.We observe that, in our simulations, the time-averaged electric field smoothly converges to a stationary steady-state, effectively resolving the multi-valuedness problem of the local neoclassical theory.A universal property of the global electric field observed in our simulations is that its radial derivative has an approximately parabolic profile extending from the inflection point to the zero of the electric field, making it possible to define the width of the transition region as the distance between the roots of the parabolic fit of E ′ r at its minimum.
We have also explored how the location of the transition layer, its radial structure, and the maximum value of the electric field gradient depend on the electron-ion temperature ratio, the normalised ion gyroradius and the ion collisionality.Most of the results can be understood in terms of a simple scaling law (equation ( 28)) derived under the assumption that the ion root is realised when the ions are in the √ ν-regime and the electrons are in the 1/ν-regime.The resulting criterion correctly predicts the variation of the location of the transition region as well as the magnetic geometry according to the dimensionless plasma parameters .
Future areas of research include the simulation of microinstabilities in plasmas with a root transition and a selfconsistent neoclassical electric field profile, and verifying the behaviour of the transition layer in the low ion gyroradius limit.An unambiguous comparison between the theory and experimental measurements of the radial electric field in W7-X is underway.
The index i is equal to either 1 or 2, and the weighting factors are h 1 = 1 and h 2 = K.If D 11 ∼ (KT) ζ then the integrals in equation (21) imply δ 12 = ζ, so that δ 12 characterises the temperature dependence of L 11 .

Figure 1 .
Figure 1.Benchmarking of the global, neoclassical radial electric field calculated with EUTERPE against the local neoclassical code Neotransp for the W7-X standard configuration.(a) the input profiles for the simulation (note that the electron and ion densities are equal).(b) the resulting electric field and its radial derivative: the black line shows the time-averaged EUTERPE solution with a window size of 0.5τ i .The red and blue lines are the electron-and ion-root solutions from Neotransp.The vertical lines from left to right mark the locations of the maximum, inflection point and the zero of the global electric field, and are denoted further as rmax, r i , and r 0 , accordingly.

Figure 2 .
Figure 2. Comparison of the ion (left) and electron (right) bootstrap current (top), energy flux (middle), and particle flux (bottom) for the W7-X standard configuration.The black lines represent the time-averaged quantities obtained from EUTERPE.Orange, dashed lines show the 1 − σ error bars, and the blue and red lines are the electron-and ion-root solutions from Neotransp.

Figure 3 .
Figure 3. Simulations of the radial electric field (with a reduced number of markers) for the standard, low-mirror and high-mirror configurations of W7-X.Top plots show the electric field, and the bottom ones its derivative.Vertical lines indicate the maximum of the electric field, its inflection point, and zero, accordingly.

Figure 4 .
Figure 4. Parabolic fit at the minimum of E ′ r for the example global field shown in figure 1(b).The vertical lines from left to right mark the locations of rmax, r i , and r 0 , accordingly.The orange, solid line denotes the best parabola fit at the minimum of E ′ r , i.e. at the point r i .For this example, the zeros of the parabola appear in the vicinity of rmax and r 0 .This observation lead to the definition of the submaximal transitions.

Figure 5 .
Figure 5. Dimensionless plasma parameters for the experimental profiles used from the W7-X standard configuration.(a) the radial dependence of the electron-ion temperature ratio τ and the critical ratio τ * calculated from equation (28).The intersection of these lines provides the radial position r * .(b) the radial dependence of the normalised ion gyroradius and the ion collisionality for the same profiles.

Figure 6 .
Figure 6.Dependence of various radial locations related to the electric field on the considered dimensionless plasma parameters for W7-X standard configuration.The blue lines (top) show the r 0 location, where Er = 0.The green lines (bottom) show the rmax location, i.e. the maximum of Er.The orange lines show the r i location, that is the inflection point of Er.The red lines indicate r * , the solution to τ = γτ * .

Figure 7 .
Figure 7. Dependence of various definitions of the transition layer width on the considered dimensionless plasma parameters in W7-X standard configuration.The blue lines (labelled w) are calculated from the distance between the roots of the parabolic fit of E ′ r at the inflection point of Er.The green lines (labelled r 0 − rmax) are the distance between the zero and the maximum of the electric field.The orange and red lines are to be interpreted similarly.

Figure 8 .
Figure 8. Dependence of the minimum of E ′ r * in the transition layer on the considered dimnsionless plasma parameters in W7-X standard configuration.Note that the vertical axes have different ranges, and that no root transition appears for ατ ⩽ 0.5 and therefore E ′ r * is not plotted for these cases.

Figure 9 .
Figure 9. Examples of electric field transitions in W7-X standard configuration that cannot be explained by local neoclassical transport theory.(a) a situation (with ατ = 0.5) where no electron-root is present in the global solution even though it appears in local calculation.(b) it is possible that the zero-crossing of the electric field in the global simulation can appear outside of the local electron-ion-root overlap region.This plot corresponds to αρ * = 1.5, but similar behaviour can also be seen for αρ * = 1.25 (not shown).

Table 4 .
Comparison of several quantities of interest for different magnetic configurations of W7-X.All the quantities are given in centimetres, except for E ′ r * , which is measured in m −1 .Configuration 2(r i − rmax) 2(r 0 − r i ) r 0 − rmax w E

Table 1 .
Comparison between the local and global formalism.The table lists the terms entering the particles' equations of motion ( Ṙ0 ∥ ), and the source terms ( Ṙ1 and v1 ∥ ), appearing on the RHS of equation

Table 3 .
Convergence tests for several quantities of interest in simulations with Np,e = 10 9 and 10 8 electron-macroparticles.All the quantities are given in centimetres, except for E ′ * r , which is measured in m −1 .
It appears that increasing τ has diminishing importance beyond a certain threshold and thus decreasing ν * i becomes essential for increasing |E ′ r * |.The role of ρ * i appears to be less important for turbulence suppression.The simulated transition layers are approximately 3 − 6 cm wide, and whether the * i (which in this work is equivalent to small n) and large τ values.