Free-boundary tokamak equilibria

Physica 11X (1982) 103-131 North-Holland Publishing Company

FREE-BOUNDARY

103

TOKAMAK

EQUILIBRIA

J.M. AKKERMANS Association Euratom-FOM, The Netherlands

FOM-Instituut

voor Plasmafysica

“Rijnhuizen’:

P.O. Box 7, 3430 AA Nieuwegein,

Received 4 May 1982 Revised 8 July 1982

Free-boundary equilibria of high-p tokamak plasmas with surface currents are investigated both analytically and numerically in the limit of a small inverse aspect ratio. The method described to solve this non-linear problem exploits fast Fourier transforms and conformal mapping techniques, and uses an iterative procedure on the plasma boundary. The plasma may be surrounded either by a vacuum region or by a layer of force-free currents. The influence of the cross-sectional shaping of the conducting shell is studied. Numerical and analytical results are given for circular, elliptic, D-shaped, and racetrack-like wall shapes. Force-free currents enhance the coupling between the shape of the wall and that of the plasma. It is concluded that for given outward displacement and horizontal width of the plasma, both force-free currents and wall shaping, in particular the triangularity, have a favourable effect on the poloidal beta. The problem of the existence and the nature of the equilibrium limit for the poloidal beta is solved.

1. Introduction In this paper we discuss free-boundary equilibria of high-beta tokamaks. In the free-boundary analysis of equilibria one asks for the cross-sectional shape of a plasma freely moving inside a given conducting shell. The paper investigates the general case of toroidal plasmas with surface currents that are surrounded by a region carrying force-free currents and are bounded by a conducting wall of arbitrary shape. These two characteristics, wall shaping and a force-free instead of a vacuum region between the plasma and the wall, are expected to be favourable for obtaining a high value of beta, the parameter that measures the efficiency of magnetic plasma confinement. The interest in studying toroidal plasmas and optimizing their characteristics with respect to magnetohydrodynamic equilibrium and stability stems from the success of tokamak research in the seventies. Two different approaches are encountered. In the fixed-boundary approach the plasma cross section is prescribed. This approach is not self-consistent because the wall may intersect flux surfaces. From the mathematical point of view it leads to a linear but ill-posed problem. An approach closer to the experimental situation is to self-consistently determine the plasma shape resulting from the interaction with a conducting wall. This more complicated free-boundary problem is well-posed, but intrinsically non-linear. Within the fixed-boundary approach, Gajewski [l] and Kerner et al. [2] studied optimization of the plasma cross section with respect to the equilibrium beta in a straight geometry. Toroidal studies employing a skin-current model and the high-beta tokamak ordering were carried out by Freidberg et al. [3,4] and D’Ippolito et al. [5]. In ref. [S] the effects of surrounding the plasma by force-free currents were also considered. The analysis of force-free currents is of interest for the investigation of screw-pinches, experimentally studied at Rijnhuizen [6], since such currents are commonly produced during pinch discharges. The application of the fixed-boundary model of a high-p tokamak with force-free currents to the screw-pinch has been considered by Rem [7]. A variety of methods has been developed for solving the free-boundary problem. Methods similar to the hodograph method borrowed from fluid dynamics are applied to different types of free-boundary 0378-4363/82/0000-0000/$02.75

@ 1982 North-Holland

104

J.M. Akkermans

/ Free-boundary tokamak equilibria

problems in refs. [8-lo]. Grad et al. [8] study the case where the plasma is confined by a vertical field. Shercliff [9] investigates the deformation of the plasma due to parallel line currents, whereas Vabishchevich and Degtyarev [lo] study plasmas surrounded by a conducting shell. Stevens [ll] proposes an iteration procedure on the boundary of a plasma that is confined by applying a vertical field. Suzuki [12] and Qing and Zhou [13] develop iterative procedures for the solution of the Grad-Shafranov equation while Steffen [14] uses an integral equation method. These studies are restricted to the low-beta regime. Within a skin-current model, Goedbloed [15] has investigated high-beta tokamak plasmas confined within a conducting shell with circular cross-section. A very general, but very time-consuming numerical method for computing three-dimensional free-boundary equilibria has been developed by Bauer et al. [16]. A review of numerical codes for the computation of ideal magnetohydrodynamic equilibria has been given by Lackner [17]. None of these studies considers force-free currents. This paper extends the approach of ref. [15] to the general case of a force-free outer region and arbitrary wall shapes. It is the aim to obtain a method which is both fast and accurate. A fast procedure is needed to permit scans of the parameter space without undue effort. A high accuracy is of importance because equilibrium calculations are an essential part of the stability calculations. The overall accuracy of the stability studies is limited by the accuracy of the equilibrium calculation. A specific feature of the method adopted is the extensive use of conformal mapping methods and of Fourier transforms. Another specific feature of the present approach is that it allows an analytical study of the properties of the free-boundary plasma equilibrium. In section 2 the equations for the model are derived. In section 3 analytical approximations for the equilibrium plasma shape and the poloidal beta are presented and compared with numerical results for the case of a circular wall and a vacuum region, which serves as the reference case. Non-circular walls and a vacuum outer region are the subject of section 4, whereas force-free regions bounded by a circular wall are discussed in section 5. Section 6 considers the combined effects of wall shaping and force-free currents outside the plasma. Finally, in section 7 some remarks are made on the maximum poloidal-beta limit for the free-boundary equilibrium.

2. The model We consider axisymmetric toroidal systems. The core consists of a dense plasma of constant pressure, confined by currents flowing on the plasma surface only. The dense plasma is surrounded either by a pressureless plasma region carrying force-free currents or simply by a vacuum region, which in turn is bounded by a perfectly conducting shell with an up-down symmetric poloidal cross section. A schematic of the geometry is given in fig. 1. Next, we introduce dimensionless coordinates in the poloidal plane by normalizing the horizontal wall radius to unity (see fig. 2a). The poloidal cross-section of the wall is labelled C,, and that of the plasma Co. Furthermore, we have two important geometric plasma parameters, viz. the shift of the plasma centre A and the plasma half-width a. In the following analysis we assume these two parameters to be given, and we will show that fixing these parameters suffices to determine the problem (apart from some scaling constants). The free-boundary problem can now be formulated as follows: given the wall shape C, and the plasma parameters A and a, determine the plasma shape Co and other physical variables of interest such that pressure balance is established. 2.1. Maxwell According toroidal:

equations to this approach,

in the plasma

region

we have a vacuum

magnetic

field that is purely

J.M. Akkermans

/ Free-boundary

tokamak

105

equilibria

R

Fig. 1. (R, z, 4) coordinates used to describe the geometry of the model.

B = (RdR)&e4, In the force-free

region

Bo=B(R the magnetic

= R,).

(1)

field also has a poloidal

component.

This region

is described

by

div h = 0 ,

(24

curlti=

(2b)

&,

where hats denote quantities in the outer region. For simplicity, we assume (Y to be constant. problem is further simplified by making the assumption of a small inverse aspect ratio:

and the additional

assumption

on the force-free

The

constant

ab-E.

WI

This latter condition precludes the possibility of toroidal field reversal in the force-free region. Furthermore, it ensures that the poloidal variation of the toroidal field component is approximately the same as in eq. (1): &,&‘, = B+/B, = [lWith regard V&=0,

I%,-&(R

&xl-l,

to the poloidal

field component,

= R,). it follows

(4) that to leading

order

in E

V*.&,-&~~.-~_

Here, V is the two-dimensional operator, V* = (8/6’y, -~Sl&‘x), and

(5) nabla operator, V = (a/&x, a/ay), V* is the two-dimensional b, is a dimensionless two-dimensional poloidal field vector,

curl nor-

106

J.M. Akkermans

malized

/ Free-boundary

tokamak

equilibria

at Co:

lip= lipl
(6)

brackets

indicate

averaging dl = If(S)

@,)=#&,dl/(jdl,

the curve C, being given in polar constant cy by a constant r, defined

q* = a&o

f/

dl

I+ = ae&/(BJ

over the contour

+ I]“’

C,:

d6,

coordinates by r = f(G). according to

1

(7) For convenience,

we have

replaced

.

the

@b)

As the parameter q* measures the total toroidal current in the plasma, the force-free constant r is proportional to the ratio of the force-free and plasma current densities. This is discussed in more detail in section 5. Deriving the poloidal field from a flux function + : 6, = -V*$, we find from eq. (5) that the problem reduces to solving a two-dimensional Poisson equation in the force-free region: A+==,

(94 on C, ,

Pb)

& = I& on Co .

(9c)

I/J = & = constant $=O,

The index v indicates the normal derivative. The inhomogeneous part of the Poisson equation represents the strength of the force-free currents; the vacuum case is simply given by r’ = 0. From the Maxwell equations we have derived a problem, yielding a relationship between plasma shape and magnetic field around the plasma, that is very similar to problems in two-dimensional potential theory. This suggests to treat the present problem by one of the standard tools used in this theory, viz. conformal mapping. 2.2. Conformal

mapping

Let us first split the poloidal problem I(, = c#l+ PI/P For the particular

so that Al/lP= 1.

)

flux function

A4=0.

solution,

we choose

into a harmonic

part 4 and a particular

solution

I,V’of the

(104

J.M. Akkermans

Next, we apply Green’s second Specializing to the plasma boundary,

Mro) = - f [Ch,

WW

/ Free-boundary

theorem we find

tokamak

107

equilibria

to 4 in the force-free

region

bounded

by CO and

C,.

(11)

- W-0, GM-@] dlh+ f G,(ro, O$(ri) dli , Cl

G

where the indices 0 and 1 refer to the curves CO and C,, respectively. G represents Green’s function, AG = S(r - r’), G = 0 on Ci. However, the Green function is not known for a general boundary C,. This is where the conformal mapping techniques enter. We map the physical plane (the z-plane) onto a computational l-plane, such that the wall boundary C, is mapped onto a circle, and the plasma boundary is centred with respect to the origin at the horizontal axis. The existence of such a mapping of the region bounded by C, is guaranteed by Riemann’s mapping theorem. This mapping, represented by

(12)

z=&& m=O can be thought to be effected by two successive mappings. Ci into a circle while leaving the origin unchanged: 2 = 2

&P

The first one, which is numerical,

transforms

.

(134

m=l

Next, an explicit 2’ = (l+

Moebius

transformation

is used to centre

the plasma

S)/(l + sg) .

(13b)

The purpose of the conformal mapping (13a) is to exploit the known properties of Green’s function for the unit circle. The centring of the plasma by means of the Moebius transformation (13b) has a different aim, viz. to employ Fourier analysis, even in the case of high-beta plasmas that are highly shifted outward. A scheme is given in fig. 2. The known parameter pair (A, a) is replaced by a pair (8,~) that can be computed from A and a by the mappings (13a) and (13b). In the numerical computations S and y are

Fig. 2. Conformal putational L-plane.

mapping. (a) Physical plane. Dimensionless coordinates, z = x + iy; (b) intermediate Because of the up-down symmetry only the upper part is shown.

Y-plane;

(c) com-

108

J.M. Akkermans

I Free-boundary

tokamak

equilibria

then treated as the basic parameters. In general, the coefficients +,,, of the mapping (13a) have to be computed numerically. A method to do this based on fast Fourier transforms has been developed by Henrici [18], and has been extensively discussed by Goedbloed [19] in connection with problems of magnetohydrodynamics. The coefficients pu, of eq. (12) are obtained from &, 6 and y. Since Laplace’s equation is invariant under conformal mapping, the two mappings lead to an equation for the harmonic function $J in the computational plane that has the same structure as eq. (11). But now C, has become the unit circle, so that the Green function is known:

G@,P’) =

1

2pp’ cos(8 - e’) + p’2 & ln[1p2- -2pp’ cos(8 - e’) + p’p”

1R

(14)

’

where (p, 19)are polar coordinates in the computational l-plane: { = p exp(i8). Using the properties of this function and eliminating the auxiliary function 4 in favour of Ic/and I+Q~ (see eq. (lOa)), the analogue of eq. (11) in the l-plane becomes

= p{ j

G&h P+bp@i) dli - t+P@o) - j

[G&o,

pi,)@‘@;) - G@,, ,$,)~,?py@(,)] dlb).

(15)

CO

Cl

Here, dl and v denote the arc length and the normal derivative in the computational plane. Eq. (15) is rewritten in a form that is more convenient for numerical and analytical treatment:

$

[G(8, e’)+ c]v(e’) dfY d (R(e) + f G(e, ef)Q(ef) de! - f L(e, eyqey =~($vode)/

de’} .

(16)

Here, the following definitions have been introduced. Describing the plasma contour Co by p = p. = g(O), and denoting by h = ]dz/dl] the scale factor of the conformal mapping, we define d=

ho(e)[g*(e) + g2(ey de,

(174

C= h/d,

(1nJ)

u(e)/f v(e) de = 6p(e)ho(e)[g2+ g2]Td.

(I7c)

Furthermore.

L(e, et) =

Gv(po,

p;)[g*(el) + g*(ey

_ (1 2_ITp2)g’2(1 + g”) - gg’(1 + g”) cos 4 + gg’(1 - g”) sin 4 I ’ [ (1 - 2gg’ cos 4 + g*g’*)(gz - 2gg’ cos + + g”)

(18)

where g = g(8), g’= g(0’) and 4 = 0 - 8’. The kernel G is obtained from eq. (14) by specializing to and p’ = p;, = g(8’).

p = p. = g(8)

J.M. Akketmans

Considering

now the particular

I Free-boundary

solution,

tokamak

109

equilibria

we have

(194 o(e) = @(po)[g2 + g2y2 m

=

-gPlg

+

1C kpU2kgZk(e) +t 2

n=l

k=l

&i + 4

2 (~2+ 2k)pU,+kjhkg”+2k(e) cos ne,

2 2 pn+kpkg”(e) COSfle -@(e)

n=l

W’b)

k=O

.

(I9c)

k=O

In eq. (19~) we have employed the properties of the Poisson kernel to calculate the contour integral. Inserting eqs. (17)-(19) into eq. (15) we obtain eq. (16). This equation is suitable for numerical treatment, which will be discussed in section 2.4. It gives a relationship between the plasma shape function g(8) and the computational quantity v(O) which is directly related to the poloidal magnetic field at the plasma edge. All quantities occurring in eq. (16) are expressed in terms of known conformal mapping quantities and the plasma shape g(8). Notice that only contour integrals along the plasma boundary appear. In the case of a vacuum outer region, eq. (16) simplifies considerably because the right-hand side vanishes in that case. 2.3. Pressure

balance

A second relation between condition that at equilibrium boundary:

plasma shape and poloidal the kinetic and magnetic

field at the plasma edge is provided pressure should match across the

by the plasma

(20) The jump in the fields is generated by the currents flowing on the plasma the poloidal variation of the poloidal field along the plasma boundary BP = B,[P + (1 - &/B;)(l where

we have introduced

surface.

+ EX)-~] “2 ,

the plasma

P-c, together magnetic

eq. (4) we find

(21)

beta

p = 2plB;. In this paper

From

(22) we employ

the high-beta

tokamak

approximation,

expressed

by (23)

with the additional conditions of eq. (3). In this ordering, plasma confinement results from the well of the toroidal field (in contrast to the low-beta regime, p - E’, where the poloidal field is

110

J.M. Akkermans

mainly responsible written as

for plasma confinement).

/ Free-boundary

tokamak

equilibria

To leading order in E, the poloidal field may now be

8, = &[ 1 - I&B: + p + 2p&x]“2 .

(24)

The term 2fl.e~ shows the influence of the toroidal geometry effects as a result of which the plasma is driven outward. At p = 0 the poloidal field is constant along Co. Increasing p implies an increasing variation of 8,, until a limit is reached when BP = 0 for x = A - a. To describe the amount of variation of the poloidal field along Co, we define

so that 0 G k2 s 1. Small values of k2 correspond to the low-beta regime, and values close to 1 describe the high-beta regime. The value k2 = 1 corresponds to the high-beta limit just mentioned. Rewriting eq. (24) in terms of k2, and changing to the dimensionless normalized poloidal field on Co, &, we obtain (see eqs. (6) and (7)) & = F(8)/(F)

(264

,

F(8) = [l + +k2(x - A - ~)/a]“~ .

The parameter

(26b)

k2 can also be related to the poloidal beta, & = 2p/(@,

according to

a&, = tk2/(F)‘. Combination

(27)

of eqs. (20), (22) and (8b) yields the scaling law for high-beta tokamaks

PI (28)

uap, = pq * 2/us ,

All of these quantities are expressed here in terms of the coordinates of the physical plane. Transforming to the l-plane with the aid of the mappings (13a) and (13b) and using definition (17~) the final result for ~(0) is (29) The poloidal beta becomes a~& = bk2

[f

!r0[g2+~2]1iZde/$v(t?)d0]2.

Eq. (29) is the basic result from pressure balance considerations plasma shape and the poloidal field at the plasma boundary. 2.4. Outline of the numerical

(30) that relates g(8) and u(e), i.e., the

procedure

Eqs. (16) and (29) are the two fundamental equations of our free-boundary problem. For a given plasma shape each of these equations yields a solution for the poloidal magnetic field at the plasma boundary. However, in general the two solutions for the poloidal field will not coincide. This is the way

J.M. Akkermans

/ Free-boundary

tokamak

111

equilibtia

in which the non-linear character of the free-boundary problem is manifested. Numerically, the solution of the problem can be obtained by the following three-step procedure: (i) for a given plasma boundary (i.e., g(8)), find the poloidal field at the plasma boundary (i.e., v(0)) from eq. (16); (ii) do the same from eq. (29); (iii) iterate on g(B) until the two solutions obtained for u(8) coincide (within some prescribed accuracy). This procedure for the solution of the free-boundary problem will be discussed below. First step. Eq. (16) is the reformulation of the Maxwell equations in the force-free outer region, leading to a Poisson (or in the vacuum case, Laplace) equation with Dirichlet-type boundary conditions. Hence, the problem is well-posed, in contrast to the fixed-boundary problem considered in ref. [5]. The application of conformal mapping has provided the possibility to exploit the techniques of fast Fourier transforms. Therefore, let us write down the Fourier expansion of the functions g(8) and v(0): m

g(e) =

&To +

x

a, cos

2

urn cosme,

mo,

urn

III=1

v(e)=&+ and let us define Km,

2

n-

ff

m=l

for any kernel

U,

=-$f

=$ f

g(O)cos m0 de,

@lb)

v(e)cosmede,

K(f3, 0’) Fourier

K(e, e’) cos mecos ne’ de

(314

coefficients

K,,,, such that

de’.

(32)

Because of the up-down symmetry of the problem we only need the cosine part of the expansions. Applying the operator (2/7r) $ de cos me to eq. (16) and using the above definitions, we obtain ma1

for

(33) Here, we have used L,,,,, = 2&,. The symbols p,,, q,,, r, denote the Fourier coefficients of the functions P, Q, R, and O,,= 2v,Ju0. For a given plasma boundary g(B) the right-hand side of this equation is known. Hence, eq. (33) can be formulated as a matrix equation Gu=f,

(34)

with known G and f. The solution of this problem can be found by means of standard linear equation solvers. Of course, in the numerical solution one has to truncate the number of Fourier harmonics. The code FREAK, written for the present free-boundary problem, allows for a maximum of 40 harmonics. The m = 0 component gives an expression for the unknown constant C measuring the magnetic flux between the plasma and the wall:

87~

= - 2 n=l

GO,% -

(1-

rpqold)Goo

+ 2rp(2r0

- p&d + 2&’

2

(G,,“q,, - Lo,,p,,)/d

.

(35)

n=l

Upon solving eq. (34) for the coefficients ii,, it is straightforward to obtain C from eq. (35). Note that the first step provides the poloidal field function u(O) apart from a normalization factor 42: this factor

112

J.M. Akkermans

/ Free-boundary

tokamak

equilibria

is found in the second step of the procedure. A final note concerns the kernels and their Fourier expansions. The kernel L is non-singular and therefore does not constitute a numerical problem. However, the kernel G is singular for 8 + 8’ (4 + 0). This difficulty is circumvented by writing

G(8, f3’) = -& ln[

1+J(4w3

2y2(1- cos 4) l-2y*cos4-y4

(364 (=I

The kernel J is non-singular. The Fourier can be found analytically. Then

of the term of eq. (36a) containing

the singularity,

(n = 01,

J mo >

Gm,=

coefficients

$ (1-

Jm -

(37)

y2”‘)&, ,

(n >O).

This completes the discussion of the first step of the procedure. Second step. We write the pressure balance equation (30) in the form v’(e) =

A(e)+ k%(e) .

(38)

For a given plasma boundary g(8), the functions A(8) and B(8) are known will be discussed below). Fourier expansion of eq. (38) yields

(the reason

for squaring

v(0)

d: = a,,, + k2b,,, , where a,, b,,, are the Fourier coefficients obtained in the second step. The quantities the first step, we can write d’, = &,,,

d’, =i

,

of A(B) and B(8), and dZ symbolizes the solution for v’(e) k* and vo/2 are still unknown. However, from the solution of

wo

u’*(e)cosmede,

where c, are known coefficients obtained imposing (throughout the procedure)

from i$,, in the first step. Now, k2 and vo/2 may be obtained

by

d: = d:’ ,

(41)

k2 = -(aocl - alco)l(bocl - hco),

(424

vo/2= [(boa1 -

Wb)

d:, = df

and

so that

blao)l(bocl - hco)YR .

In general, for m 32, we have d’, # dU,. Third step. This part consists of an iteration

procedure

on g(8) such that finally

df, = d!L for all m

J.M. Akkermans / Free-boundary tokamak equilibria

-1.0

1.0

0.0

Fig. 3. Plasma shape at high beta; a = 0.52, A = 0.41, a& physical plane.

-1.0

= 0.76. Circular wall, vacuum

1.0

0.0

region;

(a) computational

plane; (b)

within a prescribed accuracy. Typically, in the program FREAK the allowed mean quadratic error in the field is set equal to the value lo-‘. This iteration procedure can be schematically described as follows: (i) Determine the conformal mapping coefficients y, 6, &,, and p-L, from the wall shape and the known parameter pair (A, a). (ii) Guess a plasma boundary g(0) or u,,, (see below), starting with a circle, g(B) = y, in the first iteration. (iii) Determine d’, and d: as outlined before. (iv) Terminate if d’, = & within the prescribed accuracy for all m. If not, guess a new plasma boundary and go back to (iii). Guessing the plasma boundary in the (N + 1)th iteration is done according to the expression &‘+I)

= ur)

+

(y,,,(d; - df,,) ,

(m 32).

(43)

Since the coefficients d, are related to the square of the magnetic field at the plasma edge, the choice (43) implies that the plasma boundary is distorted corresponding to the pressure imbalance exerted at the boundary. Thereby it is assumed that to first approximation the mth harmonic of the plasma distortion responds to the mth harmonic of the pressure imbalance only. This displays the physical picture behind the iteration procedure. The procedure described is due to Goedbloed [15]. In his approach, CQ is a free parameter to be determined by trial and error. This is a time-consuming procedure since the optimum value of a turns out to be a function of m as well as A and a. In the present work we use the expression which is obtained from the analytical solution found by linearizing eqs. (16) and (29) with respect to the distortion of the plasma shape:

om = 2y[(m

(1- y’“) - 1+ T’y)+ (m + 1- Py)y2”].

Starting with a circle (a,,, = 0 for rn L l), this ensures that all coefficients order in S at the second iteration. This will be demonstrated in section analytic justification for the assumption underlying eq. (43), viz. that the shape essentially couples to the m th harmonic of the pressure imbalance. The scheme outlined here results in a successful iteration procedure. vacuum region, convergence is obtained up to a@, = 0.8 for the maximum

a,,, will be correct to leading 3. Moreover, we will give an m th harmonic of the plasma For a circular wall and a number of 40 harmonics. An

114

J.M. Akkermans

/ Free-boundary

tokamak

equilibria

analysis of the harmonic content shows that the limitation of our non-linear procedure is due to the cut-off of the number of harmonics. A typical case showing what kind of plasma boundaries are obtained at very high betas is given in fig. 3. Further results will be discussed in the remainder of this paper.

3. Circular

wall

In this section we discuss the case which was numerically investigated by and its comparison with the numerical provides a reference case for discussing currents surrounding the plasma. The Maxwell and pressure-balance

of a circular wall with a vacuum region surrounding the plasma, Goedbloed [15]. Here, we concentrate on an analytical approach results. The case of a circular wall and a vacuum outer region the more complicated effects of wall shaping and of force-free equations

(16) and (29) here simplify

to

[G(t), Or)+ C]u(O’) de’ = 0,

(45)

1’

1 * (1 - Sy)(y - 6g2 - (1 - Sy)g cos e) l/T v(0) = ho(B)[g2 + g2]‘j2[ 1 - T k y (1 + 2sg cos 8 + 62g2) where

ho(e) is the scale factor h,,(8) = ]dz/d&

of the Moebius

transformation

= (1 - a2)/(1 + 26g cos 8 + a2g2).

at the plasma

(46) boundary (47)

This mapping centres the plasma: t = A 2 a +cJ = _’ y, where y = g(0) = g(r). From this condition we find the relationship between the known input parameters A and a and the new parameters 6 and y,

A = 6(1-

r’)/(l-

S2y2),

(484

a = y(l-

S2)/(1 - 6*y2),

(48b)

or, inversely,

6 = (I- a2 + A2)/2A - [((l - a2 + A2)/2A)2 - 1]1’2,

y = (I+

a2-

A2)/2a - [((l+

a2- A2)/2a)*-

lll’.

(48~)

(48d)

The normalized poloidal field at the plasma edge & is obtained from eq. (17c), and the poloidal beta is found from eq. (30). In order to proceed with an analytic approach, we linearize eqs. (45) and (46) with respect to the deviation of the plasma boundary from a circular shape. First, we observe that in the case of zero beta (k2 = 0 in eq. (46)) the solution for the plasma cross-section is a circle with A = S = 0. (In view of the toroidal geometry it may seem strange that the solution for the plasma shape could be a circle that is centred with respect to the tube. This fact derives from the expansion in the inverse aspect ratio, where the outward shift is a higher-order effect in E at /3 = 0.) The poloidal field then is constant along the circumference of the plasma. Accordingly, in the case of zero beta

J.M. Akkermans

/ Free-boundary

tokamak

11.5

equilibria

ud2 = vo/2 = y ) cr, = v, =

0,

(494 (m 5 1).

Wb)

Finite-beta effects will result in an outward displacement of the plasma and a deviation from a circular shape. This suggests to consider the plasma deformation, the outward shift, the plasma beta, and the angular variation of the poloidal field as small parameters in a power series expansion. Defining s=ui3,/2-

y,

(50a)

urn cos me,

GObI

00

u(e) =

2 VI=1

we consider s and ~(0) as the small parameters with respect to the plasma shape. Note that the m = 0 and m = 1 Fourier components of the plasma shape are not independent of the coefficients for higher m.

m

u1=

-

x

@lb)

oinl+1,

m=l

which follows from the property g(0) = g(r) = y. Consequently, s = -Z”,=, u2,,,. In addition to s and u(e), in eqs. (45) and (46) we have as small parameters v,(m Z=l), k* and 6. With the aid of the foregoing considerations, to first order in the plasma distortion the Green function for the unit circle can be represented as G(8, e’) = H(8,

et) + Lute)+ u(e’)]K(e, t-l’)

(524

)

where

H(e,

et)

G

& ln[ 1 _ 2(ud2)2(1 - cos4)

]

2(uo/2)2 cos 4 + (uo/2)4

qe, e’) = The Fourier

1 - y4 2y2 cos 4 + y4).

47ry(l-

representation

HO, = 4 ln(u&)

Ken = 2/y * 80, ,

(52b)

’

of the kernels

* SO”,

(52.c)

H and K is diagonal:

1 H,,=-m[l-(uo/2)2m].Smn,

K,,,, = yzm-l - S,,,,,

Next, we use these results

in order

(m >O),

(53a)

(m >O).

to solve the integral

Wb) equation

(45). Applying

the operator

116

J.M. Akkermans

f

1 Free-boundary

tokamak

equilibria

dOcos m6

to eq. (45) we find, after some algebra, as solution of the integral equation, to leading order in CT,,,, d:, = ~$2, d’, = v&,/2,

(54) (55)

(m 2 11,

where

Furthermore. (57) In first approximation, we have a linear relationship between the perturbations of the magnetic field at the plasma edge and the plasma shape. The factor m occurring in eq. (56) stems from the singularity in the Green function (cf. eq. (37)). In the pressure balance equation (46) the small parameters are 6, k2 and the deviation of the plasma boundary from a circle. Squaring v(O) and applying the operator

to eq. (46), we find to leading order in S and o,,,: d;I = 2y(l-

+k2)+ 4ys,

d”, = 2yom + 2(-l)“-‘(m

(58) + l)S”-‘ym+‘(mk2/8

- Sy) ,

(m 3 1).

(59)

The first term in eq. (59) gives the contribution of the plasma distortion to the pressure-balance. The second term, which is the most important one, contains the contribution of the outward displacement of the plasma. The condition (41) now yields expressions for the unknowns vO/2 and k2: v0/2 = ~(1 - tk2) + s , k2 = 8Sy + 8~1y/(l-

ew

r’),

where cl is determined according to eq. (51b). The solution for the shape coefficients a, is now obtained from the requirement (55) and (59). This gives

(61)

d’, = d!!, and eqs.

J.M. Akkermans

/ Free-boundary

m i (m’- l)Smym+l(l - y2” (m-l)+(m+l)y*” )7

%=(-I)

117

rokamak equilibria

(62)

(m 22),

to leading order in 6. This is our basic analytical result for the case of a circular wall surrounding the plasma. All other physical variables of interest can now be expressed y by means of eqs. (60) (61) and (62). By employing eqs. (4%~) and (48d) one physical quantities in terms of the original input parameters A and a. This confirms whole problem is determined by fixing these two parameters (apart from the scaling

and vacuum fields in terms of S and may express these the claim that the constants 6, E, and

Bo). It is seen from eq. (62) that the individual plasma shape coefficients a, are of order 8” (i.e., A”’ in physical space). Having neglected terms of order o’(0) in the derivation, this implies that the full solution g(t9) for the plasma shape is correct up to third order in S (or A). This solution can be written as g(8) = y - aZ(i - cos 28) with u2 and u3 Now, what putational and its inverse C(z)

u3(~~~

8 -

cos

38))

(63)

given by eq. (62). about the plasma shape in physical space. 7 The transformation between the comthe physical plane is effected by the Moebius transformation z(C) = (L + a)/(1 + SY) and = (z - a)/(1 - Sz). It follows that in the physical z-plane the plasma cross section is given

by f(8)

= a + A cos 6 - (CJ*+ A2/4a)[l

- cos 219]- (Us + Ao2/a(l-

which is correct to third order in A. To interpret this result, plasma centre (x = A, y = 0). Eq. (64) then transforms to (c3 + AaaJ(l

f(T) = a - a*[1 - cos 27]-

a’))[cos

it is instructive

6 - cos 361, to shift the origin

- u’))[cos r - cos 371,

2Au =

1 _

The second-order 6&y)=

>I1

A2

u2

l+

(I-

expression

l-262 (I_

u2)2

+

for the poloidal

-

Before comparing the above numerics are in order.

field at the plasma

A*(1- 2~~)cos

:_“a,2

cos

results

6

+

(l-u’)2

with the numerical

plasma shifted -2uzlu. axis.) central

w

.

u2)2

to the

(65)

which may be compared with the result (63) for the computational plane. From eq. (65) we can draw the following conclusions with respect to the free-boundary shape. To first order in the plasma shift the plasma cross section is represented by a circle outward over a distance A. To second order in A, it is a shifted ellipse with ellipticity e = (Throughout this paper, the ellipticity is defined as e = (vertical axis - horizontal axis)/horizontal To third order we have a D-shape like plasma shifted outward, with the flat side turned to the axis. For the poloidal beta we find

@P

WI

edge reads

26

(67)

. computations,

some remarks

with regard

to the

J.M. Akkermans

118

/ Free-boundary tokamak equilibtia

Table I Comparison of numerical and analytical results for the plasma Fourier coefficients cm. Circular wall, a = 0.4, A = 0.3

r’ = 0 (a&& = 0.338)

r’ = 1 (as&, = 0.506)

m

cm (num.)

a, (anal.)

Us (num.)

a,

0

0.97209 -0.00937 -0.03320 0.00955 -0.00143 -0.OOoO9 0.00017

0.97268 -0.00883 -0.03284 0.00845 -0.00182 0.00037 -0.00007

1.01654 -0.01780 -0.05280 0.01755 -0.00426 0.00059 0.00021

0.98896 -0.01246 -0.04001 0.01185 -0.00274 0.00058 -0.00012

1 2 3 4 5 6

The analytical results derived justify some aspects of the First, eqs. (55), (56), and (59) justify, as a first approximation, harmonics of the plasma boundary and the poloidal field, the function has a diagonal representation to leading order in the a$ = 0, in the second iteration a:) coincides with eq. (62) if

anI = 2y[(m

shape

(l- Y2") - 1) + (m + l)y2”]

(anal.)

numerical iteration scheme (see eq. (43)). the separate consideration of different main reason being the fact that the Green plasma distortion. Second, if we start with we use in eq. (43)

(68)

’

which is expression (44) for r’ = 0. This choice ensures that all individual coefficients will have the correct value to leading order in 6 in the second iteration. Although this does not guaranfee a converging iteration scheme, it turns out to work well for the whole parameter space also at later steps of the iteration. A limit to the niethod outlined here is presented by the cut-off for the number of harmonics. Approaching the high-beta limit (A’ = A/(1 - a)+ 1 or k*+ 1, see section 7), the plasma becomes more and more distorted and the number of harmonics needed to accurately describe the plasma boundary quickly increases. With a cut-off at 40 harmonics the maximum values reached are approximately A’ = 0.90 and k* = 0.99. A comparison of analytical and numerical results for the plasma shape is presented in table I and fig. 4 for a typical case. For clarity, in fig. 4 the third-order result is not shown. The maximum deviation

-1 .o

0.0

Fig. 4. Comparison of numerical (full line), first-order (dot-dashed cross section for a circular wall; a = 0.4, A = 0.3, r’ = 0.

1.0 line) and second-order

(dashed

line)

results for the plasma

J.M. Akkermans

f Free-boundary

tokamak

equilibria

119

0.6

0.4

0.2

0.0

0

Fig. 5. Numerical (full line), first/second-order Circular wall, vacuum outer region.

(dot-dashed line) and third-order

(dashed line) results for the poloidal beta.

between the numerical and the third-order analytic result is less than 1% in this case, showing that there is a good global agreement. Considering the separate plasma coefficients, we see from table I that there is a good agreement for low values of the Fourier order m. For increasing m, the relative error increases. This is due to the neglect in the analytic approximation of cross-terms of the plasma shift and the lower-order shape coefficients. Finally, we remark that although for fat plasmas (i.e., plasmas having a large filling factor with respect to the tube) the shape becomes more circular, the relative error in the shape coefficients depends only weakly on the plasma thickness. In fig. 5 we compare analytic and numerical results for the poloidal beta as a function of the relative plasma shift A’ = A/(1 - a). Good agreement is obtained even for quite high values of A’, particularly for fat plasmas. The analytic third-order result indicates a possible cross-over of the curves for thin and fat plasmas. This cross-over is indeed found numerically.

4. Wall shaping Let us now investigate the case of non-circular region. Numerically, we have studied wall shapes X,,H = cosk

+ Ci sin x + C, sin 2x),

YWtil= C, sin x .

walls, with the plasma being surrounded represented by the expression

by a vacuum

(69) (70)

This includes circles, ellipses (Ci = C, = 0, C, # l), D-shapes (C, # 0, C, = 0) and racetracks (C, = 0, C’, f 0). Generally speaking, the coefficients q& of the mapping (13a) that conformally transforms regions with boundaries of this type onto a circular region, can be found only by numerical means.

120

J.M. Akkermans

/ Free-boundary

tokamak

equilibria

However, we can say something about the expected influence of the wall shape on the plasma shape and the poloidal beta by analyzing walls that are slightly deformed from a circular shape. This problem we address first. We represent such slightly non-circular shapes by m

n)7mcos m6,

rw,tr(fi) = qO/2 + 2

(71)

m=l

where that

we treat

the n7, (m 2 1) as small

rIoi = 1 - 2 7)&,

771=

-

ItI=,

To leading

order

4 m+1=

77,

2 m=l

41

al>,

=

From

the property

rWal,(0)= r,,,l(rr) = 1 it follows

r/Zm+l.

in the wall distortion, (m

parameters.

(72)

the conformal

mapping

(13a) is known

to be

770/2 .

(73)

Now we may perform a similar analysis as in the foregoing section. Note that the n’s are small parameters that are independent of the parameter 6. The solution of the integral eq. (45) is given by eq. (55) as in the circular case. The wall shape explicitly enters the pressure balance equation (29). To leading order in n,,, its solution is dE = &(circle) where

dt(circle)

+ 2(m + l)n,~~+~

for the plasma

a, =

(m 2 1) 9

is given by eq. (59). Here,

(m + l)[n,

I?#+1

we have neglected

cross-terms

of 6 and nlm. Further,

shape

(75) coefficients

we obtain

+ (-l)“_‘(m - l)s”]y”+‘(l(m - 1) + (m + l)y2”

The second term is the circular contribution from the wall shape and behaves like em--qnY

(74)

.

dT: = (no/2)&(circle) Solving

,

72”)

7

(ma:).

of eq. (62). The first term contains

.

(76) the leading

contribution

(77)

Consequently, to this order, a wall deformation of Fourier order m induces a change in the plasma shape in the same Fourier order. This change is smaller if the plasma is farther away from the wall. Regarding the influence of a wall distortion n1 on the plasma shape coefficient a,, where 1# m, we infer from eq. (29) that the leading contribution is of order G“v, with ]k f I] = m. Accordingly, these contributions are due to finite-beta effects. The next question concerns the influence of the wall shape on the poloidal beta. To leading order in the plasma shift and the wall distortion we find

@p

=

%‘(a

-

?1)[

I-

3772 +

m$l

12mY2-3

,

(784

J.M. Akkermans

1 Free-boundary

tokamak

-1.o

-1.o

equilibria

121

0.0

0.0

1.0

1.0

Fig. 6. Comparison of the plasma cross section and the poloidal beta for different wall shapes at a = 0.52, A = 0.41, r’ = 0. See fig. 3b for the circle: (a) Ellipse (Cl = 0.0, Cz = 0.0, C, = 1.5), a.&, = 0.82; (b) racetrack (Cl = 0.0, C2 = -0.5, C3 = 1.5), aa& = 0.75; (c) D-shape (C, = 0.5, CZ = 0.0, C, = 1.5), aa&, = 0.89; (d) inverted D-shape (C, = -0.5, C2 = 0.0, C, = 1.5), a~& = 0.67.

where 77 =

7)1-

R/(1 - 73.

(78b)

Thus, even and odd components produce different effects. If we have a wall shape with odd Fourier components, in the case of zero beta the plasma is not centred with respect to the tube. For instance, for a D-shaped wall weighted to the left the condition a& = 0 corresponds to a situation in which the plasma centre is also shifted to the left. This is understandable because the occurrence of odd components implies that there is no left-right symmetry of the wall. Even components, on the other hand, only enter the expression for the poloidal beta in a combination with the plasma shift. For elliptical walls we obtain a nice and simple expression in terms of the physical parameters A and a: a&?, =

3

[l+ ie(l+

a’)],

(79)

122

J.M. Akkermans

ae?

I

/ Free-boundary

tokamak

equilibria

P 0.6

0.0 0,. 0

0.2

0.4

0.6 d

0.8

1.0

A’=A/(l-a1

Fig. 7. Poloidal beta versus the relative plasma shift for a D-shape (C, = 0.5, C2 = 0.0, C, = 1.5) for a = 0.2 and a = 0.7. Full lines: circle, r’ = 0; dot-dashed lines: D-shape, r’ = 0; dashed lines: D-shape, f’ = 0.5.

which holds to leading order in the ellipticity e of the wall. It is seen that at a@, = 0 the plasma is centred with respect to the tube. For given A and a, elongation of the wall has a favourable influence on the poloidal beta. The wall distortion does not show up to leading order in the expression for the poloidal field 6, at the plasma edge. For the same plasma shift and width, the plasma shape is shown for different wall shapes in fig. 6. The corresponding result for a circular wall is found in fig. 3. A parameter scan of the poloidal beta versus the plasma shift at different plasma widths for elliptical walls is given in ref. [6]. Fig. 7 shows a similar scan for a JET/INTOR-like D-shape, for a thin (a = 0.2) and a fat (a = 0.7) plasma (dot-dashed lines). At low values of the plasma shift, the general trends derived above regarding the influence of the wall on the plasma shape and the poloidal beta are confirmed by the numerical results. For high values of the plasma shift, we see from fig. 6 that for all wall shapes the outer side of the plasma takes on the shape of the wall. The inner side of the plasma becomes increasingly flatter when the plasma shift is increased. The poloidal beta is increased from UE/~~= 0.76 for a circular wall to US&, = 0.82 for an ellipse with e = 0.5. D-shaping leads to a further increase to a.#, = 0.89. The favourable influence of wall shaping, in particular the triangularity, on the poloidal beta is further demonstrated in fig. 7. The improvement is larger for fat plasmas because their shapes correspond more closely to the shape of the wall. The fact that in fig. 7 the curves for the D-shape do not pass through the origin confirms the conclusion derived above that for a wall shape containing odd Fourier components the plasma is not centred with regard to the tube at zero beta.

5. Force-free currents In order to separately study the effects of force-free ourselves to circular walls. The force-free constant

currents surrounding the plasma we first restrict r’ occurring in the integral equation (16) is

J.M. Akkermans

/ Free-boundary

tokamak

123

equilibria

considered a free parameter. To be compatible with the high-beta tokamak ordering, however, it must be of order one with respect to the inverse aspect ratio. To obtain an approximate analytical solution for the plasma shape, we first note that, to leading order in the plasma distortion, L(8, 0’) = rK(0, 0’) in eq. (16), where K is defined in eq. (52~). Furthermore, for a circular wall we have p. = S, ,u, = (-l)m-18m-1(1 - S’) (m 2 l), which is to be inserted in eq. (19) for the functions related to the particular solution of Poisson’s equation. Solving eqs. (16) and (29) for the plasma shape coefficients, we find, to leading order in S, a,

=

l)S”ym+’ (1 - 72”) + r’S”y”[$r(m

(-l)“_‘{(m’-

- m(1 - r’) + y2(1-

r’“)]}[(m

+ l)(l - y2)(1 - 7”“)

- 1+ T’y) + (m + 1 - rly)y*“]-’

.

(80)

For r’ = 0, this equation reduces to eq. (62) for the case of a vacuum outer region. From eq. (80) it is derived that employing the expression (44) in the numerical iteration scheme (43) yields values for urn that are correct to leading order in S in the second step of the iteration. For the poloidal beta we obtain r’)] )

U&P, = S[2y + rl(lor in terms of the physical a&& = &

(814

quantities

[2a + r’(1 - a’)] .

@lb)

This expression is valid to first order in A. For a simple physical interpretation of this result, force-free constants r and r’ = T/u. From the definition ref. [5], we obtain

we first discuss the physical meaning of the @a) for r, which is identical to the one used in

for the case of a low-beta plasma and a circular wall. Here, j,+w and j,+,Qe are the toroidal current densities in the force-free region and the plasma, respectively. r is also related to the safety factor 4 in the outer region, since 4 =

q*(r/ay[i

+ $(9/d

-

q-1 ,

(83)

in the same approximation as used above. Thus, for r = 2, i.e., for equal force-free and plasma current densities, the magnetic field in the outer region has a constant pitch. According to ref. [S], r is the appropriate parameter for describing the effects of force-free currents in magnetohydrodynamic stability. In describing the equilibrium properties, however, it is the ratio of the total toroidal currents, not that of the densities, which is important. From eq. (82) we obtain

4, FFII+,p,pP = r’(1 - u2)/2u . Eq. (84) allows a simple

physical

w interpretation

of the expression

(81) for the poloidal

beta, since

(85)

124

J.M. Akkermans

I

I Free-boundary

tokamak

equilibria

I

I

I

0.4

0.6

0.8

0.8 a@

P 0.6 t

-I D.U

0.2

j

1.0

A’=A/(l-a)

Fig. 8. Influence of force-free currents on the poloidal beta for a circular wall; a = 0.2 and 0.7. Full lines: r’ = 0; dot-dashed lines: r’ = 0.5: dashed lines: r’ = 1.

This implies that, for low values of the plasma shift, the poloidal beta is proportional to the ratio of the toroidal currents flowing in the force-free region and on the plasma surface. If these currents are of the same magnitude, the poloidal beta is approximately doubled. This favourable influence of the force-free currents on the equilibrium poloidal beta is also clear from the numerical results shown in fig. 8. It is seen, however, that for higher values of the plasma shift the increase in the poloidal beta due to the force-free currents becomes relatively smaller. The reason for this will be discussed in section 7, where it is argued that, for A’+ 1, the limiting value for the poloidal beta does not depend on the amount of force-free currents. From table I it follows that, for a fixed plasma shift, force-free currents give rise to a stronger distortion of the plasma shape than in the vacuum case. This can also be deduced from eq. @Cl),which is seen to be in reasonable agreement with the numerical results at intermediate values of the poloidal beta. For m 9 1 the plasma shape coefficients behave as urn - P(l + 14,FF/14,Ppe) to leading order in S. At high values of the shift, the inner side of the plasma is more flattened and the outer side follows the wall shape more closely than in the case of vacuum fields. If the poloidal beta, instead of the shift, is kept constant the plasma is more circular and the outward displacement is smaller as compared to the vacuum case.

6. The general case In this section we consider the combined effects of wall shaping and of force-free currents surrounding the plasma on the equilibrium. An approximate analytical solution for the plasma shape is found along the same lines as discussed in the previous sections. To leading order in the shift S and in

J.M. Akkermans

the deviation nm of the wall from a circular (29) is given by

o;,

=

uI

m

+

/ Free-boundary

shape,

cm+lh,~m+‘(l -

y2”)+ r177mr”[m(l(m-l+T’y)+(m+l-r’y)y2”

rokamak equilibria

the solution

125

of the equilibrium

y2) - y2(1 - r’“)]

r')l

(16) and

(86)

2

where a; is given by eq. (80). The first term, & contains the leading contribution displacement, i.e., the finite-beta effects. The second term is the main contribution when neglecting finite-beta effects. For m Z+ 1 it behaves as

a;, - %Yrn[Y+ ml-

equations

due to the outward from the wall shape

(87)

9

which is to be compared with the result (77) for the vacuum case. Furthermore, we can derive from eq. (86) that Sgn(&,,,/X’) = Sgn(n,) for all m. Consequently, force-free currents enhance the coupling between the distortion of the wall from a circular shape and that of the plasma. This can be physically understood in terms of the interaction between the toroidal currents that are flowing parallel in the force-free region and on the plasma. The enhancement is roughly proportional to the ratio I +,rrw/14,ppof the force-free and plasma currents. This general conclusion from the analysis is confirmed by the numerical results (see fig. 9). As is shown in fig. 9a, for a D-shaped wall the plasma becomes much more D-shaped when force-free currents are present. Also the effects of finite beta are seen. They lead to a plasma that is flattened at the inner side. As already pointed out in section 5, force-free currents enhance this flattening. In fig. 9b we show the results for a racetrack wall similar to that of SPICA-II, a screw-pinch experiment presently under construction at Rijnhuizen [6]. In such experiments force-free currents are routinely produced.

L

-1

.o

0.0

1.0

-1 .o

I

0.0

1.0

Fig. 9. Effects of force-free currents on the plasma shape. (a) D-shape (C, = 0.5, C2 = 0.0, C, = 1.5) for a = 0.7 and A = 0.05. Solid line: r’ = 0, ae& = 0.23. Dashed line: r’ = 0.5, a&& = 0.31. (b) Race-track (Cl = 0.0, C2 = -0.35, C, = 2.33) for a = 0.24 and A = 0.24. Solid line: f’ = 0, a&& = 0.17. Dashed line: r’ = 1.0, ae& = 0.46.

126

J.M. Akkermans

I Free-boundary tokamak equilibria

The wall being strongly elongated, the force-free currents induce a strong elongation of the plasma. In the vacuum case, the plasma elongation is about 1.2. In the force-free case, corresponding to a force-free current that is roughly twice the plasma current (cf. eq. (84)) the plasma elongation is increased to 3.6. This conclusion agrees with that of Hoekzema [20] who performed a free-boundary analysis for the SPICA experiment. For the poloidal beta the analysis leads to a very simple result. To leading order in the wall distortion and the plasma shift, we obtain a&P, = a&&[1 + rl(1 - r2)/2rl ,

cw

where a@, is given by eq. (78). In first approximation the expression for the poloidal beta is factorized into a part describing the effects of the wall shape, and a part containing only the effects of the force-free currents. With the aid of eq. (84) we find a simple physical interpretation for the expression for the poloidal beta:

This equation is valid for small values of the plasma shift. We conclude that, for any wall shape, force-free currents have a favourable influence on the poloidal beta. This is demonstrated in fig. 7, where the numerical results for the poloidal beta are given as a function of the relative plasma shift, for a JET/INTOR-like D-shaped wall. A case is shown with and without force-free currents (dashed and dot-dashed lines respectively). Also the reference case, a circular wall and a vacuum region, is shown (solid lines). Compared to the reference case, the combined effects of wall shaping and force-free currents lead, for given plasma shift and width, to a substantial increase in the poloidal beta. Fig. 7 may be compared with fig. 8, which gives the effects of the force-free currents only (dot-dashed curves). 7. The high-beta limit In this section we address the problem of the poloidal-beta limit. This is the limit d’+ 1, i.e., the situation in which the plasma displacement is so large that the gap between the plasma and the wall is

-1 .o

0.0

1.0

Fig. 10. Limiting plasma shape for A ‘+ 1. D-shape (Cl = 0.5, Cz = 0.0, C3 = 1.5). a = 0.7 and a&&, = 1.26. This is the limiting configuration of the situation depicted in fig. 9a.

J.M. Akkermans

/ Free-boundary

tokamak

equilibria

127

closed. The limiting plasma shape is already suggested by figs. 3, 6, and 9a. At the inner side the plasma becomes completely straight, whereas the outer side assumes the shape of the wall (see fig. 10). That this is indeed the limiting equilibrium solution can be checked by insertion in the basic equilibrium equations. Force-free currents do not change this limit configuration as may be seen from physical considerations. At the flat side the plasma current is purely poloidal, while at the outer side the toroidal plasma current is completely compensated by equal and opposite currents flowing in the wall, as a consequence of Amp&e’s law. Accordingly, in the limit there is no interaction between the plasma currents and the (toroidal) current in the force-free region. The poloidal magnetic field &, at the outer plasma edge assumes a maximum value 2[a~&,,,~]“~ at the horizontal axis. It monotonically decreases along the circumference, until a zero is reached at the sharp edges of the plasma boundary. At the inner side of the plasma boundary & equals zero. A similar statement applies to the poloidal field along the wall. In view of this limiting situation, it is not surprising that the numerical method described in section 2 breaks down at high values of,the shift (in practice, at A ’ = 0.90). The limit A’ = 1 cannot be reached by an accurate continuously increasing A’, since the sharp edges of the plasma cross section prevent representation in terms of a finite Fourier expansion. Also the conformal mapping method breaks down, because in the limit S + 1 the Moebius transformation (13b) becomes singular. Hence, the numerical method proposed is fundamentally limited by the cut-off for the number of Fourier harmonics. Given the limiting plasma configuration discussed above, however, we can compute the limit for the poloidal beta directly from eq. (27). From the observation that the poloidal field becomes zero at the inner side of the plasma edge it follows that k* = 1 at A’ = 1 (cf. eq. (25)). The numerical results for the maximum poloidal beta computed from eq. (27) are shown in fig. 11. Results are given for different wall shapes in order to give an impression of the effects of the elongation and the triangularity of the wall. For some wall shapes it is possible to obtain analytical expressions for the maximum poloidal beta.

1.8

I

I

I

I

1.6 acP

t

P

1.4

1.2

-+.

0.6

-k._

0.X.0.0

_.L I

I

I

I

0.2

0.4

0.6

0.8

--‘T-”

1.0

Fig. 11. Maximum poloida beta as a function of plasma width for the wall shapes depicted in fig. 6. Solid line: circle. Dashed line: ellipse (Cl = 0.0, c2 = 0.0, C3 = 1.5). Dot-dashed line: D-shape (CI = 0.5, C2 = 0.0, C3 = 1.5). Cross-dashed line: inverted D-shape (CI = -0.5, c2 = 0.0, C3 = 1.5). Dotted line: Racetrack (C, = 0.0, Cz = -0.5, C3 = 1.5).

J.M. Akkermans

128

/ Free-boundary

tokamak

equilibria

1.8

1.6

0.8

0.6 t

0

-k._ .a-:0.0

__.4-I

I

I

I

0.2

0.4

0.6

0.8

-.

1.0

Fig. 12. Effects on the maximum poloidal beta due to different parametrizations of the wall shape in the case of a racetrack. Full line: racetrack with elongation E = 2.33. Long-dashed line: numerical parametrization (C, = 0.0, C2 = -0.35, Cs = 2.33). Short-dashed line: racetrack with E = 1.5. Dotted line: numerical parametrization (C, = 0.0, C2 = -0.5, C, = 1.5).

For elliptic

walls we find

ae& max= aa

[l+ e($- a)]-)+ e(&+ fia(la))]E(a)-

1

[l+

[l +te]

arcsinG

(W

(1- a)[1 + e(&+ &a)]K(a) I 2 ’

to leading order in e. However, eq. (90) also gives a very good approximation for high values of the ellipticity. K(a) and E(a) are the complete elliptic integrals of the first and second kind, respectively. Eq. (90) may be considered to be the high-beta analog of the low-beta expression (79). It is seen that the value of the maximum poloidal beta is not very sensitive to the ehipticity. The influence of elongating the wall may also be studied by considering a rectangular wall. For a rectangular wall with vertical elongation E we have

a+Lax

=

1 i

-

1

2 .

EjJ4a 3

(91)

1

Again we see that elongation in itself does not have a very strong influence. By increasing E from zero up to infinity the maximum poloidal beta is increased from 9/16 to 1. The effects of D-shaping of the wall might be studied in a simple approximation by considering a triangular wall pointing outward. In this case

a%,

max

where

= $[l + sin aI2 ,

a is the angle

between

(92) the horizontal

axis and the outer

side of the wall;

CYis related

to the

J.M. Akkermans

/ Free-boundary

rokamak equilibria

129

elongation according to sin (Y= E/d/4+ E*. We note that eqs. (91) and (92) coincide with the expressions found from a fixed-boundary analysis [5], as they should. We conclude that the combination of triangularity and a moderate elongation has a favourable influence on the poloidal-beta limit. Similar conclusions may be drawn from the numerical results shown in fig. 11. A few remarks are in order. The value of the maximum permissible poloidal beta in the limit a + 0 (where the plasma is a point plasma situated at the wall) appears to be strongly dependent on the way the wall shape has been parametrized. This is shown in fig. 12, where we compare the results for the numerical approximation (70) to the racetrack shape with the exact values given by the expression

a,$,, max=

I

(E - 1) + d/a(l

(E - 1) + 2[E(a)

-

a) + arcsind/a

- (1 - a)K(a)]/X&

* I

’

(93)

The short-dashed and dotted lines give the results for the exact racetrack shape and its numerical approximation, respectively, for the case E = 1.5. The solid and long-dashed lines give the results for the SPICA racetrack which has elongation E = 2.33. The wall shapes themselves are shown in figs. 6b and 9b, respectively. If the wall is completely straight near the horizontal axis of symmetry, as is the case for a rectangular wall or a racetrack, one has a.$, + 1 for a + 0. Small deviations as occur in the numerical parametrization of the wall, lead to large changes in the poloidal beta for a +O. However, this limit is of no physical significance. For higher values of the plasma width the numerical results are seen to be in good agreement with the exact values of eq. (93). The other numerical wall shapes considered, including the inverted D-shape, are more or less elliptical near the horizontal axis. For a +O, the plasma only “sees” this part of the wall and so cannot distinguish between the different wall shapes. Noting that according to eq. (90) the poloidal beta does not depend on the ellipticity for a + 0, this explains that the poloidal beta has approximately the same value at vanishing plasma width (see fig. 11). Considering now the case a# 0, it is easily seen by inspection that for thin plasmas the limiting plasma cross-sections for a circular or elliptic outer shell have stronger triangularity than for a racetrack. Hence, one would expect that the poloidal-beta limit for a racetrack is smaller than for a circular or an elliptic wall in the case of thin plasmas. For fat plasmas, however, the situation is reversed. These conclusions are confirmed by the numerical results shown in fig. 11. The physical interpretation of the above results is not straightforward, however. It may be argued that consideration of the poloidal-beta limit is not significant from a physical point of view, because the limiting situation will not occur in practice. On the other hand, the tendencies discussed above with respect to the effects of wall shaping agree with the conclusions arrived at in the previous sections regarding the low-beta regime. These conclusions are not drastically changed if they are not formulated in terms of the plasma width but, instead, are related to another parameter, e.g., the volume compression ratio. Consideration of the poloidal-beta limit is also of interest with respect to the force-free currents. The effects of these currents can simply be summarized by saying that if force-free currents are present the plasma is closer to the high-beta limit, with respect to both the plasma shape and the poloidal beta, at moderate values of the shift. The equilibrium limitations on the plasma beta itself follow from the scaling law (28) for high-beta tokamaks. The equilibrium beta depends both on the poloidal beta and on the total current flowing in the plasma. The limit on the plasma current is to be determined from stability considerations. 8. Summary A method to obtain free-boundary been described. Within the framework

equilibria of high-beta of a skin-current model

tokamaks with force-free currents has and the high-p tokamak ordering p - E,

130

J.M. Akkermans

/ Free-boundary tokamak equilibria

a Poisson equation for the poloidal flux is derived, subject to the boundary conditions that the plasma edge as well as the wall are surfaces of constant flux. Applying Green’s theorem, Poisson’s equation is transformed to an integral equation for the poloidal magnetic field at the plasma edge. A second equation for the poloidal field is provided by the requirement of kinetic and magnetic pressure balance across the plasma boundary. Non-circular cross sections of the wall are investigated by mapping the region interior to the wall onto the unit disk by means of a conformal mapping. Next, the plasma is centred by means of a Moebius transformation so that the techniques of fast Fourier transformations can be used in solving the equilibrium equations. The solution of these equations is found through an iterative procedure on the Fourier harmonics of the plasma boundary. In each iteration they are changed in response to those of the pressure imbalance until pressure balance is established. The response coefficients are obtained from the analytical solution found by linearizing the equilibrium equations with respect to the plasma shape. This provides a fast and accurate scheme that converges even for highly non-circular plasma and wall shapes. Apart from being helpful with respect to the numerics, the analytical solution of the linearized equations is useful in explaining the characteristics of free-boundary equilibria. A variety of wall shapes could be studied due to the speed of our method. The free-boundary plasma shape corresponds to the shape of the wall at the outer side, whereas at the inner side it becomes increasingly flatter with increasing outward shift of the plasma (i.e., with increasing poloidal beta). In the limiting situation, when the poloidal field becomes zero at the inside and a maximum for the poloidal beta is reached, the inner side of the plasma becomes completely straight, while the outer side coincides with the shape of the wall. Cross-sectional shaping of the wall favourably influences the poloidal beta. This applies to D-shaping in particular. Cross-sectional shaping of the plasma, the discussion of which was initiated in refs. [21-241, has recently also been investigated within the framework of the flux-conserving tokamak (FCT) concept [25-271. The idea here is to obtain high-P equilibria through a sequence of equilibria produced by increasing the pressure in an arbitrary way while holding the safety factor 9 fixed. The conclusions with respect to the optimum plasma shape are similar to those reached in this paper: D-shaped cross sections with mild elongation being most favourable. Such shapes also seem to be optimal with respect to stability [28]. Finally, we remark that screw-pinch experiments are of interest in relation to the FCI concept, since in screw-pinches high-0 equilibria are obtained with a given (flat) q-profile. In this paper the limit of a small inverse aspect ratio has been considered. Rebhan and Salat [29] investigated plasma-shape optimization in the framework of the surface-current model for general values of the aspect ratio. However, their studies are restricted to the low-beta regime and do not take into account the effects of a conducting wall. Surrounding the plasma by a layer of force-free currents, as occurs in screw-pinch discharges and as theoretically analyzed in this paper, enhances the coupling between the shape of the wall and that of the plasma. In this way the cross-sectional shaping of the plasma is facilitated. The resulting plasma shape is closer to the aforementioned limiting configuration, as compared to the case where the plasma is surrounded by a vacuum region. For any wall shape, force-free currents considerably increase the poloidal beta. It has been demonstrated that the effects of force-free currents on the free-boundary equilibrium and the poloidal beta are quantitatively related to the ratio of the total toroidal currents flowing in the force-free outer region and the plasma.

Acknowledgements The author is indebted to Dr. Ir. J.P. Goedbloed, useful discussions in the course of this investigation.

Prof. Dr. F. Engelmann,

and Dr. J. Rem for many

J.M. Akketmam

1 Free-boundary

rokamak equilibtia

131

This work was performed as part of the research programme of the association agreement of Euratom and the “Stichting voor Fundamenteel Onderzoek der Materie” (FOM) with financial support from the “Nederlandse Organisatie voor Zuiver-Wetenschappelijk Onderzoek” (ZWO) and EURATOM.

References [l] [2] [3] [4] [5] [6]

[7] [8] [9] [lo] [ll] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28]

[29]

R. Gajewski, Phys. Fluids 15 (1972) 70. W. Kemer, D. Pftrsch, and H. Tasso, Nucl. Fusion 12 (1972) 433. J.P. Freidberg and F.A. Haas, Phys. Fluids 16 (1973) 1909; 17 (1974) 440. J.P. Freidberg and W. Grossmann, Phys. Fluids 18 (1975) 1494. D.A. D’Ippolito, J.P. Freidberg, J.P. Goedbloed and J. Rem, Phys. Fluids 21 (1978) 1600. R.M.O. Galvao, J.P. Goedbloed, J. Rem, J.M. Akkermans, C. Bobeldijk, E.J.M. van Heesch, J.A. Hoekzema, A.F.G. van der Meer, D. Oepts and A.A.M. Omens, in Plasma Physics and Controlled Nuclear Fusion Research 1980 (IAEA, Vienna, 1981) Vol. II, p. 325. J. Rem, Rijnhuizen Report RR 81-130 (1981). H. Grad, A. Kadish and D.C. Stevens, Commun. Pure Appl. Math. 27 (1974) 39. J.A. Shercliff, J. Plasma Phys. 21 (1979) 347. P.N. Vabishchevich and L.M. Degtyarev, Dokl. Akad. Nauk 247 (1979) 1342 [Sov. Phys. Dokl. 24 (1979) 6121. D.C. Stevens, Phys. Fluids 17 (1974) 222. Y. Suzuki, Nucl. Fusion 14 (1974) 345. Qing Cheng-rui and Zhou Yu-mei, Acta Phys. Sin. 29 (1980) 106 [Chin. Phys. 1 (1981) 111. B. Steffen, J. Comput. Phys. 27 (1978) 180. J.P. Goedbloed, Phys. Fluids 25 (1982) 852. F. Bauer, 0. Betancourt and P. Garabedian, A Computational Method in Plasma Physics (Springer Verlag, New York, 1978). K. Lackner, Cornput. Phys. Commun. 12 (1976) 33. P. Henrici, SIAM Rev. 21 (1979) 481. J.P. Goedbloed, Comput. Phys. Commun. 24 (1981) 311. J.A. Hoekzema, in Pulsed High Beta Plasmas, D.E. Evans, ed. (Pergamon Press, Oxford, 1976) p. 535; see also: Rijnhuizen Report RR 78-112 (1978). L.S. Solovev, V.D. Shafranov and E.I. Yurchenko, in Plasma Physics and Controlled Nuclear Fusion Research, I%9 (IAEA, Vienna, 1969) Vol. I, p. 173. G. Lava], H. Luc, E.K. Maschke, C. Mercier and R. Pellat, in Plasma Physics and Controlled Nuclear Fusion Research 1971 (IAEA, Vienna, 1971) Vol. II, p. 507. J. Callen, B. Coppi, R. Dagazian, R. Gajewski and D. Sigmar, in Plasma Physics and Controlled Nuclear Fusion Research 1971 (IAEA, Vienna, 1971) Vol. II, p. 451. L.A. Artsimovich and V.D. Shafranov, ZhETF Pis. Red. 15 (1972) 72 [JETP Letters 15 (1972) 511. J.F. Clarke and D.J. Sigmar, Phys. Rev. Lett. 38 (1977) 70. R.A. Dory and Y.-K.M. Peng, Nucl. Fusion 17 (1977) 21. G. Bateman, MHD Instabilities, (MIT Press, Cambridge, Mass., 1978) ch. 8. A.M.M. Todd, M.S. Chance, J.M. Greene, R.C. Grimm, J.L. Johnson and J. Manickam, Phys. Rev. Lett. 38 (1977) 826; G. Bateman and Y.-K.M. Peng, Phys. Rev. Lett. 38 (1977) 829; L.A. Charlton, R.A. Dory, Y.-K.M. Peng, D.J. Strickler, S.J. Lynch, D.K. Lee, R. Gruber and F. Troyon, Phys. Rev. Lett. 43 (1979) 1395. E. Rebhan and A. Salat, Nucl. Fusion 18 (1978) 1639 and 20 (1980) 839.