SCIDAC Center for Extended MHD Modeling

EXTENDED MHD MODELING:
BACKGROUND, STATUS AND VISION
Dalton D. Schnack
Center for Energy and Space Science
Science Applications International Corp.
San Diego, CAOVERVIEW• The Extended MHD model• The computational challenges:– Extreme separation of time scales– Extreme separation of spatial scales– Extreme anisotropy– Importance of geometry, boundary conditions– Causality: can’t parallelize over time!
At least as challenging as hydrodynamic turbulence!
• Present computational approaches:– Implicit time differencing– Specialized spatial grids• Status of present models• Vision for integrated modelingDEFINITIONS• Hydrodynamics - A mathematical model that describes themotion of a continuous, isotropic fluid• Magnetohydrodynamics (MHD) - A mathematical model that
describes the motion of a continuous, electricallyconducting fluid in a magnetic field– Hydrodynamics and Maxwell equations coupled throughLorentz body force and Ohm’s law• Ideal MHD - the fluid has infinite electrical conductivity(zero resistivity)• Resistive MHD - The fluid has finite conductivity andresistivity• Extended MHD - additional effects of electron dynamicsand/or non-Maxwellian species

MODERN TOKAMAKS ARE RICH IN MHD ACTIVITY
Example: DIII-D shot 86144

3/2 NTM
2/1 wall locking

1/1 sawteeth

Disruption•Sawtoothing discharge•3/2 NTM triggered at 2250 msec
•2/1 locks to the wall

MODELING REQUIREMENTS• Slow evolution
Nonlinear fluid model required• Plasma shaping
Realistic geometry required• High temperature
Large “Reynolds’ numbers”• Low collisionality
Extensions to resistive MHD required• Strong magnetic field
Highly anisotropic transport required• Resistive wall
Non-ideal boundary conditions required

APPROACHES• Quasi-equilibrium:Ñp = J´ B + constraints- “Magnetohydrostatics” (MHS?)- Eliminates all waves- Basis for 1-1/2 dimensional transport models- Extension to 3-D?
• Time dependent- Solve 2-fluid equations- Retain all normal modes- Focus of present SciDAC efforts

FLUID MODELS••

Kinetic models of plasmas based on distribution function foreach charge species
Satisfies kinetic equationdfa= åC [fa ,f b ]dtbfa (x,v,t ) - six dimensions plus time

- computationally impractical for time scales of interest••••

Fluid models derived by taking successive velocity moments ofkinetic equation– Reduce dimensionality by 3
Hierarchy of equations for n, v, p, P, q, …….
Equations truncated by closure relations– Express high order moments in terms of low order moments– Capture kinetic effects in these moments
Result is Extended MHD

2-FLUID MODEL• Maxwell (no displacement current):¶B= -Ñ ´ E ,¶tÑ ´ B = m 0J ,• Momentum, energy, and continuity for each species (a = e, i):æ ¶va
öma na ç+ va × Ñva ÷ = -Ñ ×Pa + qa na (E + va ´ B) + å Rab + Smaè ¶tøb3¶pa+ va × Ñpa = - pa Ñ × va - Pa : Ñv a - Ñ × qa + Qa2¶t¶na= -Ñ × (na v a ) + San¶t• Current and quasi-neutrality:
Ja = na qa va ,n = ne = Zni

SINGLE FLUID FORM
• Add electron and ion momentum equations:æ ¶v
örç + v ×Ñv÷ = -Ñ × P¢ + J ´ Bè ¶tø• Subtract electron and ion momentum equations (Ohm’s law):
E=- v´B123 +Ideal MHDhJ{+Resistive MHD

1 1- nJ´B
ne1+n 441442
3
Hall Effecté¶Jù11Ñ × (Pe¢ - nPi¢ ) ++ Ñ × (vJ+ Jv)úê2ë
ûne(1+) 24444
0 w pe (1+ n ) ¶t
144n4
43 e1444442444443
Diamagnetic Effects
and Closures

Electron Inertia

All effects beyond resistivity constitute Extended MHD

COMPUTATIONAL CHALLENGES•Extreme separation of time scales– Realistic “Reynolds’ numbers”– Implicit methods•

Extreme separation of spatial scales– Important physics occurs in internal boundary layers– Small dissipation cannot be ignored– Requires grid packing or adaption•Extreme anisotropy– Special direction determined by magnetic field
Accurate treatment of B × Ñ operator is important- Requires specialized gridding

SEPARATION OF TIME SCALES<t{At{S<<t{evolSound transit time

Affvén transit time<<MHD evolution time

=Lundquist number: StRtA~ 10 8tA
N<<<<t{RResistive diffusion time>> 1Explicit time step impractical:Dt<DxtAL»

Require implicit methodst evol

IMPLICIT METHODS• Partially implicit methods– Treat fastest time scales implicitly– Time step still limited by waves• Semi-implicit methods– Treat linearized ideal MHD operator implicitly– Time step limited by advection– Many iterations• Fully implicit methods– Newton-Krylov treatment of full nonlinear equations– Arbitrary time step
– Still a research projectLINEAR SOLVER REQUIREMENTS•

10Extremely large condition number : > 10 !!

– Specialized pre-conditioners– Anisotropy•

Ideal MHD is self-adjoint– Symmetric matrices– CG•

Advection and some 2-fluid effects (whistler waves) are not selfadjoint– Need for efficient non-symmetric solvers•

Everything must be efficient and scalable in parallel•

Should interface easily with F90

SEPARATION OF SPATIAL SCALES• Important dynamics occurs in internal boundary layers– Structure is determined by plasma resistivity or otherdissipation– Small dissipation cannot be ignored• Long wavelength along magnetic field• Extremely localized across magnetic field:d /L~S-a<< 1 for S >> 1

• It is these long, thin structures that evolve nonlinearly on
the slow evolutionary time scale

EXTREME ANISOTROPY• Magnetic field locally defines special direction in space• Important dynamics are extended along field direction, very
narrow across it• Propagation of normal modes (waves) depends strongly on
local field direction• Transport (heat and momentum flux) is also highly
anisotropic==> Requires accurate treatment of operator B × Ñ

Inaccuracies lead to “spectral pollution” andanomalous perpendicular transport

GRIDDING AND SPATIAL REPRESENTATION• Spatial stiffness and anisotropy require special gridding– Toroidal and poloidal dimensions treated differently• Toroidal (f, primarily along field)– Long wavelengths, periodicity => FFTs (finitedifferences also used)• Poloidal plane (R,Z)– Fine structure across field direction– Grids aligned with flux surfaces (~ field lines)– Unstructured triangular grids– Extreme packing near internal boundary layers• Finite elements– High order elements essential for resolving anisotropies• Dynamic mesh adaption in research phase

POLOIDAL GRIDS

DIII-D poloidal cross-section withflux aligned grid (NIMROD)

Circular poloidal cross-section withtriangles and grid packing (M3D)

Poloidal grids from SciDAC development projects

BEYOND RESISTIVITY - EXTENDED MHD• 2-fluid effects– Whistler waves (Hall term) require implicit advance with
non-symmetric solver– Electron inertia treated implicitly– Diamagnetic rotation may cause accuracy, stabilityproblems• Kinetic effects - influence of non-Maxwellian populations– Analytic closures• Seek local expressions for P, q, etc.
– Particle closures• Subcycle gyrokinetic df calculation• Minority ion species - beam or a-particles

STATUS••

2 major SciDAC development projects for time-dependent models– M3D - multi-level, 3-D, parallel plasma simulation code
• Partially implicit• Toroidal geometry - suitable for stellarators• 2-fluid model• Neo-classical and particle closures– NIMROD - 3-D nonlinear extended MHD• Semi-implicit• Slab, cylindrical, or axisymmetric toroidal geometry• 2-fluid model (evolving computationally)• Neo-classical closures• Particle closures being de-bugged
Both codes have exhibited good parallel performance scalingOther algorithms are being developed in the fusion programSTATUS - RESISTIVE MHD

Sawtooth in NSTX computed by M3D

Stellarator ballooning mode computed by M3D

STATUS - RESISTIVE MHD

Secondary magnetic islands generated during
sawtooth crash in DIII-D shot 86144 by NIMRODSTATUS - EXTENDED MHD•

Effect of energetic
particle population onMHD mode•

Subcycling ofenergetic particle
module within MHD
codes

•

M3D agrees well withNOVA2 in the linear
regime•Energetic particles arebeing incorporated
into NIMRODSTATUS - EXTENDED MHD
Neo-classical tearing modes with NIMROD
using analytic closure

Computed NTM stability
diagram

Sawtooth and secondaryislands

Secondary island growth
Varying SNEXT STEP - INTEGATED MODELING••••

Non-local kinetic physics, MHD, and profile evolution are all interrelated– Kinetic physics determines transport coefficients– Transport coefficients affect profile evolution– Profile evolution can destabilize of MHD modes– Kinetic physics can affect nonlinear MHD evolution (NTMs,TAEs)– MHD relaxation affects profile evolution– Profiles affect kinetic physics
Effects of kinetic (sub grid scale) physics must be synthesized
into MHD models– Extensions to Ohm’s law (2-fluid models)– Subcycling/code coupling– Theoretical models (closures), possibly heuristic
Effects of MHD must be synthesized into transport models
Predictions must be validated with experimental data

VISION: VDE EVOLUTION
EQUILIBRIUM
RECONSTRUCTION

DESIGN

RF CODE:
SOURCESLINEARSTABILITY
CODE

TRANSPORTCODE
GK CODE:
TRANSPORT
COEFFICIENTSCOMPARISONSTABLE?
EDGE
CODE

STRESSANALYSIS
CODE

EXPERIMENTAL
DATA

CURRENTS IN SHELLUNSTABLE?

NONLINEAR
EXTENDED MHD
CODESIMULATION DATAVISION: SAWTOOTH CYCLEEQUILIBRIUM
RECONSTRUCTION

GK CODE:TRANSPORT
COEFFICIENTSEXPERIMENTAL
DATATRANSPORTMODULERELAXEDPROFILESSTABLE?
COMPARISONLINEARSTABILITY
CODE WITH
ENERGETIC PARTICLESUNSTABLE?

ENERGETICPARTICLE
MODULENONLINEAR
EXTENDED MHD
CODESIMULATION DATAENABLING COMPUTER SCIENCE TECHNOLGIES•••••

Largest, fastest computers!
– But intermediate computational resources often neglected,and…
– The computers will never be large or fast enough!
Algorithms– Parallel linear algebra– Gridding, adaptive and otherwiseData structure and storage– Adequate storage devices– Common treatment of experimental and simulation data
– Common tools for data analysis
Communication and networking– Fast data transfer between simulation site and storage site
– Efficient worldwide access to data
– Collaborative tools– Dealing with firewalls
Advanced graphics and animation

SUMMARY•

Predictive simulation capability has 3 components– Code and algorithm development– Tightly coupled theoretical effort– Validation of models by comparison with experiment• Integration required for:– Coupling algorithms for disparate physical problems– Theoretical synthesis of results from different models– Efficient communication and data manipulation• Progress is being made with Extended MHD– Integration of energetic ion modules into 3-D MHD– Computationally tractable closures– Resistive wall modules
Need to bring a broader range of algorithms and codes to bear foroverall fusion problem