Magnetic reconnection and m = 1 oscillations in current carrying plasmas

ANNALS OF PHYSIC~

112, 443-476 (1978)

Magnetic Reconnection and m = 1 Oscillations in C u r r e n t Carrying Plasmas G. ARA, B. BASU, AND B. CoPPI Massachusetts Institute of Technology, Cambridge, Massachusetts 02139 AND G. LAVAL,* M. N. ROSENBLUTH, AND B. V. WADDELLf Institute for Advanced Study, Princeton, New Jersey 08540 Received August 1, 1977

The problem of magnetic field reconnection and the onset of instabilities resulting from it is analyzed for current carrying plasmas in high temperature regimes of thermonuclear interest. In order to compare our theoretical analysis of modes having azimuthal wave number m = I with present experimental observations and to predict those in future experiments producing high electron temperatures, we have developed a model that includes simultaneously the effects of finite electrical resistivity, ion gyroradius, electron drift wave frequency and ion-ion collisions. The viscous term arising from the last effect prevents the considered modes from developing a radial structure with typical "wavelength" shorter than the ion gyroradius and makes the moment equation description adequate for the regimes of interest. The relevant growth rates are found to be strongly decreasing functions of the electron temperature and, as the temperature increases, tend to achieve a form similar to that previously derived for higher m modes which can be found only in the presence of finite electrical resistivity. The relationship of these modes with the m = 1 ideal MHD internal kink mode is studied and the effects of toroidal geometry are briefly discussed.

I. INTRODUCTION A current carrying plasma c o l u m n in which r e c o n n e c t i o n o f the magnetic field lines is allowed by finite plasma resistivity is subject to c u r r e n t driven instabilities associated with it. Instabilities of this type have, in fact, been p r o p o s e d to explain solar flares [1] a n d have been identified as p l a y i n g a role in the d i s r u p t i o n [2] a n d minidisruptions, identified as sawtooth oscillations [3] o f soft X-ray emission, i n magnetically confined plasmas. Both in the case o f solar flare models a n d o f the experimentally observed minidisruptions, a current driven m o d e is excited following * Permanent Address: Laboratoire PMI, Ecole Polytechnique, 91128 Palaiseau Cedex, France. t Present Address: Oak Ridge National Laboratory, Oak Ridge, Tennessee 37830. 443 0003-4916/78/1122-0443$05.00/0 Copyright © 1978 by Academic Press, Inc. All rights of reproduction in any form reserved.

444

ARA ET AL.

a thermal instability that raises both the electron temperature and the current density in the central region of the plasma column. In the present paper we mostly consider a cylindrical plasma column of finite length and limit our attention to modes with both azimuthal and longitudinal wave numbers equal to unity. Modes with these wave numbers can be found on the basis of the ideal M H D approximation, where no reconnection of magnetic field lines is allowed, and are called internal kink modes [4]. However, the experimental observations, as well as the theoretical indication [5] that the internal kink modes are stable in toroidal configurations under realistic conditions has led us to investigate the properties of modes which owe their existence to the reconnection of the magnetic field lines and lead to a redistribution of the current density, that is necessary in order to explain the sawtooth oscillations mentioned above. Our presentation is organized as follows: In Section II, we describe the main features of the confinement configuration we consider and analyze the limits within which the ideal M H D approximation can be applied. In Section III, we study the effects that the introduction of finite electrical resistivity has. Reconnection of the magnetic field lines is produced and, as a result, magnetic islands are formed. The growth rate of the internal kink mode analyzed in Section II is significantly increased near the condition for ideal M H D marginal stability while a new mode that is labeled as "reconnecting" mode is found in the limit where the ideal M H D approximation would only predict a singular radial plasma displacement, i.e., with locally infinite amplitude. The results of the analysis of the nonlinear evolution of the considered instability are also discussed. In Section IV, we introduce the set of moment equations that include all the nonideal M H D terms besides electrical resistivity that need to be considered. The main effects to be taken into account are, in fact, those of finite ion gyroradius, finite electron drift wave frequency and ion-ion collisions. In Section V, the equations for the plasma perturbations that can take place in the plasma layer where magnetic reconnection occurs are given for high temperature regimes of thermonuclear interest. In Section VI, we produce a solution of the problem neglecting for simplicity the effects of ion-ion collisions. We find that the growth rates of both the resistive internal kink mode [6] and of the reconnecting mode [7] are drastically reduced. However, the spatial structure of these modes tends to develop radial oscillations with scale distances of the order of or shorter than the ion gyroradius. In Section VII, we introduce the effects of ion-ion collisions in addition to those mentioned in the previous section. One of the effects is to nearly eliminate the fine radial structure of the modes found in the previous section and, thus, to restore the validity of the adopted moment equations and of the results obtained from them. The growth rate is further reduced while a (real) frequency of oscillation is found even in the frame of reference where no radial electric field exists. The phase velocity of the resulting modes is in the direction of the electron diamagnetic velocity in agreement with the experimental observations.

MAGNETIC RECONNECTION AND m =

1

MODES

445

In Section VIII, we discuss our results and compare them with those that are presently available from the experimental observations. The presented theory can, in fact, explain the magnitude of the measured frequency of oscillation, the observed direction of the phase velocity and mode radial profile. Other features, such as the lengthening of the period of sawtooth oscillations and the reduction o f their amplitude as the plasma density increases, are shown to be consistent with our considerations. In fact, m = 1 modes are observed on top of these oscillations just before they reach their maximum. We also point out that, in view of the sharp decrease of the growth rates of the considered modes as the electron temperature is increased, we expect that magnetically confined plasmas to be produced in future experiments with higher electron temperatures than those realized so far will be considerably less affected by them.

II. IDEAL M H D THEORY We consider at first a current carrying cylindrical plasma column with B = Bo(r) eo + l~(r) e z ,

(11.1)

representing the azimuthal and axial components of the magnetic field. The length o f the cylinder is L = 2~rR, and we shall refer to this configuration in order to simulate an axisymmetric toroidal one having a relatively large aspect ratio with major radius R. This simulation has a number of shortcomings that we shall attmept to identify as we proceed. We assume that the current distribution is such that q(r) = r B J R B o ,

(II.2)

is smaller than unity for r < r0. Here q(r) : 27r/~(r), ~(r) is the rotational transform. We refer in particular to low-fl plasmas in which Bz ~ constant over the plasma column, Bo 2 ~ B~ 2, and thus the variation o f q(r) is only related to that of the current density J~ :

c d 4rrr dr (rBo).

(11.3)

We consider perturbations from the equilibrium of the form = ~(r) exp(--icot + imO q- ik~z),

(11.4)

and limit our attention to modes for which ks = 1/R and m = --1. If we indicate the relevant propagation vector by k : (m/r) eo -+- k,e~, we have k • B = (mBo/r + k~Bz) : (q -- 1) Bo/r,

and k • B = 0 for r -----ro where q(ro) :

i.

(II.5)

446

ARA ET AL.

If we describe the considered perturbations by the set of ideal M H D equations, an instability can be found, the so called internal kink mode [4], with its growth rate expressed by

5' : (VAo/ro) An,

(11.6)

Bo d In q V~o -- (4~p)a/2 d In r '

(II.7)

where 5' = Im to,

evaluated at r = r 0 , and Art is a parameter of order

AH "~ (ro/g) 2,

(II.8)

such that the ideal M H D mode can be found for AH ~ 0, and AH ---- 0 corresponds to the ideal M H D marginal stability. We recall that one of the most restrictive equations of the ideal M H D approximation is the so-called "frozen-in-law," E + (1/c) V × B ---- 0

(II.9)

which compels the fluid to move with the magnetic field lines and prevents their reconnection. Given the order of magnitude of An, the results obtained in cylindrical geometry cannot be strictly applied to the toroidal configuration that is simulated. In fact, it was first indicated in Ref. [5] that the ideal M H D internal kink is stable under realistic conditions when evaluated for a toroidal configuration. In carrying out the stability analysis we distinguish two asymptotic regions: an inner layer of width 3r 0 < r 0 , centered around r ---- ro, where the complete dynamical equations are used; and, an outer region described by the ideal M H D equations in which inertia is neglected. The solution in the inner layer, where the appropriate inertial term is retained, is asymptotically matched to that in the outer region and thus related to the external value of the perturbed potential energy. We assume that the value of this is obtained by a perturbative method involving a small parameter ca, which can be either the ratio of the transverse to longitudinal wavelengths, the inverse aspect ratio, or a measure of noncircularity of the cross section. In all cases, the relevant calculations have to be carried out to second order in ca, the system being marginally stable at the zeroth and first orders. Then we assume that the p e r t u r b e d potential energy can be written in the form ~W :

3 W 0 --~ E1 3 W 1 -~ ,12 3 W 2

with ~W0 : J'g "(Hog)dr, ~Wa : f g "(Hag) dr, and 3W2 ---- J'g "(H2g)dr, where g is the perturbed displacement, Ho, Ha, and H~ are Hermitian operators. The lowestorder contribution 3 Wo reduces to the well-known cylindrical form 3Wo ---- ~ r3(k • B0)~

dr,

(II.I0)

MAGNETIC RECONNECTION AND m = 1 MODES

447

where a is the radius o f the p l a s m a column. The m i n i m u m value 8 Wmin. o f the p o t e n t i a l energy is f o u n d to be given b y

~l.(Hlgo) dr -~- f go.(H2~o) dr],

3Wmin = e12 I f

(II.1l)

where T0 is a p u r e r a d i a l displacement, which is c o n s t a n t for 0 < r < r o a n d is zero for r > ro • The first-order d i s p l a c e m e n t has to be a solution o f H0gl + Hlgo :

0,

(II.12)

a n d the s e c o n d - o r d e r s o l u t i o n is given b y H 0 ~ + HI~I + H ~ o = 0.

(11.13)

I n t e g r a t i n g Eq. (II. 13) f r o m r = 0 to r = r < r o we find

r3(k • Bo) 2 where we have t a k e n ~:o = ~:~o : we find

(~)

~E12s% ['Vmin '

--

(II.14)

constant for r < r o . O n the o t h e r hand, for r > ro,

d~2~. __ dr

const rS(k

01.15)

• B0) 2 "

Therefore, the o u t e r solution correct to second o r d e r in E1 c a n be written as

= ~0 + ~ ? ~ , where, for r < ro, d~2

~:0 = ~:~o = const;

dr

~ Wmin

~r3(k • B0) 2 ¢1z '

(II.16)

and, for r > r0, ~o :

0;

d~z _

dr

SWmin(r : ro) ~ j a ( k . Bo) 2 ~a~ "

(II.17)

H e r e we have used the n o t a t i o n ~:~ = ~: a n d chosen the c o n s t a n t o n the r i g h t - h a n d side o f Eq. (II.15) to be a suitable f o r m for j o i n i n g the solutions across r = ro. N e a r the inner layer, the outside ~: is then given b y

~(X) = ~ = --c/x, 595/112/2-x5

C/X,

for

X --~ 0 _ ,

(II.18)

for

x --+ 0 + ,

(II.19)

448

ARA ET AL.

where x -~ (r -- ro)/ro, e = 8Wmin(r = ro)/[~roS{(k • B0)'}2], (k "Bo)' ~ (d/dr) × (k • Bo) evaluated at r = ro, and we have used the fact that k • B 0 = 0 at r = r o . Choosing the solution of the inner layer equation, where the inertia term is included, in such a way that ~ ~ ~:o~ as x ~ -- ~ and ~ --~ 0 as x ~ ~ and then equating d~/dx obtained from Eq. (ILl8) with d~/dx obtained from the inner solution for x - + --oo yields the growth rate given in Eq. (II.6) where now An can be given explicitly in terms of 8 Wmin by ~tH =

--'rr ~ Wmln(r =

ro)/(~ooroBoq')Zr=~0 •

(II.2o)

We notice that the instability occurs if 8 W m i n < 0 for given profiles of pressure and current density. We observe that in toroidal geometry we cannot neglect the contribution of the kinetic energy resulting from the displacement along the toroidal field. This contribution, a m o n g others that are due to toroidal effects, is included in the relevant toroidal expression [5] for 8 Wmin.

I I I . RESISTIVE THEORY

We consider briefly the effects of finite electrical resistivity simply by replacing the "frozen-in-law," Eq. (II.9), with the generalized Ohm's law

E + (1/c)(V × 13) = ~J

(IliA)

where ~ is the electrical resistivity. The main effect of finite plasma resistivity is to produce reconnection of the magnetic field lines around the surface r = r0 at which q(r) = 1 and, as a result, magnetic islands are formed. In addition, (a) a significant increase in the growth rate of the ideal M H D internal kink mode is realized in the range of parameters where this is predicted to exist, i.e., f o r h n ~>0, (b) the existence of a new unstable mode, which we call the "reconnecting" mode, is made possible in the regime An < 0 where the ideal M H D theory does not predict any instability. Both the internal kink mode [6] and the reconnecting mode [7] are driven by the magnetic energy in the poloidal field that is converted into plasma thermal and flow energy through an inner layer centered around r = r 0 . A . Resistive Internal K i n k M o d e

When Eq. (III.1) is used, instead of Eq. (11.9), the equations which describe the dynamics of the plasma within the inner layer can be written as [6] ~t2 d2~/dx z = x d~b/dx 2,

(III.2)

~b =

(III.3)

x~: + (e/A) d2~b/dx 2,

MAGNETIC RECONNECTION AND

m = 1

MODES

449

where x ~ (r - - ro)/ro, ~ is the radial displacement, 4, ~- iBlfl(roF' ), BI~ is the radial c o m p o n e n t of the perturbed magnetic field, F ' =- (Bo/r)(dq/dr)I~=~ ° , A ~ - - i o J z u , =-- "ru/'rn , "cu ~ to~ V Ao is a typical Alfv6n time, "rR ~ (4rrro2/~c 2) being the resistive diffusion time. F o r the derivation of Eqs. (Ili.2) and (111.3) see Section V. We introduce the small p a r a m e t e r 3 < 1, such that x ~-~ 3. Then the "resistive" ordering for which all terms in Eqs. (Ili.2) and (111.3) are comparable, and which describes the resistive internal kink mode, is (III.4) F r o m this ordering we immediately find a ~ 3 ~ El/3. We look for solutions of Eqs. (III.2) and (III.3) that can m a t c h smoothly with the relevant M H D solutions outside the inner layer. I f we reduce Eqs. (II1.2) and (III.3) to a fourth order equation in {:, there are four independent solutions with different asymptotic behavior as I x I --~ oo. One solution is ~: = constant, two others behave as exp[4-x2/(4EA) '/2] and the last solution is such that d~/dx ~ (const)/x ~ as I x / --+ ~ . We are interested in the last solution which we write in the f o r m ~: = ½~:o~+ ~:ooo(x), where the function ~:oad is odd in x. Notice that the f o r m o f the solution outside the inner layer given in Section II implies that 1

dE

~-,o

~® dx

--

1

_

An

(1II.5)

7r x 2 "

Thus we impose that the inner solution satisfy the condition x~ d

1

- -77" a n ,

--

2 d x (In ~odd) ~

(I11.6)

for X --+ - - c~ as ~odd ~ {:~o/2. It is possible to find the appropriate analytical solution of Eqs. (III.2) and (III.3) by introducing the even function: x(x)

=-

x d4,1dx - - 4, = a S d ~ l a x + xo~ ,

where X~o = const. Then ~ = -1-2

dx(x - - x~)

and

4, = - - X - - x

In addition ~:® ~ 2 J'2 (dx/dx) dx/x and d~/dx ~ Xo~ = (21H/Tr)

dx(dx/dx)/x.

- - X ~ / X L Thus Eq. (IIL6) reduces to

.f~ dx

dx

dx

x

(III.7)

and Eqs. (III.2) and (III.3) to

E,)l ( d 2 x \

dx 2

2 d~_xI _ ( x 2 + A2) X = - - x ' x ~ .

x

~.Tt, /

(II1.8)

450

ARA ET AL.

We next introduce k ~ x13, and put 3a = EA, then Eq. (III.8) becomes, d2x

2 dX

d~~

~ d~

(~2 + ,~z/2) X :

--x2x~o ,

(III.9)

where A --~ ~/E113, The solution of Eq. (111.9) can be presented in the form of an integral representation (see Appendix A) as XIx~ :

1 -- --2~3/2 fol dt(1 -- t)(a31214)-514 (1 + t) -d31214)-51a exp(-- ½t&z)

(III.10)

provided Re ~3/2 > 1. The eigenvalue equation can be obtained by evaluating the integral in Eq. (III.7), which gives

+ 5)/4)

~H ~- )tH/Ell3.

(III.11)

We can examine some of the limiting cases. In the M H D marginal stability case, i.e., An = 0, Eq. (III.11) yields h

=

E 1/'3.

(II1.12)

The corresponding ~:, consistent with the boundary conditions, is = ½~o[1 -- erf(k/21/2)].

(111.13)

Thus we see that the finite resistivity destabilizes the otherwise stable ideal M H D internal kink mode. For 0 < ~n < 1, Eq. (III.11) gives ~_ 1 + ~Ai4~li 2,

and, in the ideal MHD limit (E ~ 0), i.e., AH >~ 1, we recover the results of Section II, namely, A = AH, X : X~x2/( x2 + h2) (III.14) ~: : ½~:~[1 -- (2/~r) arctan(x/h)] B. Reconnecting M o d e [7]

We consider the possibility that An can be negative in toroidal configurations and introduce two small parameters: 3 < 1, such that x ~ 3, and 3/] hH[ < 1. The asymptotic ordering by which the lowest-order forms of Eqs. (III.2) a n d (III.3) describe the reconnecting mode is (1/qS) d2~/dx 2 ~ ([ ~HI ~)-L

(III.15)

MAGNETIC RECONNECTION AND

m = 1

MODES

451

Equation (III.15) implies that we can write ~b = ~bo ~ ~bl, where ~b0 = constant, and ~b~(x)/¢o ~-~ 8/I An I in the sense that (d~bddx)/4J o ~.~ ] AH ]-x. By the term by term inspection of Eqs. (III.2) and (III.3) we find, to lowest order, ¢/¢o ~ 3-1

(III. 16)

e/h ,~ [ h n l 8,

(IXI.17)

k 2 ~ 3z/I ass I.

(III. 1 8)

Combining Eqs. (II1.17) and (III.18) we obtain the scalings of a and 3 in terms o f and l assl as it ~-, ~SlSII Ass I~i5,

,~I511 Ass I*15.

a ,-~

/'/,'

Then, 8/I ass I < 1 implies I ass I > O/3 and a < O/~.

soft x - ray

:

1

I

,

lJk,

J I

1

2 ms/div FIG. 1. Sawtooth behavior of so•X-ray emission from the hot central core of toroidal discharges. Taken from Ref. 30. Reprinted with permission of the authors, T. F. R. Group, Nucl. Fusion 17, 1283 (1977). We can now rewrite Eqs. (III.2) and (III.3), keeping only the lowest-order terms, as dx 2

--X

dx 2 E

d2¢1

(III.19) (III.20)

The subscript --1 on ~ has been used to indicate that ~--1/~0 ~ 3-1" F r o m Eqs. (III. 19) and (III.20) we obtain the following inhomogeneous equation for ~:-1 A2 d2~--1

dx~

(a/,) x~#_~ = (al,) Gx.

(III.21)

452

ARA ET AL.

We have to find the solution to Eq. (111.21) that matches properly with the solutions valid outside the inner layer. The form of the solution outside the inner layer implies that as we approach the inner layer from either side, i.e., as I x [ ~ I(r -- ro)/ro I -+ O, = ~:~[1 + (2tH/rr) X-l],

as

x --~ 0_,

= ~(An/'rr) x -1,

as

x -+ 0+.

For [ x I -~ 0, the perturbed radial magnetic field outside the inner layer, ¢, is related to ~ as

4, = - x ~ , and, hence, ¢ = --~o[()tH/~r) + :'Co],

as

x ~ 0_,

¢ = --()tH/rr) ~®,

as

x -+ 0+.

The discontinuity in ( d ¢ / d x ) / ¢ across x = 0 is, therefore, given by

A' - (d¢ldx)l¢ l.-.o+ -- (d¢ldx)/¢ I~o_ --------~/a,,.

(III.22)

The appropriate matching condition for the solutions of the inner layer equations is that (d~bJdx)/¢s 1 ~ 8 = ( d ¢ / d x ) / ¢ ].~.o+

and

(dGIdx)lG 1~->-~ = (d¢oldx)l¢ I~*o_, where ~b and ¢~ are outside and inner layer solutions respectively. Combining the two conditions we get x=+½~ = A' = [(dCsldx)/~bs]~=_~,a Using the asymptotic expression Cs = ¢0 + order, as

f+~ (d~Cddx~) dx

¢1 we

-

-

rrl)tn .

(111.23)

can rewrite Eq. (III,23), to lowest

--(~1~.) ¢0,

(111.24)

where we have extended the limits of integration to infinity since ¢1 is nonzero only within x ~ 8. Using Eq. (III.20), Eq. (III.24) becomes '~oo

(a/c)

f~

oo

[1 + x&ll¢o] dx = --~la,~.

(Ill.25)

MAGNETIC RECONNECTION AND m = 1

MODES

We introduce the dimensionless quantities ~ ~ x / 8 , 84 ---- EA. Then, Eq. (III.21) becomes

453

Y --= (~/¢0) ~:-1, and set

d 2 y / d ~ 2 - - ~ y = ~.

(1II.26)

We can give the explicit solution of Eq. (Ill.26) by an integral representation of the form /t

1

Y(~) = -- ½~ J0 dt(1 - - t2) - 1 / ' exp(-- ½tk~).

(1II.27)

The matching condition, Eq. (III.25), in terms of Y becomes d,,i

(111.28)

oo

Substituting Eq. (II1.27) into Eq. (Ill.28) and using the definition 84 = cA we obtain, after integration, the eigenvalue equation 1 = --(2x/~,~u/zr)[F(3/4)]2 as/%-W',

0II.29)

which yields, for ha < 0, i.e., A' > 0,

A = A0~3/5 I AH I-4/5,

(III.30)

where, A o =- [21/2{['(3/4)}2/~r] -4/5 ~ 1.37. Using the definition 84 = cA, we obtain = A~/',2/51

A,, 1-1/5.

(III.31)

In Fig. 2, we have presented the radial profile of ~-1 • We notice that the profile is similar to that of a resistive tearing mode [8, 9] that is found for m >~ 2. We also present, in Fig. 3, the corresponding profiles for the (resistive) internal kink mode Y

FIG. 2. Radial profile of the lowest order radial displacement (f-l) for the pure resistive reconnecting mode: Y[--=~/~o)~:-1]vs ~(~- x/g).

454

ARA ET AL.

that is found for An ) 0 in order to illuminate its difference with the resistive reconnecting mode. In order to connect the results we have just presented with those of the nonlinear analysis reported in Refs. [10-13], we introduce the vector potential A such that B -----V × A and recall that 3 A / O t = - - c ( E + V~)

(III.32)

where ¢ is the scalar potential. I f eh indicates the unit vector in the direction defined by the considered helical perturbation, we have (III.33)

e~ = e~ - - ( k ~ r / m ) eo ,

and,

eh • V q ~ =

0.

I-

ro~----d

&-Go

i\bi 1 o

v'Qq ..... ro

r

FIG. 3. Radial profiles of the radial displacement (D for the resistive internal kink mode: (a) )tn = 0; (b) ~n > 0. The vertical dotted lines indicate the approximate boundaries of the inner layer. We introduce the magnetic flux 7 t through the helical ribbon defined by the axis of the cylinder and the helix intersecting the point (r, 0, z) such that 7 j ~ eh • A. Then, we have ~7"t/8t = --c(% "E).

(III.34)

Notice that the flux function 71 introduced here is proportional to the function ~b introduced earlier in this section. In fact, 7"t(x) = B o ( r o ) ( d q / d x )

[~=o ¢(x).

I f C is an arbitrary vector we have % • (C X B) = - - C • VW. Hence, from Eqs. (III.1) and (III.34) we obtain dttt /dt =

- - ~lc(en • J)

(III.35)

rn = 1

MAGNETIC RECONNECTION AND

MODES

455

and since ej, ~ e.(k.r ~ m) in our case, we have dT-' / dt =

--

,~c]z

011.36)

,

where

d~P/dt =- OTt/~t + V • V7j. From Ampere's law, the current density J~ is given by Jz =

4~" Vx2~v +

m

'

(III.37)

where 1 O(O) lO S V ± 2 - - r O ~ r~rr + r 2 002" Then, Eq. (III.36) becomes dW -- */c2 (V±2W +

dt

47r

2k~B~ ] -----m--i"

(Ill.38)

We note that, for Bo ~ B~ and k~r ~ m, V~ may be neglected and we can take V • V = 0. Consequently, we introduce the velocity stream function S such that V = VS × e~.

(111.39)

Then, applying the operator e. • V x to the momentum balance equation

p d/dt V = --Vp + (J × B)/e,

(III.40)

and retaining only the lowest order terms in the expansion parameter k~r, we obtain

4rrp d/dt V±2S : --e~ • (V±~ x V±V±2~),

(III.41)

where p is the plasma mass density, d/dt denotes the convective time derivative, and V Ldenotes spatial derivatives perpendicular to z. The linearized forms of Eqs. (III.38) and (III.41) can be used instead of Eqs. (III.2) and 0II.3) in order to describe the considered modes.

C. Nonlinear Analysis The nonlinear evolution of the internal kink mode has been studied numerically in Ref. [12] using Eqs. (III.38) and (III.41), and restricting to perturbations having a fixed helical symmetry. This implies that the considered perturbations from a cylindrical configuration are functions of r,t, and the helical coordinate r =~ mO + k~z only.

456

ARA ET AL.

The principal results of Ref. [12] is that the m = 1 island growth is very rapid even in the highly nonlinear phase. In fact, the island continues to grow exponentially at more or less the linear growth rate until the current is completely flattened within the q = 1 singular surface. This result confirms the heuristic argument of Kadomtsev [I0]. For numerical integration, an alternating-direction implicit scheme was used to advance Eq. (III.38) in time while Eq. (III.41) was advanced explicitly. Although the experimental value of c 1 is approximately 5 x 108, the code was run for E-1 = 5 x 10~ in order to avoid numerical problems associated with a small inner layer width. Nevertheless, since the mode growth rate is proportional to Ewa, the discrepancy between the two values of ¢ should cause no qualitative difference in the results. The safety factor profile was chosen so that q at the origin was 0.9, q at the limiter was 3.4, and the radius of the q = 1 "singular" surface was 0.2a, where a is the minor radius. For the parameters specified in the preceding paragraph, the logarithm (base ten) of the fluid kinetic energy is plotted as a function of t[~rR in Fig. 4. Initially, the system is given a m ----- 1 perturbation such that the m a x i m u m width of the magnetic island is 4 x 10-~a. The linear growth rate calculated from the initial slope of the curve in Fig. 4 agrees with the value obtained from Eq. (III.12) to within 16 %. In summary, the numerical results show that nonlinearly the kinetic energy grows exponentially at approximately the linear growth rate until a m a x i m u m of 1 eV per particle is reached. (Kadomtsev's heuristic theory [10] predicts that about 3 eV per particle should be released by the magnetic field). By the time the kinetic energy is

I

i

l

i

i

t ~-7

£ -8

I

f

I

TIME ([O-BTR)

FXG. 4. Logarithm (base 10) of the fluid kinetic energy as a function of time.

MAGNETIC

RECONNECTION

AND

m

~

1 MODES

457

maximum, the axial current density has flattened inside the singular surface and the safety factor at the plasma center has increased to approximately unity. After reaching a maximum, the kinetic energy decays by a factor of 5 at an average rate of about 1/16 of the linear growth rate. In Fig. 5, the current density is plotted as a function ofr along the line going through the x-point and the center of the magnetic island. At time t = 1.92 × 10-3"rR , 2

50 - -

I

i

,

~

i

,

i

i

i

i

I

I

I

I

I

I

I

I

l

I

I

I

I

I

I

1

I

I

I

f

I

I

I

I

I

I

i

200 150 I O0

50 0

200

1501 f0C .5C C

~ _ _ J

20C

1.0C .5C C

I

I

20C 1,5C

10C

.5C -II20

F I G . 5.

-080

-040

000 RADIUS/rw

0.40

Ol80

Axial current density as a function of r at selected times (normalized to the skin time).

the current begins to flatten in the center of the island while a skin current develops at the x-point. When the kinetic energy reaches the maximum, the current is flat through most of the plasma interior and the skin current is quite large. Then, as the kinetic energy decreases, the skin current disappears and the current through the plasma interior remains fairly flat. The magnitude of the current at the plasma center corresponds to a q of approximately unity. (Notice that the total current in Fig. 5 does not remain constant because in the code the electric field at the wall, rather than the total current, is held fixed). The time required for the flattening to occur is essen-

458

ARA ET AL.

tially the time for the kinetic energy to reach the maximum value, i.e., the time for the magnetic island width to increase to approximately 2r0. In the results presented in Ref. [12], the time evolution of the electron temperature was not included. Recent results show that the electron temperature profile also becomes flat inside the q = 1 singular surface in a way that is similar to the evolution of the current density. This process occurs because the heat, which diffuses rapidly along the magnetic field lines, flows from the plasma core to the region around the singular surface as the magnetic island grows. Furthermore, if the evolution of the electron temperature is included in the code, the flux contours become circular once the kinetic energy reaches a maximum, in contrast to the complicated structures presented in Ref. [12] and in agreement with Kadomtsev's heuristic argument [10]. Finally, it should be emphasized that the preceding result can be affected by finite gyroradius and toroidal corrections, which have not been incorporated into the present numerical scheme.

IV, T w o SPECIES TRANSPORT EQUATIONS

We adopt the following set of equations [14] to describe the dynamics of a collisional plasma comprising of electrons and ions: an~/at q- V • (n~V,) = 0,

(IV.I)

m~n, d/dt 1~ V~ ~- - - V p , - - V . rr -~ e~n,(E q- V × B/c) q- R~,

(IV.2)

~n~ d/dt I, T~ q- p~V • V~ = 0,

c(V X E) =

--~B/~t,

e(V X 13) ----4zrJ,

(IV.3) (IV.4)

(IV.5)

where p . = n~T~,

and

d/dt I~ ~ O/St + V . • V.

In these equations, n~ is the particle density, V, is the fluid velocity, m~ is the mass, p~ is the scalar pressure, e~ is the electric charge, T~ is the temperature (measured in energy units) rr~ is the stress tensor and R~ represents the collisional momentum transfer for a given species ~ [~ = e(i) for electrons (ions)]. The other quantities have their usual meanings. Adding the electron and ion equations of motion and neglecting the electroninertia and electron-stress tensor, we obtain the plasma equation of motion p[~Vi/~t -k (Vi • V) Vii ---- --Vp -- V • rri -1- (J × B)/c,

(IV.6)

where p ~ men~ -k m i n i ~- min(n~ = n i ~ n ) is the plasma mass density, p = p, + p, is the total pressure, J = en(Vi - - V~) is the total current density. We have used the relation Ri = - - R , . Furthermore, neglecting electron-inertia and electron-stress

MAGNETIC RECONNECTION AND m z

1 MODES

459

tensor (valid when electron-ion collisions are sufficiently frequent), the electron equation of motion yields the generalized Ohm's law (since it relates the current to the electric field): ~7,,J, 4- ~j_J. = E 4- ( V e x B ) / c 4 - V p ~ / e n 4 - ( O . 7 1 / e ) V , , T ~

OV.7)

where, we have used R, = en(~l~J , 4- ~ / . J ± ) - 0.71nV,~T~, ~7's are electrical resistivities and the symbols II and _L indicate the component parallel or perpendicular to the magnetic field, respectively. Equations (IV.l) and 0V.3)-(IV.7) constitute a closed set of equations suitable for studying low-frequency and long-wavelength phenomena including the various non-MHD effects, such as finite electrical resistivity, finite ion gyroradius, finite drift wave frequency and ion-ion collisions.

V. PERTURBED LINEARIZED EQUATIONS

We consider a low-fl plasma and refer to an equilibrium configuration in which the equilibrium quantities depend only on r. Since we intend to consider the effects of finite ion gyroradius and relevant diamagnetic velocity, we discuss briefly the equilibrium state as described by the equations - - d p d d r - - en[E~ + (1/c)(V~oB~ - - V~Bo)] = O,

(v.1)

- - d p d d r 4- en[E,r 4- (1/c)(VioB~ - - Vi~Bo)] -~ O,

iv.z)

where Er is the equilibrium radial electric field and the equilibrium magnetic field is B 0 ---- eoBo(r) 4- e~Bz(r), B o ~ B~. I n order to simulate the equilibrium of a toroidal configuration, usually considered for the neoclassical transport theory, we may take [15] Viz =

( c / B o ) [ E r - - (1/He) d p i / d r ] ,

(V.3)

Vio = O,

(V.4)

Vez = Viz - - J~/ne,

(v.5) (v.6)

Veo - - --(e/neB~) dpddr - - e(Er/B~).

We can carry out our analysis in the frame of reference in which E,. = 0, i.e., the frame of reference moving relative to the laboratory frame of reference with a velocity cEr/Bo in the z direction. In this frame of reference, the equilibrium velocities do not have the E x B0 drift components and the frequency oJ is replaced by the Doppler shifted frequency (o'

=

o~ -

k~(cE/So).

(V.7)

460

ARA

ET AL.

Thus, except for the appearance of an additional (real) frequency of oscillation the stability analysis remains unaffected by the inclusion of an equilibrium radial electric field. Hereafter, we shall use the subscripts 0 and 1 to indicate equilibrium and perturbed quantities respectively. As indicated in Sections II and III, all the non-MHD effects that we consider here, as well as the ion inertia, play an important role only within a "singular" layer around r ~ r0, while outside this layer the plasma dynamics is determined by the ideal M H D equations with the ion inertia being neglected. In this section we shall derive the relevant equations which are valid within this singular layer.

A. Singular Layer Equations We first recognize that in the ion stress tensor rri there are contributions from the finite ion gyroradius [16] as well as contributions from the ion-ion collisions. It can be shown [17] that in the linearized momentum balance equation, Eq. (IV.6), the contribution from the part of the ion stress tensor due to the finite gyroradius cancels exactly the part of the inertial term arising from the equilibrium velocity. The components of the linearized Eq. (IV.6) are then, using Eq. (IV.5),

--ito'pVil r : --(d/dr)[pl -~- ( B o + t*±(d2/dr2) Vilr,

--ito'pVilo

=

• B~)/47r]

+ [iF(r)B~, -- 2BoB~o/r]/4rr

(v.8)

--im[p~ + (Bo " B~)/47r]/r + [iF(r)Bao + BoB~r/r

+ B~(d/dr) ao]/4~r + tz±(d2/d#) Nil 0 ,

(V.9)

Here F(r) ~ k • B o = (mBo/r + k~B~), oo' is the Doppler shifted frequency as introduced before, and/zz is the transverse collisional viscosity coefficient given by /z. = (3/10) n~Tdg2i2-ri, where g2i is the ion gyrofrequency and ri is the ion-ion collision time. In deriving Eqs. (V,8) and (V.9) we have neglected BI~ (valid for low/3) and Vil~ (valid ifk=r ~ m), and used V • V1 ~ 0 recalling that V • V1 ~ OVlr/~r for Bo ~ B~ and k~r ~ rn. We have also made the standard assumption that within the singular layer of small thickness dZfl/dr~>~fa/rZ where f l stands for the perturbed quantities. Eliminating [pl + (B0 • B1)/47r] from Eqs. (V.8) and (V.9), using V • V1 ~ 0, V • B1 = 0, and neglecting the variation of the equilibrium quantities within the singular layer, we obtain for rn = --1 modes,

47rp (to, _ i tx± dZ ) d2 d~ P ~dr ~r-~ Vnr = --F(r) -d~ BIr. We introduce the radial displacement variable ~: given by

Nil r =

(V.IO)

--i(to' -- to,i)~,

MAGNETIC RECONNECTION AND rn = 1 MODES

where co,i ~ - - [ c k z / ( n e J B o ) ] ( d p i o / d r ) -~---[c/(neBoro)](dpio/dr diamagnetic frequency. Then Eq. (V.10) becomes

CO*i) (COt __

4rrip(oJ' i

d2

i P

dr 2

)

d2£

--

461

) at r = to, is the ion

r(r)

d2 Yfi

B~,

(V.11)

•

From the radial component of the linearized form of Eq. OVA), using the linearized forms of Eq. (IV.5) and Eq. (IV.7), we can derive the following equation valid inside the singular layer for kzr ~ 1 : i(co' - - t3,~) Bar = - - i F ( r )

lies ,

d2 - - [.71c/(eBor)] F(r) Tel -- (c2v,,/4=) ~ B1,

(V.12)

where oS.e ~ co,e + [.7lc/(eBor)](dTeo/dr), co.e =-- [c/(neBor)](@~o/dr) being the electron drift wave frequency. We notice that the transverse resistivity % does not play any role. In fact, resistive instabilities are physically due to the appearance o f a longitudinal electric field E, = ~ ~jJ , , which decouples the motion of the lines o f force from that of the plasma and thus facilitate the relase of magnetic free energy. From the linearized electron energy balance equation, assuming incompressibility, we find (V.13)

Tel = --i[(aTeo/dr)/(co'-- co,e)] Vel~,

and the linearized form of the equation of continuity yields n~l

--

--i[(dn~oldr)l(co' - - co,a)] V<~a,,

o~ : e, i.

(v.14)

The quasineutrality condition then gives

Vel,. = -- i(co' -- CO,e) ~,

(V. 15)

using the definition Vil, ~ --i(co' -- co,i)~:. Combining Eqs. (V.12), (V.13), and (V.15) we finally obtain BI~ = iF(r) ~ + i[(c%/,,/47r)/(co' -- dJ,e)] d2/dr z B l r

We now introduce the dimensionless singular layer variable x ~--

.

(r -- ro)/r o

(V.16) and define

A = --ico'rH,

~i ~

--O),eTH

,

--O)~iTH

,

~b ~ iBlr/(r aF/dr),=~o ,

where 7n is the typical Alfven time introduced in Section III and oD,,, co,~ are evaluated at r = r0, r0 being the radial position of the singluar surface. We recall

462

ARA ET AL.

that F(r = r0) = 0, and then expanding Eqs. (V.11) and (V.16) about r = ro we can cast them into the dimensionless forms (A -- iAi)(A~" -- D~"')

=

(V.17)

x~b",

¢ = --xse + [e/(A -- i,~,)1 tb".

(V.18)

The resistivity parameter E has already been introduced in Section III and D ~ (tz±~'n)/(pro~) ~ (pi/ro) ~ vii"rn represents the diffusion coefficient for the transport of transverse ion momentum that is due to ion-ion collisons, vi~ is the collision frequency, and pi is the ion gyroradius. In Eqs. (V.17) and (V.18), and henceforth, primes on the functions denote differentiation with respect to x. The fourth derivative with respect to x has been denoted by the superscript "'. The singular layer equations derived in this section are the required generalizations of Eqs. (111.2) and 0II.3) when the relevant "kinetic" effects such as finite-ion gyroradius, finite-electron (ion)-drift wave frequency, and ion-ion collisions are included. In the following sections we shall study their effects on the stability features of the pure resistive modes discussed in Section lII.

VI. FINITE ION GYRORADIUS AND DRIFr WAVE FREQUENCY EFFECTS

The solution of the coupled singular layer equations has to be matched properly with the ideal M H D solutions valid outside the singular layer according to the prescription given in Section IlL This procedure will give the relevant eigenvalue equation. In this section we consider the eigenvalue proglem in the limiting case where D = 0 but the other "kinetic" effects are retained. We shall treat the cases: (A) Resistive internal kink mode, which is predicted to exist for An ~ 0; and (B) Reconnencting mode, which appears in the regime An < 0. A . Resistive Internal K i n k M o d e [18-21]

The asymptotic orderings by which the singular layer equations describe the internal kink mode have already been mentioned in Section III for the pure resistive case. The required generalizations for the discussion of this section are 4,If ~

8; , l ( a

-

ia,) ~

82 ~

x(a -

iaO,

(VIA)

where 3 measures the characteristic width of the singular layer and has to be determined. Following the analysis of Section II1, we introduce the generating function X defined as x ( x ) ~- x ¢ / -

4, = 2~(2, -

i;~3 /i' + x~o

MAGNETIC RECONNECTION AND

m = 1 MODES

463

where X~o is a constant. Then, Eq. (V.17) with D = 0 is automatically satisfied and Eq. (V. l 8) becomes

X" -- 2X'/X -- (A + B x 2) X + BX~ x~ = 0,

(VI.2)

where A = (2t -- i)t,)/E, B = (2t -- iA,)/[eA(2t -- iAi)]. The condition for matching the singular layer solutions properly with the outside M H D solutions remains the same as in Eq. (III.7), namely, X~o = (2/rr) A~

dx X'/x.

(VI.3)

We define ~ ~ x/3, and set ~, = ,~(a -

~,)/(~ -/~).

(VI.4)

Then, Eq. (VI.2) becomes

dZx d~ 2

2 dX ~ d~

( S¢2 4- Az/2)

X =

--x2x~o,

(VI.5)

where, A -- [Ji(A -- i~,)(~ -- i,~i)]l/z, _--- t/O/3

Equations (VI.5) and (III.9) are exactly similiar and we can immediately write the solution as

X/Xo~ = 1 -- TAa/~ fo: dt(1 -- t)(A~/21~)-514(1 + t)-('43a/4)LSl 4 exp(-- ½t&2), (VI.6) provided Re A a / 2 > 1. The eigenvalue equation, derived from Eq. (VI.3) using Eq. (VI.4), is [~(;~ _ i~,)1:/2

=

~n Agz4 1([-4 8/2 -- 11/4) 8 F([A3/~ + 5]/4)'

(VI.7)

where we have defined ~n = 2tn/e:/z. Notice that for )tdi) = 0, Eq. (VI.7) reduces to Eq. (III.11). The eigenvalue equation has to be solved numerically in order to determine the complete dependence of complex ~ on ~,. It is, however, possible to examine analy.tically some of the interesting limiting cases. 595/Iiz]2-I6

464

ARA ET AL.

(i) In the ideal M H D marginal stability case for which ?tH-----0, X~o = 0, Eq. (VI.7) yields A = 1, that is, t ( I -- il,)(1

-- il,) =

1.

(VI.8)

Integrating by parts once, Eq. (VI.6) becomes

X = --X~( A3/2 + 1)/(A3/2 -- 1) -- 2x~[A3/2/(A 3/2 -- 1)] ddt {(1 -t- t) -(A31~/a)-S/4 exp(-- ½t~2)).

× f : dt (1 -- t~tA3/~_x)/4,

(VI.9)

In the limiting case of Xo~ = 0, A = 1, Eq. (VI.9) yields the nontrivial solution X(:¢) = -- (a/2 ~/z) exp(-- ½~z),

(VI. 10)

where = lim [X~/(Aa/2- 1)] Xoo--~O

is an undetermined constant that can be fixed by the boundary condition on ~:. The corresponding ~, consistent with boundary condition, is ~(~) = ½~®[1 -- erf(k/21/2)].

(VIA 1)

In the higher temperature regimes where 1, > 1, Eq. (VI.8) has a root with positive real part (meaning instability) given by

I ,'~ (ZdZ,) t~-~[1 q- i(T,/T,)(T,/T, -- 1) t73],

(VI.12)

where we have taken It -------(TdT~) I,. The other two roots of Eq. (VI.8), I ___ it~ and I _~ iI~ are stable. We first notice that the growth rate of the resistive internal kink mode is considerably reduced by the electron drift wave frequency, the reduction factor being (T,/T~) I~ 2 [Note that I = 1 in the pure resistive case]. Also, for T, > T~, the mode acquires a frequency of oscillation in our frame of reference with phase velocity in the direction of the electron diamagnetic velocity. This feature is preserved for all values of t , . For example, if t , < 1 the unstable root is given by I ~ 1 -]- iI,(1 -- T,/T,)/3,

(VI.13)

for t, = --(T,/T,) 1~. (ii) In the case where 0 < t n < 1, the eigenvalue equation becomes t ( f -- il,)(t

-- il,) ~

[1 q- 1n/{rr1(t

-- il,)}1/2] 2.

(VI.14)

MAGNETIC RECONNECTION AND

m ----- 1

MODES

465

(iii) In the other limiting case for which E --~ 0 and ,~n, ~, ~,, )~i are considerably larger than unity, Eq. (VI.7) reduces to

a(a - / h i ) _~ a,,~,

(vi.15)

where we have used the asymptotic property of the gamma functions, namely, 1-'(z)/l'(z + a) --+ z -a for z >~ 1. The unstable root of Eq. (VI. 15) is ,~ _~ ;~.[1 -

;~7/(8;~,,D] +

i.~,12,

(VI.16)

for An > [ hi [. Once again, we observe that the ideal M H D growth rate (h ~ AH) is reduced by the finite drift wave frequency effect. However, the mode now has a phase velocity in the direction of the ion diamagnetic velocity unlike the case An = 0. Therefore, there has to be a value of ~n for which Im h = 0 and the characteristic phase velocity changes its direction. If we do the analysis for the case '~H > 1 more carefully keeping more terms in the asymptotic expansion for the gamma function we find that indeed Im h vanishes when

(w.17)

~/t '~ 2.4 1(2)~, -- ~)/~, ]~/a.

We solved Eq. (VI.7) numerically and the behavior of Im )~ as a function of ,~n is shown in Fig. 6 for the particular choice of,~ ---- --(TdTe ) ~ . From the foregoing analysis we can conclude that the growth rate of the resistive internal kink mode is reduced by the drift wave frequency for all values of An >/0. The mode can be suppressed most effectively in relatively high temperature regimes and for plasma parameters such that 0 ~< ;~n "~ l, although it cannot be ~ompletely stabilized.

o

I

FIG. 6. Imaginary part of ~ for the resistiveinternal kink mode as a function of ~,vcorresponding to ~, = 1.0 and various values of TdT~.

466

A R A ET AL.

Another effect of the finite drift wave frequency is to introduce fine structure (spatial oscillations) in the radial profile of this mode (see Fig, 7 and compare it with Fig. 3). Using Eqs. (VI.4) and (VI.8) we find that, for An = 0, 3 is given by t~2

=

e/(,~ - - i/~e),

and the generating function X h a s oscillations within a gaussian envelope [see Eq. (VI.10)]. The width of the envelope, 6R, is given by 3R2/r02 m E/(Re h),

(VI.18)

and the typical "wavelength" of the oscillations, 3t, is given by 312/ro~ ~ ~ / I m ( h - i?Q[.

I

t,I

(VI.19)

i

I

ro a

.t

I I I I I I _

......

0

r

FIG. 7. Radial profile of the resistive internal kink mode when IH = 0, ~, ¢ 0. The vertical dotted lines indicate theapproximate boundaries of the singular (inner) layer. It should be noted that when ;t~ --~ 0, Im A ~ 0 and, hence, 3t -+ o% i.e., the oscillations disappear. In the limiting case ~ >~ 1, we find

rO:;V(TdT.), implying that

p,(v/vAo)(ro/r.),

(vi.20)

and, correspondingly, ~,Ip~ ~

( V j VAo)(Te/Ti)a/2 (ro/r,) ~-3/~,

(VI.21)

where pi is the ion gyroradius, V8 ---- (Tdm31/2 is the ion-sound speed, VAo is the Alfv6n speed, and r,~ = I dln n/dr [-1. Since ~, ,~ Ta,/2n1/3, the right-hand side of Eq. (VI.21) goes like (TdT~) x/2 TJ/4. Thus in the high electron temperature regimes,

MAGNETIC RECONNECTION AND

m

=

1 MODES

467

8s can be less than p, and the considered fluid description will break down. At this point, however, the ion-ion collisional viscosity effect represented by D in Eq. (V. 17) can no longer be ignored and comes to the rescue.

B. Reconnecting Mode [7] The features of the reconnecting mode change substantially as we consider the high temperature limit where the classical resistivity becomes so small that the resistive growth rate given in Eq. (III.30) becomes comparable with the electron drift wave frequency. In this limit, Eqs. (III.19) and (III.20) are modified as: A(A -- iA,) £"1 = x¢;',

(VI.22)

¢o = - x £ : 1 + [ 4 0 - iz,)] ¢2.

(vi.23)

The equation for £-1 becomes A(A -- iA3 £:1 -- [(A -- iA~)/e] x2£_1 = [(Z -- iA,)/E] CoX,

(VI.24)

and the matching condition is now given by [(A -- iA,)/e] f_+~ [1 + x£-1/¢0] dx = --Tr/AH.

(VI.25)

From Eqs. (VI.22) and (VI.23) we observe that all terms are of the same order if, A(A -- it,) ~ ~/I An l, E/(A -- iA,) ,-~ 3 [AH 1. The solution of Eq. (VI.24) can be given in the form of an integral representation as

£-1(x)/¢o

=

- ½<~x

s:

dt

(1 -- t2) -~14 exp(-- ½etx2),

(VI.26)

where a2(A) = (a -- iA~)/[eA(A -- iai)].

(VI.27)

After performing the integration in Eq. (V1.25), we obtain the eigenvalue equation, for Ass < 0, [A(A --

iA3] TM (A

iA,P/a = ~/'1 Ass I -~ A05/;

(VI.28)

provided Re cr > 0, where A 0 has been defined after Eq. (III.30). We recognize by inspection of Eq. (VI.26) that Re a > 0 corresponds to well-localized normal modes. It is convenient to introduce the following normalized variables:

A =- At AH 1415IAo,

A,=_ A, I M ?/5IAo, A,~A,

I A~ l'/5/Ao,

468

~ A ET AL.

Then Eq. (VI.28) becomes

A(A - ~Ai)(A - ~A~)~ = ~.

(VI.29)

At low temperatures such that Ae(i) can be ignored, we recover the pure resistive growth rate, namely, z{ = e3/5, and the mode decays away from the singular surface. As the temperature increases (low collisionality regimes) such that A, >~ E3/5, the reconnecting mode acquires a frequency o f oscillation in the moving frame of reference, which is in the sense of, and comparable with, the electron drift wave frequency, and the growth rate is reduced. As the electron temperature is increased, Re cr decreases and passes through zero for Re./I > 0. This occurs, for example, when A~ ~ 1.57E3/5, Re./I ___ .55~3/5, and Im ,/I ~ ¢ a / 5 ~ .64A~. These estimates can be obtained from Eqs. (VI.27) and (VI.29) if we assume A, = - A d Z When A, > A, ~ 1.57¢a/5, Im ,/I ~ A~ (i.e., Re oY -+ 6,~) and Re a becomes negative for Re A > 0, i.e., the unstable mode ceases to be spatially localized. The stability of systems in such a situation where no localized normal modes exist has been considered in Refs. [22] and [23]. Following their analysis it can be shown that the reconnecting mode remains unstable with reduced growth rate even at high temperatures where the mode acquires an important electrostatic component and the relevant unstable solutions belong to a continuum of nonlocalized modes which are of "fluid-drift" type with Re o J ' ~ c5,~. By superimposition of these nonlocalized solutions it is possible to construct a convective wave packet which propagates outwards from the unstable plasma region undergoing spatial amplifications before entering a new plasma region where it becomes damped. In Fig. 8 we have presented the radial profiles of ~:-1 for the particular case when Re cr = 0, i.e., A~ = Ac • The characteristic width of the spatial oscillations in ~:-z, which are due to the finite drift wave frequency, is given by 8~osc/ro ~ ~

I(Im ~)1-1,

(VI.30)

and it decreases with temperature. Numerical estimates indicate that the temperature for which Re cr = 0 is already below the values that have been achieved in existing experiments such as Alcator and that the corresponding values of 8ose is somewhat larger than the ion gyroradius. As higher temperature regimes are approached, 3ose becomes of the order of or smaller than the ion gyroradius. At this point the effects due to ion-ion collisional viscosity have to be included. We recall at this point that the asymptotic connection of the inner solution with the outer one is made for

and the mode of interest does not become electrostatic to the extent that ¢ ,-' 4(~, - ~&) ¢".

MAGNETIC RECONNECTION AND m ---- 1 MODES

469

Therefore, in order that the mode does not become electrostatic we require that IS

-

-

iA, I >~ ,/[ A, I.

(VI.31)

We refer to the eigenvalue equation, Eq. (VI.28), and observe that in the high temperature regime

ReY

(a)

ImY

AAA

(b)

AAA

VVV

FIG. 8. Radial profiles o f the reconnecting m o d e at a high t e m p e r a t u r e t h a t c o r r e s p o n d s to R e ~ = 0, i.e., A, = Ao : (a) R e Y[------(SAb0)~_d vs ~(------x/8); (b) I m Y v s ~. W e have scaled 8 to t h e s a m e distance as in the pure resistive case for illustration purpose.

and so the requirement given in Eq. (VI.31) can be satisfied if A~z ~ 1/I A~ [.

(VI.33)

It is easy to verify that this is true for a realistic choice of the values of the relevant physical parameters.

VII. EFFECTS OF ION-ION COLLISIONS We have shown in Section VI that in the high electron temperature regimes both the resistive internal kink mode and the reconnecting mode acquire spatially modulated radial profiles which are due to the finite electron drift wave frequency. The typical "wavelength" of these spatial oscillations tends to become of the order of, or smaller

470

ARA ET AL.

than, the ion gyroradius and, as indicated earlier, the ion-ion collisional viscosity effect tends to become important. From Eq. (V. 17) we observe that this occurs when (ro=/a,2) D ~> Re A.

(VII.l)

For the resistive internal kink modes in the high temperature regimes, if we use for Re A and 81 the expressions given in Eqs. (VI. 12) and (VI.21) respectively, Eq. (VII. 1) implies that El/a/,~e5 ~ vii(ro/VAo)(VAo/Vs)2(Ti/Te)2(rn/ro) 2.

(VIi.2)

It is possible to verify that this inequality is satisfied for values of the relevant parameters that are close to those of existing experiments such as Alcator. Similar estimates can be obtained for the reconnecting mode. In particular, it is found that, for the reconnecting mode," the colhSlonal .... " Viscosity effect tends to become important at higher temperatures because the growth rate decreases more slowly as the temperature increases and the "wavelength" of the oscillations remains larger than the ion gyroradius over a longer range of temperatures. An exact analytical treatment of the eigenvalue problem including the effects of ion-ion collisions is a difficult task because of the complicated nature of the singular layer equations. Instead, we shall present here approximate analytical methods of studying the stability feature s in the high density and high electron temperature regimes. The eigenvalue problem has been solved numerically by Ara et aL in Ref. [24]. The predictions of the approximate analysis presented here agree well with the numerical results both qualitatively and quantitatively. We write the singular layer equations in the form a(~" - - D , ~ " )

(VII.3) (VII.4)

= x~b",

= --x~ + b¢",

where, a -~ A()t - - iAi) ,

b -~ 4 0 D,

-

iA~),

(vii.s)

~ D/A,

and consider the special case A~ : 0. In order to make use of a variational method we Fourier analyze Eqs. (VII.3) and (VII.4), and then eliminating ~b(k), obtain dZO dk 2

a(1 + D,k2)(1 4- b k 2 ) (I + bk2) 2

a --

where,

~(k)

k =

(1 + bk~)~/~ ~ t ~ '

3b

• (~),

(vii.6)

MAGNETIC RECONNECTION AND rn = 1 MODES

471

~:(k), ~b(k) being the Fourier transforms of ~¢(x) and ~b(x) respectively. In terms of the new variable y ~-- b~/~k, Eq. (VII.6) becomes d2qb dy 2

~(1 + D * * y2)(1 + yZ)a _ 3 tp(y), (1 + y~)2

(VII .7)

where D** ~ D , / b , o~ ~ a/b. F r o m Eq. (VII.7) we can construct the variational form f+-2 dy [qb2/(1 ~- yZ),, _ (dg)/dy)2]

~x =

(VII.8)

f-+o~ dy (1 -/D**y2)(1 - / y 2 ) ~2 • For D** = 0, and AH = 0, c~ = 1 and qb(y) is given by ~o(Y) = [const/(1 + y2)~/2] e x p ( l y 2 ) . Using (Po as a trial function and evaluating the relevant integrals we obtain c~ ~ (1 -/ ½D**) -x, that is, A(A -- iA~)(A -- iAi) ~

1 -- 1/3(~

_

i/~)2(~

_

i~),

(VlI.9)

where we have defined ~ ~ A/el/3, ~ ( i ) ~ A,(i)/el/a and /3 ~ D/E. Consider the quantity

/3=_ /~±~'n prJ

3 (~., mi~ 1/2 - - Y6

~n~ /

~ '

where fir ~ (8rrnT,)/Bo 2. For typical tokamak parameters, T , / T i = 2 and /3~ = 5 × 10 -3, and, hence, /) ~ 1/10. Since /3 is very small, Eq. (VII.9) with the equality sign is a good approximation to the required eigenvalue equation. In order to find the unstable root of the eigenvalue equation in the limiting case of~ > 1,~,and/3 < 1, w e s e t

= ~o + ~ , where ] ~o/,~ ] ~ to lowest order,

1 ,~o/A~ I ~

I ~/'~o I, and we assume Ai = --/~e for simplicity. Then, Re Ao = 1/,~( 1 q-d), Im ~o = d~d(1 + d),

where d ~ / ) / 2 . To next order, we find Re ,~ ---- (Re ,~o)[1 -- 2d2/(1 q- d)2],

(VII.10)

472

ARA ET AL.

and Im ~ = Im ~0 + [d/~,(1 q- d)][(Re ~0) ~ - ( I m A0)~].

(vn.ll)

We, therefore, can conclude that, in the present devices where/5 ~ 1, viscous stabilization does not occur although the growth rate is reduced. However, the present analysis indicates the possible importance of collisional viscosity as a stabilizing mechanism in future experiments with physical parameters that correspond to /3 > 1. For a hydrogen plasma with TJT~ = 2, /3 > 1 corresponds to/3, > 6 %. For such large values o f / 3 or/3~, however, the present analysis is inadequate. The additional effect of collisional viscosity is to induce a small but finite frequency of oscillation even for E, = 0, which corresponds to to' = to. Furthermore, according to the numerical analysis reported in Ref. [24], the most significant effect of the collisional viscosity is to broaden or nearly eliminate the spatial oscillations in the radial profile of ~:, as a result of which the typical "wavelength" of the oscillations remains effectively larger than the ion gyroradius. This reassures the validity of the fluid description adopted here. Similar conclusions are also derived [25] for the reconnecting mode which, as indicated earlier, acquires a frequency of oscillation to' ~ dJ,, in the high density and high electron temperature regimes of furture interest. To summarize, while the growth rates of the m = 1 modes are sharply reduced by the "kinetic" effects considered in this paper, the frequency of oscillation of these modes can be expressed as to = t 3 . J ( ~ H , f)) -q- cE~/(roB~)

(VII.12)

where F is a function of ~H a n d / ) . In particular, F ~ 1 for '~H ~ -- 1 (reconnecting mode), and F ~ ½/) [see Eq. (VII.11)] for ~hr = 0 (resistive internal kink mode) a n d / ) ~ 1. VIII. CONCLUSIONS A number of observations can be made in the conclusion of the present paper: (a) In a toroidal configuration the ideal M H D internal kink mode can be made stable by the effects of finite toroidicity and the mode that can be excited in its place is associated with the finite electrical plasma resistivity. In the limit of relatively small values of resistivity, this mode that we call "reconnecting" has a radial profile that is peaked around the surface r = r 0 where q(r) = 1 and is consistent with the available experimental observations [3]. (b) In the limit of relatively high temperatures, the effects of finite ion gyroradius, electron drift wave frequency, and ion-ion collisions modify significantly the radial mode profile and reduce drastically its growth rate. This leads us to believe that in the higher temperature regimes to be achieved by future magnetic confinement experiments plasmas will be more immune from the effects of these modes.

MAGNETIC RECONNECTION AND rn ~

1 MODES

473

(c) Presently, m = 1 oscillations are observed to arise in the last phase of the sawtooth oscillations that are excited in toroidal plasmas (see Fig. 1) as well as at the onset of the disruptive instabilities [26]. (d) Following Ref. [3], we consider the rise phase of the sawtooth, for r ~< r0, to result from a thermal instability that tends to produce an increasingly peaked temperature and current density profile. Given the relatively small growth rates predicted by the analysis given in Sections VI and VII, that applies to existing experiments, a relatively steep current density profile is needed to trigger the relevant instability. (e) The onset of the rn ---- 1 mode is identified by that of a precursor oscillation preceeding the sharp drop in the sawtooth. The observed phase velocity of this precursor oscillation is in the direction of the electron diamagnetic velocity. According to the analysis developed in Sections VI and VII, the observed frequency of oscillation can be explained by the one that is found in a frame of reference where there is no equilibrium radial electric field and is related to the electron drift wave frequency for ,~/ < 0, /) :/: 0. Another possibility is to assume that the ions are electrostatically confined by a radial electric field so that the relevant mode would always have (even for An ---- 0 and b ---- 0) a finite frequency due to the resulting Doppler shift. (f) The abrupt drop in the central plasma temperature, corresponding to that of the sawtooth, can be interpreted either as a nonlinear effect of the m = 1 fluctuation after this has reached a critical amplitude level or as a result of a change in the sign of,~n from negative to positive. In the latter case, we would have first a weakly unstable reconnecting mode and then an almost purely growing mode with a radial profile similar to that of the ideal M H D internal kink mode. (g) There is experimental evidence that the amplitude o f the sawtooth oscillations tends to increase as the density is increased while their length tends to increase. Furthermore, in the Alcator experiments, the temperature profile is observed to become narrower as the density increases. Thus we may attribute the lengthening of sawtooth oscillations to the decrease of the growth rate with density for the relevant thermal instability [27, 28]. We can describe this with a relatively simple model equation for the electron thermal energy balance at the center of the plasma column, 3 eT, E"~ [T*] 3/2 vznT,. 2 n ~ - ~70 \To! --

(VIII.l)

Here 70 = ~cz(T, = To), To is a constant reference temperature, uL represents the rate of thermal energy loss from the center of the plasma column and E~ is the applied longitudinal electric field. It is easy to verify from Eq. (VIII. 1) that the growth rate o f the resulting thermal instability and the value o f the temperature at which the instability can become explosive decrease as n increases for reasonable models o f vL.

474

.ARA ET AL.

(h) Notice that an explicit expression for the growth rate of the resistive internal kink mode when An = 0 and T~ ~ T i , obtained from Eq. (VI.12), is __ Y ~

(r.]

2

1

~s z

T~I2 '

where ~,~ = ~qjhcZ/4rr is the magnetic diffusion coefficient, ~n ~-cTJ(eBo) is the so-called Bohm diffusion coefficient, r, = - - ( d l n n / d r ) -x and V~o is the Alfv6n velocity derived with the poloidal magnetic field as indicated in Eq. (11.7). In the case of the "reconnecting" mode the corresponding expression for ~ is

t[ VAo][ r.

][ R

]at /a

1

This estimate can be obtained if we recall from Section VI that in the high temperature regimes the growth rate of the "reconnecting" mode is given by Re A ~ Eli )tH I ~/~ )t~ i~. We see that in either case ), tends to decrease sharply when the electron temperature increases. Thus we can expect that magnetically confined plasmas with higher electron temperatures than those realized so far will be considerably less affected by the onset of these modes than those produced in existing experimental devices.

APPENDIX

A: SOLUTION OF THE EQUATION FOR X

Consider Eq. (Ill.9). The solution of this equation by the method of integral representation has been given in the Appendix of Ref. [21]. Here we present an alternative method of arriving at the same integral representation. In terms of the variable z ~ ~2 Eq. (111.9) becomes d2x dz 2

1

dx

1

2z

dz

4z

(z -P- ja/~) X =

1 4

X~.

(A.1)

Substituting

X/X~o =

a~ n

tZ) exp(-- ~z),

where L," are the generalized Laguerre Polynomials satisfying [29] d2L~ ~, dL~ ~ z ~ + (~ -t- 1 - - z) ~ + nL~ ~ = 0 ,

(A.2)

yields a.z-ae-Zl~L;Zl2[n

-t- (~31~ _

1)/4] =

¼.

(A.3)

MAGNETIC RECONNECTION AND m ~--- 1 MODES

475

Next, multiply by e-Z/2z-1/2L~ 3/a and use the orthogonality relation [29]

fo~ e-~z~Lm~Ln ~ dz = ~mn (n +n! ~)!

(A.4)

to find

z-a/2e-~12L~ a/2 dz.

an (n --n!3/2)! (~3/2 -- 1 + 4n) =

(A.5)

Employing L~~ = L~+1 -- L~+_~and the relation [29]

fo°~ e-y~t~L, ~ dt = (n +n! a)! (~ -- 1)~ ),-(~+,+1)

(A.6)

to evaluate the integral in Eq. (A.5) gives

X/Xoo : 2 -1/~e-z/2 ~, (--1)" ( 4 n -

1)L~a/2(z)/(~3/2 -- 1 -k 4n).

(A.7)

n=0

Next, writing (~3/2 _

1 -k 4n)-1 = fo~° dfl e -Bta~/~-l+4~),

and using the relation [29] ) ' Ln°~(z)x n = ( l - - x ) -cz-I

exp /{ \ ~x z, ;J~

n=0

(A.8) J L /

results in an integral representation for X [introducing y -----exp(--4/3)]: 1 z 1--y X/Xo~ = 1 -- 2-5/2,~a/2 fo dy y(~/~-5)/4 (1 -k y)1/2 exp [-- ~ ( ~ ) ] .

(A.9)

as may be verified by direct substitution into Eq. (III.9). The substitution t = (1 -- y)/(1 -k y) transforms Eq. (A.9) into Eq. (III.10).

ACKNOWLEDGMENTS This work was sponsored in part by the U.S. Energy Research and Development Administration, and by the Energy Research and Development Administration under contract with the U n i o n Carbide Corporation.

REFERENCES 1. B. CovvI AND A. FRIEDLAND, Astrophys. J. 169 (1971), 379. 2. S. V. MIRNOV AND I. V. SEMINOV,in "Plasma Physics and Controlled Nuclear Fusion Research," Vol. II, International Atomic Energy Agency, Vienna, 1971, p. 401.

476 3. 4. 5. 6. 7.

8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19. 20.

21. 22. 23. 24, 25. 26. 27. 28, 29. 30,

ARA ET AL. S. yon GOELER,W. STODIEK,AND N. SAUTHOFF,Phys. Rev. Lett. 33 (1974), 1201. V. D. SHAFRANOV,Zh. Tekh. Fiz. 40 (1970), 421 [Soy. Phys. Tech. Phys. 15 (1970), 175]. M. N. BUSSAC,R. PELLAT,D. EDERY, AND J. L. SOULE,Phys. Rev. Lett. 35 (1975), 1638. B. COPPI, R. GALVAO, R. PELLAT, M. N. ROSENBLUTH,AND P. H. RUTHERFORD,M.I.T. Report PRR-7520, Cambridge, Mass., 1975; Fizika Plazmy 6 (1976), 961. Previous accounts of the reconnecting mode have been given in: B. BASUAND B. COVvI, Bull. Amer. Phys. Soc. 21 (1976), 1090; B. BASUAND B. COVPI, in "Plasma Physics and Controlled Nuclear Fusion Research," Paper CN-35/B13, International Atomic Energy Agency, Vienna, 1976; B. BAsra AND B. CovvI, M. I. T. Report PRR-76/39, Cambridge, Mass., 1976. H. P. FURTH,J. KILLEEN, AND M. N. ROSENBLUTH,Phys Fluids 6 (1963), 459. B. Covet, J. M. GREENE,ANO J. L. JOHNSON,NucL Fusion 6 (1966), 101. B. B. KADOMTSEV,Fizika Plazmy I (1975), 710 [Soy. J. Plasma Phys. 1 (1975), 389]. A. F. DANILOV,YU. N. DNESTROVSKn,D. P. KOSTOMAROV,AND A. M. POPOV, Fizika Plazmy 2 (1976), 167. B. V. WADDELL,M. N. ROSENBLUTH,D. A. MONTICELLO,AND R. B. WHITe, NucL Fusion 16 (1976), 528. A. SYKESAND J. A. WESsoN, Phys. Rev. Lett. 37 (1976), 140. S. I. BRACINSKII,in "Reviews of Plasma Physics," Vol. 1 (M. A. Leontovich, Ed.), Consultants Bureau, New York, 1965, p. 214. B. Cocci AND J. REM, Phys. Fluids 17 (1974), 184. K. V. ROBERTSAND J. B. TAYLOR,Phys. Rev. Lett. 8 (1962), 197. B. CoPPt, Phys. Rev. Lett. 12 (1964), 417. B. BASUAND B. COVVI,"Proceedings of Annual Meeting on Theoretical Aspects of Controlled Thermonuclear Research," Paper 1A-l, University of Wisconsin, Madison, Wisconsin, 1976. G. LAVAL,M. N. ROS~NBLUTnANO B. V. WADOELL,Ref. [18], Paper 1A-11. M. N. BUSSAC,D. EOERY, R. PELLAr, AND J. L. SOtrLE, in "Plasma Physics and Controlled Nuclear Fusion Research," Paper CN-35/B3, International Atomic Energy Agency, Vienna, 1976. B. BASUAND B. CovPI, M.I.T. Report PRR-76/38, Cambridge', Mass., 1976. B. Covet, G. LAVAL,R. PELLAT,AND M. N. ROSE~,mLUTrI,NucL Fusion 6 (1966), 261. P. H. RUTHERFORDAND H. P. FURTH,Princeton University Plasma Physics Laboratory, Report MATT-872, Princeton, N. J., 1971. G. ARA AND B. CoPPI, Bull. Amer. Phys. Soc. 21 (1976), 1090. G. ARA, B. BASU, AND B. CoPPI, M.I.T. Report PRR-77/7, Cambridge, Mass., 1977. S. V. MmNov ANDI. V. SEMINOV,in "Plasma Physics and Controlled Nuclear Fusion Research," Paper CN-35/A9, International Atomic Energy Agency, Vienna, 1976. B. CoPPI, in "Lecture at the Varenna Symposium on Plasma Heating," Sept., 1976. O. MANLEY, Private Communication at the Varenna Symposium on Plasma Heating, Sept., 1976. W. MAGNUS, F. OBERHETTINGER,AND R. P. SONI, "Formulas and Theorems for the Special Functions of Mathematical Physics," Springer-Verlag, New York, 1966, p. 239. T. F. R. GROUP, NucL Fusion 17 (1977), 1283.