135
Journal of Nuclear Materials 128 & 129 (1984) 135-140
AN AXIALLY AVERAGED-RADIAL TRANSPORT MODEL OF TOKAMAK EDGE PLASMAS A.K. PRINJA
and R.W. CONN
School of Engineering and Applied Science, University of Ca&knia,
Los Angeles, Los Angeles, CA 90024, USA
Key words: radial transport, scrape-off layer, pump limiter, divertor, recycling A two-zone axially averaged-radial transport model for edge plasmas is described that incorporates parallel electron and ion conduction, Iocalixed recycling, parallel efectron pressure gradient effects and sheath losses. Results for high recycling show that the radial electron temperature profile is determined by parallel electron conduction over short radial distances ( - 3 cm). At larger radius where T, has fallen appreciably, convective transport becomes equally important. The downstream density and ion temperature profiles are very flat over the region where electron conduction dominates. This is seen to result from a sharply decaying velocity profile that follows the radial electron temperature. A one-dimensional analytical recycling model shows that at high neutral pumping rates, the plasma density at the plate, nia, scales linearly with the unperturbed background density, ni,. When ionization dominates ni,/nio - exp(n& while in the intermediate regime ni,/ni, - exp( - nio). Such behavior is qualitatively in accord with experimental observations.
1. Introduction
Considerable experimental [l] and theoretical [2] effort is currently being focused on underst~~g the host of fundamental processes pertinent and unique to the boundary region in both diverted and limited plasmas. Lately, sophisiticated two-dimensional scrape-off layer models [3,4] of transport on open flux surfaces have begun to supplant the earlier one-dimensional radial or axial transport ]5,6] models, but have met with limited success for routine application. In conventional edge models [S,ll], emphasis has been on radial transport in order to provide a self-consistent link between central and boundary plasmas. Parallel losses to a target were uniformly described by the sheath model. However, this assumption breaks down in the high density, localized recycling regimes [7] which are likely to be typical of future machines. In such situations, one-dimensional models simulating flow along field lines [6,8] have demonstrated that gradients are indeed gentle over the bulk of the system but adjacent to the target (divertor or pumped limiter neutralizer plates) in regions approximately 3-5 ionization mean free paths long, recycling induces large variations in local density and temperature. Under such circumstances, finite conduction along field lines, thermal forces [9] and potentials [lo] begin to assume a significance that cannot be neglected for a detailed understanding of the edge plasma. In order to include such phenomena in a radial transport model without resorting to a fully two-dimensional model, we have developed an axially-averaged model that improves and extends conventional scrapeoff layer models. The axial region is divided into an upstre~ gradient free wne and a downstream recycling dominated zone where the sheath mode1 is applied. The two-dimensional transport equations [9] are then averaged in each zone, radial transport is neglected 0022-3115/84/$03.00 QJElsevier Science Publishers B.V. bosh-Holland Physics Pubfishing Division)
in the downstream region compared to parallel transport there and the flow velocity is prescribed as a fraction of the local sound speed [8]. The model is described in detail in section 3 and numerical results are presented and discussed in section 4. As high density low temperature edge plasmas [7,13] have become an attractive mode of operation in tokamaks, codes have focused on the physics of recycling that drive such modes [3,6-81. In order to illustrate the essential features of recycling plasmas, we present first in section 2 an analytical nonlinear model for density amplification as the background plasmas undergo a transition from transparency to opaqueness to neutrals.
2. Analytical reeyefing model A simple analytical model illustrates the physics of localized recycling in the presence of neutral gas pumping. The model adequately demonstrates the qualitative flow characteristics adjacent to divertor and pump limiter neutralizer piates observed in recent experiments 112,131 and simulations [3,6,14,15]. A simplified but physically meaningful approach begins by assuming constant temperature along the field lines (and TS = Ti) and uniform removal of neutral gas atoms at a constant rate. The plasma and neutrals are modelled as counterstrcaming fhtids interacting through electron impact ionization at a rate (uu),. In one-dimensional slab geometry the steady state continuity equations for the plasma density, ni, and neutral density, no, reduce to
q,$ - (oo),nino Vv%-
= 0,
(ou),n,n,-~==o,
where r0 is the neutral pumping time constant. V,, is the 2. PLASMA EDGE MODELLING
136
A.K. Prinja, R. W.
Conn / An axially averaged-radial
plasma flow velocity and V. the neutral fluid velocity projected along the field lines. The boundary conditions are: at z = 0;
PQ(O) = “i0
(3)
at z = a;
nO(a)K3’,=Pni(a)~I
(4)
transport model
10” Slope = 2.29
$
E 2 E
where a is the system size and p( I 1) allows for neutral absorption at the plate. Taking p = 1 corresponds to assuming 100% recycling. Details of the model will be presented elsewhere. Here we focus on the recycling coefficient, R , defined as R =ni(a)/ni(0)
(5)
for a few limiting cases of the parameter E given by (z-s-
70
Lo
7i
X0’
Thus, c -X 1 corresponds to strong neutral pumping and c B- 1 to high recycling. We are also interested in the intermediate range c - 1. For c =z 1, regular perturbation theory to first order yields, for a/L, > 1 R = 1 +p$
I
= 1 +p~~n~~(au),
(7)
10”
10’2 UPSTREAM
ION DENSITY.
10’3 n, (cm-3)
Fig. 1. Density amplification at the target as a function of unperturbed upstream density from analytical recycling model. From fig. 1 we arrive at the following scaling:
or
PIi(u) = ?I$ + C&J”,
“ia = “i0 +pTO(UU),n~O.
where B - 1, a 2 1 and c * 1, and depends on the temperature and pumping time constant. At T, = 30 eV, typical of the plasma edge in a tokamak or the halo plasma in a tandem mirror, we find that a = 1.29 and fi = 1.25. At very high densities (c 39 l), the solution breaks down but an asymptotic evaluation for p = 1, a/L, x- 1 and In R s c yields
(8)
This solution holds for 7. < l/ (ou),n, and shows that for a pumping-dominated situation the plate density increases linearly with the upstream density with the first order non-linearity appearing as the sqaure of the density. The theory may be extended to the range E - 1 by iteration, obtaining f
lnR=~E#nR), I
R=exp(~o(cu),nio) (9
where E2( x) is the exponential integral of second order [16]. For large R (strictly, large In R), and asymptotic expansion yields R = exp[ ( *o/Ti)1’2] = exp[ ( r(eu)enio)1’2]
(10)
Fig. 1 shows ni( a) as a function of rtio for three different values of electron temperatures and 7. = 5 x lo-’ s, obtained by solving eq. (9). The transition from linear to non-linear behavior is now evident. Note that in this model the neutral density at the plate scales exactly as the plasma density and that the neutral pressure will scale with no/~o. The interesting feature is that the model is a reasonable reflection of what has been observed experimentally. In fact, as expected and as is evident from fig. 1, the transition occurs at T~/T~- 1 or Lo/X0 - 1, i.e., when the ionization mean free path becomes of the same order as the pumping scale length. We note here that very clear evidence of such behavior has recently been seen in the PISCES experiment [17].
(11)
(12)
Thus, In R increases linearly with plasma density at large upstream plasma density, niO, as opposed to scaling as the square root of niO, as found in the intermediate range [see eq. (lo)]. For no pumping (c + 00) we reproduce the well known result R = l/(1
-p).
In reality, we do not expect R to increase so drastically with background density because the electron temperature will collapse, preventing further significant ionization of neutrals. The inclusion of temperature effects of course necessitates a numerical approach such as described in the next section. Nevertheless, the salient features of the recycling physics are readily reproducible by simple analytical means even though the results are neither obvious nor trivial.
3. Axially averaged radial transport model Consider the two-dimensional two f&id transport equation, as given by Braginskii [9] for a uniform mag-
137
A.K. Prinja, R. W. Corm / An axiall,, averaged-radial transport model
netic field: Continuify
ar,, !!!+2arax +,,=S,(X,Z,t),
04)
where Fl = - D,an/ax,
05)
r,, = ny,.
06)
Ion energy
a(i + a4Li at
ax
:
aqlli_ _ aZ
v
anr, a~ + Qic + SEi
11
9
Fig. 2. Schematic of edge plasma showing dowstream zones, sources and transport.
(17)
where zi = $nTi + +mnV,t,
(18)
q,,i=(:Ti+jmYI:)nVI,-K,,i~,
(19)
aT
qli=
5
-nXli~-~T,Dl~.
an (20)
Electron energy
ac -1+~+t.!&+$_
at
Qie+ SET,
(21)
upstream
and
gradient free region 1. The parameter L,, is the axial system size and L, is the interface between the two regions. Flow characteristics are symmetric about z = 0. Further, because the dimension of region 2 namely, L,, - L,, is the order of a few neutral ionization mean free paths, we neglect radial transport down stream compared to parallel transport there. However, detailed cross field structure is retained in the upstream region which in turn couples to the core plasma of the tokamak. We now introduce the following averaging operators. In Region 1; (27)
(22)
In Region 2: (23) (24)
where
Also
A = L,, - L, * L,,. (25)
(29)
Applying these operators successively to eqs. (14) through (29, we arrive at the following set of one-dimensional equations: Region I or Upstream Equations:
(26)
Here ci,c represents the total energy of the ions (electrons) and the kinetic energy component of the electron energy has been neglected [eq. (22)]. The quantity, qli,c, is the perpendicular ion (electron) heat flux and q,,i,e the corresponding parallel heat flux. The parameter Q, is the rethermalixation term and the Si represents the particle, momentum and energy sources and sinks. A reduced set of one-dimensional equations will now be obtained for n, Ti, and T,. As will be discussed and explained shortly the momentum equation is not selfconsistently solved for u,,. Based on rigorous parallel flow models and experimental observations, we divide the axial domain into two zones, the upstream region or region 1, and the down stream region, region 2, as shown in fig. 2. These are distinguished by the recycling-driven gradients in the region 2 and a relatively
1
aCei> a(qli> y+p+G41,i(LT)
ax
at =
aw at =
+ (Qic) + (SEi),
-(&W/az) I
ah) ~ ax
+
(31)
~q,,&+)
(vlW”Jaz> - (Qie> + (SE,).
(32)
Region 2 or Downstream Equations: ~+~[r,,(L,,)-r,,(LT)]
~+~[q,,i(L,,)-~lli(LT)]
(33)
=%x,0,
=
-
-YF= -
az
+Q,+s,,,
(34) 2. PLASMA EDGE MODELLING
138
A. K. Prinja, R. W. Conn / An axially averaged-radial transport model
the location of the sheath to L, + +A. This yeilds 2n(T,+
T,)L,,=mnV,f+n(,
(42)
(35) In order to “close” the system of equations [eqs. (30) through (35)], we note three points. First, we clearly need to express averages such as (c), Gier etc., in terms of the averages of the four basic fluid quantities, n, V,,, T and T,. In order to effect this we use the approximation
(fb, Ti, T,, 5,)) =f((‘$ Pi), CT,),(v,,))
(36)
and do likewise for region 2 averages. Alternatively one could specify axial profiles and integrate explicitly but x- and t-invariant axial profiles lead to inconsistenties at the interface when averaged. &. (36) is quite reasonable in region 1 where gradients are gentle, but is probably inaccurate in region 2. Other averaging schemes are being investigated. Second, and perhaps most importantly, the flow velocity in our model is prescribed. Strictly speaking, the momentum equation should be averaged and solved in conjunction with eqs. (30) through (35). Although eventually this will be done, for now we prescribe the flow as a fraction of the local sound speed, i.e., (37) where Vllj = (V) or v, etc., and Mj is the average Mach number in each region. This model is based on results of more rigorous high recycling simulations [8]. Finally, we note in eqs. (30) through (35) the appearance of fluxes computed at the interface L, and the boundary L,,. The system is only closed once. these have been expressed in terms of the average density and temperature. Since there must be a uniform transition of flow from z < L, to z > LT, we make the reasonable approximation of replacing the interface values by their averages in region 1, i.e., n(Lr)=(n);
K(Lr)=(Te);
Ti(Lr)=(Ti).
(38)
In conjunction with eq. (36), this completely specifies the losses from region 1 and sources into region 2. The fluxes at L,, representing losses from region 2 are specified using the sheath model:
+,I) = e&I(LII)~
(39)
4,,i(LII)=2yir(LII)Ti(LII)’
(40)
4,&J
(41)
= 2YJ(qwII)~
where yc = 2.9 and yi = 1 and q,(L,,) is given by eq. (37) with M = 1, that is, sound speed at the plate. It remains to relate the density and temperature at the plate to the adjacent averages. This is effected by momentum balance between the sheath and region 2. We integrate the source-free momentum equation without viscosity from
where we have used sonic flow at z = L,,. Further using eq. (37) we obtain 2n(T,+
q).,,=
(1-t M@(E+
q).
(43)
It is necessary now to make some assumptions about n (L,,) or T,,i( L,,) and obtain the other self-consistently. We assume T,,8(L[I)=ae,iTe,r~
(49
i.e., that the temperature in the sheath is a fraction of the downstream average. The density at the plate is (1+ M;)E(?;,+ n(L,,)
=
pi)
.
2(aC~+aiQ
(45)
The parallel heat conduction fluxes for ions and electrons at L,, are expressed simply as 4lli.e
=
(Z,c- CT.,)) -Kfli,e
q/2
’
(46)
where Klli,= is the parallel conductivity [9] evaluated at the average of (Ti,,) and z,e. Here we have used the fact that the average temperatures are defined at the “ mesh points”, z=L,/2 and z=L,+iA, and the gradient has been linearly approximated between these two points. The term proportional to the electron pressure gradient has not yet been included in the numerical model but will be included in the future. Its neglect is justified in region 1; in region 2 its contribution becomes significant when a large gradient in electron pressure develops. Since parallel electron thermal conductivity tends to maintain a constant temperature, a large density gradient and hence high recycling is likely to be the condition for its importance, although the flow velocity will be highly subsonic in this case. At present we also use a very simple model of recycling by setting the recycling coefficient at the plate to be R = l/(1 - p), as discussed in section 2. Self-consistent neutral calculations will be included subsequently.
4. Results and discussion Numerical results have been obtained for the high recycling mode, for which we have set p = 0.85, M, = 0.1 and M2 = 0.3. The latter values are typical for the ASDEX divertor operation [8]. For this low Mach nurnber flow ion kinetic energy transport can be justifiably neglected compared to the internal energy. The boundary conditions used in the upstream region for radial trans-
A.K. Prinja, R. W. Cvnn / An axial& averaged-radial transport model
port are: atx=O;
(n)(O, 1) = 8 x 1012 cmm3,
-Upstream Temperatures
( qLc)(O, t) = 0.9 W cm-*,
-
t) =O.l W cm-l,
(qii)(O,
at x = a, (first wall/liner); gradient boundary condition. These are values typical of a 500 kW plasma in a machine the size of PLT and it is assumed that 10% of the power flowing into the edge is carried by ions. We have set D, = IO4 = xl0 and xii = 8 x 103. Fig. 3 displays the upstream and downstream radial density profiles and fig. 4 the corresponding electron and ion temperatures. The low Mach number flow in region 1 results in an almost linear upstream density profile, but roughly corresponds to an exponential with an e-folding length of 5 cm. The radial variation of the upstream ion temperature is also not exponential but is fairly linear up to about 3 cm into the scrape-off layer. The reason is that ion energy transport is dominated by convection which, in region 1, is at low Mach number (about 0.1). The relatively high xii then causes a flatter temperature profile radially. The electron temperature, on the other hand, is dominated by parallel conduction, and the large sheath losses result in a sharper than exponential radial drop out to about 3 cm. Beyond this point, it follows the ion temperature. Here we see an interesting behavior in that as the electron temperature decays radially, the energy transport mechanism shifts from being dominated by conduction to being largely convective flow. Plotted on
.
\
Upstream Plasma
Denritl
-----Downstream Plasma
Densit]
\ \ \ \ \
J 0
0.8
DISTANCE
1.8 2.4 INTO
3.2
4.0 4.8
5.6
SCRAPE-OFF
6.4
7.2 8.0
LAYERx
km)
Fig. 3. Upstream and downstream radial density profiles for high recycling and low Mach number flow.
139
0
0.8
DISTANCE
1.8
2.4
INTO
3.2
- --Downstream Temperatures
4.0
4.8
SCRAPE-OFF
5.8
8.4
7.2
LAYER.x
8.0 (cm)
Fig. 4. Upstream and downstream radial electron and ion temperature profiles for high recycling and low Mach number flow.
a semi-log scale, we find that the initial radial decay has
an exponential decay length of 1.65 cm, while further out the tail decays at 7.2 cm. Perhaps the most fascinating result is the structure of the downstream profiles. Except for the electron temperature, these profiles do not resemble the upstream profiles. Over the first 3 cm radially, the density profile is practically flat and the ion temperature actually exhibits a broad maximum and then becomes flat at about 2 eV. This happens in a region where the electron temperature falls very sharply (e-folding at 1.5 cm), and is the cause for the peculiar ion temperature and density behavior. Because there is a large difference between the electron and ion temperature over this radial extent, the flow speed is determined by T, and thus also decays sharply with radius. This means that the convective losses, which dictate both the particle and ion energy balance, are significantly reduced causing the profiles to flatten. This precise situation might not arise at a higher Mach number in region 2 but it seems that the sharply decaying electron temperature will always occur as long as conduction dominates transport along field lines. This will definitely influence !&, Ti profiles through convective transport. We also note that radiative cooling of electrons in region 2 through collisions with neutral atoms and molecules will also affect the local electron temperature. More detailed investigations along these lines are underway. 2. PLASMA EDGE MODELLING
140
A. K. Prinja, R. W. Conn / An axially averaged-radial
Acknowledgment
This research has been supported by the United States Department of Energy, Office of Fusion Energy, under contract DE-AM03-76SFOOO34, P.A. DE-AT0380ER52062.
References [1] F. Wagner et al., Phys. Rev. Lett. 49 (1982) 1408. 12) M. Petravic, these Proceedings. [3] M. Petravic et al., Phys. Rev. Lett. 48 (1982) 326. [4] B.J. Brahms. Proc. 11th European Conf. Controlled Fusion and Plasma Physics, Vol. 7 (1983). [5] H.C. Howe, J. Nuci. Mater. 93 & 94 (1980) 17. (61 D. Post et ai., J. Nuci. Mater. 121 (1984)
transpon model
171 M. Shimada et al., J. Nuci. Mater. 121 (1984) PI W. Schneider and K. Lackner, Intern. Conf. on Plasma Physics, Goteborg (1982). 191 S.I. Braginski, Review of Plasma Physics, Ed. M.A. Leontovich (Consultants Bureau, New York, 1965) WI P.J. Harbour and J.G. Morgan, Cuiham Report, CLMR234 (December 1982). 1111 W.A. Houiberg et al., Proc. ANS 1st Meeting on Advances in Mathematical Methods for the Solution of Nuclear Engineering Problems, Munich, FRG, II, 425 (April, 1981). 1121R. Budny et al., J. Nuci. Mater. 121 (1984) in press. 1131S. Taimadge et al., J. Nuci. Mater. 111 & 112 (1982) 274. 1141A.K. Prinja et al., J. Nuci. Mater. 111 & 112 (1982) 279. I151 R.W. Conn et al., J. Nuci. Mater. 121 (1984) 1161 M. Abramowitz and I.A. Stegun, eds., Handbook of Mathematical Functions (Dover Pubi., 1972). P71 D. Goebei and R.W. Conn, these Proceedings.