TN-04_Numerical11DMomentBasedDriftKineticPeriodicBoundaryConditions
===================================================================

.. meta::
   :description: technical note

   :keywords: Report,2047357-TN-02,,M2.2,Numerical,study,of,1+1D,,moment-based,drift,kinetic,models,with,periodic,boundary,conditions,M.,Barnes1,,F.,I.,Parra1,,M.,R.,Hardman1,and,J.,Omotani2,1,Rudolf,Peierls,Centre,for,Theoretical,Physics,,University,of,Oxford,,Clarendon,Laboratory,,Parks,Road,,Oxford,OX1,3PU,,United,Kingdom,2,Culham,Centre,for,Fusion,Energy,,Culham,Science,Centre,,Abingdon,,Oxon,,OX14,3DB,,United,Kingdom,E-mail:,michael.barnes@physics.ox.ac.uk,1.,Introduction,We,expect,that,one,of,the,biggest,challenges,in,numerically,solving,drift,kinetic,equations,in,the,plasma,edge,is,treating,the,motion,of,electrons,along,the,magnetic,field.,Because,the,electrons,are,light,,they,move,rapidly,along,the,field,,placing,a,severe,stability,restriction,on,the,step,size,for,explicit,time,advance,schemes.,Unfortunately,,an,implicit,treatment,is,not,straightforward,due,to,an,implicit,dependence,of,the,electrostatic,potential,on,the,charged,particle,distribution,functions.,One,of,the,main,aims,of,our,research,is,to,develop,and,test,a,novel,analytical,model,and,associated,numerical,algorithm,for,relaxing,this,restriction.,As,a,first,step,towards,this,goal,,we,developed,a,new,code,in,the,programming,language,Julia,to,simulate,a,simple,model,for,parallel,dynamics,(described,in,our,Feb,2021,report,[1]),without,the,novel,moment-based,approach,that,we,intend,to,explore.,We,have,now,extended,the,code,to,simulate,a,modified,set,of,equations,in,which,the,density,is,removed,from,the,particle,distribution,function,and,is,evolved,separately,using,the,continuity,equation.,In,this,report,we,give,an,overview,of,this,moment-based,approach,and,discuss,the,numerical,issues,around,its,implementation.,2.,Model,equations,A,detailed,derivation,of,the,model,we,consider,is,provided,in,our,Jan,2021,report,[2].,Here,we,reproduce,a,brief,overview,of,the,model,from,our,Feb,2021,report,[1],for,the,Reader’s,convenience.,The,model,we,consider,consists,of,a,single,ion,species,of,charge,e,,a,single,neutral,species,,and,an,electron,species,modelled,as,having,a,Boltzmann,response,,all,immersed,in,a,straight,,uniform,magnetic,field,in,the,z,direction.,We,allow,for,charge,exchange,collisions,between,ions,and,neutrals,but,do,not,account,for,intra-species,collisions.,Finally,,we,assume,that,the,plasma,is,homogeneous,in,the,plane,perpendicular,to,the,magnetic,field.,With,these,assumptions,,our,model,system,of,1,2,(1),(2),(3),(4),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,equations,is,and,.,=,∂fi,∂t,+,v,(cid:107),∂fi,∂z,−,∂fn,∂t,+,v,(cid:107),e,mi,∂fn,∂z,∂φ,∂z,∂fi,∂v,(cid:107),=,Rin,(nnfi,−,−,nifn),,,Rin,(nifn,−,nnfi),,,=,−,∞,(cid:90),−∞,ns(z,,t),=,dv,(cid:107),fs(z,,v,,,t),,(cid:107),ni,=,Ne,exp,,,eφ,Te,(cid:19),(cid:18),⊥,v,⊥,dϑdv,and,v,Fs,the,marginalized,particle,distribution,function,for,species,s,,with,fs,the,components,of,the,particle,velocity,parallel,and,perpendicular,to,the,v,(cid:82),(cid:107),magnetic,field,,respectively,,ϑ,the,gyro-angle,,mi,the,ion,mass,,t,the,time,,φ,the,electrostatic,potential,,and,Rin,the,charge,exchange,collision,frequency.,⊥,,,with,For,our,boundary,conditions,,we,impose,periodicity,on,fs,in,both,z,and,v,(cid:107),periods,Lz,and,Lv(cid:107),,respectively.,There,is,also,the,option,to,impose,zero,boundary,conditions,on,z,and,v,at,the,upwind,boundary,of,the,domain.,As,fs,should,go,to,zero,(cid:107),,,imposition,of,zero,boundary,conditions,and,periodic,boundary,conditions,at,v,should,be,equivalent,as,long,as,Lv(cid:107),is,sufficiently,large.,Note,that,with,either,choice,of,boundary,conditions,,the,line-averaged,density,We,normalize,Eqs.,(1)-(4),by,defining,Lz,0,dz,ns,is,conserved.,(cid:107),→,±∞,(cid:82),˜fs,.,=,fs,,,,,,,vth,i√π,Ne,vth,i,Lz,z,Lz,v,(cid:107),vth,i,ns,Ne,eφ,Te,,,,,,,.,=,t,˜t,.,=,˜z,.,=,˜v,(cid:107),.,=,.,=,˜ns,˜φ,and,with,vth,i,.,=,˜Rin,.,=,Rin,NeLz,vth,i,2Te/mi.,In,terms,of,these,normalised,quantities,,Eqs,(1)-(4),become,(cid:112),∂,˜fi,∂˜t,∂,˜fi,∂,˜z,−,∂,˜φ,∂,˜z,1,2,+,˜v,(cid:107),∂,˜fi,∂˜v,(cid:107),=,˜Rin,−,˜nn,˜fi,−,(cid:16),˜ni,˜fn,,,(cid:17),(5),(6),(7),(8),(9),(10),(11),(12),(cid:18),∂nn,∂t,+,v,(cid:107),∂nn,∂z,(cid:19),,,eφ,Te,(cid:19),=,0,,ni,=,Ne,exp,∂ns,∂t,+,(cid:18),∂nsus,∂z,∞,dv,(cid:107),gsv,.,(cid:107),(cid:90),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,∂,˜fn,∂˜t,and,+,˜v,(cid:107),=,∂,˜fn,∂,˜z,e,˜φ,=,˜ni,=,˜Rin,−,1,√π,˜ni,(cid:16),∞,˜fi,,,(cid:17),˜nn,˜fn,−,˜fi,,d˜v,(cid:107),(cid:90),−∞,˜nn,=,1,√π,∞,(cid:90),−∞,˜fn.,d˜v,(cid:107),3,(13),(14),(15),2.1.,Moment,approach:,density,We,now,define,the,modified,distribution,function,gs,terms,of,gs,,the,system,of,equations,given,by,Eqs.,(1)-(4),becomes,.,=,fs/ns,so,that,gs,=,1.,In,dv,(cid:107),(cid:82),∂gi,∂t,+,v,(cid:107),∂gi,∂z,−,e,mi,∂φ,∂z,∂gi,∂v,(cid:107),(cid:19),ni,(cid:18),+,gi,∂ni,∂t,+,v,(cid:107),∂ni,∂z,(cid:19),=,Rinninn,(gi,−,−,gn),,,∂gn,∂t,+,v,∂gn,∂z,(cid:107),nn,(cid:18),+,gn,(cid:19),(cid:18),=,Rinninn,(gn,−,−,gi),,,and,us,=,−∞,Note,that,the,1D,continuity,equation,(19),has,replaced,the,moment,equation,(3),as,a,means,of,computing,the,density,for,each,species.,Substituting,the,continuity,equation,(19),into,the,drift,kinetic,equations,(16),and,(17),gives,∂gi,∂t,+,v,(cid:107),∂gi,∂z,−,e,mi,∂φ,∂z,∂gi,∂v,(cid:107),=,Rinnn,(gi,−,−,gn),+,gi,∂ui,∂z,−,(cid:18),v,(cid:107),−,ui,∂,ln,ni,∂z,(21),(cid:19),∂gn,∂t,+,v,(cid:107),∂gn,∂z,=,Rinni,(gn,−,−,gi),+,gn,∂un,∂z,−,v,(cid:107),−,un,(cid:18),We,normalize,Eqs.,(18)-(22),by,using,Eqs.,(6)-(11),and,by,further,defining,(cid:0),(cid:1),,,(22),(cid:0),(cid:1),∂,ln,nn,∂z,(cid:19),.,=,gsvth,i√π,˜gs,.,=,˜us,us,vth,i,.,(23),(24),and,and,In,terms,of,these,normalised,quantities,,Eqs,(18)-(22),become,∂˜gi,∂˜t,+,˜v,(cid:107),∂˜gi,∂,˜z,−,∂,˜φ,∂,˜z,1,2,∂˜gi,∂˜v,(cid:107),=,˜Rin˜nn,(˜gi,−,−,˜gn),+,˜gi,∂,˜ui,∂,˜z,−,(cid:18),˜v,(cid:107),−,˜ui,(cid:0),(cid:1),∂,ln,˜ni,∂,˜z,(cid:19),,,(25),(16),(17),(18),(19),(20),4,(26),(27),(28),(29),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,∂˜gn,∂˜t,+,˜v,(cid:107),∂˜gn,∂,˜z,=,−,∂,˜un,∂,˜z,−,(cid:18),=,0,,˜v,(cid:107),−,˜un,(cid:0),(cid:1),∂,ln,˜nn,∂,˜z,,,(cid:19),˜Rin˜ni,(˜gn,−,∂,˜ns,∂˜t,˜gi),+,˜gn,+,∂,˜ns,˜us,∂,˜z,e,˜φ,=,˜ni,,and,˜us,=,1,√π,∞,d˜v,(cid:107),˜gn˜v,(cid:107),.,(cid:90),−∞,The,above,form,for,the,equations,is,appealing,because,it,maintains,the,form,of,an,advection,equation,with,the,only,modification,being,the,addition,of,source,terms.,However,,it,can,pose,challenges,for,numerical,conservation,of,quantities,such,as,the,0th,velocity,moment,of,gs.,This,is,because,parts,of,the,source,terms,must,cancel,upon,velocity,space,integration,with,some,of,the,advective,terms.,To,ease,the,task,of,preserving,conservation,properties,numerically,,the,equations,can,be,manipulated,into,the,following,form,in,which,such,cancellations,can,be,built,into,the,discretisation:,∂˜gi,∂˜t,+,˜v,(cid:107),˜ni,∂,˜ni˜gi,∂,˜z,−,∂,˜φ,∂,˜z,1,2,∂˜gi,∂˜v,(cid:107),=,˜Rin˜nn,(˜gi,−,−,˜gn),+,˜gi,˜ni,∂,˜ni,˜ui,∂,˜z,,,∂˜gn,∂˜t,+,˜v,(cid:107),˜nn,∂,˜nn˜gn,∂,˜z,=,∂,˜ns,∂˜t,˜gi),+,˜gn,˜nn,∂,˜nn,˜un,∂,˜z,,,−,˜Rin˜ni,(˜gn,−,∂,˜ns,˜us,∂,˜z,e,˜φ,=,˜ni,,+,=,0,,and,˜us,=,1,√π,∞,(cid:90),−∞,d˜v,(cid:107),˜gn˜v,(cid:107),.,3.,Numerical,implementation,(30),(31),(32),(33),(34),The,algorithms,described,in,this,Section,have,been,implemented,in,the,code,,written,in,the,Julia,programming,language,,currently,available,on,GitHub,at,https://github.,com/mabarnes/moment_kinetics.,3.1.,Time,advance,We,evolve,Eqs.,(12)-(15),or,Eqs.,(30)-(34),using,a,member,of,the,family,of,Strong,Stability,Preserving,(SSP),Runge-Kutta,(RK),schemes;,see,,e.g.,,[3,,4,,5].,Current,SSPRK,options,implemented,in,the,code,are,SSPRK1,(forward,Euler),,SSPRK2,(Heun’s,method),SSPRK3,(Shu-Osher,method),and,four-stage,SSPRK3.,The,user,can,also,specify,the,use,of,‘flip-flop’,Lie,operator,splitting.,Operator,splitting,limits,the,time,advance,scheme,to,second,order,accuracy,in,step,size,,but,could,be,useful,for,separately,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,5,treating,different,pieces,of,physics.,Here,we,describe,the,current,default,option,,which,is,the,four-stage,SSPRK3,method,without,operator,splitting.,For,convenience,of,notation,,we,express,the,normalised,drift,kinetic,equations,for,the,ions,and,neutrals,in,the,vector,form,∂f,∂t,=,G[f,],,(35),with,f,the,solution,vector,containing,the,evolved,quantities:,f,=,if,Eqs.,(12)-,(15),are,solved,,and,f,=,(˜gi,,˜gn,,˜ni,,˜nn)T,if,Eqs.,(30)-(34),are,solved.,The,operator,G,accounts,for,parallel,streaming,,parallel,acceleration,(for,the,ions),and,charge,exchange,collisions,,as,well,as,for,the,divergence,of,the,particle,flux,if,the,density,is,separately,evolved.,The,four-stage,SSPRK3,method,for,advancing,this,system,of,equations,is,3rd,order,accurate,in,time,step,size,∆t,,with,a,Courant,number,of,two.,It,is,given,by,(cid:16),(cid:17),T,˜fi,,˜fn,f,(1),=,f,(2),=,f,(3),=,f,n+1,=,1,2,1,2,2,3,1,2,f,n,+,1,2,1,f,(1),+,2,1,(cid:0),f,(2),+,6,1,f,(3),+,2,f,n,+,(f,n,+,∆tG,[f,n]),,,f,(1),+,∆tG,f,(1),,,(cid:2),(cid:3)(cid:1),f,(2),+,∆tG,1,6,(cid:0),f,(3),+,∆tG,f,(3),(cid:2),,,(cid:3)(cid:1),(cid:2),(cid:3)(cid:1),f,(2),,,(36),where,the,superscript,n,denotes,the,time,level.,(cid:0),We,have,tested,our,implementation,of,SSP,RK2,and,4-stage,SSP,RK3,by,calculating,the,rms,error,in,the,distribution,function,after,it,is,advected,in,one,dimension,with,constant,advection,speed,˜v,=,1,for,10,transits,of,the,domain,with,length,L,=,1:,.,=,(cid:15)rms,(cid:118),(cid:117),(cid:117),(cid:116),1,Nz,Nz,j=1,(cid:88),fi(zj,,t,=,10),|,−,2.,fi(zj,,t,=,0),|,(37),The,results,when,paired,with,a,finite,difference,discretisation,(third,order,upwind),are,given,in,Fig.,1.,For,Chebyshev,pseudospectral,on,a,single,element,with,4-stage,SSPRK3,,see,Fig.,2.,3.2.,Spatial,discretisation,There,are,two,discretisation,schemes,implemented,in,the,code:,finite,differences,and,Chebyshev,(pseudo)spectral,elements.,The,user,can,choose,at,run-time,which,scheme,to,use,for,each,of,the,z,and,v,coordinates.,(cid:107),3.2.1.,Finite,difference,discretisation.,For,the,finite,difference,discretisation,,the,corresponding,coordinate,grid,is,uniform,on,the,domain,[,L/2,,L/2],,with,L,the,coordinate,box,length.,The,default,method,employed,for,derivatives,is,3rd,order,upwind,−,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,6,Figure,1.,RMS,error,as,a,function,of,time,step,size,∆t,and,varying,values,for,Nz,for,both,SSP,RK2,(solid,lines,and,circles),and,for,4-stage,SSP,RK3,(dashed,lines,and,squares).,Due,to,the,CFL,restriction,that,ties,temporal,resolution,to,spatial,resolution,,the,range,in,∆t,over,which,time,domain,errors,dominate,is,limited,for,RK2,and,is,effectively,non-existent,for,RK3.,differences,,though,1st,and,2nd,order,schemes,are,also,available,as,options.,For,an,overview,of,upwind,differences,and,a,discussion,of,the,merits,of,the,different,upwind,schemes,,see,,e.g.,[6].,The,associated,integration,weights,used,for,field-line,averages,in,z,integration,required,for,obtaining,fields/moments,are,obtained,using,and/or,for,the,v,(cid:107),the,composite,Simpson’s,rule,(sometimes,referred,to,as,composite,Simpson’s,1/3,rule):,L,0,(cid:90),dx,f,(x),h,3,≈,(N,−,1)/2,j=1,(cid:88),(f,(x2j,1),+,4f,(x2j),+,f,(x2j+1)),,,−,(38),where,N,is,the,number,of,grid,points,in,the,coordinate,x,,and,h,=,L/(N,1),is,the,uniform,grid,spacing.,The,composite,rule,(38),is,only,applicable,for,N,odd,,so,it,is,supplemented,at,the,boundary,by,Simpson’s,3/8,rule,when,N,is,even.,−,3.2.2.,Chebyshev,spectral,elements.,When,using,Chebyshev,spectral,elements,,the,corresponding,coordinate,grid,is,the,Gauss-Chebyshev-Lobatto,grid,on,each,element.,For,a,description,of,Chebyshev-Gauss,quadrature,,see,,e.g.,[7].,Inclusion,of,the,endpoints,1×10−91×10−81×10−71×10−61×10−50.00010.0010.010.000100.001000.01000(cid:15)rms∆t50100200400800∝∆t2,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,7,Figure,2.,RMS,error,as,a,function,of,time,step,size,∆t,for,4-stage,SSP,RK3,with,a,Chebyshev,pseudospectral,discretisation.,within,each,element,facilitates,enforcement,of,continuity,at,element,boundaries,,and,the,use,of,Chebyshev,polynomials,as,a,basis,enables,the,use,of,Fast,Fourier,Transforms.,In,our,code,,these,transforms,are,done,using,the,widely-used,FFTW,library,[8].,The,associated,integration,weights,used,for,field-line,averages,in,z,and/or,for,the,v,(cid:107),integration,required,for,obtaining,fields/moments,are,obtained,using,Clenshaw-Curtis,quadrature,rules,[9].,Clenshaw-Curtis,quadrature,is,convenient,,as,it,allows,for,the,use,of,endpoints,in,the,integration,domain,(which,is,dictated,by,the,use,of,a,Gauss-,Chebyshev-Lobatto,grid),while,still,exactly,integrating,polynomials,up,to,degree,N,1,,with,N,the,number,of,points,within,the,element.,−,A,1D,advection,test,demonstrating,the,spectral,accuracy,of,the,Chebyshev,scheme,on,a,single,element,is,given,in,Fig.,3,,where,the,rms,error,is,given,by,Eq.,(37).,The,maximum,stable,time,step,subject,to,the,CFL,restriction,is,plotted,as,a,function,of,the,number,of,z,grid,points,on,a,single,element,in,Fig.,4,and,as,a,function,of,the,number,of,elements,Nelem,with,Nz,=,9,fixed,in,Fig.,5.,Slight,deviations,from,the,expected,scalings,are,likely,due,to,the,numerical,dissipation,that,is,introduced,by,the,use,of,the,derivative,from,the,upwind,element,at,the,overlapping,point,at,element,boundaries,and,at,the,boundary,of,the,periodic,domain.,1×10−91×10−81×10−71×10−61×10−50.00010.0010.01(cid:15)rms∆tspectral∝∆t3,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,8,Figure,3.,RMS,error,as,a,function,of,the,number,of,grid,points,Nz.,The,time,advance,scheme,used,is,4-stage,SSP,RK3.,3.3.,conservation,properties,g,=,1.,The,field-line-averaged,density,,We,consider,two,approaches,to,ensuring,that,these,properties,are,preserved,by,the,numerical,scheme:,carefully,chosen,discretisation,of,the,equations,in,conservative,form,and,conserving,corrections,applied,at,the,end,of,each,time,level.,dz,n,,should,be,conserved,,as,should,dv,(cid:107),(cid:82),(cid:82),3.3.1.,conservative,differencing,To,conserve,particle,number,,the,discretisation,must,satisfy,0,=,∂,∂t,Nelem,j=1,(cid:88),(cid:90)Ωj,dzj,nj,=,−,dzj,∂Γj,∂zj,(cid:90)Ωj,Nelem,Nz,(39),wj,k,∂Γj,∂zj,,,(cid:19)k,(cid:18),=,−,j=1,(cid:88),(cid:88)k=1,where,Γ,=,nu,is,the,particle,flux,,Nelem,is,the,number,of,elements,,Nz,is,the,number,of,grid,points,per,element,,Ωj,is,the,domain,in,z,corresponding,to,the,jth,element,,the,superscript,j,denotes,evaluation,on,the,jth,element,,and,wj,k,is,the,integration,weight,1×10−101×10−81×10−60.00010.0116810121416(cid:15)rmsNzdata12(4Nz)Nz,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,9,Figure,4.,Maximum,stable,time,step,subject,to,the,CFL,restriction,as,a,function,of,the,number,of,Gauss-Chebyshev-Lobatto,grid,points.,The,max,stable,time,step,scales,a,bit,more,weakly,than,1/N,2,z,,,as,expected.,associated,with,the,kth,grid,point,within,the,jth,element.,For,expression,(39),to,be,satisfied,,there,are,three,requirements,that,must,be,met.,First,,the,contributions,to,the,integral,from,interior,points,on,each,element,must,cancel,with,one,another.,This,requires,a,centred,difference,scheme,if,finite,differences,are,used.,As,shown,in,Appendix,A,,it,is,automatically,satisfied,if,the,flux,derivative,is,obtained,using,the,Chebyshev,pseudo-spectral,method,,and,the,integral,is,computed,using,Clenshaw-Curtis,quadrature.,Second,,the,fluxes,must,be,continuous,across,element,boundaries.,For,this,to,be,guaranteed,,the,derivative,of,the,flux,itself,must,be,single-valued,at,the,element,boundaries.,This,forces,one,to,choose,the,flux,derivative,at,the,boundary,between,elements,to,be,the,average,of,the,flux,derivatives,on,each,element.,Finally,,the,fluxes,must,vanish,or,cancel,one,another,at,the,boundaries,of,the,simulation,domain.,This,forces,one,to,use,either,zero,boundary,conditions,or,to,again,use,the,average,of,the,flux,derivatives,on,each,element,with,the,assumption,of,periodicity.,Next,we,consider,how,to,ensure,g,=,1.,The,terms,involving,z,derivatives,dv,(cid:107),of,g,on,each,side,of,the,drift,kinetic,equation,will,cancel,as,long,as,the,differencing,is,done,consistently,for,each,term,,as,well,as,using,g,and,u,at,the,same,time,level.,The,charge,exchange,term,should,also,vanish,as,long,as,gi,and,gn,are,used,at,the,same,time,(cid:82),0.0010.010.1110∆tCFLNzdata10N−9/5z,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,10,Figure,5.,Maximum,stable,time,step,subject,to,the,CFL,restriction,as,a,function,of,the,number,of,elements,with,the,number,of,grid,points,per,element,fixed,at,9.,The,minimum,grid,spacing,scales,inversely,with,the,number,of,elements,,leading,to,a,max,stable,time,step,that,is,inversely,proportional,to,the,number,of,elements,Nelem.,(cid:107),(cid:107),−,g(v,(cid:107),(∂g/∂v,(cid:107),level.,The,only,non-trivial,term,then,is,the,parallel,acceleration.,For,this,term,to,vanish,integration,,the,fundamental,theorem,of,calculus,must,hold,numerically;,i.e.,,upon,v,,min,are,the,maximum,and,dv,,max,and,v,,min),=,0,,where,v,,max),),=,g(v,(cid:107),(cid:107),(cid:107),(cid:107),included,in,the,grid.,For,the,last,equality,,the,values,of,g,at,the,minimum,values,of,v,(cid:82),(cid:107),must,be,the,same.,This,could,be,achieved,either,with,periodic,BCs,or,boundaries,in,v,with,zero,BCs.,However,,there,is,the,complication,that,the,drift,kinetic,equation,drives,non-periodic,solutions,that,only,decay,to,zero,at,v,.,To,enforce,periodicity,g,,we,set,the,boundary,values,to,be,the,average,of,the,nominal,without,modifying,values,at,the,endpoints,in,v,(cid:107),Taken,together,,these,requirements,remove,the,possibility,of,using,numerical,dissipation,via,,e.g.,,upwind,differencing,or,upwind,fluxes,at,element,boundaries,,to,improve,numerical,stability.,This,is,problematic,for,advection-dominated,problems,like,the,one,we,are,considering,,as,in,the,absence,of,some,form,of,numerical,dissipation,the,scheme,is,unstable.,(cid:107),→,±∞,dv,(cid:107),(cid:82),.,0.010.1110∆tCFLNelemdata1/(5Ne),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,11,3.3.2.,conserving,corrections,The,second,approach,to,ensuring,numerically,the,desired,conservation,properties,is,to,correct,the,numerical,solutions,for,n,and,g,at,the,end,of,each,time,step.,For,the,density,,one,can,set,nm+1,=,ˆnm+1,+,nm,1,(cid:18),where,ˆnm+1,is,the,updated,solution,(at,time,level,m,+,1),to,the,continuity,equation,nm),=,0.,before,applying,any,conserving,correction.,This,guarantees,that,Note,that,the,superscripts,here,refer,to,the,time,level,,not,the,element,index.,The,additional,error,in,the,density,introduced,by,this,correction,is,−,(cid:82),(cid:82),dz,(nm+1,(40),−,(cid:19),(cid:82),,,dz,ˆnm+1,dz,nm,nm,1,(cid:18),dz,ˆnm+1,dz,nm,−,(cid:82),(cid:82),(cid:19),=,nm,=,nm,1,(cid:32),−,(cid:82),dz,(cid:15)m,dz,nm,dz,=,(cid:0),(cid:82),O,nm+1,exact,+,(cid:15)m,dz,nm,((cid:15)m),,(cid:33),(cid:1),(41),(cid:82),(cid:82),where,(cid:15)m,is,the,error,due,to,numerical,discretisation,,and,nm+1,in,the,limit,(cid:15)m,=,0.,exact,is,the,solution,for,ˆnm+1,A,similar,technique,can,be,applied,to,conserve,g.,In,particular,,we,set,dv,(cid:107),gm+1,=,ˆgm+1,+,gm,1,(cid:18),where,ˆgm+1,is,the,updated,solution,to,the,drift,kinetic,equation,before,applying,any,gm,=,1.,Again,,conserving,correction.,This,ensures,(δm),,where,δm,is,the,the,additional,error,in,g,associated,with,this,correction,is,(cid:82),discretisation,error.,gm+1,=,1,,provided,ˆgm+1,dv,(cid:107),dv,(cid:107),(42),O,−,(cid:19),(cid:90),(cid:82),,,(cid:82),dv,(cid:107),This,approach,is,simple,,does,not,change,the,order,of,accuracy,of,the,discretisation,scheme,and,allows,for,the,use,of,numerical,dissipation,to,improve,numerical,stability,properties.,Results,showing,its,efficacy,are,given,in,Sec.,4,4.,Numerical,results,To,benchmark,our,numerical,the,moment-based,approach,implementation,of,encapsulated,in,Eqs.,(30)-(34),,we,compare,our,simulation,results,with,the,analytical,benchmarks,developed,in,[2],and,with,the,numerical,results,obtained,by,directly,solving,the,kinetic,system,corresponding,to,Eqs.,(12)-(15).,The,results,reported,here,were,obtained,using,the,conservative,differencing,detailed,in,Sec.,3.3.1,for,the,continuity,equation,and,for,the,source,terms,,while,the,conserving,correction,given,by,Eq.,42,is,applied,to,ensure,that,gs,=,1.,We,have,initialised,the,distribution,functions,for,the,ions,and,neutrals,to,be,of,the,dv,(cid:107),(cid:82),form,˜gs,=,1/2,Te,T,s,(cid:19),(cid:18),exp,˜v2,(cid:107),−,(cid:18),Te,T,s,(cid:19),,,(43),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,12,with,ns,=,ns,+,δns,and,Ts,=,T,s,+,δTs,,and,an,overline,denoting,a,field,line,average.,The,piece,of,the,density,that,varies,along,z,,δns,,is,chosen,to,be,small,compared,to,ns,(δns/ns,=,0.001),so,that,the,system,of,equations,can,be,linearised,to,a,good,approximation.,This,facilitates,comparisons,with,the,linear,analytical,theory,presented,in,[2].,For,all,cases,shown,here,,ni,=,nn,=,Ne/2,,T,i,=,T,n,=,Te,and,mi,=,mn.,The,charge,exchange,collision,frequency,is,varied,,and,damping,rates,and,frequencies,are,extracted,by,considering,the,time,evolution,of,the,spatially-varying,component,of,the,electrostatic,potential,,δφ.,In,particular,,a,least-squares,fit,for,δφ(t)/δφ(t0),is,done,for,ϕ),to,each,simulation,to,a,function,of,the,form,exp(,obtain,the,damping,rate,γ,,frequency,ω,and,phase,ϕ.,The,results,are,given,in,Fig.,6.,There,is,good,agreement,across,a,wide,range,of,charge,exchange,collision,frequencies,,both,for,the,damping,of,finite,frequency,modes,(corresponding,to,the,solid,lines),and,to,a,zero,frequency,mode,that,appears,at,larger,collisionalities,(dashed-dotted,lines).,The,minor,discrepancy,between,the,analytical,and,numerical,damping,rates,that,is,apparent,for,the,case,with,normalised,charge,exchange,collision,frequency,near,0.7,is,due,to,the,simultaneous,presence,of,both,modes,with,similar,damping,rates.,ϕ)/,cos(ωt0,−,t0)),cos(ωt,γ(t,−,−,−,Figure,6.,Normalized,damping,rate,and,real,frequency,as,a,function,of,the,charge,exchange,collision,frequency.,Note,that,the,normalising,vth,employed,here,differs,from,the,vth,i,employed,in,the,text:,It,is,chosen,to,be,vth,=,2Ti/mi,to,facilitate,comparison,with,the,analytical,results,obtained,in,[2].,(cid:112),In,Figure,7,we,show,the,difference,in,conservation,properties,between,cases,for,which,the,mixed,conservative,differencing,and,conservative,corrections,indicated,at,the,beginning,of,the,Section,are,employed,and,those,for,which,no,conserving,correction,is,applied.,With,the,conservative,implementation,,both,the,particle,number,and,g,are,conserved,to,machine,precision,,regardless,of,numerical,resolution.,dv,(cid:107),(cid:82),5.,Future,plans,Now,that,we,have,a,proof-of-concept,implementation,of,the,moment-based,approach,to,solving,the,1+1D,kinetic,problem,,we,plan,to,extend,the,treatment,to,include,separate,00.511.522.500.511.52ω/kkvth(ni+nn)Rin/(cid:12)(cid:12)(cid:12)kkvth(cid:12)(cid:12)(cid:12)analyticalkineticmoment-based−1.4−1.2−1−0.8−0.6−0.4−0.200.511.52γ/kkvth(ni+nn)Rin/(cid:12)(cid:12)(cid:12)kkvth(cid:12)(cid:12)(cid:12)analyticalkineticmoment-based,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,13,Figure,7.,Time,evolutions,of,the,conservation,error,for,dv(cid:107)g,(left),and,for,the,particle,number,(right).,The,velocity,moment,of,g,is,guaranteed,to,be,unity,in,the,blue,case,because,the,conserving,corrections,described,in,Sec.,3.3,are,applied.,Conservative,differencing,in,the,continuity,equation,and,the,use,of,average,fluxes,at,element,boundaries,ensures,particle,conservation,for,both,simulations.,Both,simulations,were,0.3,and,used,very,low,resolution,to,run,for,T,i,=,T,n,=,Te,,normalised,Rin,enhance,the,non-conservation:,a,Chebyshev,pseudo-spectral,discretisation,was,used,with,Nz,=,5,on,one,element,and,Nv,=,9,on,five,elements.,≈,(cid:82),evolution,of,the,parallel,flow,and,parallel,pressure.,Inclusion,of,a,separately,evolved,parallel,flow,should,enable,us,to,include,kinetic,electrons,in,our,model,and,test,various,numerical,approaches,to,obtaining,the,electrostatic,potential,in,an,efficient,way.,Appendix,A.,Conservative,differencing,of,the,continuity,equation,We,start,with,the,continuity,equation,,∂n,∂t,+,∂Γ,∂z,=,0,,(A.1),and,we,expand,the,particle,flux,Γ,in,Chebyshev,polynomials,on,each,of,Nelem,elements,so,that,N,1,−,.,=,Γj,i,ˆΓj,nTn(zi),(A.2),is,the,particle,flux,at,the,ith,grid,point,within,the,jth,element,,and,Tn(zi),is,a,Chebyshev,polynomial,of,the,first,kind.,The,derivative,appearing,in,the,continuity,equation,is,then,n=0,(cid:88),∂Γj,∂z,=,N,1,−,n=0,(cid:88),ˆΓj,n,∂Tn(z),∂z,=,N,2,−,n=0,(cid:88),dj,nTn(z),,with,dj,n,related,to,ˆΓj,n,via,dj,n,−,1,=,1,cn,−,1,(cid:16),2nˆΓj,n,+,dj,n+1,,,(cid:17),(A.3),(A.4),Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,14,where,dj,N,1,=,dj,N,=,0,and,cn,=,2,if,n,=,0,and,cn,=,1,for,n,>,0.,−,The,field-line-averaged,density,should,be,conserved,according,to,the,continuity,equation.,This,can,be,achieved,by,carefully,choosing,the,discretisation,to,be,of,a,conservative,form.,In,particular,,the,discretisation,must,satisfy,0,=,∂,∂t,(cid:90)Ωj,dz,n,=,−,(cid:90)Ωj,Nz,=,−,dzj,∂Γj,∂zj,wj,i,∂Γj,∂zj,(A.5),(cid:19)i,(cid:18),Γj,1,,=,Γj,i=1,(cid:88),Nz,−,,,and,the,equality,in,the,last,where,the,integration,weights,are,given,by,the,set,line,holds,as,long,as,the,integration,scheme,is,chosen,to,exactly,integrate,polynomials,1.,The,Clenshaw-Curtis,quadrature,satisfies,this,requirement.,of,degree,less,than,Nz,−,However,,there,is,one,more,issue,to,consider:,continuity,of,the,density,across,element,boundaries.,If,we,enforce,continuity,,the,derivative,of,the,flux,must,be,single-,valued,at,the,element,boundaries.,Let,the,flux,derivative,at,the,element,boundaries,be,wh,i,},{,∂Γj,b,∂zj,=,α,∂Γj,Nz,∂zj,+,(1,α),−,∂Γj+1,1,∂zj+1,,,(A.6),with,α,[0,,1].,∈,The,time,rate,of,change,of,the,field-line-averaged,density,is,thus,∂,∂t,(cid:90),Ne,Ne,dz,n,=,dzjnj,=,j=1,(cid:90)Ωj,(cid:88),1,Ne,−,j=1,(cid:32),(cid:88),(cid:0),+,j=1,(cid:88),+,wj+1,1,wj,Nz,(cid:1),Γj,Nz,−,Γj,1,(cid:0),∂Γj,b,∂zj,−,(cid:1),wj,Nz,∂Γj,Nz,∂zj,−,wj+1,1,∂Γj+1,1,∂zj+1,.,(cid:33),(A.7),Because,Γ,is,single-valued,at,element,boundaries,,each,of,the,interior,terms,in,the,first,sum,on,the,RHS,above,cancel,,leaving,only,the,contribution,from,j,=,1,and,j,=,Nelem.,The,final,sum,is,then,the,deviation,from,exact,particle,conservation.,This,can,be,made,to,vanish,by,choosing,α,in,Eq.,(A.6),to,be,α,=,wj,Nz,+,wj+1,1,,,wj,Nz,(A.8),which,evaluates,to,1/2,if,the,grids,are,identical,in,neighbouring,elements.,The,same,treatment,must,also,be,applied,at,the,exterior,boundaries,to,ensure,exact,numerical,conservation.,Note,that,this,is,possible,with,periodic,BCs,but,not,with,a,zero,BC.,Numerical,study,of,moment-based,drift,kinetic,model,with,periodic,BCs,15,[1],M.,Barnes,,F.,I.,Parra,,and,M.,R.,Hardman.,Numerical,study,of,1d,drift,kinetic,models,with,periodic,boundary,conditions.,Excalibur/Neptune,Report,,2:2047357–TN–01–02,M2.1,,2021.,[2],F.,I.,Parra,,M.,Barnes,,and,M.,R.,Hardman.,1d,drift,kinetic,models,with,periodic,boundary,conditions.,Excalibur/Neptune,Report,,1:2047357–TN–01–02,M1.1,,2021.,[3],C.-W.,Shu,and,S.,Osher.,Efficient,implementation,of,essentially,non-oscillator,shock-capturing,schemes.,J.,Comp.,Phys.,,77:439–471,,1988.,[4],S.,Gottlieb,and,C.-W.,Shu.,Total,variation,diminishing,runge-kutta,methods.,Mathematics,of,Computation,,67:73–85,,1998.,[5],S.,Gottlieb,,C.-W.,Shu,,and,E.,Tadmor.,Strong,stability-preserving,high-order,time,discretization,methods.,SIAM,Rev.,,43:89,,2001.,[6],D.,R.,Durran.,Numerical,methods,for,wave,equations,in,geophysical,fluid,dynamics.,Springer,,1999.,[7],M.,Abramowitz,and,I.,A.,Stegun.,Handbook,of,Mathematical,Functions,with,Formulas,,Graphs,,and,Mathematical,Tables.,Dover,,New,York,,1972.,[8],Matteo,Frigo,and,Steven,G.,Johnson.,The,design,and,implementation,of,FFTW3.,Proceedings,of,the,IEEE,,93(2):216–231,,2005.,Special,issue,on,“Program,Generation,,Optimization,,and,Platform,Adaptation”.,[9],C.,W.,Clenshaw,and,A.,R,Curtis.,A,method,for,numerical,integration,on,an,automatic,computer.,Numerische,Mathematik,,2:197,,1960.

:pdfembed:`src:_static/TN-04_Numerical11DMomentBasedDriftKineticPeriodicBoundaryConditions.pdf, height:1600, width:1100, align:middle`