Analysis of equilibrium and turbulent fluxes across the separatrix in a gyrokinetic simulation
The SOL width is a parameter of paramount importance in modern tokamaks as it controls the power density deposited at the divertor plates, critical for plasma-facing material survivability. An understanding of the parameters controlling it has consequently long been sought (Connor et al. 1999 NF 39 2). Prior to Chang et al. (2017 NF 57 11), studies of the tokamak edge have been mostly confined to reduced fluid models and simplified geometries, leaving out important pieces of physics. Here, we analyze the results of a DIII-D simulation performed with the full-f gyrokinetic code XGC1 which includes both turbulence and neoclassical effects in realistic divertor geometry. More specifically, we calculate the particle and heat fluxes along the separatrix, discriminating between equilibrium and turbulent contributions. We find that the density SOL width is impacted almost exclusively by the turbulent electron flux. In this simulation, the level of edge turbulence is regulated by a mechanism we are only beginning to understand: -drifts and ion X-point losses at the top and bottom of the machine, along with ion banana orbits at the low field side (LFS), result in a complex poloidal potential structure at the separatrix which is the cause of the drift pattern that we observe. Turbulence is being suppressed by the shear flows that this potential generates. At the same time, turbulence, along with increased edge collisionality and electron inertia, can influence the shape of the potential structure by making the electrons non-adiabatic. Moreover, being the only means through which the electrons can lose confinement, it needs to be in a balance with the original direct ion orbit losses to maintain charge neutrality.
Studies of the plasma edge are indispensable for a variety of tokamak operation and performance issues such as the density and power scrape-off-layer (SOL) widthschang2017gyrokinetic (); meier2017drifts (); meier2016analysis (); myra2011reduced (); russell2015modeling (); pankin2015kinetic (); rozhansky2018structure (); reiser2017drift (), intrinsic rotation and momentum transportstoltzfus2012tokamak (); stoltzfus2012transport (); muller2011experimental (); loizu2014intrinsic (); groebner2009intrinsic (); chang2008spontaneous (); seo2014intrinsic (), edge flowschankin2007possible (); hoshino2007numerical (); kirnev2005edge2d (); labombard2004transport (); pigarov2008simulation (); pitts2005edge ()and the L-H transitionaydemir2012pfirsch (); shaing1989bifurcation (); chang2017fast ().
Out of all the above mentioned important and inter-related facets of edge physics, the SOL width is crucial in tokamak operation as it governs the power density delivered at the divertor plates. Indeed, a viable fusion reactor would need to have plasma-facing materials that can withstand large heat loads without the need of regular replacement. However, an ITER exhaust power of , concentrated on a speculated SOL width would produce a whopping of parallel power density, placing stringent constraints on the lifetime of the divertor materialsgoldston2015theoretical (). It is therefore of utmost importance to be able to predict the SOL width for a given discharge based on its parameters but also to gain an understanding of the physical mechanisms that control its size. To this end, Eicheich2013scaling () has constructed scaling relations for the power SOL width, , based on multiple regressions of relevant parameters from various experimental discharges. He has arrived at a scaling law of the form , where is the poloidal magnetic field, based on which the SOL width of ITER is inferred to be . Goldstongoldston2011heuristic () has provided theoretical justification for the above scaling by devising a heuristic model that predicts the SOL width. In summary, this model rests on the assumption that the perpendicular and curvature drift ion flows are balanced exactly by parallel flows with one half of the flux being directed to the divertors at a speed of and the other half returning upstream via the Pfirsch-Schlütter current. This mechanism is responsible for setting the density SOL width, . The next assumption is that this density channel fills with turbulent electrons and empties through Spitzer-Härm conductivity, establishing the power SOL width. The heuristic model correctly captures the inverse relation between and proposed by Eich and can thus, with reasonable accuracy, reproduce the experimental results.
Until recently, simulations of the edge have relied upon reduced fluid models in simplified geometries omitting essential physics. This changed with the introduction of XGC1chang2009compressed (); ku2009full (); ku2016new (); churchill2017pedestal (); chang2017fast () which is a full-f, Particle-In-Cell, 5D gyrokinetic code using realistic X-point geometry, optimized to simulate the tokamak edge. The code has gyrokinetic ions and offers an option for drift-kinetic, fluid or adiabatic electrons. Employing it in a study of ITER SOL width, Chang et.alchang2017gyrokinetic () discovered that the simulated width is six times larger than the one predicted by the Eich fit or the Goldston heuristic model, casting doubt on whether the heretofore assumed scaling relations continue to be valid at the ITER regime. The explanation that they provide regarding this disparity has to do with the importance of edge electron turbulenced2011convective (); myra2015turbulent (); myra2016theory () in the new regime. More specifically, they claim that the electrostatic-‘blobby’ edge turbulence is setting the scale of the electron channel. The size of a blob or a streamer (or any other coherent turbulent structure) is, at least partially, related to the width of the edge shearing layerfurno2008experimental (); bisai2005formation (), wherein the blobs are being generated. In present day tokamaks, this shearing layer width is relatively small and the density channel is dominated by ion neoclassical flows inside the SOL. On ITER however, due to the much larger size of the shearing layer, the radial extent of blobs and streamers is bigger and turbulence dominates the density width.
From the above discussion, it is clear that important figures of merit such as the SOL width are the result of a non-linear interaction between particle orbit effects (e.g., neoclassical effects and ion X-point losseschang2002x (); ku2004property ()) with the electrostatic-‘blobby’ edge turbulence. The interaction is through the radial and poloidal edge electric field, established by the ion losses and the sheared radial flows that it generates which can influence both the turbulence and the ion flows. Such interplay can often yield unexpected results and thus merits further study with state-of-the-art codes which might illuminate certain aspects of this complex process. Here, we analyze the results of an XGC1 simulation of DIII-D from the viewpoint of equilibrium and turbulent particle and heat fluxes across the last closed flux surface (LCFS). Within this framework, we attempt to understand the relationship between the direct ion orbit losses and the electrostatic potential that drives these fluxes as well as the effect of the generated shear on the turbulence. Furthermore, we try to measure the impact of turbulence on the electron transport and thereby, test the second assumption of the heuristic model which states that only turbulent electrons are responsible for filling the density SOL channel.
The organization of the paper is as follows: In Section II we give a brief description of the simulation, it’s geometry and some relevant parameters. In Section III we provide the definitions of the measured fluxes and the methodology of extracting them from the XGC1 data. We present the equilibrium and turbulent particle and heat fluxes across the separatrix and interpret the observed patterns using the notions of ion X-losses, ion banana drifts, non-adiabaticity, quasineutrality and shear flows. In Section IV we perform numerical integrations to assess the size of the different fluxes and extricate the contribution of each to the cross-separatrix transport and hence, to the SOL width. In Section V we summarize and draw conclusions.
In this section we are going to provide a description of the simulation domain and parameters. The simulation is of the DIII-Dluxon2002design () discharge # 144981, which is an H-mode heated by neutral beam injection. It is initialized with experimental profiles taken at time . The simulation inputs include experimental kinetic profiles of electron density and temperature ( and ), ion temperature (), and magnetic equilibrium, from a kinetic EFIT magnetic reconstruction. In Fig. 1 we see a poloidal cross-section of the simulation domain along with the directions of the unit vectors of the cylindrical coordinate system used in this paper. The origin of the coordinate system in the poloidal cross-section will be taken to be the point . The geometry is lower single null but with a secondary “virtual” upper X-point (outside of the simulation domain). The magnetic field is in the negative direction making the ions -drift towards the X-point and the electrons towards the top. The black dashed line denotes the separatrix. The yellow lines are closed and open flux surfaces and the black and red continuous lines are showing the poloidal projection of the banana orbits (black for ascending-red for descending parts of the orbit). Those have been plotted using contours of toroidal angular momentum. The closed dashed black line close to the center is a potato orbit. The orientation of the coordinate system and the toroidal magnetic field direction are also shown. In the rest of the paper, when we refer to the poloidal angle , of a separatrix point, we will mean the angle formed between the line passing through that point and the line passing through the low-field side midplane with both lines passing through the origin. Positive angles are in the counter-clockwise direction.
The total simulation time is with a time-step of . Although this time is not enough for the core turbulence to saturate it is adequate for the saturation of the edge turbulence. The poloidal magnetic field is (measured at the outboard midplane, in the positive direction) and the toroidal current is . Equilibrium plasma density at the edge was measured to be with a Greenwald density of . Sound speed and poloidal flow were estimated from the data to be and respectively. Safety factor close to the edge was found to be giving a connection length of about . This connection length, along with an electron thermal velocity of results in a transit frequency of . The collision times were calculated at and . The ions and electrons are weakly collisional thus, the adiabatic invariant is conserved over many bounce times.
Given that in the next section, the trapped particle orbits are going to play a fundamental role in our interpretation of the flux patterns, here, we give a few relevant numerical estimates regarding those. The fraction of trapped particles at the edge was estimated from the relation with being the aspect ratio. The bounce times of deeply trapped particles were calculated to be and . Technically, the bounce time of a deeply trapped particle whose banana orbit turning point is located at the high field side (HFS) midplane is infinite. The above bounce times are simple estimates found by dividing the length of a poloidal circuit with the magnetic drift speed of each particle and they are presented here to make the point that almost all of the ions have enough time, within the simulation, to complete a banana orbit. With , we calculate and , we can deduce that at the edge, the electrons are in the plateau regime () whereas the ions are in the banana transport regime ().
Iii The Fluxes
iii.1 Flux definitions
One of the goals of the present paper is to explore the relationship between equilibrium and turbulent processes and the associated fluxes across the separatrix. For dynamical quantities such as the plasma density , or velocity we define:
where denotes an average in either time or toroidal planes , or both. The time averages considered here are taken over a time interval late in the simulation where a quasi-steady turbulent state has been achieved.
Employing Eq. (1), the fundamental relationship for the fluxes is:
where the cross terms and vanish due to the vanishing of and .
These “fluxes” are not technically transport fluxes but rather local density weighted flows. Although collisions are fully accounted for in the simulation, the rapid electron transit time implies that the radial excursions of electrons due to either magnetic or drifts are small and can cancel without causing net transport as is well known from neoclassical theory. In this paper, we will continue to refer to them as fluxes for brevity. As will be seen, they provide a useful diagnostic for understanding basic properties of both the background plasma and the resulting turbulence. Often in the literature, we find the first piece of the rhs of Eq. (2) to be called the equilibrium flux whereas the second piece is known as the turbulent flux. In this work, we employ the same definition for the turbulent flux however, for convenience in analysis the definition of the equilibrium flux employed in this paper is slightly modified. Our equilibrium flux is given by . These two definitions of the equilibrium flux are not exactly equivalent. They are connected by the relation: . However, independent calculation of the difference , has shown that this “ringing” term, associated with temporal fluctuations of the toroidally averaged fields, is three orders of magnitude smaller than the two equilibrium fluxes indicating that, for the level of numerical accuracy assumed in this work, we can take the two definitions to be equivalent.
The definitions for the heat fluxes are similar, with the replacement of the density by the pressure of each species. Because the simulation is quasineutral, we can not distinguish between the densities or the net particle fluxes of the two species. In general both magnetic and ExB drifts contribute to the local particle fluxes crossing the separatrix at a particular poloidal position. In principle, the electron and ion ExB drift fluxes can be different due to the larger ion gyroradius. This orbit averaging effect results in the suppression of the ion fluxes by a factor . Close to the separatrix, and from power spectra of equilibrium and turbulent fields, , which gives a suppression of about for the equilibrium and about for the turbulent fluxes. In the case of the heat fluxes though, a clear distinction between ions and electrons can be made, based on their different temperatures. Finally, we should note that, given the available diagnostics, we lack the capability of calculating the turbulent polarization current and its flux. Nevertheless, in later sections, we estimate its magnitude relative to the rest of the fluxes.
iii.2 Flux calculation
In this section, we describe our method for calculating the various fluxes, using the XGC1 data, which illustrates the limitations of the particular data set and can be useful to practitioners in the field.
The calculation of radial drift velocities from the density and potential fields, whose values we know on the computational nodes, involves derivatives along the and directions. However, a plot of time and toroidally averaged and on the separatrix as the one in Fig. 3 shows that, despite the fact that a characteristic trend can be immediately identified, the values are noisy, especially at the low field side (LFS). Therefore, the strategy of interpolating the computational node values along horizontal and vertical lines and then taking derivatives along those lines to construct the vector velocity can give highly inaccurate results particularly for a sub-dominant component, such as the radial component of the velocity. For that reason, we resorted to expressing the required derivatives as derivatives with respect to the poloidal angle . This can be accomplished by writing the magnetic field in terms of the coordinate system . In this coordinate system, the magnetic field is . Using the formula for the inverse Jacobian, we can turn the expression for the radial component of the velocity, , into
For all equilibrium quantities, since the system has axisymmetry. For turbulent quantities, we can compute such derivatives using the assumption that turbulent structures are field aligned. Demanding we arrive at the equation
which can be readily used. The above expressions allow us to work with derivatives in the poloidal direction for the calculation of the velocity. This procedure has the advantage that we only need to do a single interpolation and smoothing of and around the LCFS and then take the derivative. The results thus obtained are much less noisy.
iii.3 The Particle Fluxes
First, we will describe the fluxes qualitatively and attempt to give an intuitive explanation for their shapes. Afterwards, we will seek to arrive at relevant quantitative conclusions based on their relative sizes.
iii.3.1 Equilibrium Flux
We start with the spatial distribution of the equilibrium, particle flux normal to, and hence crossing, the separatrix, which is common between the two species, as discussed above.
In Fig. 2 we observe that the magnitude of the particle flux (arrow length and color indicate flux intensity) is not uniform and exhibits a strong poloidal variation, not all of which can be explained by the difference in field strength between the inboard and outboard sides. Certainly, the flux is suppressed at the HFS compared to the LFS. However, at the LFS, the flux seems to be concentrated at two lobes, one outward, starting near the X-point and extending up to a little below the midplane, and one inward, which starts above the midplane and extends near the top of the machine. At the top, near the top X-point, the flux becomes again slightly outward, before it turns again to slightly inward at the HFS all the way to the bottom X-point.
To understand the poloidal pattern we just described, we turn to the distributions of density and potential on the LCFS which we show in Fig. 3. There, we plot the two terms of the adiabatic relation so that we get a measure of the electron departure from adiabaticity at various locations along the LCFS. Moreover, the term has the same structure as the potential itself, given that the electrons are almost isothermal on the flux surfaces. From this plot, we immediately recognize the points where the global minimum () and global maximum () are located, which are the top and bottom X-points respectively. Although the exact location of these two extreme points is a result of a complex interplay between Pfirsch-Schlütter flows and collisional effects aydemir2012pfirsch (); aydemir2009intrinsic (), it follows from the direction of the magnetic drifts: ions drift downwards and electrons upwards therefore, the potential arranges itself so that it attracts ions from the core at the top and electrons at the bottom. In both cases, the extremization of the potential, follows an analogous increase/decrease of density. In the case of the bottom X-point, an important factor that would lead to a density build-up, is the X-point ion loss. The X-point losschang2002x (); ku2004property () of hot ions leads to a lower ion temperature at the LCFS. It should be noted that, in our simulation, the well from the X-loss has already been established. Therefore, the loss cone has moved to higher energies, resulting in hot ion loss. Although the scalar ion pressure is not constant along a field linechurchill2017total (), approximate pressure balance is consistent with the lower temperature to be offset by a density increase. These two extremal points of potential establish an electric field in the HFS region. It is this electric field which in turn creates the flow which is responsible for the inward flux pinch we see in Fig. 2 on the HFS as well as the dominant outward flux on the LFS. The reason why this flux is not as high as the ones we see on the LFS, even if the electric field is stronger, is the mitigating effect of the strong magnetic field.
In addition to the global minimum and the global maximum of , we also have two local ones. The positions of those, correlate with the two bulges of opposite direction that we see on the LFS in Fig. 2. Indeed, the directions of the fluxes are completely in line with the electric field directions established by the alternating pairs of minimums and maximums of . However, in contrast to the two global extremal points of potential, it is not immediately clear what is responsible for the two local ones in the LFS. Here, we propose that this potential structure is produced by the orbits of trapped ions which, as we saw in the previous section, constitute roughly of the total ion population near the edge. A trapped ion moving on its banana orbit exits the separatrix at some point below the LFS midplane, leaving behind it a negative charge. If the electrons were completely adiabatic, negative charge would be instantly neutralized by electrons moving along the field lines to prevent the potential build-up. But, as we can readily verify from Fig. 3, there is a significant departure from adiabaticity at the LFS. Therefore, a trapped ion leaving the separatrix will create a negative potential in order to attract more ions. Conversely, the region above the LFS midplane is the region where the trapped ions re-enter the LCFS (at least those that didn’t get entrained in the parallel SOL flow) generating a positive charge. The non-adiabaticity of the electrons, forces the potential to increase, in order to attract them at that location. Thus neoclassical ion orbits combined with X-point losses are the dominant mechanism responsible for the structure of the poloidal electric field, and hence the equilibrium radial particle flux in this simulation.
The non-adiabaticity of the electrons at the separatrix, seems to be crucial for the exact shape of the edge potential. This non-adiabaticity at the edge can be understood if we compare the time rates of all terms in the electron momentum equation. More specifically, the inertial term, i.e., the transit frequency, was given in Section II as . Although large, it is comparable to the electron collision frequency (also given in II) and the turbulent frequency, which is discussed in the next section and it is found to be . Therefore, both collisions and the inertial term are important in the electron momentum equation leading to the observed departure from adiabaticity.
iii.3.2 Banana Tip distribution
The subject of neoclassical direct ion orbit losses and their contribution to the edge flows and electric fields has been treated substantially, both computationally and theoreticallydegrassie2015thermal (); battaglia2014kinetic (); miyamoto1996direct (); stacey2015distribution (); stacey2011effect (); stacey2016recent (); stacey2013interpretation (); wilks2017calculation (); wilks2016improvements (); hahn2005wall (). Here, we offer a new argument to support our claim that the potential structure at the LFS is due to the excursions of trapped ions on their banana orbits. Assuming that the trapped particles, with orbits within a banana width from the edge, will create either a charge excess or a charge hole at the point of inflection (“banana tip”), we give an estimate for the angular distribution of those banana tips at the LFS.
To get the fraction of trapped particles with banana tips in the range , where is the angle that the magnetic field is and zero is the angle at (LFS midplane), we perform the following integral:
where we have changed coordinates from velocity space to the , space, with being the energy, the adiabatic invariant. Here, denotes the sign of the particles parallel velocity.
Plugging in the angular formula for , , valid for circular flux surfaces, we turn this into the cumulative angular distribution for banana tips:
Interpreting as we get Therefore,
The poloidal distribution of the flux caused by these particle excursions is the product of the tip distribution with the actual magnetic flux. In Fig. 4 we can see this flux compared to the fluxes. We point out that indeed, the “trapped” flux, tracks qualitatively the equilibrium flux. Two things should be noted: the modelling of the tip distribution is rather simple as it assumes circular flux surfaces and completely ignores the effect of the X-point and the ion losses there. Also, a qualitative similarity between the fluxes is all we can hope for. Indeed, if the equilibrium flux pattern at the LFS is the result of the potential structure due to the trapped particle excursions, as we claim, then, the relationship between their relative sizes is far from obvious.
iii.3.3 Turbulent Flux
We proceed with the description of the turbulent flux which we show in Figs. 4 and 5. There, we observe that particle turbulence is entirely confined to the LFS and has a ballooning-like shape which peaks near the midplane, gets interrupted above it up to roughly , picks up again after that and recedes off at the top of the machine.
A probable explanation for such a shape is provided by Fig. 6. There, we overlay the plots of the turbulent flux and the shearing rate on the LCFS. We observe that as the shearing rate exceeds a certain value (roughly at ) the turbulent flux drops sharply. When the drops below this value, the turbulence raises again. Indeed, we have measured the main frequency of the turbulence in this discharge to be very close to . We can imagine therefore, that in the absence of shear, the turbulence would have created a ballooning structure in the LFS, peaking at the midplane and tapering off at the top, where the curvature drive is minimal. However, the presence of shear at exactly the right level, suppresses it locally creating the resulting patternburrell1997effects (); hahm1995flow (); terry2000suppression (); biglari1990influence ().
Recall from the previous discussion that the electrostatic potential pattern, and hence the shear, appears to be controlled by neoclassical drift-orbit effects. It is notable that the system has arranged itself so that this equilibrium shear is of sufficient magnitude to influence the turbulence. We have verified that turbulence generated Reynolds stresses play a negligible role in establishing these sheared flows, which is not unexpected for this H-mode phase of the dischargemuller2011experimental ().
In addition to shear, it is possible that the observed turbulence pattern may also be related to the existence of curvature driven modes which peak at the top of torusabdoul2017generalised (); halpern2013ideal (). Poloidal flux surface expansion, i.e., the distance between adjacent flux surfaces which is increasing towards the top of the machine, may also play a role. First, it increases the outward velocity and hence the of filamentary turbulent structures in order to maintain their field-aligned charactergalassi2017drive (). Second, flux surface expansion reduces shear more quickly than the density or temperature gradients, since the former are proportional to a second radial derivative of the potential and therefore scale more strongly than the gradients which are proportional to a single radial derivative of the profiles.
Lastly, given the turbulent frequency, we can evaluate the relative size of the ion polarization flux compared to the turbulent flux. Their ratio is approximately , which means that our inability to calculate the polarization current does not threaten the credibility of our results, at the present level of precision.
iii.4 The Heat Fluxes
In this section we present the qualitative features of the heat fluxes. Contrary to the particle fluxes where we couldn’t distinguish between the two species (apart from the, almost negligible, gyroaveraging effect of the ion orbits), for the heat fluxes we have a clear separation due to the different temperatures, as shown on Fig. 7. Therefore, in what follows, the calculation of the fluxes still assumes quasineutrality but each of the species has its own temperature. In calculating the heat fluxes, we employ Eq. (2) where in the place of velocity we substitute the pressure , with denoting the species. Therefore, the heat flux is given by the formula
where again, we call the first term of the rhs the “equilibrium” and the second the “turbulent” heat flux. The velocity here is the one.
In Figs. 8-9 we present the equilibrium and turbulent fluxes of electrons and ions respectively, along the LCFS. The shapes of both species equilibrium flux is similar except for the scaling effect of the temperature difference. This can be easily understood from Fig. 7 where we plot the temperatures as a function of the poloidal angle. There, we see that the electrons are very close to being isothermal. Therefore, the qualitative features between the particle equilibrium and electron equilibrium heat fluxes are common.
The ion temperature on the other hand, does not remain constant along the separatrix. Specifically, it drops near the two X-points indicating that, close to these areas, hotter ions have exited the plasma. High energy ions with the right pitch angle will be in the ion loss conechang2002x (); ku2004property (). An interesting feature of the ion temperature is that the average temperatures below and above the outboard midplane are not the same: the region below the midplane has an average temperature which is higher than the average temperature above. We can understand this difference in the context of the trapped ion banana orbits that we talked about in the particle fluxes section. The ions exiting the flux surface below the midplane are on the inward side of their banana orbit. They travel along a region of higher temperature therefore the temperature locally increases. Equivalently, when they re-enter above the midplane, they are coming from outside the LCFS where the temperature has dropped. Indicatively, a typical value for a banana orbit width near the edge is and the temperature drop across this distance is about , accounting in part for the observed average temperature difference. Moreover, the re-entering ions could have originated from further away in the SOL than a typical banana orbit width, making the above estimate just a lower bound for the real temperature difference. Despite the above observations, the qualitative features of the equilibrium ion heat flux are very similar to those of the equilibrium particle and electron heat fluxes.
The real difference between the two species heat flux behavior comes from the comparison between the relative sizes of equilibrium and turbulent heat fluxes. In the case of the electrons, turbulence dominates; neoclassical electron transport is negligible because of the small electron banana width. Thus, the local equilibrium electron flux shown in Fig. 8 will be further reduced by cancellation from bounce-orbit averaging. Turbulence is the sole mechanism for net electron particle loss across the separatrix. It must compensate for the net loss of ions in order to maintain charge balance. Turbulent electron energy loss follows.
In contrast, for the ions, it seems that the equilibrium process, set up by the neoclassical mechanism of trapped particle orbits, is indeed the dominant mode of ion loss, in this simulation, both for heat and particles. Both species LFS turbulent ballooning shapes are being interrupted in the high shear region. Although the ion turbulent heat flux never recovers and drops to zero, the electron turbulent heat flux picks up again after the initial dip. Moreover, in the region above the midplane, the electron turbulent heat flux is almost twice the size of the equivalent ion one, despite the big difference of the absolute temperature of the two species, with the ions being significantly hotter. This is a clear indication that, in the case of electrons, turbulence is absolutely necessary for the plasma to retain charge balance, whereas in the case of ions, the neoclassical processes are more effective.
In Figs. 10-11 we present the effective diffusivities which are defined as , and we separate them again into equilibrium and turbulent contributions. The equilibrium diffusivities are presumably driven by neoclassical ion physics as discussed previously, but because of ambipolarity, equilibrium and turbulent diffusivities are fundamentally coupled. A first observation is that both equilibrium and turbulent diffusivities are of a similar order of magnitude, consistent with this coupling.
In Ref. 2017APS..DPPTP11079, , methods developed for application to pedestal turbulence provide information about the modes present, based on the relative sizes of the transport channels. It is interesting to apply those concepts to our situation for qualitative insight. Here we find comparable turbulent diffusivities near the midplane in particle, electron heat and ion heat channels. This is suggestive of an underlying mode with an MHD character, possibly resistive ballooning or Kelvin Helmholtz in this electrostatic simulation. Near the top of the torus the electron diffusivity is dominant, which may suggest a mode driven primarily by the electron temperature gradient, or a mode with an electron drift character. Indeed we find our turbulence has a frequency very close to the electron diamagnetic drift frequency.
A final observation is that, based on the scale of the above diffusivities, the time scale for further profile evolution is much larger than the total time of the simulation. We have indeed verified that the profiles remain constant throughout, which is a clear indication of achievement of a quasi-steady turbulent state.
Iv Quantitative Results
iv.1 Separatrix-crossing Fluxes
Here, we calculate the integrated currents and flux surface averages of the above fluxes and reach some conclusions regarding their relative importance in establishing the SOL width.
First, we start with some definitions. The fluxes of the previous section are normal to the separatrix components of the total flux, i.e., . For a quantity , the flux surface average of it is defined as:
For the calculation of such integrals, we must know the Jacobian factor the inverse of which is given by the well known formula:
with being the distance between and the point on the flux surface. In the second step of the above calculation we have replaced the usual formula , which holds in a circular geometry, with , which is an approximation for a general shape flux surface, provided the integrand is non-singular and single valued, that minimizes the integration error. Using this formula, the flux surface average (FSA) of a scalar is given by
If we perform a flux-surface average on the continuity equation, we arrive at the expression:
where and . Therefore, the physically relevant quantity is , which is calculated using the explicit formula:
We can also define the surface integral of the flux across the LCFS. This should be interpreted as the particle or heat current, i.e., the number of particles per second that cross the separatrix due to the various fluxes. This will be denoted as and calculated using
with similar expressions for the heat currents.
In Table 1, we present the integrals over the whole separatrix. We have also kept for reference the numbers for flux from the magnetic drifts. In Tables 2-3 we give the results for the integrated heat fluxes.
Here, we also define an averaging factor, , which is defined as and shows the degree of “cancellation” of each particular flux. If we assume that the density and potential are constant on a flux surface, then the FSA’s of both the equilibrium and the magnetic fluxes can be shown to be identically zero. Close to the edge of the plasma, as we saw in the previous section, the constancy of those quantities on the separatrix is violated hence, we don’t expect the to be exactly zero. Indeed, as we see in Table 4, the averaging factor for both the magnetic drift flux and the equilibrium flux are very small but nonzero, indicating that the bulk of these fluxes are returning into the main plasma whereas, the of the turbulent flux is larger. ’s of the heat fluxes are similar.
|Ion Heat Fluxes|
|Electron Heat Fluxes|
A few comments about the relative sizes of those integrals are in order: In the particle fluxes and ion heat flux case, the size of the turbulent current is roughly three times larger than the equilibrium one. In the case of the electron heat fluxes though, the turbulent current is almost ten times larger than the equilibrium one indicating the importance of turbulence for the electrons, supporting the main argument of the previous section.
iv.2 Particle Balance in the SOL
In this sub-section, we explore particle balance in the SOL by comparing the flux of electrons exiting across the separatrix with the exhaust flux in the SOL, which is dominated by the parallel flows. The competition between these cross-field and parallel transport processes is what sets the density SOL width.
Regarding the separate contribution of each flux to the setting of the SOL width the ratio of the currents doesn’t tell the full story. Despite the fact that the turbulent flux is fully localized at the outboard midplane, contributing solely to the LFS SOL, the equilibrium fluxes have an important inward component at the HFS. Moreover, as we have mentioned in the previous sections, the particle equilibrium flux cannot take electrons outside of the separatrix. Therefore, we start with the assumption that the current from the turbulent particle flux is the total electron current exiting the separatrix. Now, we take a Gaussian box whose one side is the part of LCFS starting at the point where the turbulent particle flux begins and ending where it ends. The top and bottom of this box are horizontal lines, starting at these two locations and extending outwards, deep inside the SOL. The extent of these lines is left to to be determined by the behavior of the normal fluxes on them. The right hand side of the box is unspecified but we assume that no flux leaks from this side of the box to the wall.
In Figs. 12-13 we have plotted the normal component of the particle fluxes at the bottom and top caps of the Gaussian box respectively. The normal component of the particle flux includes contributions from the flux, the magentic drift flux and, importantly, the parallel flows. We are not concerned with distinguishing the contributions from equilibrium and turbulent therefore, we just present the total. An important note here is that, since we do not distinguish between equilibrium and turbulence, in the calculation of the velocity we take . We believe that the error introduced by such an approximation is minimal since, the contribution to the normal fluxes is only important in a very narrow region close to the separatrix. There, the poloidal magnetic field is almost zero hence, . Generally, the parallel flux contribution is dominant everywhere and the only place where the other fluxes seem to be relevant is at a very narrow channel close to the LCFS. In both figures, the negative direction is downwards, which means that for the top cap, negative fluxes are entering the box whereas for the bottom cap negative fluxes are leaving.
In the case of the top cap, it seems that there is a well defined point where all fluxes drop to zero, around . This seems to be a reasonable place to truncate the integration of the fluxes in order to find the total current crossing the top of the box. The total current thus found is , entering the box.
For the bottom cap, things are not so clear cut. Although the fluxes coming from and magnetic drifts do indeed drop to zero, the parallel flux seems to remain constant. This had to be expected since there is a significant degree of ionization happening in the SOL. It is very hard to precisely locate a point to truncate the integration, meaning to find where the SOL ends and what is left after that is simply ionization. We believe that a sensible point to stop our integration is at (shown in figure 12). The total current found from this choice is , exiting the box. If we extend the integration further out, this number will increase.
As we said above, the total electron current crossing into the box from the separatrix, has to be the one coming from the turbulent flux given in Table 1 as . If we compare this number to the total current found from the top and bottom of the boxes, we see that it accounts for . A direct conclusion that can be drawn from such a comparison is that indeed, as presumed in the Goldston model, the density channel is filled up almost entirely with turbulent electron flux. Of course, had we taken the limit of integration of the bottom cap further out, this percentage would have dropped however, it is reasonable to assume that this result would be severely contaminated by atomic processes that take place outside the LCFS and have nothing to do with turbulence and flows that cross the separatrix.
Here, we need to note that even though we can safely claim that the electron contribution to the SOL channel is almost entirely turbulent, a similar statement can not be made for the ions. From the data that we have available, we have no way of isolating the turbulent and the neoclassical contributions since turbulence, in the ions’ case, exists on top of a neoclassical background.
If we repeat the previous calculation for electron heat fluxes, we find both the top and bottom caps of the box to have well defined limits where all fluxes go to zero. However, the comparison of the total electron heat current exiting the Gaussian box from the caps to the electron heat current entering from the separatrix, shows that the latter is more than twice the former. A plausible explanation for this result is that heat is dissipated in this volume due to various atomic processes, predominantly ionization. A detailed account of the relative strengths of such processes is beyond the scope of this paper.
V Summary and Conclusions
In this paper we have analyzed the results of an XGC1 DIII-D simulation focusing on the fluxes exiting the separatrix. We reconstructed the flux distinguishing between the equilibrium and the turbulent part. The poloidal patterns of the two particle fluxes were presented. For the equilibrium flux, we traced it’s poloidal shape to the LCFS potential. We interpreted this potential structure as the result of ion drifts. Specifically, we could clearly distinguish the impact of the -drifts and X-point losses on the potential’s global maximum and minimum at the bottom X-point and top of the machine, respectively. For the poloidal potential variation at the LFS, we argued that it’s due to the trapped ion orbits close to the edge that exit and re-enter the confinement. We constructed a simplified model for the angular distribution function of the banana orbits’ inflection points and showed that the resulting “trapped-particle” flux follows the shape of the equilibrium flux, to support our argument. Moreover, we noted that LFS edge electrons exhibit a significant departure from adiabaticity. We provided as reasons for this departure the fact that the turbulent frequency and the collision rate are comparable to the transit frequency. The turbulent flux pattern was also given. It appears to be strongly ballooning but with a significant region of suppression. We contended that this suppression is coming from the ion-neoclassically generated shear which seems to be at the right level to affect the turbulence exactly at that region.
We presented the heat fluxes for ions and electrons at the separatrix. From the difference in size between the turbulent heat fluxes of the two species it is clear that turbulence is more important for the electrons since it is their only way of losing confinement. Ions are exiting the separatrix through neoclassical mechanisms, therefore, turbulent electron heat flux has to increase in order to maintain quasineutrality. The electrons were found to be isothermal whereas the ion temperature exhibits significant drops near the X-points, in accordance with the predictions of X-loss theory. The small temperature difference between the upper and lower halfs of the LFS can be, in part, explained by ion banana orbits exiting from a hotter and re-entering from a cooler part of the machine. The equilibrium and turbulent diffusivities of those fluxes were calculated, providing insight regarding the modes present. Numerical integrations were performed to calculate the flux-surface-average fluxes and the outgoing particle and heat currents from the separatrix. Lastly, we confirmed that anomalous electron diffusion accounts for almost all of the electrons filling the SOL particle channel.
Given the importance of the predictability of the SOL width for modern tokamaks, more state-of-the-art gyrokinetic simulations need to be analyzed to elucidate the complex physical processes that give rise to it. The interaction between the ion drifts and the turbulence seems to be critical for setting the size of the channel and the study of the fluxes and flows at the edge provides illuminating insights into the nature of the system dynamics. In future publications we will examine simulations from different tokamaks and look into the characteristics of the turbulent “blobs” from this and other simulations. These efforts will be directed towards obtaining a better understanding of the interaction of neoclassical and turbulent physics, and the extent to which the present results are applicable to other plasma and device regimes.
I. Keramidas Charidakos would like to acknowledge useful help with the data from Jugal Chowdhury at the early phase of this project. We acknowledge computing resources on Titan at OLCF through the 2015 INCITE and the 2016 ALCC awards. Work supported by the U.S. Department of Energy Office of Science, Office of Fusion Energy Sciences under Award Number DE-FG02-97ER54392 and by subcontract SO15882-C with PPPL under the U.S. Department of Energy HBPS SciDAC project.
-  JW Connor, GF Counsell, SK Erents, SJ Fielding, B LaBombard, and K Morel. Comparison of theoretical models for scrape-off layer widths with data from compass-d, jet and alcator c-mod. Nuclear fusion, 39(2):169, 1999.
-  RJ Goldston. Heuristic drift-based model of the power scrape-off width in low-gas-puff h-mode tokamaks. Nuclear Fusion, 52(1):013009, 2011.
-  Choong Seock Chang, S Ku, Alberto Loarte, Vassili Parail, Florian Koechl, Michele Romanelli, Rajesh Maingi, J-W Ahn, T Gray, J Hughes, et al. Gyrokinetic projection of the divertor heat-flux width from present tokamaks to iter. Nuclear Fusion, 57(11):116023, 2017.
-  ET Meier, RJ Goldston, EG Kaveeva, MA Makowski, S Mordijck, VA Rozhansky, I Yu Senichenkov, and SP Voskoboynikov. Drifts, currents, and power scrape-off width in solps-iter modeling of diii-d. Nuclear Materials and Energy, 12:973–977, 2017.
-  ET Meier, RJ Goldston, EG Kaveeva, MA Makowski, S Mordijck, VA Rozhansky, I Yu Senichenkov, and SP Voskoboynikov. Analysis of drift effects on the tokamak power scrape-off width using solps-iter. Plasma Physics and Controlled Fusion, 58(12):125012, 2016.
-  JR Myra, DA Russell, DA DâIppolito, J-W Ahn, Rajesh Maingi, RJ Maqueda, DP Lundberg, DP Stotler, SJ Zweben, J Boedo, et al. Reduced model simulations of the scrape-off-layer heat-flux width and comparison with experiment. Physics of Plasmas, 18(1):012305, 2011.
-  David A Russell, Daniel A D’Ippolito, James R Myra, John M Canik, Travis K Gray, and Stewart J Zweben. Modeling the effect of lithium-induced pedestal profiles on scrape-off-layer turbulence and the heat flux width. Physics of Plasmas, 22(9):092311, 2015.
-  AY Pankin, T Rafiq, AH Kritz, GY Park, CS Chang, D Brunner, RJ Groebner, JW Hughes, B LaBombard, JL Terry, et al. Kinetic modeling of divertor heat load fluxes in the alcator c-mod and diii-d tokamaks. Physics of Plasmas, 22(9):092511, 2015.
-  V Rozhansky, E Kaveeva, I Senichenkov, and E Vekshina. Structure of the classical scrape-off layer of a tokamak. Plasma Physics and Controlled Fusion, 60(3):035001, 2018.
-  D Reiser and T Eich. Drift-based scrape-off particle width in x-point geometry. Nuclear Fusion, 57(4):046011, 2017.
-  T Stoltzfus-Dueck. Tokamak-edge toroidal rotation due to inhomogeneous transport and geodesic curvature. Physics of Plasmas, 19(5):055908, 2012.
-  T Stoltzfus-Dueck. Transport-driven toroidal rotation in the tokamak edge. Physical review letters, 108(6):065002, 2012.
-  SH Müller, JA Boedo, KH Burrell, RA Moyer, DL Rudakov, WM Solomon, et al. Experimental investigation of the role of fluid turbulent stresses and edge plasma flows for intrinsic rotation generation in diii-d h-mode plasmas. Physical review letters, 106(11):115001, 2011.
-  J Loizu, P Ricci, FD Halpern, S Jolliet, and A Mosetto. Intrinsic toroidal rotation in the scrape-off layer of tokamaks. Physics of Plasmas, 21(6):062309, 2014.
-  RJ Groebner, KH Burrell, WM Solomon, et al. Intrinsic toroidal velocity near the edge of diii-d h-mode plasmas. Nuclear Fusion, 49(8):085020, 2009.
-  CS Chang and S Ku. Spontaneous rotation sources in a quiescent tokamak edge plasma. Physics of Plasmas, 15(6):062510, 2008.
-  Janghoon Seo, CS Chang, S Ku, JM Kwon, W Choe, and Stefan H Müller. Intrinsic momentum generation by a combined neoclassical and turbulence mechanism in diverted diii-d plasma edge. Physics of Plasmas, 21(9):092501, 2014.
-  AV Chankin, DP Coster, N Asakura, G Corrigan, SK Erents, W Fundamenski, HW Müller, RA Pitts, PC Stangeby, and M Wischmeier. A possible role of radial electric field in driving parallel ion flow in scrape-off layer of divertor tokamaks. Nuclear Fusion, 47(8):762, 2007.
-  K Hoshino, A Hatayama, N Asakura, H Kawashima, R Schneider, and D Coster. Numerical analysis of the sol/divertor plasma flow with the effect of drifts. Journal of nuclear materials, 363:539–543, 2007.
-  GS Kirnev, G Corrigan, D Coster, SK Erents, W Fundamenski, GF Matthews, and RA Pitts. Edge2d code simulations of sol flows and in–out divertor asymmetries in jet. Journal of nuclear materials, 337:271–275, 2005.
-  B LaBombard, JE Rice, AE Hubbard, JW Hughes, M Greenwald, J Irby, Y Lin, B Lipschultz, ES Marmar, CS Pitcher, et al. Transport-driven scrape-off-layer flows and the boundary conditions imposed at the magnetic separatrix in a tokamak plasma. Nuclear fusion, 44(10):1047, 2004.
-  A Yu Pigarov, SI Krasheninnikov, B LaBombard, and TD Rognlien. Simulation of parallel sol flows with uedge. Contributions to Plasma Physics, 48(1-3):82–88, 2008.
-  RA Pitts, P Andrew, X Bonnin, AV Chankin, Y Corre, G Corrigan, D Coster, I Duran, T Eich, SK Erents, et al. Edge and divertor physics with reversed toroidal field in jet. Journal of Nuclear Materials, 337:146–153, 2005.
-  AY Aydemir. Pfirsch–schlüter current-driven edge electric fields and their effect on the l–h transition power threshold. Nuclear Fusion, 52(6):063026, 2012.
-  Ker-Chung Shaing and EC Crume Jr. Bifurcation theory of poloidal rotation in tokamaks: A model for l-h transition. Physical Review Letters, 63(21):2369, 1989.
-  CS Chang, S Ku, GR Tynan, R Hager, RM Churchill, I Cziegler, M Greenwald, AE Hubbard, and JW Hughes. Fast low-to-high confinement mode bifurcation dynamics in a tokamak edge plasma gyrokinetic simulation. Physical review letters, 118(17):175001, 2017.
-  RJ Goldston. Theoretical aspects and practical implications of the heuristic drift sol model. Journal of Nuclear Materials, 463:397–400, 2015.
-  Thomas Eich, AW Leonard, RA Pitts, W Fundamenski, RJ Goldston, TK Gray, A Herrmann, A Kirk, A Kallenbach, O Kardaun, et al. Scaling of the tokamak near the scrape-off layer h-mode power width and implications for iter. Nuclear fusion, 53(9):093031, 2013.
-  CS Chang, S Ku, PH Diamond, Z Lin, S Parker, TS Hahm, and N Samatova. Compressed ion temperature gradient turbulence in diverted tokamak edge. Physics of Plasmas, 16(5):056108, 2009.
-  S Ku, CS Chang, and PH Diamond. Full-f gyrokinetic particle simulation of centrally heated global itg turbulence from magnetic axis to edge pedestal top in a realistic tokamak geometry. Nuclear Fusion, 49(11):115021, 2009.
-  S Ku, Robert Hager, Choong-Seock Chang, JM Kwon, and Scott E Parker. A new hybrid-lagrangian numerical scheme for gyrokinetic simulation of tokamak edge plasma. Journal of Computational Physics, 315:467–475, 2016.
-  RM Churchill, CS Chang, S Ku, and J Dominski. Pedestal and edge electrostatic turbulence characteristics from an xgc1 gyrokinetic simulation. Plasma Physics and Controlled Fusion, 59(10):105014, 2017.
-  DA DâIppolito, JR Myra, and SJ Zweben. Convective transport by intermittent blob-filaments: Comparison of theory and experiment. Physics of Plasmas, 18(6):060501, 2011.
-  JR Myra, DA D’Ippolito, and DA Russell. Turbulent transport regimes and the scrape-off layer heat flux width. Physics of Plasmas, 22(4):042516, 2015.
-  JR Myra, DA Russell, and SJ Zweben. Theory based scaling of edge turbulence and implications for the scrape-off layer width. Physics of Plasmas, 23(11):112502, 2016.
-  I Furno, B Labit, M Podestà, A Fasoli, SH Müller, FM Poli, P Ricci, C Theiler, S Brunner, A Diallo, et al. Experimental observation of the blob-generation mechanism from interchange waves in a plasma. Physical review letters, 100(5):055004, 2008.
-  N Bisai, A Das, S Deshpande, R Jha, P Kaw, A Sen, and R Singh. Formation of a density blob and its dynamics in the edge and the scrape-off layer of a tokamak plasma. Physics of plasmas, 12(10):102515, 2005.
-  CS Chang, Seunghoe Kue, and H Weitzner. X-transport: A baseline nonambipolar transport in a diverted tokamak plasma edge. Physics of Plasmas, 9(9):3884–3892, 2002.
-  Seunghoe Ku, Hoyul Baek, and CS Chang. Property of an x-point generated velocity-space hole in a diverted tokamak plasma edge. Physics of plasmas, 11(12):5626–5633, 2004.
-  James L Luxon. A design retrospective of the diii-d tokamak. Nuclear Fusion, 42(5):614, 2002.
-  AY Aydemir. An intrinsic source of radial electric field and edge flows in tokamaks. Nuclear Fusion, 49(6):065001, 2009.
-  Randy M Churchill, John M Canik, CS Chang, R Hager, Anthony W Leonard, Rajesh Maingi, R Nazikian, and Daren P Stotler. Total fluid pressure imbalance in the scrape-off layer of tokamak plasmas. Nuclear Fusion, 57(4):046029, 2017.
-  John S deGrassie, Jose A Boedo, and Brian A Grierson. Thermal ion orbit loss and radial electric field in diii-d. Physics of Plasmas, 22(8):080701, 2015.
-  DJ Battaglia, KH Burrell, CS Chang, S Ku, JS Degrassie, and BA Grierson. Kinetic neoclassical transport in the h-mode pedestal. Physics of Plasmas, 21(7):072508, 2014.
-  Kenro Miyamoto. Direct ion orbit loss near the plasma edge of a divertor tokamak in the presence of a radial electric field. Nuclear fusion, 36(7):927, 1996.
-  Weston M Stacey and Matthew T Schumann. The distribution of ion orbit loss fluxes of ions and energy from the plasma edge across the last closed flux surface into the scrape-off layer. Physics of Plasmas, 22(4):042504, 2015.
-  Weston M Stacey. The effect of ion orbit loss and x-loss on the interpretation of ion energy and particle transport in the diii-d edge plasma. Physics of Plasmas, 18(10):102504, 2011.
-  WM Stacey. Recent developments in plasma edge theory. Contributions to Plasma Physics, 56(6-8):495–503, 2016.
-  WM Stacey, M-H Sayer, J-P Floyd, and RJ Groebner. Interpretation of changes in diffusive and non-diffusive transport in the edge plasma during pedestal buildup following a low-high transition in diii-d. Physics of Plasmas, 20(1):012509, 2013.
-  TM Wilks, WM Stacey, and TE Evans. Calculation of the radial electric field from a modified ohm’s law. Physics of Plasmas, 24(1):012505, 2017.
-  TM Wilks and WM Stacey. Improvements to an ion orbit loss calculation in the tokamak edge. Physics of Plasmas, 23(12):122505, 2016.
-  SH Hahn, Seunghoe Ku, and CS Chang. Wall intersection of ion orbits induced by fast transport of pedestal plasma over an electrostatic potential hill in a tokamak plasma edge. Physics of plasmas, 12(10):102501, 2005.
-  KH Burrell. Effects of e b velocity shear and magnetic shear on turbulence and transport in magnetic confinement devices. Physics of Plasmas, 4(5):1499–1518, 1997.
-  TS Hahm and KH Burrell. Flow shear induced fluctuation suppression in finite aspect ratio shaped tokamak plasma. Physics of Plasmas, 2(5):1648–1651, 1995.
-  PW Terry. Suppression of turbulence and transport by sheared flow. Reviews of Modern Physics, 72(1):109, 2000.
-  H Biglari, PH Diamond, and PW Terry. Influence of sheared poloidal rotation on edge turbulence. Physics of Fluids B: Plasma Physics, 2(1):1–4, 1990.
-  PA Abdoul, David Dickinson, CM Roach, and Howard Read Wilson. Generalised ballooning theory of two-dimensional tokamak modes. Plasma Physics and Controlled Fusion, 60(2):025011, 2017.
-  Federico D Halpern, Sebastien Jolliet, Joaquim Loizu, Annamaria Mosetto, and Paolo Ricci. Ideal ballooning modes in the tokamak scrape-off layer. Physics of Plasmas, 20(5):052306, 2013.
-  Davide Galassi, Patrick Tamain, Hugo Bufferand, Guido Ciraolo, Ph Ghendrih, C Baudoin, Clothilde Colin, Nicolas Fedorczak, N Nace, and Eric Serre. Drive of parallel flows by turbulence and large-scale e b transverse transport in divertor geometry. Nuclear Fusion, 57(3):036029, 2017.
-  M. Kotschenreuther, X. Liu, D. Hatch, L. Zheng, S. Mahajan, A. Diallo, R. Groebner, A. Hubbard, J. Hughes, C. Maggi, S. Saarelma, and JET Contributors. Gyrokinetic analysis of pedestal transport. In APS Meeting Abstracts, page TP11.079, October 2017.