Diffusion and Transport in Axisymmetric Geometry Stephen C. Jardin - - PowerPoint PPT Presentation

MHD Simulations for Fusion Applications Lecture 2

Diffusion and Transport in Axisymmetric Geometry

Stephen C. Jardin Princeton Plasma Physics Laboratory

CEMRACS ‘10 Marseille, France July 20, 2010

1

10-10 10-2 104 100

SEC.

CURRENT DIFFUSION

10-8 10-6 10-4 102 ωLH

• 1

Ωci

• 1

τA Ωce

• 1

ISLAND GROWTH ENERGY CONFINEMENT SAWTOOTH CRASH TURBULENCE ELECTRON TRANSIT

(a) RF codes (b) Micro- turbulence codes (c) Extended- MHD codes (d) Transport Codes

These 4 areas address different timescales and are normally studied using different codes

Transport codes solve the same equations as the Extended MHD codes, but in 2D rather than 3D

3

(d) Transport Codes R Z φ

φ ∂ = ∂

axisymmetry First we will consider a model problem to better understand the timescales in 2D, then will proceed to a general approach

R Z φ Plasma Cross Section (toroidal current in) A tokamak needs an externally generated “vertical field” for equilibrium. A purely vertical field will produce a nearly circular cross-section plasma.

4

Tokamak Equilibrium Basics: Need for a vertical field

We first consider a very crude “rigid plasma” model to better understand where the slow resistive time scale comes from

R Z φ R Z φ Good Curvature

• Stable to vertical mode
• Oblate plasma
• low beta limits

Bad Curvature

• Unstable to vertical mode
• Elongated plasma
• higher beta limits

5

In actual tokamak experiments, external poloidal field is not purely straight but has some curvature to it

R Z φ R Z φ Good Curvature

• Stable to vertical mode
• Oblate plasma
• low beta limits

Bad Curvature

• Unstable to vertical mode
• Elongated plasma
• higher beta limits

6

In actual tokamak experiments, external field is not purely straight but has some curvature to it

R Z φ R φ Z External poloidal field with curvature can be thought

• f as a superposition of vertical and radial field.

If plasma column is displaced upward, the force J x B = IP x BRext will accelerate it further upward. Same for downward. Alfven wave time scale: very fast! Vertical Instability

7

R φ Z Plasma with current IP Describe the plasma as a rigid body

• f mass m with Z position ZP.

Assume time dependence eiωt ZP Equation of motion: Circuit equation for wall:

8

A nearby conductor will produce eddy currents which act to stabilize Conducting wall with dipolar current IC

• +

inertia conductor external field inductance resistance plasma coupling

R φ Z Plasma with current IP Conducting wall with dipolar current IC ZP Three roots:

9

Introduce plasma velocity VP = i ωZP to get a 3x3 matrix eigenvalue equation for ω of standard form +

Three roots: These are high frequency (~10-7 sec) stable oscillations that are slowly damped by the wall resistivity This is the unstable mode. Very slow (~ 10-1 sec), and independent of plasma mass.

10

With only passive conductor, still an unstable root but much smaller. Not on Alfven wave time scale but on L/R timescale of conductor.

R φ Z Plasma with current IP Conducting Wall with current IC ZP Three roots: This “rigid” mode is easily stabilized by adding a pair of feedback coils of

• pposite sign, and applying a voltage

proportional to the plasma displacement

• or its time integral or time derivative (PID)

11

Total stability is obtained by adding an active feedback system which only needs to act on this slower timescale.

To model this “vertical instability” in realistic geometry, and take the non- rigid motion of the plasma into account, we take advantage of the fact that the unstable mode does not depend on the plasma mass (or inertia), and the stable modes are very high frequency and very low amplitude.

2

( ) ( ) 3 3 2 2 ( , )

i i i ij j i i i i i j P

n n t t nM p t p p p J t d L I M I J G R R dR R I V dt

φ

η η

≠

∂ + ∇ = ∂ ∂ = −∇× ∂ ∂ +

• ∇

+ ∇ = × ∂ + × = ∂ ⎛ ⎞ + ∇ + = − ∇ + ⎜ ⎟ ∂ ⎝ ⎠ = ∇× ⎡ ⎤ + + + = ⎢ ⎥ ⎣ ⎦

∑ ∫

V B E V V V J B E V B J q V V J B i i i

We start with the basic MHD + circuit equations and apply a “resistive timescale ordering” Introduce small parameter 1 ε

~ ~ ~ ~

i

V R t η ε ∂ ∂ V E ∼ ∼

12

To model this “vertical instability” in realistic geometry, and taking the non-rigid motion of the plasma into account, we take advantage of the fact that the unstable mode does not depend on the plasma mass (or inertia), and the stable modes are very high frequency and low amplitude.

2

( ) ( ) 3 3 2 2 ( , )

i i i ij j i i i i i j P

n n t t nM p t p p p J t d L I M I J G R R dR R I V dt

φ

η η

≠

∂ + ∇ = ∂ ∂ = −∇× ∂ ∂ +

• ∇

+ ∇ = × ∂ + × = ∂ ⎛ ⎞ + ∇ + = − ∇ + ⎜ ⎟ ∂ ⎝ ⎠ = ∇× ⎡ ⎤ + + + = ⎢ ⎥ ⎣ ⎦

∑ ∫

V B E V V V J B E V B J q V V J B i i i

We start with the basic MHD + circuit equations and apply a “resistive timescale ordering” Introduce small parameter 1 ε

~ ~ ~ ~

i

V R t η ε ∂ ∂ V E ∼ ∼

2

ε

All equations pick up a factor of , in all terms, which cancels out, except in the momentum equation, where the inertial terms are multiplied by .

ε

2

ε

13

To model this “vertical instability” in realistic geometry, and taking the non-rigid motion of the plasma into account, we take advantage of the fact that the unstable mode does not depend on the plasma mass (or inertia), and the stable modes are very high frequency and low amplitude.

2

( ) 3 3 2 2 ( , ) ( )

i i ij j i i i i i P i j

n n t t p p p p J t d L I M I J G R nM R dR R I V dt t

φ

η η

≠

∂ + ∇ = ∂ ∂ = −∇× ∂ ∇ = × + × = ∂ ⎛ ⎞ + ∇ + = − ∇ + ⎜ ⎟ ∂ ⎝ ⎠ = ∇× ⎡ ⎤ + + + = ⎢ ⎥ ⎣ ∂ +

• ∂

⎦ ∇ +

∑ ∫

V B E J B E V V B J q V V V V J B i i i

We start with the basic MHD + circuit equations and apply a “resistive timescale ordering” Introduce small parameter 1 ε

~ ~ ~ ~

i

V R t η ε ∂ ∂ V E ∼ ∼

2

ε

This allows us to drop the inertial terms in the momentum equation, and replace it with the equilibrium equation. Huge simplification…. removes Alfven timescale

14

To model this “vertical instability” in realistic geometry, and taking the non-rigid motion of the plasma into account, we take advantage of the fact that the unstable mode does not depend on the plasma mass (or inertia), and the stable modes are very high frequency and low amplitude.

2

( ) 3 3 2 2 ( , )

i i ij j i i i i i j P

n n t t p p p p J t d L I M I J G R R dR R I V dt

φ

η η

≠

∂ + ∇ = ∂ ∂ = −∇× ∂ ∇ = × + × = ∂ ⎛ ⎞ + ∇ + = − ∇ + ⎜ ⎟ ∂ ⎝ ⎠ = ∇× ⎡ ⎤ + + + = ⎢ ⎥ ⎣ ⎦

∑ ∫

V B E J B E V B J q V V J B i i i

This is the set of equations we solve to simulate control of the plasma position and shape. There are 3 production codes that solve these nonlinear equations in 2D and are used to design and test control strategies.

• TSC

(PPPL)

• DINA (Russia)
• CORSICA (LLNL)

15

16

(1) (2) t η ∂ = −∇× ∂ + × = B E E V B J

φ R Z Consider first the vector magnetic field equation The most general form for an equilibrium axisymmetric magnetic field is:

( )

(3) g φ φ = ∇ ×∇Ψ + Ψ ∇ B

(4) g t t φ φ ∂Ψ ∂ ∇ ×∇ + ∇ = −∇× ∂ ∂ E

Substitute (3) into (1): Take dot product of (4) with

φ ∇

[ ]

2

1 (5) g R t φ φ ∂ = −∇ ∇× = ∇ ∇ × ∂ E E i i

Noting that , the remaining part of (4) becomes

t t φ φ ∂Ψ ∂Ψ ∇ ×∇ = −∇× ∇ ∂ ∂

2

(6) ( ) R C t t φ ∂Ψ = ∇ + ∂ Ei

Constant can be taken to vanish to match boundary condition that Ψ=0 at R=0

2 2

is "flux function" g is toroidal field function =0 1 0, R φ φ Ψ ∇ ∂ = ∇ = ∂ B i

(2) η + × = E V B J

φ R Z Recall:

[ ]

2

1 (5) g R t φ ∂ = ∇ ∇ × ∂ E i

2

(6) R t φ ∂Ψ = ∇ ∂ Ei

( )

*

g g φ φ μ φ φ = ∇ ×∇Ψ + Ψ ∇ = ∇× = Δ Ψ∇ + ∇ ×∇ B J B

Use (2) to eliminate E from (5) and (6):

( ) ( )

2 2 2 2

1 1 (7) g g R g t t R R R g η φ φ φ η μ φ ⎡ ⎤ ∂ = ∂ ⎡ ⎤ = ∇ −∇ × × ∇ − + ∇ ∇ ×∇Ψ + ∇ ⎢ ⎥ ∂ ⎣ ⎦ + ∇ × ⎣ ⎦ ∂ V B J V V i i i

( )

2 * 2

(8) R t t R φ η φ η μ ∂Ψ = ∂Ψ = − − × ∇ + ∇ + Δ Ψ ∂ ∂ ∇Ψ V V B J i i i

* 2 2

R R− Δ Ψ ≡ ∇ ∇Ψ i

φ R Z

( )

g φ φ μ = ∇ ×∇Ψ + Ψ ∇ = ∇× B J B

( )

2 * 2 2 2

( ) 3 3 2 2 1 And, the equilibrium constraint: n n t p p p t t g g R g t R R p η η μ η φ φ μ ∂ + ∇ = ∂ ∂ ⎛ ⎞ + ∇ + + ∇ = ⎜ ⎟ ∂ ⎝ ⎠ ∂Ψ + ∇Ψ = Δ Ψ ∂ ⎡ ⎤ ∂ + ∇ − ∇ ∇ ×∇Ψ − ∇ = ⎢ ⎥ ∂ ⎣ ⎦ ∇ = × V q V V J V V V J B i i i i i i

Summary of scalar equations: Note:

• The pressure and magnetic field variables obey separate time

advancement equations, yet they must always satisfy the equilibrium constraint

• Each of the equations contains the velocity variable V, yet there is no

equation to advance V.

• The heat flux vector q is very anisotropic, much larger parallel to the

magnetic field than perpendicular to it.

19

Because of the anisotropy of the heat conduction, we want to transform to a moving coordinate system aligned with the magnetic flux surfaces. At any given time, we will define the non-

• rthogonal flux coordinate system:

( )

( ), ( ), ( ) ψ θ φ x x x We also have the inverse representation: ( , , , ) t ψ θ φ x This has the associated volume element: and Jacobian:

[ ]

1

d d d d Jd d d J ψ θ φ τ ψ θ φ ψ θ φ ψ θ φ

−

= ≡ ≡ ∇ ×∇ ∇ ∇ ×∇ ∇ i i We next define the coordinate velocity at a particular location as: ( , , ) ψ θ φ

, , C

t ψ θ φ ∂ = ∂ x u

, , , , , , C

t t t t t

ψ θ φ ψ θ φ ψ θ φ

α α α ∂ ∂ ∂ ∂ ∂ ∂ = + ⇒ = − ∇ ∂ ∂ ∂ ∂ ∂ ∂

x x

x u x i i For any scalar function α: Also, one can verify the relation for the time derivative of the Jacobian:

, , C

J J t ψ θ φ ∂ = ∇ ∂ u i

( )

ψ ψ ψ = Ψ ∇ = Bi

20

, , C

t t ψ θ φ ∂ ∂ = − ∇ ∂ ∂

x

u i

[ ]

1 , , C

J J J t ψ θ φ ψ θ φ

−

∂ ≡ ∇ ×∇ ∇ = ∇ ∂ u i i The fluid velocity that appears in the MHD equations is divided into two parts:

C R

= + V u u

Velocity of the coordinates Velocity relative to the moving coordinates actual fluid velocity

( ) ( )

, , , ,

( ) ( ) ( )

C C R C R R

n n t n n n n t nJ J n t

ψ θ φ ψ θ φ

∂ + ∇ = ∂ ∂ − ∇ + + ∇ + ∇ + = ∂ ∂ + ∇ = ∂

x

V u u u u u u i i i i i

Transform the continuity equation to the moving frame: Time dependent coordinate transformation

( , ), ( , ), t t ψ θ φ x x

, , C

t ψ θ φ ∂ = ∂ x u

, , C

t t ψ θ φ ∂ ∂ = − ∇ ∂ ∂

x

u i

[ ]

1 , , C

J J J t ψ θ φ ψ θ φ

−

∂ ≡ ∇ ×∇ ∇ = ∇ ∂ u i i

( ) ( )

( )

3/5 3/5 2/5 2 * 2 2 2

( ) ( ) 2 5 1

R R R R R

nJ J n t p J J p Jp t t J g g J g t R R R η η μ η φ φ μ

−

∂ + ∇ = ∂ ∂ ⎡ ⎤ + ∇ + ∇ − = ⎣ ⎦ ∂ ∂Ψ + ∇Ψ = Δ Ψ ∂ ⎡ ⎤ ∂ ⎛ ⎞ + ∇ − ∇ ∇ ×∇Ψ − ∇ = ⎢ ⎥ ⎜ ⎟ ∂ ⎝ ⎠ ⎣ ⎦ u u q J u u u i i i i i i

Apply the same technique to all the time-advance equations: These equations are now in the moving frame. Because these are all conservation equations, only the relative velocity appears in the equations!

, ,

t t

ψ θ φ

∂ ∂ → ∂ ∂

( , ), ( , ), t t ψ θ φ x x

( ) ( ) ( )

R R

nJ J n nV nV t t ψ ψ ∂ ∂ ∂ ′ ′ ⎡ ⎤ + ∇ = → + ∇ = ⎣ ⎦ ∂ ∂ ∂ u u i i

Next, we integrate the equations in the poloidal angle around the flux surface. Use the property, that for any vector A

( ) ( ) ( ) ( )

2 2 2

2 2 2 J J J J d J d J V

π π π

ψ θ ψ θ π θ π ψ θ π θ ψ θ ψ ψ ∂ ∂ ∇ = ∇ + ∇ ∂ ∂ ∂ ∂ ∇ = ∇ + ∇ ∂ ∂ ∂ ′ = ∇ ⎡ ⎤ ⎣ ⎦ ∂

∫ ∫ ∫

A A A A A A A i i i i i i i Here, we have defined the differential volume and surface average:

2 2

2 ( ) 2 dV V J d a Jad d V

π π

π ψ π θ θ ψ ′ ≡ = = ′

∫ ∫

Apply to continuity equation:

( , ), ( , ), t t ψ θ φ x x

( )

( )

( )

3/5 3/5 2/5 2 * 2 2 2

(* ( ) 2 5 ) 1

R R R R

nV nV t p V p V p V V t t gV R gV R V g t R ψ ψ ψ ψ η ψ ψ η μ η ψ ψ ψ μ

− − −

∂ ∂ ′ ′ ⎡ ⎤ + ∇ = ⎣ ⎦ ∂ ∂ ⎡ ⎤ ∂ ∂ ∂ ′ ′ ′ ′ ⎡ ⎤ + ∇ + ∇ − = ⎢ ⎥ ⎣ ⎦ ∂ ∂ ∂ ⎣ ⎦ ∂Ψ + ∇Ψ = Δ Ψ ∂ ⎡ ⎤ ∂ ∂ ′ ′ ′ + ∇ − ∇ ∇ = ⎢ ⎥ ∂ ∂ ⎣ ⎦ u u q J u u i i i i i i

After integrating all the equations over a flux surface: These are now 1-dimensional equations for the surface averages. We can use one of these equations to eliminate the relative velocity from the others. Note that Equation (*) is for the derivative of the toroidal magnetic flux inside a flux surface

( , ), ( , ), t t ψ θ φ x x

24

( ) [ ]

2

1 2 1 2 1 2 d d d d

φ π φ

ψ φ τ π φ τ π φ π

= =

Φ = ∇ = ∇ = =

∫ ∫ ∫ ∫

B B B S B S i i i i

φ

Toroidal flux is the integral of the toroidal magnetic field over the area inside a surface

( )

2 2 2

(*) 1

R

gV R gV R V g t R η ψ ψ ψ μ

− −

⎡ ⎤ ∂ ∂ ′ ′ ′ + ∇ − ∇ ∇ = ⎢ ⎥ ∂ ∂ ⎣ ⎦ u i i

Equation (*) is for the derivative of the toroidal magnetic flux within a flux surface. If we choose the relative velocity so as to make the time derivative vanish, we can identify the flux coordinate as the toroidal flux.

R

u

( )

g φ φ = ∇ ×∇Ψ + Ψ ∇ B

( )

( )

2 2 2

1 2 d d d JgR d gV R

ψ π ψ

ψ φ τ θ ψ π ψ

− −

Φ = ∇ = ′ =

∫ ∫ ∫ ∫

Bi

ψ → Φ const. ψ =

25

* *

( )

R R

t f φ ψ η φ ψ μ η ψ μ ∂Ψ ∇ ×∇ = ∂ ⎡ ⎤ ∇ ×∇ ∇Ψ − Δ Ψ = ⎢ ⎥ ⎣ ⎦ ∇Ψ − Δ Ψ = u u i i i i There is a constraint on the relative velocity in that the coordinate ψ must remain a flux coordinate as time involves Presently undetermined We now determine this function by requiring that the flux coordinate ψ be the toroidal flux inside a surface

2 2 2 2 * 2 2 * 2 2 2 2

1 1 ( ) 1 1 ( )

R

gV R V g R gV R f gV R V g R g gV R V g gR R f gV R R η ψ ψ μ η η ψ ψ μ μ η η ψ μ μ η ψ μ

− − − − − −

′ ′ ∇ − ∇ ∇ = ′ ′ ′ + Δ Ψ − ∇ ∇ = ′ ′ ∇ ∇Ψ Δ Ψ − ∇ ∇ = − = − ′ u i i i i i

26

*

( )

R

f η ψ μ ∇Ψ − Δ Ψ = u i

2 2

1 ( ) g gR f R η ψ μ

−

∇ ∇Ψ = − i

2 * 2

1 (**)

R

g gR R η μ

−

⎡ ⎤ ∇ ∇Ψ ⎢ ⎥ ⎢ ⎥ ⇒ ∇Ψ = Δ Ψ − ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ u i i Recall: Put these together: Relative velocity is determined by the equations themselves once we identify the physical meaning of the flux coordinate. This says that if the toroidal magnetic flux within a magnetic surface is used as the flux coordinate, this Equation (**) gives the fluid velocity relative to the flux surface velocity. We can therefore eliminate this velocity from each of the previous equations for the density, pressure, poloidal flux.

( ) ( )

2 2

define: 1 1 2

R L

g gR V R ψ η ψ π μ

−

Γ ≡ ∇Ψ ∇ ∇Ψ ≡ u i i

ψ → Φ

[ ]

( )

( )

3/5 3/5 2/5 2

( ) 2 5 1 2

L

nV nV t p V p V p V V t V t ψ η π

−

∂ ∂ ′ ′ + Γ = ∂ ∂Φ ∂ ∂ ∂ ⎡ ⎤ ′ ′ ′ ′ ⎡ ⎤ + Γ + ∇ − = ⎣ ⎦ ⎢ ⎥ ∂ ∂Φ ∂Φ ⎣ ⎦ ∂Ψ = ∂ q J i ( )

( )

( )

3/5 3/5 2/5 2 * 2 2 2

( ) 2 5 1

R R R R

nV nV t p V p V p V V t t gV R gV R V g t R ψ ψ ψ ψ η ψ ψ η μ η ψ ψ ψ μ

− − −

∂ ∂ ′ ′ ⎡ ⎤ + ∇ = ⎣ ⎦ ∂ ∂ ⎡ ⎤ ∂ ∂ ∂ ′ ′ ′ ′ ⎡ ⎤ + ∇ + ∇ − = ⎢ ⎥ ⎣ ⎦ ∂ ∂ ∂ ⎣ ⎦ ∂Ψ + ∇Ψ = Δ Ψ ∂ ⎡ ⎤ ∂ ∂ ′ ′ ′ + ∇ − ∇ ∇ = ⎢ ⎥ ∂ ∂ ⎣ ⎦ u u q J u u i i i i i i

2 * 2

1

R

g gR R η μ

−

⎡ ⎤ ∇ ∇Ψ ⎢ ⎥ ⎢ ⎥ ⇒ ∇Ψ = Δ Ψ − ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ u i i

( ) ( )

2 2

1 1 2

R L

g gR V R η ψ ψ π μ

−

∇ ∇Ψ Γ ≡ ∇Ψ ≡ u i i

ψ → Φ

28

( )

1 2 d d ψ φ τ π Φ = ∇ =

∫ ∫

B B S i i

φ

Just as the toroidal flux is the integral of the toroidal magnetic field over the area inside a surface, we can define the poloidal magnetic flux in a similar way. Recall:

( )

1 2 1 2 2

PF

d Jd d d ψ θ τ π φ θ ψ θ φ π π Ψ = ∇ = ∇ ×∇Ψ ∇ = Ψ

∫ ∫

Bi i

Similarly:

Rotational Transform:

( )

1 ( )

PF

d d q ι ψ ψ Ψ ≡ = Φ

[ ]

( )

2/3

( ) 3 2

L L

N N t Q K V V t dV t d σ ι

−

∂ ∂ ′ ′ + Γ = ∂ ∂Φ ∂ ∂ ∂ ′ + = ∂ ∂Φ ∂Φ ∂ = ∂ Φ

( )

2 * 2

1 g gR R η μ

−

⎡ ⎤ ∇ ∇Ψ ⎢ ⎥ ⎢ ⎥ ⇒ Γ Φ = Δ Ψ − ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ i

( )

2 2 2 2 2

1 2 (2 ) 5 2

L

g gR V R V K q R Q V q p πη ψ μ π μ

−

∇ ∇Ψ ≡ ∇Φ ′ ≡ ⎡ ⎤ ′ ≡ ∇Φ + Γ ⎢ ⎥ ⎣ ⎦ i i

( )

5/3 PF

N nV pV d d σ ι ψ ′ ′ ≡ ′ ≡ Ψ ≡ Φ

These are the evolution equations for the 1D Adiabatic Variables. Note: time derivatives are zero if there is no dissipation :

p ∇ = × J B

Also, This is solved in a way that the adiabatic variables stay fixed during solution.

Grad-Hogan Method

• H. Grad and J. Hogan, PRL, 24 1337 (1970)

Just needs heat fluxes to close equations!

( )

η = = q

p ∇ = × J B

This must be solved in a way that the adiabatic variables stay fixed during solution.

Grad-Hogan Method

• H. Grad and J. Hogan, PRL, 24 1337 (1970)

( )

*

g g φ φ μ φ φ = ∇ ×∇Ψ + Ψ ∇ = ∇× = Δ Ψ∇ + ∇ ×∇ B J B

* 2

dp dg R g d d Δ Ψ + + = Ψ Ψ

Grad-Shafranov Equation Must express p and g in terms of adiabatic variables:

( ) ( ) ( )

5/3 2

q g p V dV d R σ

− −

Ψ ′ Ψ = Ψ = Ψ

Alternate advancing the adiabatic variables in time, and then solving the equilibrium equation holding the adiabatic variables fixed.

Time sequence of using the TSC code to model the evolution of a highly elongated plasma in the TCV tokamak. At each instant of time, the vacuum vessel is providing stabilization on the fast (ideal MHD) time scale. The external coils are both feedback stabilizing the plasma and providing shaping fields as they slowly elongate it to fill the entire vessel. In this case, there were 4 PID feedback systems corresponding to:

• Vertical position
• Radial position
• Elongation
• Squareness

Marcus, Jardin, and Hofmann, PRL, 55 2289 (1985) 31

32

Codes can also accurately model the current drive action

• f the OH coils.

Simulation of flattop phase of a basic tokamak discharge. (a) At start of flattop, OH coil has current in same direction as plasma current (b) Flux in plasma uniformly increases due to resistive

• dissipation. OH and

Vertical field coils adjust boundary values so flux gradient in plasma remains almost unchanged. (c) At end of flattop, OH coil has current in opposite direction as plasma current. Vertical field coil OH coil

φ = ∇ ×∇Ψ

P

B

33

PF5 OH PF1AL PF1AU PF3L PF3U

Simulation of NSTX discharge evolution

As a validation exercise, we have simulated the evolution of a NSTX discharge using the experimental values of the coil currents as the preprogrammed currents. To control the plasma in the simulation, several feedback systems need to be added to the coil groups. The “goodness”

• f the simulation is measured by how

small the current in these feedback systems is to still match other measured quantities (such as the flux in flux loops). In general, we find that if we can match the plasma density and temperature evolution, then we can predict the plasma current evolution very accurately.

( ) ( ) ( )

i PP FB

I t I t I t = +

IOH vs time IPF3U vs time IPF3L vs time IPF5 vs time IPF1AU vs time IPF1AL vs time Simulation IOH has feedback added to match experimental plasma current IP Simulation IPF3U and IPF3L have vertical stability feedback added Simulation IPF5 have radial feedback system added experiment simulation

PF5 OH PF1AL PF1AU PF3L PF3U Red are simulation flux loop data and green are experimental data. Origin of each curve is approximate position of flux loop around machine. Excellent agreement!

Summary

• Resistive time-scale dynamics arise both from dissipation in plasma

(resistivity, thermal conductivity) and from resistance in nearby conductors

• On this timescale, plasma inertia can be ignored
• Equations can be written in a moving flux-coordinate system and averaged
• ver the flux surfaces to give a set of 1D equations + time
• Dissipation coefficients only give relative motion of particles, energies,

poloidal, and toroidal fluxes, so one can be chosen as reference

• Choosing the amount of toroidal flux within a surface as the reference flux

coordinate is the most natural for tokamaks

• These equations evolve the 1D adiabatic variables in time. We also must

solve the equilibrium equation in a way that the adiabatic variables stay fixed during the solution….. Grad-Hogan method

36

37