J. Plasma
J.
Plasma Fusion
Fusion Res.
Res. SERIES,SERIES,Vol.
Vol.6 6(2004)(2004)152–155
000–000

Analysis of the MHD Instability Driving Mechanism
in 3D Heliotron and Quasiaxisymmetric Systems
COOPER W. Anthony, NARUSHIMA Yoshiro1, WATANABE Kiyomasa1, YAMAZAKI Kozo1,
SUZUKI Chihiro1 and OKAMURA Shoichi1
Centre de Recherches en Physique des Plasmas, EPFL, Lausanne, Switzerland
1
National Institute for Fusion Science, Toki 509-5292, Japan(Received: 9 December 2003 / Accepted: 22 March 2004)

Abstract
A module for the three dimensional (3D) ideal magnetohydrodynamic (MHD) code TERPSICHORE has beenformulated and implemented to investigate the driving and stabilising mechanisms associated with global and local
ideal MHD modes in 3D magnetic confinement systems. The energy principle that describes the MHD stability
behaviour is expressed as δ WP = (1/2) ∫ ∫ ∫ d 3x (C 2 – D|ξs|2 ), where C 2 is the stabilising element due to the compression
and bending of magnetic field lines, whereas the instability driving term is D|ξs|2 which can be separated into twocomponents. One of these components is proportional to the pressure gradient responsible for local ballooning andMercier modes. The second component is proportional to the parallel current density j·B/B2. This term is responsible
for global internal and external kink modes. However, in finite β plasmas, the surface varying component of j·B/B2 is
proportional to the pressure gradient. Although this could complicate the distinction between ballooning and surface
varying parallel current driven kinks, it is found that the perturbed energy structures each generate are different incharacter. In addition, the average and surface varying contributions of j·B/B2 to δ WP can be evaluated separately.
Applications to a currentless 10-period Heliotron demonstrates that the ballooning-interchange mechanism is dominant
in destabilising ideal MHD modes. In a 2-period quasiaxisymmetric device with finite bootstrap current, the kinkmode mechanism driven by the surface varying component of j·B/B2 is dominant near marginal stability.
Keywords:
ideal MHD energy principle, stellarators, quasi-axisymmetry (QAS), instability drive, Pfirsch-Schlüter current,
3D equilibrium and stability code

1. Introduction

currents. A rough identification of Pfirsch-Schlüter driven kinks
and kink modes driven by BC, OH or ECCD is undertaken
by suppressing the contribution of one or the other in δ WJ.
Applications to a 10-period currentless Heliotron [1] and a 2period quasiaxisymmetric (QAS) [2] stellarator reactor withfinite bootstrap-like current are investigated.

Ideal magnetohydrodynamic (MHD) instabilities impose
limits on the β value that can be achieved in magneticallyconfined plasmas. The bending and compression of field lines
balance the interaction of the pressure gradient with themagnetic field line curvature that drives ballooninginterchange modes and the parallel current density that drives
internal and external kink modes. A diagnostic code isdeveloped to analyse the driving and stabilising terms of theinternal plasma potential energy δ WP. The driving terms canbe separated into pressure gradient ballooning-interchangeδ WD and parallel current density δ WJ contributions. The
Pfirsch-Schlüter current is related to the flux surface varying
component of the parallel current density and is proportional
to the pressure gradient. A finite flux surface average componentof the parallel current density is a consequence of bootstrap(BC), ohmic (OH) or electron cyclotron current drive (ECCD)

2. The components of the energy principle
The incompressible ideal MHD stability problem isreduced to the solution of the equation δ WP + δ WV – ω 2δ WK= 0, where δ WP is the internal plasma potential energy, δ WV
is the energy in the vacuum region surrounding the plasma,–ω 2δ WK is the kinetic energy and ω 2 is the eigenvalue. InBoozer magnetic coordinates [3], expressions for δ WP, δ WV
and δ WK described in Ref. [4] have been implemented in the3D TERPSICHORE stability code [5]. The driving terms ofδ WP can be separated into a pressure gradient p′ component©2004 by The Japan Society of Plasma
Science and Nuclear Fusion Research

Corresponding author’s e-mail: wilfred.cooper@epfl.ch1152

Cooper W.A. et al., Analysis of the MHD Instability Driving Mechanism in 3D Heliotron and Quasiaxisymmetric Systems12δ WD = −×+∫ ∫2πLs10ds0m/n = 3/2 structure. As anticipated, the stabilising δ WC 2 term
is positive everywhere. We have also verified that thereconstructed δ WP profile reproduces exactly the profile
calculated spectrally in the TERPSICHORE code. An analysisof the full 3D distribution of the different components of theEnergy Principle is also instructive. We present the2π∫ dθ p′(s)dφ0g p′(s ) ∂ g J ( s )ψ ′′( s ) − I ( s )Φ ′′( s )s 2+−(ξ )∂sB2
B22π2πL1 1ds s dφdθ p′(s ) Bs2 g B⋅∇(ξ s )2002 0B∫ ∫∫and a parallel current density j·B/B2 component

δ WJ = −∫ ∫δW 2Cg B2 ⎛ j ⋅ B ⎞s 2
⎟ + ψ ′( s )Φ ′′( s ) − Φ ′( s )ψ ′′( s ) (ξ )⎜| ∇s |2 ⎝ B 2 ⎠Cg B ⋅∇(ξ s )20

δ WJ
δ WD

The stabilising contribution then is extracted by evaluatingδ WC2 = δ WP – δ WD – δ WJ and must be verified to be positivedefinite everywhere. The flux surface averaged component ofj·B/B2 is given, as calculated in a Boozer magnetic coordinate
reference frame by the expression [J(s)I′(s) – I(s)J′(s)]/[ψ ′(s)J(s) – Φ ′(s)I(s)] = ( j·B/B2)00, where J(s) and I(s) are
the toroidal and poloidal current flux functions, respectivelyand Φ (s) and ψ (s) are the toroidal and poloidal magnetic flux
functions, respectively. Primes (′ ) indicate derivatives withrespect to the radial variable s. By suppressing the flux
surface varying component of j·B/B2 associated with thePfirsch-Schlüter current, we can estimate the impact ofbootstrap, ohmic or ECCD currents to the kink driving termby evaluating

δ WJ 0 = −PP∫

I (s ) gsθ
J ( s ) gs φ+2| ∇s | g | ∇s |2 g

δWJ

2πLsD∫ ∫×4

∫

2π
1 1
dsdφdθ j ⋅ 2B+
002 0B×x 10

2π
2π
1 1Lj⋅Bds s dφdθ 200B
2 0

δW ,δW ,δW ,δW 2

×550Fig. 1

0.1

0.2

0.3

0.4

0.5
S

0.6

0.7

0.8

0.9

1The stabilising (δWC ), the internal plasma potential
energy ( δW P), the kink (δ WJ) and the ballooninginterchange (δWD) profiles in a currentless 10-period
Heliotron. The curves labelled from top to bottom areδWC , δWP , δWJ and δWD, respectively.
222π
2π1 1Ls⎛ j⋅B⎞
dsdφdθ ⎜ 2 ⎟0002⎝ B ⎠ 00

∫ ∫∫g B2 ⎛ j ⋅ B ⎞s 2
⎟ + ψ ′( s )Φ ′′( s ) − Φ ′( s )ψ ′′( s ) (ξ )⎜| ∇s |2 ⎝ B 2 ⎠ 00

2π2π1Lj⋅Bdθ ⎛⎜ 2 ⎞⎟+ 1 ds s dφ002 0⎝ B ⎠ 00

∫ ∫×∫I ( s ) gsθ
J ( s ) gs φ+2| ∇s | g | ∇s |2 gg B ⋅ ∇(ξ s )2 .

The difference between δ WJ and δ WJ0 yields an estimate of
the impact of the Pfirsch-Schlüter current as a driving
mechanism for kink modes.

3. Application to currentless 10-period
Heliotron
We analyse the different contributions to the Energy
Principle in a current-free 10-period Heliotron [6]. The
profiles for δ WP, δ WC 2, δ WD and δ WJ are shown in Fig. 1 fora case with β ~– 2%. The ballooning-interchange drivingmechanism δ WD dominates over the parallel current driving
mechanism δ WJ for the most unstable mode of the system, aFig. 2a The structure of the integrands of the functionals δWD(top) and δWJ (bottom) for a 10-period currentless
Heliotron on a toroidal cross section at the beginningof a period (φ = 0). The scale for δWD is –2.6 × 10–3 <δWD < 7 × 10–4 and the scale for δWJ is –1.4 × 10–3 < δWJ< 7 × 10–4.

2153

Cooper W.A. et al., Analysis of the MHD Instability Driving Mechanism in 3D Heliotron and Quasiaxisymmetric Systems

distributions of integrands of the functionals δ WD and δ WJ inFigs. 2a–c on cross sections at the beginning, at one quarterand at midperiod. The most negative regions of δ WD
concentrate at the tips of the elliptic cross sections
corresponding to the locations where the magnetic field line
curvature is most destabilising. This characterises the mode
structure as ballooning. The most negative regions of δ WJ
concentrate at the inner edge of the torus and are a factor 2smaller that δ WD. Though the Pfirsch-Schlüter current isformally proportional to p′, the δ WJ structure it drives does
not display ballooning features. The observations madesuggest where to physically position fluctuation detectors toidentify the type of mode. In practice, heliotron systems like
LHD [1] display finite currents due to neutral beam and1.5x 105

δ WC2

δ WJ , δ WP , δ WJ0 , δ WD , δ WC21Fig. 2b The structure of the integrands of the functionals δWD(top) and δWJ (bottom) for a 10-period currentless
Heliotron on a toroidal cross section at one quarter ofa period (φ = π /20). The scale for δWD is –2.6 × 10–3 <δWD < 7 × 10–4 and the scale for δWJ is –1.4 × 10–3 < δWJ< 7 × 10–4.

0.5
δWD0

δ WJ0
δ WP

0. 5

1

1. 5
δW

J

20

Fig. 3

0.1

0.2

0.3

0.4

0.5
S

0.6

0.7

0.8

0.9

1

The profile of the different components that contributeto the energy principle δWP (the curves shown from top
to bottom are δWC , δWD , δWJ 0 and δWJ, respectively) ina 2-period QAS reactor with finite bootstrap-likecurrent.
2x 10

5I=3.78MA
δ WD
0I=5.39MA
δ WJ , δ WJ0 , δ WD

δ WJ0δ WJ

100Fig. 2c The structure of the integrands of the functionals δWD(top) and δWJ (bottom) for a 10-period currentless
Heliotron on a toroidal cross section at midperiod (φ =π /10). The scale for δWD is –2.6 × 10–3 < δWD < 7 × 10–4
and the scale for δWJ is –1.4 × 10–3 < δWJ < 7 × 10–4.

Fig. 4

3154

0.1

0.2

0.3

0.4

0.5
S

0.6

0.7

0.8

0.9

1The profiles of δWD , δWJ and δWJ0 in a 2-period QAS
stellarator reactor with toroidal currents of 3.78 MA and5.39 MA. The curves shown from top to bottom
correspond to 1) (δWD , 5.39 MA), 2) (δWD, 3.78 MA), 3)(δWJ0 , 3.78 MA), 4) (δWJ , 3.78 MA), 5) (δWJ0 , 5.39 MA)
and 6) (δWJ, 5.39 MA), respectively.

Cooper W.A. et al., Analysis of the MHD Instability Driving Mechanism in 3D Heliotron and Quasiaxisymmetric Systems

bootstrap effects. The impact of these currents will be
addressed in future work.

bootstrap, ohmic or ECCD currents. Test applications havebeen performed on a 10-period current-free Heliotron and on
a 2-period QAS stellarator reactor with a fixed modelbootstrap current profile and current magnitude that isadjusted arbitrarily. In the currentless Heliotron, theballooning-interchange driving mechanism dominates that of
the Pfirsch-Schlüter current. Furthermore, the ballooninginterchange destabilising energy concentrates in the region of
weak magnetic field line curvature at the tips of the elliptic
flux surfaces. The energy in the sub dominant Pfirsch-Schlüter
current contribution localises at the inside edge of the torus.
In the finite current QAS system investigated, the ballooninginterchange contribution is globally stabilising. The PfirschSchlüter current contribution constitutes the main destabilisingmechanism near marginal stability. For higher toroidal
currents, the bootstrap current contribution to destabilise thekink mode becomes more important and can be quantified
through the evaluation of the functional δ WJ0 . The different
locations within the plasma where the unstable contributionof the kink driving term δ WJ and the ballooning-interchange
driving term δ WD concentrate can provide a very usefulindicator of the type of unstable mode that can be detectedfrom an experimental viewpoint.

4. Application to 2-period QAS withfinitetoroidal current
We investigate a sequence of equilibria at β ~– 4.25%computed with VMEC [7] with a fixed hollow bootstrap-liketoroidal current profile that model a 2-period QAS stellarator
reactor in which the magnitude of the current is varied
arbitrarily. The formulation is valid also for peaked currents
from neutral beam injection or ECCD in as much as thecurrent profile and magnitude are prescribed as input to theVMEC code which computes selfconsistent equilibria. This
procedure, in turn, affects the δ WJ and δ WJ0 distributions
when the stability is investigated using TERPSICHORE.
However, we have not examined peaked current profiles inthis paper. The bootstrap current causes the rotationaltransform profile to cross the critical ι = 1/2 resonant surface
whent the toroidal current approaches 3.8 MA and an external
m/n = 2/1 kink becomes destabilised. The profiles of δ WP,δ WC 2, δ WD, δ WJ and δ WJ0 are shown in Fig. 3 for a case with
3.78 MA toroidal current. The configuration is weakly
unstable and driven primarily by the Pfirsch-Schlüter current
because δ WJ is dominantly negative and an order ofmagnitude larger than δ WJ0. Surprisingly, on average theballooning-interchange drive δ WD is positive, thus stabilising.
Locally though, it is negative in regions of destabilisingmagnetic field line curvature. In Fig. 4, we display the profilesof δ WD, δ WJ and δ WJ0 for the weakly unstable 3.78 MA
current case and for a more strongly unstable 5.39 MA current
case. It is seen that for higher currents that the contributionof the BC to the kink driving mechanism becomes more
important and almost comparable in magnitude to that of thePfirsch-Schlüter current.

Acknowledgements
The authors are thankful to Dr. S.P. Hirshman for theuse of VMEC. This work was supported by the Fonds
National Suisse pour la Recherche Scientifique and byEuratom. The computational results were performed on theNEC-SX5 platform at the Centro Svizzero di Calcolo
Scientifico in Manno, Switzerland.

References[1] A. Iiyoshi, A. Fujiwara, O. Motojima, N. Ohyabu andK. Yamazaki, Fusion Technol. 17, 169 (1990).
[2] S. Okamura et al., Nucl. Fusion 41, 1865 (2001).
[3] A.H. Boozer, Phys. Fluids, 26 496 (1983).
[4] W.A. Cooper, Plasma Phys. Control. Fusion 34 1011
(1992).
[5] D.V. Anderson, W.A. Cooper, R. Gruber, S. Merazzi and
U. Schwenn, Scientific Computing on Supercomputers
II, Plenum Press, New York and London 159 (1990).
[6] Y. Nakamura, T. Matsumoto, M. Wakatani et al., J.
Comput. Phys. 128 43 (1996).
[7] S.P. Hirshman and O. Betancourt, J. Comput. Phys. 96
99 (1991).

5. Summary and conclusions
A diagnostic routine has been developed to analyse thestabilising and destabilising contributions to the ideal MHD
Enery Principle using the eigenfunctions computed with the3D TERPSICHORE code. The driving term can be separated
into contributions from the parallel current density, which isresponsible for kink modes, and the pressure gradient, whichis associated with ballooning and interchange modes. The
parallel current density contribution can be roughly separated
into a Pfirsch-Schlüter term (which is proportional to thepressure gradient but does not generate a ballooning structure)
and a flux surface average term associated with finite4155