Table Of ContentRadially local approximation of the drift kinetic equation
H. Sugama,1,2 S. Matsuoka,3 S. Satake,1,2 and R. Kanno1,2
1)National Institute for Fusion Science, Toki 509-5292, Japan
2)Department of Fusion Science, SOKENDAI (The Graduate University for Advanced Studies), Toki 509-5292,
Japan
3)Japan Atomic Energy Agency, 178-4, Wakashiba, Kashiwa 277-0871, Japan
(Dated: 24 March 2016)
Anovelradiallylocalapproximationofthe driftkinetic equationispresented. The newdriftkinetic equation
6
that includes both E B and tangential magnetic drift terms is written in the conservative form and it
1
×
0 has favorable properties for numerical simulation that any additional terms for particle and energy sources
2 are unnecessary for obtaining stationary solutions under the radially local approximation. These solutions
satisfy the intrinsic ambipolarity condition for neoclassical particle fluxes in the presence of quasisymmetry
r
a ofthe magnetic fieldstrength. Also,anotherradiallylocaldrift kinetic equationis presented,fromwhichthe
M
positivedefinitenessofentropyproductionduetoneoclassicaltransportandOnsagersymmetryofneoclassical
transport coefficients are derived while it sacrifices the ambipolarity condition for neoclassical particle fluxes
3
in axisymmetric and quasi-symmetric systems.
2
PACS numbers: 52.25.Dg, 52.25.Fi, 52.25.Xz, 52.55.Hc
]
h
p
- I. INTRODUCTION difficultyinobtainingthe stationarysolutionofthe local
m
drift kinetic equation. However, the new local drift ki-
s Effects of neoclassical transport1–3 on plasma confine- neticequation,whichiswrittenintheconservativeform,
a
has favorable properties for numerical simulation such
l ment are more significant in stellarator and heliotron
p plasmas than in tokamak plasmas because, in the for- thatanyadditionaltermsforparticleandenergysources
.
s mer, radial drift motions of trapped particles in heli- are unnecessary for obtaining stationary solutions. In
c cal ripples enhance particle and heat transport due to addition, it satisfies the intrinsic ambipolarity condition
i
s nonaxisymmetry of the magnetic configuration.4–6 Con- for neoclassical particle fluxes in axisymmetric systems
y ventionalcalculationsofneoclassicaltransportfluxesare as well as in quasi-symmetric helical systems.16,17 The
h
present work also treats interesting issues regarding the
p doneapplyingradiallylocalapproximationtosolvingthe
[ driftkineticequation,inwhichvd f areoftenneglected entropy production rate and Onsager symmetry18,19 for
as a small term of higher order·i∇n the normalized gy- neoclassical transport equations resulting from the new
2 roradius parameter δ ρ/L. (Here, v , f, ρ, and L local drift kinetic model.
v ∼ d
represent the guiding center drift velocity, the deviation Therestofthispaperisorganizedasfollows. InSec.II,
1
9 of the guiding center distribution function from the lo- we consider the full drift kinetic model based on Little-
4 cal Maxwellian equilibrium distribution, the gyroradius, john’s guiding-center equations20 without radially local
3 andtheequilibriumscalelength,respectively.) However, approximation. Particle,energy,andparallelmomentum
0 in stellarator and heliotron plasmas, this v f term balance equations are derived from the full drift kinetic
d
. ·∇
1 is known to be influential on the resultant neoclassical equation. These balance equations are flux-surface aver-
0 transport because it significantly changes orbits of par- agedtoconfirmthattheycontainthesecond-orderterms
6 ticles trapped in helical ripples. Therefore, at least, the in δ, which represent neoclassical transport across flux
1 E Bdriftpartv f inv f hasbeenkeptinmost surfaces. Also,expandingthedistributionfunctionabout
E d
: × ·∇ ·∇
v studies of neoclassical transport in helical systems.7–15 the local Maxwellian, we rewrite the drift kinetic equa-
i Recently, it was shown by Matsuoka et al.13 that the tion to explicitly show that the thermodynamic forces
X
neoclassical transport is significantly influenced by re- defined by the backgrounddensity and temperature gra-
r
a taining the magnetic drift tangential to flux surfaces in dients and the parallel electric field cause the deviation
v f for the magnetic configurationofLHD especially f from the local Maxwellian. In Sec. II, a new drift
d
·∇
whentheradialelectricfieldisweak. However,aspointed kinetic model is constructed by applying radially local
by Landreman et al.,14 stationary solutions of the drift approximation to Littlejohn’s guiding-center equations
kinetic equation with radially local approximation used with keeping E B and tangential magnetic drift ve-
×
require additional artificial sources (or sinks) of parti- locities. The new local drift kinetic equation for f is
cles and energy when the above-mentioned drift terms showntobe compatible withthe stationarysolutionand
are retained. In this paper, a novel radially local drift togiveintrinsicambipolarparticlefluxesforaxisymmet-
kinetic equation,whichincludes bothE Bandtangen- ric and quasi-symmetric systems. In Sec. IV, we present
×
tial magnetic drift motions, is presented. The radially another radially local drift kinetic equation, from which
local guiding center motion equations do not satisfy the the positive definiteness of entropy production due to
conservationlawofthephase-spacevolumewhilethefull neoclassicaltransportandOnsagersymmetryofneoclas-
guidingcentermotionequationsdo. Thisfactcausesthe sical transport coefficients are derived although this lo-
2
caldrift kinetic equationno longerguaranteesrigorously which can be proved by using Eqs. (3) and (4).
the intrinsic ambipolarity of neoclassical particle fluxes Thedriftkineticequationforthedistributionfunction
for axisymmetric and quasi-symmetric systems. Finally, F(X,U,µ,t) is given by
conclusions are given in Sec. V.
∂ ∂
+X˙ +U˙ F(X,U,µ,t)=C(F)+ (6)
∂t ·∇ ∂U S
(cid:18) (cid:19)
II. FULL DRIFT KINETIC MODEL
where the total time derivative is denoted by˙= d/dt.
A. Drift kinetic model based on Littlejohn’s In the right-hand side of Eq. (6), C(F) is the collision
guiding-center equations term and the additional term is given to represent ex-
S
ternalparticle,momentum,and/orenergysourcesifany.
We denote the guiding-centervariablesby (X,U,ξ,µ), Here, is considered to be of the second order in δ. We
S
canalso treateffects ofturbulent fluctuations by Eq.(6)
where X represents the position vector of the guiding
if we regard the second-order additional term as the
center, U the parallel velocity, ξ the gyrophase defined
S
ensemble average of the product of fluctuation parts in
by the azimuthal angle of the gyroradius vector around
the electromagnetic fields and the distribution function
themagneticfieldline,andµthemagneticmoment. The
asshowninRefs.21and22wherethenotation isused
Lagrangianfortheguiding-centermotionisgivenbyLit-
D
tlejohn20 as instead of to represent the term including the effects
S
ofturbulentfluctuations. Using Eq.(5), the drift kinetic
e mc
L= A+mUb X˙ + µξ˙ H, (1) equation can be rewritten in the conservative form as
c · e −
where the Ham(cid:16)iltonian H is(cid:17)given by ∂(DF) + (DFX˙ )+ ∂(DFU˙) =D[C(F)+ ]. (7)
∂t ∇· ∂U S
1
H = mU2+µB+eΦ. (2)
2
B. Particle, energy, and parallel momentum balance
Here, Φ denotes the electrostatic potential. Using
equations
Eqs. (1) and (2), the guiding-center motion equations
are derived as
Multiplying Eq. (7) with an arbitrary function
dX c
=V Ub+ b (mU2b b+µ B eE∗), (t,X,U,µ) which is independent of the gyrophase ξ
dt gc ≡ eB∗ × ·∇ ∇ − A
k and taking its velocity-space integral, the balance equa-
dU 1 tion for the density variable d3vF in the X-space is
= b (µ B eE)+Ub b V , A
gc derived as
dt −m · ∇ − ·∇ ·
R
dξ ∂
=Ω, d3vF + d3vF X˙
dt ∂t A ∇· A
dµ (cid:18)Z (cid:19) (cid:18)Z (cid:19)
dt =0, (3) = d3v F ˙+[C(F)+ ] , (8)
A S A
whereΩ=eB/(mc), =∂/∂X,E Φ c−1∂A/∂t, Z (cid:16) (cid:17)
B A, E∗ ∇ Φ c−1∂A∗≡/∂−t,∇B∗− A∗, where
B∗≡ ∇B∗×b, andA≡∗−∇A+−(mc/e)Ub are us≡ed,∇an×dthe d ∂ ∂
k ≡ · ≡ ˙ = A = A +X˙ +U˙ A, (9)
guiding-center drift velocity Vgc is defined by the right- A dt ∂t ·∇A ∂U
hand side of the equation for dX/dt. On the right-hand
side of the equation for dU/dt in Eq. (3), the last term and the velocity-space integral is denoted by d3v =
Ub b V is smallerthanother terms by the orderof 2π dU dµD for gyrophase-independent integrands.
·∇ · gc R
δ =ρ/Lwhere ρ andL representthe gyroradiusandthe For the case of = 1, Eq. (8) reduces to the time-
gradient scale length given by L B/ B Φ/ Φ. evoRlutionRequatioAn for the density d3vF,
∼ |∇ |∼ |∇ |
TheJacobianfortheguiding-centervariablesiswritten
as ∂ d3vF + d3vFX˙R = d3v . (10)
∂(x,v) B∗ ∂t(cid:18)Z (cid:19) ∇·(cid:18)Z (cid:19) Z S
k
D =det = (4)
∂(X,U,ξ,µ) m Inderiving Eq.(10), the conservationlaw, d3vC(F)=
(cid:20) (cid:21)
0,isused. However,itisnotedthat,ifweusethecollision
where x and v denote the particle position vector and R
operator obtained by the transformation from the par-
the velocity vector, respectively. Then, the conservation
ticle coordinates to the guiding-center coordinates with
of the phase-space volume d3xd3v = Dd3XdUdξdµ is
finite-gyroradiuseffects taken into account, the velocity-
represented by
spaceintegral d3vC(F)doesnotvanishbutitbecomes
∂D ∂(DU˙) the opposite sign of the divergence of the classical parti-
+ (DX˙ )+ =0, (5) cle flux as showRn in Refs. 23–25. Here and hereafter, we
∂t ∇· ∂U
3
assumethattheexpressionofC(F)isthesameasthatof and rewrite Eq. (16) as
the Landau collision operator given in the particle coor-
dinates for simplicity so that d3vC(F)= 0 is satisfied ∂ (nmu )+b ( P) =neE +F + d3v mU. (18)
and the classical transport is neglected. ∂t k · ∇· k k S
R Z
We next consider the energy =H [see Eq. (2)] as
E A Here, the density n, the parallel flow velocity u , the
in Eq. (8) and obtain the energy balance equation, k
pressure tensor P, and the parallel friction force F are
k
∂ defined by
d3vF + d3vF X˙
∂t E ∇· E
(cid:18)Z (cid:19) (cid:18)Z (cid:19) n= d3vF
= d3v F ˙+[C(F)+ ] , (11) Z
E S E
Z (cid:16) (cid:17) nu = d3vFU
k
where the total time derivative of the energy is written
Z
as P=P +π
CGL 2
˙ = dE PCGL = d3vF mU2bb+µB(I bb)
−
E dt
Z
=e∂Φ(X,t) +µ∂B(X,t) e∂A∗(X,t) X˙ .(12) π2 = d3vFm(cid:0)U(X˙⊥b+bX˙ ⊥) (cid:1)
∂t ∂t − c ∂t ·
Z
We easilyseefromEq.(12)that ˙ =0forthe stationary F = d3vC(F)mU, (19)
k
E
electromagnetic field. When we use the kinetic energy, Z
whereX˙ X˙ (X˙ b)b. Note thatthe pressuretensor
W = 1mU2+µB = eΦ, (13) P consis⊥ts≡of th−e Ch·ew-Goldbeger-Low (CGL) tensor26
2 E −
P and the viscosity tensor π of the second order
CGL 2
another form of the energy balance equation is given by in δ, where π2 satisfies π2 : I = π2 : bb = 0 and the
deviation of F from the local Maxwellian distribution is
∂
d3vFW + d3vFWX˙ considered to be of (δ).
∂t ∇· ItiswellknownthOat,if weusethe originalBoltzmann
(cid:18)Z (cid:19) (cid:18)Z (cid:19)
kinetic equation instead of the drift kinetic equation in
= d3v FW˙ +[C(F)+ ]W , (14)
S Eq. (7), we can derive the momentum balance equation,
Z (cid:16) (cid:17)
where the total time derivative of the kinetic energy is ∂ u
(nmu)+ P=ne E+ B +F+ d3v mv,
written as ∂t ∇· c × S
dW d dΦ (cid:16) (cid:17) Z (20)
W˙ = = E e wheretheBoltzmannkineticequationisassumedtoalso
dt dt − dt
contain the source term . In Eq. (20), the particle flow
=µ∂B(X,t) +eE∗ X˙. (15) nu, the pressure tensorSP, and the friction force F are
∂t ·
defined by nu = d3vF, P = d3vFmvv, and F =
The parallel momentum balance equation is derived d3vC(F)mv, where, exactly speaking, F = F(x,v,t)
R R
from Eq. (8) with =mU as representsthe particle distributionfunction givenby the
A R
solution of the Boltzmann kinetic equation and it has a
∂
d3vFmU + d3vFmUX˙ gyrophase dependence that is not included in the solu-
∂t(cid:18)Z (cid:19) ∇·(cid:18)Z (cid:19) tion of the drift kinetic equation. Comparing Eqs. (18)
and(20), weseethatEq.(18)coincideswiththeparallel
= d3v FmU˙ +[C(F)+ ]mU . (16)
S component of the exact momentum balance equation in
Z (cid:16) (cid:17) Eq.(20)exceptthattheformercontainsnmu ∂b/∂tand
We now use ·
the non-CGL viscosity tensor expressed differently from
the one in the latter.
d3vFmUX˙ d3vFmU˙
∇· − We now consider general toroidal configurations, for
(cid:18)Z (cid:19) Z
which the magnetic field is written in terms of the flux
= d3vFmU2b + d3vFb (µ B eE) coordinates (s,θ,ζ) as
∇· · ∇ −
(cid:18)Z (cid:19) Z
B=ψ′ s θ+χ′ ζ s, (21)
+ d3vFmUX˙ d3vFmUX˙ (b )b ∇ ×∇ ∇ ×∇
⊥ ⊥
∇· − · ·∇
(cid:18)Z (cid:19) Z whereθ andζ representthe poloidalandtoroidalangles,
=b d3vF mU2bb+µB(I bb) respectively, and s is an arbitrary label of a flux surface.
· ∇· − The poloidal and toroidalfluxes within a flux surface la-
(cid:18) (cid:20)Z
(cid:0) beled by s are givenby 2πψ(s) and 2πχ(s), respectively.
+mU(X˙ b+bX˙ ) eE d3vF, (17)
⊥ ⊥ − k The derivative with respect to s is denoted by ′ = d/ds
(cid:17)i(cid:17) Z
4
so that ψ′ = dψ/ds and χ′ = dχ/ds. Taking the flux- C. Drift kinetic equation expressed in terms of flux
surface average of the covariant toroidal component of coordinates
Eq.(20)andmaking the summationoverspecies, weob-
tain the expression for the radial current as6 Using the flux coordinates (s,θ,ζ), the drift kinetic
equation, Eq. (6), is rewritten as
χ′ ∂
e n us = m n u + ( P )
c ah a ai a∂th a aζi h ∇· a ζi ∂ ∂ ∂ ∂ ∂
a a (cid:20) +s˙ +θ˙ +ζ˙ +U˙ F(s,θ,ζ,U,µ,t)
X X ∂t ∂s ∂θ ∂ζ ∂U
(cid:18) (cid:19)
d3v amavζ , (22) =C(F)+ (26)
− S S
(cid:28)Z (cid:29)(cid:21)
where
wherethesuperscriptsandthesubscriptζ representthe
covariant radial component and contravariant toroidal d
[s˙,θ˙,ζ˙]= [s,θ,ζ]
component given by taking the inner products with s dt
∇
and ∂x/∂ζ, respectively, and the subscript a is used to ∂
explicitly show the particle species. Using the symmetry = +X˙ [s(X,t),θ(X,t),ζ(X,t)]. (27)
∂t ·∇
property of the pressure tensor P, we can show that, for (cid:18) (cid:19)
axisymmetric toroidal systems, InEq.(27),thefunctionss(X,t),θ(X,t),andζ(X,t)are
definedbytheinverseofX=X(s,θ,ζ,t),wheretisgen-
( P) = 1 ∂ (V′ Ps ), (23) erally included as a parameter. Denoting the Jacobian
h ∇· ζi V′∂s h ζi for the flux coordinates (s,θ,ζ) by
waxhiesyremPmζset=ric∇tso·roPid·a∂lxs/y∂stζe.mIsn,1a7xwiseymhamveetric and quasi- √g =det ∂(∂s(,Xθ,)ζ) = [ s ( 1θ ζ)], (28)
(cid:20) (cid:21) ∇ · ∇ ×∇
( PCGL)ζ =0. (24) the conservationlaw of the phase-space volume, Eq. (5),
h ∇· i
and the conservative form of the drift kinetic equation,
Then, using Eqs. (22)–(24) and P = PCGL + π2, we Eq. (7), are rewritten as
find that, even for axisymmetric toroidal systems in the
stationary state (∂/∂t = 0) with = 0, the surface- ∂(√gD) ∂(√gDs˙) ∂(√gDθ˙) ∂(√gDζ˙)
S + + +
averaged radial current does not vanish exactly due to ∂t ∂s ∂θ ∂ζ
the second-order viscosity tensor π as shown by
2 ∂(√gDU˙)
+ =0, (29)
χ′ 1 ∂ ∂U
e n us = [V′ (π )s ]. (25)
c ah a ai V′∂s h a2 ζi and
a a
X X
However, it is shown in Ref. 17 that (π )s is a small ∂(√gDF) ∂(√gDFs˙) ∂(√gDFθ˙) ∂(√gDFζ˙)
h a2 ζi + + +
quantityof (δ3)inaxisymmetricsystemswithup-down ∂t ∂s ∂θ ∂ζ
O
symmetry(aswellasinquasi-axisymmetricsystemswith ∂(√gDFU˙)
stellarator symmetry) where all terms in the toroidal + =√gD[C(F)+ ], (30)
∂U S
momentum balance equation given from Eq. (22) van-
ish up to (δ2). The same argument as above can be respectively.
done for oOther quasi-symmetric systems such as quasi- For an arbitrary function (s,θ,ζ,U,µ,t) which is in-
A
poloidally-symmetric and quasi-helically-symmetric sys- dependentofthegyrophaseξ,thephase-spaceintegralis
tems if stellarator symmetry holds. On the other hand, written as
in axisymmetric systems without up-down symmetry,
hla(πriat2y)sζcion=diOti(oδn3) is neotnguuasran=tee0d.isTnhoetn,authtoemaamtibciaploly- 2πZ d3XZ dUZ dµDA
satisfied on the secaonadhoardeariin δ because of the third- =2π ds dθ dζ√g dU dµD
A
order radial parPticle fluxes (c/e χ′V′)∂[V′ (π )s ]/∂s Z I I Z Z
a h a2 ζi
driven by the second-order shear viscosity tensor com- = dsV′ d3v , (31)
ponents (π )s [here, it is useful to formally regard the A
a2 ζ Z (cid:28)Z (cid:29)
electric charge as the (δ−1) quantity27 so that the ra- where
O
dial current due to the third-order radial particle flux is
1
immediately found to be of the second order]. However, = dθ dζ√g (32)
eveninthis axisymmetricbutup-downasymmetriccase, h···i V′ ···
I I
thesecond-orderradialneoclassicalparticlefluxesdriven
represents the flux-surface average and
by the CGL tensors still automatically satisfy the am-
bipolarity condition for the radial current up to the first dV
order.1–3 V′ = ds = dθ dζ√g (33)
I I
5
denotestheradialderivativeofthevolumeV(s)enclosed where the zeroth-order density n and temperature T
0 0
within a flux surface labeled by s. We now integrate areflux-surfacefunctionsindependentof(θ,ζ),andtheir
Eq. (30) with respect to the coordinates (θ,ζ,U,µ) to time dependence follows the transport ordering, ∂/∂t =
obtain (δ2). The parallel electric field E is given by
k
O
∂ ∂
V′ d3v F + V′ d3v F s˙ 1∂A
∂t(cid:18) (cid:28)Z A(cid:29)(cid:19) ∂s(cid:18) (cid:28)Z A (cid:29)(cid:19) Ek =−b· ∇Φ+ c ∂t
(cid:18) (cid:19)
=V′ d3v F ˙+[C(F)+ ] , (34) BEk e BEk
(cid:28)Z (cid:16) A S A(cid:17)(cid:29) =BhhB2ii +(cid:18)Ek−BhhB2ii(cid:19), (41)
where
where Φ=Φ Φ . We now define the first-orderdistri-
˙ = dA = ∂A +s˙∂A +θ˙∂A +ζ˙∂A +U˙ ∂A. (35) bution f by −h i
A dt ∂t ∂s ∂θ ∂ζ ∂U
e
e l BE
The time-evolution equation for the surface-averaged F =F 1+ dl E Bh ki +f, (42)
density d3vF isderivedfromEq.(34)withA=1as 0" T0 Z (cid:18) k− hB2i (cid:19)#
∂ (cid:10)RV′ (cid:11)d3v F + ∂ V′ d3v F s˙ where ldl represents the integral along the magnetic
∂t ∂s field line. Then, substituting Eqs. (41) and (42) into
(cid:18) (cid:28)Z (cid:29)(cid:19) (cid:18) (cid:28)Z (cid:29)(cid:19) R
Eq. (26) yields
=V′ d3v . (36)
S df F W 5 eUB
(cid:28)Z (cid:29) = 0 Vs X +X + X
For the casesof A=W, Eq. (34) reduces to the surface- dt T0 (cid:26) gc(cid:20) 1 2(cid:18)T0 − 2(cid:19)(cid:21) hB2i1/2 E(cid:27)
averagedenergy balance equation, +CL(f)+ (δ2), (43)
O
∂ ∂
V′ d3v FW + V′ d3v FW s˙ where the thermodynamic forces are defined by
∂t ∂s
(cid:18) (cid:28)Z (cid:29)(cid:19) (cid:18) (cid:28)Z (cid:29)(cid:19) 1 ∂p ∂Φ ∂T BE
X = 0 e , X = 0, X = h ki ,
=V′ d3v FW˙ +[C(F)+ ]W , (37) 1 −n ∂s − ∂s 2 − ∂s E B2 1/2
0
S h i
(cid:28)Z (cid:16) (cid:17)(cid:29) (44)
where W˙ is given by Eq. (15). In Eqs. (36) and (37), and CL(f) represents the linearized collision operator.
d3v F s˙ and d3v FWs˙ represent the radial neo- Note that all terms explicitly shown on the right-hand
classical transport fluxes of particles and energy, respec- side of Eq. (43) are of the first order in δ. Using the
t(cid:10)iRvely, whic(cid:11)h are(cid:10)rRegardedas o(cid:11)f (δ2) assuming that the transport ordering ∂/∂t= (δ2) and f = (δ), the left-
deviation of F from the local MOaxwellian is of (δ) (see hand side of Eq. (43) is wriOtten as O
Sec.II.D).Theradialtransportfluxesof (δ2)aOreconsis- df ∂ ∂ ∂ ∂
tentwiththe so-calledtransportorderinOg1 whichimplies dt = Vgsc∂s +Vgθc∂θ +Vgζc∂ζ +U˙ ∂U f(s,θ,ζ,U,µ)
∂/∂t= (δ2) in Eqs. (36) and (37). (cid:18) (cid:19)
O + (δ3)
O
1 ∂( fVs) ∂( fVθ) ∂( fVζ) ∂( fU˙)
D. Expansion about a local Maxwellian distribution = D gc + D gc + D gc + D
∂s ∂θ ∂ζ ∂U
!
D
The zeroth-order solution F of the drift kinetic equa- + (δ3), (45)
0
O
tion, Eq. (26), is given by the local Maxwellian,
whereVs =X˙ s,Vθ =X˙ θ, Vζ =X˙ ζ,and =
m 3/2 W √gD. Sgicnce Vs·∇= (gδc), the·r∇adialgdcrift te·r∇m Vs∂fD/∂s
F =n exp gc O gc
0 0 in Eq. (45) is of the second order in δ and this gives
2πT −T
(cid:18) 0(cid:19) (cid:18) 0(cid:19) rise to global or finite-orbit-width effects on neoclassical
3/2
m eΦ transport.
=n exp E − , (38)
0
2πT − T
(cid:18) 0(cid:19) (cid:18) 0 (cid:19)
which annihilates the collision term,
III. RADIALLY LOCAL APPROXIMATION
C(F )=0. (39)
0
The total time derivative of F is written as Undertheradiallylocalapproximationmadehere,the
0
guiding center equations are written as
dF dlnn dlnT mv2 3 1 dW
0 0 0
=F0 + dX
dt (cid:20) dt dt (cid:18)2T0 − 2(cid:19)− T0 dt (cid:21) dt =Vg(rcl) ≡Ub+(Vg(rcl))⊥,
∂lnn e ∂ Φ ∂lnT
=F X˙ s 0 + h i + 0 dU µ
0 ·∇ ∂s T ∂s ∂s = b B+Ub b V(rl)
(cid:26) (cid:20) 0 dt −m ·∇ ·∇ · gc
W 3 eUEk + (δ2), (40) dµ = 1(V(rl)) (mU2b b+µ B), (46)
×(cid:18)T0 − 2(cid:19)(cid:21)− T0 (cid:27) O dt −B gc ⊥· ·∇ ∇
6
wherethesecond-orderpart Φ c−1∂A/∂toftheelec- We should note that influences of the magnetic and
−∇ −
tricfieldEisneglectedandtheguidingcenterdriftveloc- E B drift motions aresignificant mainly for precession
×
ityinEq.(3)isreplacedbyV(rl)ewhichhasnoradialcom- driftorbitsoftrappedparticles,andthatparticlesinthe
gc
ponent: V(rl) s=0. The component of V(rl) perpen- region,Λ<Λ0, are passing ones whose orbits almostco-
gc gc
·∇ incide with field lines. Therefore, even if the functional
dicular to the magnetic field is denoted by by (V(rl)) .
gc ⊥ form of α(Λ) and the value of Λ are changed, the ar-
We later impose the condition, (Vg(rcl))⊥(µ = 0) = 0, in tificial reduction factor α(Λ) for Λ0 < Λ0 is expected to
ordertoderiveappropriatebalanceequationsofparticles, cause little change in resultant passing particles’ orbits
energy, and parallel momentum [see Eqs. (54), (57), and exceptthatthelimitingcondition,lim (V(rl)) =0,
Λ→+0 gc ⊥
(61)] by removing improper sources and/or sinks at the
is rigorously satisfied. However, this insensitivity to the
boundary µ=0 in the velocity-space integral domain.
form of α(Λ) remains a future subject to be verified by
InEq.(46),the magneticmomentµ isallowedtovary
numerical simulations.
in time such that conservation of the kinetic energy of
It also should be mentioned that the radially local ap-
the particle W =mU2/2+µB,
proximation described by Eqs. (46) and (48) is indepen-
dW dU dµ dent of what poloidal and toroidal angles are chosen for
=mU +B +µV(rl) B =0, (47)
dt dt dt gc ·∇ the flux coordinates. This is a favorable property that
is lost in Ref. 13. We see that the radially local guid-
issatisfied. Itmightappearthatthe energy =W+eΦ
E ingcenterequationsgivenbyEqs.(46)andthe Jacobian
should be conserved instead of W. However, using Φ
≃ D = B∗/m for the phase-space coordinates (X,U,ξ,µ)
Φ , we find that the difference eΦ between and W is k
h i E [see Eq. (4)] do not satisfy the conservation law of the
approximately constant along the radially local guiding
phase-space volume as shown in Eq. (5). This violation
center orbit and accordingly the conservation of W is
ofthephase-space-volumeconservationoccursevenifB∗
reasonable under the radially local approximation. k
We now define (V(rl)) by removing the radial com- is used instead of B in the denominator on the right-
ponent from (Vgc)⊥ agcs ⊥ hand side of Eq.(48) to define (Vg(rcl))⊥. In the next sec-
tion, we consider another Jacobian in order to recover
(V ) s
(V(rl)) =α(Λ) (V ) gc ⊥·∇ s the conservation law although, in this section, a sim-
gc ⊥ gc ⊥− s2 ∇ pler approximate Jacobian D B/m is used. Also, we
(cid:18) |∇ | (cid:19) 0 ≡
c s hereafteremploy(X,W,U,ξ)asphase-spacecoordinates.
=α(Λ)eB(b×∇s) ∇s2 ·[mvk2b·∇b+µ∇B] Then,fromthe JacobianD0 B/m for(X,U,ξ,µ)with
(cid:18)|∇ | µ = (W 1mU2)/B, the J≡acobian for (X,W,U,ξ) is
dΦ − 2
+e derived as 1/m, which is constant in the phase space.
ds
(cid:19) Using V(rl), dU/dt, and dW/dt=0 given by Eqs. (46)
c mU2 s gc
=α(Λ) (b s) +µ ∇ B with (48) under the radially local approximation, the
eB ×∇ B s2 ·∇ driftkineticequationforthefirst-orderdistributionfunc-
(cid:20)(cid:18) (cid:19)|∇ |
4π dP dΦ tion f(X,W,U) in the stationary state is written as
+mU2 +e , (48)
B2 ds ds
(cid:21)
where Bk∗ and E∗ in the definition of Vgc given by ∇·(fVg(rcl))+ ∂∂U fddUt
Eq. (3) are replaced with their lowest-orderparts B and (cid:18) (cid:19)
d−u(cdeΦd/tdos)s∇atsis,fyretshpeecctoinvedliyt,ioann(dVtg(hrcle))f⊥ac(µto=r α0()Λ=) 0is.iHnterroe-, = FT00 (cid:26)Vgsc(cid:20)X1+X2(cid:18)WT0 − 52(cid:19)(cid:21)+ hBeU2iB1/2XE(cid:27)
theratioofthemagneticmomentµtothekineticenergy +CL(f). (51)
W = mU2/2+µB is used to define the dimensionless
parameter,Λ µB /W,where B is the maximum
≡ max max Theradialcomponentoftheguidingcenterdriftvelocity
value of B on the flux surface. This parameter Λ is a
Vs on the right-hand side of Eq. (51) is given by
measureforclassifyingtheguidingcentermotionintoei- gc
ther passingortrappedorbit. As Λincreasesfrom0 and
approaches to 1, the orbit changes from the passing to c 1
Vs = [ s (b B)] mU2+W . (52)
the trapped one. Then, we assume that gc eB2 ∇ · ×∇ 2
(cid:18) (cid:19)
lim α(Λ)=0 (49)
Λ→+0 In deriving Eq. (52) from the guiding center drift ve-
locity given in Eq. (3), only the lowest-order terms in
whileα(Λ)=1exceptforaninterval,0 Λ<Λ ,where
0
≤ δ is retained and the formula, s [b (b )b] =
Λ ( 1)is asmallpositive constantvalue. Forexample,
0 ∇ · × · ∇
≪ s (b B)/B, obtained from the MHD equilibrium
α(Λ) is defined by
∇cond·itio×n ∇[n (s)T (s)] = (4π)−1( B) B is used.
0 0
∇ ∇× ×
sin(πΛ/2Λ ) (Λ<Λ ) The fact that the Jacobian is constant is used in deriv-
α(Λ)= 0 0 (50)
1 (Λ Λ ). ing Eq. (51) which is rewritten by using the flux surface
0
(cid:26) ≥
7
coordinates (s,θ,ζ) as wheretheparallelandperpendicularheatfluxesaregiven
by
1 ∂ ∂
(√gfV(rl) θ)+ (√gfV(rl) ζ) 5
√g ∂θ gc ·∇ ∂ζ gc ·∇ q = d3vf W T U,
(cid:20) (cid:21) k − 2
∂ dU Z (cid:18) (cid:19)
+ ∂U f dt q⊥1 = 5p0cX2 s b,
(cid:18) (cid:19) 2 eB ∇ ×
= FT00 (cid:26)Vgsc(cid:20)X1+X2(cid:18)WT0 − 25(cid:19)(cid:21)+ hBeU2iB1/2XE(cid:27) q(⊥r2l) =Z d3vf(cid:18)W − 25T(cid:19)(Vg(rcl))⊥, (58)
+CL(f). (53)
and the collisional heat generation is defined by
Here, we should note that, in Eq. (53), partial deriva-
Q= d3vCL(f)W. (59)
tives of the first-order distribution function f are taken
onlywithrespecttothethreevariables(θ,ζ,U)andthat Z
the radial coordinate s and the kinetic energy W enter In Eq. (57), q(rl) is the second-order flux like Γ(rl) in
⊥2 ⊥2
f(s,θ,ζ,W,U) as constant parameters. Eq.(54). Taking the flux surface averageof Eq. (57), we
Takingthevelocity-spaceintegralofEq.(51)yieldsthe obtain
continuity equation in the stationary state,
Q =0, (60)
h i
(Γ b+Γ +Γ(rl))=0. (54)
∇· k ⊥1 ⊥2 which represents the collisional heat exchange balance
that needs to be satisfied in the stationary state. Un-
TheparallelandperpendicularparticlefluxesinEq.(54)
equal temperatures T = T can occur in the case of
are defined by a0 6 b0
m /m 1 or m /m 1 where the characteristic
a b a b
≪ ≫
time of the collisional thermal equilibration between the
Γ =n u = d3vfU, species a and b is much longer than the 90◦ scattering
k 0 k
Z times due to like-species collisions characterized by τaa
n cX
Γ⊥1 =n0u⊥1 = 0 1 s b, and τbb. Then, Cab(fa0,fb0) does not vanish even for
eB ∇ × the local Maxwellian distribution functions f and f
a0 b0
Γ(rl) =n u(rl) = d3vf(V(rl)) , (55) given by Eq. (38) and it describes the above-mentioned
⊥2 0 ⊥2 gc ⊥ slow collisional thermal equilibration although the lin-
Z
earized collision operator CL used for the Eq. (53) does
where the velocity-space integral is written in terms of not include this equilibrium part of the collision term.
the variables W and U by However, the heat generation Q , which is defined by
ab
Eq.(59)withthelinearizedoperatorCL forcollisionsbe-
ab
2π +∞ +√2W/m tweendifferentspeciesaandb,generallyremainsnonzero
d3v = dW dU. (56) (even for the case of T = T ). Therefore, Eq. (60),
Z m Z0 Z−√2W/m which is rewritten as aQ0 b0 Q = 0 (recall
h ai ≡ b6=ah abi
Q 0), is considered to be the physically reasonable
ThediamagneticflowΓ⊥1 =n0u⊥1 andtheparallelflow coanadi≡tionthatshouldbesatisfiedPinthemulti-speciessta-
Γ are of the first order in δ =ρ/L while Γ(rl) =n u(rl)
k ⊥2 0 ⊥2 tionarystateoftheradiallylocalmodelwithoutrequiring
is of the second order. The collisional particle conserva-
additional heat source or sink.
tion law, d3vCL(f) = 0, is used to obtain Eq. (54).
Multiplying Eq.(51) with mU andtaking its velocity-
Also, it should be noted that the boundary condition,
spaceintegralgivethe parallelmomentumbalanceequa-
(Vg(rcl))⊥(µR=0)=(Vg(rcl))⊥(U = 2W/m)=0, is used tion,
±
for deriving Eq. (54) as well as the energy and parallel
p
momentum balance equations [see Eqs. (57) and (61)] BE
b [ p + (π +π(rl))]=n eBh ki +F (61)
from Eq. (51). We find that the flux surface average · ∇ 1 ∇· 1 2 0 B2 k
of the left-hand side of Eq. (54) automatically vanishes h i
so that no particle source is required for obtaining the wherethefirst-orderpressurep1andtheviscositytensors
stationary solution. Thus, the radially local approxima- π and π(rl) are defined by
1 2
tion presented here has self-consistency with neglecting
2
the radial transport that causes variationin the surface- p = d3vfW,
1
averagedparticles’ number [see Eq. (36)]. 3
Z
Next, we multiply Eq.(51) with (W 5T/2)and take 1
its velocity-space integral to derive − π1 = d3vf(mU2−µB) bb− 3I ,
Z (cid:18) (cid:19)
(q b+q +q(rl))=Q, (57) π(2rl) = d3vfmU (Vg(rcl))⊥b+b(Vg(rcl))⊥ , (62)
∇· k ⊥1 ⊥2 Z (cid:16) (cid:17)
8
and the parallel friction force is given by (Aw)+A s ( B)/B2 = (Aw), which leadto
∇· ∇ · ∇× ∇·
Eqs. (67) and (68), respectively. Then, using Eqs. (66)–
F = d3vCL(f)mU. (63) (68), we also have
k
Z ∂x
( π(rl)) =0. (69)
The first-orderviscosity tensor π1 is written in the form ∂ζ · ∇· 2
of the traceless part of the CGL pressure tensor as π = (cid:28) (cid:29)
1
(p p )(bb 1I),wherep andp representtheparallel It is found from Eq. (65) and Eqs. (67)–(69) that the
k− ⊥ −3 k ⊥
and perpendicular pressures, respectively. On the other second-orderviscositytensorπ(rl)intheradiallylocalap-
2
hand, the second-order viscosity tensor π(rl), which is proximation cannot contribute to the neoclassical trans-
2
given by the correlation between the parallel velocity U port like the first-order pressure p and viscosity tensor
1
and the perpendicular drift velocity (Vg(rcl))⊥ , can not π1.
be written in the CGL form. We now multiply Eq. (61) We now use Eqs. (61) and (68) to rewrite the expres-
withthemagnetic-fieldstrengthB andtakeitsmagnetic- sion of the radial neoclassical particle flux in Eq. (65)
surface averageto derive as
c ∂x
hB·[∇·(π1+π(2rl))]i=n0ehBEki+hBFki, (64) Γncl = eχ′ (cid:28)∂ζ ·(∇p1+∇·π1)](cid:29)
whichisusedlatertoderiveanalternativeexpressionfor cBζ n ehBEki + Fk . (70)
the neoclassical particle flux. − eχ′ 0 B2 B
(cid:18) h i (cid:28) (cid:29)(cid:19)
The radial neoclassical particle flux is written as
Here, the first surface-averaged part on the right-hand
of Eq. (70) represents the nonaxisymmetric part of the
Γncl = d3vfVgsc neoclassical particle flux,6,19
(cid:28)Z (cid:29)
c s c ∂x
= ∇ [b ( p + π )] Γna = ( p + π )] . (71)
e B · × ∇ 1 ∇· 1 eχ′ ∂ζ · ∇ 1 ∇· 1
(cid:28) (cid:29) (cid:28) (cid:29)
c ∂x
= ( p + π ) Then, we find from Eqs. (70) and (71) that the radial
eχ′ ∂ζ · ∇ 1 ∇· 1 electric current is written as
(cid:28) (cid:29)
cB b
− eχζ′ B ·(∇p1+∇·π1) , (65) eaΓnacl = eaΓnaa, (72)
(cid:28) (cid:29) a a
X X
where Vs =V s is given by Eq. (52). Derivation of
Eq. (65)gcuses thgce·f∇ollowing formula, where the quasineutrality ana0ea = 0 and the colli-
sional momentum conservation F =0 are used.
P a ka
Foraxisymmetricandquasi-axisymmetricsystems,we
s b ∂x B
χ′∇ × = ζb, (66) have P
B ∂ζ − B
∂p ∂x
andthe Boozercoordinates(s,θ,ζ),28 forwhichthe con- ∂ζ1 = ∂ζ ·(∇·π1) =0, (73)
travariantpoloidalandtoroidalcomponents, B andB , (cid:28) (cid:29) (cid:28) (cid:29)
θ ζ
ofthemagneticfieldBareflux-surfacefunctions. Wesee from which
fromEq.(65)thattheneoclassicalparticlefluxiscaused
Γna =0, (74)
by the spatial gradients of the first-order pressure and
viscosity tensor. It can be shown that the second-order and
viscosity tensor π(rl) defined in Eq. (62) satisfies
2 cB BE F
Γncl = ζ n eh ki + k (75)
s −eχ′ 0 B2 B
∇ b ( π(rl)) =0 (67) (cid:18) h i (cid:28) (cid:29)(cid:19)
B · × ∇· 2
(cid:28) h i(cid:29) are derived. Here, the quasi-axisymmetry means that
the magnetic field strength B = B is independent of
and | |
the toroidal angle ζ. For derivation of Eq. (73), the ζ-
b independence of B and √g = (4π2)−1(dV/ds) B2 /B2
( π(rl)) =0. (68) h i
B · ∇· 2 in the Boozer coordinates10 is used [see also Eq. (25) in
(cid:28) (cid:29)
Ref. 17]. It is confirmed from Eqs. (72) and (74) that
In deriving Eqs. (67) and (68), it is convenient to write the solution f of the radially local drift kinetic equation
π(rl) = A(Bw+wB) with w (b s)/B. Then, we shown in Eq. (51) or (53) gives the neoclassical particle
fin2d w ( π(rl))=w2B A+≡A(B×∇w w w w fluxes which satisfy the ambipolarity condition,
B) = w·2B∇· 2A A s (·∇ w) = ·∇(A· s− w·)∇and·
(b/B)( π·(∇rl))=− ∇(Aw· )∇+×A(w B∇b· b∇B×w)/B= eaΓnacl =0, (76)
· ∇· 2 ∇· ·∇ · − ·∇ · Xa
9
automatically in axisymmetric and quasi-axisymmetric Here, CT(f ) C (f ,F ) and CF(f ) C (F ,f )
ab a ≡ ab a b0 ab b ≡ ab a0 b
systems. The intrinsic ambipolarity can be proved in are called test- and field-particle collision operators, re-
the same way for all other quasi-symmetric systems spectively,andtheysatisfytheadjointnessrelations,30,31
such as quasi-poloidally-symmetric and quasi-helically-
f g
symmetric systems. d3v a CT(g )= d3v a CT(f ),
F ab a F ab a
It is seen from Eqs. (22)–(25) that the radial current Z a0 Z a0
is closely related to the toroidal viscosity or the radial T d3v fa CF(g )=T d3v gb CF(f ), (79)
transport of the toroidal momentum. As remarked after a F ab b b F ba a
Z a0 Z b0
Eq. (25), the ambipolarity condition is not guaranteed
and Boltzmann’s H-theorem,30,31
on the second order in δ for axisymmetric systems with-
out up-down symmetry (as well as quasi-axisymmetric f
systems without stellarator symmetry) because of the Ta d3vFa [CaTb(fa)+CaFb(fb)]
component(π )s ofthesecond-ordernon-CGLviscosity Z a0
a2 ζ f
tensor. We also note that the radialneoclassicalparticle +T d3v b [CT(f )+CF(f )] 0. (80)
b F ba b ba a ≤
flux defined by Eq. (65) is the second-order flux driven Z b0
bythefirst-orderCGLtensorwhichbecomesadominant Strictly speaking, the adjointness relations and the H-
part for nonaxisymmetric systems although it does not theorem are rigorously satisfied by the linearized Lan-
containthe third-orderflux due to the second-orderten- dau collision operator only for the case of T = T al-
a b
sor. Therefore, Eq. (76) should be interpreted to imply thoughtheyarestillapproximatelyvalidevenforT =T
a b
that the intrinsic ambipolarity condition for the axisym- when (m /m )1/2 or (m /m )1/2(1 T /T ) are s6mall
a b b a b a
metric case can be correctly treated only up to the first enough.30 −
order by the present radially local approximation. On Since the drift kinetic equations for different particle
the other hand, the second-order neoclassical radial flux speciesarecoupledwitheachotherduetothefieldparti-
h(πa2)sζi of the toroidal momentum in the axisymmetric clecollisionoperators,fa dependsnotonlyonthermody-
butup-downasymmetriccasecanalsobeevaluatedusing namic forces (X ,X ,X ) but also on those for b=a,
a1 a2 E
the solution f of the radially local drift kinetic equation (X ,X ). Accordingly,wefindthatΓncl,qncl,andJ6 in
b1 b2 a a E
even without resort to the radially global model. This Eq.(77)arerelatedtothethermodynamicforcesthrough
can be done by substituting the solution f into the for- theneoclassicaltransportequationswhicharewrittenas
mula for the toroidalmomentum transportflux givenby
Eq. (18) in Ref. 17. [It is confirmed from Eqs. (11) and Γncl = (L11X +L12X )+L1 X ,
a ab b1 ab b2 aE E
(13)inRef.29that,ifusingthedefinitionofπ givenby
2 Xb
Eq. (19) in the present work to evaluate h(πa2)sζi, only a qancl/Ta = (L2a1bXb1+L2a2bXb2)+L2aEXE,
partoftheresultfromEq.(18)inRef.17isreproduced.]
b
X
J = (L1 X +L2 X )+L X . (81)
E Eb b1 Eb b2 EE E
b
IV. ENTROPY PRODUCTION RATE AND ONSAGER X
SYMMETRY ASSOCIATED WITH NEOCLASSICAL Here, the neoclassical transport coefficients
TRANSPORT EQUATIONS (L11,L12, ) are regarded as functions of the variables
ab ab ···
[E ( dΦ/ds), s B, s (b b)] which charac-
s
≡ − ∇ ·∇ ∇ · ·∇
TheneoclassicalradialparticlefluxΓncl,heatfluxqncl, terize the perpendicular guiding center velocity (V(rl))
gc ⊥
and parallelelectric currentJ = BJ / B2 1/2 are de- defined in Eq. (48).
E k
h i h i
fined in terms of the solution f of Eq. (51) or (53) by We here examine the Onsager symmetry of the neo-
classical transport coefficients. In order to prove the
Onsager symmetry, the adjointness relations written in
Γncl = d3vf Vs ,
a a gca Eq. (79) and the phase-space-volume conservation along
(cid:28)Z (cid:29)
the collisionless guiding center orbit are required as
5
qancl = d3vfaVgsca W − 2Ta , shown in Refs. 19 and 29. However, in the radially
(cid:28)Z (cid:18) (cid:19)(cid:29) local model based on Eq. (51), the latter condition
JE =hB2i−1/2 ea d3vfaU , (77) ∇·Vg(rcl)+∂(dU/dt)/∂U =0isbrokensothattheOnsager
a (cid:28)Z (cid:29) symmetry is not satisfied.
X
As noted before Eq. (51), the Jacobian for the phase-
where the subscript a denotes the particle species. The space coordinates(X,W,U,ξ) is givenby 1/m. Here, we
linearizedcollisionoperatorinEq.(51)forthespeciesais consider a modified Jacobian,
definedintermsofthebilinearoperatorC forcollisions
ab
between the species a and b by D =[1+d (X,W,U)]/m, (82)
W ∗
[C (f ,F )+C (F ,f )]. (78) whichdiffersfromtheonementionedabovebythecorrec-
ab a b0 ab a0 b
tiontermd of (δ)[seeEq.(84)below]. Thistermd is
Xb ∗ O ∗
10
determined by assuming that the Jacobian D satisfies The collision term CL(f) in Eq. (51) is replaced by
W
the conservation law of the phase-space volume element CL(f ) in Eq. (88) where the deviation of CL(f ) from
∗ ∗
written as CL(f) is of (δ2) and it is neglected. It is shown from
O
Eq. (83) that the differential operator satisfies the an-
∂ dU V
(D V(rl))+ D =0, (83) tisymmetry relation,
∇· W gc ∂U W dt
(cid:18) (cid:19)
d3vα β = d3vβ α , (90)
where Vg(rcl) and dU/dt are given by Eqs. (46) and (48). V − V
(cid:28)Z (cid:29) (cid:28)Z (cid:29)
We can rewrite Eq. (83) as
where α and β are arbitrary smooth functions on the
∂ ∂ dU ∂ phase space.
V(rl) θ +V(rl) ζ + lnD
gc ·∇ ∂θ gc ·∇ ∂ζ dt ∂U W Replacingf withf∗inEq.(77),wecandefinemodified
(cid:18) (cid:19) transport fluxes, Γncl, qncl, and J , the values of which
∂ ∂ dU ∂ ∗a ∗a ∗E
= V(rl) θ +V(rl) ζ + ln(1+d ) agreewiththoseofΓncl,qncl,andJ ,respectively,tothe
gc ·∇ ∂θ gc ·∇ ∂ζ dt ∂U ∗ a a E
(cid:18) (cid:19) lowest order in δ because f∗ = f[1+ (δ)]. Then, sub-
= V(rl) ∂ dU . (84) stituting the solution f∗ of Eq. (88) inOto the definitions
−∇· gc − ∂U dt of the modified transport fluxes, we can derive the neo-
(cid:18) (cid:19)
classical transport equations relating (Γncl,qncl,J ) to
Noting that the last line of Eq. (84) is of (δ), we can ∗a ∗a ∗E
O (Xb1,Xb2,XE). Thesetransportequationstakethesame
take the correction term d as a small quantity of (δ).
The left-hand side of Eq. ∗(84) represents the derivOative formsasthoseinEq.(81),andweuse(L1∗1ab,L1∗2ab,···)to
represent the modified transport coefficients which cor-
of lnD along the radially local guiding center orbit la-
W respond to (L11,L12, ) in Eq. (81), respectively. It is
beledbytheconstantparameters(s,W). Assumingthat ab ab ···
shown in the same way as in Sec. III that no additional
the guiding center orbit ergodically covers the (θ,ζ,U)
sources and/or sinks are required to obtain stationary
space, D = (1+d )/m is determined by Eq. (84) ex-
W ∗ particle and energy balances from Eq. (88)
cept for a factor that is an arbitrary function of (s,W).
We now multiply Eq. (88) for particle species a with
InordertouniquelyspecifyD =(1+d )/m,weimpose
W ∗ T f /F and take its velocity-space integral, flux-
a ∗a a0
another constraint,
surface average, and summation over species. Then, we
obtain
+√2W/m
dUd∗ =0, (85) T ΓnclX +qnclX +J X
*Z−√2W/m + a ∗a a1 ∗a a2 ∗E E
a
X(cid:0) (cid:1)
f
where representstheflux-surfaceaveragedefinedin = T d3v ∗a CT(f )+CF(f ) 0,(91)
Eq.(32h)·.··Oiwing to the condition in Eq. (85), d∗ is given − a,b a(cid:28)Z Fa0 ab ∗a ab ∗b (cid:29)≥
as a small correction and it does not affect the surface- X (cid:2) (cid:3)
averagedvelocityintegralofthe equilibriumdistribution where the inequality is due to the H-theorem given in
function F as shown by Eq. (80). Equation (91) means that the neoclassical
0
transport process is subject to the second law of ther-
modynamics: the summation of products between the
d3v(1+d )F = d3vF =n , (86)
∗ 0 0 0 transport fluxes and forces equals the entropy produc-
(cid:28)Z (cid:29) (cid:28)Z (cid:29)
tion rate expressed in terms of the linearized collision
where d3v is given by Eq. (56) and F0 is the local operator, which is positive definite.
MaxwelliandefinedinEq.(38)withtheequilibriumden- Since the differential operator and the linearized
sity n0Rand temperature T0 given as flux-surface func- collision operator CL satisfy the aVntisymmetry relation
tions. in Eq. (90) and the adjointness relations in Eq. (79),
We next define another distribution function f∗ by respectively, we can use the same procedures as in
Ref. 29 to prove that the modified transport coefficients
f∗ ≡ 1+fd∗, (87) (wLr1∗it1atbe,nL1∗a2asb,···) obey the Onsager symmetry relations
and use Eq. (83) to rewrite Eq. (51) in terms of f∗ as Li∗jab(β)=Lj∗iba(−β) (i,j =1,2),
f = F0 Vs X +X W 5 + eUB X Li∗aE(β)=−Li∗Ea(−β) (i=1,2),
V ∗ T0 (cid:26) gc(cid:20) 1 2(cid:18)T0 − 2(cid:19)(cid:21) hB2i1/2 E(cid:27) L∗EE(β)=L∗EE(−β), (92)
+CL(f∗), (88) where β [Es, s B, s (b b)] represent the
≡ ∇ · ∇ ∇ · · ∇
variables associated with the perpendicular guiding cen-
where the differential operator is defined by
V ter velocity (V(rl)) as explained after Eq. (81). Note
gc ⊥
dU ∂ that the change from β to β corresponds to turning
V ≡(1+d∗) Vg(rcl)·∇+ dt ∂U . (89) (Vg(rcl))⊥ in the opposite direc−tion.
(cid:18) (cid:19)
