Multiscale numerical study of the ELM-induced sputtering on the tungsten castellated divertor gap

Material sputtering and erosion are key issue in fusion science. In the steady state, tungsten sputtering can be maintained at a low level. However, tungsten is seriously sputtered by the hot plasma from edge localized mode (ELM) bursts. The time evolution and spatial distribution of sputtering on the divertor gap are still not clear. This unclearness influences the understanding of ELM-induced sputtering and the explanation of experimental phenomenon. In this study, the evolution of tungsten flux generated by ELM burst is obtained via hybrid multiscale simulation and a double-peak of tungsten flux is found. The first peak is produced by the ELM electrons, which provides great sheath potential, and the second peak is generated by the energetic ions from ELM. The castellated divertor gap (CDG) can influence the tungsten sputtering distribution due to the change of the impacting angle. The sputtering and erosion on CDG are more serious than that on the mono-block top, but the tungsten source from the gap is very insignificant because the gap size is small. The simulation results are compared with EAST experiments and a qualitative consistency is obtained.


Introduction
Material erosion and sputtering on plasma-facing materials (PFMs) are important issues in plasma-wall interaction 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.
(PWI), which can cause significant PFM damage and plasma contamination [1].Consequently, tungsten (W) is selected as the divertor material for fusion devices, such as EAST, JET, and ITER [2].Advantageously, tungsten has a high melting point, high sputtering threshold, low physical sputtering yield, high chemical stability, and low tritium inventory [3,4].Using a radiative divertor and other techniques, the plasma temperature can be reduced to below 10 eV in the steady state [5,6], which is significantly smaller than the sputtering threshold energy (209 eV for deuterium).Even with the plasma temperature below 10 eV, tungsten can be sputtered by impurities that will inevitably be present in ITER and future devices.However, the energy of the ions produced by the edge localized mode (ELM) can easily reach and even exceed the sputtering threshold energy, yielding a considerable amount of tungsten particles [7].Even in the detachment state, the energetic particles can cause serious tungsten sputtering due to the ELM 'burn-through' which is observed in the JET ITERlike wall (ILW) [8].Thus, the ELM burst in H-mode plasmas plays an important role in tungsten erosion and damage [7,9,10].Furthermore, as tungsten is a high-Z material and can cause strong radiation and energy losses, the tolerance of tungsten impurity is considerably lower than that of other impurities [11].
Various experimental and modeling studies on ELMinduced tungsten sputtering have been performed to understand the mechanism behind the sputtering process.For example, the relations between ELM-resolved divertor erosion and ELM frequency, power loss, and plasma parameters have been investigated in the JET [12].Furthermore, the tungsten impurity sputtering due to type I ELMs by pedestal electron cyclotron emission and target Langmuir probe measurements has been estimated in JET [10].The spatial distribution of tungsten erosion on the strike point tile has been investigated in ASDEX Upgrade [13].A systematical experimental observation was performed in EAST to investigate the influence of ELM frequency, pedestal temperature, collision, etc., on tungsten sputtering [14].From the theoretical and modeling aspects, an analytical kinetic model for particle recycling and tungsten erosion in the divertor during ELMs was developed and validated in DIII-D [15,16] and JET ILW [17].Carbon and tungsten sputtering in ELM burst was modeled using a one-dimensional (1D) particle-in-cell (PIC) code, and the consistency between simulation and experiment was realized [18].This implies the efficiency of PIC modeling for ELM burst was proven.The geometry of the divertor target is essential during PWI, and Dai et al performed a predictive study on it [19].
The details of transient PWI during ELM burst have not been fully understood.For example, Dai and Wang reported that the heat and particle flux on a target vary with time and that a double-peak time evolution of heat flux exists due to the velocity difference between the electrons and ions from ELM burst [20][21][22].However, previous studies have not considered this temporal evolution.Moreover, 0.5 mm-sized castellated divertor gaps (CDGs) have been reported to significantly influence the heat flux distribution [22], which has also not been considered.The lack of understanding about the details of the transient PWI process makes the explanation of the experiment phenomena descriptive and qualitative.To obtain a more precise explanation, two critical questions need to be answered.Question one: how do the double-peak particle and heat flux temporal evolution influence tungsten sputtering?Furthermore, what is the temporal evolution of tungsten flux during ELM burst?Question two: what is the influence of CDG on tungsten sputtering in the background of the doublepeak flux?
This study aims to answer the two questions.For this, a hybrid multiscale PIC simulation is performed, which is introduced in section 2. The ELM-induced tungsten sputtering flux is elaborately discussed in section 3. Section 4 presents the conclusions and the answers to the two questions.

Method
The questions are not easy to answer in both experiments and simulations.The castellated divertor is shown in figure 1.The plasma flux (along the toroidal direction) is denoted.The CDG is 0.5 mm width and is perpendicular to the plasma flux.The length of mono-blocks in the toroidal direction is 12 mm.The duration of the ELM burst is basically less than 1 ms [14].The time resolution and space accuracy of diagnostics on the divertor target (Langmuir probes and impurity spectrum, for example) are not sufficient for observing the details of ELM-induced tungsten sputtering processes (microsecond and sub-millimeter) in EAST.Moreover, few numerical models can include both ELM transport and PWI processes due to the significant temporal and spatial differences.Thus, the main challenge of simulating the sputtering process during ELM bursts is the significant temporal and spatial differences between ELM transport and PWI.The energetic particles from ELM bursts move along the magnetic field line from the upstream (source) to the target.The magnetic connection length between the source and the target is over 100 m.Furthermore, the time of reflecting an integral ELM burst should be close to or longer than 1 ms.In contrast, the physics of plasma sheath and the gyro-motion of particles must be considered during PWI modeling.To obtain accurate and reliable results, the time step and grid size should be of the order of the ps and µm in the fusion plasma environment.In other words, nearly one billion calculation steps and nearly one million grids are required for a 1D simulation to reflect a relatively realistic ELM-induced PFM sputtering process.In most PWI simulations, 1D calculation is not sufficient because the influence of the gap needs to be considered, which requires two-dimensional (2D) or even three-dimensional simulation.This is almost an impossible task under the current calculation capacities.Therefore, multiscale modeling has been proposed for the ELM-induced sputtering process.The simulation is divided into two main steps.In the first step (the macroscopic simulation), the ELM plasma transport is modeled by 1D PIC simulation to reflect the electron and ion evolution in the SOL region.In the second step (the microscopic simulation), a 2D PIC simulation is performed to study the effect of plasma sheath and the particles trajectory.This approach has been elaborately introduced in Niu et al [22].
The geometry of the CDG in the simulation is shown in figure 2(a).In figure 2(a), points A, B, C, D, and E are defined to describe the particle flux on CDG.The cardinal coordinates of the points A, B, C, D, and E are (x, z) = (−1, 0), (0, 0), (0, −1), (0.5, −1), and (0.5, 0), respectively (denoted by the X sign).The vertical surface facing the plasma flux is the wetted surface (WS).The mono-block horizontal surface between A and B is the horizontal top (HT).The gap bottom (GB) is between C and D. The vertical surface that is back to the plasma flux between D and E is described as the shadowed surface (SS).The region at z ∼ 0 and 0 < x < 0.5 mm is defined as the gap entrance (GE) because particles enter the gap through this region.For convenience, the coordinate of the CDG position is described by the curve distance L, which denotes the distance from the point A along CDG.The curve coordinates of the points A, B, C, D, and E are 0, 1, 2, 2.5, and 3.5 mm, respectively.The vertical top (VT) is a short segment on the top part of the wetted surafce in the range of (1, 1 + δL) mm in the curve coordiante L. The space interval of the distance resolution of the flux on CDG is δL = 0.0035 mm.
To reflect the effect of ELM burst on the particle flux, two-component approximation is applied in the PWI modeling.The particle source is assumed to comprise hot and cold particles.Hot particles denote the particles from ELM burst with high temperatures.Cold particles denote particles from the background plasma with low temperatures.Based on this approximation, the initial velocity of the injected particles is modeled through the double Gaussian distribution, namely, the velocity distribution functionf Here, f h and f c represent the hot and cold components, respectively.The temperatures (T h and T c ) and weights (ω h and ω c ) of the hot and cold particles are obtained via least-square fitting from the 1D macroscopic simulation.
The sputtering yield as a function of impacting energy and angle is calculated with the aid of the TRIM code [23,24].The range of impacting energy is from 200 eV to 20 keV.The impacting angle is defined as the angle between the particle velocity and the normal direction of the hitting surface, as the figure 2(b).The impacting angle is from 0 to 89.9 • with interval 10 • .

Plasma transport from ELM source to the outer target
Due to the considerable difference between the masses of electrons and ions, the transit time of parallel transport from the source edge (s = L − = 40 m) to the outer target (s = 0) is different.Here the L − is the boundary of source where particles are injected.s is the curve length of magnetic field line.The target locates at s = 0. Compared to cold ions, hot electrons from the ELM move faster and reach the target earlier.Figure 3 displays the particle position and kinetic energy at different times.Six time points are selected to represent the different phases during ELM burst.The horizontal axis is the connection length from the source boundary at s = L − to the outer target (s = 0).The vertical coordinate represents the kinetic energy of particles, which is , where v ∥ and v ⊥ are the parallel and perpendicular velocities, respectively.To clearly depict the energy evolution, the particles are differentiated by the energy 0.1 keV.The yellow and violet dots represent electrons and ions whose E i > 0.1 keV, respectively.On the other hand, the electrons and ions with E i < 0.1 keV are represented with green and cyan points, respectively.The figure also shows the electron and ion kinetic energy profiles (red and blue lines, respectively).
Figure 3(a) displays the particle energy before ELM burst (t = 0).The energy of almost all particles is below 0.1 keV.The electron and ion temperatures at the source (s ∼ L − ) are 50 eV.The electron and ion temperatures at the target are T e ≈ T i ≈ 20 eV.In figure 3, the flux direction is from right to the left.
At t = 0.05 ms (figure 3(b)), the high-energy electrons have reached the target but quite a few high-energy ions are in the SOL region (0 < s < L − ).Additionally, the energy of background electrons significantly increases.Furthermore, a significant energy gap exists between the high-and low-energy electrons, which is obviously caused by the different energies of ELM and background electrons.Some electrons with kinetic energy of close to 0.1 keV are also present, which may be produced by the Coulomb collision between electrons.Thus, the electron energy profile in the region s ∼ 30 m is enhanced, while the E i profile does not significantly change.
In figure 3(c) (t = 0.125 ms), high-energy electrons (yellow dots) are present in the entire SOL region, while highenergy ions (violet dots) are mainly present near the upstream (s > L − /2).In other words, ample high-energy electrons have arrived at the target, while only a few high-energy ions have reached the target.The different E e and E i profiles reveal this difference between electrons and ions.The electron energy near the target is significantly increased to about 0.4 keV, which is significantly greater than the energy in the source region.This is because the energetic electrons have reached the target.The ion energy near the source is greater than  In figure 3(d), the energy of high-energy ions is significantly greater than that of high-energy electrons.This is because while the faster electrons have reached the target, the slower electrons are still in the SOL region.This tendency can also be seen in the energy profiles.E e is almost constant while a significant peak of E i exists at the s ∼ 30 m. E i dramatically increases near the target due to the strong electric sheath produced by high-energy electrons.The energy gaps of both electrons and ions disappear.At t = 0.625 ms, the high-energy electrons disappear and only low-energy ions are present in the SOL region.Furthermore, the temperature of low-energy electrons decreases.Both the E e and E i profiles become almost constant and E e < E i .At t = 1.25 ms, the energy of highenergy ions decreases and only several high-energy ions exist in the SOL region.Additionally, the energies of low-energy electrons and ions decrease to below those at t = 0.625 ms.
The difference in the parallel transport of electrons and ions leads to the different evolution of the plasma temperature at the sheath entrance.The sheath entrance is found when the quasineutral condition is broken (ρ ̸ = 0).Figure 4(a) displays the temporal evolution of the electron (red) and ion (blue) temperatures at the sheath entrance.The electron temperature T e dramatically increases when ELM burst commences.It reaches the maximum value of T e,m = 0.4 keV at about t = 0.15 ms and then quickly decreases.The ion temperature T i starts to gradually increase at t = 0.125 ms, reaches its peak value T i,m = 0.2 keV at t = 0.4 ms, and then decreases to 20 eV (the temperature before ELM).
Seven characteristics time points (t 0 -t 6 ) are considered to conveniently describe the energy evolution.The time point t 1 represents the time when T e reaches T e,m , i.e.T e (t 1 ) = T e,m .When t = t 3 , T e (t 3 ) = T i (t 3 ).The time point t 2 is located between t 1 and t 3 .At t 4 , T i reaches the maximum value T i,m , i.e.T i (t 4 ) = T i,m .T i (t = t 5 ) = 0.5T i,m .The impacting energy is calculated to determine the time when the ELM-induced sputtering stops, as shown in figure 4(b).K i is the estimated impacting energy of ions, i.e.
The time point t 6 is defined as K i (t = t 6 ) = Eth, where Eth is the threshold energy of physical sputtering.The period between t 0 and t 3 is called the electron-dominating phase (EDP) because T e > T i .On the other hand, the phase between t 3 and t 6 is called the ion-dominating phase (IDP) because T e < T i .The period after t 6 is called the post-ELM phase (PEP).In PEP, since the impacting energy is less than the threshold energy, no sputtered particles are generated in principle.Thus, PEP is not analyzed herein.

Effect of plasma sheath on particle trajectory and velocity
As T e and T i change with time, the electrostatic potential, which is closely related to them, dramatically changes with time.Figure 5 compares the electrostatic potential at t = t 1 (EDP) and t = t 4 (IDP).The electric field points to the material surface near the mono-block, as denoted by the black arrows.Two main differences exist between the electric fields in EDP and IDP.The first difference is that the potential drop Φ w in EDP is significantly greater than that in IDP because it is mainly determined by T e .The second difference is the potential profile at GE.The electrostatic field gradually decreases at GE (namely, E z < 0 at GE, as the red arrows show in figure 5(a)) in EDP.In contrast, a potential peak is present at GE in IDP, as shown in figure 5(b).This peak causes the z-component of the electric field E z to change its sign when z > 0 and z < 0.More specifically, E z (z > 0) > 0 (as shown by the blue arrows) and E z (z < 0) < 0 (as shown by the red arrows).The potential peak becomes an energy barrier and changes the trajectory of the charged particles, as illustrated in figure 6.The peak stems from the difference in the Larmor radius of electrons and ions [25].Electrons are shadowed by SS and cannot reach GE due to the small Larmor radius.In contrast, ions have a large Larmor radius and can circumvent SS and enter the gap.Thus, a positive charged space forms at GE.The positive charging effect becomes more significant as the Larmor radius of ions increases.The comparison of the ion trajectories with initial parallel velocity v ∥ = 116 km s −1 and perpendicular velocity v ⊥ = 100 km s −1 in E(t 1 ) (EDP) and E(t 4 ) (IDP) are shown in figures 6(a) and (c), respectively.The figure also shows the velocity evolution with time.In EDP, the z-component of velocity (v z ) gradually decreases with time in E(t 1 ) in (b), implying that E z is almost constant.Due to the negative E z , ions impact the lower part of CDG (a).The velocity in the x-direction v x slowly increases before t < 13 µs and dramatically decreases after 13 µs.This is because when ions are close to the gap surface, the negative E x strongly attracts them.On the other hand, the abstract value of v z does not increase after 13 µs in E(t 4 ) because the E z at GE is positive, and ions are repelled by a positive E z .Thus, ions impact HT in IDP (c).

Incident deuterium flux distribution on CDGs.
Based on the two-component assumption, the deuterium (D) flux is composed of hot and cold components.By applying the leastsquare method, the temperature and weight of hot and cold ions can be obtained, which are listed in table 1.
The weight ω h denotes the percentage of hot ions in the total flux, and the percentage of cold ions is ω c = 1 − ω h .At t = t 1 , ω h = 0, implying that no hot ions have reached the target.Due     to the Coulomb collsition between electrons and ions, the temperature of cold ions increases to 0.03 keV.When t = t 2 , the faster particles of the hot ions reach the target.Hence, the temperature of hot ions is as high as 0.8 keV.However, the ratio of hot ions is very small (only 1.4%).Thus, the average temperature at the sheath entrance is T is = 0.06 keV.With time, more hot ions reach the target and the percentage increases, leading to the increase of T is .In contrast, the temperature of hot ions decreases.When t = t 4 , the percentage of hot ions reaching the maximum value is ω h = 72%, and T h = 0.24 keV and T c = 0.08 keV.The temperature of cold ions increases with ω h and T h .The could bedue to Coulomb collision.h,D (t = t 1 ) ∼ 0. So, the mechanism of the first peak is the increase of the electrostatic potential, which accelerates ions.Moreover, ions can reach deeper parts of the gap.This is because the electric field at GE drags ions to the gap, which is illustated in figure 6(a).The Γ D distribution at t = t 2 is similar to that at t 1 and Γ D, HT and Γ D, VT are 0.6 × 10 24 m −2 s −1 and 0.5 × 10 24 m −2 s −1 , respectively.
When t ⩾ t 3 (i.e., the transit time of EDP and IDP), the flux distribution changes.The cold ion flux drops significantly to a low value (∼0.2 × 10 24 m −2 s −1 ) and the hot ion flux dramatically increases.The values on HT are 0.5 × 10 24 , 0.7 × 10 24 , and 0.4 × 10 24 m −2 s −1 at t 3 , t 4 , and t 5 , respectively.While the hot ion flux are 0.1 × 10 24 , 0.5 × 10 24 , and 0.2 × 10 24 m −2 s −1 at t 3 , t 4 , and t 5 , respectively.Thus, the second flux peak at t 4 is mainly caused by hot ions.The total deuterium flux Γ D at VT is significantly greater than that on HT.The values of Γ D, VT are 10 24 , 1.8 × 10 24 , and 10 24 m −2 s −1 at t 3 , t 4 , and t 5 , respectively.The hot ion flux on the VT is approximately equal to Γ D, VT .In other words, the high flux on the VT is mainly caused by the hot ions.The decay length λ w is ∼ 0.7 mm.This value is greater than that at t 0 because of the greater Larmor radius of ions compared to electrons.In EDP, the contributions of hot ions is quite minor.In contrast, the percentage of hot ions is close to or higher than 50% in IDP.
Two main differences exist between the deuterium flux Γ D in EDP and IDP on CDG.First, the D flux on the top surface in EDP is close to or even greater than that on VT.In contrast, Γ D, HT is significantly less than Γ D, VT in IDP.Second, the decay length λ w in EDP is close to the gap depth (1 mm), and in IDP, λ w ∼ 0.7 mm.

Impacting velocity and angle distribution on CDG.
The impacting energy and angle play important roles in the physical sputtering, as the figure 8 shows.The horizontal axis is the impacting energy in keV while the vertical axis is the impacting angle in degree.The black contour lines with number represent the sputtering yield.The red, purple, blue, and green lines represent the evolution of impacting energy and angle of hot ions on HT, hot ions on WS, cold ions on HT, and cold ions on WS, respectively.The O signs denote the start point, and the arrows display the evolution directions.The yellow rectangle shows the region that the impacting energy is less than the sputtering threshold.
The impacting energy of cold ions (green & blue lines) at t = t 0 is 0.08 keV, which is less than the sputtering threshold energy Eth = 0.21 keV.Thus, the majority particles before ELM burst cannot produce sputtering atoms.In present simulation, the maximum impacting energy is 0.31 keV, so some high energy ions can cause sputtering.Therefore, at t = t 0 the sputtering exists but is weak.The impacting energy of cold ion jumps to 1.53 keV at t = t 1 when serious sputtering is produced.After that, the impacting energy drops monotonously and returns to 0.08 keV at the end of ELM burst.
The impacting energy of cold ions on the HT and WS is the same because the impacting energy is determined by the potential drop and irrelevant with the landing point.By contrast, due to the difference in the glancing angles, the The evolution of impacting energy and angle of hot and cold ions on HT and WS and the dependence of the sputtering yield Y on the impacting energy and angle.The horizontal axis is the impacting energy in keV while the vertical axis is the impacting angle in degree.The black contour lines with number represent the sputtering yield.The red, purple, blue, and green lines represent the evolution of impacting energy and angle of hot ions on HT, hot ions on WS, cold ions on HT, and cold ions on WS, respectively.The O signs denote the start point, and the arrows display the evolution directions.The yellow rectangle shows the region that the impacting energy is less than the sputtering threshold.most-possible impacting angle of cold ions on WS and HT is ∼5 • and 55 • at t = t 0 , respectively.In EDP, the ions are greatly accelerated by the strong electric field, and the impacting angle on WS increases while that on HT decreases.The sputtering yield increases from zero to around 0.02.In IDP, the impacting angle returns to the initial values.Although the impacting angle significantly changes, the sputtering yield is mainly influenced by the impacting velocity because it is not sensitive to the impacting angle in this region (∂Y/∂θ ∼ 0).
The impacting energy of hot ions at t = t 2 is very high (∼4 keV) and then drops monotonously with time, as the red and purple lines show.Due to the difference in the impacting angles (θ(t 0 ) = 52 • on HT and θ(t 0 ) = 28 • on WS), the sputtering yields are different: Y h,WS ≈ 0.021 and Y h,HT ≈ 0.028.The sputtering yield on WS is smaller than that on HT.The impacting energy decreases with time, but the impacting angle on HT slightly increases meanwhile the value of θ on WS decreases.Both Y h,WS and Y h,HT decreases with time.Y h,WS drops to 0.018 and 0.012 at t = t 3 and t = t 4 , respectively.Meanwhile Y h,HT ≈ 0.024 at t = t 3 and Y h,HT ≈ 0.015 at t = t 4 .The decrease speeds of Y h,WS and Y h,HT are similar.The minimum value is 0.17 keV, which is smaller than E th .This means that the hot ions cannot produce sputtering at the end of ELM burst.

Sputtered tungsten flux evolution during the ELM burst.
Figure 9 shows the temporal and spatial evolutions of tungsten flux Γ W (t, L).The horizontal and vertical axes represent time and distance, respectively.HT, WS, GB, and SS are denoted in sequence on the right.The tungsten distribution is obtained by the cubic interpretation of the time points from t 0 to t 6 .The tungsten flux is zero between t 6 and t = 3 ms (PEP).EDP, IDP, and PEP are also denoted.The tungsten flux evolution is similar with the deuterium impacting flux (figure 7(c)).There are two peaks on both HT and VT.The peaks are at t 1 and t 4 .
The sputtering flux at the VT is extremely high at t 4 .
The sputtered tungsten source distribution in the pre-ELM period (t 0 ) is almost zero because the impacting velocity is less than the threshold velocity, as illustrated preciously in figure 8.In EDP t = t 1 , the tungsten flux on the top surface is Γ W, HT = 1.3 × 10 22 m −2 s −1 , while Γ W, VT decreases to 1.1 × 10 22 m −2 s −1 and exponentially decreases to the GB.The decay length of tungsten sputtering is λ w = 1 mm.This distribution is consistent with the deuterium flux distribution at the same time.At t = t 2 , Γ W, HT ≈ Γ W, VT ≈ 1.0 × 10 22 m −2 s −1 .The sputtered tungsten flux on HT is higher or close to that on VT between t 0 and t 2 (i.e., EDP).At  EDP, IDP, and PEP.At the first peak (t = t 1 ), the tungsten flux generated by cold ions is equal to the total tungsten flux, namely, Γ c,W (t 1 ) = Γ W (t 1 ).This means the all the sputtering particles are generated by cold ions.The case is the same at t = t 2 .After t = t 3 , the contribution of cold ions to Γ W significantly decreases.However, a long tail of Γ W caused by the cold ions exists until t 6 .At about t = t 3 , the contribution of hot ions on Γ W quickly becomes dominant and yields the second peak.At t 3 , the tungsten flux Γ c,W = 0.08 × 10 22 m −2 s −1 and Γ h,W = 0.2 × 10 22 m −2 s −1 on HT.The percentage of the contribution of hot ions is 71%.Meanwhile, Γ c,W = 0.2 × 10 22 m −2 s −1 (24%) and Γ h,W = 0.5 × 10 22 m −2 s −1 (76%) on WS.After that, the contribution from cold ions drops and the contribution of hot ions increases.The contributions of hot ions are 94% and 90% at t 4 and t 5 on HT, respectively, and are 93% and 90% at t 4 and t 5 on VT, respectively.This means the contribution of hot ions is dominant in the IDP.Based on these results, the mechanism of the two sputtering peaks is clear.The first peak is generated by the cold ions that are accelerated by the strong electric field.Hot ions play a very minor role in the first peak.The second peak stems from the energetic ions from ELM burst.
The fluence in EDP and IDP is defined as the integral from t 0 to t 3 and from t 3 to t 6 , respectively, Figure 11 presents the of the tungsten sputtering fluences F W of hot and cold ions in EDP and IDP on HT and WS.F W induced by cold ions in the EDP and IDP are labelled by red and green, respectively while tungsten fluences by hot ions in the EDP and IDP are labeled by blue and purple, respectively.The contributions on HT and WS are significantly different.On the HT, the first sputtering peak stemming from the cold ions in EDP contributes to 43% (F c,EDP ) of the tungsten sputtering.The sputtering tail generated by the cold ions in IDP (F c,IDP ) is small (4%).The highest contribution is from the hot ions in IDP (50%) (F h,IDP ).The contribution of hot ions in EDP (F c,EDP ) is negligible (1%) as very few hot ions are present in EDP.
The contribution of hot ions in IDP (F h,IDP ) on WS is dominant (60%).In contrast, the contribution of cold ions is relatively small (5% in IDP and 32% in EDP).However, in deeper parts of the gap, the contributions of hot and cold ions in IDP decays very quickly due to the electrostatic sheath potential characteristics.Thus, tungsten sputtering by cold ions in EDP becomes important on the lower part of the CDG.
The total erosion on the HT and VT is 5.9 × 10 18 (m 2 shot) −1 and 9.8 × 10 18 (m 2 shot) −1 , respectively.Namely, the erosion rate on the VT is much serious than that on the HT.The erosion rate (the erosion depth per shot) is 0.16 nm shot −1 on the VT.This value is greater than previous estimation [26].By assuming that the ELM frequency is 50 Hz, the time of production a 1 mm damage on VT is only 1.2 × 10 5 s (∼35 d).Thus, the erosion effect on the gap edge is a serious issue.The tungsten source S W (in 10 20 s −1 ) is obtained to compare the contributions of tungsten from HT and WS.The flow is obtained by integrating Γ W from point A to B and B to E: The factor 2 mm comes from the distance between strike point and the landing position of the flux [22].The factor 12 is introduced because the mono-block dimension in the toroidal direction is 12 mm and the length between A and B is 1 mm.The factor 2π R CDG /12.5mm is the number of times the mono-block and CDG repeat in the toroidal direction.Here, an implicit assumption is made that the ELM flux on the entire toroidal circumference is the same.It is approximately true in a relatively long time and many ELM bursts are triggered.The contributions by hot and cold ions on HT and CDG are: The contributions of hot and cold ions to the tungsten source on HT and CDG are shown in figure 12.As the size of HT is considerably greater than the CDG width, the contribution from HT is dominant.The cold ions contribute to about 46% of the tungsten production and hot ions contribute to 50%.In contrast, the contributions of cold and hot ions on CDG are only 1% and 2%, respectively.This means that the global contribution of the CDG is quite minor, although the physical sputtering on the VT is serious.

Influence of the peak temperature on the sputtering yield
The peak temperature T peak during ELM burst determines the amount of sputtered tungsten particles.In this study, five different peak temperatures in ELM burst are selected T peak = 0.6, 0.8, 1.0, 1.2, and 1.6 keV.The tungsten flux evolution on HT and VT are shown in figures 13(a) and (b), respectively.The tungsten flux evolution is similar under different peak temperatures.Two peaks of Γ W exist on both HT and VT.On HT, the first peak is greater than the second peak.In contrast, on VT, the second peak is significantly greater than the first peak.This means that the tungsten generation process is quite similar under different peak temperatures.The benchmark of the numerical modeling of ELMinduced sputtering is performed.As the temporal detail of tungsten flux evolution is not provided by EAST, the effective sputtering yield is compared.The experimental effective sputtering yield Y EES is defined from the ratio of the gross tungsten flux (by the tungsten spectrum) to the saturation ion flux (from Langmuir probes) with shoot number #69634, 63 884, 63 885, and 63 887 [14].The numerical effective sputtering yield (NES) Y NES is defined in a similar way: Please note that Y NES is measured by gross impurity production.Y EES and Y NES are shown in figure 14 with red and blue points, respectively.The experimental effective sputtering yield Y EES is a function of T ped , which is the pedestal temperature in the steady state, as shown in the sub-figure (i) in figure 14.The pedestal temperature is around several hundred eV.The numerical effective sputtering yield (Y NES ) is a function of T peak , which is the peak value of T SOL (the plasma temperature in the SOL region), as shown in the sub-figure (ii) in figure 14. T SOL is around 50 eV in the steady state.No available experimental data for T peak from EAST at present but the estimated value of T peak is around several hundred eV.The power of EAST is lower than the power of ITER, so we chose relatively higher values of T peak for future experiments.As the sputtering flux on the VT is significantly greater than Both the Y NES, VT and Y NES, HT obey the parabolic relation with T peak , which implies that the sputtering becomes more serious with T peak .The study results also provide a possible explanation for the unexplained phenomenon that some experimental effective sputtering yields are greater with lower pedestal temperature.For example, Y EES = 0.089% with T ped = 0.37 keV in #63 885 and Y EES = 0.074% with T ped = 0.46 keV in #69634.This could be because the ELM filament lands on HT in #69634, producing a low tungsten sputtering yield.Meanwhile, part of the ELM filament touches the vertical edge of CDG and generates more tungsten impurities.It is because the data (Y EES = 0.089% with T ped = 0.37 keV in #63885, the ∆ point) are close to the fitting line of Y NEP on VT (the square points).

Influence of re-deposition and self-sputtering
The effects of re-deposition and self-sputtering should also be considered in the ELM-induced sputtering, as shown in figure 15  The initial velocity of sputtered tungsten atoms is randomly selected in Maxwellian distribution.The velocity direction is uniformly selected.The averaged energy of sputtered W atoms is ĒW = 12 eV from TRIM calculation.ĒW does not strongly depend on the impacting energy and angle of deuterium ions, so ĒW is the same for different impacting energy and angle.The ionization of tungsten atoms is modelled by using Monte Carlo method.The mean-free-path of ionization (W atom to W+) is λ ∼ 1 mm.The ionization transit time is τW,ionz ∼ λ/ √ ĒW /m W = 0.6 µs.Tungsten ions are traced by using the 4th order Runge-Kutta method.The transit time of ionization to redeposition is τW, redep ∼ 7.5 µs.Thus, the transit time of sputtering to redeposition is τW = τW,ionz + τW, redep ∼ 8 µs.Therefore, the sputtering and redeposition occur almost simultaneously because τtotal is significantly smaller than the duration of ELM burst (∼0.3 ms [22]).
There are two peaks of redeposition flux Γ redep on HT, which is similar to the temporal evolution of tungsten sputtering flux Γ W .The reason of this similarity is that the transmit time τW is short compared to the flux evolution time.The first and second peaks of Γ redep are 1.33 × 10 22 m −2 s −1 and 0.59 × 10 22 m −2 s −1 , respectively.These values are close to the peaks of Γ W (1.35 × 10 22 m −2 s −1 and 0.56 × 10 22 m −2 s −1 , respectively).This means the redeposition rate on HT R redep (HT) = Γ redep (HT) /Γ W (HT) is high.Γ redep at t 4 (0.59 × 10 22 m −2 s −1 ) is greater than the counterpart Γ W (0.56 × 10 22 m −2 s −1 ).The reason is that the sputtered atoms from VT deposit on HT.
The total redeposition rate on CDG are R redep (CDG) = 95%, 94%, 88%, and 86% at t = t 1 , t 3 , t 4 , and t 6 , respectively.This result is basically consistent with the results from DIII-D [27].The high redeposition rate (95%) is due to the higher electrostatic field in the EDP.In IDP, the electrostatic field becomes weaker and R redep (CDG) drops.There is a significant difference between the net erosion Γ net = Γ W − Γ redep and Γ W due to the high redeposition rate.The peak net erosion on HT is Γ net (HT) ∼ 0.16× 10 22 m −2 s −1 in the EDP (at t = t 1 ), which means the erosion is reduced about 1 − Γ net /Γ W ∼ 88% by redeposition.VT is eroded while SS is deposited.
The self-sputtering of W ions Γ self is shown in figure 15(c).The most intense self-sputtering is on HT.The reason is that the redeposition flux on HT is high.In EDP, the self-sputtering W flux on SS and VT is similar.By contrast, Γ self (SS) is significantly greater than Γ self (VT) in IDP.The reason is that in EDP, the sputtering on VT is not very strong and the contribution of the sputtering flux from VT on Γ redep (SS) is not very high.In IDP, the sputtering on VT is strong, which enhances Γ redep (SS) and Γ self (SS) rises.
There are two peaks of Γ self .The first peak is significantly greater than the second one on the CDG.The high selfsputtering yield in the first peak is caused by the high impacting energy of W ions K W .The averaged K W are 304 eV, 274 eV, 191 eV, 104 eV, 97 eV and 95 eV at t = t 1 , t 2 , t 3 , t 4 , t 5 , and t 6 , respectively.

Conclusion and summary
In this study, multiscale simulation was performed and the following conclusions were obtained.
Character of ELM transport and its influence.The ELM transport was modeled via macroscopic 1D PIC simulation and the energy evolution of electrons and ions was investigated.Due to the difference in the velocities of the electrons and ions from ELM burst, hot electrons arrive at the target earlier than hot ions.In this phase, the electron temperature is higher than the ion temperature.Thus, this phase is called EDP.When hot ions reach the target, most hot electrons have already been absorbed by the target.Hence, the ion temperature is greater than that of electrons, and this phase is called IDP.The ELM transport influences the temporal evolution of the plasma temperature and density at the sheath entrance, which are critical factors for investigating the impacting velocity and angle.
Sheath evolution and its influence on sputtering.The electrostatic sheath potential was simulated in the microscopic step via 2D PIC simulation.The results show that the sheath potential profiles are different in EDP and IDP, which significantly affects the dynamics of ions.The impacting velocity is determined by the initial velocity (the velocity of ions entering the sheath entrance) and the potential drop.In EDP, ions are strongly accelerated by the sheath electric field because of the high electron temperature.In IDP, the effect of potential drop is less important because the initial velocity is large.Moreover, the sheath potential influences the impacting angle, and the influence is different on the mono-block HT and VT.
The answer to question one (How does ELM transport influence the sputtering process and evolution as the velocity of electrons and ions from ELM burst is different?)ELM transport results in a double-peak distribution of the sputtered tungsten flux.The first peak is in EDP and is produced by the cold ions accelerated by the strong electrostatic sheath potential.The electrostatic sheath potential is generated by the high-temperature electrons from ELM burst.The second peak is directly caused by the hot ions from ELM burst.
The answer to the question two (What is the influence of CDG on tungsten sputtering?)The physical sputtering on the WS of CDG is different from that on the mono-block HT due to the difference in the impacting angle.In EDP, the sputtered tungsten production on CDG is less than that on HT and the number of tungsten particles from WS is dramatically greater than that from the HT in IDP.The tungsten fluence from CDG is about twice of that from the top surface, implying that the erosion on CDG is much more serious than that on HT.However, integration over the entire mono-block shows that the total amount of tungsten production from CDG sputtering is only 3% of the overall tungsten production because the sputtering on CDG is very local and the area of the HT is much greater than that of CDG.
A simple application of this study.The sputtering yield under different ELM peak temperatures (0.6, 0.8, 1.0, 1.2, and 1.6 keV) was also investigated.The study results are comparable to the EAST experimental data and show a qualitative consistency.This study provides a possible explanation of an experimental phenomenon.The influences of redeposition and self-sputtering are also considered.The redeposition rate on the HT is high (∼90%) and the erosion on the HT surface is largely reduced by redeposition.By contrast, the redeposition on the vertical surface is relatively low and the net erosion is significant.Meanwhile, the SS is net deposited.The impurity production by self-sputtering is about 13% of the amount of tungsten particles sputtered by deuterium ions.Thus, the selfsputtering is not the major mechanism responsible for ELMinduced sputtering.

Figure 2 .
Figure 2. (a) The CDG geometry in the simulation.The plasma flux is denoted by the red arrow with glancing angle Ψ = 5 degree.The cardinal coordinates of the points A, B, C, D, and E are (x, z) = (−1, 0), (0, 0), (0, −1), (0.5, −1), and (0.5, 0), respectively.The horizontal top (HT), wetted surface (WS), gap bottom (GB), shadowed surface (SS) corresponds to the segment AB, BC, CD, and DE, respectively.The small segment at the top of WS is defined as vertical top (VT).(b) the impacting angle of two sample particles.The arrows represent the velocity direction and the dash lines are the normal directions.

Figure 3 .
Figure 3. Particle positions and energy distributions at different times.The green and cyan dots represent low-energy electrons and ions, respectively.The yellow and violet dots represent high-energy electrons and ions, respectively.The red and blue lines indicate the kinetic energy profiles of electrons and ions, respectively.

Figure 4 .
Figure 4. (a) Time evolution of electron and ion temperatures at the sheath entrance and (b) the estimated impacting energy of the ions.

Figure 5 .
Figure 5. Electrostatic field directions at (a) t = t 1 (EDP) and (b) t = t 4 (IDP).The electrostatic potential profiles are also shown in the corresponding sub-figures.The arrows denote the direction of the electric field E.

Figure 6 .
Figure 6.Trajectories and velocity components of ions with parallel velocity of 116 km s −1 and perpendicular velocity of 100 km s −1 in the electrostatic field at t = t 1 ((a) and (b)) and t = t 4 ((c) and (d)).In figure (a) and (c), the blue line represents the coordinate of the sampled ion.Figure (b) and (d) show the x, y, and z component of velocity.
Figure 6.Trajectories and velocity components of ions with parallel velocity of 116 km s −1 and perpendicular velocity of 100 km s −1 in the electrostatic field at t = t 1 ((a) and (b)) and t = t 4 ((c) and (d)).In figure (a) and (c), the blue line represents the coordinate of the sampled ion.Figure (b) and (d) show the x, y, and z component of velocity.

Figure 7 .
Figure 7.The deuterium flux distribution on the CDG as a function of time and space.The unit of deuterium flux is 10 24 m −2 s −1 .The sub-figures (a), (b), and (c) represent the flux of cold ion flux, hot ion flux, and total flux, respectively.The horizontal axis is time in ms.Vertical axis represents the curve distance from point A on the CDG.The geometry of the CDG is shown in the attached figure and the points A to E corresponds to the ones on the vertical axis.The surfaces of the CDG are also labeled on the right.

Figure 7
figure and the points A to E corresponds to the ones on the vertical axis.The surfaces of the CDG are also denoted on the right.At t = t 0 , the total deuterium flux (Γ D ) is equal to the flux of cold ion (Γ c,D ) and Γ h,D is zero.Γ D on HT is 0.1 × 10 24 m −2 s −1 , while that on the VT is Γ D, VT = 0.3 × 10 24 m −2 s −1 .Γ D, VT > Γ D, HT because the glancing angle on VT is ψ B,HT = 5 • while that on WS is ψ B,WS = 1 − ψ B,HT = 85 • .The Γ D on WS exponentially decreases as the depth increases, and Γ D ∼ 0 at the depth z ∼ −0.3 mm.The depth at which Γ D (z) ∼ 0 is defined as the decay length λ w .There are two flux peaks of total D flux on the CDG.The first one appears at t = t 1 and the second one is at t = t 4 .When t = t 1 = 0.15 ms (EDP), the total D flux on HT significantly increases to 0.8 × 10 24 m −2 s −1 and the flux is 0.6 × 10 24 m −2 s −1 on the VT.The cold ion flux is almost equal to total flux, namely, Γ c,D (t = t 1 ) ≈ Γ D (t = t 1 ).Meanwhile

Figure 8 .
Figure 8.The evolution of impacting energy and angle of hot and cold ions on HT and WS and the dependence of the sputtering yield Y on the impacting energy and angle.The horizontal axis is the impacting energy in keV while the vertical axis is the impacting angle in degree.The black contour lines with number represent the sputtering yield.The red, purple, blue, and green lines represent the evolution of impacting energy and angle of hot ions on HT, hot ions on WS, cold ions on HT, and cold ions on WS, respectively.The O signs denote the start point, and the arrows display the evolution directions.The yellow rectangle shows the region that the impacting energy is less than the sputtering threshold.

Figure 9 .
Figure 9. Temporal and spatial evolutions of the tungsten flux during ELM burst.The color scale shows the tungsten flux in units of /m 2 s.The horizontal and vertical axes represent time and distance, respectively.HT, WS, GB, and SS are denoted in sequence on the right.EDP, IDP, and PEP are also denoted.
t 3 = 0.3 ms, the relation between Γ W, HT and Γ W, VT changes: Γ W, HT = 0.3 × 10 22 m −2 s −1 and Γ W, VT = 0.7 × 10 22 m −2 s −1 .Γ W, VT is significantly greater than Γ W, HT .This relation does not change in the entire IDP (t 4 , t 5 , and t 6 ).Γ W, VT is 1.6 × 10 22 and 0.3 × 10 22 m −2 s −1 at t 4 = 0.42 ms and t 5 = 0.75 ms, respectively.Moreover, Γ W, HT is 1.4 × 10 22 and 0.2 × 10 22 m −2 s −1 at t 4 and t 5 , respectively.The numbers on the color scale denote the tungsten flux values at the counterpart time numbers (t 1 -t 5 ).The upper and lower rows represent the tungsten flux on VT and HT, respectively.The contributions of hot and cold ions on the tungsten flux on HT and VT are shown in figures 10(a) and (b), respectively.The red and blue areas represent the values of the tungsten impurity generated by the hot and cold impacting deuterium particles, respectively.The black line is the total amount of sputtered tungsten.The dashed lines separate the

Figure 10 .
Figure 10.Comparison of the temporal evolution of the flux of sputtered tungsten due to hot and cold ions on (a) HT and (b) VT.

Figure 11 .
Figure 11.Fluence of tungsten impurity on CDGs due to hot and cold ions in EDP and IDP.

Figure 12 .
Figure 12.Tungsten source S as a function of time and the contributions of hot and cold ions on HT and CDG on the total source.

Figure 13 .
Figure 13.Tungsten flux on the (a) MB horizontal top (HT) and (b) vertical top (VT) of the gap edge with different peak temperatures.

Figure 14 .
Figure 14.Effective tungsten sputtering yield on HT (squares) and VT (circles) under different peak temperatures.The experimental and numerical effective sputtering yield are denoted by red and blue points, respectively.The X, +, ∆, and ∇ points show the experimental data in the # 69634, 63 884, 63885, and 63 887.The circle and square points show the numerical effective sputtering yield on the VT and HT, respectively.

Table 1 .
Weights and temperatures of hot and cold ions.