Reduced model for H-mode sustainment in unfavorable ∇B drift configuration in ASDEX Upgrade

A recently developed reduced model of H-mode sustainment based on interchange-drift-Alfvén turbulence description in the vicinity of the separatrix matching experimental observations in ASDEX Upgrade has been extended to experiments with the unfavorable ∇B drift. The combination with the theory of the magnetic-shear-induced Reynolds stress offers a possibility to quantitatively explain the phenomena. The extension of the Reynolds stress estimate in the reduced model via the magnetic shear contribution is able to reproduce the strong asymmetry in the access conditions depending on the ion ∇B drift orientation in agreement with experimental observations. The Reynolds stress profile asymmetry predicted by the magnetic shear model is further extended by comparison with GRILLIX and GENE-X simulations matched with comparable experiments in realistic X-point geometry. The predictions of the radial electric field well depth and its difference between the favorable and unfavorable configurations at the same heating power from the extended model also show consistency with experimental measurements.


Introduction
Access to the high confinement regime (H-mode) [1] in tokamaks is today considered one of the key ingredients for developing a fusion reactor.However, to this day several aspects of the physics of accessing H-mode are still not fully understood.
One of them is the dependence on the X-point position relative to the main plasma, known for more than 30 years since the times of ASDEX [1] and reproduced in many other devices such as ASDEX Upgrade [2], Alcator C-Mod [3] and DIII-D [4].In the so-called favorable ∇B drift configuration associated with the ion ∇B drift pointing towards the active X-point, ample empirical evidence has led to the development of scaling laws for the threshold power necessary to access the H-mode.
However, in the unfavorable configuration (ion ∇B drift pointing away from the active X-point) the power threshold is typically observed to be at least 2 times larger than in the favorable configuration [5,6,7].Furthermore, the radial electric field E r shear at the plasma edge has been observed to be significantly weaker in the unfavorable configuration with respect to the favorable configuration at the same heating power and density, see [8] and references therein.. Understanding this asymmetry in H-mode access with respect to the X-point orientation is significant both for testing models of H-mode access as well as assessing the operational range of alternative confinement regimes such as the I-mode, which is typically accessed in the unfavorable configuration due to the H-mode being less accessible [9].
Attempts at explaining the behavior often focus on the experimental observation that the radial electric field E r profile at the edge of the plasma before the transition to H-mode is very different between the two configurations.Several theories and physical mechanisms for explaining this behavior have been proposed, among others mechanisms dependent on drift directions such as SOL flows setting the boundary conditions for E r [10], neutral penetration [11], or mechanisms based on the X-pointinduced asymmetry such as ion orbit losses [12], and acceleration due to different Reynolds stresses [13].The magnetic-geometry-induced Reynolds stress theory was able to recover the difference in E r profiles in Tore Supra [14], and more recently an extended version of the model similarly recovered E r profiles in WEST [15].
In this article a somewhat different approach based on the separatrix operational space (SepOS) model [16] is taken, rather focusing on an established H-mode with a typical E r ∼ ∇p i /(en i ) profile [17], but with turbulence suppression rates facilitated by the Reynolds stress [18].The Reynolds stress in turn is assumed to depend on the X-point configuration along the lines of Reynolds stress symmetry breaking as in [13] and [19].Therefore, this article focuses conceptually rather on H-mode sustainment and thereby the H-L back-transition, rather than on the H-mode access and the L-H transition.The possible impact of any hysteresis effects will be discussed later.
The article is structured as follows: In section 2 the empirical observations in ASDEX Upgrade are interpreted in the framework of the SepOS model and differences between the favorable and unfavorable configurations are empirically characterized by a single factor.In order to find a possible explanation for this factor, the Reynolds stress X-point dependence mechanism from [13] is briefly reviewed and recast in a more intuitive form in order be included in the SepOS in section 3. To quantify and better understand the Reynolds asymmetry, state-of-the-art full-f simulations are analyzed in section 4, extending the preliminary analysis in [20].Finally, the insight gained from simulations is combined with the SepOS model to consistently explain the empirical factor.As a consistency check, a proof-of-principle reconstruction of differing L-mode E r profiles from [8] is attempted in section 6.
The reader may notice that the article is long and dense with equations, with many plus and minus signs.We beg for the reader's understanding and patience, since surely one must admit that a solution to a problem eluding understanding for the last 30 years is unlikely to be trivial.

Empirical ∇B difference from the perspective of the SepOS model
Conceptually, the H-mode turbulence suppression criterion derived in [18] and employed in the SepOS model [16] is given by the following balance of the energy transfer from turbulence to the shear flow (which suppresses the turbulence) on the left-hand-side and the turbulent energy production on the right-hand-side: ϕ,pe ⟨u y ⟩∂ x ⟨ũ x ũy ⟩ > γ e (E tk + E te ) + γ i E ti (1) with the Reynolds stress ⟨ũ x ũy ⟩ and the mean shear flow ⟨u y ⟩ on the left-hand-side, representing together the so-called Reynolds work -transfer of turbulent energy to the flow.The x and y coordinates represent the locally field-aligned radial and bi-normal (to both the magnetic field and x, i.e. close to poloidal) directions, respectively.The transfer is made less effective by low adiabaticity represented by the cross-phase δ ϕ,pe between the potential and pressure.The turbulent energy production is approximated by the characteristic quasi-linear growth rates of electron γ e and ion γ i modes and the kinetic E tk and free energies E te , E ti .Through further quasi-linear approximations and by normalizing the equations to the kinetic energy scale E tk the criterion in normalized units can be understood as a relation between an effective turbulence energy suppression rate γ suppr,eff (stabilizing term) and an effective turbulence energy production rate γ prod,eff (destabilizing term) α c =κ 1.2 (1 + 1.5δ) with the normalized collisionality turbulence parameter α t (closely related to the cylindrical safety factor q cyl and the typical edge electron collisionality ν * ∝ q cyl Rn e /T 2 e ) describing the cross-phase of potential and pressure perturbations (i.e.drift-wave or interchange transport), the typical normalized electrostatic wavenumber scale k 2 ES ∝ β (previously called k EM ‡), the critical normalized pressure gradient α c capturing shaping dependencies (of elongation κ and triangularity δ) and the normalized curvature ω B , set by the major radius R normalized by the electron pressure gradient length λ pe , describing interchange growth.The ion to electron pressure gradients ratio τ i Λ pi (assumed to be close to 1) represents the mean flow approximated by the diamagnetic flow ⟨u y ⟩ ≈ ∇p i /neB (which normalizes to τ i Λ pi ≈ 1) on the left-hand-side and approximates an ITG-like growth rate on the right-hand-side.The approximation of the growth rate by ∝ √ ω B is consistent with curvature in the vicinity of the outer midplane (OMP).For details the reader is referred to Ref. [16].
Note that the effective suppression and production rates are by construction rather close to a shearing rate ∂ x u y and linear growth rates, respectively, reminiscent of the more commonly applied shear suppression criteria ∂ x u y ∼ γ [21].However, here these effective rates are not exactly the same thing due to being constructed from energy transfer arguments.
For the present discussion it suffices to point out that all the terms in (2) are functions of measurable quantities: Firstly the equilibrium parameters such as the cylindrical safety factor q cyl , major radius R and shaping (average elongation and triangularity), and secondly of kinetic quantities: the electron temperature T e and density n e and the edge pressure gradient length λ pe := −p e /∇ r p e .The three kinetic ‡ the original name in [16] comes from the derivation of a typical wavenumber scale below which electro-magnetic (EM) effects become significant.However, subsequently a factor 2 larger wavenumber scale is used for the quasilinar approximation of electro-static (ES) turbulence behavior (i.e. 2 times larger wavenumbers are expected to be far from EM effects, hence should be dominantly ES), therefore, the EM subscript is somewhat misleading and is replaced with the more accurate ES quantities are estimated at the very edge of the confined plasma in the vicinity of the separatrix in ASDEX Upgrade from Spitzer-Härm power balance analysis and Thomson Scattering profiles for stationary plasma discharge windows of ∼ 200 − 400 ms as explained in detail in Ref. [16].
For each of the windows where this estimation is done the ratio of γ suppr,eff /γ prod,eff as the criterion for (2) for H-mode turbulence suppression can be thus evaluated and compared with the actual identified confinement regime.This is shown in the first row of figure 1 for several combinations of magnetic fields and plasma currents and for the favorable and unfavorable drift cases.The range of these parameters of the databases is summarized in table 1.All the discharges considered here are in deuterium, typically shortly after boronization with no additional impurity seeding, hence a low effective charge Z eff ≈ 1.24 ± 0.13 (part of the displayed errorbars) is assumed as in [16].While Z eff is not actually known at the separatrix location, the comparison with a core Z eff estimate shows no systematic difference between the favorable and unfavorable configuration datasets presented here and agrees with the assumed Z eff range.For the unfavorable cases both lower and upper single null (LSN and USN) discharges with corresponding unfavorable field orientations are used.The most important observation is that the ratio γ suppr,eff /γ prod,eff at the horizontal boundary separating L-modes and H-modes remains constant across a wide density range, as well as a wider range of magnetic fields and currents as summarized in table 1, effectively spanning the whole operational range of ASDEX Upgrade.For the favorable case the ratio is 1 by construction, since (2) was constructed in such a way so as to describe the data.In other words, this is a data-driven approach which attempts to capture the scaling of the data through physics-informed parameters.But the fact that a constant, albeit higher, ratio is also found for the unfavorable cases is a non-trivial result.This shows that not only these normalized variables manage to capture the various (and even non-monotonic) scalings of the L-H boundary with current, magnetic field, density, etc., but additionally this specific ratio of the terms apparently contains the physics missing to explain the difference between the favorable and unfavorable cases.
The ratio near the L-H boundary for the unfavorable case appears to be about 2.5.This means that for a given effective production rate in the unfavorable case one needs a 2.5 times higher suppression rate compared to an equivalent favorable case for H-mode to exist.This can be also understood as the suppression rate being only 1/2.5 = 0.4 = 40% as efficient as in the equivalent favorable case with the same production rate.The main The first row shows the ratio of the two sides of (2), i.e. the γ suppr,eff /γ prod,eff ratio.The second row shows the database points in the (n e , T e ) space normalized to a standard scenario and the blue H-mode suppression line corresponding to (2) for the standard scenario.The third row shows the corresponding power crossing the separatrix P SOL estimated from the normalized temperatures.Note that the normalization for some discharges (lower toroidal field scenarios) results in a P SOL not actually achievable in ASDEX Upgrade.
question now becomes: What does this factor mean and why does it have this value?
Before this question is tackled, it is illustrative to show how this larger factor in the unfavorable case impacts the power threshold.To this end, the criterion (2) can be interpreted also in the (n e , T e ) operational space which is related to the power entering the SOL P SOL ∝ T 3.5  e .For a given equilibrium (i.e.fixed q, R, shaping) and for each n e value (2) can be solved as a non-linear equation for T e , given a relation for the gradient length λ pe (T e , n e , . . .).For these cases the H-mode ELMy scaling λ pe ∝ ρ s,pol (1 + C • α 1.9 t ) is used as in [16] which describes the widening of the gradient lengths at higher collisionalities.The solution at each density T e,LH (n e ) for the standard AUG scenario and shape with magnetic field B t = 2.5 T and current I p = 0.8 MA then becomes the blue line, similar to an L-H power threshold, in the second row of figure 1.Since the experimental T e points come from a wider range of current and magnetic fields and shapes, they are shown normalized by the ratio of T e,LH (n e ) for the standard I p , B t and their scenario T e,LH (n e , 2.5 T, 0.8 MA)/T e,LH (n e , B t , I p ).The blue line again well separates L-mode and H-modes.Most importantly, the line is significantly higher (at the same density) in the unfavorable case than in the favorable case.
Through the formula P SOL ∝ T 3.5 e λ q (α t ) the normalized temperature translates to a significant increase in the L-H power threshold as shown in the third row.Interestingly, due to the widening of the gradient and power widths dependent on α t , the P SOL is about 2.5 times larger at higher densities in the unfavorable configuration, though at lower densities (less widening) it could be even more than 4 times higher compared to the favorable case at the same density.Therefore, it qualitatively reproduces this difference in power threshold increase in the low and high density branches reported in [8].
One possible explanation for the factor 2.5 is that either the left or right-hand side of (2) is actually missing some factor.For instance, if the stabilizing term on the lefthand side of (2) had an additional factor of 0.4 (i.e. was only 40 % efficient) in the unfavorable and 1 in the favorable cases, then (2) could describe both configurations consistently.An obvious candidate for such an explanation is the pressure gradient ratio τ i Λ pi which was only assumed to be close to 1 in typical discharges in the favorable case.However, dedicated measurements of this ratio in the vicinity of the separatrix show that it is very unlikely that the ratio would be systematically τ i Λ pi ≈ 0.4 only in the unfavorable cases.
Therefore, a different mechanism for changing the terms was sought.
In the next two sections 3 and 4 an explanation for the factor 2.5 is sought by extending the model to have the energy transfer through the Reynolds stress dependent on the X-point geometry based on ideas from [13,19].However, the result of these sections is essentially just a single numerical factor and several general observations, and so the impatient reader may jump to section 5 to see how these are used in explaining the 40% factor.

Coupling with the shear-induced Reynolds stress model
Because the Reynolds stress enters in the stabilizing term of (1), it follows that the theory of the impact of the X-point asymmetry on the Reynolds stress as described in Ref. [13] is an obvious candidate for extending the SepOS model.In the following the pieces of the two theories [16] and [13,14] relevant for the coupling are reviewed and the theories are shown to be ideally suited for such a coupling.
As seen in ( 1), the terms are approximated mainly through characteristic wavenumbers and linear growth rates under the assumption of k x = k y = k ⊥ set to k ES which as will be illustrated later in section 5 is generally consistent with the shear suppression criterion.However, for the following considerations based on eddy tilting in [13,14] it is necessary to distinguish k y and k x .
In particular, the Reynolds stress is approximated in [16,14,13] through because the turbulent velocity fluctuation amplitudes are approximated by typical scales of fluctuating E × B velocities ũx ≈ ik y φ and ũy ≈ −ik x φ with the fluctuating plasma potential φ.The not-so-trivial derivation of the minus sign on the right-hand side of ( 3) is discussed in more detail in section Appendix A. The magnitude of the magnetic field B is "hidden" by the normalization of wavenumbers to the hybrid (sound) gyroradius ρ s = √ m i T e /eB, i.e. k x is in fact k x ρ s in SI units.The actual amplitude of the plasma potential fluctuation energy φ2 in (2) cancels through the normalization by the kinetic energy, therefore, only the description of the Reynolds stress through wavenumber scales is sought for the coupling.The averaging ⟨. . .⟩ brackets represent an ensemble average over turbulent fluctuations and also a flux surface average, because the Reynolds stress can poloidally vary on a flux surface, as will be shown later, and only the "net" suppression effect is of interest.
The poloidal profile of the Reynolds stress can be deduced from the relationship with the structure tilt: At a given poloidal location on a flux surface the covariance of the turbulent velocity fluctuations resulting in the Reynolds stress is not only due to the amplitude (variance) of the fluctuations, but also due to the correlation of the velocity fluctuations.The simplest model for correlated velocity fluctuations is a linear correlation ũy ∝ ũx .Locally (at a given poloidal location), such a co-linearity can be related to the average tilt of stretched (sheared) structures as illustrated in figure 2. The proportionality constant α u can be related to the corresponding wavenumbers through with α k the ratio of the wavenumbers.The ratios can be related to the geometric angles with the hat accent α as illustrated in figure 2. Note that the ratio of the wavenumbers α k is opposite (or in geometric angle π/2 shifted) with respect to α u , because it reflects the angle of the corresponding wavevector (in a plane-wave-like approximation), rather than the geometric tilt of the ellipse-like structure.Although much of the theoretical considerations in [13,14] are based on α k , in the following α u will be used instead, because it more naturally follows the intuitive geometric structure tilt mental image in real space.With this formalism the local Reynolds stress can be expressed as where ⟨. . .⟩ t represents the local (at a given poloidal location) ensemble-only average §.
The original assumption of k x = k y in [16] is, therefore, equivalent to setting the average tilt to αu ∼ −45 • as in the right image in figure 2. This would be consistent with a strong flow shear near the transition to H-mode.However, here the theory of [14] comes in with a more elaborate relation for α k also dependent on the magnetic shear.In this theory the impact of both the flow shear as well as the magnetic shear along the field line on ballooning-like structures elongated along the field line (small k ∥ ) originating dominantly around the outer midplane (OMP) at poloidal angle θ = 0 is considered, as illustrated in figure 3.An initial radial wavenumber spectrum with zero mean is assumed, so that on average (in ensemble sense) such structures will experience a tilt along the field line on the flux surface § the subscript t is motivated by the fact that it is typically estimated through a sufficiently long time average expressed through the formula [14] (recast here in terms of α u ) with the typical eddy life time over which the structure experiences the flow shears τ , the mean shearing rate ∂ x u y , and the normalized magnetic shear parameter ŝ = ∂ ln(q)/∂ ln(x) = (x/q)∂ x q.For convenience the two tilt components are split into the flow-shear-induced tilt α (u) u and the magnetic-shear-induced tilt α (ŝ) u (which incorporates the minus sign).
The minus sign in α u reflects the fact that a positive magnetic shear over a counterclockwise poloidal angle θ displacement tilts the structure clockwise.Note that in Ref. [13] ∂ x u y is denoted as V ′ E and the relation ( 6) is between k r and k ⊥ .Indeed, a certain equivalence between the poloidal and bi-normal wavenumbers is assumed, which at the very edge of the plasma with a sufficiently high safety factor is likely reasonable.
The formula (6) remains the same for any helicity (because ŝ does not depend on the sign of q) and magnetic field orientation (incorporated in the flow shear ∂ x u y sign).Broadly speaking, in a typical tokamak ŝ ≈ 2 > 0 holds at the very edge of the confined plasma near the separatrix [23].Therefore, in a given magnetic equilibrium shape α (ŝ) u is expected to remain the same.On the other hand, the sign of the flow shearing rate and thus also α u can be changed in the experiment by the toroidal field direction, and this can, therefore, either act favorably or unfavorably relative to the magnetic shear.
However, the favorable tilting superposition above the OMP shown in figure 3 would become unfavorable below the OMP.Therefore, a certain poloidal asymmetry is required for a non-zero net effect.
To this end [13] introduced another abstraction convenient for the coupling; The separation of the poloidal dependence of the potential fluctuation energy into a function F ϕ representing only its poloidal dependence through φ2 (θ) = F 2 ϕ (θ) φ2 (θ = 0).This separation facilitates the cancellation of φ2 in the approximations leading to (2).In [13] k y was assumed to be poloidally constant in a weak magnetic flux expansion limit, thus only the poloidal dependence of the potential fluctuation energy was considered.For greater generality here the actual radial velocity fluctuation energy ũ2 x ≈ k 2 y φ2 will be considered and the assumption of a poloidally constant k y will be further investigated in the following sections.Therefore, the F function will be defined here as Combining this abstraction with (5) results in where the proportionality constant α RS between the Reynolds stress and the average radial velocity fluctuations energy ⟨ũ 2 x ⟩ is the average (now both in ensemble and flux surface sense) structure tilt weighted by the poloidal distribution of velocity energy F 2 .The average tilt is again separated into components due to the magnetic shear α RS .The form of (8) with ⟨ũ 2 x ⟩ factored out is particularly convenient for the normalization of (1) by the total kinetic energy at the OMP as discussed later in section 5, because then the not well known φ2 cancels, highlighting how these theories are well suited for coupling.
However, it is important to note that the approximation of E tk through φ2 in (2) also assumes k y = k y (θ = 0) is constant on the flux surface, therefore, it will likely underestimate or overestimate ⟨ũ 2 x ⟩ if k y = const is not the case.Thus, for the purposes of the inclusion of α RS in (2) it is necessary to rescale it to that approximation as α RS,0 in order to give the appropriate Reynolds stress estimate by As will be shown in more detail in section 5, ultimately α RS,0 scales the lefthand-side term of (2), according to the balance of α To determine what values of α RS,0 are realistic and whether this model works, higher fidelity simulations will be investigated in the following section 4.
One might be inclined to also deduce that a seemingly second point of coupling would be the mean poloidal flow ⟨u y ⟩ in the suppression criterion (1), because its modification by the Reynolds stress is described by the poloidal momentum balance equation developed in [13] and slightly modified in [15].
However, here the two theories employ somewhat different assumptions due to being applied in different circumstances.In [16] the turbulence is assumed to be suppressed (in average sense) in H-mode conditions described by (1), and thus the mean flow is approximated by the diamagnetic ion flow ⟨u y ⟩ ≈ (∇p i )/enB owing to the observation that in H-mode in both the favorable and unfavorable configurations, and in the absence of significant momentum input, the E r minimum (well depth) is typically observed to be close to this value set by neoclassical physics [17,24].
Therefore, in this model the Reynolds stress acts only as a "conduit" which enables the transfer of any free turbulent energy (represented on the right-hand side of (1)) generated due to steep pressure gradients back to the mean flow (set by the pressure gradient), thereby maintaining suppression.

Equilibrium Reynolds stress asymmetry model
To validate a model for α RS one would ideally use a full poloidal profile of the Reynolds stress measured in an experiment.While several measurements of the Reynolds stress in various devices and using different diagnostic techniques exist, most have been done only in the vicinity of the OMP or did not cover a sufficiently dense poloidal profile to enable a flux-surface average estimate.
Therefore, one can attempt to validate against the next best thing, which is a stateof-the-art turbulence simulation well matched with experimental observations.For this particular study the best option at the time of writing appeared to be the GRILLIX full-f fluid simulation of edge plasma in realistic X-point equilibrium geometry, validated and well matched against an AUG L-mode discharge thanks to the inclusion of a simple neutrals diffusion model [25].Similarly, the GRILLIX simulations of both favorable and unfavorable TCV L-mode discharges in the TCV-X21 dataset [26] will be analyzed to investigate the difference between favorable and unfavorable configurations.It is however fair to note that in particular the AUG GRILLIX simulation has not fully matched the experimental radial electric field E r profile [27], therefore, the flow patterns discussed below could in principle be somewhat different in reality.
In the following quantities from the GRILLIX simulations will be shown averaged over the converged portion of a given simulation, i.e. they will be shown with the ensemble approximation average ⟨. . .⟩ t applied.However, for simplicity this bracket operator will be omitted in the following, except when its application is not obvious.
Since F 2 is assumed to be approximately an even function of θ (symmetric about the OMP) and α (ŝ) u approximately odd in θ, it is evident that for a non-zero Reynolds stress originating from the magnetic shear component the average ⟨F 2 α (ŝ) u ⟩ must be non-zero either due to a poloidal asymmetry of F 2 or α (ŝ) u .In [14] u was assumed to be a simple odd function of θ with a constant ŝ on the flux surface, and thus the asymmetry had to come from F 2 .To this end F 2 was approximated by a poloidally periodic Gaussian centered around the OMP, representing a ballooned fluctuation drive envelope related to normal curvature, but truncated around the X-point due to resistivity effects around the X-point.Initial investigation of the GRILLIX F 2 (both potential and velocity) for the reference simulation shown in figure 4 indeed suggests that this picture is reasonable.Nevertheless, it is also evident that the profiles of ⟨ũ 2 x ⟩ and φ2 do not have the same shape, thus the assumption of a constant k y is not valid.
When it comes to calculating a flux surface average of the Reynolds stress, it is, however, necessary to not look only at the geometric angle as was done in [13].For this reason, poloidal profiles of quantities will be shown against the constant volume differential angle θ V as defined for instance in [28], because θ V is the effective angle over which the true flux surface average is performed as Therefore, poloidal profiles with respect to θ V represent the distribution∥ of contributions to the flux surface average.The reference location θ V = 0 is associated with the OMP where interchange-like structures are assumed to dominantly originate.Essentially, θ V is "stretched out" with respect to the geometric angle around locations of low poloidal field.
In figure 5 the poloidal profiles of the discussed quantities at ρ pol = 0.999 are shown.This radial location just a little inside the separatrix is chosen, because it corresponds to the model in (1), and also conveniently there is little flow shear ∂ x u y ∼ 0 which enables to study the impact of the magnetic shear in a somewhat isolated manner.
Firstly, it is evident in figure 5 panel (a) that the Reynolds stress profile displays a significant asymmetry, particularly in the peak-to-peak magnitudes.In order to investigate to what extent this asymmetry is due to asymmetries in F 2 the profiles of the potential fluctuation energy φ2 and the radial velocity fluctuation energy ũ2 x ,  corresponding to the un-normalized F 2 profile shape in the case where no poloidal variation of k y /B is neglected or included, respectively, are shown in figure 5 panel (b).Interestingly, in terms of the flux-surface angle θ V the profiles appear to be quite symmetric compared with the Reynolds stress asymmetry.Particularly in terms of their poloidal extent, the truncation by the primary X-point below the OMP as shown in figure 4 and assumed in [13] is not as pronounced anymore.The reason for this approximate symmetry in θ V particularly near the X-point regions is that the lower B p in these regions compensates through (10) for the truncation asymmetry in the realspace geometry.The small difference between the peak heights of ũ2 x would account only partly for the much larger difference in the magnitudes of the Reynolds stress peaks.
Therefore, the pronounced Reynolds stress asymmetry appears to be not just due to F 2 in this case.Instead it seems the dominant origin of the asymmetry has a different reason.Inspired by the explanation offered in [19] that it is the asymmetry of the local magnetic shear, an attempt to generalize the α  x and potential fluctuation φ2 energies (b), average effective tilt structure angle α u (c) and effective poloidal wavenumber ky at the reference radial location ρ pol = 0.999 in GRILLIX simulation of AUG discharge 36190 [25].

following.
Reduced model for H-mode sustainment in unfavorable ∇B drift configuration in ASDEX Upgrade16 The generalization of the magnetic shear of field lines is where the geometric local magnetic shear S l,g is integrated along the field line length l up to angle θ relative to θ = 0 at the OMP.It is directly related to the local magnetic shear S l defined in [29] through S l,g = sign(q)S l .Here the sign of the safety factor q carries the signs of the magnetic field vectors relative to the fixed typical geometry orientation, i.e. counter-clockwise toroidal and poloidal angles, because in [29] S l is defined with respect to the magnetic field vector orientation and will change signs depending on the field orientation, whereas for the following geometric considerations S l,g is relative to the aforementioned fixed geometry orientation.This sign cancellation by q is also illustrated in the high-aspect-ratio symmetry limit S l qR → ŝ.
In figure 5 panel (c) the effective angle tilt profile α GRILLIX u = ⟨ũ x ũy ⟩ t /ũ 2 x in the simulation is compared with the tilt angle α (S l ) u calculated through (11) with the local magnetic shear.The agreement is already quite reasonable in terms of roughly capturing the asymmetry, which is due to the local magnetic shear approaching 0 towards the active X-point coming from the OMP.Nevertheless, the agreement can be further improved by considering the poloidal variation of k y in the following.
An effective poloidal wavenumber is shown in figure 5  x /( φ/B t ) 2 based on the E × B velocity prescription.The local toroidal magnetic field strength B t is used since AUG is not a high-aspect-ratio tokamak and B t varies poloidally on a flux surface.The toroidal field is used in here instead of the total field one would expect in an E × B velocity, because GRILLIX approximates the perpendicular gradient in the poloidal plane only.
It is evident that kGRILLIX y significantly increases in the proximity of X-points relative to its value near the OMP.The decrease further beyond the X-points is not truly representative, because both ũ2 x and φ2 vanish towards 0 (the latter slower) and there is not much turbulence to speak of, hence their ratio is not entirely meaningful in those regions.While approximating ky analytically is generally difficult due to the non-linear nature of the studied turbulence regime, the proximity of the increase in ky to the X-points motivates the following argument based on the magnetic field line scales.
Let us assume that through field line helicity the poloidal wavenumber k y is related to a parallel wavenumber k ∥ .Further assuming there is a dominant parallel wavenumber, perhaps related to a critical normalized gradient α c as assumed in [16], we can relate the phase increment along a field line increment ∆l with the same increment projected onto a poloidal arclength increment ∆l pol as k ∥ ∆l ≈ k y ∆l pol .Since the ratio of the length increments is given by the field line pitch angle, the poloidal wavenumber can be approximated as k(B) where k∥ is a reference parallel wavenumber.Setting k∥ to agree with kGRILLIX y at the OMP the comparision of k(B) y with kGRILLIX y in figure 5 panel (d) shows that this approximation indeed captures the trend of increasing k y towards the X-points rather well.Naturally, there is a big disagreement once the estimation of kGRILLIX y is meaningless due to the vanishing turbulence levels.
Considering that ( 6) is derived in eikonal formulation in [13] effectively through an increment of the radial wavenumber ∆k x = k y ŝ∆θ, (11) can be generalized using the k y approximation from (12) as As can be seen in figure 5 panel (c) this refined approximation α of the magneticshear-induced tilt agrees remarkably well with the simulated effective tilt α GRILLIX u .This result points towards the magnetic shear asymmetry as being the dominant origin of the Reynolds stress asymmetry, since it accounts for the large asymmetry in α GRILLIX u , which in turn accounts for the main asymmetry in the Reynolds stress.
Altogether, the average tilt calculated from the profiles shown in figure 5 following (8) (i.e. the ratio of the Reynolds stress and radial fluctuation energy flux surface integrals) gives α RS ≈ −0.24.When the variation of the poloidal wavenumber is not taken into account and the result is rescaled according to (9), the average tilt is α RS,0 ≈ −0.35.
The next step is to investigate whether this magnetic-shear-induced asymmetry remains the same between the favorable and unfavorable configuration when only the magnetic field direction is changed.To this end a similar analysis has been applied to the TCV-X21 favorable and unfavorable cases and the results are displayed in figure 6.Additionally, the analysis was performed also for the corresponding simulation of the favorable case with GENE-X [30].For these cases the radial reference location a little deeper inside ρ pol = 0.99 is shown in order to to highlight the impact of the finite flow shear.
It is clearly seen in figure 6 panel (a) that the Reynolds stress profiles significantly differ between the favorable and unfavorable GRILLIX cases.However, their shape is rather similar, the main difference appears to be mainly of an offset-like nature.The favorable GENE-X case is rather different, though still has a sine-like dependence.
Firstly, the asymmetry of the of the potential and radial velocity energy fluctuations are investigated to determine if they are the dominant source of the Reynolds stress asymmetry.Clearly, the asymmetry in the potential and radial velocity fluctuation envelopes in panels (b) and (d), respectively, between the two GRILLIX cases is quite significant.But this is not due to an X-point truncation (the X-point location does not change) as assumed in [13], but rather a "skewing" of the whole envelope.This might suggest that additional effects also play a role in the Reynolds stress asymmetry.
In [13] it was proposed that the envelope could poloidally move due to the Bloch shift x energies (b), average effective tilt structure angle α u (c), potential fluctuation energy φ2 (d), and effective poloidal wavenumber ky (e) at the reference radial location ρ pol = 0.99 in GRILLIX simulations of favorable and unfavorable TCV-X21 cases [26] as well as the GENE-X favorable case [30].
of the interchange growth rate to the poloidal location where the wavenumbers are the smallest, i.e. to the location with α u = 0.However, the peaks of φ2 in figure 6 panel (d) go in the opposite direction, rather following the flow direction.
To further investigate this additional φ2 asymmetry, a similar analysis was also done for the corresponding favorable GENE-X simulation [30].There the asymmetry of φ2 was found to be significantly smaller.Furthermore, the envelope is even more ballooned, likely due to the presence of trapped electron mode turbulence (TEM).Inclusion of a Landau fluid closure into GRILLIX [31] to get a more realistic heat flux model, i.e. taking it closer to GENE-X, also leads to more symmetric profiles.This suggests that the large envelope asymmetry may not be entirely realistic and could be down to the limitations of the Braginskii model used in those GRILLIX simulations.
Therefore, since the φ2 asymmetry is likely unrealistic, we turn our attention again to the estimated average tilt α u in figure 6 panel (c).The general asymmetry shape of α GRILLIX,GENE−X u besides the offset appears to be the same in both the favorable and unfavorable cases.The α (S l ,B) u formula again captures the α GRILLIX,GENE−X u asymmetry shape quite well, up to a numerical factor of 0.7.This overestimation (in absolute value) by α (S l ,B) u can be attributed to the effective correlation ρ corr = ⟨ũ x ũy ⟩ t / ũ2 x ũ2 y being significantly below 1 especially around the OMP, i.e. to the departure from the near-perfect linear correlation of velocities assumed in (4).This may be due to to the turbulence regime being different compared to the AUG case and/or could be a consequence of the smaller size of the TCV geometry, where the ratio of the turbulence structures size and the equilibrium geometry scales may be larger than in AUG.
Nevertheless, the main asymmetry in the shape of the Reynolds stress profile in all cases shape again appears to be mainly due to the magnetic shear asymmetry, as was the case with the AUG simulation.
The offset-like difference seen most clearly above the OMP between α GRILLIX u in the favorable and unfavorable GRILLIX cases is related to the finite flow shear and consistent with (4), where α (u) u is the offset which changes sign according to the shear flow orientation.The offset between the two cases is also very small near the active X-point, which can be understood by the local shearing rate vanishing ∂ x u y → 0 due to strong flux expansion, but the lifetime τ not approaching infinity at the same time since even in this region turbulent structures have a finite average lifetime, i.e. α (u) u → 0 in this location.
Looking at the average tilts, one finds α RS ≈ −0.27 and α RS,0 ≈ −0.85 for the favorable case, and α RS ≈ 0.07 and α RS,0 ≈ 0.26 for the unfavorable case, respectively.Assuming that the magnetic α For completeness let us state that for the favorable GENE-X simulation α RS ≈ −0.11 and α RS,0 ≈ −0.31 is found.Since the average tilt profile α GENE−X u in figure 6 panel (c) crosses 0 close to θ V = 0, it appears it is less affected by the flow shear component.Therefore, it can be viewed a being more representative of mainly the magnetic shear component, similarly as appeared to be the case in the AUG simulation discussed before.
Therefore, altogether both the AUG GRILLIX and GENE-X TCV cases, as well as the estimate from the difference of the two favorable/unfavorable GRILLIX TCV cases, appear to point to the average tilt due to the magnetic shear at the very edge of the plasma being about α (ŝ) RS,0 ≈ −0.3 in a LSN configuration.

Combined separatrix H-mode sustainment model
For the final combination of ( 1) and ( 8) it is necessary to approximate τ ∂ x u y .Since (1) aims to describe H-mode suppression in the vicinity of the OMP, it is reasonable to assume τ ∂ x u y ∼ 1 as in typical turbulence shear-suppression criteria 1/τ ∼ ∂ x u y .However, this is likely to be fulfilled only around the OMP, because further away towards the X-points the E × B flow (assuming ⟨ϕ⟩ t is close to a flux function) and the corresponding velocity shear is weaker due to flux expansion, as was already suggested by the analysis of figure 6.Therefore, likely |α On the other hand, because around the OMP likely |τ

is a reasonable assumption in the state of H-mode suppression and so is the original approximation of
for the other terms in (1).Specifically, the ∂ x gradient of the Reynolds stress in (1) giving the Reynolds force was in [16] ¶ approximated by the typical turbulence scale ∂ x ≈ k ES .Since this operator acts already at the OMP on the average Reynolds stress, this approximation is retained.Similarly, the total turbulent kinetic energy at the OMP remains , which is then equivalent to the radial velocity energy.
In principle this radial velocity energy ⟨ũ 2 x ⟩ t (θ = 0) does not have to be equal to the flux-surface average ⟨ũ 2 x ⟩, but the analysis in section 4 suggests that they are sufficiently close, therefore, the original approximation is also retained.
The mean flow is again approximated as the ion diamagnetic velocity in the electron diamagnetic (hence no negative sign) direction.To capture its poloidal orientation along θ, as shown in figure 3, dependent on the magnetic field orientation, the approximation is generalized as where the model normalization (for details see [16]) is done by the cold-ion sound speed c s = T e /m i and the drift scale ρ s = √ T e m i /eB, and the positive sign of the toroidal magnetic field, sign(B t ) > 0, corresponds to the counter-clockwise orientation as viewed from above, i.e. along the toroidal angle φ shown in figure 3. ¶ In [16] the actual approximation used the unnecessary Fourier analysis ansatz with the imaginary unit ∂ x ≈ ik ES which then necessitated another imaginary unit in the mean velocity ⟨u y ⟩ = −iτ i Λ pi .Through the full treatment of the velocity orientation in ( 15) and ( 14) this argument is no longer required.
This approximations implies that the toroidal rotation component is assumed to have no significant radial shear, because ⟨u y ⟩ approximates the shearing rate as will be explained in the following.The approximation could also possibly somewhat overestimate the mean flow in L-mode conditions.Nevertheless, even with this caveat the criterion ( 1) is still able to separate L-mode and H-mode conditions sufficiently well as seen in figure 1.
In order to get the signs of the suppression term dependent on the field and flow orientations right, it is necessary to revisit the suppression energy transfer term on the left-hand-side of (1).The term ⟨u y ⟩∂ x ⟨ũ x ũy ⟩ describing the net transfer of energy to the flow is often called the "Reynolds work".Although "Reynolds work" is the customary name, it is actually a transfer of energy, so it corresponds in terms of units to a power (per mass).This term is an approximation of the actual suppression transfer of only turbulence energy to the flow ⟨ũ x ũy ⟩∂ x ⟨u y ⟩ as derived in [18].The approximation holds under the assumption that the local sources of flow momentum are negligible ∂ x (⟨ũ x ũy ⟩⟨u y ⟩) ≈ 0 as shown in the following partial integration Most importantly, the last term has a minus sign when the flow orientation is fully taken into account.For the approximation of the radial derivative of the Reynolds stress again the typical turbulence scale is used, because the Reynolds stress in the outer shear layer is the largest (in absolute value) near the separatrix with the largest shear and lowest (in absolute value) deeper inside towards the E r well with vanishing shear.Altogether, the normalized stabilizing term + on the left-hand side of ( 2) is finally approximated as which is nearly the same as in (2) except for the factor α RS,0,B := sign(B t )α RS,0 .This factor now captures the difference between the favorable and unfavorable ∇B drift directions.This is illustrated in table 2 for the typical experimental cases in AUG.For its construction it is important to recall that α (u) RS,0 depends on the orientation of the magnetic field as α (u) RS,0 ∝ ∂ x ⟨u y ⟩ ∝ sign(B t ), because at the outer shear layer the flow is expected to be (in absolute value) greater deeper inside, and approaching 0 near the separatrix, while α (ŝ) RS,0 is determined by the X-point location in the equilibrium.While α (u) RS,0 is not actually known, it can be inferred in the following from the fact that α RS,0,B ≈ 1 holds robustly for the favorable cases: Assuming α (ŝ) RS,0 ≈ −0.3 estimated from GRILLIX and GENE-X in section 4 for the LSN configuration is sufficiently general, this means the shear flow component is approximately α (u) RS,0 ≈ −0.7 + except for the non-adiabaticity correction 1/(1 + δ 2 ϕ,pe ).Since it is significant only at very high densities, it will not be written explicitly in the following in the interest of clarity in this case where in LSN the favorable magnetic field orientation is B t < 0 This is consistent with the expectation that it will be, in absolute value, somewhat lower than 1.For the unfavorable LSN case it is assumed that α (ŝ) RS,0 remains the same, but the shear flow now changes sign as does, therefore, α (u) RS,0 ≈ +0.7.Its magnitude is assumed to be the same at the point of H-mode turbulence suppression, which is consistent with experimental observations [8].This gives a factor of α RS,0,B ≈ 0.4, which is in very good agreement with the experimental observations in section 2. In other words, if the suppression term is by a factor 0.4 = 1/2.5 weaker, the original stabilizing term in (2) would indeed need to be 2.5 times larger.
For the USN cases the shear-induced tilt component is expected to change sign by symmetry but retains its value since most USN AUG equilibria are quite comparable to up-down mirrored LSN equilibria, i.e. α (ŝ) RS,0 ≈ +0.3.The resulting value of the α RS,0,B factor then shows which B t orientation is the favorable or unfavorable one in this configuration as shown in table 2.

X-point sign(B t )
∇B drift α Table 2. Comparison of the favorable vs. unfavorable configuration from the perspective of the global fields orientation with the average tilt factor.
Altogether, the extended model with α RS,0,B correctly describes the favorable and unfavorable behavior consistently with the usual ∇B drift orientation, though it does not have a direct physical connection with the ∇B drift.It also does not depend on the sign of the helicity, consistently with the ∇B drift behavior.

Consistency with E r profile measurements
In order to attempt to recover the difference in radial electric fields in L-mode in the different configurations reported in [8], one possibility is to attempt to use the momentum balance equation from [13] or even its combination with the spectral filament model in [15].However, those models have several free parameters for which firstprinciples-based estimates are difficult to obtain, such as an ad-hoc poloidal momentum diffusivity or the width of the tanh(x)-like gate function simulating the transition to SOL conditions.Furthermore, the α (ŝ) RS estimates obtained from the simulations are almost three times lower than those used in [15].These various parameters can have a significant impact on the final profile.Other parameters such as the average radial velocity fluctuation energy scale ⟨u 2 x ⟩, which are not directly measured, can be estimated, for instance using the sheared spectral filament model as in [15], but that model is based Reduced model for H-mode sustainment in unfavorable ∇B drift configuration in ASDEX Upgrade23 on interchange physics only and likely is not entirely appropriate for lower collisionalities with significant drift-wave presence.
Therefore, an attempt to recover the observed E r difference is made directly based on the SepOS framework with known quantities.The effective bi-normal flow velocity related to E r can be obtained from the suppression term on the left-hand side of (1) if it is normalized by the "Reynolds force" (per mass) ∂ x ⟨ũ x ũy ⟩ (neglecting the nonadiabaticity correction for clarity again).In terms of dimensions analysis the suppression energy transfer rate -i.e.power (per mass) divided by a force (per mass) indeed gives a velocity v = P/F .Expressed through the updated suppression criterion (16) this gives In order to understand to which location this mean flow ⟨u y ⟩ corresponds, it is illustrative to rewrite (17) as which expresses the implicit assumptions involved: Firstly, ⟨u y ⟩ actually approximates the true shearing rate ∂ x ⟨u y ⟩ in (15) which implicitly assumes that ⟨u y ⟩ ≈ 0 near the boundary.Secondly, the radial gradient length is assumed to again be comparable to the turbulence scale ∂ r ∝ k ES as before.But the prefactor α RS,0,B effectively makes the gradient scale similar near the L-H boundary in both the favorable and unfavorable configurations, because in the unfavorable case at the higher temperature near the L-H boundary k ES ∝ √ β e is larger, but the smaller α RS,0,B compensates that.This roughly corresponds to the E r well width being roughly the same in both considerations, in agreement with experimental observations [8].Finally, the suppression rate γ suppr,eff is expected to approximate (at least in terms of its scaling) the "average" shearing rate.From this rough formula it is evident that ⟨u y ⟩ actually corresponds to a radial location deeper inside.Since γ suppr,eff is assumed to be an approximation of the "average" shearing rate across the outer shear layer, it follows that ⟨u y ⟩ should correspond roughly to the E r minimum.
Further assuming that the background toroidal flow component is either small or is not significantly sheared (consistent with observations in [8]) and can thus be subtracted as an offset, the electric field (without the toroidal component) can be estimated as E r ≈ u y B t , and removing the normalization, the (negative) E r minimum estimate at the OMP from the SepOS model can be expressed as where, additionally, the flux expansion f x ∼ 1.8 accounts for the narrow gradient length at the OMP and the R geo /(R geo + a geo ) rescales to the toroidal field B t to the OMP relative to the reference B t used in the normalization.It may seem strange to use the separatrix temperature T e,sep for the estimation of the E r minimim deeper inside.One might expect some extrapolation deeper inside.But it is there in fact simply because T e,sep is used for the normalization.In cases with significant, but insignificantly sheared, toroidal rotation (19) represents rather −min(E OMP r + v i × B).The E r well minimum is in experiments typically taken as a proxy for the shearing rate in cases where toroidal rotation is insignificant [17].In other words, here the SepOS approximation goes the other way than in the experiment.
To simulate a power scan from L-mode towards H-mode at fixed density and global parameters, the separatrix temperature T e,sep is virtually scanned towards the expected temperature at the L-H boundary T e,LH (n e ) as a solution of (2) (i.e. the blue line in figure 1).To estimate the E r minimum in L-mode following (19) it is necessary to first approximate γ suppr,eff in such conditions.This can be done by approximating the ratio of the suppression and production rates γ suppr,eff /γ prod,eff by the proximity to the L-H boundary as γ suppr,eff γ prod,eff ≈ T e,sep T e,LH 1.5 (20) where the ratio T e,sep /T e,LH expresses the proximity to the L-H boundary.The exponent 1.5 is motivated by the fact that the rates ratio scales as γ suppr,eff /γ prod,eff ∝ T 1.5 e,sep across the wider database shown in figure 1.
To this end for each temperature point up to T e,LH (n e ) and the other known parameters the effective production rate γ prod,eff is estimated following the approximations on the right-hand side of (2).However, in L-mode a different scaling for λ pe than the one used for H-mode is required.Therefore, the L-mode scaling λ pe = 17.3α 0.298 t mm from [32] will be used.Due to the large scatter of measured λ pe relative to the scaling, the prediction will be shown also with errorband representing λ pe ± 4 mm which is representative of the scatter.
The combination of ( 20) and (19) gives the simulated −min(E OMP r ) evolution from L-mode towards H-mode shown in figure 7 and figure 8 as dashed lines for the favorable and unfavorable configurations.The T e,LH (n e ) points in either configuration are represented by star symbols.The results of these calculations are also compared with the actual E r minimum measurements reported in [8] for the medium density, USN and low density, LSN cases in figure 7 and figure 8, respectively.Since in the LSN cases the toroidal rotation component appears to be significant, the NEOART estimate of this component v i,neo × B with little shear is subtracted and the errorbars appropriately enlarged.In the USN cases the toroidal rotation component appears to be insignificant, therefore, the actual measured E r values are shown in figure 7.
The measurements agree with the prediction within uncertainties quite well, given that the estimations above were expected to capture the scaling at most.Nevertheless, it is clear that at the same heating power, i.e. the same T e,sep , the E r minimum is predicted to be shallower in the unfavorable configuration than in the favorable one, consistent with the experimental trends.
The gradual development of min(E r ) dashed line may seem at odds with  .Predicted (negative) radial electric field field E r minimum at the outer midplane for the favorable and unfavorable configurations from the SepOS models vs. experimental measurements (circle points) reported in [8] for the medium density USN case.The star symbol shows the expected E r at the L-H boundary.
experimental observations where min(E r ) typically develops more abruptly as the heating power is increased towards the transition to H-mode [8].However, it is important to note that here the figures are shown with respect to the separatrix temperature T e,sep .If they were shown with respect to P SOL ∝ T 7/2 e,sep , the rise would indeed appear to be much faster.
The predicted E r near the L-H boundary (star symbols) is, however, quite similar, again in agreement with the general experimental observation that the L-H transition appears to occur at roughly the same shearing rate or E r minimum (as the proxy of the shear with negligible or subtracted, insignificantly sheared toroidal rotation consistent with observations in [8]) in either configuration.The predicted values near the L-H boundary are close to the typically experimentally observed E r minimum [8].The calculations would suggest that the E r near the L-H boundary is slightly lower in the unfavorable case.However, this may be due to the employed λ pe scaling which was derived only for favorable L-mode cases which did not include higher temperatures  .Predicted (negative) radial electric field field E r minimum at the outer midplane for the favorable and unfavorable configurations from the SepOS models vs. measurements (with the toroidal component v i,neo × B estimated by NEOART subtracted) reported in [8] for the low density LSN case.The star symbol shows the expected E r minimum at the L-H boundary.
> 80 eV and thus may be simply extrapolating too far outside its domain of validity.Generally, the employed λ pe scaling was developed to capture the general L-mode trends and did not specifically target L-modes close to the L-H transition.
The prediction of the absolute value of |min(E r )| (neglecting toroidal rotation) or equivalently γ suppr,eff near the L-H boundary generally appears not to be entirely constant with density, possibly slightly decreasing towards higher densities.While this may be partly due to the imperfect scaling, it is broadly consistent with measurements reported in [8,33].
Regardless of all the aforementioned uncertainties and imperfect assumptions, the main result is that the model robustly reproduces the general qualitative trend of min(E r ) being deeper in the favorable configuration relative to the unfavorable case at the same heating power.

Discussion
Since the criterion 1 conceptually focuses on the H-L back-transition, it is worth considering how far from the L-H transition it may be.Although the L-H and H-L transitions can be conceivably different in certain aspects, such as in the heating power hysteresis, the typical difference is comparable to the uncertainty of the power crossing the separatrix estimate employed here.Therefore, the analysis presented herein cannot reliably distinguish the hysteresis effect anyway.Furthermore, past studies in AUG have shown that local edge quantities such as the pedestal top pressure (as a proxy for the flow shear driven by the pressure gradient) or the radial electric field do not exhibit the hysteresis seen in the power [24,34,35].Therefore, as the study in this article using separatrix kinetic quantities to approximate the flow shearing of turbulence is conceptually similar, the hysteresis effects can be expected to be limited.Finally, one might argue that for a reactor knowledge of requirements for reaching a sustainable, "locked-in" H-mode is of ultimate importance.
The original SepOS model [16], includes the tilt by k x = k y .Structures are tilted with k x = k y ∂ x u y τ which with the shear suppression criterion results in k x = k y .The shear flow can interact effectively with the turbulence only if the inverse eddy turn-over time or vorticity Ω = ∇ × ũ is approximately equal to the shearing rate ∂ x u y .This means ∂ x ũy ∼ ∂ x u y the radial derivatives of the binormal velocity are of similar order.Unlike for density and temperature fluctuations, the fluctuations in binormal velocity are not smaller than the background velocities, both being on the order of diamagnetic velocity.This is due to the fact that the fluctuations are small but their gradients are not, thus ∇ φ ≈ ∇p ≈ ∇p.As these are eddies, ∂ x ũy ∼ ∂ y ũx ∼ k 2 y φ.When the inverse eddy turn-over time is similar to the shearing rate ∂ x ũy ∼ ∂ x u y , approximating ⟨ũ x ũy ⟩∂ x u y by u y ∂ x ⟨ũ x ũy ⟩ seems reasonable.
Nevertheless, these considerations about the fluctuation amplitude should not be interpreted as leading to a zonal flow of the same amplitude as the mean flow.This is because although the fluctuation amplitude of individual turbulent structures may be large, the average Reynolds work contribution to the total flow is still expected to be only a small fraction relative to the mean flow dominated by the mean pressure gradient.
In the SepOS model, the characteristic wavenumber is not determined by the linear growth rates, i.e. it is not approximated by the dominant mode.The main principle of linear growth rate determination by determining its maximum is not applied.Therefore, it is not a linear stability calculation.Intrinsically, small scales are first suppressed by the shear flow [36,37,38].So one should look for the largest scale to be suppressed and that should be the transition from electrostatic to electromagnetic turbulence k ES .If this scale is suppressed by the shear flow, all wavenumbers above it should also be suppressed, because typical turbulence spectra decay towards higher wavenumbers, and at the same time the transfer term ∂ x ⟨ũ x ũy ⟩ ∝ k 3 is stronger.Turbulence on smaller wavenumbers is electromagnetic and then responsible for the remaining transport in the H-mode.The growth rates enter as effective nonlinear growth rates, with only the cross phase being approximated quasi-linearly.This leads to similarities with the calculations for linear growth rates.
The impact of the Reynolds stress in the end appears to not be described exclusively by either the theory in [13] or [19].Rather it is a combination of the best parts of both theories.Specifically, the theory of [13] seems to have correctly recognized the importance of the competition between the magnetic and flow shear, but attributed the origin of the Reynolds stress asymmetry to the fluctuation envelope which does not appear to agree with the presented GRILLIX and GENE-X simulations.On the other hand, [19] recognized the impact of the local magnetic shear asymmetry, but did not recognize the impact of the flow shear and instead used an ad-hoc ansatz based on the ∇B drift.
It is worth stressing that neither of the GRILLIX or GENE-X simulations analyzed actually simulates H-mode-like conditions yet, nor are the matches with the experimental radial electric field (or potential) profiles particularly good.Therefore, it is possible that a better matched simulation of H-mode conditions might reveal a behavior not observed in the current analysis.Nevertheless, the very good match of the general equilibriumbased arguments with the simulation results across two completely different machines and for both favorable and unfavorable configurations offers hope that the conclusions regarding the average tilt asymmetry shape being dominated by the local magnetic shear asymmetry may be sufficiently general.Preliminary analysis of on-going AUG H-mode GRILLIX simulations [39] indeed suggests that the conclusions about the average tilt components following (8) hold in that case as well.
The conclusions are expected to hold when the very edge of the plasma is in a turbulence regime dominated by drift-wave-interchange turbulence susceptible to magnetic shear effects, and the particle source is sufficiently localized relatively far away from this edge layer -in the sense it has limited impact on any turbulence velocity energy asymmetries.Contrarily, the arguments above of course cannot account for differences in symmetric equilibria, such as in limited circular plasmas considered in the experimental comparison in [14], where indeed only an asymmetry in the velocity distribution (possibly due to a particle source determined by the limiter location) can account for differences between favorable and unfavorable configurations.Also in double null diverted configurations the magnetic symmetry implies that any asymmetries in behavior will likely have a different origin, possibly the impact of a different E r field in the SOL due to different neutrals behavior.Such fine details related to a balance of particle sources and sinks could of course be very machine-dependent, possibly giving rise to often contradictory results observed in the community.
The presented results may appear to suggest that the difference between the favorable and unfavorable configuration, or even the L-H transition physics generally, is only due to the very edge of the plasma just inside the separatrix, in the outer shear layer at the pedestal foot.This may seem to be at odds with the observations that the shear suppression occurs first around the inner shear layer deeper inside [40] and more generally with the whole pedestal region build-up during the L-H transition.However, this is not in contradiction, because the conditions at the outer shear layer in the vicinity of the separatrix can be seen rather as a necessary, but not sufficient, condition for the establishment of a fully turbulence-suppressed pedestal.The choice of this position in the vicinity of the separatrix is rather a convenient balance between several motivations: Firstly, it enables a direct connection to power exhaust considerations, which makes this reduced model particularly well-suited for reactor design considerations.Secondly, in this region the very simple drift-fluid model basis is likely still applicable and exploitable for such semi-analytical, quasi-linear estimations.
The general physical model of the Reynolds work energy transfer being different and contributing less to turbulence suppression by shear flows in the unfavorable configuration is qualitatively similar to the observations in Ref. [7] where the energy transfer through the Reynolds stress to zonal flows in the unfavorable configuration was found to be weaker due to simultaneous transfer to Geodesic acoustic modes (GAM).The overall model for the turbulence suppression in H-mode by the transfer of turbulence energy to the shear flow is conceptually similar to related models using a similar transfer term [41,42].Some of the differences with respect to these models is the explicit quasilinear approximation of the production terms in (2) and the approximation of the mean flow through only the diamagnetic component in (14).Specifically, in the above the generated zonal flow contribution to the mean flow is considered to be negligible relative to the diamagnetic component due to the low turbulence level in H-mode.The impact of the mean flow shear on the production terms also isn't taken directly into account.The extension of the model towards these considerations may offer opportunities to describe phenomena such as limit cycle oscillations as was done in [41].

Conclusion and summary
The combination of the SepOS model [16] and a combination of ideas about the magnetic-shear-induced Reynolds stress symmetry breaking from [13] and [19] provides a model which can quantitatively match the experimental observations of the increased power threshold in the unfavorable ion ∇B drift configuration in ASDEX Upgrade.The main mechanism for increasing the power threshold is the weakening of the Reynolds stress which acts as a "conduit" for transferring turbulence energy to the mean flow, thereby maintaining turbulence suppression.The Reynolds stress is weakened due to the flow shear acting on average against the local magnetic shear in unfavorable configuration.This mechanism has been quantified using the average structure tilt.The total average tilt is determined mainly by the combination of the toroidal magnetic field direction determining the mean flow shear direction and the local magnetic shear asymmetry induced by the single null X-point geometry.
The quantification of the average local magnetic shear effect in terms of the average structure tilt in several GRILLIX and GENE-X simulations leads to a common estimate which agrees well with the observed favorable/unfavorable difference describing the larger ASDEX Upgrade database.Interestingly, the model describes the difference between the favorable and unfavorable configurations consistently with the usual ∇B drift orientation, though for physics reasons not directly related to the ∇B drift.
Additionally, the combined model can also qualitatively explain the difference in the radial electric field E r between the favorable and unfavorable configurations at the same heating power, as seen experimentally.The predicted E r minimum values appear to even quantitatively match quite well the experimentally observed values, though due to the uncertainties involved a claim of a full quantitative match of E r would be premature.
Future efforts will focus on the comparison of this model with upcoming H-mode GRILLIX simulations, as well as with other simulation codes, in order to determine how general the picture established here is.Similarly, the establishment and analysis of a larger, multi-machine database of separatrix operational space quantities in the framework of an on-going ITPA activity is expected to eventually yield also validation of this model on a wider scale from the experimental side.

Figure 1 .
Figure 1.H-mode suppression criterion compared with the AUG separatrix quantities database for the favorable and unfavorable in the left and right columns, respectively.The first row shows the ratio of the two sides of (2), i.e. the γ suppr,eff /γ prod,eff ratio.The second row shows the database points in the (n e , T e ) space normalized to a standard scenario and the blue H-mode suppression line corresponding to (2) for the standard scenario.The third row shows the corresponding power crossing the separatrix P SOL estimated from the normalized temperatures.Note that the normalization for some discharges (lower toroidal field scenarios) results in a P SOL not actually achievable in ASDEX Upgrade.

Figure 2 .
Figure2.Sketch of how the shape of an E × B eddy (vortex) around a potential perturbation φ determines its contribution to the Reynolds stress ⟨ũ x ũy ⟩ t .With no co-linearity of velocities (left), the average covariance of velocities over all its regions is 0, while with co-linearity (right) -specifically anti-correlation -the regions with one (ũ x ũy ) sign dominate.Sketch inspiration from[22].The blue arrows indicate the corresponding wavenumber vector ⃗ k and wavelength 1/k scales in a plane-wave-like approximation.The angles with the hat accent α correspond to the tilt ratios through α = tan(α)

Figure 3 .
Figure 3. Schematic of tilting (represented by α u ) of a structure originating at the outboard midplane due to the flow shear ∂ x u y over a characteristic time τ (blue), and the magnetic shear ŝ along field lines (orange).The displayed flux surface corresponds to a radial location just inside the separatrix, i.e. the outer shear layer.The displayed orientation of the magnetic field B t < 0 and the plasma current I p > 0 (negative helicity q < 0) relative to the counter-clockwise toroidal angle φ corresponds to the most common orientation in AUG.
determined by the combination of the X-point-induced asymmetry and flow shear orientation, respectively.
The weighting by 1/B p essentially represents flux expansion, because the flux surface average is actually a volume average over an infinitesimally radially narrow (differential) volume around the flux surface.

Figure 4 .
Figure 4. Radial velocity fluctuation energy ũ2 x and potential fluctuation energy φ2 profile in ASDEX Upgrade GRILLIX simulation of L-mode just a little inside the separatrix, shown with respect to the geometric angle θ geo .Locations in the proximity of the active or the secondary X-points are shown by dashed lines, corresponding to local minima of the poloidal magnetic field B p .
expression is made in the Reduced model for H-mode sustainment in unfavorable ∇B drift configuration in ASDEX Upgrade15

Figure 5 .
Figure 5. Poloidal profiles of the Reynolds stress ⟨ũ x ũy ⟩ t (a), the radial fluctuation velocity ũ2x and potential fluctuation φ2 energies (b), average effective tilt structure angle α u (c) and effective poloidal wavenumber ky at the reference radial location ρ pol = 0.999 in GRILLIX simulation of AUG discharge 36190[25].
panel (d), which is evaluated from the two profiles in panel (b) as kGRILLIX y = ũ2

Figure 6 .
Figure 6.Poloidal profiles of the Reynolds stress ⟨ũ x ũy ⟩ (a), the radial fluctuation velocity ũ2x energies (b), average effective tilt structure angle α u (c), potential fluctuation energy φ2 (d), and effective poloidal wavenumber ky (e) at the reference radial location ρ pol = 0.99 in GRILLIX simulations of favorable and unfavorable TCV-X21 cases[26] as well as the GENE-X favorable case[30].

Figure 7
Figure 7. Predicted (negative) radial electric field field E r minimum at the outer midplane for the favorable and unfavorable configurations from the SepOS models vs. experimental measurements (circle points) reported in[8] for the medium density USN case.The star symbol shows the expected E r at the L-H boundary.

Figure 8
Figure 8. Predicted (negative) radial electric field field E r minimum at the outer midplane for the favorable and unfavorable configurations from the SepOS models vs. measurements (with the toroidal component v i,neo × B estimated by NEOART subtracted) reported in[8] for the low density LSN case.The star symbol shows the expected E r minimum at the L-H boundary.

Table 1 .
Range of global equilibrium parameters of the favorable and unfavorable separatrix operational space databases