Status Of Gyro-Landau-fluid Development In BOUT

Status of gyro-Landau-fluiddevelopment in BOUT C. H. Ma1, 2X. Q. Xu2, P. W. Xi1, 2, T. Y. Xia2,3, T.F. Tang1,4, A. Dimits2, M. V. Umansky2, S. S. Kim5 , I. Joseph2, P Snyder61FusionSimulation Center, School of Physics, Peking University, Beijing, ChinaLivermore National Laboratory, Livermore, CA 94550, USA3Institute of Plasma Physics, Chinese Academy of Sciences, Hefei, China4Dalian University of Technology, Dalian, China5WCI Center for Fusion Theory, NFRI, Korea6General Atomics, San Diego, CA 92186, USA2LawrencePresented at 2015 BOUT mini-workshopDecember 16, 2015 Livermore, CaliforniaLLNL-PRES-6801851

Outline Introduction Code structure Normalization Implementation User interface Applications Benchmarks Physics results Summary2

3 1 GLF simulations in BOUT Turbulent transport constrains thepedestal gradient in the edge; The linear KBM physics in EPED1successfully predicts the H-modepedestal height and width; 3 1 gyro-Landau-fluid (GLF)model is implemented in BOUT to include the KBM turbulenceeffects in nonlinear ELMsimulations.1.P. B. Snyder , et al., NF(2011)3

3 1 Gyro-Landau-fluid model We utilize the gyrofluid model1 developed by P. Snyder and G.Hammett; 3 1 model: (𝑛𝑖 , 𝑣 𝑖 , 𝑝 𝑖 , 𝑝 𝑖 , 𝜛, 𝐴 , 𝑝 𝑒 , 𝑝 𝑒 ) Full FLR effects (Padé approximation): 𝑘 𝜌𝑖 1; Parallel Landau damping: non-local transport; Non-Fourier methods Both in collisionless2 and weakly-collisional3 limits Toroidal resonance; Non-isotropic response (𝑝 𝑝 ) In long-wavelength limit (𝑘 𝜌𝑖 1) and isotropic assumption(𝑝 𝑝 ), the set of equations is reduced to 6-field Landaufluid model2 with gyro-viscosity.1.2.3.4.P. B. Snyder and G. W. Hammett, PoP (2001)A. M. Dimits, et al., PoP (2014)M.V. Umansky, et al., J. Nucl. Mater. (2015)T. Y. Xia, X. Q. Xu and P. W. Xi, Nucl. Fusion (2013)4

Gyrofluid equations are derived by momentshierarchy from gyrokinetic equations𝑓 𝒙, 𝒗Kinetic model: 6D 𝑓𝑒 𝒗𝑓 𝑣 𝑬 𝒗 𝑩 𝑓 0 𝑡𝑚Average over gyro-motion which𝜇 is adiabatically conserved𝑓 𝑅, 𝑣 , 𝜇Gyrokinetic model: 5D 𝑓 𝑓 𝑹 𝑓 𝑣 0 𝑡 𝑣 moments𝑛 𝑅 , 𝑢 𝑅 ,𝑝 𝑅 , 𝑝 𝑅 ,𝑞 𝑅 , 𝑞 𝑅 , Gyrofluid model: 3D𝑛 𝑓𝑑 3 𝑣𝑛𝑢 𝑓𝑣 𝑑 3 𝑣𝑝 𝑚 𝑓 𝑣 𝑢 2 𝑑 3 𝑣𝑝 𝑚 𝑓𝐵𝜇𝑑 3 𝑣𝑞 3𝑚𝑣𝑡2 𝑛0 𝑢 𝑚 𝑓𝑣 3 𝑑 3 𝑣𝑞 𝑚𝑣𝑡2 𝑛0 𝑢 𝑚 𝑓𝐵𝜇𝑣 𝑑 3 𝑣5

Full set of ion equationsin 3 1 GLF modelFLR effectContinuityCompressionLandau dampingToroidal closure6

Vorticity formulation is used with fullelectron response in 3 1 GLF modelFLR effectContinuityCompressionLandau dampingEnergy fluxDrift Alfven WavePoisson equation: Have better numerical property than 𝑛𝑒 equation1. P.W. Xi, X.Q. Xu, et al., Nucl. Fusion (2013)7

Carefully chosen closures areessential to match kinetic effectsGyro averaging and Padé approximation:Ion Landau closures:Electron Landau closures:Toroidal closures:Energy fluxToroidal closure 2 (Imaginary)Toroidal closure 3 (Real)8

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary9

Initialization (physics init)Physics initRead optionsEquilibriumNormalizationCoefficients Equilibrium cases (P0,N0, T0): 1: 𝜂𝑖 scan profiles; 2: cyclone case; 4: tanh functionprofiles; 5: Self-consistentbootstrap current grid; 6: Real geometry withexperimental profiles.Save profilesEnd10

Time evolution (physics run)Physics runCalculate other fieldsGet electric fieldGyro-averaged quantitiesClosuresEvolving equations 𝑛𝑒 , 𝑇 𝑖 , 𝑇 𝑖 , 𝑇 𝑒 , 𝑇 𝑒 𝐽 , 𝑉 𝑒 Real space: 𝑛𝑖 , 𝑢 𝑖 Electrostatic potential: 𝜙 Φ, 𝐴 Parallel Heat flux: 𝑞 𝑖 , 𝑞 𝑖 , 𝑞 𝑒 , 𝑞 𝑒 Toroidal resonance𝑑𝑑𝑑𝑑 Ion: 𝑑𝑡 𝑛𝑖 , 𝑑𝑡 𝑢 𝑖 , 𝑑𝑡 𝑃 𝑖 , 𝑑𝑡 𝑃 𝑖𝑑𝑑𝑑𝑑 Electron: 𝑑𝑡 𝑈, 𝑑𝑡 𝐴 , 𝑑𝑡 𝑃 𝑒 , 𝑑𝑡 𝑃 𝑒End11

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary12

Basic normalized quantities Normalization parameters: 𝐿, 𝑇, 𝑁, 𝐵𝑉 𝐿 𝑇 , 𝑉 2 𝑉𝐴2 𝐵2 𝜇0 𝑚𝑖 𝑁 , Ω 𝑒𝐵 𝑚𝑖 , 𝐶𝑛𝑜𝑟 Ω𝑇 Define Evolving variablesother important variables13

Normalized ion equations14

Normalized electron equations15

Normalized vorticity equationsWhere:16

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary17

Gyro-average operator Gyro-average operator with Pade approximation: Multiply the denominator in the operator:1 𝜌𝑖2 2 Φ 𝜙 This equation is solved using Laplacian inversion in the code:18

Modified Laplacian operators The modified Laplacian operators can be expressed as thesubtract of gyro-average operators:19

Poisson equation The gyro-kinetic Poisson equation is Define Poisson equation becomes With Pade approximation In the code20

Poisson equationwith adiabatic response The adiabatic response is used in the electrostatic simulations; The gyro-kinetic Poisson equation with adiabatic response is1/2𝑛0 𝜙 𝑛0𝑛 Γ01/2 1 Γ0 𝜙 Γ0 𝑛𝑖 𝑏 0 𝑇𝑖 𝑇𝑒0 𝑇𝑖0𝑇𝑖0 𝑏 After normalization and Pade approximation, the equationbecomes𝑇0 𝑛𝑖22 2𝜌𝑖 𝜙 𝜙 1 𝜌𝑖 𝐶𝑛𝑜𝑟 𝑛0 In the code21

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary22

Options in the input fileOptionDescriptionelectrostaticSolve Poisson equation with adiabatic response,instead of solving the electron equations.eHallElectron drift wave termsGyroaverageGyro-average and FLR effect termsFLR effectContinuityCompressible termsCompressionParallel viscosityIsotropic (default: false)Average the parallel and perpendicular pressureLandau damping iLandau damping terms for ionsLandau damping wcoll iCollisional terms in ion Landau dampingLandau damping eLandau damping terms for electronsLandau damping wcoll eCollisional terms in electron Landau dampingEnergy fluxEnergy flux terms in toroidal closuresToroidal closure2The imaginary part of 𝜔𝑑 terms in toroidal closuresToroidal closure3The real part of 𝜔𝑑 terms in toroidal closures23

Options in the input file cont. curv model: Controls the implementation ofcurvature term (𝑖𝜔𝑑 ) 1: 𝑖𝜔𝑑 2: 𝑖𝜔𝑑 3: 𝑖𝜔𝑑 𝑇0𝑏 𝜅 ;𝑒𝐵𝑇0𝑏 𝐵 ;𝑒𝐵2𝑇01𝑏 𝜅 𝑏 𝐵 ;2𝑒𝐵B𝑇0𝑩′𝑏 𝜅′ , where 𝜅 3 𝑒𝐵𝐵 4: 𝑖𝜔𝑑 𝜇0 𝑃0 the magnetic curvature calculated in the code; 5: 𝑖𝜔𝑑 𝑇02𝑒𝐵𝐵22is1B𝑏 𝜅′ 𝑏 𝐵 .24

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary25

Our 3 1 GLF code is benchmarked with other gyrokinetic and gyro-fluid code in ES ITG simulations With adiabatic electron response,the 3 1 GLF results are in goodagreement with other gyrofluidcode (GLF 3 0) and gyorokineticcode (FULL).𝜙1. C.H. Ma, X.Q. Xu, et al., PoP (2015)Adiabatic response:1/2𝑛0 𝜙 𝑛0 𝜙𝑛0 Γ01/2 1 Γ0 𝜙 Γ0 𝑛𝑖 𝑏𝑇𝑇𝑒0𝑇𝑖0𝑇𝑖0 𝑏 𝑖 26

The Landau closure with collisions hasimplemented and tested in BOUTN0 1014cm3N0 1013cm3 The nonlocal heat flux hascloser results for the strongercollision case; The implementation inBOUT is well benchmarkedwith Maxim’s results1.1. M.V. Umansky, et al., J. Nucl. Mater. (2015) 27

The linear growth rates of 3 1 and 6-field model agreewell in lower mode numbers28

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary29

The global beta-scan with a series of selfconsistent equilibrium Jet like equilibrium; The temperature profile is fixed; The density and pressure profilesincrease when 𝛽 increases. Weakly-collisional case, 𝜈𝑒 0.1 𝜂𝑖 0.685, the same for all cases.30

Kinetic physics has stabilizing effects onballooning modes This is kinetic ballooning modes(KBM) because there is no instabilitywithout curvature drive; The real frequency is around thetheoretical prediction; Since 𝜂𝑖 1, threshold of KBM andIBM is about the same.Same ion diamagneticstabilizing effects!1. C.H. Ma, X.Q. Xu, et al., PoP (2015)31

The relative energy loss increases withincreasing beta When beta increases, therelative energy loss of initialcrash from all channelsincrease; Convective energy loss isdominant because of thelarge density height; Energy loss from ion is largerthan the loss from electron; Electron perturbation isdamped by the Landaudamping effect, which islarger than ions by a factor of𝑚𝑖 /𝑚𝑒 .32

KBM is unstable below IBM thresholdwhen temperature gradient is large Concentric circular magneticsurfaces without shift; The unstable threshold arethe same for the constanttemperature case (𝜂𝑖 0).; The KBM is unstable underideal ballooning modethreshold when 𝜂𝑖 2. There is no second stableregion in these casesbecause the Shafranov shifteffects is missed; i 0 i 233

The shift increases with beta in ourequilibrium The shift and maximumelongation increase withbeta in the equilibrium.Radial profile of Elongation34

Second stable region of KBM is foundin the self consistent beta scan Beta scan in a series of selfconsistent grids; The linear growth rate for theballooning mode peaks at𝛽 0.66%; The KBM is unstable when𝛽𝑐 0.4% and 𝛼𝑐 1.8; The second stable region ofKBM is observed when𝛽 0.9% and 𝛼 3.75; The growth rate of idealballooning mode is larger thanKBM in this case.Parameters:𝑛 20, (𝑘𝜃 𝜌𝑖 0.11)𝑞 2.0, 𝑠 2.70𝑅 3𝑚, 𝑎 2𝑚𝜖𝑛 1/18.6, 𝜂𝑖 𝜂𝑒 3.135

Initial GLF simulations in X-point geometryk s0.280.560.841.121.4080100DIII-Dshot 144977 low currentshot 144981 high currentgrowthrate(/ )0.3n 300.20.10.0204060mode number n The GLF simulations show reliable mode structures The DIII-D magnetic and plasma profiles are ideal P-B mode stable The mode is at the outside midplane side, driven by the bad curvature The Landau damping and toroidal closures are turned off in theseinitial simulations36

The Landau damping closures have similarimpact on the ELM size as flux-limited heat fluxLandau damping:𝑞 𝑗 𝑛0𝑖𝑘 𝑘𝐵 𝑇𝑗8𝑣𝑇𝑗, 𝑗 𝑖, 𝑒0.5𝜋𝑘 𝜆𝑗Flux limiting:𝜅𝑆𝐻 𝜅𝐹𝑆𝑞 𝑗 0 𝑇𝑗𝜅𝑆𝐻 𝜅𝐹𝑆 Nonlinear simulation shows that the energy loss of am ELMare similar with Landau damping closure or flux-limited heatflux in 6-field Landau-fluid simulations.1. C.H. Ma, X.Q. Xu, et al., PoP (2015)37

Outline Introduction Code structure Normalization Implementation Inputs Applications Benchmarks Physics results Summary38

Summary The 3 1 GLF model is implemented in the BOUT framework for pedestal turbulence and transport; The 3 1 GLF model is well benchmarked with othergyrokinetic, gyrofluid and two-fluid codes in bothelectromagnetic and electrostatic regimes; The energy loss of an ELM increases with beta nearthe first stable region; The second stable region of ballooning mode isfound in our self-consistent beta scan with the globalequilibrium.39

Install the 3 1 module Get the bout glf git repository:git clone rs/bout glf/www/git/bout glf.git Switch to the modomegad branch:git checkout bout modomegad Compile the code (on cori):./configure --with-netcdf /global/u2/c/chma/cori/local -with-fftw /global/u2/c/chma/cori/localmakecd examples/glfkbm3-1make40

Running example Launch the linear run job:qsub bout cori debug.sh Get the linear growth rate:python growthrate.py -v -f P data Example output:Growth rate from linear fit from -20 to -1 is: 0.149107741

Dec 16, 2015 · The Landau damping closures have similar impact on the ELM size as flux-limited heat flux Nonlinear simulation shows that the energy loss of am ELM are similar with Landau damping closure or flux-limited heat flux in 6-field Landau-fluid simulations. 1. C.H. Ma,