TN-13-2_NumericalReducedModelCoupling2D2VDriftKineticIons2D3VKineticNeutralsHelicalMagne
========================================================================================

.. meta::
   :description: technical note

   :keywords: Report,2047357-TN-13,Numerical,study,of,a,reduced,model,coupling,2D+2V,drift-kinetic,ions,and,2D+3V,kinetic,neutrals,in,a,helical,magnetic,field,with,wall,boundaries,M.,R.,Hardman1,,J.,Omotani3,,M.,Barnes1,,and,F.,I.,Parra3,1,Rudolf,Peierls,Centre,for,Theoretical,Physics,,University,of,Oxford,,Clarendon,Laboratory,,Parks,Road,,Oxford,OX1,3PU,,United,Kingdom,2,Princeton,Plasma,Physics,Laboratory,,P.O.,Box,451,,Princeton,,New,Jersey,08540,,United,States,3,Culham,Centre,for,Fusion,Energy,,Culham,Science,Centre,,Abingdon,,Oxon,,OX14,3DB,,United,Kingdom,E-mail:,michael.hardman@physics.ox.ac.uk,1.,Introduction,The,aim,of,this,report,is,to,describe,an,initial,investigation,into,the,challenges,associated,with,developing,2D,edge,plasma,models,with,kinetic,ion,and,neutral,species.,We,focus,on,two,fundamental,aspects,of,the,problem:,we,consider,the,implementation,of,model,ion-neutral,collision,operators,,and,the,implementation,of,simple,wall,boundary,conditions.,There,are,other,challenges,to,developing,numerical,simulations,in,the,plasma,edge,,for,example,,a,comprehensive,model,must,also,deal,with,the,magnetic,flux,coordinate,singularity,at,the,separatrix,between,the,core,plasma,and,the,scrape,off,layer,,as,well,as,the,inherently,global,nature,of,the,problem,in,both,the,radial,and,poloidal,directions.,For,a,full,review,of,the,physics,challenges,,please,see,[1].,Previous,reports,in,this,series,have,examined,the,1D,problem,of,neutrals,and,ions,interacting,on,a,single,magnetic,field,line,,neglecting,cross-field,physics,[2,,3,,4,,5].,Both,the,‘standard’,Eulerian,drift-kinetic,approach,and,the,novel,‘moment,kinetics’,approach,proposed,in,[6,,7,,8,,9],were,implemented,and,tested,[3,,4,,5].,It,was,demonstrated,that,the,‘moment,kinetics’,approach,is,a,viable,scheme,for,solving,the,system,numerically,,although,some,problems,were,encountered,with,wall,boundary,conditions,[5].,Analytical,tests,of,the,implemented,code,were,developed,and,automated,,and,an,initial,physics,study,was,carried,out,looking,at,the,difference,between,the,ion,and,neutral,distribution,functions,in,steady,state,in,the,case,with,wall,boundary,conditions.,The,major,differences,between,the,2D,and,1D,models,are,the,following.,Cross-,field,drifts,now,enter,into,the,ion,drift,kinetic,equation,,and,the,difference,between,charged,particle,and,neutral,particle,trajectories,manifest,themselves,in,qualitatively,different,velocity,space,structure,for,the,two,different,types,of,particle,species.,Whilst,1,2,DRIFT-KINETIC,IONS,COUPLING,TO,KINETIC,NEUTRALS,2,neutral,particles,are,only,accelerated,by,interparticle,collisions,,charged,particles,are,also,constrained,to,perform,helical,cyclotron,orbits,around,magnetic,field,lines,by,the,Lorentz,force.,This,means,that,whilst,neutral,particles,must,be,described,with,three,velocity,coordinates,,we,can,describe,the,motion,of,charged,particles,with,just,two,velocity,coordinates,,by,averaging,over,the,rapid,cyclotron,motion.,This,fact,has,consequences,for,the,implementation,of,any,kinetic,model.,For,simplicity,in,the,initial,development,of,the,implementation,of,the,2D,model,,in,this,report,we,restrict,our,attention,to,the,‘standard’,drift-kinetic,approach.,We,take,as,our,starting,point,the,2D,model,with,wall,boundary,conditions,proposed,in,[7].,We,present,the,2D,model,in,the,next,section.,In,section,3,,we,discuss,the,differing,velocity-space,grids,of,charged,and,neutral,particles,,and,the,consequent,challenges,for,the,simulation.,In,sections,4,and,5,we,consider,the,wall,boundary,conditions,for,ions,and,neutrals,,and,electrons,,respectively.,In,section,6,,we,write,down,the,normalised,system,of,equations,that,are,implemented,in,the,2D,code.,In,section,7,,we,briefly,describe,the,numerical,algorithms,underpinning,the,implementation,of,the,model.,In,section,8,,we,describe,tests,that,verify,the,implementation,of,the,2D,model,using,the,method,of,manufactured,solutions.,Finally,,in,section,9,,we,provide,suggestions,for,future,development.,2.,Drift-kinetic,ions,coupling,to,kinetic,neutrals,Here,we,provide,a,brief,overview,of,the,2D,kinetic,model,that,was,derived,in,detail,in,[7].,The,system,that,we,consider,consists,of,a,single,ion,species,of,charge,e,and,mass,mi,,a,single,neutral,species,of,mass,mn,=,mi,,and,an,electron,species,that,is,modelled,as,having,a,Boltzmann,response.,The,magnetic,field,is,taken,to,be,helical,,with,the,form,B,=,Bzˆz,+,Bζ,ˆζ,,(1),where,(r,,z,,ζ),are,cylindrical,coordinates,,and,Bz,and,Bζ,are,constant,in,z,,r,,and,ζ.,See,figure,1.,Helical,field,lines,approximate,the,geometry,of,the,open-field-line,scrape-,off,layer,region,of,plasma,that,is,confined,in,an,axisymmetric,toroidal,magnetic,field.,With,this,interpretation,we,can,identify,z,as,the,vertical,coordinate,along,the,axis,of,symmetry,,r,as,the,major,radial,coordinate,and,ζ,as,the,toroidal,angle.,We,include,a,simple,model,charge,exchange,collision,operator,acting,between,ions,and,neutrals,,and,a,simple,model,ionisation,collision,operator,involving,ions,,electrons,and,neutrals.,We,do,not,account,for,charged-species,collisions,that,couple,ions,to,ions,,or,electrons,to,ions.,To,be,consistent,with,the,lack,of,radial,variation,in,the,magnetic,field,,we,impose,periodic,boundary,conditions,in,r,,but,we,accommodate,either,wall,boundary,conditions,or,periodic,boundary,conditions,in,z.,We,assume,that,the,plasma,equilibrium,is,independent,of,ζ.,The,model,kinetic,equations,are,as,follows.,We,have,a,drift-kinetic,equation,for,2,DRIFT-KINETIC,IONS,COUPLING,TO,KINETIC,NEUTRALS,3,Figure,1:,The,form,of,the,magnetic,field,and,the,left-handed,coordinates,(r,,z,,ζ),in,terms,of,the,right-handed,cartesian,(x,,y,,z).,Taken,from,[7].,To,make,a,connection,to,the,geometry,of,the,scrape-off,layer,–,the,open,field,lines,surrounding,the,toroidal,flux,surfaces,of,the,core,plasma,–,we,should,think,of,r,as,the,major,radial,coordinate,,z,as,the,vertical,axial,coordinate,,and,ζ,as,the,toroidal,angle.,ions,∂Fi,∂t,(cid:18),+,bzv∥,−,Er,B,(cid:19),∂Fi,∂z,+,+,Ez,B,∂Fi,∂r,∂Fi,∂v∥,−,Rin,(nnFi,−,ni,⟨Fn⟩),+,Rionne,⟨Fn⟩,+,Si,,ebzEz,mi,=,(2),(cid:113),where,Fs,is,the,particle,distribution,function,for,species,s,,ns,is,the,particle,number,density,for,species,s,,v∥,=,v,·,b,and,v⊥,=,|v,−,v∥b|,are,the,components,of,the,particle,velocity,v,that,are,parallel,and,perpendicular,to,the,magnetic,field,direction,b,=,B/B,,ζ,.,The,operator,⟨·⟩,=,(cid:82),π,respectively,,with,B,=,−π,dϑ/2π,is,the,gyroaverage,,ϑ,is,the,gyro-angle,,t,is,the,time,,bz,=,Bz/B,,ϕ,is,the,electrostatic,potential,,Rin,and,Rion,are,the,charge,exchange,and,ionization,collision,frequency,factors,,respectively,,and,Si,is,a,source,function,that,may,be,used,to,inject,particles,,momentum,,and,heat,,or,to,facilitate,a,test,via,the,method,of,manufactured,solutions.,In,addition,to,equation,(2),we,must,solve,the,2D-3V,kinetic,equation,for,neutrals:,z,+,B2,B2,∂Fn,∂t,+,vz,∂Fn,∂z,+,vr,∂Fn,∂r,=,−Rin,(niFn,−,nnFi),−,RionneFn,+,Sn,,(3),where,Sn,is,the,neutral,source,function.,We,use,velocity,coordinates,aligned,with,ˆz,,ˆr,,and,ˆζ,to,describe,the,advection,of,the,neutral,species;,here,,we,use,the,coordinates,vz,=,v,·,ˆz,,vr,=,v,·,ˆr,,and,vζ,=,v,·,ˆζ.,3,EVALUATING,THE,COLLISION,OPERATOR,TERMS,4,The,ion,density,ni,is,computed,by,integrating,over,the,ion,distribution,function:,ni(z,,r,,t),=,2π,(cid:90),∞,(cid:90),∞,dv∥,−∞,0,dv⊥v⊥Fi(z,,r,,v⊥,,v∥,,t).,The,same,is,true,for,the,neutral,density,nn:,nn(z,,r,,t),=,(cid:90),∞,(cid:90),∞,(cid:90),∞,dvr,dvζ,−∞,−∞,−∞,dvzFn(z,,r,,vζ,,vr,,vz,,t).,(4),(5),The,electrostatic,potential,is,computed,by,enforcing,quasineutrality,ni,=,ne,and,using,a,Boltzmann,response,for,electrons,,i.e.,,ni,=,ne,=,Ne,exp,(cid:19),,,(cid:18),eϕ,Te,(6),where,the,constant,Ne,is,either,taken,to,be,a,reference,density,,or,calculated,through,a,simple,electron,sheath,model.,We,obtain,the,electric,fields,by,differentiation,of,ϕ:,Ez,=,−,∂ϕ,∂z,,,and,Er,=,−,∂ϕ,∂r,.,(7),3.,Evaluating,the,collision,operator,terms,There,are,two,major,differences,between,the,2D,and,1D,models.,The,first,is,that,the,E,×,B,drift,now,appears,,causing,spatial,advection,across,field,lines,radially.,The,second,difference,is,that,the,charged,particles,and,neutral,particles,have,very,different,velocity,space,dynamics.,A,large,component,of,this,project,is,the,investigation,of,how,to,efficiently,couple,neutrals,to,charged,particles,in,a,kinetic,code.,Charged,particles,and,neutral,particles,couple,in,two,distinct,ways:,through,the,wall,boundary,condition,,and,through,inter-particle,collisions.,In,this,section,,we,discuss,a,problem,related,to,the,fast,implementation,of,the,collision,operators.,Charged,particles,in,strong,magnetic,fields,undergo,a,helical,motion:,particles,are,free,to,travel,along,the,magnetic,field,line;,but,the,v,×B,component,of,the,Lorentz,force,causes,particles,to,execute,a,rapid,cyclotron,(gyro-),motion,in,the,plane,perpendicular,to,B,,with,a,gyration,radius,ρi,much,smaller,than,either,the,radial,or,axial,length,scales,Lr,and,Lz,,respectively.,Therefore,,it,is,natural,to,use,the,field-aligned,velocity,coordinates,(v∥,,v⊥,,ϑ),,with,which,we,can,write,v,=,v∥b,+,v⊥(cos,ϑ,ˆr,−,sin,ϑ,b,×,ˆr).,(8),This,yields,a,great,simplification,in,the,description,of,the,charged,particle,species,because,in,strong,magnetic,fields,dv∥,dt,cs,Lz,where,the,sound,speed,is,cs,=,(cid:112)2Te/mi,and,the,cyclotron,frequency,is,Ωi,=,eB/mi.,Because,we,are,not,typically,interested,in,motions,on,frequencies,comparable,to,Ωi,,we,can,conveniently,average,over,the,ϑ,motion,,resulting,in,the,drift-kinetic,equation,(2).,dv⊥,dt,dϑ,dt,∼,Ωi,,(9),≪,∼,∼,3,EVALUATING,THE,COLLISION,OPERATOR,TERMS,5,None,of,these,simplifications,apply,to,the,neutral,species,,which,experience,no,forces,due,to,the,electromagnetic,fields.,In,the,absence,of,collisions,,neutral,particles,execute,rectilinear,motion,,freely,crossing,magnetic,field,lines.,It,is,therefore,natural,to,consider,a,neutral,velocity,space,grid,that,is,aligned,with,the,spatial,coordinates,,i.e.,,v,=,vzˆz,+,vrˆr,+,vζ,ˆζ,(10),Such,a,grid,may,be,useful,if,the,kinetic,neutral,code,must,interface,with,legacy,software,that,uses,a,particular,spatially,aligned,coordinate,system.,There,are,,however,,major,disadvantages,to,using,different,velocity,space,coordinates,for,the,charged,particle,and,neutral,particle,species.,The,simple,reason,for,this,is,that,there,is,a,nontrivial,transformation,between,(vz,,vr,,vζ),and,(v∥,,v⊥,,ϑ).,Combined,with,the,collision,operator,terms,in,equations,(2),and,(3),,this,nontrivial,transformation,requires,a,3D,interpolation,of,distribution,function,data,at,every,time,integration,step.,Interpolation,between,3D,grids,is,necessarily,slow,and,likely,to,be,a,source,of,error,unless,an,advanced,interpolation,scheme,is,used,(such,as,a,spectral,interpolation,or,cubic,spline).,We,now,write,the,velocity,coordinate,transformation,for,later,use.,Using,equations,(8),and,(10),,one,can,show,that,and,vz(v∥,,v⊥,,ϑ),=,v∥bz,−,v⊥,sin,ϑbζ,vr(v∥,,v⊥,,ϑ),=,v⊥,cos,ϑ,vζ(v∥,,v⊥,,ϑ),=,v∥bζ,+,v⊥,sin,ϑbz,,,v∥(vz,,vr,,vζ),=,vzbz,+,vζbζ,v⊥(vz,,vr,,vζ),=,(cid:113),r,+,(vζbz,−,vzbζ)2,v2,.,(11),(12),To,evaluate,the,cross-species,collision,operator,terms,in,equation,(2),,we,require,the,neutral,distribution,function,Fn(vz,,vr,,vζ),at,every,(z,,r),to,be,gyroaveraged,and,given,on,the,(v∥,,v⊥),grid.,This,is,achieved,by,writing,⟨Fn⟩,(v∥,,v⊥),=,(cid:90),π,−π,Fn(vz(v∥,,v⊥,,ϑ),,vr(v∥,,v⊥,,ϑ),,vζ(v∥,,v⊥,,ϑ)),dϑ,2π,.,(13),Similarly,,to,evaluate,the,charge,exchange,collision,operator,terms,in,equation,(3),we,require,the,charged,particle,distribution,function,Fi(v∥,,v⊥),at,every,(z,,r),to,be,interpolated,,i.e.,,Fi(vz,,vr,,vζ),=,Fi(v∥(vz,,vr,,vζ),,v⊥(vz,,vr,,vζ)).,(14),We,carry,out,the,necessary,interpolation,using,a,readily,available,Julia,package:,Interpolations.jl,[10].,For,simplicity,and,speed,of,evaluation,we,use,a,linear,interpolation,suitable,for,non-uniformly,spaced,data.,4,BOUNDARY,CONDITIONS,ON,IONS,AND,NEUTRAL,PARTICLES,6,A,simple,strategy,to,avoid,the,cost,of,interpolation,would,be,to,choose,to,use,the,same,velocity,space,grid,for,both,charged,and,neutral,particles.,In,the,standard,drift-,kinetic,approach,where,the,moments,are,not,evolved,separately,,this,strategy,would,only,introduce,the,need,to,evaluate,vz,,vr,,and,vζ,via,the,transform,(11).,We,have,not,attempted,to,implement,this,approach,numerically,due,to,time,constraints,,although,we,highlight,this,approach,as,a,possible,future,optimization,of,the,model.,4.,Boundary,conditions,on,ions,and,neutral,particles,We,must,specify,boundary,conditions,on,the,evolved,particle,distribution,function,in,the,coordinates,where,differential,operators,appear:,z,,r,and,v∥.,No,boundary,conditions,are,required,on,v⊥,,vz,,vr,,or,vζ,,as,no,advection,occurs,in,these,coordinates,,and,no,Fokker-Planck,collision,operator,currently,appears,in,the,model.,The,physical,boundary,condition,in,v∥,is,that,Fi(z,,r,,v∥,→,±∞,,v⊥,,t),→,0.,(15),This,condition,may,easily,be,imposed,by,either,forcing,zero,incoming,particles,in,the,v∥,advection,,or,by,enforcing,that,Fi,is,periodic,in,v∥,at,the,ends,of,the,v∥,domain.,This,is,equivalent,to,the,zero,boundary,condition,if,the,v∥,domain,is,sufficiently,large.,For,the,boundaries,in,r,,one,can,imagine,several,physically,interesting,scenarios:,a,periodic,boundary,condition,,in,which,case,the,model,represents,something,like,a,local,‘flux,tube’;,a,Dirichlet,boundary,condition,as,in,a,‘global’,model;,and,a,limited,configuration,with,conditions,that,model,direct,contact,to,the,vacuum,vessel.,In,this,report,we,only,consider,the,first,possibility.,Finally,,we,consider,two,types,of,boundary,condition,in,z;,periodic,boundary,conditions,,and,‘wall’,boundary,conditions,that,model,the,presence,of,the,‘sheaths’,at,the,ends,of,open,field,lines,on,the,vacuum,vessel,or,divertor,plates.,The,periodic,boundary,condition,is,implemented,just,as,a,testing,exercise:,we,do,not,implement,the,vorticity,equation,that,is,necessary,to,calculate,the,potential,in,closed-field-line,geometry,[8].,Only,the,wall,boundary,condition,option,is,physically,consistent.,To,impose,the,wall,boundary,condition,,we,assume,that,the,entrance,to,the,sheath,occurs,at,the,ends,of,the,domain,in,z.,We,also,assume,that,ions,that,exit,the,simulation,domain,continue,on,to,the,wall,,where,they,recombine.,No,ions,enter,the,domain,from,the,walls,,giving,a,zero,incoming,boundary,condition,for,the,ions:,Fi(z,=,−Lz/2,,v∥,>,Er/Bz,,v⊥,,t),=,0,Fi(z,=,Lz/2,,v∥,<,Er/Bz,,v⊥,,t),=,0.,(16),Note,that,in,the,case,that,Er,=,0,then,the,condition,(16),reduces,to,the,equivalent,1D,condition.,The,boundaries,of,the,domain,(corresponding,to,the,sheath,entrances),are,taken,to,be,at,z,=,−Lz/2,and,z,=,Lz/2.,Neutrals,that,leave,the,domain,are,assumed,to,hit,the,wall,and,thermalise,at,the,temperature,of,the,wall,,Tw.,Ions,that,recombine,at,the,wall,also,re-enter,as,neutrals.,The,resulting,boundary,condition,on,the,neutrals,5,A,SIMPLE,ELECTRON,SHEATH,MODEL,is,Fn(z,=,−Lz/2,,vz,>,0,,vr,,vζ,,t),=,(cid:0)Γi,−Lz/2,+,Γn,−Lz/2,Fn(z,=,Lz/2,,vz,<,0,,vr,,vζ,,t),=,(cid:0)Γi,Lz/2,+,Γn,Lz/2,(cid:1),FKw,(cid:1),FKw,(cid:16),vz,,(cid:16),(cid:113),vz,,(cid:113),r,+,v2,v2,ζ,(cid:17),r,+,v2,v2,ζ,,,where,FKw(vz,,vt),(cid:112)v2,is,the,Knudsen,cosine,distribution,[11],,and,(cid:19)2,.,=,3,π,(cid:18),mi,2Tw,|vz|,z,+,v2,t,(cid:18),exp,−,(cid:19),mi,(v2,z,+,v2,t,),2Tw,(cid:17),,,7,(17),(18),Γi,−Lz/2,.,=,2π,(cid:90),Er/B,(cid:90),∞,dv∥,dv⊥v⊥,(cid:12),(cid:12),(cid:12),(cid:12),bzv∥,−,Er,B,(cid:12),(cid:12),(cid:12),(cid:12),Fi(z,=,−Lz/2,,r,,v∥,,v⊥,,t),(19),−∞,(cid:90),0,0,(cid:90),∞,dvz,Γn,−Lz/2,.,=,(cid:90),∞,dvr,dvζ,|vz|,Fn(z,=,−Lz/2,,r,,vz,,vr,,vζ,,t),(20),and,−∞,−∞,−∞,Γi,Lz/2,.,=,2π,(cid:90),∞,(cid:90),∞,dv∥,dv⊥v⊥,Γn,Lz/2,.,=,Er/B,(cid:90),∞,dvz,0,(cid:90),∞,dvr,(cid:90),∞,0,−∞,−∞,(cid:12),(cid:12),(cid:12),(cid:12),bzv∥,−,Er,B,(cid:12),(cid:12),(cid:12),(cid:12),Fi(z,=,Lz/2,,v∥,,v⊥,,t),(21),dvζ,|vz|,Fn(z,=,Lz/2,,vz,,vr,,vζ,,t),(22),are,the,combined,fluxes,of,neutrals,and,ions,towards,the,walls,at,z,=,−Lz/2,and,z,=,Lz/2,,respectively.,In,the,1D,model,implementation,[5],,where,Er,=,0,,we,note,that,the,wall,boundary,condition,is,imposed,on,the,marginalised,distributions,functions,fi,=,2π,(cid:82),dv⊥v⊥Fi,and,fn,=,(cid:82),dvζdvrFn.,Whilst,the,1D,and,2D,implementations,are,able,to,share,almost,all,subroutines,related,to,the,imposition,of,the,wall,boundary,condition,,we,note,that,this,slight,difference,between,the,1D,and,2D,implementations,means,that,the,marginalised,distribution,of,neutral,particles,fn,must,be,proportional,to,the,marginalised,Knudsen,distribution,fKw,=,2π,(cid:82),dv⊥v⊥FKw,at,the,wall.,Please,see,sections,2,and,2.1,of,[5],for,further,details.,5.,A,simple,electron,sheath,model,In,addition,to,the,boundary,conditions,for,the,kinetic,species,described,in,the,last,section,,we,have,implemented,a,simple,electron,sheath,model,that,is,compatible,with,a,Boltzmann,response,in,the,bulk,plasma,[12].,This,sheath,model,allows,us,to,impose,a,wall,potential,ϕW,,,in,place,of,a,fixed-in-time,constant,density,Ne,in,equation,(6).,The,model,is,obtained,as,follows.,First,,we,demand,that,there,is,zero,net,current,to,the,wall,plate,at,the,boundaries,in,z,,i.e.,,the,current,parallel,to,the,magnetic,field,line,vanishes:,J∥,=,J∥,i,+,J∥,e,=,0,at,z,=,±Lz/2,(23),5,A,SIMPLE,ELECTRON,SHEATH,MODEL,8,where,J∥,e,and,J∥,i,are,the,electron,and,ion,contributions,to,the,current,,respectively.,Second,,we,assume,that,the,electron,distribution,function,is,a,Maxwellian,close,to,the,entrance,of,the,sheath,,i.e.,,Fe,≈,ne,π3/2,(cid:18),me,2Te,(cid:19)3/2,(cid:32),exp,−,me(v2,∥,+,v2,⊥),2Te,(cid:33),.,(24),We,assume,that,the,deviation,in,Fe,from,the,Maxwellian,is,due,to,the,high,velocity,tail,of,electrons,that,are,absorbed,by,the,wall.,The,electron,distribution,function,changes,with,z,in,the,sheath,according,to,the,illustration,in,figure,2.,In,particular,,we,note,that,Figure,2:,An,illustration,of,the,assumed,form,of,the,electron,distribution,function,Fe,as,a,function,of,v∥,near,the,wall,at,z,=,−Lz/2.,There,is,a,high,velocity,tail,of,electrons,that,reach,the,wall,despite,the,repelling,wall,potential,,as,as,a,result,,are,lost,to,the,distribution.,The,cut,off,velocity,vcut,drops,to,zero,at,the,wall,itself,,yielding,the,half-,Maxwellian,that,is,used,to,compute,J∥,e,for,(25).,This,argument,is,developed,in,full,in,[12].,any,electrons,that,have,reached,the,wall,are,unable,to,return.,Third,,we,assume,that,J∥,e,,Te,and,Ne,are,constant,across,the,sheath,,and,that,the,Boltzmann,response,(6),continues,to,hold.,Hence,,by,integrating,the,distribution,function,of,electrons,(24),at,the,wall,at,z,=,−Lz/2,,we,obtain,a,relationship,between,J∥,e,and,Ne,,Te,and,ϕW,:,J∥,e,=,ene,(cid:114),Te,2πme,=,eNe,exp,(cid:18),eϕW,Te,(cid:19),(cid:114),Te,2πme,Finally,,we,rearrange,equation,(25),for,Ne,,with,the,result,(cid:19),J∥,i,ecs,(cid:114)4πme,mi,eϕW,Te,Ne,=,−,exp,−,(cid:18),.,Equation,(26),determines,the,constant,Ne,to,be,used,in,equation,(6).,(25),(26),In,the,Bulk,PlasmaIn,the,sheathAt,the,wall,6,NORMALISATION,OF,THE,SYSTEM,OF,MODEL,EQUATIONS,9,norm.,variable,˜t,definition,t(cref/Lref),˜z,˜r,˜v∥,˜v⊥,˜vz,˜vr,˜vζ,(cid:101)Ne,˜ns,(cid:101)ϕ,(cid:101)Er,(cid:101)Ez,˜Rin,˜Rion,(cid:101)Fs,(cid:101)Ss,ρ∗,bz,z/Lref,r/Lref,v∥/cref,v⊥/cref,vz/cref,vr/cref,vζ/cref,Ne/nref,ns/nref,eϕ/Tref,eLrefEr/Tref,eLrefEz/Tref,Rin(nrefLref/cref),Rion(nrefLref/cref),refπ3/2/nref),Fs(c3,refπ3/2Lref/nrefcref),cref/LrefΩref,Ss(c3,Bz/B,ref.,quantity,definition,Lref,Tref,nref,cref,mi,Bref,Ωref,ref.,length,(m),ref.,temperature,(KeV),ref.,density,(m−3),(cid:112)2Tref/mi,(ms−1),ion,mass,(kg),ref.,B,(T),eBref/mi,(s−1),Table,1:,Definitions,for,normalised,and,reference,quantities,used,in,the,report.,6.,Normalisation,of,the,system,of,model,equations,We,now,define,the,normalisations,used,when,solving,numerically,the,model,equations,In,terms,of,the,(2)-(7).,We,define,normalised,and,reference,quantities,in,Table,1.,normalised,variables,,the,model,system,equations,takes,the,following,form.,Equation,(2),takes,the,form,∂,(cid:101)Fi,∂˜t,(cid:16),+,bz,˜v∥,−,ρ∗,2,(cid:101)Er,(cid:17),∂,(cid:101)Fi,∂,˜z,+,ρ∗,2,−,˜Rin,(cid:101)Ez,(cid:16),∂,(cid:101)Fi,∂˜r,+,∂,(cid:101)Fi,∂˜v∥,=,bz,(cid:101)Ez,2,(cid:68),(cid:101)nn,(cid:101)Fi,−,(cid:101)ni,(cid:69)(cid:17),+,˜Rion(cid:101)ne,(cid:69),(cid:68),(cid:101)Fn,+,(cid:101)Si,,(cid:101)Fn,(27),where,ρ∗,=,cref/LrefΩref,is,a,small,parameter,that,measures,the,size,of,the,ion,gyro-,orbits,compared,to,the,macroscopic,size,of,the,system.,Note,that,bz,is,also,a,small,6,NORMALISATION,OF,THE,SYSTEM,OF,MODEL,EQUATIONS,10,parameter,in,the,formal,ordering,of,the,2D,model,[7,,8].,Equation,(3),takes,the,form,(28),(29),(30),(31),(32),∂,(cid:101)Fn,∂˜t,+,˜vz,∂,(cid:101)Fn,∂,˜z,+,˜vr,∂,(cid:101)Fn,∂˜r,(cid:16),=,−,˜Rin,(cid:17),(cid:101)ni,(cid:101)Fn,−,(cid:101)nn,(cid:101)Fi,−,˜Rion(cid:101)ne,(cid:101)Fn,+,(cid:101)Sn,,and,quasineutrality,,equation,(6),,takes,the,form,(cid:101)ni,=,(cid:101)ne,=,(cid:101)Ne,exp,(cid:32),(cid:33),.,(cid:101)ϕ,(cid:101)Te,The,normalised,electric,fields,are,defined,by,(cid:101)Ez,=,−,∂,(cid:101)ϕ,∂,˜z,,,and,(cid:101)Er,=,−,∂,(cid:101)ϕ,∂˜r,.,and,the,densities,have,the,following,definitions,(cid:101)ni,=,1,√,π,(cid:90),∞,(cid:90),∞,d˜v∥,−∞,0,2d˜v⊥˜v⊥,(cid:101)Fi.,and,nn,=,(cid:90),∞,(cid:90),∞,(cid:90),∞,d˜vr,d˜vζ,d˜vz,(cid:101)Fn.,1,π3/2,−∞,The,boundary,conditions,(15),and,(16),carry,over,trivially.,The,wall,boundary,−∞,−∞,condition,on,the,neutrals,takes,the,form,(cid:16),(cid:101)Fn(˜z,=,−Lz/2Lref,,˜vz,>,0,,˜vr,,˜vζ,,˜t),=,(cid:101)Γi,−Lz/2,+,(cid:101)Γn,−Lz/2,(cid:101)Fn(˜z,=,Lz/2Lref,,˜vz,<,0,,˜vr,,˜vζ,,˜t),=,(cid:16),(cid:101)Γi,−Lz/2,+,(cid:101)Γn,−Lz/2,(cid:17),(cid:17),(cid:16),(cid:16),˜vz,,˜vz,,(cid:101)FKw,(cid:101)FKw,(cid:113),r,+,˜v2,˜v2,ζ,(cid:113),r,+,˜v2,˜v2,ζ,(cid:17),(cid:17),,,,,(33),where,|˜vz|,z,+,˜v2,t,is,the,normalised,Knudsen,cosine,distribution,[11],,and,the,normalised,particle,fluxes,are,z,+,˜v2,t,),(cid:101)Tw,(cid:101)FKw,=,(cid:112)˜v2,(34),exp,(˜v2,−,(cid:18),(cid:19),√,3,π,(cid:101)T,2,w,(cid:101)Γi,−Lz/2,.,=,1,√,π,(cid:101)Γn,−Lz/2,.,=,1,π3/2,(cid:90),ρ∗,(cid:101)Er/2,(cid:90),∞,d˜v∥,2d˜v⊥˜v⊥,−∞,(cid:90),0,0,(cid:90),∞,d˜vz,(cid:90),∞,d˜vr,−∞,−∞,−∞,(cid:12),(cid:12),(cid:12)bz,˜v∥,−,ρ∗,2,(cid:12),(cid:12),(cid:101)Fi(˜z,=,−Lz/2Lref,,˜r,,˜v∥,,˜v⊥,,˜t),(cid:12),(cid:101)Er,d˜vζ,|˜vz|,(cid:101)Fn(˜z,=,−Lz/2Lref,,˜r,,˜vz,,˜vr,,˜vζ,,˜t),(35),(36),and,(cid:101)Γi,Lz/2,.,=,1,√,π,(cid:90),∞,(cid:90),∞,d˜v∥,2d˜v⊥˜v⊥,(cid:12),(cid:12),(cid:12)bz,˜v∥,−,ρ∗,2,(cid:12),(cid:12),(cid:101)Fi(˜z,=,Lz/2Lref,,˜r,,˜v∥,,˜v⊥,,˜t),(cid:12),(cid:101)Er,(37),ρ∗,(cid:101)Er/2,(cid:90),∞,0,(cid:90),∞,d˜vz,(cid:101)Γn,Lz/2,.,=,1,π3/2,0,−∞,−∞,(cid:90),∞,d˜vr,d˜vζ,|˜vz|,(cid:101)Fn(˜z,=,Lz/2Lref,,˜r,,˜vz,,˜vr,,˜vζ,,˜t),(38),Finally,,the,electron,sheath,model,for,Ne,,equation,(26),,takes,the,normalised,form,(cid:101)Ne,=,−,(cid:114),4πme,mi,(cid:32),exp,−,(cid:33),(cid:101)ϕw,(cid:101)Te,J∥,i,enrefcs,.,(39),7,NUMERICAL,IMPLEMENTATION,11,7.,Numerical,implementation,The,algorithms,used,to,implement,the,2D,model,are,identical,to,those,used,to,implement,the,1D,model,described,in,previous,reports.,To,avoid,repetition,,we,simply,list,these,algorithms,here,,with,further,details,left,to,previous,reports,[2,,3,,4,,5].,We,evolve,the,normalised,systems,of,equations,with,a,time-advance,scheme.,To,maximise,the,efficiency,of,the,algorithm,with,regards,to,numerical,stability,,we,utilise,a,Strong,Stability,Preserving,(SSP),Runge-Kutta,(RK),scheme,[13,,14,,15].,By,default,,we,use,the,four-stage,SSPRK3,method,,although,two,and,three,stage,SSPRK,schemes,are,options,in,the,code.,The,implementation,of,the,time-stepping,algorithm,was,tested,in,the,1D-1V,code,(figures,1,and,2,of,[3]).,For,the,spatial,and,velocity,discretisation,both,finite,difference,and,Chebyshev,spectral,methods,are,implemented.,For,the,finite,difference,option,,a,uniformly,spaced,grid,is,used.,Derivatives,are,evaluated,with,the,third,order,upwind,scheme,[16],and,integration,is,carried,out,using,the,composite,Simpson’s,rule,[3].,For,the,Chebyshev,spectral,elements,option,,the,coordinate,grid,is,the,Gauss-Chebyshev-Lobatto,grid,on,each,element,[17].,Inclusion,of,the,endpoints,within,each,element,facilitates,enforcement,of,continuity,at,element,boundaries,,and,the,use,of,Chebyshev,polynomials,as,a,basis,In,our,code,,these,transforms,are,done,enables,the,use,of,Fast,Fourier,Transforms.,using,the,widely-used,FFTW,library,[18].,The,associated,integration,weights,used,for,integration,are,obtained,using,Clenshaw-Curtis,quadrature,rules,[19].,Clenshaw-Curtis,quadrature,allows,for,the,use,of,endpoints,in,the,integration,domain,while,still,exactly,integrating,polynomials,up,to,degree,Ngrid,−,1,,with,Ngrid,the,number,of,points,within,the,element.,The,implementation,of,the,spatial,discretisation,was,tested,in,the,1D-1V,code,(figures,3,,4,and,5,of,[3]).,The,scale,of,the,2D-3V,model,,having,5,dimensions,for,the,neutral,species,,and,4,for,the,ion,species,,requires,a,larger,memory,requirement,than,can,be,supported,by,a,single,core,of,a,typical,computer.,To,reach,large,enough,resolutions,to,achieve,the,results,presented,in,this,report,,it,was,necessary,to,implement,shared,memory,OpenMPI,support,[20].,The,previous,implementation,of,the,1D-1V,code,[2,,3,,4,,5],allowed,for,the,development,of,a,suite,of,automatic,functional,tests.,The,tests,come,in,two,classes.,In,the,first,class,,the,automated,1D-1V,tests,check,the,numerical,results,for,sound,wave,damping,of,small,amplitude,fluctuations,against,analytically,obtained,results,[21,,4].,In,the,second,class,of,test,,results,from,converged,‘nonlinear’,simulations,with,finite,fluctuation,amplitudes,are,stored,and,tested,against,the,output,of,the,code,after,a,new,revision.,The,current,2D-3V,version,of,the,code,passes,all,previous,1D-1V,automated,tests,,suggesting,that,the,time-stepping,and,spatial,discretisation,algorithms,have,been,preserved,in,the,generalisation,of,the,code,to,the,2D-3V,case.,All,revisions,of,the,‘moment,kinetics’,code,are,written,in,the,Julia,programming,language,,and,are,currently,available,on,GitHub,at,https://github.com/mabarnes/,moment_kinetics.,The,latest,2D-3V,code,is,held,in,the,branch,https://github.com/,8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS12,mabarnes/moment_kinetics/tree/radial-vperp-standard-DKE-with-neutrals.,8.,Testing,the,numerical,implementation,with,the,method,of,manufactured,solutions,The,normalised,system,of,equations,(27)-(32),is,difficult,to,solve,analytically,,especially,with,the,wall,boundary,conditions,described,in,sections,4-5.,This,presents,a,challenge,for,testing,the,numerical,implementation,of,the,2D,model.,To,overcome,this,challenge,the,method,of,we,use,a,standard,technique,from,the,fluid,mechanics,community:,manufactured,solutions.,The,great,power,of,the,method,of,manufactured,solutions,is,that,tests,derived,from,this,method,can,be,tuned,to,focus,on,individual,operators,,or,the,whole,code,,as,desired.,In,the,following,,we,present,tests,which,can,be,used,separately,to,cover,the,time,advance,algorithm,,the,spatial,advection,,the,velocity,advection,,the,collisions,,periodic,boundary,conditions,,and,the,wall,boundary,conditions.,New,tests,can,be,developed,as,features,are,added,to,the,model.,i,n,and,(cid:101)F,M,S,We,use,the,following,method,to,test,the,2D,model,equations,with,manufactured,solutions,[22].,First,,we,propose,manufactured,,target,,solutions,(cid:101)F,M,S,that,we,can,write,down,in,closed,form,with,known,functions.,The,manufactured,solutions,are,only,restricted,to,obey,the,boundary,conditions,that,are,imposed,on,the,system,of,equations.,Second,,we,compute,manufactured,sources,(cid:101)Si,and,(cid:101)Sn,by,inserting,the,manufactured,solutions,(cid:101)F,M,S,into,equations,(27),and,(28),and,evaluating,the,necessary,derivatives,and,velocity,integrals,analytically.,Third,,we,perform,a,numerical,simulation,using,the,manufactured,sources.,We,run,the,simulation,for,a,fixed,time,˜t,=,O(1),(to,allow,numerical,errors,to,accrue),for,varying,resolutions.,We,check,that,the,numerical,solution,is,consistent,with,the,target,manufactured,solution,by,confirming,that,the,numerical,error,reduces,towards,zero,as,the,resolution,increases.,To,measure,the,error,between,the,numerical,solution,and,the,target,manufactured,solution,,we,introduce,the,following,measures.,We,define,an,error,on,the,ion,and,neutral,species,density,and,(cid:101)F,M,S,n,i,(cid:115)(cid:80),ϵ((cid:101)ns),=,i,j,|(cid:101)ns(zi,,rj),−,(cid:101)nM,S,s,(zi,,rj)|2,,,(40),NrNz,s,where,(cid:101)nM,S,is,the,target,manufactured,density,of,species,s,,and,Nr,and,Nz,are,the,total,number,of,points,in,the,r,and,z,grids,,respectively.,We,define,an,error,on,the,ion,distribution,function,(cid:118),(cid:117),(cid:117),(cid:116),ϵ(,(cid:101)Fi),=,(cid:80),i,j,k,l,|,(cid:101)Fi(v∥i,,v⊥j,,zk,,rl),−,(cid:101)F,M,S,i,Nv∥Nv⊥NrNz,(v∥i,,v⊥j,,zk,,rl)|2,,,(41),where,Nv∥,and,Nv⊥,are,the,total,number,of,points,in,the,v∥,and,v⊥,grids,,respectively.,Finally,,we,define,an,error,on,the,neutral,distribution,function,(cid:118),(cid:117),(cid:117),(cid:116),ϵ(,(cid:101)Fn),=,(cid:80),i,j,k,l,m,|,(cid:101)Fn(vzi,,vrj,,vζk,,zl,,rm),−,(cid:101)F,M,S,n,Nvz,NvrNvζ,NrNz,(vzi,,vrj,,vζk,,zl,,rm)|2,,,(42),8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS13,where,Nvz,,,Nvr,,and,Nvζ,are,the,total,number,of,points,in,the,vz,,vr,and,vζ,grids,,respectively.,The,errors,ϵ((cid:101)ns),,ϵ(,(cid:101)Fi),,and,ϵ(,(cid:101)Fn),measure,the,average,error,per,point,in,the,array,of,the,numerical,solution.,As,noted,previously,[22],,Julia,has,support,for,elements,of,symbolic,algebra,via,the,Symbolics.jl,package,[23,,24].,This,allows,us,to,partially,automate,the,process,of,calculating,the,manufactured,sources.,Currently,,the,Symbolics.jl,package,supports,symbolic,differentiation,,allowing,for,a,great,reduction,in,the,effort,of,developing,manufactured,solutions,tests.,Unfortunately,,symbolic,integration,is,not,yet,supported,by,the,Symbolics.jl,package,,so,the,target,manufactured,solutions,must,be,chosen,to,be,sufficiently,simple,to,integrate,by,hand.,8.1.,Tests,with,periodic,boundary,conditions,in,r,and,z,We,first,describe,a,manufactured,solutions,test,that,is,compatible,with,periodic,boundary,conditions,in,r,and,z.,We,choose,target,manufactured,solutions,that,depend,on,time,and,space,in,a,nontrivial,manner,,allowing,for,tests,of,the,timestepping,algorithm,and,the,advection,in,r,,z,,and,v∥.,We,choose,the,ion,density,to,be,of,the,form,(cid:101)ni,=,3,2,+,sin(2π˜t),10,(cid:18),sin,(cid:19),(cid:18),2πr,Lr,+,sin,(cid:18),2πz,Lz,(cid:19)(cid:19),.,We,choose,the,neutral,density,to,take,the,form,(cid:101)nn,=,3,2,+,sin(2π˜t),10,(cid:18),cos,(cid:19),(cid:18)2πr,Lr,+,cos,(cid:18)2πz,Lz,(cid:19)(cid:19),.,(43),(44),We,specify,the,ion,and,neutral,distributions,as,Maxwellians,with,constant,temperatures,and,no,flow,,i.e.,,and,(cid:101)Fn,=,(cid:101)nn,exp,(cid:0)−˜v2,z,−,˜v2,r,−,˜v2,ζ,(cid:101)Fi,=,(cid:101)ni,exp,(cid:0)−˜v2,∥,−,˜v2,⊥,(cid:1),,,(cid:1),.,(45),(46),The,gyroaverage,of,Fn,in,equation,(27),is,computed,analytically,by,noting,that,˜v2,∥,+,˜v2,ζ,.,This,test,could,easily,be,extended,to,include,a,specified,temperature,or,mean,velocity,moment,in,addition,to,the,density,moment.,⊥,=,˜v2,r,+,˜v2,z,+,˜v2,Having,specified,target,solutions,,we,must,calculate,(cid:101)Si,and,(cid:101)Sn,for,specific,input,physics,parameters,–,the,target,solutions,can,be,used,to,test,cases,with,or,without,collisions,,and,with,any,specific,choice,of,the,helical,magnetic,field.,We,first,focus,on,a,case,with,bz,=,0.5,,ρ∗,=,1.0,,(cid:101)Te,=,(cid:101)Ne,=,1.0,and,˜Rin,=,˜Rion,=,0.,We,must,also,describe,the,numerical,input,parameters.,We,use,spectral-element,coordinate,grids,for,all,dimensions,[3].,These,grids,are,characterised,by,Ngrid,,the,number,of,points,in,each,spectral,element,,and,Nelement,,the,number,of,elements.,The,total,number,of,grid,points,is,N,total,grid,=,(Ngrid−1)Nelement+1.,We,identified,that,the,biggest,source,of,numerical,error,comes,from,the,velocity,grid,discretization.,In,the,following,,we,present,tests,where,we,8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS14,Figure,3:,Plot,of,the,measures,of,the,error,on,the,ion,and,neutral,density,,and,the,ion,and,neutral,distribution,function,,defined,in,equations,,(40),,(41),,and,(42),,respectively.,The,simulations,use,the,manufactured,solutions,(43)-(46),,and,retain,all,but,the,collision,terms,in,equations,(27)-(32).,The,scan,is,in,the,number,of,elements,Nelement,in,the,spectral-element,scheme,in,the,velocity,space,coordinates,v∥,,v⊥,,vz,,vr,,and,vζ.,The,relatively,small,errors,observed,for,Nelement,=,16,suggests,that,the,implementation,is,performing,well.,have,increased,the,number,of,elements,Nelement,in,the,velocity,space,grids,v∥,,v⊥,,vz,,vr,,and,vζ,,holding,fixed,Ngrid,=,5,for,all,coordinates,,and,holding,fixed,Nelement,=,2,in,the,coordinates,r,and,z.,We,take,the,maximum,velocity,vmax/cref,=,6.0,on,all,grids.,We,simulate,for,160,time,steps,with,∆˜t,=,0.002.,The,results,of,these,simulations,are,shown,in,figure,3.,The,numerical,error,in,all,fields,is,either,small,,or,reduces,for,increasing,Nelement.,At,the,largest,Nelement,=,16,,the,error,in,the,ion,and,electron,densities,appears,to,be,converged.,The,residual,error,is,likely,due,to,the,finite,size,of,∆˜t,,see,fig,2,in,[3].,The,input,files,for,these,simulations,are,provided,in,Appendix,A.1.,In,a,second,test,,we,now,introduce,the,charge-exchange,and,ionisation,collision,operators.,We,use,the,input,parameters,described,above,,with,˜Rin,=,˜Rion,=,1.0,,and,we,again,scan,in,the,velocity,space,resolution.,The,results,of,the,simulations,appear,in,figure,4.,We,again,see,convergence,with,increasing,Nelement,,but,at,a,slower,rate,than,in,the,case,without,collisions.,We,can,attribute,the,slower,convergence,to,the,presence,of,the,linear,interpolation,in,the,collision,operators,,as,discussed,in,section,3.,The,linear,interpolation,is,likely,to,be,more,inaccurate,for,larger,values,of,Ngrid,,where,the,spacing,of,the,grid,points,in,each,element,become,more,nonuniform.,The,input,files,for,these,simulations,are,provided,in,Appendix,A.2.,8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS15,Figure,4:,Plot,of,the,measures,of,the,error,on,the,ion,and,neutral,density,,and,the,ion,and,neutral,distribution,function,,defined,in,equations,,(40),,(41),,and,(42),,respectively.,The,simulations,use,the,manufactured,solutions,(43)-(46),,and,retain,all,terms,in,equations,(27)-(32).,The,parameter,Nelement,is,in,the,number,of,spectral,elements,in,the,velocity,space,coordinates,v∥,,v⊥,,vz,,vr,,and,vζ.,The,size,of,the,errors,observed,for,Nelement,=,16,suggests,that,the,implementation,is,performing,well.,8.2.,Tests,with,wall,and,sheath,boundary,conditions,in,z,We,now,describe,tests,for,the,implementation,of,the,wall,and,sheath,boundary,conditions,described,in,sections,4-5.,For,simplicity,,we,consider,periodic,boundary,conditions,in,r.,We,must,construct,a,target,distribution,function,that,satisfies,the,necessary,boundary,conditions.,The,ion,distribution,function,(cid:101)Fi,must,satisfy,the,condition,(16).,We,can,achieve,this,by,constructing,a,distribution,function,from,the,velocity,coordinate,v∥,=,˜v∥,−,ρ∗,(cid:101)Er/2bz.,We,choose,(cid:34),(cid:101)Fi,=,H,(cid:0)v∥,(cid:1),v2,∥,(cid:19),(cid:18),1,2,+,z,Lz,n+(r),+,H,(cid:0)−v∥,(cid:1),v2,∥,(cid:18),1,2,(cid:19),(cid:18)1,2,−,+,(cid:19),(cid:19),z,Lz,z,Lz,n−(r),(cid:35),n0(r),exp,(cid:0)−v2,∥,−,˜v2,⊥,(cid:1),,,(47),+,(cid:18),1,2,−,z,Lz,with,H,(x),the,Heaviside,function,,taking,the,values,1,for,x,>,0,and,0,for,x,<,0.,In,the,Julia,implementation,H,(0),=,1/2.,Note,that,the,forward,going,part,of,the,distribution,vanishes,at,z,=,−Lz/2,and,the,backward,going,part,vanishes,at,z,=,Lz/2,–,ensuring,that,condition,(16),holds.,In,principle,,we,could,specify,any,functions,for,n+(r),,n−(r),8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS16,Figure,5:,Plot,showing,the,error,measures,on,the,ion,and,neutral,densities,,and,the,ion,and,neutral,distribution,functions,for,a,test,case,with,wall,boundary,conditions,in,z,,no,radial,coordinate,(meaning,(cid:101)Er,=,vr,=,0),,and,the,full,ion,and,neutral,velocity,spaces.,The,target,manufactured,solutions,(47)-(53),are,used.,Electrons,are,taken,to,be,Boltzmann,with,(cid:101)Ne,=,1.0,,and,no,collisions,are,included,in,the,simulation.,Convergence,of,the,errors,is,observed,for,increasing,number,Nelement,of,spectral,elements,in,the,velocity,space,coordinates.,and,n0(r),that,satisfy,the,radial,boundary,conditions.,We,take,n+,=,n−,=,n0,,with,n0,=,1,+,1,20,sin,(cid:18)2πr,Lr,(cid:19),.,(48),We,choose,the,radially,varying,component,of,n0,to,be,relatively,small,to,ensure,that,any,E,×,B,drifts,are,easy,to,resolve,with,vmax/cref,=,6.0,–,note,that,radial,variation,in,n0,induces,radial,variation,in,ni,,and,hence,,in,ϕ,,thus,generating,an,Er,̸=,0.,This,choice,is,necessary,,because,the,manufactured,solution,is,proportional,to,a,shifted,Gaussian,in,v∥,,with,the,shift,proportional,to,Er/B.,If,Er/B,is,set,too,large,,then,the,bulk,of,the,Gaussian,is,off,grid.,The,Julia,Symbolics.jl,package,does,not,currently,support,symbolic,or,numerical,integration.,To,specify,the,ion,density,,needed,to,compute,the,potential,via,quasineutrality,,we,compute,the,velocity,integral,of,equation,(47),analytically,,with,the,result,(cid:101)ni,=,n+(r),4,(cid:18),1,2,+,z,Lz,(cid:19),+,n−(r),4,(cid:18),1,2,−,z,Lz,(cid:19),+,n0(r),(cid:18),1,2,+,(cid:19),(cid:18)1,2,z,Lz,−,(cid:19),.,z,Lz,(49),To,write,down,the,neutral,target,distribution,function,for,the,MMS,test,,we,also,require,8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS17,Figure,6:,Plot,showing,the,error,measures,on,the,ion,and,neutral,densities,,and,the,ion,and,neutral,distribution,functions,for,a,test,case,with,wall,boundary,conditions,in,z,,no,radial,coordinate,(meaning,(cid:101)Er,=,vr,=,0),,and,the,full,ion,and,neutral,velocity,spaces.,The,target,manufactured,solutions,(47)-(53),are,used.,Electrons,are,taken,to,be,Boltzmann,with,(cid:101)Ne,determined,by,the,sheath,through,equation,(39),,and,no,collisions,are,included,in,the,simulation.,Convergence,of,the,errors,is,observed,for,increasing,number,Nelement,of,spectral,elements,in,the,velocity,space,coordinates.,the,particle,fluxes,(cid:101)Γi,−Lz/2,and,(cid:101)Γi,Lz/2.,We,find,that,(cid:101)Γ−Lz/2,=,bzn−(r),√,2,π,,,and,(cid:101)ΓLz/2,=,bzn+(r),√,2,π,.,(50),For,the,neutral,particles,we,can,specify,any,target,distribution,function,that,satisfies,the,condition,(33).,We,choose,(cid:34),(cid:101)Fn,=H,(˜vz),(cid:101)Γi,−Lz/2,(cid:19)2,(cid:18)1,2,−,z,Lz,(cid:35),+,(cid:101)Γn,(cid:101)FKw,(cid:16),˜vz,,(cid:113),r,+,˜v2,˜v2,ζ,(cid:17),(cid:34),+,H,(−˜vz),(cid:101)Γi,Lz/2,(cid:19)2,(cid:18)1,2,+,z,Lz,(cid:35),+,(cid:101)Γn,(cid:101)FKw,(cid:16),˜vz,,(cid:113),r,+,˜v2,˜v2,ζ,.,(cid:17),(51),The,neutral,particle,flux,(cid:101)Γn,is,a,constant.,That,we,can,identify,this,constant,with,the,particle,flux,is,demonstrated,by,noting,that,1,π3/2,(cid:90),0,(cid:90),∞,(cid:90),∞,d˜vr,d˜vz,−∞,−∞,−∞,d˜vζ,|˜vz|,(cid:101)FKw,(cid:16),˜vz,,(cid:113),r,+,˜v2,˜v2,ζ,(cid:17),=,1.,(52),8,TESTING,THE,NUMERICAL,IMPLEMENTATION,WITH,THE,METHOD,OF,MANUFACTURED,SOLUTIONS18,We,must,also,evaluate,the,neutral,density.,We,find,that,(cid:101)nn,=,(cid:32),√,π,3,4,(cid:101)T,1/2,w,(cid:101)Γi,−Lz/2,(cid:19)2,(cid:18),1,2,−,z,Lz,+,(cid:101)Γi,Lz/2,(cid:19)2,(cid:18),1,2,+,z,Lz,(cid:33),+,2(cid:101)Γn,.,(53),To,perform,these,integrals,,we,used,the,identities,√,(cid:90),∞,0,(cid:90),∞,0,√,u,x2,+,u2,exp,−,du,=,(cid:18),(cid:19),u2,v2,π,2,v,erfc,(cid:17),(cid:16),x,v,exp,(cid:19),,,(cid:18)x2,v2,x,erfc(x),dx,=,1,4,,,and,(cid:90),∞,0,x2,erfc(x),dx,=,1,√,3,π,.,(54),The,above,target,manufactured,solutions,can,be,used,to,test,the,implementation,of,the,wall,and,sheath,boundary,conditions.,We,first,describe,a,test,for,the,wall,boundary,condition,in,z,in,isolation,from,the,radial,coordinate,and,sheath,boundary,condition.,For,our,physics,parameters,we,take,(cid:101)Ne,=,1.0,(obviating,the,need,for,an,electron,sheath),,bz,=,1.0,,(cid:101)Te,=,(cid:101)Tw,=,1.0,and,˜Rin,=,˜Rion,=,0.,For,the,numerical,input,we,take,∆˜t,=,0.001,,200,time,steps,,Ngrid,=,1,and,Nelement,=,1,for,r,(so,that,no,radial,advection,is,carried,out),,and,Ngrid,=,6,for,all,other,dimensions.,We,take,Nelement,=,2,for,z,,and,we,vary,Nelement,for,the,velocity,space,coordinates.,The,resulting,error,measures,from,these,simulations,are,provided,in,figure,5.,We,note,the,rapid,convergence,of,the,four,error,measures,as,Nelement,increases.,We,also,note,that,it,is,important,for,the,integration,scheme,that,an,even,number,of,Nelement,is,used.,This,is,to,make,sure,that,the,discontinuity,in,the,neutral,distribution,(cid:101)FKw,appears,at,an,element,boundary.,If,an,odd,number,of,elements,were,utilised,,no,convergence,would,be,observed.,We,can,use,an,identical,test,setup,to,check,implementation,of,the,electron,sheath,boundary,condition.,To,evaluate,the,ion,current,at,the,wall,we,use,the,fluxes,evaluated,in,equation,(50).,Instead,of,specifying,(cid:101)Ne,=,1.0,as,a,constant,,we,calculate,(cid:101)Ne,via,equation,(39).,The,resulting,numerical,errors,from,this,series,of,simulations,is,shown,in,figure,6.,We,see,almost,(but,not,quite),identical,results,as,in,figure,5.,As,expected,from,the,reports,on,the,1D,wall,boundary,condition,in,the,1V,case,[5],,the,implementation,handles,well,the,wall,boundary,condition,in,the,3V,case,,providing,that,the,velocity,space,resolution,is,high,enough.,The,input,files,for,the,simulations,in,the,1D-3V,wall,tests,are,provided,in,Appendix,A.3,and,Appendix,A.4.,Finally,,we,describe,the,results,of,testing,the,wall,boundary,condition,in,a,case,with,both,z,and,r,spatial,dimensions.,To,avoid,complication,in,evaluating,the,electrostatic,potential,,in,this,test,we,revert,to,a,simple,Boltzmann,response,for,electrons,with,(cid:101)Ne,=,1.0.,We,take,Ngrid,=,6,and,Nelement,=,2,for,the,r,dimension,,and,we,take,ρ∗,=,1.0.,Otherwise,,we,use,the,physics,and,numerical,input,parameters,described,above.,We,again,vary,the,number,of,velocity,space,spectral,elements.,We,plot,the,resulting,error,measures,from,the,simulations,in,figure,7.,In,contrast,to,the,other,tests,presented,in,this,report,,figure,7,does,not,show,good,convergence,as,Nelement,increases.,In,fact,,more,detailed,investigations,reveal,that,the,error,results,from,deviations,between,the,numerical,and,symbolic,solutions,at,the,wall,boundary,,see,figures,8,and,9.,This,9,CONCLUSIONS,AND,FUTURE,PLANS,19,suggests,that,the,numerical,integration,scheme,is,not,handling,well,the,discontinuities,in,the,ion,distribution,function,that,now,depend,on,(cid:101)Er.,It,is,worth,noting,that,the,electric,fields,show,larger,errors,than,the,ion,density,,indicating,that,problems,with,the,spatial,differentiation,may,also,contribute,to,the,issue,highlighted,here.,Further,investigations,are,required,to,determine,an,adequate,solution,to,this,problem.,The,input,files,for,the,simulations,in,the,2D-3V,wall,tests,are,provided,in,Appendix,A.5.,Figure,7:,Plot,showing,the,error,measures,on,the,ion,and,neutral,densities,,and,the,ion,and,neutral,distribution,functions,for,a,test,case,with,wall,boundary,conditions,in,z,,periodic,boundary,conditions,in,r,,and,the,full,ion,and,neutral,velocity,spaces.,Equations,(27)-(32),are,solved,with,a,Boltzmann,electron,response,and,(cid:101)Ne,=,1.0.,Collisions,are,removed,by,setting,˜Rin,=,˜Rion,=,0.0.,The,target,manufactured,solutions,(47)-(53),are,used.,Convergence,of,the,errors,is,not,observed,for,increasing,number,Nelement,of,spectral,elements,in,the,velocity,space,coordinates,,suggesting,an,underlying,problem,with,the,implementation,of,the,numerical,scheme.,9.,Conclusions,and,future,plans,The,work,presented,in,this,report,consists,of,two,parts.,First,,the,implementation,of,,and,description,of,,a,2D,(drift-),kinetic,model,that,includes,both,ion,and,neutral,species,and,wall,boundary,conditions.,Second,,the,development,of,manufactured,solutions,tests,that,can,be,used,to,check,the,accuracy,of,the,implementation.,To,our,knowledge,,the,model,presented,here,is,unique,in,the,plasma,physics,community,in,supporting,drift-kinetic,ions,,fully,kinetic,neutrals,,and,wall,boundary,conditions,for,both,species.,The,code,has,the,potential,to,be,developed,into,a,leading,9,CONCLUSIONS,AND,FUTURE,PLANS,20,(a),(c),(b),(d),Figure,8:,Plots,showing,moments,at,the,last,timestep,of,the,simulation,with,Nelement,=,16,in,figure,7.,We,plot,(a),the,symbolic,,target,ion,density,ni,,(b),the,numerically,computed,ion,density,ni,,(c),the,symbolic,,target,neutral,density,nn,,(d),the,numerically,computed,neutral,density,nn,,Note,that,the,scales,on,each,plot,differ,,and,that,the,error,in,concentrated,near,the,boundaries,in,z.,physics,simulation,software,for,modelling,plasmas,in,the,tokamak,edge.,Some,of,the,features,that,would,be,desirable,in,generalising,the,model,are,available,in,other,codes,,such,as,GKEYLL,[25],and,COGENT,[26],,making,it,plausible,that,our,model,could,be,developed,in,a,relatively,short,timeframe.,Important,features,to,be,developed,in,a,generalisation,of,the,model,could,include,the,following:,a,fluid,or,kinetic,electron,model;,a,vorticity,equation,for,computing,the,electric,field,in,the,closed,field,line,regions;,nontrivial,tokamak,geometry,,including,a,separatrix;,more,complicated,collision,models,including,ion-ion,collisions;,more,accurate,sheath,models;,and,electromagnetic,fluctuations.,The,manufactured,solution,tests,that,we,have,developed,here,were,critical,in,the,rapid,implementation,of,the,code,,allowing,us,to,quickly,check,the,correct,implementation,of,new,features,as,they,were,added.,For,example,,in,the,results,shown,in,this,report,,the,tests,have,revealed,the,difficulties,that,the,numerical,integration,9,CONCLUSIONS,AND,FUTURE,PLANS,21,(a),(c),(b),(d),Figure,9:,Plots,showing,the,fields,at,the,last,timestep,of,the,simulation,with,Nelement,=,16,in,figure,7.,We,plot,(a),the,symbolic,,target,Er,,(b),the,numerically,computed,Er,,and,finally,(c),the,symbolic,,target,Ez,and,(d),the,numerically,computed,Ez.,Note,that,the,scales,on,each,plot,differ,,and,that,the,error,in,concentrated,near,the,boundaries,in,z,,and,are,larger,in,the,electric,fields,than,the,moments,shown,in,figure,9.,algorithm,has,with,dealing,with,the,wall,boundary,condition,–,the,wide,coverage,of,the,implemented,tests,means,that,we,should,be,able,to,easily,pinpoint,the,area,of,the,algorithm,that,must,be,addressed.,The,results,from,the,manufactured,solutions,tests,are,not,only,useful,in,determining,if,there,is,a,problem,with,the,algorithm:,the,tests,also,give,information,about,the,resolution,needed,to,resolve,different,types,of,position,space,or,velocity,space,features.,The,main,priority,for,future,development,of,the,2D,implementation,should,be,to,track,down,and,fix,the,numerical,issue,that,leads,to,poor,convergence,with,velocity,space,resolution,when,wall,boundary,conditions,are,imposed,,see,figure,7,and,the,corresponding,discussion.,Preliminary,evidence,suggests,that,the,error,is,due,to,the,treatment,of,the,discontinuity,of,the,ion,distribution,function,at,z,=,±Lz/2,and,v∥,−,Er/Bz,=,0,that,is,introduced,by,the,wall,boundary,condition.,However,,the,dominant,errors,are,observed,in,the,electric,fields,and,the,ion,and,neutral,densities,,APPENDIX,A,SUPPORTING,DOCUMENTATION,FOR,THE,MANUFACTURED,SOLUTIONS,TESTS22,rather,than,in,the,distribution,function,itself.,Determining,the,source,of,these,numerical,errors,is,the,object,of,future,work.,The,presence,of,the,automated,functional,tests,of,the,1D-1V,code,and,the,manufactured,solutions,tests,in,the,2D-3V,code,will,make,this,software,attractive,to,scientists,in,the,plasma,physics,community,that,desire,a,robust,simulation,tool,with,automated,self-testing.,With,this,in,mind,,we,have,made,a,start,in,automating,the,manufactured,solutions,tests.,In,the,ideal,case,,a,manufactured,solution,test,would,be,run,every,time,a,revision,was,made,to,the,code,,and,the,result,would,be,returned,to,the,developer.,This,is,already,possible,with,the,1D-1V,functional,tests.,However,,in,the,case,of,the,manufactured,solutions,tests,on,the,2D-3V,model,,significant,computing,resources,are,required,across,multiple,cores,to,perform,a,2D-3V,model,test.,In,the,future,,computing,frameworks,will,be,needed,in,order,to,facilitate,rapid,automated,testing.,Appendix,A.,Supporting,documentation,for,the,manufactured,solutions,tests,simulations,used,to,The,report,were,generated,by,the,branch,https://github.com/mabarnes/moment_kinetics/tree/,radial-vperp-standard-DKE-with-neutrals,,with,the,latest,commit,at,the,time,of,writing,being,19ba012d11c5a0fb2c55d45547fa070ec27f30e7.,the,data,presented,in,this,create,In,this,appendix,we,give,URL,links,to,the,input,files,used,to,generate,the,simulation,data.,To,run,a,simulation,use,the,following,command,$,julia,-O3,--project,run_moment_kinetics.jl,input.toml,To,post,process,the,simulation,data,and,to,generate,the,plots,in,this,report,run,the,following,command:,$,julia,-O3,--project,run_MMS_test.jl,Appendix,A.1.,2D-3V,periodic,boundary,condition,tests,without,collisions,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_nel_r_2_z_2_vpa_16_vperp_16.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_nel_r_2_z_2_vpa_8_vperp_8.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_nel_r_2_z_2_vpa_4_vperp_4.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_nel_r_2_z_2_vpa_2_vperp_2.toml,APPENDIX,A,SUPPORTING,DOCUMENTATION,FOR,THE,MANUFACTURED,SOLUTIONS,TESTS23,Appendix,A.2.,2D-3V,periodic,boundary,condition,tests,with,collisions,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_cxiz_nel_r_2_z_2_vpa_16_vperp_16.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_cxiz_nel_r_2_z_2_vpa_8_vperp_8.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_cxiz_nel_r_2_z_2_vpa_4_vperp_4.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-sound-wave_cheb_cxiz_nel_r_2_z_2_vpa_2_vperp_2.toml,Appendix,A.3.,1D-3V,Wall,boundary,condition,test,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_1_z_2_vpa_16_vperp_16.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_1_z_2_vpa_8_vperp_8.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_1_z_2_vpa_4_vperp_4.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_1_z_2_vpa_2_vperp_2.toml,Appendix,A.4.,1D-3V,Wall,boundary,condition,and,electron,sheath,test,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals-with-sheath_nel_r_1_z_2_vpa_16_vperp_16.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals-with-sheath_nel_r_1_z_2_vpa_8_vperp_8.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals-with-sheath_nel_r_1_z_2_vpa_4_vperp_4.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals-with-sheath_nel_r_1_z_2_vpa_2_vperp_2.toml,Appendix,A.5.,2D-3V,Wall,boundary,condition,test,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_2_z_2_vpa_16_vperp_16.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_2_z_2_vpa_12_vperp_12.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,runs/2D-wall_cheb-with-neutrals_nel_r_2_z_2_vpa_8_vperp_8.toml,https://github.com/mabarnes/moment_kinetics/blob/radial-vperp-standard-DKE-with-neutrals/,APPENDIX,A,SUPPORTING,DOCUMENTATION,FOR,THE,MANUFACTURED,SOLUTIONS,TESTS24,runs/2D-wall_cheb-with-neutrals_nel_r_2_z_2_vpa_4_vperp_4.toml,[1],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–03–01,[2],Barnes,M,,Parra,F,I,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–02–02,[3],Barnes,M,,Parra,F,I,,Hardman,M,R,and,Omotani,J,2021,Excalibur/Neptune,Report,2047357–,TN–04,[4],Barnes,M,,Parra,F,I,,Hardman,M,R,and,Omotani,J,2021,Excalibur/Neptune,Report,2047357–,TN–06,[5],Barnes,M,,Parra,F,I,,Hardman,M,R,and,Omotani,J,2021,Excalibur/Neptune,Report,2047357–,TN–08,[6],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–05–01,[7],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–07–01,[8],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–09–01,[9],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–11–01,[10],URL,https://github.com/JuliaMath/Interpolations.jl,[11],Knudsen,M,1916,Annal.,Phys.,353,1113,[12],Stangeby,P,C,2000,The,Plasma,Boundary,of,Magnetic,Fusion,Devices,(Bristol,and,Philadelphia:,Institute,of,Physics,Publishing),[13],Shu,C,W,and,Osher,S,1988,J.,Comp.,Phys.,77:439–471,[14],Gottlieb,S,and,Shu,C,W,1998,Mathematics,of,Computation,67:73–85,[15],Gottlieb,S,,Shu,C,W,and,Tadmor,E,2001,SIAM,Rev.,43:89,[16],Durran,D,R,1999,Numerical,methods,for,wave,equations,in,geophysical,fluid,dynamics,(Springer),[17],Abramowitz,M,and,Stegun,I,1972,Handbook,of,Mathematical,Functions,with,Formulas,,Graphs,,and,Mathematical,Tables,(New,York:,Dover),[18],Frigo,M,and,Johnson,S,G,2005,Special,issue,on,Program,Generation,,Optimization,,and,Platform,Adaptation,Proceedings,of,the,IEEE,93(2):216–231,[19],Clenshaw,C,W,and,Curtis,A,R,1960,Numerische,Mathematik,2:197,[20],Byrne,S,,Wilcox,L,C,and,Churavy,V,2021,JuliaCon,Proceedings,1(1),,68,[21],Parra,F,I,,Barnes,M,and,Hardman,M,R,2021,Excalibur/Neptune,Report,2047357–TN–01–02,[22],Barnes,M,2022,Excalibur/Neptune,Report,2047357–TN–12,[23],Gowda,S,,Ma,Y,,Cheli,A,,Gwozdz,M,,Shah,V,B,,Edelman,A,and,Rackauckas,C,2021,arXiv,preprint,arXiv:2105.03949,[24],URL,https://symbolics.juliasymbolics.org/dev/,[25],Hakim,A,,Mandell,N,,Bernard,T,,Francisquez,M,,Hammett,G,and,Shi,E,L,2020,Phys.,Plasmas,27,042304,[26],Dorf,M,,Dorr,M,,Hittinger,J,,Cohen,R,and,Rognlien,T,2016,Phys.,Plasmas,23,056102

:pdfembed:`src:_static/TN-13-2_NumericalReducedModelCoupling2D2VDriftKineticIons2D3VKineticNeutralsHelicalMagne.pdf, height:1600, width:1100, align:middle`