Pedestal particle balance studies in JET-ILW H-mode plasmas

JET-ILW type I ELMy H-modes at 2.5 MA/2.8 T with constant NBI heating (23 MW) and gas fuelling rate were performed, utilising edge localised mode (ELM) pacing by vertical kicks and plasma shaping (triangularity, δ ) as tools to disentangle the effects of ELMs,


Introduction
In an H-mode plasma [1] the confinement is improved due to the formation of a narrow transport barrier at the edge.Inside the transport barrier, the energy and particle transport are reduced, and a steep pressure gradient is formed, which gives rise to a pressure pedestal.The heightened edge pressure also leads to improved global confinement due to core profile stiffness [2][3][4][5].However, the steep pressure gradient can trigger edge localised modes (ELMs) [6], followed by a transient loss of energy and particles which are deposited onto the plasma facing components.Understanding pedestal physics is essential to predict and optimise plasma performance in current and future tokamak experiments.
The H-mode pedestal is governed by at least three interacting processes: pedestal stability, transport and sources.The stability of the pedestal in a type I ELMy regime [6] is generally described by MHD theory [7].The high pressure gradient and the edge current density at the edge drive Peeling-Ballooning (P-B) modes unstable [8,9], triggering an ELM.P-B stability provides an ultimate limit on the pedestal pressure at the onset of the ELM, but the inter-ELM turbulent [10][11][12] and neoclassical [13,14] transport, along with the sources, will determine the relative contributions from the pedestal density and temperature.The exact role of the edge particle source and pedestal transport in setting the density pedestal is still an open question [15].Reduced pedestal models [16][17][18][19][20] have proven successful in predicting the pedestal electron pressure for a wide range of plasma scenarios, but they lack a first principle based, predictive model for the edge density.
In this study, we investigate the balance of particle sources and transport in the latter part of the ELM cycle where the pedestal recovery is typically slow compared to the transient crash and fast recovery phase [21,22].The evolution of the pedestal particle content is evaluated from detailed analysis of high resolution profile reflectometry data.The edge particle source is estimated by means of interpretative edge-scrape-off-layer (edge-SOL) fluid transport simulations using the EDGE2D-EIRENE code [23][24][25].Comparison of these results allows an assessment of the contributions of sources and transport to edge particle balance.The role of MHD stability in the pedestal formation is also investigated.
The paper is organised as follows.Section 2 describes the analysed JET-ILW H-mode experiments and presents the comparison of pedestal density and temperature profiles.The analysis of the linear ideal MHD stability of the investigated pedestals is discussed in section 3. The interpretative EDGE2D-EIRENE simulations to estimate the edge particle source is presented in section 4. The results on the balance of the source and transport during the inter-ELM pedestal evolution is discussed in section 5, followed by the conclusions in section 6.

Experiments
In this work, we investigate JET-ILW type I ELMy H-mode experiments that utilised plasma shaping (triangularity) and ELM pacing by vertical kicks as tools to disentangle the effects of edge stability, ELMs, and inter-ELM particle transport on the pedestal particle balance.It has been observed in several tokamak experiments, that the pedestal confinement typically improves with increasing triangularity [22,[26][27][28][29][30][31][32][33][34][35].In JET and ASDEX Upgrade, the increase of pedestal pressure with triangularity, at otherwise fixed engineering parameters, is due to higher pedestal density being obtained for similar pedestal temperature [22, 26-29, 31, 32].It has also been reported from these two devices that the ELM frequency is reduced at higher triangularities [22,26,30,31].One proposed explanation for the beneficial effect of triangularity is the improved stability of P-B modes before the ELM event (pre-ELM) [8,[36][37][38].Increased triangularity decouples peeling and ballooning modes opening up second stability access leading to a higher stability boundary.
The application of a fast and controlled vertical plasma motion (known as vertical kicks in JET) for frequency control of ELMs has been successfully demonstrated on several tokamaks [39][40][41][42].This method on JET relies on the vertical stabilization control system to drive the fast vertical plasma motion.Vertical kicks can trigger ELMs by introducing an intermittent, local perturbation of the current density close to the separatrix.If the current perturbation is large enough to destabilise the peeling component of the edge instability responsible for the ELM onset, an ELM is triggered before the pedestal would naturally reach the P-B stability limit [39].
Comparative type I ELMy H-modes were performed at 2.5 MA plasma current and 2.8 T toroidal magnetic field with constant NBI (Neutral Beam Injection) power of 23 MW, no Ion Cyclotron Resonance Heating (ICRH) in the main heating phase 11 and constant deuterium gas fuelling rate at 4.4 × 10 22 electron s −1 .A lower triangularity (δ = δ lower /2 + δ upper /2 = 0.2) and a higher triangularity (δ = 0.3) pulse with natural ELM frequency and a higher triangularity (δ = 0.3) pulse with 36 Hz vertical kicks for ELM pacing were selected for detailed comparative analysis.The main waveforms of these three pulses are shown in figure 1. Strike point sweeping with 4 Hz frequency was introduced in these pulses for operational reasons, to reduce the heat load on the tungsten coated divertor tiles.This introduces oscillations in some of the measured plasma parameters.
The pedestal electron density (n e ) and temperature (T e ) profiles from Thomson scattering (TS) [43] and the ion temperature (T i ) profiles from the edge charge exchange spectroscopy (CXRS) system [44,45] measuring Ne and C impurities are shown in figure 2. These figures show all of the TS and CXRS profiles from the steady phase of the discharge, i.e. including both ELMs and inter-ELM phases.Besides the raw data points, figure 2 shows the modified tanh (mtanh) [46] fitted profiles as well.Figure 2(c) also shows the raw ion temperature data on the top of the mtanh fit of the electron temperature profiles (same as in figure 2(b) for comparison of T e and T i .The ion and electron temperatures are similar at the pedestal top within the measurement uncertainties of the edge CXRS system, but T i in the edge transport barrier and at the separatrix cannot be resolved.
As shown in figure 2, the pedestal pressure significantly improves with increasing triangularity mostly due to an increase in pedestal density and the change in electron and ion temperatures is small.Also, the natural ELM frequency decreases from 39 Hz at δ = 0.2 to 25 Hz at δ = 0.3.When f ELM is increased by vertical kicks at δ = 0.3 to 37 Hz (close to the natural f ELM of the δ = 0.2 pulse), the pedestal density and thus the pressure is somewhat decreased, but both still significantly higher than in the δ = 0.2 case.This observation is discussed in view of the pedestal stability in the next section.The pedestal kinetic profiles for the three representative discharges examined in this paper.These figures show all of the TS and CXRS profiles from the steady phase of the discharge, i.e. including both ELMs and inter-ELM phases.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kick in magenta.(a) Electron density data from TS and corresponding mtanh fits.(b) Electron temperature data from TS and corresponding mtanh fits.(c) Ion temperature data from CXRS.The underlying solid lines are the mtanh fits of the TS electron temperature, same as in figure (b).This is to confirm that Te ≈ T i at the pedestal top.(d) Electron pressure data from TS and corresponding mtanh fits.(The kinetic profiles are radially aligned so that the separatrix temperature is 100 eV.).

Pedestal MHD stability
Vertical kicks can help to investigate whether physics mechanisms other than improved P-B stability play a role in the observed increase in pedestal performance at higher triangularity.With the introduction of vertical kicks, the ELMs are triggered before the pedestal would naturally reach the MHD stability limit.In the δ = 0.3 plasma with vertical kicks, the pedestal pressure is degraded (by ≈25%) with respect to the δ = 0.3 discharge with natural f ELM .Despite this artificial degradation of the pedestal MHD stability, the pedestal pressure is still significantly higher (by ≈30%) in the δ = 0.3 kicked pulse compared to that of the δ = 0.2 counterpart.This implies that increased plasma triangularity may also affect inter-ELM transport, and thus lead to increased pedestal confinement.
Ideal MHD stability of the pedestal was investigated with the HELENA fixed boundary equilibrium [47] and the ELITE ideal MHD stability codes [8,9].The j − α pedestal stability diagram for the three representative discharges is shown in figure 3, where α is the normalised pressure gradient 12 and j is the normalised current density.The HELENA equilibrium is self-consistent with respect to the current profile that was evaluated using Redl's bootstrap current formula [49] and assuming a fully diffused Ohmic current with neoclassical resistivity.The pre-ELM (80%-97%) 13 mtanh fitted TS data was used as input kinetic profiles, assuming T e = T i (consistent with charge exchange measurements at the pedestal top) and line averaged Z eff with Be as a single impurity.The kinetic profiles are radially aligned so that the separatrix temperature is 100 eV.
The P-B stability boundary (solid lines in figure 3) was evaluated using γ MHD > 0.5 × ω dia as stability criterion, where γ MHD is the linear growth rate and ω dia is the diamagnetic frequency.The stability boundary is very similar for all three 12 The normalised pressure gradient is defined as in [48] ∂Ψ , where V is the volume enclosed by the flux surface, R 0 is the geometric centre of the plasma and Ψ is the poloidal flux. 13All ELM cycles are normalised to a relative time scale, where 100% corresponds to the onset of the ELM crash and 0% corresponds to the ELM crash of the preceding ELM cycle.The pre-ELM phase is defined as the 80%-97% of the ELM cycle.discharges and no clear stabilisation effect from higher δ can be seen.This is because the increase in triangularity from 0.2 to 0.3 is too small and other parameters of the experimental equilibria (such as β N , β pol , volume, etc) affecting the stability boundary more.In a separate study (not shown here), a small stabilising effect of triangularity was recovered when all other parameter differences between discharges were removed and only the shaping effect was investigated by increasing δ from 0.2 to 0.3.
The stars in figure 3 show the operational point of the pedestal as obtained in the experiment with their associated error bars derived from the uncertainty of the mtanh fit parameters.The operational points for the pulses with natural ELM frequency are close to the P-B boundary thus the ELM trigger is consistent with the P-B paradigm.The operational point is further from the stability boundary for the δ = 0.3 pulse where vertical kicks were applied to increase f ELM .This is expected as the vertical kick introduces a current perturbation at the edge which triggers the ELM in conditions when the pedestal would otherwise still be stable to P-B modes.The P-B stability analysis used the equilibrium and kinetic profiles from this stable period preceding the triggered ELMs.
It was found that the pressure pedestal is wider in both δ = 0.3 pulses than in the δ = 0.2 case.This is shown in figure 4, where the edge pressure gradients are compared for the three pedestals.The wider pedestal at δ = 0.3 allows for improved P-B stability and thus leads to higher pedestal pressure.This suggests that the improved pedestal pressure at higher δ is not necessarily a sole pre-ELM P-B stability effect, but a change in inter-ELM transport could possibly play an important role in setting the higher density and thus higher pressure when the triangularity is increased.

Estimate of the edge particle source
The role of particle transport in the inter-ELM phase is studied in section 5 in more detail where the pedestal particle balance is discussed.In this section, we present the interpretative EDGE2D-EIRENE [23][24][25] simulations which were performed to estimate the edge particle source to examine the balance of source and transport in the pedestal.EDGE2D is a 2D fluid code with realistic geometry of the SOL and divertor region, which is coupled to EIRENE, a Monte Carlo code used to calculate the neutral particle distribution.

Simulation set-up and interpretative transport coefficients
We ran EDGE2D-EIRENE in interpretative mode, where the perpendicular transport coefficient of electron particle diffusion D ⊥ (Γ e,⊥ = D ⊥ ∇ ⊥ n e ), electron and ion heat transport χ e,i,⊥ (q e,i,⊥ = −n e,i χ e,i,⊥ ∇ ⊥ T e,i ) and the divertor pump albedo were iterated until the solution fitted the measured inter-ELM (40-80% fraction of the ELM cycle) upstream n e and T e profiles (measured by TS and the fast Li-beam emission spectroscopy system [50]).The ion and electron heat transport coefficient profiles (χ e,⊥ and χ i,⊥ ) were assumed The mtanh fit of the pre-ELM (80%-97%) electron pressure profile (a) and its gradient (b) for the three discharges analysed.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kick in magenta.The error bars are derived from the uncertainty of the mtanh fit parameters.Note that while the pressure gradient is similar for the δ = 0.2 and the δ = 0.3 with kicks cases, α as shown in figure 3 is different.This is because the normalisation of the pressure gradient (α) also includes the plasma volume that is smaller for the δ = 0.3 plasmas.to be the same, due to lack of constraints to distinguish between them.
The separate contribution of ELM and inter-ELM transport was not studied with EDGE2D-EIRENE.The simulations were run until convergence, thus time independent solutions were obtained.The building up of the plasma stored energy between ELMs (dW/dt), that is equivalent to the ELM power loss (P ELM ), was not taken into account in the timeindependent simulations.Thus, the input power in EDGE2D-EIRENE was set to the power crossing the separatrix inter-ELM (P inter−ELM = P sep − P ELM ).
Cross-field drifts and a particle pinch term were not introduced in these simulations.The edge particle pinch may have an important role in the particle transport [51][52][53], but in time independent simulations the experimental n e profile could be reproduced with different variations of the diffusion coefficient and the pinch velocity, due to the lack of constraints.Thus, the particle diffusion coefficient here is regarded as an effective parameter describing the transport needed to exhaust the particle source from the plasma.
The grid for EDGE2D is generated such that it is aligned to the plasma equilibrium reconstructed using kinetic EFIT [54].The grid extends to ∼10-15 cm inside the separatrix towards the core plasma.The gas fuelling rate and location was set in accordance with the experiment.The recycling coefficient was set to 1 on the walls and the divertor targets 14 .The experimental NBI fuelling rate was estimated using the PEN-CIL code [55], although this was found to be small compared to the neutral flux crossing the separatrix.The divertor configuration around the outer strike point cannot be modelled precisely in EDGE2D-EIRENE.The wall structure had to be slightly modified around the outer strike-point so that the grid does not cross wall surfaces, as described in [56] in more detail.

EDGE2D-EIRENE solutions
As described above, the perpendicular transport coefficients were iterated until the experimental upstream kinetic profiles could be matched.However, it has been found that multiple sets of transport coefficients can reproduce the kinetic profiles within the experimental uncertainty without simultaneously constraining the main chamber recycling.This is shown in figure 5, where the transport coefficients and the corresponding upstream n e and T e profiles are shown for the δ = 0.2 discharge.For these four simulations, all the input parameters were kept the same except for the perpendicular transport coefficients.It is shown in figure 5, that the particle diffusion coefficient changes by an order of magnitude in the pedestal and the SOL, while the heat transport coefficient changes by a factor of 3 in the pedestal.Regardless of this significant variation in the perpendicular transport, the upstream profiles remain almost unchanged.It is important to note that the divertor target quantities did not vary significantly either in this scan.The variation in j sat , n e and T e at the outer target was between ±30%.For the inner target, the highest D ⊥ case is significantly colder than the other three cases.j sat , n e and T e for the three lower D ⊥ cases, however, are also within ±30% variation.
The above behaviour can be understood by examining the role of main chamber recycling in fuelling the plasma.The recycling flux is determined by ion flux arriving to the edge of the EDGE2D plasma grid.The recycling coefficient was set to 1 in these simulations (except for the pumping surfaces in the divertor), thus all particles are recycled as neutrals.By increasing the particle diffusion coefficient, the ion flux arriving at the edge of the plasma and at the wall surfaces is also increased, resulting in a higher recycling flux.The increased recycling flux leads to a higher neutral density and such increased particle source and it is this increased particle source in the plasma that compensates the effect of increased particle transport.In this way, multiple sets of particle diffusion coefficients can reproduce very similar upstream kinetic profiles.
It is therefore necessary to also constrain the main chamber particle source, so as to extract the unique set of D ⊥ and χ ⊥ profiles that matches the investigated pedestal.Figure 6 shows the ionisation source inside the separatrix as a function of D ⊥,ETB for the four simulations presented in figure 5.It can be seen that for these cases, the pedestal source increase approximately linearly with D ⊥,ETB resulting in similar upstream density and temperature profiles (as shown in figure 5) despite different set of transport coefficients.This also means that additional simulation constraints that carry information on the neutral density in the main chamber are needed to estimate the pedestal particle source.

Main chamber Dα line radiation constraint in EDGE2D-EIRENE
In order to constrain the pedestal particle source, we compared the emitted deuterium Balmer-α (D α ) radiation by the main chamber plasma in the EDGE2D-EIRENE simulation to that of the experimental measurements.The line-of-sight (LOS) of the corresponding spectrometer crosses the plasma in the main chamber, but not in the divertor as shown in figure 7(c).In this way, the measured line intensity carries direct information about the ionisation source in the SOL and the pedestal, but it does not integrate through the orders of magnitude brighter divertor emission.However, the D α emission from the divertor is reflected by the JET-ILW metallic walls [57], hence it is captured along the horizontal LOS making significant contribution to the measured line intensity.Thus reflections need to be taken into account in the analysis by forward modelling.
To estimate the effect of reflections on the measured signal, the experimental D α brightness in the divertor was taken from imaging spectroscopy measurements [58].Imaging spectroscopy provided a time-averaged (during the inter-ELM phases over the steady phase of the discharge), 2D tomographic reconstruction of the D α line intensity distribution in the divertor.The contribution of the reflected light from the divertor to the horizontal D α line intensity measurement was estimated using the CHERAB code [59,60] as illustrated in figure 7(c).CHERAB includes definitions of the full JET-ILW machine geometry and estimates of the metallic tile spectral bidirectional reflectance distribution functions.It uses the RAYSECT [61] raytracing engine to provide estimates for the reflected light intensity.We used the experimental divertor D α emission as input for CHERAB, assuming zero emission from the main chamber.In this way, the background reflected emission of the measured D α intensity on the horizontal LOS can be approximated.
The measured D α light intensity as function of time is shown in grey in figure 7. The difference of the measured signal and the background reflected emission provides an estimate for the line intensity emitted by the plasma in the main chamber (the SOL and the pedestal).These estimated main chamber D α line intensities are shown in figures 7(a) and (b) in black for the steady phase of the δ = 0.2 and δ = 0.3 natural f ELM discharges after removing the spikes caused by ELMs.The obtained signal is affected by large oscillations, in part due to the applied strike point sweeping, which increases the uncertainty on this method to constrain the neutral density in the plasma.Nonetheless, this experimental comparison provides a sanity check on the overall level of main chamber D α emission in EDGE2D-EIRENE simulations and thus excludes some of the solutions obtained with unrealistic level of particle source in the pedestal and SOL.
To compare with the estimated experimental main chamber D α line intensities, synthetic diagnostic data was produced from the EDGE2D-EIRENE solutions by evaluating (using ADAS atomic data [62]) and integrating the D α emission along the diagnostic LOS.This is shown in figures 7(a) and (b) with the coloured lines.For the simulations using higher D ⊥ values, the synthetic D α emission is significantly larger than the experimental estimate for both cases.Thus, the simulations flagged as D ⊥,ETB = 0.03 m 2 s −1 and 0.06 m 2 s −1 for the δ = 0.2 case, and D ⊥,ETB = 0.04 m 2 s −1 and 0.06 m 2 s −1 for the δ = 0.3 case are the ones deemed experimentally feasible solutions.
Since the measured D α emission along the horizontal LOS is dominated by reflections from the divertor, unfortunately it is not possible to select a unique solution from the EDGE2D/EIRENE simulations and thus fully constrain the evaluation of the edge particle source.Nevertheless, in the next section, we discuss the balance of the source and transport during the inter-ELM pedestal evolution in view of the above, experimentally feasible EDGE2D-EIRENE solutions with the consideration of the large uncertainty on the pedestal particle source.

The balance of the source and transport during the inter-ELM pedestal evolution
The main aim of this investigation is to characterise the balance of sources and transport that sets the density pedestal.In current tokamak experiments, the dominant particle fuelling in the pedestal comes from the ionisation of neutral atoms penetrating into the confined region.The dominant channels of particle transport are ELM losses, neoclassical and turbulent transport.Other loss channels, for example ion orbit losses and ripple losses, are assumed to be small.Our aim is to estimate the contribution of these different channels.For this, the transient nature of ELMs needs to be considered as well.Figure 8(a) shows the time evolution of the particle content (offset to zero at the ELM onset) evaluated from profile reflectometry during an ELM cycle for the low δ = 0.2 and δ = 0.3 natural f ELM pulses.At the ELM crash, significant amount of particles is lost from the plasma.But most of the lost particles are quickly recovered.The start of the ELM crash is denoted by the solid black line, while the end of the 'fast' recovery phase is indicated with the black vertical dashed line in figure 8(a).Experimental evidence [63] suggests that the fast recovery may be driven by an increased recycling flux as a result of the increased particle flux arriving to the divertor targets and limiter as a result of the ELM crash.A sign of this increase in particle source can be seen in the spectroscopy data in figure 8(b), where the D α light intensity in the main chamber is shown.As discussed in section 4, the quantitative interpretation of this signal is challenging, but it shows that qualitatively the ionisation source in the main chamber is increased following the ELM crash.
The 'fast' recovery phase is followed by a slower and longer phase (to the right of the horizontal dashed line).In the slower recovery phase, the main chamber D α light intensity is approximately constant.We focus on the slow recovery phase to study the pedestal particle balance assuming that the rate of change ( The question we seek to answer is to see whether the pedestal recovery is dominated by the source or the transport.To answer this question, we are combining together the data analysis procedures and simulation tools presented in previous sections.The workflow is summarised in a flowchart in figure 9. dN/dt is evaluated from the experimentally measured evolution of the particle content using high resolution profile reflectometry as explained earlier in this section.The source is taken from the experimentally relevant EDGE2D-EIRENE simulations from section 4.This requires the introduction of a simulation constrain for the main chamber D α emission in EDGE2D-EIRENE.The main chamber D α intensity measurement is largely affected by reflections from the divertor.The reflections are accounted for by using the CHERAB code which takes the experimental divertor D α emission distribution as input.The magnitude of pedestal transport is then calculated as a difference between the source and dN/dt.The transport channel is then further divided into contributions from neoclassical and turbulent transport15 by estimating the neoclassical particle flux using the NEO drift-kinetic code [64,65].The NEO code takes a plasma equilibrium and corresponding kinetic profiles as input and calculates the neoclassical particle and heat flux profiles as a result.The peak of the particle flux profiles in the pedestal is used in this analysis.Note that the turbulent particle flux is not modelled in this work, it is simply deduced from equation (1), assuming that the total transport (T) is the sum of the neoclassical and turbulent fluxes.
The bar charts in figure 10 compare the magnitude of the source with that of the transport for the two investigated plasmas (δ = 0.2 in figure 10(a) and δ = 0.3 in figure 10(b)).We show two results for both plasmas.These represent the estimated particle source from the 'experimentally relevant' EDGE2D-EIRENE simulations as shown in figure 7. We regard the two solutions (for each pulse) as a lower and upper bound estimate for the particle source, and the corresponding transport flux is calculated respectively using the experimentally measured dN/dt.
For the δ = 0.3 pulse (figure 10(b)), it is clear that both the source and transport terms are relevant and the difference of these two large terms results in a small dN/dt.The contribution of neoclassical transport is small, D NC ≈ 0.004 m 2 s −1 for the δ = 0.2 case and D NC ≈ 0.006 m 2 s −1 for the δ = 0.3 case.This suggests that the role of turbulent transport is important in the evolution of the pre-ELM pedestal and it cannot be neglected.In the case of the δ = 0.2 discharge, the picture is less clear.Due to the uncertainties on the main chamber D α radiation constraint, the relative difference in the upper and lower bounds is greater than that of the δ = 0.3 pulse.In the D ⊥,ETB = 0.03 m 2 s −1 case the source seems to be dominant over transport, but in the D ⊥,ETB = 0.06 m 2 s −1 case the source and transport terms are comparable.In summary, we conclude that the role of turbulent particle transport, in the pulses analysed, cannot be neglected from the pedestal particle balance.But, further improvements required both in spectroscopic measurements and edge transport simulations to be able to quantitatively compare the pedestal particle balance between different pedestals.

Discussion and conclusions
In the present paper, the pedestal particle balance and the role of pedestal MHD stability have been investigated in JET-ILW type I ELMy H-mode plasmas at different triangularities and ELM frequency.We have observed that the pedestal pressure is higher at higher δ even when the pedestal stability is artificially degraded by introducing vertical kicks, suggesting that improved P-B stability is not solely responsible for improved pedestal performance.This is likely due to the increased width of the pedestal at higher δ that allows for higher pedestal pressure at similar pressure gradient.The physical mechanism responsible for the wider pedestal has not been investigated in this work, but a change in inter-ELM pedestal transport as a result of enhanced shaping offers a candidate explanation.Unfortunately the workflow presented in this work cannot provide quantitative comparative evidence of transport between the pedestals with different triangularities as the estimate on the pedestal particle source is too uncertain.
The literature contains several gyrokinetic studies pertaining to the impact of plasma shaping on core turbulence [66][67][68], but we have not found any systematic study investigating the effect of triangularity on pedestal turbulent transport.In reference [66], Belli mentions that nonlinear gyrokinetic simulations of core turbulence capture some of the shaping effects found experimentally, but they do not completely explain the degree of this dependence on triangularity.It may be that much of the experimentally observed strong triangularity dependence comes from the edge turbulence, which sets the boundary conditions for core transport.
In current tokamaks, the characterisation of the neutral source is essential to understand the formation of the density pedestal structure, as also suggested by several studies [15,[69][70][71].However, due to the lack of direct measurement of the neutral density, unravelling the role of the edge particle source and transport in setting the density pedestal structure is highly challenging.In this work, we used edge-SOL transport simulations (EDGE2D-EIRENE), together with various plasma measurements for the estimation of the edge particle source.Comparing the estimated pedestal particle source with the experimentally inferred evolution of the pedestal particle content inter-ELM, we found that the role of turbulent particle transport, in general, cannot be neglected from the pedestal particle balance.Cases exist where a large source term is balanced by a large transport term resulting in a relatively slow recovery of the pedestal.This implies that for detailed pedestal prediction, the properties of both the particle source and transport need to be characterised.
The workflow presented in this study offers a way to constrain the overall magnitude of the pedestal particle source.However, the separate contribution of divertor and main chamber recycling, and external gas injection to pedestal fuelling has not been investigated in detail and requires further analysis.In EDGE2D-EIRENE, the grid on which the fluid equations are solved (EDGE2D grid) does not extend to the limiters.Thus, the simulation of main chamber recycling requires assumptions on the incidence ion flux.We used a projection algorithm that results in all the ion flux arriving at the edge of the plasma grid being recycled in the main chamber.In the simulations presented in this work, the source of neutrals crossing the separatrix was very much dominated by main chamber recycling.Neutrals recycled in the divertor were found to contribute much less to the pedestal fuelling even in the simulations with low particle transport coefficients.The choice of the projection algorithm possibly contributes to this behaviour.In reality, some of the ion flux arriving to the edge of the EDGE2D grid may arrive at the divertor target before it would hit the main chamber limiters.A different projection algorithm that redirects some fraction of the ion flux to recycle at the divertor target (for e.g.exponential fall off with an exponentially decreasing fraction going to the target) would be expected to affect the contribution from the divertor recycling to pedestal fuelling.But due to the lack of experimental constraints, the choice of the projection algorithm is somewhat arbitrary.Further analysis focusing on the ratio of main chamber and divertor recycling would require the edge transport simulation codes that allow for extended plasma grids, covering the whole vessel and thus is outside the scope of this paper.In our simulations, the total amount of pedestal fuelling is experimentally constrained, but the balance between divertor and main chamber recycling is not studied.It is also important to note that the inner divertor in the simulations would be colder and more opaque if cross-field drifts were turned on, thus contributing significantly more to the fuelling of the pedestal [72].
Recently, exceptional progress has been made in pedestal gyrokinetic studies [12,[73][74][75][76][77][78][79].These analyses have identified various micro-instabilities to be potentially responsible for plasma transport in the H-mode pedestal.However, the focus was mostly on the understanding of the pedestal heat transport.Our work has shown that in general, inter-ELM turbulent particle transport plays an important role in setting the density pedestal.Thus, the study of micro-instabilities responsible for particle transport at the edge is a worthwhile line of future work towards a fully predictive pedestal model.Interpretative edge-SOL modelling, such as those presented in this paper, or more direct information about the particle source in the form of neutral density measurements where possible [80][81][82][83] could provide information for gyrokinetic simulations to include a realistic particle source term and such advance the understanding of the density pedestal formation.

Figure 1 .
Figure 1.The main plasma parameters of the three representative discharges examined in this paper.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kicks in magenta.The grey area indicates the analysed time interval.(a) The total plasma current, Ip and the toroidal magnetic field, Bt.(b) The NBI and ICRH heating power.(c) The total radiated power.(d) The core and edge line averaged density.(e) The normalised plasma β, β N .(f) The core electron temperature.(g) The external D 2 fuelling gas rate.(h) The average plasma triangularity, δ = δ lower /2 + δupper/2.

Figure 2 .
Figure 2. The pedestal kinetic profiles for the three representative discharges examined in this paper.These figures show all of the TS and CXRS profiles from the steady phase of the discharge, i.e. including both ELMs and inter-ELM phases.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kick in magenta.(a) Electron density data from TS and corresponding mtanh fits.(b) Electron temperature data from TS and corresponding mtanh fits.(c) Ion temperature data from CXRS.The underlying solid lines are the mtanh fits of the TS electron temperature, same as in figure (b).This is to confirm that Te ≈ T i at the pedestal top.(d) Electron pressure data from TS and corresponding mtanh fits.(The kinetic profiles are radially aligned so that the separatrix temperature is 100 eV.).

Figure 3 .
Figure 3. Linear MHD pedestal stability diagram for the three discharges analysed.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kick in magenta.The operational points corresponding to the pre-ELM phase of the discharges are indicated with the stars and the respective error bars as a function of the normalised pressure gradient (αmax) and the normalised current density (⟨j edge ⟩max/j).The solid lines show the P-B stability boundary.

Figure 4 .
Figure 4.The mtanh fit of the pre-ELM (80%-97%) electron pressure profile (a) and its gradient (b) for the three discharges analysed.The δ = 0.2 plasma is shown in orange, the δ = 0.3 in green and the δ = 0.3 with kick in magenta.The error bars are derived from the uncertainty of the mtanh fit parameters.Note that while the pressure gradient is similar for the δ = 0.2 and the δ = 0.3 with kicks cases, α as shown in figure3is different.This is because the normalisation of the pressure gradient (α) also includes the plasma volume that is smaller for the δ = 0.3 plasmas.

Figure 5 .
Figure 5. Inter-ELM TS profiles for ne and Te (in grey) and Li-BES data for ne in the SOL (in black) in the steady phase of the δ = 0.2 discharge.The upstream ne and Te profiles of the interpretative EDGE2D-EIRENE simulations with different set of D ⊥ coefficients are shown in colour.The simulations are labelled with the corresponding D ⊥ value inside the edge transport barrier, D ⊥,ETB .

Figure 6 .
Figure 6.The particle source integrated over the plasma volume inside the separatrix is shown as a function of D ⊥,ETB inside the edge transport barrier (ETB) for the four interpretative EDGE2D-EIRENE simulations with different set of D ⊥ coefficients (#96448, δ = 0.2).

Figure 7 .
Figure 7.The measured Dα light intensity as function of time is shown in grey for the δ = 0.2 and δ = 0.3 natural f ELM discharges.The signals after the removal of ELMs and offsetting the reflected light are shown in black.The synthetic Dα line intensities for the corresponding main chamber LOS are shown in colour from the EDGE2D-EIRENE solutions using different set of transport coefficients.

Figure 8 .
Figure 8.(a) Time evolution of the particle content during an ELM cycle (offset to 0 at the ELM onset) and (b) time evolution of the Dα light intensity in the main chamber during an ELM cycle for the low δ = 0.2 and δ = 0.3 natural f ELM pulses.

Figure 9 .
Figure 9.The analysis workflow describing the characterisation of the particle balance.The orange boxes refer to experimental data analysis, while the green boxes indicate simulations.The turbulent transport (red box) is an output in the sense that it is calculated from the other terms.

Figure 10 .
Figure 10.Comparison of the source and transport in the pedestal for the 'experimentally relevant' EDGE2D-EIRENE solutions for the two investigated plasmas.(a) δ = 0.2 in orange, (b) δ = 0.3 in green.