Multiconstrained Variational Problems in Magnetohydrodynamics, II: Slow Evolution

Alexander Lifschitz  Bruce Turkington  Joel Spruck

Alexander Eydeland Department of Mathematics University of Massachusetts Amherst, MA 01003 August 9, 1990

1 Introduction 2 Slow evolution equations

We begin by recalling the full set of equations governing a plasma-vacuum system confined in a toroidal device such as a tokomak. Under the usual assumptions of ideal magnetohydrodynamics, the equations valid in the plasma region are

\frac{Dρ}{Dt} + ρ∇ · V = q_{M} (2.1)

\frac{Dζ}{Dt} + ζ∇ · V = q_{S} (2.2)

ρ\frac{DV}{Dt} + ∇p − J × B = Q (2.3)

\frac{∂B}{∂t} − ∇ × (V × B) = 0 (2.4)

∇ × B = J, ∇ · B = 0 ; (2.5) and the equations valid in the vacuum region are ∇ × B = \overline{J}, ∇ · B = 0 . (2.6)

Here ρ, ζ, p, V, B and J denote the mass density, entropy density (per unit volume), pressure, velocity, magnetic field and current density, respectively. The notation D/Dt: = ∂/∂t + V · ∇ is used for the convective derivative. We suppose that in the plasma there are distributed sources of mass, entropy and momentum production given by q_{M}, q_{S} and Q, respectively; and that in the vacuum there are (external field) coils carrying a given current density \overline{J}. We assume here and throughout that all equations are expressed in nondimensional variables.

It is worthwhile to note that in contrast to usual practice we use the entropy per unit volume ζ rather than the entropy per unit mass s(= ζ/ρ). This choice has two virtues: first, it permits the parallel development of the conservation laws (2.1) and (2.2); and second, it guarantees the convexity of the internal energy (per unit volume)

U(ρ, ζ) = e^{ζ/ρ}ρ^{γ}/(γ − 1) (1 < γ < ∞) (2.7)

as a function of ρ, ζ, where γ is the adiabatic index (ratio of specific heats). The corresponding equation of state in this notation takes the form

p = e^{ζ/ρ}ρ^{γ} . (2.8)

Most of our subsequent development actually applies to the general case in which U(ρ, ζ) is taken to be any smooth and strictly convex function, and the equation of state is derived from the generalized thermodynamic relation p = ρ ∂U/∂ρ+ζ ∂U/∂ζ− U (see [?]). However, we shall restrict our discussion to the familiar case (2.8) of an ideal gas for the sake of definiteness.

By means of the first law of thermodynamics, the energy production q_{E} can be related to q_{M} and q_{S}. Differentiating (2.7) gives

dU = \frac{(γ − 1)p}{ρ} [(γ − \frac{ζ}{ρ}) dρ + dζ] ,

which clearly means that

q_{E} = \frac{(γ − 1)p}{ρ} [(γ − \frac{ζ}{ρ})q_{M} + q_{S}] . (2.9)

Consequently, in practical applications q_{M} and q_{E} (say) may be considered as prescribed data instead of q_{M} and q_{S}.

The above equations hold in a toroidal region D. On the fixed boundary ∂D, which is assumed to be a perfectly conducting shell, the normal component of B vanishes. On the free boundary surface S, the plasma-vacuum interface conditions hold — namely, p = 0, n · B = 0 on S and n × B is continuous across S, where n is the unit normal on S. The fixed and free boundaries are therefore magnetic surfaces at every instant of time.

Let us suppose for the moment that the external sources of mass, entropy and momentum do not exist and that the external current density does not vary in time. Then the above evolutionary problem has a class of static (V = 0) equilibrium (∂/∂t = 0) solutions. In the plasma region these solutions satisfy the standard equilibrium equations

J × B = ∇p, ∇ × B = J, ∇ · B = 0 , (2.10)

while in the vacuum region they satisfy (2.6); on the fixed and free boundaries they fulfill the conditions stated above. It is important to notice that in equilibrium the density ρ and entropy ζ do not have a precise meaning, since only the pressure p, which is the combination (2.8), enters in the governing equations (2.10) and hence can be determined.

Now let us suppose that (in appropriate non-dimensional units) the external sources are small and vary slowly in time, and that the external currents vary slowly in time. Then it is possible to derive approximate equations describing the slow evolution of a plasma-vacuum system that is almost in equilibrium at every instant of time. In order to obtain these equations we introduce a dimensionless parameter ϵ << 1 and scale the unknowns as follows:

ρ = ρ̂, ζ = ζ̂, p = p̂, V = ϵV̂, B = B̂ , (2.11)

where the caretted unknowns depend on the scaled variables

x̂ = x, t̂ = ϵ t . (2.12)

Such a scaling is based on the assumption that the given external sources and currents can be expressed in the form

q_{M} = ϵq̂_{M}(x, ϵt), q_{S} = ϵq̂_{S}(x, ϵt), Q = ϵ^{2}Q̂(x, ϵt), J = Ĵ(x, ϵt) . (2.13)

Upon substituting of these expressions into the governing equations (2.1)–(2.6) and the associated fixed and free boundary conditions, and after dropping the carets, we find that all of those equations remain unchanged except for the force balance equation which assumes the form

ϵ^{2}ρ \frac{DV}{Dt} + ∇p − J × B = ϵ^{2}Q . (2.14)

Now formally neglecting the O(ϵ^{2}) terms in (2.14), we therefore conclude that the adiabatically slow evolution of the plasma-vacuum system caused by the presence of external sources and currents is governed by a set of reduced equations: the conservation laws (2.1), (2.2), the equilibrium equations (2.11) in the plasma and (2.6) in the vacuum, and the flux-freezing (induction) equation (2.4). In addition, we find that the fixed and free boundary conditions are unchanged.

3 Axisymmetric solutions and their constraints

We henceforth assume that the plasma-vacuum system is axisymmetric. The toroidal region D = {x = (r, ϕ, z) : (r, z) ∈ Ω, 0 ≤ ϕ < 2π} is then defined by its cross-section Ω in the usual cylindrical coordinates. The invariance of the system with respect to the toroidal angle ϕ implies that the slow evolution equations derived in §2 can be replaced by a simpler set of equations. For the sake of clarity in these considerations we assume that the velocity field is purely poloidal,

V = (V^{r}, 0, V^{z}) . (3.1)

(A modification which allows for toroidal flow of arbitrary magnitude can be made using the method that follows.) The axisymmetric slow evolution equations govern the unknowns ρ(r, z, t), ζ(r, z, t), p(r, z, t), Ψ(r, z, t) and f(r, z, t), where the magnetic field and current density are written in terms of flux functions Ψ and f according to

B = ∇Ψ × ∇ϕ + f∇ϕ (3.2) J = ∇f × ∇ϕ + (LΨ)∇ϕ , (3.3)

where

L: = −r \frac{∂}{∂r} (\frac{1}{r} \frac{∂}{∂r}) − \frac{∂^{2}}{∂z^{2}} .

In the plasma region these equations are

\frac{Dρ}{Dt} + ρ∇ · V = q_{M} (3.4)

\frac{Dζ}{Dt} + ζ∇ · V = q_{S} (3.5)

− r^{−2}f∇f + r^{−2}(LΨ)∇Ψ = ∇p, ∇f × ∇Ψ = 0 (3.6)

\frac{DΨ}{Dt} = 0 (3.7)

\frac{D}{Dt} (r^{−2}f) + r^{−2}f∇ · V = 0 ; (3.8)

and in the vacuum region they are

∇f = 0, LΨ = r\overline{J}_{ϕ} . (3.9)

Here ∇ = (∂/∂r, 0, ∂/∂z), and so ∇a×∇b can be identified with ∂(a, b)/∂(r, z). The external current density \overline{J} = \overline{J}_{ϕ}∇ϕ is prescribed to be purely toroidal; for instance, it may be realized as a finite collection of elementary (axisymmetric) current coils. The fixed and free boundary conditions can be formulated under axisymmetry as follows. The normalized boundary condition on the shell is taken to be

Ψ = 0 on ∂Ω . (3.10)

The plasma-vacuum interface conditions are expressible as

(3.11) p = 0, Ψ = σ_{0} on S, ∇Ψ, f are continuous across S;

the flux constant σ_{0} (say) then determines the magnetic surface S. The slow evolution problem for an axisymmetric plasma-vacuum system is thus fully described.

As explained in Part I, the (equilibrium) force balance equations (3.6) reduce to the Grad-Shafranov equation

LΨ = f(Ψ, t)f′(Ψ, t) + r^{2}p′(Ψ, t), f = f(Ψ, t), p = p(Ψ, t) , (3.12)

where prime denotes differentiation with respect to Ψ. Consequently, (3.12) can replace (3.6) in the above set of reduced equations. The profile functions f(Ψ, t) and p(Ψ, t) are called “surface quantities,” being constant on the magnetic surfaces {Ψ = σ} (σ ≥ σ_{0}).

By virtue of (3.7), each magnetic surface {Ψ = σ} moves with the plasma flow V so that the poloidal flux σ (between the magnetic surface and the shell) is conserved. Moreover, the toroidal flux, mass and entropy within each flux tube are constrained by the motion; these quantities are given by, respectively,

∫_{{Ψ>σ}} r^{−1}f drdz, ∫_{{Ψ>σ}} rρ drdz, ∫_{{Ψ>σ}} rζ drdz (σ ≥ σ_{0}) . (3.13)

The general classes of constraints of motion governed by axisymmetric ideal magnetohydrodynamics consist of functionals of the form

(3.14)  C_{F} = ∫_{Ω} r^{−1}fΦ(Ψ) drdz
(3.15)  C_{M} = ∫_{Ω} rρΦ(Ψ) drdz
(3.16)  C_{S} = ∫_{Ω} rζΦ(Ψ) drdz ,

where Φ(s) is any real function (with suitable regularity properties) that is supported in the interval σ_{0} < s < +∞. As is easily seen, constraints on these functionals are implied by corresponding constraints on the σ-parametrized families listed in (3.13); indeed, the obvious identity

∫_{Ω} aΦ(Ψ) drdz = ∫_{σ_{0}}^{+∞} Φ′(σ)dσ ∫_{{Ψ>σ}} a drdz ,

valid for any integrable function a = a(r, z), can be used to represent these functionals in terms of the quantities (3.13). The significance of the classes of constraints C_{F} , C_{M} and C_{S} (as Φ runs through a suitable class of functions) derives from the fact that each member of these classes can be evolved in time without (explicit) reference to the velocity field V. Specifically, the following identities hold:

(3.17)  \frac{dC_{F}}{dt} = 0
(3.18)  \frac{dC_{M}}{dt} = ∫_{Ω} rq_{M}Φ(Ψ) drdz
(3.19)  \frac{dC_{S}}{dt} = ∫_{Ω} rq_{S}Φ(Ψ) drdz .

Of course, these identities represent the conservation of toroidal flux, mass and entropy within the magnetic surfaces. Their verification is a straightforward consequence of the equations (3.4), (3.5), (3.7), (3.8). For example, the calculation needed to prove (3.18) is as follows:

\frac{dC_{M}}{dt} = ∫_{Ω} r [\frac{∂ρ}{∂t}Φ(Ψ) + ρΦ′(Ψ) \frac{∂Ψ}{∂t}] drdz
= ∫_{Ω} r [q_{M}Φ(Ψ) + ρΦ′(Ψ) (V · ∇Ψ + \frac{∂Ψ}{∂t})] drdz
= ∫_{Ω} rq_{M}Φ(Ψ) drdz .

In addition to the above constraints of motion within the plasma, there is another constraint of motion — namely,

F_{0} = ∫_{Ω} r^{−1}f drdz , (3.20)

which expresses conservation of the total (plasma plus vacuum) toroidal flux. In order to verify the corresponding identity

\frac{dF_{0}}{dt} = 0 , (3.21)

it is necessary to recall that the Maxwell equation

\frac{∂B}{∂t} + ∇ × E = 0

holds in D, while the tangential components of E vanish on ∂D. Then (3.21) follows immediately by applying Stokes’ formula.

4 Variational method for relaxed problems

The slow evolution equations for an axisymmetric plasma-vacuum system as given in §3 are degenerate in the sense that they do not include the convective derivative DV/Dt. This set of equations is therefore underdetermined with respect to the evolution of the velocity field V. Nevertheless, the velocity component V^{⊥} normal to the magnetic surfaces {Ψ = σ} in the plasma, can be determined from the poloidal flux convection equation (3.7); namely, V^{⊥} can be defined by

V^{⊥} = −|∇Ψ|^{−1} \frac{∂Ψ}{∂t} , (4.1)

thus making (3.7) valid pointwise everywhere that the normal n = |∇Ψ|^{−1}∇Ψ is itself defined. The conservation laws (3.8), (3.4) and (3.5) for toroidal flux, mass and entropy, on the other hand, have ambiguous meaning since each of them also involves the tangential velocity to the magnetic surfaces. A natural way to rectify this degeneracy is to relax the requirement that (3.8), (3.4) and (3.5) hold at every point in the plasma to the weaker requirement that toroidal flux, mass and entropy be conserved within every magnetic surface. By relaxing these equations in such a manner we achieve two goals. First, we obtain a self-consistent formulation of the slow evolution problem (involving only V^{⊥} = n · V). Second, we arrive at a problem which at each instant of time has the same variational structure as the equilibrium problem (P_{∞}) considered in Part I. Thus, we can invoke the method developed in Part I to solve the relaxed slow evolution problem.

In order to formulate the above mentioned relaxation we introduce the functionals

(4.2)  F_{σ} = ∫_{Ω} r^{−1}f (Ψ − σ)_{+} drdz
(4.3)  M_{σ} = ∫_{Ω} rρ (Ψ − σ)_{+} drdz
(4.4)  S_{σ} = ∫_{Ω} rζ (Ψ − σ)_{+} drdz

parameterized by the flux variable σ which runs through the range of Ψ in the plasma. (Henceforth we write simply σ ≥ σ_{0}, understanding that F_{σ}, M_{σ} and S_{σ} vanish when σ exceeds max Ψ.) Of course, these functionals are just the constraints of motion C_{F} , C_{M} and C_{S} corresponding to the particular choice Φ(s) = (s − σ)_{+} : = max(s − σ, 0). Moreover, they are identical with (minus) the σ-antiderivatives of the classical quantities displayed in (3.13). In terms of these constraints of motion we can express the relaxation of the conservation laws (3.8), (3.4) and (3.5), respectively, as the equations

(4.5)  \frac{dF_{σ}}{dt} = 0
(4.6)  \frac{dM_{σ}}{dt} = ∫_{Ω} rq_{M} (Ψ − σ)_{+} drdz
(4.7)  \frac{dS_{σ}}{dt} = ∫_{Ω} rq_{S} (Ψ − σ)_{+} drdz

Also, on account of the free boundary we impose the additional equation

\frac{dF_{0}}{dt} = 0 , (4.8)

where F_{0} is the constraint of motion defined in (3.20). The precise formulation of the relaxed slow evolution problem can now be stated: the plasma region is governed by the (equilibrium) force balance equations (3.12), the magnetic surface convection equation (3.7), the integral conservation laws (4.5)–(4.8), and the equation of state (2.8); the vacuum region is governed by the field equations (3.9); and, the fixed and free boundary conditions are given by (3.10) and (3.11), respectively.

The above relaxation from (pointwise) conservation laws to (integral) constraint of motion can be interpreted in terms of averaging over magnetic surfaces. In the notation adopted by Grad et al., let

⟨a⟩: = ∫_{\{Ψ=σ\}} \frac{ra}{|∇Ψ|}dℓ

denote the surface average of a function a = a(r, z). Then we claim that (4.5)–(4.7) are equivalent to the equations

\frac{∂}{∂t}⟨r^{−2}f⟩ = 0, \frac{∂}{∂t}⟨ρ⟩ = ⟨q_{M}⟩, \frac{∂}{∂t}⟨ζ⟩ = ⟨q_{S}⟩ ,

respectively, where the surface averaged quantities involved depend upon σ and t. Indeed, recalling the relevant discussion given in §1 of Part I, we find that these equations are identical with the evolution equations for F_{σ}, M_{σ}, S_{σ} after applying ∂^{2}/∂σ^{2} to each of the identities (4.5)–(4.7) in σ. Consequently, we conclude that the relaxed slow evolution problem is identical with the surface averaged Grad-Hogan equations. We prefer however to replace the concept of surface averaging by the equivalent concept of constraints of motion in order to expose the natural variational structure of the problem.

The justification on physical grounds for averaging over magnetic surfaces can be summarized as follows. The plasma region is foliated by toroidal magnetic surfaces on which the helical field lines wind so that, at least generically, each field line is dense in its associated surface. Therefore, transport processes (such as *****) beyond the scope of the governing equations under consideration effectively enforce the postulated averaging in reality. On this basis the relaxed (averaged) equations constitute a realistic model of the adiabatically slow evolution under study. On the same grounds we shall assume that the external sources are prescribed as surface quantities q_{M} = q_{M}(Ψ, t) and q_{S} = q_{S}(Ψ, t), although this simplification is not strictly necessary.

We now proceed to give a variational method of solving the relaxed slow evolution problem. The prescribed data for this problem consist of the external sources q_{M}(s, t) and q_{S}(s, t), and the flux function \overline{ψ}(r, z, t) for the poloidal field induced by the external current density \overline{J}_{ϕ}(r, z, t) according to

L\overline{ϕ} = r \overline{J}_{ϕ} in Ω, \overline{ψ} = 0 on ∂Ω . (4.9)

Also, appropriate initial data are to be prescribed, but this issue is deferred for now. The total flux function Ψ can be split into the sum Ψ = ψ + \overline{ψ}, where ψ is obviously the flux function induced by current density supported in the plasma region. In what follows, ψ is used as an unknown rather than Ψ, since \overline{ψ} is prescribed. In the present context the variational principle which characterizes the (equilibrium) force balance at every instant of time is a simple extension of the variational principle developed in Part I. The total energy

E(ψ, f, ρ, ζ) = \frac{1}{2} ∫_{Ω}[r^{−1}|∇ψ|^{2} + r^{−1}f^{2} + r U(ρ, ζ)] drdz (4.10)

serves as the objective functional; the constraint functionals are supplied by the constraints of motion F_{0}, F_{σ}, M_{σ}, S_{σ} (σ ≥ σ_{0}) defined above, which now can be considered as functionals of (ψ, f, ρ, ζ) depending explicitly on t through \overline{ψ}(t). At every instant of time t we let (ψ^{∗}, f^{∗}, ρ^{∗}, ζ^{∗}) denote the solution of the constrained minimization problem

E(ψ, f, ρ, ζ) → min subject to
F_{0}(f) = F_{0}^{∗}, F_{σ}(ψ, f; t) = F_{σ}^{∗}
M_{σ}(ψ, ρ; t) = M_{σ}^{∗}, S_{σ}(ψ, ζ; t) = S_{σ}^{∗} (σ ≥ σ_{0})

(4.11)

corresponding to (instantaneous) constraint values F_{0}^{∗}, F_{σ}^{∗}, M_{σ}^{∗}, S_{σ}^{∗}. In doing so we assume that the minimizer for (4.11) is unique (at least locally along a trajectory in the solution space), even though a general uniqueness theorem is not available. Since the validity of this assumption can be verified computationally, we can consider the formal solution map

(F_{0}^{∗}, F_{σ}^{∗}, M_{σ}^{∗}, S_{σ}^{∗}) ⟶ (ψ^{∗}, f^{∗}, ρ^{∗}, σ^{∗}) (4.12)

to be well defined for the purposes of discussion. With this map in hand we are able to pose the relaxed slow evolution problem as a (σ-parametrized) family of ordinary differential equations:

\{ \begin{matrix}\frac{dF_{0}^{∗}}{dt} & = & 0, \frac{dF_{σ}^{∗}}{dt} = 0 \\ \frac{dM_{σ}^{∗}}{dt} & = & _{∫Ω}rq_{M} (ψ^{∗} + \overline{ψ}, t) (ψ^{∗} + \overline{ψ} − σ)_{+} drdz \\ \frac{dS_{σ}^{∗}}{dt} & = & _{∫Ω}rq_{S} (ψ^{∗} + \overline{ψ}, t) (ψ^{∗} + \overline{ψ} − σ)_{+} drdz , \end{matrix} (4.13)

where σ ≥ σ_{0} (actually σ runs through the range of Ψ = ψ + \overline{ψ}, which is invariant in time). In principle, the (slowly evolving) solution trajectory (ψ^{∗}(t), f^{∗}(t), ρ^{∗}(t), ζ^{∗}(t)) can be advanced in time according to these equations.

The formal variational principle (4.11) is very similar to (P_{∞}) introduced in Part I, which also is a minimization problem with (continuously) infinite families of constraints. As we stress in Part I, we prefer to formulate a more relaxed version of this kind of variational problem in which a finite family of constraints is imposed. The resulting multiconstrained variational problem is then amenable to standard methods of analysis and computation. For this purpose we follow exactly the construction given in Part I. Let σ_{0} < σ_{1} < ⋯ < σ_{n} < +∞ be a partition of the interval σ_{0} ≤ σ < +∞, and let Δσ_{i}: = σ_{i} − σ_{i−1}. Relative to this partition, let

Φ_{i}(s) : = \frac{1}{Δσ_{i}} ∫_{σ_{i=1}}^{σ_{i}} (s − σ)_{+} dσ − \frac{1}{Δσ_{i+1}} ∫_{σ_{i}}^{σ_{i+1}} (s − σ)_{+} dσ (4.14)

be a finite family of “basis functions” (i = 1, …, n). These particular functions have the property that their s-derivatives Φ′_{i}(s) are the usual finite element functions. Now recalling the general classes of constraints of motion given in (3.14)–(3.16), we define the functionals

(4.15)  F_{i}(ψ, f; t) = ∫_{Ω} r^{−1}fΦ_{i}(ψ + \overline{ψ}) drdz
(4.16)  M_{i}(ψ, ρ; t) = ∫_{Ω} rρΦ_{i}(ψ + \overline{ψ}) drdz
(4.17)  S_{i}(ψ, ζ; t) = ∫_{Ω} rζΦ_{i}(ψ + \overline{ψ}) drdz ,

where the explicit dependence on t enters through \overline{ψ}. The finitely constrained minimization problem corresponding to (4.11) can be stated as

E(ψ, f, ρ, ζ) → min subject to

F_{0}(f) = F_{0}^{∗}, F_{i}(ψ, f; t) = F_{i}^{∗}
M_{i}(ψ, ρ; t) = M_{i}^{∗}, S_{i}(ψ, ζ; t) = S_{i}^{∗} (i = 1, …, n) ,

(4.18)

where the constraint values are derived from those in (4.11) according to

F_{i}^{∗} = \frac{1}{Δσ_{i}} ∫_{σ_{i−1}}^{σ_{i}} F_{σ}^{∗} dσ − \frac{1}{Δσ_{i+1}} ∫_{σ_{i}}^{σ_{i+1}} F_{σ}^{∗} dσ

and similarly for M_{i}^{∗} and S_{i}^{∗}. In turn, there is a finite dimensional system of ordinary differential equation corresponding to (4.13). Introducing the vector of constraint values

X^{∗}: = (F_{0}^{∗}, F_{i}^{∗}, M_{i}^{∗}, S_{i}^{∗}) ∈ R^{3n+1} ,

this system can be written as

\frac{dX^{∗}}{dt} = A(X^{∗}; t) ; (4.19)

the nonlinear operator A, which involves the solution map for (4.18), has the 3n + 1 components (i = 1, …, n)

A_{0} = 0 , A_{i} = 0
(4.20)  A_{n+i} = ∫_{Ω} rq_{M} (ψ^{∗} + \overline{ψ}, t) Φ_{i} (ψ^{∗} + \overline{ψ}) drdz
A_{2n+i} = ∫_{Ω} rq_{S} (ψ^{∗} + \overline{ψ}, t) Φ_{i} (ψ^{∗} + \overline{ψ}) drdz

Here (ψ^{∗}, f^{∗}, ρ^{∗}, ζ^{∗}) is the minimizer for (4.18) corresponding to constraint values X^{∗}; again it is assumed to be uniquely determined by X^{∗}.

The above discretization of the σ-parametrized family of constraints of motion may be viewed as a specific kind of averaging between the magnetic surfaces {Ψ = σ_{i}} (i = 0, 1, …, n). Consequently, the proposed solutions governed by (4.18), (4.19) satisfy the conservation laws for r^{−2}f, ρ and σ only in the sense of a volume average over {σ_{i−1} < Ψ < σ_{i}} (i = 1, …, n). Thus, the slow evolution problem is relaxed even further. On the other hand, the force balance equations hold exactly for any finite partition {σ_{i}}. This important fact can be demonstrated by calculating the Lagrange multiplier rule for (4.18). The resulting equations are (omitting stars)

(4.21)  Lψ = fλ · Φ′ (ψ + \overline{ψ}) + r^{2}[ρμ · Φ′(ψ + \overline{ψ}) + ζν · Φ′ (ψ + \overline{ψ})]
(4.22)  f = λ_{0} + λ · Φ (ψ + \overline{ψ})
(4.23)  \frac{∂U}{∂ρ} (ρ, ζ) = μ · Φ (ψ + \overline{ψ})
(4.24)  \frac{∂U}{∂ζ} (ρ, ζ) = ν · Φ (ψ + \overline{ψ}) ,

where λ_{0}, λ_{i}, μ_{i}, ν_{i} (i = 1, …, n) are the multipliers associated with the solution (ψ, f, ρ, ζ); here the notation λ · Φ = ∑_{i=1}^{n}λ_{i}Φ_{i} is used. In fact, these variational equations can be reduced to the Grad-Shafranov equation (3.12). The poloidal current profile f(Ψ) is defined by (4.22). The profiles ρ(Ψ) and ζ(Ψ) are determined by solving (4.23) and (4.24), which is possible because of the strict convexity of U(ρ, ζ); then the pressure profile p(Ψ) is defined and satisfies

p′(Ψ) = (ρ∂U/∂ρ + ζ∂U/∂ζ − U)′
= ρ(∂U/∂ρ)′ + ζ(∂U/∂ζ)′
= ρμ · Φ′(Ψ) + ζν · Φ′(Ψ) .

From this relation the equivalence of (4.21) with (3.12) is immediate. For the purpose of reference the specific profiles involved here are recorded:

ρ(Ψ) = \{ν · Φ(Ψ) exp (\frac{μ · Φ(Ψ)}{ν · Φ(Ψ)} − γ)\}^{1/(γ−1)}

ζ(Ψ) = \{γ − \frac{μ · Φ(Ψ)}{ν · Φ(Ψ)}\} ρ(Ψ) .

The variational formulation of the relaxed slow evolution problem described above has the virtue that it applies in a general setting. However, we now give a simplified version which is less involved and hence more suited to numerical simulation. Namely, we assume that the plasma is isentropic and homentropic, meaning that

ζ = s_{0}ρ and q_{S} = s_{0}q_{M} (4.25)

for a constant s_{0} (independent of x and t); here we also ensure that the sources are consistent with the assumption. The equation of state for such a plasma is polytropic: p = e^{s_{0}}ρ^{γ}. As a consequence, the two constraints of motion M_{σ} and S_{σ} can be amalgamated into one:

\frac{d}{dt} ∫_{Ω} rp^{1/γ} (Ψ − σ)_{+}drdz = e^{s_{0}/γ} ∫_{Ω} rq_{M} (Ψ − σ)_{+}drdz . (4.26)

It is then straightforward to derive a corresponding version of the relaxed slow evolution equations (4.13). For this purpose, we recall the unknown

g: = (\frac{2p}{γ − 1})^{\frac{1}{2}} , (4.27)

introduced for technical convenience in Part I, and we define the functional

G_{σ}(ψ, g) = ∫_{Ω} rg^{2/γ}(ψ + \overline{ψ} − σ)_{+}drdz (4.28)

which (up to a constant factor) is identical with the constraint of motion occurring in (4.26). Then we can pose the variational problem governing the solution at every instant of time in the form

(P_{∞}) \{ \begin{matrix}E(ψ, f, g) → min subject to \\ F_{0}(f) = F_{0}^{∗}, F_{σ}(ψ, f) = F_{σ}^{∗}, G_{σ}(ψ, g) = G_{σ}^{∗} (σ ≥ σ_{0}) , \end{matrix}

with the quadratic energy functional

E(ψ, f, g) = \frac{1}{2} ∫_{Ω} [r^{−1}|∇ψ|^{2} + r^{−1}f^{2} + rg^{2}] drdz . (4.29)

Similarly, we can state the multiconstrained variational problem corresponding to (P_{∞}) in the form

(P_{n}) \{ \begin{matrix}E(ψ, f, g) → min subject to \\ F_{0}(f) = F_{0}^{∗}, F_{i}(ψ, f) = F_{i}^{∗}, G_{i}(ψ, g) = G_{i}^{∗} (i = 1, …, n) , \end{matrix}

employing the obvious notation. Of course, (P_{∞}) and (P_{n}) are precisely the same (at a given instant of time) as the equilibrium variational problems treated in Part I. Therefore, the globally convergent iterative algorithm developed in §3 and §4 of Part I can be invoked in the present context to define the solution map (for fixed t)

X^{∗} = (F_{0}^{∗}, F_{i}^{∗}, G_{i}^{∗}) ⟶ (ψ^{∗}, f^{∗}, g^{∗}) .

As a result, this simplified version of relaxed slow evolution problem is completely solved in principle, since the constraint values X^{∗} are advanced in time according to the obvious modification of the system (4.19).

It remains to discuss the initial conditions for the system (4.19) or its modifications. In most realistic situations, the external sources or currents that drive the slow evolution are applied to a known equilibrium configuration of interest. For instance this is the case when either plasma heating on adiabatic compression is initiated. Therefore, it is reasonable to assume that the initial constraint values are derived from a given equilibrium solution as found in Part I.

References

[1] G. Bateman, MHD Instabilities, MIT press, Cambridge, Massachusetts, 1978.

[2] J. Blum, Numerical Simulation and Optimal Control in Plasma Physics, with Applications to Tokomaks, Wiley/Gauthiers-Villars, Chichester/New York, 1989.

[3] R.A. Dory and Y.-K. M. Peng, High pressure flux-preserving tokomak equilibria, Nucl. Fusion 17 (1977), 21-31.

[4] A. Eydeland, J. Spruck and B. Turkington, Multiconstrained variational problems of nonlinear eigenvalue type: new formulations and algorithms, Math. of Comp., in press.

[5] J. P. Freidberg, Ideal Magnetohydrodynamics, Plenum Press, New York, 1987.

[6] H. Grad, Magnetic confinement fusion energy research, Proceeding of Symposia in Applied Math. 21 (1977), 3-40.

[7] H. Grad, Survey of 1\frac{1}{2} D transport codes, Courant Institute of Mathematical Sciences Report MF-93, New York University, 1978 (unpublished).

[8] H. Grad, P.N. Hu and D.C. Stevens, Adiabatic evolution of plasma equilibria, Proc. Nat. Acad. Sci. 72 (1975), 3789-3793.

[9] H. Grad, P.N. Hu, D.C. Stevens and E. Turkel, TITLE?, Plasma Phys. Contr. Nucl. Fusion Res. 2 (1976), 355-365.

[10] H. Grad and H. Rubin, Hydromagnetic equilibria and force-free fields, Proc. of 2nd U. N. International Conf. on the Peaceful Uses of Atomic Energy, U.N., Geneva (1958), vol. 31, 190-197.

[11] M. D. Kruskal and R. M. Kulsrud, Equilibrium of a magnetically confined plasma in a toroid, Physics of Fluids 1(1958), 265–274.

[12] A. E. Lifschitz, Magnetohydrodynamics and Spectral Theory, Kluwer Academic Publishers, Dordrecht, 1989.

[13] L. Woltjer, Hydromagnetic equilibrium,III: axisymmetric incompressible media, Astrophysical J. 130 (1959), 400-404.

[14] ———-, Hydromagnetic equilibrium,IV: axisymmetric compressible media, Astrophysical J. 130 (1959), 405-413.