How fluctuation intensity flux drives SOL expansion

Predictions of heat load widths λq based on particle orbits alone are very pessimistic. This paper shows that pedestal peeling-ballooning (P-B) magnetohydrodynamic (MHD) turbulence broadens the stable scrape-off layer (SOL) by the transport, or spreading, of fluctuation energy from the pedestal. λq is seen to increase with Γε , the fluctuation energy density flux. We elucidate the fundamental physics of the spreading process. Γε increases with pressure fluctuation correlation length. P-B turbulence is seen to be especially effective at spreading, on account of its large effective mixing length. Spreading is shown to be a multiscale process, which is enhanced by the synergy of large and small-scale modes. Pressure fluctuation skewness correlates well with the spreading flux–with the zero crossing of skewness and Γε spatially coincident–suggesting the role of coherent fluctuation structures and the presence of intermittency in λq broadening. λq∼Bp−1 scaling persists for the broadened SOL. We show that the spreading flux increases for increasing pedestal pressure gradient ∇P0 and for decreasing pedestal collisionality υped∗ . This trend is due to the dominance of peeling modes for large ∇P0 and low υped∗ . Ultimately, we see that a state of weak MHD turbulence, as for small ELMs, is very attractive for heat load management. Our findings have transformative implications for future fusion reactor designs and call for experimental investigations to validate the observed trends.

Entrainment, the tendency of a patch of turbulence to expand into a neighboring laminar region, is a well-known phenomenon.A classic example of entrainment is the turbulent 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.wake [1], which expands or 'spreads' behind a moving object such as a boat, and expands downstream (distance ∼x) with width w ∼ x 1/3 .Entrainment is a mixing process, where a turbulent region invades and mixes into a laminar region.Thus, entrainment is naturally characterized by a mixing length which is considered as the distance traversed by a fluid mass before it loses its individuality by mixing with neighboring masses [2], related to the characteristic length scale of the turbulence.Entrainment occurs in wakes, boundary layers and jets, and in similar phenomena in plasma physics, where it is called turbulence spreading [3][4][5][6][7][8][9][10][11].Since the interface between laminar and turbulent regions are likely to be multifractal, mean field theory approaches are limited, and descriptions of entrainment involving jets [12], avalanching [13] and directed percolation [14] have been developed.
In this paper, we show how the magnetohydrodynamics (MHD) turbulent core (specifically the pedestal) of a plasma confinement device can entrain an otherwise laminar boundary (i.e.scrape-off layer (SOL)), thus naturally broadening or thickening the latter.This is of great interest, since such a broadened boundary (SOL) enables a broader distribution of the heat load on the plasma facing components, specifically on tokamak divertor plates.Significantly, this work represents the first simulation/modeling endeavor to incorporate the convective transport of fluctuation energy, which distinguish it from the predominantly diffusive mechanism employed in previous published models.
Improved confinement in H-mode leads to ExB shear suppression of SOL turbulence.As a consequence, heat load width collapses to unacceptably narrow sizes.This trend is captured by pessimistic scalings predicted by the heuristic drift-based model (HD model) [15] and ITPA multi-machine experimental (Eich) scaling law [16], in which the SOL width λ q is inversely proportional to the poloidal magnetic field B p in low-gas-puff H-mode tokamak plasmas.Extrapolation to ITER leads to pessimistic results, with λ q ⩽ 1 mm.
Recently, both experiments and simulations show that in small/grassy edge localized modes (ELMs) regimes, the ELM size remains small but the SOL width λ q is broadened, while high plasma confinement persists, in contrast to the case of type-I ELMs [17][18][19][20][21][22][23].Here ELMs stand for the edge-localized modes, being characterized by quasi-periodic relaxation events occurring at the edge pedestal H-mode plasmas.Experimental observations from AUG and JET show the broadening of λ q for H-mode plasma with high separatrix density or collisionality [24][25][26].The mechanisms for the λ q broadening in the high density/ collisionality regimes are: (1) parallel confinement time in the SOL is enhanced due to large thermal resistivity at high collisionality [27,28]; (2) perpendicular transport is enhanced by the onset of additional instabilities in the SOL, such as resistive MHD modes, etc [29].In addition, simulations by the XGC1 and BOUT++ codes predict a significant broadening of λ q for ITER [30,31], on account of enhanced radial transport.BOUT++ transport simulations confirm that a transition from magnetic drift to turbulent transport regime occurs as the turbulent cross-field diffusivity increases [31][32][33].Even though both XGC and BOUT++ simulations claimed victory for ITER, they did not identify or elucidate the physical mechanism underpinning SOL broadening.In this paper, we show that in the small/grassy ELM regime, the SOL width is significantly broadened by turbulence spreading from the pedestal to the SOL.Along the way, we see that the peak divertor heat flux is much smaller than that for type-I ELMs, and fundamental physics of the spreading process is elucidated.
EAST experiments show that there is a strong relationship between the divertor particle flux width λ js and upstream density fluctuation level S fluc near the separatrix, as shown by figure 1.Here the data are the results of statistical average calculation based on multiple discharges with the similar parameters as listed in this letter.In the ELM free regime, the density fluctuation level is low, and λ js remains steady (blue squares).This is comparable with that for inter-Large ELM regime (black crosses) and is also consistent with the prediction made by the HD model.Upon entering the grassy ELM regime, the density fluctuation increases, and the divertor particle flux width increases concurrently.This suggests that turbulence may be at work broadening the SOL.However, this requires the transport of turbulence energy into the SOL.The question remains though, what is the physics mechanism underpinning the SOL width broadening driven by turbulence?Thus, we turn our attention to turbulence spreading.In this letter, we use BOUT++ [34] turbulence nonlinear simulations to investigate the impact of turbulent transport on heat flux width broadening, for the specific four EAST experimental discharges with shots #103 751, #103 745, #103 748 and #090 949 [19,35].We show that turbulence spreading is responsible for heat load broadening, and that SOL broadening is enhanced for steep pedestal pressure gradient and low collisionality.
Physics model-To study turbulence spreading from the pedestal to the SOL, a simple model is used to identify the key physical variables for the SOL width broadening.Let us consider only the radial convective evolution of the pressure equation: Here the damping and cross-field coupling are neglected.Fluctuations are defined by taking the difference to the The evolution of the turbulence energy can be thus derived as follows: The first term on the right-hand side in equation ( 4) represents energy transfer from the mean field to the turbulence.This term accounts for the local drive of turbulence in the SOL.The second term corresponds to the divergence of a flux Γ e = v r p2 , which signifies the transport of the fluctuation intensity from the pedestal to the SOL, effectively transfer it from unstable to stable regions.This term characterizes the spreading of turbulence, while the flux itself Γ p , is specifically recognized as the fluctuation energy density flux Here c s is the sound speed, p 0 is the equilibrium pressure.To compare these two processes, we introduce the production ratio R a = Γ e /P sol , which is the ratio of spreading Γ e = ´SOL sep − ∂ ∂r v r p2 dr = v r p2 | sep from pedestal to the SOL to the locally initiated turbulence production P sol = ´SOL sep − v r p ∂⟨p⟩ ∂r dr in the SOL.R a characterizes the nature of the origin of SOL turbulence.When R a > 1, the SOL width broadening is mainly induced by the turbulence spreading from the pedestal to the SOL.When R a < 1, the SOL width broadening is due to locally produced SOL turbulence.If R a ∼ 1, then effects of edge turbulence spreading, and local SOL production make comparable contributions.In this letter, we will focus on the impact of turbulence spreading from pedestal to the SOL on the heat flux width broadening.
The BOUT++ six-field two-fluid turbulence code [34] is used to calculate the turbulence spreading and SOL width for EAST experimental equilibria with low poloidal magnetic field B p = 0.16T for shot #090 949 (small ELM) and high B p = 0.21T for shot #103 745 (small ELM), #103 748(small ELM) and #103 751 (large ELM).First, the flux surface averaged pressure fluctuation intensity (p/p) 2 is calculated to investigate the relation between the edge turbulence and the SOL width λ q .Figure 2(a) shows the simulation results for the time evolution of pressure fluctuation intensity (p/p) 2 at the last closed flux surface (LCFS), 2(b) shows the divertor heat flux width calculated accordingly, for these four different EAST discharges.The ideal peeling-ballooning (P-B) mode is stable inside the pedestal for shot #103 748 and there is no ELM crash, with flat SOL profiles [35].The fluctuation intensity driven by drift modes is saturated at the lowest level in the nonlinear phase, as shown by the black curve.The diverter heat flux width is ∼4 mm, as shown by the black star on the figure 2(b).The simulated λ q is close to the value obtained from the Eich scaling law with λ Eich q ∼ 4.1 mm [16] and calculated by the HD model with λ HD q ∼ 4.5 mm [15].For shot #103 745 and #090 949, small ELMs are triggered by marginal P-B instability inside the pedestal [19,35].For shot 103 751, the pedestal is unstable to the P-B mode, with large linear growth rate.A strong avalanche process occurs during the nonlinear phase, leading to a large ELM.The fluctuation intensity increases from the ELM-free to the small, and further to the large, ELM regime.The divertor heat flux width λ q increases as the fluctuation intensity at the LCFS increases, as shown by figure 2(b).The resulting λ q is larger than that calculated from the Eich scaling law.These simulation results show a trend consistent with experimental results (see figure 1).The BOUT++ simulated particle flux λ js in the small/grassy ELM regime is comparable with EAST experimental measurements by probes [19].It is worth noting that the simulated results for large ELMs correspond to the ELM crash phase, commonly referred to as the intra-large ELM phase.The measuring λ js during the intra-large ELM phase in experiments can be challenging.To investigate how the turbulence generated in the pedestal spreads to the SOL, along with the correlation between the turbulence spreading and SOL width broadening, the fluctuation energy density flux Γ ε = c 2 s ⟨ v r (p/p 0 ) 2 ⟩ is examined.The flux surface averaged and time averaged Γ ε at LCFS in the nonlinear phase is calculated using the BOUT++ turbulence code.Figure 3(a) shows 3D plot of the heat flux width λ q vs poloidal magnetic field B p and fluctuation energy density flux Γ ε , while figures 3(b1) and (b2) are the corresponding 2D plots.If there is no ELM, the fluctuation level saturates at a low level, and the divertor heat flux width follows the Eich scaling law [16], as shown by the blue star and blue circle on the black solid curve of figure 3(b1).As the separatrix fluctuation energy density flux increases, the heat flux width is broadened significantly, for both high B p (starts in figure 3(b1)) and low B p (circles and bullets in figure 3(b1)).Figure 3(b2) shows the relation between the heat flux width λ q and fluctuation energy density flux Γ ε scanned from the ELM-free to the ELM regimes.λ q is significantly broadened by turbulence spreading from the ELM-free (blue circle for high B p and blue star for low B p ) to the small-ELM (from red bullet to green circle in figure 3(b2)) and further to the large ELM (green bullet) regime.P-B turbulence is especially effective at spreading, on account of its large effective mixing length as shown in figure 4. In the ELM free regime, λ q is small due to the small micro-turbulence mixing length l c .In the small ELM regime, λ q increases compared to ELM-free cases, likely due to the synergistic effects of the weak P-B and drift modes.The blue shaded region of figure 4, where significant scatter is observed, may be viewed as a cross-over regime between drift and turbulence dominated regime.As the transition from the small to the large ELM  leading to smaller λ q in a small ELM regime, as demonstrated by the transition from the hollow black star to hollow yellow star in figure 3. Furthermore, based on equation ( 4), and subsequent discussion, it is evident that both the spreading from pedestal to the SOL and local SOL turbulence production contribute to the broadening of the SOL width.BOUT++ simulations clearly demonstrate that the locally driven SOL turbulence, by steepening of SOL profiles, results in an increased broadening of the heat flux width.This effect is illustrated in figure 3, where the transition is observed from the blue circle to the blue bullet for shot #103 748, and from the red bullet to the red circle for shot #103 745.Overall, in the ELM-free regime, λ q follows the Eich scaling law.However, in the ELM regime, λ q is predominantly determined by the outward turbulence spreading at the separatrix.λ q is broadened for increasing edge fluctuation energy density flux Γ ε [36].λ q ∼ B −1 p scaling persists for the broadened SOL.
Drift-Alfvén instability (DAI) and weak P-B instability play synergistic roles in the onset of ELM and turbulence spreading in the small ELM regime.There is no ELM triggered by DAI alone.λ q follows the Eich scaling law, as shown by the blue star in figure 3(b1) for shot #90 949.The small ELM is triggered by P-B mode for cases both with (red star in figure 3(b)) and without DAI (black star in figure 3(b)).Without DAI, the ELM crash driven by weak low-n peeling mode in the pedestal occurs rapidly (figure 5(a)).The introduction of DAI delays the onset of the ELM crash (figure 5(b)) due to the mode-mode coupling [37,38].Following the ELM crash, the combination and interaction of large-scale P-B and small-scale drift-Alfvén modes enhance turbulence spreading.
Figure 5 illustrates the intermittency resulting from a weak low-n peeling turbulence, characterized by the occurrence of sporadic events over time.The skewness of pressure fluctuations, S p , shows strong correlation with the spreading flux [39], suggesting the role of coherent fluctuation structures and the presence of intermittency in λ q broadening, as depicted by figure 6.The spatial coincidence between the zero crossing of Γ ε and the location denoted by S p ∼ 0 indicates the origin of turbulence, and its alignment with flux transport process.Positive S p represents outward turbulence spreading with positive Γ ε , while negative S p represents inward turbulence spreading with negative Γ ε .S p near the separatrix is larger for the case with DAI, as compared to the case without DAI.Γ ε increases as S p increases.The time-averaged Γ ε in the nonlinear saturated phase for with DAI is 2.7 × 10 9 m 3 s −3 , ∼2.3 time larger than that without DAI (1.16 × 10 9 m 3 s −3 ).This leads to increased heat flux width broadening, as shown from the red star to black star in figure 3(b).
Having established that the heat flux width is broadened by enhanced turbulence spreading from the pedestal to the SOL, we now investigate how the pedestal parameters affect the fluctuation energy density flux Γ ε .Scans of peak pedestal pressure gradient ∇P 0 and pedestal collisionality υ * ped are performed [35].Figure 7 shows the 3D plot of fluctuation energy density flux Γ ε vs pressure peak gradient ∇P 0 and pedestal collisionality υ * ped .Here the black curves are ∇P 0 scans with low collisionality υ * ped = 0.108 (solid curve) and high collisionality υ * ped = 1 (dashed curve); while the red curves are υ * ped scans with small ∇P 0 ∼ 200 kPa m −1 (solid curve) and large ∇P 0 ∼ 400 kPa m −1 (dashed curve).The turbulence spreading from the pedestal to the SOL strongly depends on the pedestal plasma parameters ∇P 0 and υ * ped .The fluctuation energy density flux at LCFS increases as ∇P 0 increases and υ * ped decreases.This leads to the question of what is the physics mechanism for the drop of the fluctuation energy density flux Γ ε as the pedestal collisionality υ * ped increases?BOUT++ nonlinear simulations indicate that turbulence spreading from pedestal to SOL depends on the radial mode structure from linear to nonlinear phase, since the latter sets the mixing length for transport.For high collisionality plasmas, the mode structure is narrower, so the effective mixing length and the outward turbulence spreading decrease, leading to lower fluctuation energy density flux Γ ε .The low-n peeling mode induces more fluctuation energy flux Γ ε for low collisionality, as compared to the high-n ballooning mode for high collisionality.This trend is due to the wider radial mode structure of the peeling mode.However, strong peeling turbulence can induce a larger type-I ELM crash and larger heat load on the divertor.Thus, weak peeling turbulence-akin to MHD turbulence-such as occurs in the small ELM regime, is an attractive option.
In conclusion, improved confinement in H-mode leads to ExB shear suppression of SOL turbulence.As a consequence, the heat load width collapses to unacceptably narrow sizes and pessimistic scaling predicted by the HD model.We show that a state of weak MHD turbulence-as for small ELMs-is very attractive for SOL broadening driven by the transport of fluctuation energy from the pedestal to the SOL.The fundamental physics of the spreading process is elucidated in the letter.We demonstrate that λ q is proportional to the fluctuation energy density flux Γ ε at the last close flux surface (LCFS).Spreading is shown to be a multiscale process in the small ELM regime, which is enhanced by the synergy of large-scale P-B and smallscale drift modes.Pressure fluctuation skewness correlates well with the spreading flux, suggesting the role of coherent fluctuation structures and the presence of intermittency in λ q broadening.The outward turbulence spreading from pedestal to the SOL iptsshownto be sensitive to pedestal parameters, such as pedestal pressure gradient ∇P 0 and pedestal collisionality υ * ped .The fluctuation energy density flux Γ ε increases for increasing ∇P 0 and decreasing υ * ped .Peeling modes are seen to be especially effective for spreading on account of their large radial extent, which endows them with a large mixing length.The fluctuation energy density flux Γ ε induced by low-n peeling modes, which dominate at low collisionality, is larger than for high-n ballooning modes, which occur at high collisionality, on account of the wider mode structure of the former.This research offers new insights into SOL width broadening driven by outward turbulence spreading in the small ELM regime and identifies key parameters and suggestions for experiments prior to ITER.To mitigate Type-I ELMs in the ITER baseline scenario, it is imperative to lower the pedestal pressure height, consequently entering a state of marginal plasma instability.This state is accompanied by the presence of weak MHD turbulence, ultimately leading to the broadening of the SOL width [32].This strategic compromise emerges as the most optimal approach available to ITER for proficiently maintaining highconfinement conditions and effectively handling divertor heat exhaust.
This work was performed under the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract No. DE-AC52-07NA27344, LLNL-JRNL-845 666 and by UCSD under Contract No. DE-FG02-04ER54738.This work was also supported by the Users with Excellence Program of Hefei Science Center, CAS under grant Nos.2021HSC-UE014.The authors would like to express their gratitude for valuable discussions with X. Chu, Z. Y. Li, M. Y. Cao, T. Long, R. J. Hong, F. Khabanov, G. R. McKee, the BOUT++ team and the EAST pedestal experimental team.

Figure 1 .
Figure 1.The divertor particle flux width vs density fluctuation at outer midplane near the separatrix for EAST ELM-free (blue squares), inter-Large ELM (black crosses) and grassy ELM (red stars).

Figure 2 .
Figure 2. (a) Time evolution of fluctuation intensity at last close flux surface (LCFS); (b) heat flux width vs time averaged fluctuation intensity at LCFS in the nonlinear phase.
zonal and time average ⟨•⟩ = ´dθdψ dt as the mean quantity p = p − ⟨p⟩.For simplicity, the radial background flow is neglected (⟨v r ⟩ = 0) and the plasma turbulence flow is treated as incompressible, with ∂ ∂r v r + 1 r ∂ ∂θ v θ = 0.The evolution of the mean pressure, neglecting sources, sinks and dissipation for simplicity, is given by:

Figure 3 .
Figure 3. (a) 3D plot of heat flux width λq vs poloidal magnetic field Bp and fluctuation energy density flux Γε; (b) 2D plot of heat flux width λq vs poloidal magnetic field Bp (b1) and fluctuation energy density flux Γε (b2).The solid curves in (a) and (b1) are for the Eich scaling law [16] and the dashed curves are for the error bars.

Figure 4 .
Figure 4. Radial correlation length of pressure fluctuation lc vs. heat flux width λq.

Figure 5 .
Figure 5.Time evolution of fluctuation energy density flux for without (a) and with (b) drift-Alfvén instability.

Figure 6 .
Figure 6.Radial profiles of normalized fluctuation energy density flux Γε (blue) and skewness (red) for without (a) and with (b) drift-Alfvén instability.Here fluctuation energy density flux is normalized to the max value for each case.