TN-01_CoDesignPicMethodsOtherElementsNeptune
============================================

.. meta::
   :description: technical note

   :keywords: 2047355-TN-01-3:,Co-design,of,PIC,methods,and,other,elements,of,NEPTUNE.,Ben,F,McMillan,April,21,,2021,This,report,is,the,deliverable,associated,with,Milestone,6.,The,NEPTUNE,project,scope,includes,simulating,the,dynamics,of,the,plasma,and,neutrals,in,the,edge,region,of,the,plasma.,Depending,on,the,ﬁdelity,re-,quired,,and,the,situations,of,interest,,this,requires,solving,both,for,ﬂuid-like,motion,in,collisional,regions,,as,well,as,essentially,free-streaming,kinetic,be-,haviour,for,less,collisional,plasma,constituents,and,neutrals.,Various,computational,methods,to,solve,these,kinetic,problems,exist,,and,we,will,here,be,discussing,particle-based,Monte-Carlo,methods,,where,the,com-,putational,markers,are,used,to,follow,the,trajectory,of,a,small,subsample,of,the,physical,particles.,These,are,in,more,general,settings,referred,to,as,La-,grangian,methods,and,distinct,from,Eulerian,methods,,in,which,distributions,are,evolved,on,a,ﬁxed,grid.,Particle-In-Cell,methods,,and,Monte-Carlo,meth-,ods,designed,to,evolve,the,neutral,population,,are,quite,similar,in,concept,,and,we,will,refer,to,them,as,PIC,methods,for,the,purposes,of,this,document.,We,restrict,the,discussion,to,using,particle,methods,to,solve,kinetic,(rather,than,ﬂuid),problems.,Integrating,these,methods,with,other,components,of,the,NEPTUNE,project,may,require,careful,planning.,Interaction,of,the,proposed,particle,method,with,two,other,work,packages,was,highlighted,in,the,bidding,process,as,needing,explicit,consideration,,and,a,discussion,of,these,interactions,as,well,as,a,plan,for,going,forward,the,principal,purpose,of,this,document.,Firstly,,because,the,gyrokinetics,work,package,will,determine,the,appropri-,ate,Fokker-Planck-Maxwell,equations,that,well-magnetised,particles,follow,,it,is,critical,to,determine,whether,the,particle,methods,can,be,made,to,support,these,equations.,The,gyrokinetics,work,package,is,also,investigating,a,moment-based,numerical,scheme,to,discretise,these,equations,in,a,way,suitable,for,implicit,solution;,we,will,explain,the,relationship,between,this,scheme,and,the,moment-,based,particle,discretisation,that,is,designed,to,reduce,noise,,but,also,can,be,leveraged,as,part,of,an,implicit,method.,Secondly,,the,NEPTUNE,project,has,selected,a,spatial,representation,(for,various,plasma,quantities,and,electromagnetic,ﬁelds),based,on,ﬁnite,elements,,on,an,unstructured,grid,,with,a,view,to,using,curved,element,boundaries,in,1,order,to,conform,to,the,magnetic,ﬁeld,topology,and,the,physical,ﬁrst,wall,of,the,tokamak.,This,is,intended,to,be,provided,through,the,Nektar++,library.,The,advanced,particles,methods,implemented,in,this,exploratory,stage,are,however,based,around,simple,uniformly,spaced,regular,Cartesian,meshes.,We,outline,design,considerations,of,a,particle,solver,interfaced,with,Spectral/hp,curvilinear,ﬁnite,element,methods.,We,also,explain,what,features,of,the,representation,of,the,electromagnetic,ﬁelds,are,desirable,for,accurate,particle,tracing.,1,What,are,PIC,methods,useful,for?,The,relative,computational,advantages,of,various,numerical,schemes,for,solving,kinetic,problems,will,determine,where,each,one,is,used,,and,for,which,class,of,particles,,or,where,kinetic,solution,is,feasible,at,all.,Particle,based,methods,are,intuitively,attractive,for,solving,kinetic,problems,,but,Eulerian,methods,converge,faster,as,a,result,of,the,Monte-Carlo,sampling,used,in,PIC,methods.,PIC,methods,are,often,favored,over,Eulerian,methods,when,the,following,conditions,apply,•,Three,velocity,space,directions,need,to,be,resolved.,•,Collision,operators,are,easily,written,in,the,form,of,a,stochastic,diﬀerential,equation,,but,not,as,a,PDE.,•,It,is,diﬃcult,to,ﬁnd,an,appropriate,velocity-space,grid,to,represent,the,particle,distribution,function,on,(eg.,multiple-beam,distributions).,•,It,is,diﬃcult,to,ﬁnd,an,optimal,spatial,grid,(as,PIC,codes,are,less,severely,impacted,by,running,on,a,non-optimal,spatial,grid).,In,general,PIC,methods,tend,to,be,relatively,easy,to,implement,and,paral-,lelise,compared,to,Eulerian,methods.,This,,and,a,willingness,to,tolerate,some,statistical,noise,,mean,that,early,examples,of,complex,plasma,kinetic,codes,tend,to,be,PIC-based,(e.g.,witness,the,evolution,of,core,and,edge,gyrokinetic,codes).,PIC,codes,are,dominant,for,neutral,particle,pushing,and,for,6D,Vlasov-Maxwell,(i.e.,non,gyrokinetic,PIC),because,except,in,special,cases,6D,Eulerian,codes,are,computationally,uncompetitive.,For,example,,if,full,fast,particle,orbits,need,to,be,tracked,to,determine,when,and,where,their,large,orbits,hit,the,wall,,PIC,is,a,natural,choice,for,this,problem.,Some,general,discussion,on,particle,methods,is,given,in,ref.,[1].,In,edge,plasmas,,ﬂuid,modelling,of,at,least,some,of,the,species,and,some,of,the,regions,is,essential.,As,a,result,,Eulerian,and,Lagrangian,kinetic,models,need,to,be,coupled,to,ﬂuid,models.,To,provide,maximum,predictive,and,numerical,performance,,it,is,highly,desirable,for,the,Neptune,project,to,be,able,to,choose,where,and,when,to,use,each,kind,of,model,in,a,ﬂexible,and,correctly,coupled,way.,The,remainder,of,the,document,described,interfaces,between,the,numerical,discretisations,of,these,models.,2,2,Edge-relevant,Gyrokinetic,formalism,We,brieﬂy,introduce,the,baseline,mathematical,formulation,(independent,of,discretisation),of,the,Vlasov,equations,for,the,initial,stage,of,the,NEPTUNE,project.,3,Reformulation,of,kinetic,equations,in,terms,of,closure,to,ﬂuid,equations,The,Fokker-Planck,equation,for,the,evolution,of,a,species,distribution,function,f,(Z),,with,Z,the,phase,space,coordinate,(representing,position,and,velocity,coordinates),may,be,written,in,the,form,(cid:104),˙Zf,=,S(f,),(1),(cid:105),+,∇Z.,∂f,∂t,with,S(f,),a,generalised,source,term,(including,self,and,non-self,collisions).,For,phase-space,conserving,equations,of,motion,(∇.,˙Z,=,0),,this,may,be,written,in,the,conservative,form,∂f,∂t,+,˙Z.∇Zf,=,df,dt,=,S(f,).,(2),As,in,the,Eulerian,representation,,it,is,often,useful,to,separately,represent,the,evolution,of,a,background,distribution,parameterised,using,ﬂuid,moments,,and,the,remaining,kinetic,response.,That,is,,we,introduce,a,splitting,f,(Z,,t),=,f0(Z,,t),+,g(Z,,t).,In,standard,(non-transformed),coordinates,the,phase,space,vector,may,be,decomposed,into,velocity,and,position,as,Z,=,(R,,v).,We,then,choose,the,parameterisation,f0(Z,,t),=,n(R,,t),[2πkT,(R,,t)]3/2,exp,(cid:20),m{v,−,v0(R,,t)}2/2,kT,(R,,t),(cid:21),(3),which,is,,again,,a,shifted,Maxwellian,,where,the,time,and,space,varying,pa-,rameters,v0,,T,and,n,,which,are,the,mean,velocity,,temperature,and,density,associated,with,the,distribution,f0,,are,as,yet,not,speciﬁed.,We,require,that,the,moments,1,,v,and,v2,of,g,vanish;,this,is,both,a,condition,on,the,initial,condition,,as,well,as,the,evolution,of,g.,By,taking,moments,of,the,Fokker-Planck,equation,,one,then,obtains,three,ﬂuid,equations,for,the,time,evolution,of,f0.,For,the,moment,,we,restrict,the,analysis,to,considering,a,Lorentz-like,force,F,(with,the,property,that,energy,transfer,to,the,ﬂuid,is,only,through,acceleration).,In,this,case,the,moment,equations,may,be,written,in,a,conservative,form,(c.f.,eqs.,57-59,of,Brunner,et.,al.,,which,are,essentially,identical,but,for,speciﬁc,forms,of,force/collisions,,and,with,momentum,equation,in,convective,form),as,+,∇.[nv0],=,(cid:90),d3v,S,,∂n,∂t,3,(4),∂mnv0,∂t,+,∇.[mnv0v0,+,P1],+,∇[nkT,],−,nF0,=,(cid:90),d3v,mvS,,(5),∂[3nkT,+,mnv2,0],∂t,+,∇.,(cid:20),v0{5nkT,+,mnv2,0},+,2P1.v0,+,(cid:90),d3v,m(v,−,v0)(v,−,v0)2g,(cid:21),−2nF0.v0,=,(cid:90),d3v,mv2S.,(6),where,F0,is,the,mean,per-particle,force,,q[E,+,v0,×,B],for,a,Lorentz,force,and,(cid:90),P1,=,d3v(v,−,v0)g(v,−,v0),(7),is,the,anisotropic,part,of,the,pressure,tensor,(i.e.,the,Stress,tensor).,The,ﬂuid,viscosity,and,heat,conduction,which,would,be,found,through,an,ordering,argu-,ment,in,the,collisional,limit,are,here,expressed,in,terms,of,the,kinetic,correction,g.,In,the,last,two,moment,equations,,the,quantities,in,the,time,derivative,are,momentum,and,kinetic,energy,,but,these,may,be,expanded,and,combined,to,ﬁnd,equations,in,the,‘convective,form’,where,the,only,time-derivatives,that,appear,are,the,velocity,and,temperature,instead.,Note,that,the,kinetic,correction,(closure),term,in,the,energy,equation,would,normally,be,written,in,terms,of,a,collisional,closure,in,a,ﬂuid,scheme,,as,diﬀusive,heat,ﬂuxes,,so,switching,between,collisional,and,fully,kinetic,solutions,to,these,moment,equations,can,be,achieved,in,a,natural,way.,The,Maxwell,equation,source,term,is,a,function,of,current,,which,may,be,written,in,terms,of,the,ﬂuid,quantities,as,qnv0,,thus,the,ﬂuctuation,g,only,enters,as,a,closure,term,in,the,energy,equation.,In,the,highly,collisional,limit,,where,the,heat,ﬂuxes,are,negligible,,we,recover,a,conventional,(multi-species),set,of,plasma,ﬂuid,equations,simply,by,setting,g,→,0.,For,time,evolution,of,g,,we,require,for,consistency,dg,dt,=,−,df0,dt,+,S(g,+,f0).,(8),We,have,considered,the,Vlasov-Maxwell,system,here,,but,conceptually,gyro-,Vlasov-Maxwell,follows,the,same,derivation,path,,and,although,the,ﬂuid,equa-,tions,diﬀer,somewhat,,the,equation,for,g,in,the,symbolic,form,above,is,equivalent,for,Vlasov,and,gyroVlasov,equations.,In,1D,,with,the,simulation,aligned,along,a,constant,background,ﬁeld,,which,is,the,initial,conﬁguration,to,be,considered,,these,systems,are,identical.,Note,that,the,above,equation,is,also,valid,for,neutral,particles.,It,is,beyond,this,point,that,the,proposed,methods,in,NEPTUNE,diﬀer.,Actually,,the,ﬂuid,equations,,and,the,coupling,to,the,gyro-Maxwell,or,gyro-,Poisson,equations,is,identical,,so,the,same,ﬂuid,and,ﬁeld,solvers,may,be,used.,However,,the,representation,of,g,diﬀers.,Firstly,,and,obviously,,in,the,Eulerian,approach,,g,is,represented,using,coeﬃcients,on,a,ﬁxed,spatial,grid,,whereas,in,the,particle-based,approach,,g,is,represented,using,markers.,4,Secondly,,the,Eulerian,method,described,in,Ref.,[2],performs,a,transfor-,mation,on,velocity,space,,shifting,and,rescaling,the,local,coordinate,V,=,M,(v,,v0,,T,),such,than,F0(V,),=,f0(M,−1V,)/n,is,an,unshifted,Maxwellian,of,temperature,and,density,unity.,The,grid-based,or,spectral,Eulerian,scheme,is,then,written,in,terms,of,V;,this,allows,the,grid,to,adapt,to,regions,of,low,tem-,perature,by,reﬁning,the,grid,resolution,so,the,spacing,is,constant,in,terms,of,local,thermal,velocity.,Although,it,would,also,be,possible,to,perform,this,transformation,in,the,Particle-based,scheme,,it,is,not,helpful,,because,the,eﬀective,velocity,resolution,depends,on,where,the,markers,are,in,velocity,space,,rather,than,their,coordi-,nate,labelling.,Particle,schemes,adapt,to,velocity,space,automatically,because,the,marker,phase-space,density,remains,(in,a,statistical,sense),proportional,to,the,particle,phase-space,density,,so,markers,in,cold,regions,naturally,have,low,typical,velocities.,4,Particle-speciﬁc,implementation,of,the,moment-,based,scheme,We,brieﬂy,summarise,the,method,to,be,implemented;,for,Vlasov-Maxwell,sys-,tems,,this,was,implemented,in,Ref.,[3].,We,consider,N,markers,(computational,macroparticles),loaded,in,phase,space,with,positions,Zi(t),with,a,phase,space,density,K(Z),=,lim,N→∞,(cid:90),(cid:88),dV,δ6(Z,−,Zi(t)),N,(9),that,represents,how,likely,we,are,to,ﬁnd,a,marker,in,a,small,region,of,phase,space,near,Z.,Due,to,Liouville’s,theorm,,we,ﬁnd,that,dK/dt,=,0.,We,store,the,initial,value,Ki,=,K[Zi(0)],at,each,marker,position,and,identify,the,eﬀective,volume,of,phase,space,associated,with,each,marker,as,1/K(Z).,To,represent,the,ﬂuctuating,distribution,function,,we,deﬁne,a,weight,gi,at,each,particle,location,g(Z,,t),=,(cid:88),N,1,Ki,gi(t)δ6(Z,−,Zi(t)),with,gi(0),=,g(Zi(0),,0),and,and,if,we,evolve,gi(t),via,dZi(t),dt,=,˙Z(Zi,,t),dgi,dt,=,−,df0,dt,(cid:12),(cid:12),(cid:12),(cid:12)Z=Zi,+,S(Zi),5,(10),(11),(12),then,by,integrating,over,any,volume,K,,we,can,show,that,in,the,large-marker,limit,,this,will,satisfy,dg,dt,+,df0,dt,−,S(g,+,f0),≡,Q,=,0.,(13),K,d3vQ,=,0.,That,is,,unlike,standard,PIC,,in,the,weak,form,,that,is,limN→∞,which,trivially,satisﬁes,df,/dt,=,0,exactly,for,the,numerical,distribution,and,numerical,orbits,,the,numerical,Vlasov,equation,is,satisﬁed,only,in,the,weak,form,and,large,particle,number,limit,when,using,control-variates[4].,(cid:82),Note,that,this,the,standard,control,variates[5,,6],method.,Confusingly,,this,is,often,referred,to,as,the,‘δf,’,method,,but,‘δf,’,is,also,used,to,denote,methodologies,that,are,valid,only,when,the,ﬂuctuation,is,small.,There,is,no,such,approximation,used,here;,there,is,no,formal,restriction,on,the,size,of,g,,although,when,f,and,g,are,the,same,size,,the,numerical,advantages,may,be,lost.,The,principal,diﬀerence,between,this,control-variates,method,and,that,im-,plemented,in,codes,such,as,ORB5,and,EPOCH,is,time-evolution,of,the,back-,ground,distribution,function,f0,(however,,a,project,is,underway,to,implement,background,density,variation,in,ORB5).,In,a,Galerkin,method,,moments,of,the,particle,distribution,appear,multiplied,by,basis,functions.,For,example,,consider,a,representation,of,density,n(x),=,(cid:88),k,nkΛk(x),(14),where,nk,are,the,nodal,density,values,and,Λ,the,basis,function.,In,the,weak,form,,this,may,be,equated,to,the,velocity-space,integral,of,the,particle,represen-,tation,of,the,distribution,function,f,.,Although,this,may,straightforwardly,be,extended,to,control-variates,,to,simplify,the,presentation,,we,will,use,a,standard-,PIC,representation,,without,background,extraction,,and,represent,the,full,dis-,tribution,on,markers,using,fi,(for,f0,=,0,we,have,fi,=,gi).,We,then,have,(cid:90),(cid:90),(cid:90),(cid:90),(cid:90),dxΦ(x),dvf,=,dxΦ(x),dv,δ6(Z,−,Zi),=,dxΦ(x),nkΛk(x),(cid:88),i,fi,Ki,(cid:88),k,(15),for,all,trial,functions,Φ(x),in,the,ﬁnite,element,representation:,since,the,basis,functions,span,this,space,,it,is,suﬃcient,if,this,equality,is,satisﬁed,for,Φ(x),=,Λh∀h.,Because,,if,we,choose,trial,functions,in,the,same,space,as,the,density,representation,,we,have,Φ(x),=,(cid:80),(cid:90),h,nhΛh(x),,this,gives,us,(cid:90),∀h,,dxλh,δ(x,−,xi),=,dxΛh(x),nkΛk(x),(16),(cid:88),i,fi,Ki,(cid:88),k,or,(cid:88),∀h,,λh(xi),i,fi,Ki,(cid:88),=,Mhknk,k,(17),where,we,have,deﬁned,a,mass,matric,Mhk,representing,the,integral,on,the,RHS.,We,may,then,ﬁnd,a,FEM,representation,of,the,density,by,performing,the,sum,on,the,LHS,and,inverting,the,mass,matrix.,6,In,control-variates,form,,if,the,background,density,may,be,represented,spa-,tially,using,the,FEM,basis,functions,,this,equation,(eq,17),holds,with,fi,replaced,by,gi,and,nk,replaced,with,the,density,ﬂuctuation.,5,Computational,advantages,of,moment-based,schemes.,This,moment-based,scheme,(whether,in,Eulerian,or,Lagrangian,form),allows,a,close,connection,to,the,ﬂuid,equations,,and,allows,the,ﬂuid,and,ﬁeld,equations,to,be,partially,decoupled,from,the,kinetic,equations.,One,way,this,can,be,exploited,is,that,in,regions,where,kinetic,eﬀects,are,minimal,,the,kinetic,equations,no,longer,need,to,be,solved,at,all.,An,additional,advantage,is,that,the,coupled,Vlasov-Maxwell,system,can,be,quite,stiﬀ,in,the,sense,that,there,are,oscillations,and,relaxations,on,much,faster,timescales,than,the,dynamical,timescales,of,interest.,For,example,,Alfven,waves,propagate,at,close,to,the,speed,of,light,in,certain,edge,regions,,whereas,even,the,highest,energy,particles,typically,propagate,at,one,percent,of,this,speed.,In,the,usual,cases,,these,rapid,oscillations,are,associated,with,ﬂuid,processes,(i.e.,waves),and,the,kinetic,eﬀects,,such,as,nonlocal,thermal,conduction,,appear,on,longer,timescales.,It,is,computationally,advantageous,to,evolve,or,relax,these,waves,by,solving,low-dimensional,(3D),implicit,or,explicit,ﬂuid,equations,rather,than,(5-6D),kinetic,equations.,Even,if,the,ﬂuid,motion,needs,to,be,solved,explicitly,,this,simpliﬁcation,can,allow,a,massive,computational,speedup.,For,the,particle,method,,the,moment-based,scheme,also,allows,for,noise,reduction,by,extracting,out,the,background,distribution,and,integrating,it,an-,alytically.,That,is,,consider,the,moment,(cid:90),n,=,dV,f,=,(cid:90),dV,f0,+,(cid:90),dV,g.,(18),Because,f0,may,be,evaluated,analytically,,only,the,integral,over,g,needs,to,be,evaluated,using,a,Monte-Carlo,method;,and,the,RMS,error,of,the,Monte-,Carlo,integration,scales,like,the,amplitude,of,g,,so,if,g,is,small,,then,this,can,lead,to,a,radical,decrease,in,the,number,of,markers,needed,to,attain,a,speciﬁc,accuracy.,For,core,plasma,simulations,<,g2,>,/,<,f,2,>∼,10−4,this,allows,10,000,times,fewer,markers,to,be,used.,In,the,edge,,ﬂuctuations,are,large,,but,relatively,high,collisionality,means,that,the,departures,from,Maxwellian,are,not,necessarily,large,(e.g.,low,density,high,velocity,tails,are,present).,Maintaining,noise,reduction,over,long,simulation,timespans,requires,additional,machinery,to,maintain,good,sampling,of,phase,space,and,avoid,excessive,weight,spreading.,6,Proof,of,principle,implementation,A,1D,,two-moment,model,is,implemented,in,the,minepoch,code,that,is,on,the,github,repository.,7,In,the,test,case,setup,(an,input,ﬁle,‘singlestream.deck’,is,given,for,this,case),is,that,a,low,temperature,single-species,collisionless,neutral,species,is,given,an,initial,density,and,velocity,perturbation,of,the,form,n,=,n0[1,+,0.3,sin(4πx/Lx)],,v,=,v0[1,+,sin(6πx/Lx)].,(19),The,details,of,parameters,for,the,test,case,setup,are,described,in,detail,in,the,testing,code.,The,purpose,of,this,test,case,is,a,proof,of,principle,demonstration,that,the,control,variates,method,is,able,to,reduce,the,statistical,sampling,error,in,a,case,with,order,one,density,ﬂuctuations,(as,one,might,ﬁnd,in,the,edge):,this,is,an,advancement,over,control,variates,schemes,used,in,the,core,of,tokamaks,,which,cannot,handle,large,variations.,The,low,temperature,means,the,energy,equation,is,not,needed,,and,as,the,evolution,is,force,free,,a,simpliﬁed,momentum,evolution,equation,is,suﬃcient.,For,this,setup,,because,the,evolution,is,force-free,,the,method,of,character-,istics,may,be,used,to,ﬁnd,the,solution,analytically,(in,an,implicit,form).,We,show,the,time-evolution,of,the,density,in,the,simulation,in,ﬁgures,1,,2,and,3.,Figure,1:,Density,versus,position,at,an,intermediate,simulation,time.,At,t,=,2.6,×,10−9,s,,with,f,=,f0,+,g,in,the,control-variates,simulation,,we,have,<,g2,>,/,<,f,2,>=,1.1,×,10−3.,That,is,,despite,the,variation,of,density,by,a,factor,of,2,in,the,x,domain,,the,kinetic,correction,to,the,particle,distribution,g,is,∼,30,times,smaller,than,the,overall,distribution,function,f,.,The,squared,noise,amplitude,scales,like,<,g2,>,/N,,,so,either,one,may,use,1000,times,fewer,markers,in,the,control-variates,simulation,for,the,same,noise,,or,1000,times,lower,squared,noise,amplitude,with,the,same,number,of,markers.,8,Figure,2:,Colour,plot,of,density,versus,position,and,time,in,a,standard,(no,control,Variates),PIC,simulation.,7,Finite,elements,and,particle,methods,,gyroki-,netic,examples.,Finite,element,method,meshes,are,fairly,regularly,used,as,a,basis,for,particle-,in-cell,codes.,Two,examples,of,gyrokinetic,particle-in-cell,codes,are,the,XGC[7],code,and,ORB5[8];,the,author,of,this,document,is,a,developer,of,the,ORB5,code.,The,ORB5,code,used,a,logically,Cartesian,tensor-product,B-spline,represen-,tation,in,magnetic,coordinates.,The,elements,are,curved,in,laboratory,Cartesian,space,,and,this,allows,a,simple,structured,grid,that,nonetheless,conforms,to,the,curved,magnetic,ﬁeld,topology,of,the,tokamak,core.,This,approach,is,not,able,to,handle,the,X-point,of,scrape-oﬀ,layer,region,directly.,We,have,magnetic,coordinates,(s,,θ,,ζ),,which,are,toroidal,coordinates,with,s,∈,[0,,1],representing,the,radial,distance,from,the,magnetic,axis,,χ,the,straight-,ﬁeld,line,poloidal,angle,,and,ζ,the,(geometrical),toroidal,angle.,Taking,the,equal,spaced,case,for,simplicity,,ﬁelds,are,discretely,represented,as:,φ(s,,θ,,ζ),(cid:88),(cid:88),(cid:88),i=[0,Ni],j=[0,Nj,],k=[0,Nk],φi,j,kΛ,(cid:21),(cid:20),s,−,si,δs,Λ,(cid:21),(cid:20),θ,−,θj,δθ,Λ,(cid:21),(cid:20),ζ,−,ζk,δζ,(20),Here,,Λ(x),are,(quadratic,or,cubic),B-spline,functions,(Fig.,4),,which,serve,as,compact-support,basis,functions,for,the,ﬁnite-element,scheme.,For,quadratic,basis,functions,,the,representation,has,continuous,values,and,ﬁrst,derivatives,(i.e.,it,is,C,1),and,for,cubic,,the,second,derivatives,are,also,continuous,(it,is,C,2).,XGC,uses,a,block-structured,regular,mesh,,also,based,on,low-order,ﬁnite,elements.,The,blockwise,decomposition,allows,the,code,to,represent,,using,one,9,Figure,3:,Colour,plot,of,density,versus,position,and,time,in,a,control-Variates,PIC,simulation.,block,the,core,plasma,using,elements,that,conform,to,magnetic,ﬁeld,lines,to,capture,the,strong,anisotropy.,The,block,region,around,the,X-point,is,not,required,to,conform,to,the,magnetic,geometry,,but,the,low,poloidal,ﬁeld,in,this,region,simpliﬁes,the,handling,of,anisotropy.,Additional,blocks,in,the,scape,oﬀ,layer,allow,the,full,geometry,to,be,handled.,Particle-ﬁeld,operations,are,seen,in,a,radically,diﬀerent,way,in,these,ﬁnite-,element,based,PIC,codes,to,traditional,PIC,codes,,which,use,ﬁnite,diﬀerence,or,volume,concepts.,Traditional,PIC,(e.g.,EPOCH),considers,the,markers,to,have,a,smooth,‘cloud,of,charge’,associated,with,all,the,physical,particles,the,marker,is,representing,,leading,to,a,’shape,function’,which,may,be,used,to,deﬁne,a,continuous,particle,density,throughout,the,simulation,domain.,This,smooth,ﬁeld,may,then,be,evaluated,at,grid,points,in,a,ﬁnite,diﬀerence,scheme,,or,integrated,over,grid,volumes,in,a,ﬁnite,volume,scheme,(e.g.,EPOCH).,Finite-,element,based,PIC,schemes,,on,the,other,hand,,consider,the,markers,as,point-,particles,without,an,explicit,spatial,extent.,However,,functionally,the,smooth,ﬁeld,representation,in,a,ﬁnite,element,code,plays,a,very,similar,role,,and,the,eﬀective,charge,density,as,represented,on,the,FEM,grid,is,smooth.,All,these,codes,involve,direct,evaluation,of,the,ﬁelds,at,the,particle,posi-,tions,,which,is,logically,straightforward.,First,,one,determines,which,grid,cell,in,ORB5,the,particles,coordinates,are,evolved,in,magnetic,the,particle,is,in:,coordinates,,and,these,coordinates,can,be,individually,mapped,to,grid,indexes,in,each,direction.,The,compact,support,then,allows,an,evaluation,using,a,rela-,tively,small,number,of,basis,functions.,For,tensor-product,based,schemes,,about,two,multiplications,are,needed,per,contributing,basis,function,,so,even,though,in,the,cubic,case,64,basis,functions,contribute,in,3D,(cubic,spline,elements,span,4,cells,in,each,direction,,see,ﬁg,4),,the,per-evaluation,cost,is,not,extreme.,For,a,fully,unstructured,,curvilinear,grid,,one,can,in,principle,evaluate,the,10,Figure,4:,B-Splines,on,a,unit-interval,grid.,ﬁeld,at,arbitrary,locations,,but,this,may,not,be,straightforward,or,eﬃcient.,It,may,in,general,be,complicated,or,ineﬃcient,to,determine,which,grid,cell,a,particle,is,in,based,on,its,global,coordinate,position.,The,mapping,used,to,allow,distorted,element,positions,is,the,forward,mapping,from,element-wise,local,coordinates,to,global,Cartesian,coordinates,and,the,inverse,mapping,may,be,diﬃcult,to,compute,,so,even,if,the,grid,cell,can,be,determined,,the,local,coordinate,position,in,the,cell,may,be,annoying,to,compute.,In,this,case,,it,may,be,desirable,to,keep,track,of,the,cell-local,coordinates,of,the,particle,,and,,when,it,exits,the,cell,,to,update,which,cell,it,is,in,based,on,its,crossing,of,boundaries.,Using,a,cell-local,coordinate,scheme,complicates,somewhat,the,particle,time-,stepping,in,the,sense,that:,•,The,numerical,equations,may,have,to,be,rewritten,in,a,general,curvilinear,form.,•,The,coordinates,are,continuous,between,cells,,but,not,necessarily,smooth,,which,will,lead,to,low,order,of,time-convergence,of,timestepping,schemes,without,any,additional,measures,being,taken.,•,One,must,keep,track,of,which,cell,the,particle,is,in,and,hand,it,between,cells,during,the,timestep,,possibly,across,multiple,cell,boundaries.,Overall,,the,use,of,cell,coordinates,rather,than,global,coordinates,creates,a,much,tighter,coupling,and,more,complicated,interface,between,the,particle,and,FEM,code,elements.,But,tight,coupling,may,be,required,in,any,case,to,handle,practical,numerical,issues,such,as,parallelisation.,11,8,Timestepping,accuracy,and,continuity,A,particle,timestepping,equation,symbolically,of,the,form,dZ,dt,=,F,(Z(t),,t),(21),usually,involves,a,dependence,on,certain,ﬁelds,(density,ﬁelds,leading,to,particle,drag,,electromagnetic,ﬁelds,for,charged,particles),that,are,represented,,in,the,case,of,Nektar++,,using,the,FEM.,Consider,a,trajectory,Z(cid:48)(t),that,is,locally,analytic,(i.e.,smooth,to,arbitrary,order),and,the,reduced,equation,dZ,dt,=,F,(Z(cid:48)(t),,t),(22),If,F,(Z,,t),is,a,piecewise,polynomial,function,in,Z,,then,as,Z(cid:48)(t),crosses,element,boundaries,,F,will,be,non-smooth,in,time.,If,derivatives,up,to,the,Nth,are,continuous,(a,C,N,representation),then,a,simple,time-stepping,scheme,will,have,a,maximum,local,convergence,rate,of,order,N,+,2,on,timesteps,where,this,boundary,is,crossed,,and,N,+,2,globally,,since,most,intervals,will,not,cross,a,cell,boundary.,This,means,,for,example,,that,a,standard,fourth-order,Runge-Kutta,method,,with,4th,order,global,accuracy,,will,only,achieve,that,order,of,convergence,trac-,ing,particles,in,ﬁelds,with,at,least,C,2,smoothness.,Although,order,of,con-,vergence,is,not,a,complete,diagnostic,of,the,usefulness,of,a,scheme,,this,indi-,cates,that,high-order,timestepping,may,not,be,appropriate.,Whether,high-order,timestepping,accuracy,is,required,is,quite,problem,dependent.,Special,purpose,integrators,may,be,built,to,break,the,time,domain,up,into,regions,so,that,in,each,time,sub-interval,the,particle,does,not,cross,any,cell,boundaries,,but,this,clearly,creates,signiﬁcant,extra,complexity.,Particle,tracing,in,magnetised,plasmas,is,a,somewhat,special-purpose,prob-,lem,,where,the,evolution,equations,may,often,be,split,into,two,portions,F0,+,F1,,with,F0,representing,(normally),the,long-wavelength,background,magnetic,ﬁeld,,and,F1,associated,with,ﬁne-scale,low-amplitude,turbulent,ﬂuctuation.,Because,of,the,extreme,anisotropy,of,transport,,it,is,often,important,to,very,accurately,solve,for,the,particle,motion,along,the,background,ﬁeld,F0,,in,order,that,small,errors,in,tracing,unperturbed,particle,orbits,do,not,overwhelm,the,small,trans-,port,eﬀect,due,to,F1.,Note,that,in,a,tokamak,the,particle,drifts,(in,gyrokinetic,theory),depend,on,derivatives,of,the,magnetic,ﬁeld,(∇B);,this,means,that,even,deﬁning,these,drifts,requires,a,high,order,of,continuity,of,the,magnetic,ﬁeld,representation.,This,issue,will,presumably,also,appear,in,Eulerian,gyrokinetics.,Due,to,the,long,wavelength,smoothness,of,the,background,ﬁelds,,even,if,these,derivatives,formally,lack,continuity,,the,discontinuities,may,not,be,very,large.,The,standard,solution,in,tokamak,codes,is,to,represent,the,background,ﬁelds,in,a,way,that,allows,high-order,derivatives,to,be,taken,,perhaps,on,a,simple,Cartesian,grid,,where,high-order,continuous,representations,are,straightforward,,12,and,then,evaluate,them,in,a,setup,phase,on,the,grid,points,on,which,the,dy-,namics,will,be,solved,,which,may,be,a,semi-structured,or,unstructured,mesh.,Note,that,,depending,on,the,timestepping,scheme,chosen,,the,particle,evo-,lution,may,need,ﬁeld,evaluation,at,intermediate,time,points.,9,Particles,in,Nektar++,,existing,implementa-,tion,,and,proof,of,principle.,Particles,have,been,implemented,as,a,ﬁlter,within,the,Nektar++,framework,,and,exploited,for,the,purposes,of,considering,erosion,wear,on,surfaces[9].,This,framework,(currently),allows,for,certain,kinds,of,particle,motion,relevant,to,solving,the,Navier-Stokes,equations,,whose,motion,is,subject,to,local,ﬂuid,forces,,but,not,for,particles,to,impact,on,the,ﬂuid,motion.,Several,kinds,of,particle,motion,are,possible,,and,are,schematically,either,ﬂuid,particles,,where,the,particles,have,the,local,ﬂuid,velocity,v,and,satisfy,the,equation,dr,dt,=,v(r,,t),or,ﬁnite-sized,particles,,which,follow,equations,of,motion,dr,dt,=,V,,m,dV,dt,=,F,(V,,v,,t),(23),(24),of,particles,subject,to,forces,due,to,the,ﬂuid.,The,ﬁlter,implements,linear,multistep,methods,to,solve,these,equations;,these,are,quite,adaptable,to,plasma,problems,,but,adaption,to,use,other,kinds,of,schemes,would,not,be,excessively,complicated.,The,existence,of,the,particles,ﬁlter,in,nektar++,,as,well,as,a,Hasegawa-,Wakatani,solver,,has,allowed,us,to,set,up,a,basic,coupled,plasma-particles,test,case,,as,a,preliminary,to,exploring,the,impact,of,e.g.,curved,simulation,elements,on,the,particle,tracing,problem.,As,a,demonstration,,we,initialise,particles,on,a,uniform,grid,and,allow,them,to,be,advected,by,the,ﬂuid,motion,(ﬁg.,5).,10,Co-design,recommendations,•,The,similarity,of,the,Particle,and,Eulerian,moment-based-schemes,mean,that,a,common,codebase,should,be,possible,and,this,will,avoid,duplication,of,eﬀort.,Until,that,point,,only,very,simple,ﬂuid,solvers,are,required.,•,Particle,methods,for,gyrokinetic,systems,are,generally,easy,to,adaptable,in,terms,of,trajectory,equations,of,the,chosen,gyrokinetic,equations,,so,it,is,not,crucial,to,ﬁx,these,at,an,early,stage,of,the,project.,•,The,representation,of,background,and,perturbed,electromagnetic,ﬁeld,needs,special,care,if,particle,orbits,are,to,be,accurately,solved.,Some,13,Figure,5:,Particle,trajectories,on,the,(x,,y),plane.,Red,markers,show,initial,positions,,black,shows,positions,at,t,=,4.,precalculation,may,be,necessary,for,derivatives,of,magnetic,ﬁelds.,Espe-,cially,if,perturbed,magnetic,ﬁelds,may,be,large,careful,thought,needs,to,be,given,to,creating,a,smooth,representation,of,ﬁelds,or,indirectly,evaluating,derivatives.,•,Thought,needs,to,be,given,to,co-designing,timestepping,schemes,for,com-,bined,ﬂuid-particle,codes;,particle,methods,can,use,e.g.,linear,multistep,methods,,or,special,purpose,integrators,(Boris,algorithm),but,the,choice,depends,somewhat,on,the,diﬃculty,of,evaluating,ﬁelds,at,intermediate,timepoints.,•,A,dedicated,study,would,be,needed,to,examine,particle-tracking,in,piece-,wise,curvilinear,grids,if,the,Neptune,software,is,to,use,particle,methods.,There,is,infrastructure,already,existing,to,perform,this,work,with,Nek-,tar++,,both,for,conventional,ﬂuid,problems,or,coupled,to,plasma-speciﬁc,problems.,References,[1],NEPTUNE,authors.,M2.3.1,options,for,particle,algorithms.,Technical,re-,port,,CCFE,,2020.,[2],Michael,Barnes,Felix,I.,Parra,and,Michael,Hardman.,2047357-TN-01-2,,1d,drift,kinetic,models,with,periodic,boundary,conditions.,Technical,report,,Oxford,University,,2021.,14,[3],S.,Brunner,,E.,Valeo,,and,J.,A.,Krommes.,Collisional,delta-f,scheme,with,evolving,background,for,transport,time,scale,simulations.,Physics,of,Plas-,mas,,6(12):4504–4521,,1999.,[4],B.,F.,McMillan,and,L.,Villard.,Accuracy,of,momentum,and,gyrodensity,transport,in,global,gyrokinetic,particle-in-cell,simulations.,Physics,of,Plas-,mas,,21(5):052501,,2014.,[5],A.Y.,Aydemir.,A,uniﬁed,Monte,Carlo,interpretation,of,particle,simulations,and,applications,to,non-neutral,plasmas.,Physics,of,Plasmas,,1:822–831,,1994.,[6],Roman,Hatzky,,Trach,Minh,Tran,,Axel,Knies,,Ralf,Kleiber,,and,Simon,J.,Allfrey.,Energy,conservation,in,a,nonlinear,gyrokinetic,particle-in-cell,code,for,ion-temperature-gradient-driven,modes,in,theta-pinch,geometry.,Physics,of,Plasmas,,9(3):898–912,,2002.,[7],S.,Ku,,C.,S.,Chang,,R.,Hager,,R.,M.,Churchill,,G.,R.,Tynan,,I.,Cziegler,,M.,Greenwald,,J.,Hughes,,S.,E.,Parker,,M.,F.,Adams,,E.,D’Azevedo,,and,P.,Worley.,A,fast,low-to-high,conﬁnement,mode,bifurcation,dynamics,in,the,boundary-plasma,gyrokinetic,code,XGC1.,Physics,of,Plasmas,,25(5):056107,,2018.,[8],E.,Lanti,,N.,Ohana,,N.,Tronko,,T.,Hayward-Schneider,,A.,Bottino,,B.F.,McMillan,,A.,Mishchenko,,A.,Scheinberg,,A.,Biancalani,,P.,Angelino,,S.,Brunner,,J.,Dominski,,P.,Donnel,,C.,Gheller,,R.,Hatzky,,A.,Jocksch,,S.,Jolliet,,Z.X.,Lu,,J.P.,Martin,Collar,,I.,Novikau,,E.,Sonnendr¨ucker,,T.,Ver-,nay,,and,L.,Villard.,ORB5:,A,global,electromagnetic,gyrokinetic,code,using,the,PIC,approach,in,toroidal,geometry.,Computer,Physics,Communications,,page,107072,,2019.,[9],Manuel,F.,Mej´ıa,,Douglas,Serson,,Rodrigo,C.,Moura,,Bruno,S.,Carmo,,Jorge,Escobar-Vargas,,and,Andr´es,Gonz´alez-Mancera.,Erosion,wear,evaluation,using,Nektar+,+.,In,Spencer,J.,Sherwin,,David,Moxey,,Joaquim,Peir´o,,Peter,E.,Vincent,,and,Christoph,Schwab,,editors,,Spectral,and,High,Order,Methods,for,Partial,Diﬀerential,Equations,ICOSAHOM,2018,,pages,419–,428,,Cham,,2020.,Springer,International,Publishing.,15

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