Scaling of up-down asymmetric turbulent momentum flux with poloidal shaping mode number in tokamaks

1.2 MB·

by  Justin Ball

Scaling of up-down asymmetric turbulent
momentum flux with poloidal shaping mode number
6
1 in tokamaks
0
2

r
p Justin Ball and Felix I. Parra
A
Rudolf Peierls Centre for Theoretical Physics, University of Oxford, Oxford OX1

4 3NP, United Kingdom
1
Culham Centre for Fusion Energy, Culham Science Centre, Abingdon OX14 3DB,

United Kingdom
]
h
p E-mail: [email protected]
-
m
Abstract.
s
a Breakingtheup-downsymmetryoftokamaksremovesaconstraintlimitingintrinsic
l
p momentum transport, and hence toroidal rotation, to be small. Using gyrokinetic
. theory, we study the effect of different up-down asymmetric flux surface shapes on
s
c the turbulent transport of momentum. This is done by perturbatively expanding the
i
s gyrokinetic equation in large flux surface shaping mode number. It is found that the
y
momentumfluxgeneratedbyshapingthatlacksmirrorsymmetry(whichisnecessarily
h
p up-down asymmetric) has a power law scaling with the shaping mode number.
[ However,themomentumfluxgeneratedbymirrorsymmetricfluxsurfaceshaping(even

3 if it is up-down asymmetric) decays exponentially with large shaping mode number.
v These scalings are consistent with nonlinear local gyrokinetic simulations and indicate
0
that low mode number shaping effects (e.g. elongation, triangularity) are optimal for
6
creating rotation. Additionally it suggests that breaking the mirror symmetry of flux
5
3 surfaces may generate significantly more toroidal rotation.
0
.
1
0
PACS numbers: 52.25.Fi, 52.30.Gz, 52.35.Ra, 52.55.Fa, 52.65.Tt
6
1
:
v 1. Introduction
i
X
r Bulk toroidal rotation has been shown to be beneficial for plasma performance in
a
tokamaks. It can stabilize the resistive wall mode, which allows violation of the
Troyon limit [1, 2, 3, 4], a fundamental constraint on how much plasma pressure
can be confined with a given magnetic field [5]. Exceeding the Troyon limit directly
improves the economic viability of a tokamak power plant [6, 7]. Furthermore, a strong
gradient in toroidal rotation can reduce energy transport by shearing turbulent eddies
[8, 9, 10, 11, 12].
There are several mechanisms that currently generate rotation in tokamaks. Beams
of neutral particles [13] and radio frequency waves [14] are commonly used to heat the
Scaling of momentum flux with poloidal shaping mode number 2
plasma, but can also drive toroidal momentum. This externally injected momentum is
significant in current experiments, but is expected to diminish in larger devices [15].
Self-generated momentum transport, driven by plasma turbulence, is observed,
even in the absence of external injection. It is called “intrinsic” momentum transport,
but it is generally weak, creating rotation less than a tenth of plasma sound speed
[16, 17]. Recently however, the strength of intrinsic rotation was explained through
a symmetry of the gyrokinetic model [18, 19], a set of equations that are believed to
govern turbulence in the core of tokamaks [20]. This symmetry constrains the turbulent
transport of momentum to be on the order of ρ ρ /a 1, the ratio of the ion
∗ i
≡ (cid:28)
gyroradius to the tokamak minor radius [21, 22, 23]. Reference [21] shows that, in the
absenceofpreexistingrotation, thisconstraintholdsaslongasthetokamakfluxsurfaces
are up-down symmetric (i.e. have mirror symmetry about the midplane). Further
investigation suggests that breaking the up-down symmetry of the magnetic geometry
is a practical means to generate significant plasma rotation [24, 25, 26, 27]. Hence it
appears that up-down asymmetry is the most promising method to generate significant
intrinsic momentum transport in a reactor-scale, initially stationary plasma [28].
Subsequent work has demonstrated a new symmetry of the gyrokinetic model [29].
This symmetry means that poloidally translating all high order flux surface shaping
effectsbyasingletiltanglehaslittleeffectonthetransportpropertiesoftheequilibrium.
Thishasimportantconsequencesforintrinsicrotationgeneratedbyup-downasymmetry
because it creates a distinction between mirror symmetric tokamaks and non-mirror
symmetric tokamaks, which we will explore in depth. “Mirror symmetric” refers to
tokamaks that have flux surfaces with reflectional symmetry about some line in the
poloidal plane. When the line of symmetry is the midplane the mirror symmetric
tokamak can also be said to be “up-down symmetric.” “Non-mirror symmetric”
tokamaks have flux surfaces that do not have reflectional symmetry about any line
in the poloidal plane.
In this work we compare the intrinsic momentum transport in magnetic geometries
with different up-down asymmetric shaping effects. In section 2, we present the
electrostatic gyrokinetic model and give a generalized version of the local Miller
equilibrium, appropriate for specifying unusual up-down asymmetric configurations.
Then, we expand the gyrokinetic equation order-by-order in large shaping mode number
to compare the momentum flux generated by different types of flux surface shaping.
In doing so we will present two distinct arguments concerning the momentum flux
generatedbythelocalequilibrium. First,insection2.3,wecalculatehowthemomentum
flux scales with the shaping effect mode number given a specific set of simplified, non-
mirror symmetric geometries. This is designed to give a concrete illustration of the more
abstract and general scaling argument for non-mirror symmetric geometries presented
in section 2.4. Second, in section 2.5, we apply the symmetry presented in reference
[29] to establish the scaling of momentum flux with shaping mode number in mirror
symmetric(butstillup-downasymmetric)configurations. Theninsection3wecompare
theanalyticresultsofsection2tononlinearlocalgyrokineticsimulations. Lastly, section
Scaling of momentum flux with poloidal shaping mode number 3
4 gives a summary, a broad interpretation of the analytic scalings, and some concluding
remarks.
2. Analytic gyrokinetic analysis
Gyrokinetics [30, 31, 32, 33, 34, 35, 36, 37, 38, 39] is a theoretical framework used
to study plasma behavior with perpendicular wavenumbers comparable to the ion
gyroradius(k ρ 1)andtimescalesmuchslowerthantheparticlecyclotronfrequencies
⊥ i
∼
(ω Ω Ω ). Fundamentally, gyrokinetics relies on an expansion in ρ ρ /a 1,
i e ∗ i
(cid:28) (cid:28) ≡ (cid:28)
where ρ is the ion gyroradius and a is the tokamak minor radius. These particular
i
scales have been experimentally shown to be appropriate for modeling turbulence
[20]. In deriving the gyrokinetic equations, we expand the distribution function,
f = f +f +..., and assume the perturbation is small compared to the background
s s0 s1
(f ρ f ) [40]. Additionally, for tokamak plasmas, axisymmetry implies radially
s1 ∗ s0
∼
confined orbits and the transport timescale usually exceeds the collisional timescale. As
aresult, thelowestorderdistributionfunctionisassumedtobeMaxwellian(f = F ).
s0 Ms
Here
(cid:18) m (cid:19)3/2 (cid:18) m w2(cid:19)
s s
F n exp (1)
Ms s
≡ 2πT − 2T
s s
is the Maxwellian distribution function for species s, n is the particle density, m is
s s
the particle mass, T is the temperature, and w(cid:126) is the velocity in the frame rotating
s
with the plasma. In this work we will choose to neglect both electromagnetic effects (for
simplicity) and pre-existing rotation (because we are interested in generating rotation
in a stationary plasma).
Given these assumptions, we can change the coordinates of the Fokker-Plank and
quasineutrality equations from real-space coordinates to the guiding center position, i.e.
theaveragepositionoftheparticleasitspiralsaroundamagneticfieldline. Thenwecan
average over the gyrophase angle ϕ, i.e. the angle that determines the particle location
on its circular motion perpendicular to the magnetic field. This gives the two governing
equations of electrostatic gyrokinetics: the gyrokinetic equation and a modified version
of the quasineutrality equation.
The electrostatic gyrokinetic equation, in the absence of rotation and collisions, can
be Fourier-analyzed in the directions perpendicular to the magnetic field and written as
[21]
(cid:12) (cid:20)(cid:18) (cid:19) (cid:18) (cid:19)(cid:21)
∂hs +w ˆb (cid:126)θ ∂hs(cid:12)(cid:12) +i w2 + B µ (k v +k v ) w2 k µ0 dp h
∂t || ·∇ ∂θ (cid:12) || m ψ dsψ α dsα − || αΩ Bdψ s
w ,µ s s
||
(cid:12)
∂hs(cid:12) ZseFMs∂ φ ϕ
+a (cid:12) + φ ,h (cid:104) (cid:105) (2)
s|| ϕ s
∂w (cid:12) {(cid:104) (cid:105) }− T ∂t
|| θ,µ s
(cid:20) 1 dn (cid:18)m w2 3(cid:19) 1 dT (cid:21)
s s s
+ik φ F + = 0
α ϕ Ms
(cid:104) (cid:105) n dψ 2T − 2 T dψ
s s s
Scaling of momentum flux with poloidal shaping mode number 4
(cid:0) (cid:1)
in the k ,k ,θ,w ,µ,ϕ,t coordinate system. Here h is the Fourier-analyzed
ψ α || s
nonadiabatic portion of the distribution function, t is the time, w is the component of
||
ˆ (cid:126)
the velocity parallel to b B/B, the magnetic field unit vector, θ is the usual cylindrical
≡
poloidal angle shown in figure 1, µ m w2/2B is the magnetic moment, k is the radial
≡ s ⊥ ψ
wavenumber, k is the wavenumber within the flux surface and perpendicular to the
α
magnetic field, µ is the permeability of free space, Ω Z eB/m is the gyrofrequency,
0 s s s
≡
p is the plasma pressure, ψ is the poloidal magnetic flux, φ is the Fourier-analyzed
electrostatic potential, Z is the particle charge number, and e is the electric charge of
s
the proton. The magnetic drift coefficients are given by
(cid:12)
I(ψ)ˆ (cid:126) ∂B(cid:12)
v b θ (cid:12) (3)
dsψ
≡− Ω B ·∇ ∂θ (cid:12)
s ψ
 (cid:16) (cid:17)
(cid:12) (cid:12) ˆb (cid:126)θ (cid:126)α
1 ∂B(cid:12) ∂B(cid:12) · ∇ ×∇
vdsα  (cid:12) (cid:12) , (4)
≡− Ω ∂ψ(cid:12) − ∂θ (cid:12) B
s θ ψ
where I(ψ) = RB is the toroidal field flux function,
ζ
(cid:90) θ
α ζ dθ(cid:48)A (ψ,θ(cid:48)) (5)
α
≡ −
θα(ψ)
is a coordinate that selects a particular field line from a given flux surface,
I(ψ)
A (ψ,θ) (6)
α ≡ R2B(cid:126) (cid:126)θ
·∇
is the integrand in the definition of α, and θ (ψ) is a free function that determines the
α
field line selected by α = 0 on each flux surface. The parallel acceleration is given by
(cid:12)
µ ˆ (cid:126) ∂B(cid:12)
a b θ (cid:12) . (7)
s||
≡ −m ·∇ ∂θ (cid:12)
s ψ
The nonlinear term is
(cid:88) (cid:0) (cid:1) (cid:0) (cid:1) (cid:0) (cid:1)
φ ,h = k(cid:48) k k k(cid:48) φ k(cid:48) ,k(cid:48) h k k(cid:48) ,k k(cid:48) (8)
{(cid:104) (cid:105)ϕ s} ψ α − ψ α (cid:104) (cid:105)ϕ ψ α s ψ − ψ α − α
k(cid:48) ,k(cid:48)
ψ α
and the gyroaverage is given by
(cid:18) (cid:19)
k √2µB
⊥
... = J (...), (9)
(cid:104) (cid:105)ϕ 0 Ω √m
s s
where
(cid:114)
(cid:12) (cid:12)2 (cid:12) (cid:12)2
k = k2 (cid:12)(cid:126)ψ(cid:12) +2k k (cid:126)ψ (cid:126)α+k2 (cid:12)(cid:126)α(cid:12) (10)
⊥ ψ(cid:12)∇ (cid:12) ψ α∇ ·∇ α(cid:12)∇ (cid:12)
is the perpendicular wavenumber and J (...) is the nth order Bessel function of the
n
first kind.
Scaling of momentum flux with poloidal shaping mode number 5
The quasineutrality equation can be Fourier-analyzed and written as [21]
(cid:32) (cid:33)−1
(cid:88) Z2e2n (cid:88) 2πZ eB (cid:90)
φ = s s s dw dµ h . (11)
|| s ϕ
T m (cid:104) (cid:105)
s s
s s
Solving the gyrokinetic and quasineutrality equations for h and φ allows us to calculate
s
the electrostatic, turbulent flux of toroidal angular momentum according to [21]
(cid:42) (cid:42) (cid:43) (cid:43)
(cid:28)(cid:90) (cid:29)
(cid:16) (cid:17)
Π R d3wh˘ m R(w(cid:126) eˆ ) δE(cid:126) eˆ (12)
s s s ζ ζ
≡− · ·
ψ
∆ψ ∆t
(cid:42)
4π2i (cid:88) (cid:73) (cid:90)
= k dθJBφ(k ,k ) dw dµ h ( k , k ) (13)
V(cid:48) α ψ α || s − ψ − α
kψ,kα
(cid:18)I i kψ µB2J (k ρ )(cid:19)(cid:29)
1 ⊥ s
w J (k ρ )+ ,
|| 0 ⊥ s
× B Ω B m k ρ
s s ⊥ s ∆t
˘
where h is the non-adiabatic perturbed distribution function (before Fourier analysis),
s
δE(cid:126) is the turbulent electric field, ... (2π/V(cid:48))(cid:72)2πdθJ (...) is the flux surface
(cid:104) (cid:105)ψ ≡ 0
average, V(cid:48) 2π(cid:72)2πdθJ,
≡ 0
(cid:12) (cid:16) (cid:17)(cid:12)−1 (cid:16) (cid:17)−1
J (cid:12)(cid:126)ψ (cid:126)θ (cid:126)ζ (cid:12) = B(cid:126) (cid:126)θ (14)
(cid:12) (cid:12)
≡ ∇ · ∇ ×∇ ·∇
(cid:82)
is the Jacobian, ... ∆ψ−1 (...) is the coarse-grain average over a radial
(cid:104) (cid:105)∆ψ ≡ ∆ψ
distance ∆ψ (which is larger than the scale of the turbulence, but smaller than the
(cid:82)
scale of the device), ... ∆t−1 (...) is the coarse-grain average over a time
(cid:104) (cid:105)∆t ≡ ∆t
∆t (which is longer than the turbulent decorrelation time), and kψ (cid:126)k (cid:126)ψ =
⊥
(cid:12) (cid:12)2 ≡ · ∇
k (cid:12)(cid:126)ψ(cid:12) +k (cid:126)ψ (cid:126)α.
ψ(cid:12) (cid:12) α
∇ ∇ ·∇
We note that the following eight coefficients contain all the information about the
(cid:12) (cid:12)2 (cid:12) (cid:12)2
flux surface geometry: ˆb (cid:126)θ, B, v , v , a , (cid:12)(cid:126)ψ(cid:12) , (cid:126)ψ (cid:126)α, and (cid:12)(cid:126)α(cid:12) . In an
dsψ dsα s|| (cid:12) (cid:12) (cid:12) (cid:12)
· ∇ ∇ ∇ · ∇ ∇
(cid:126) (cid:126)
up-down symmetric tokamak, the coefficients v , a , and ψ α are necessarily
dsψ s||
(cid:12) (cid:12)2 (cid:12) (cid:12)2 ∇ · ∇
odd in θ, while ˆb (cid:126)θ, B, v , (cid:12)(cid:126)ψ(cid:12) , and (cid:12)(cid:126)α(cid:12) are even. As shown in reference
dsα (cid:12) (cid:12) (cid:12) (cid:12)
· ∇ ∇ ∇
[21], the parity of the geometric coefficients in an up-down symmetric tokamak has
important consquences for overall symmetry properties of the gyrokinetic equations.
(cid:0) (cid:1) (cid:0) (cid:1)
The equations become invariant to the k ,k ,θ,w ,µ,t k ,k , θ, w ,µ,t
ψ α || ψ α ||
→ − − −
coordinate system transformation, which is not true in up-down asymmetric devices.
(cid:0) (cid:1)
This symmetry means that, given any solution h k ,k ,θ,w ,µ,t , we can construct
s ψ α ||
(cid:0) (cid:1)
a second solution h k ,k , θ, w ,µ,t that will also satisfy the gyrokinetic
s ψ α ||
− − − −
equations. From equation (13) we see that this second solution will have a momentum
flux that cancels that of the original. These two solutions are each valid for different
initial conditions, but since the tokamak is presumed to be chaotic, both solutions will
arise within a turbulent decorrelation time (statistically speaking). This demonstrates
Scaling of momentum flux with poloidal shaping mode number 6
that, in the gyrokinetic limit, the time-averaged momentum flux must be zero in an
up-down symmetric tokamak.
In subsection 2.1 we present a local MHD equilibrium specification that is
appropriate for flux surfaces with arbitrary shaping. Then in subsection 2.2 we briefly
preface the asymptotic expansion of the gyrokinetic model in large shaping mode
number. Insubsection2.3, wefirstcalculatethegeometriccoefficientsinthelargeaspect
ratiolimitfromtheMHDequilibriumforarealistic,butsimpleexamplegeometry. Using
thisgeometryweexpandthegyrokineticequationstodeterminehowthemomentumflux
scales with the mode number of the symmetry-breaking effect. This concrete example
serves to illustrate the derivation for a general geometry without expanding in aspect
ratio, which is detailed in subsection 2.4. Lastly, in subsection 2.5 we explain why
mirror symmetric geometries are a special case and should be expected to have weak
momentum transport.
2.1. Up-down asymmetric local Miller equilibrium
We will calculate the local value of the geometric coefficients that appear in the
gyrokinetic equation by using the local Miller geometry model [41]. The Miller
equilibrium model is a way of specifying the local tokamak equilibrium in the vicinity
of a single flux surface of interest. The local equilibrium is completely described by
the shape of the flux surface of interest (labeled by r ), how this shape changes with
ψ0
r (the minor radial coordinate), and four scalar quantities. Traditionally B (the
ψ 0
on-axis toroidal magnetic field), q (the safety factor), sˆ (r /q)dq/dr (the magnetic
ψ0 ψ
≡
shear), and dp/dr (the pressure gradient) are used. Typically, a combination of vertical
ψ
elongation and positive triangularity are used to specify the flux surface shape, but in
this work we will use a completely general flux surface shape specification (similar to
that presented in reference [29]).
Since we know that flux surfaces must be periodic in poloidal angle, we are free to
Fourier analyze and express them as an infinite series of shaping modes. We will choose
to specify the shape of the flux surface of interest in polar form (see figure 1) as
(cid:32) (cid:33)
(cid:88) ∆ 1
m
r (θ) =r(r ,θ) = r 1 − cos(m(θ+θ )) , (15)
0 ψ0 ψ0 tm
− ∆ +1
m
m
where m is the shaping mode number. Note that this is a completely general Fourier
decomposition. The strength of each shaping effect is set by the parameter ∆ . If
m
only one shaping effect is present then ∆ = b/a, where b and a are the maximum and
m
minimum distance of the flux surface from the magnetic axis respectively. When m = 2,
this definition reduces to the usual elongation (typically denoted by κ). The tilt angles,
θ , control the relative strength of the sine and cosine terms for every m. Lastly, the
tm
flux surface label r determines the constant Fourier term. Note the distinction between
ψ
a (the minimum distance of a flux surface from the magnetic axis) and r (a flux surface
ψ
label that, as we will see, is defined through equation (19)).
Scaling of momentum flux with poloidal shaping mode number 7
Z
rr ((θθ))
00
θθ
●● R
R
0
Figure 1. An example flux surface of interest, r (θ), needed by equation (19) for the
0
Miller local equilibrium model. The (R,Z,ζ) coordinate system is defined such that
the toroidal angle ζ and the plasma current are coming out of the page.
Differentiating equation (15) radially, we find
(cid:12)
∂r (cid:12) (cid:88)
(cid:12) =1 δ∆ cos(m(θ+δθ )), (16)
m tm
∂r (cid:12) −
ψ ψ0 m
where
(cid:115)
(cid:18)∆ 1 2r d∆ (cid:19)2 (cid:18) ∆ 1dθ (cid:19)2
m ψ0 m m tm
δ∆ − + + mr − (17)
m ≡ ∆ +1 (∆ +1)2 dr ψ0∆ +1 dr
m m ψ m ψ
(cid:18) (cid:30)(cid:18) (cid:19)(cid:19)
1 dθ 2r d∆
tm ψ0 m
δθ θ + arctan mr (∆ 1) (∆ 1)+ (18)
tm tm ψ0 m m
≡ m − dr − ∆ +1 dr
ψ m ψ
for each m. Note that for the local equilibrium all radially varying quantities are
evaluated at r = r (or equivalently ψ = ψ ), the flux surface of interest. The change
ψ ψ0 0
in the strength, d∆ /dr , and tilt, dθ /dr , of each mode would be determined by
m ψ tm ψ
the global MHD equilibrium. This is governed by the Grad-Shafranov equation [42] and
requires the entire radial current profile. In Appendix A we derive these quantities using
a constant current profile in the limits of large aspect ratio and weak shaping. However,
in the local Miller equilibrium model the radial variation of the flux surface shape is an
input used to construct the poloidal magnetic field. After calculating the poloidal field,
the Grad-Shafranov equation is used to calculate all higher order radial derivatives and
approximate the global equilibrium.
In summary, the flux surface geometry for the Miller local equilibrium model is
completely specified by equations (15), (16), and
(cid:12)
∂r (cid:12)
r(r ,θ) =r (θ)+ (cid:12) (r r ) (19)
ψ 0 ψ ψ0
∂r (cid:12) −
ψ ψ0
Scaling of momentum flux with poloidal shaping mode number 8
R(r ,θ) =R +r(r ,θ)cos(θ) (20)
ψ 0 ψ
Z(r ,θ) =r(r ,θ)sin(θ), (21)
ψ ψ
where r , ∆ , θ , δ∆ , δθ , R , q, sˆ, dp/dr , and B (the tokamak major radius)
ψ0 m tm m tm 0 ψ 0
are inputs. We note that if r (θ) = r ( θ+θ ) and ∂r/∂r = ∂r/∂r for
0 0 − 0 ψ|ψ0,θ ψ|ψ0,−θ+θ0
any θ , then the tokamak is mirror symmetric, otherwise it is non-mirror symmetric.
0
Similarly, if r (θ) = r ( θ) and ∂r/∂r = ∂r/∂r , then the tokamak is up-
0 0 − ψ|ψ0,θ ψ|ψ0,−θ
down symmetric (as well as mirror symmetric), otherwise it is up-down asymmetric.
The full calculation of all eight geometric coefficients is shown in Appendix B, but
for brevity here we will only calculate them to lowest order in (cid:15) a/R 1 (i.e. the
0
≡ (cid:28)
inverse aspect ratio). To lowest order in aspect ratio B B and a 0, so we can
0 ||s
(cid:12) (cid:12)2 → (cid:12) (cid:12)2 →
focus on the other six (ˆb (cid:126)θ, v , v , (cid:12)(cid:126)ψ(cid:12) , (cid:126)ψ (cid:126)α, and (cid:12)(cid:126)α(cid:12) ). In this limit the
dsψ dsα (cid:12) (cid:12) (cid:12) (cid:12)
·∇ ∇ ∇ ·∇ ∇
momentum flux, given by equation (13), becomes
(cid:73)
2πiR B (cid:88) (cid:16) (cid:17)−1
0 0 ˆ (cid:126)
Π = k dθ b θ (22)
s (cid:72) (cid:16)ˆ (cid:126) (cid:17)−1 α ·∇
dθ b θ kψ,kα
·∇
(cid:90)
dw dµ w J (k ρ )φ(k ,k )h ( k , k )
|| || 0 ⊥ s ψ α s ψ α
× − −
to lowest order in (cid:15) 1. For ease of notation we will not use q, sˆ, or dp/dr as inputs
ψ
(cid:28)
to the Miller local equilibrium model. Instead, we will choose to replace q by dψ/dr
ψ
(see equation (B.3)). Also, when we expand to lowest order in aspect ratio, we will find
that we can replace both dp/dr and sˆ (derived from dI/dψ) with
ψ
(cid:32) (cid:33)
(cid:18)dψ (cid:19)−1 (cid:18)dψ (cid:19)−1 dp dI
sˆ(cid:48) 2+r µ R2 +R B (23)
≡ ψ0 dr 0 0 dr dr 0 0dψ
ψ ψ ψ
(cid:18)dψ (cid:19)−1
=2 r µ j R ,
ψ0 0 ζ 0
− dr
ψ
where I(ψ) RB is the toroidal field flux function, j is the current density in the
ζ ζ
≡
toroidal direction, and R is the major radial coordinate. We can make this replacement
because the toroidal current, which appears on the right side of the Grad-Shafranov
equation (see equation (27)), is a flux function to lowest order in aspect ratio. We note
that if the flux surfaces are exactly circular and dp/dr = 0, then sˆ(cid:48) = sˆ.
ψ
Using our geometry specification given by equations (15) through (21) and
employing
∂(cid:126)r/∂u ∂(cid:126)r/∂u
(cid:126) 2 3
u = × , (24)
1
∇ ∂(cid:126)r/∂u (∂(cid:126)r/∂u ∂(cid:126)r/∂u )
1 2 3
· ×
where (u ,u ,u ) is a cyclic permutation of (r ,θ,ζ), we can directly calculate the
1 2 3 ψ
poloidal field,
dψ
(cid:126) (cid:126) (cid:126)
B = ζ r . (25)
p ψ
dr ∇ ×∇
ψ
Scaling of momentum flux with poloidal shaping mode number 9
ˆ (cid:126)
This allows us to calculate b θ,
·∇
m ∂R
sˆ (cid:126)
v = b θ , (26)
dsψ
Z e ·∇ ∂θ
s
(cid:12) (cid:12)2
and (cid:12)(cid:126)ψ(cid:12) = R2B2 to lowest order in aspect ratio. However, (cid:126)α contains second-
(cid:12)∇ (cid:12) p ∇
order radial derivatives, which are not specified. The Miller model determines them by
ensuring that the Grad-Shafranov equation [42],
(cid:32) (cid:33)
(cid:126)
ψ dp dI
R2(cid:126) ∇ = µ j R = µ R2 I , (27)
∇· R2 0 ζ − 0 dψ − dψ
is satisfied. With considerable work (shown in Appendix B), we can use the Grad-
Shafranov equation to calculate that
∂α ∂α
(cid:126) (cid:126) (cid:126)
α = ψ + θ, (28)
∇ ∂ψ∇ ∂θ∇
where
∂α = (cid:90) θ(cid:12)(cid:12)(cid:12) dθ(cid:48)∂Aα +A (ψ,θ ) dθα (29)
α α
∂ψ − (cid:12) ∂ψ dψ
θα ψ
= (cid:90) θ(cid:12)(cid:12)(cid:12) dθ(cid:48)(cid:18)∂Aα(cid:19) +(cid:20) B0 ∂lp(cid:126)ψ (cid:126)θ(cid:48)(cid:21)θ(cid:48)=θ (30)
− (cid:12) ∂ψ R3B3∂θ(cid:48)∇ ·∇
θα ψ orthog 0 p θ(cid:48)=θα
(cid:18) (cid:19)
B ∂l dθ
0 p α
+
R B ∂θ dψ
0 p θ=θα
and
∂α B ∂l
0 p
= A (ψ,θ) = (31)
α
∂θ − −R B ∂θ
0 p
to lowest order in aspect ratio. Here
(cid:18)∂A (cid:19) B ∂l (cid:18)dψ sˆ(cid:48) 2 (cid:19)
α 0 p
= − +2κ (32)
∂ψ R2B2∂θ(cid:48) dr r R B p
orthog 0 p ψ ψ0 0 p
is the part of ∂A /∂ψ that remains if the (r ,θ,ζ) coordinate system is orthogonal,
α ψ
(cid:16) (cid:17) (cid:126)ψ (cid:18)∂l (cid:19)−3(cid:18)∂R∂2Z ∂2R∂Z(cid:19)
ˆ (cid:126)ˆ p
κp ≡− bp ·∇bp · (cid:12)(cid:12)∇(cid:126)ψ(cid:12)(cid:12) = ∂θ ∂θ ∂θ2 − ∂θ2 ∂θ (33)
(cid:12) (cid:12)
∇
is the poloidal magnetic field curvature (defined such that the inwards normal direction
ˆ (cid:126)
is positive), b B /B is the poloidal field unit vector, and l is the poloidal arc length,
p p p p
≡
defined such that
(cid:115)
∂l (cid:114)∂(cid:126)r ∂(cid:126)r (cid:18)∂R(cid:19)2 (cid:18)∂Z(cid:19)2
p
= = + . (34)
∂θ ∂θ · ∂θ ∂θ ∂θ
Scaling of momentum flux with poloidal shaping mode number 10
The form of equation (32) is useful because it does not contain any radial derivatives
(except dψ/dr which is an input to the calculation) and distinguishes the important
ψ
(cid:12) (cid:12)2
term: the poloidal curvature. This allows us to find (cid:126)ψ (cid:126)α, (cid:12)(cid:126)α(cid:12) , and
(cid:12) (cid:12)
∇ ·∇ ∇
(cid:18) (cid:19)
1 B ∂R ∂R∂α (cid:16) (cid:17)
0 (cid:126) (cid:126) (cid:126)
v = + ψ θ ζ (35)
dsα
Ω R ∂ψ ∂θ ∂ψ∇ · ∇ ×∇
s 0
to lowest order in aspect ratio.
2.2. Asymptotic expansion ordering
We know from reference [21] that, unless the up-down symmetry of the geometric
coefficients is broken, the time-averaged momentum flux will always be zero to lowest
order in ρ ρ /a 1. We will investigate the consequences of breaking the up-down
∗ i
≡ (cid:28)
symmetry using different shaping effects. To do this we will expand equations (2), (11),
and (22) in m 1 using
(cid:29)
h =h +h +h +h +... (36)
s s0 s1 s2 s3
φ =φ +φ +φ +φ +..., (37)
0 1 2 3
where the subscript indicates the order in m−1 1. This expansion separates the long
(cid:28)
spatial scale coordinate θ, from the short spatial scale coordinate
z mθ. (38)
≡
Distinguishing the variation on each scale, e.g. h (θ,z) and φ(θ,z), means that
s
(cid:12) (cid:12) (cid:12)
∂ (cid:12) ∂ (cid:12) ∂ (cid:12)
(cid:12) = (cid:12) +m (cid:12) . (39)
∂θ(cid:12) ∂θ(cid:12) ∂z(cid:12)
w ,µ z,w ,µ θ,w ,µ
|| || ||
Ultimately we will only be interested in large scale phenomena, so we will need to
average quantities in z using
1 (cid:73) π(cid:12)(cid:12)
(...) (cid:12) dz(...), (40)
≡ 2π (cid:12)
−π θ
but we must still manipulate the z-dependent portion, given by
((cid:103)...) (...) (...). (41)
≡ −
2.3. Practical non-mirror symmetric shaping in the gyrokinetic model
In this section we will expand the large aspect ratio gyrokinetic, quasineutrality, and
momentum flux equations order-by-order to determine the scaling of the momentum
flux with m 1. Hence formally we require that (cid:15) 1 for the aspect ratio expansion
(cid:29) (cid:28)
and also that (cid:15) m−1 1 for the subsidiary expansion in shaping mode number.
(cid:28) (cid:28)