CD-EXCALIBUR-FMS0072-M4c.2_ExcaliburFusionModellingSystem ========================================================= .. meta:: :description: technical note :keywords: ExCALIBUR,2-D,integrated,particle,and,continuum,model,M4c.2,Version,1.00,Abstract,The,report,describes,work,for,ExCALIBUR,project,NEPTUNE,at,Mile-,stone,M4c.2.,An,integrated,particle,and,continuum,model,has,been,de-,veloped,in,2d3v,,where,2d3v,implies,that,there,is,variation,in,two,space,dimensions,(2-D),and,3,velocity-phase,space,coordinates.,NEKTAR++,pro-,vides,the,2-D,finite,elements,to,represent,the,plasma,as,a,fluid,,whereas,3-D,velocity-phase,space,is,populated,with,particles,to,represent,the,neu-,tral,particles,and,the,complete,calculation,is,orchestrated,by,the,NESO-,PARTICLES,software,developed,under,project,NEPTUNE.,Much,of,the,new,code,uses,SYCL,to,pave,the,way,to,Exascale,in,a,performance-portable,manner.,Initial,tests,using,a,small,range,of,boundary,conditions,and,sources,demonstrate,successful,strong,coupling,between,the,fluid,and,the,new,par-,ticle,code.,Additionally,,summaries,are,provided,of,work,by,a,grantee,to,develop,a,high-dimensional,finite,element,library,intended,to,provide,a,com-,plementary,approach,to,the,use,of,particles.,UKAEA,REFERENCE,AND,APPROVAL,SHEET,Client,Reference:,UKAEA,Reference:,CD/EXCALIBUR-FMS/0072,Issue:,Date:,1.00,March,21,,2023,Project,Name:,ExCALIBUR,Fusion,Modelling,System,Prepared,By:,Name,and,Department,James,Cook,Will,Saunders,Owen,Parry,Signature,N/A,N/A,N/A,Date,March,21,,2023,March,21,,2023,March,21,,2023,BD,Reviewed,By:,Wayne,Arter,March,21,,2023,Project,Technical,Lead,1,1,Introduction,The,reach,of,this,milestone,is,to,create,a,proxyapp,that,combines,at,least,one,complex,element,from,several,aspects,of,exhaust,physics.,This,2d3v,proxyapp,is,based,on,an,advection-diffusion-,reaction,(ADR),equation,system,where,the,source,terms,arise,from,the,ionisation,of,neutral,parti-,cles,which,are,treated,as,Lagrangian,markers.,The,system,of,three,coupled,fluid,equations,for,number,density,,velocity,and,energy,,encapsulated,in,a,state,vector,(cid:126)ς,and,expressed,in,a,form,familiar,to,the,area,of,CFD,and,specifically,Finite,Volume,formulations,as,∂(cid:126)ς,∂t,+,∇,·,(F((cid:126)ς)),=,G((cid:126)ς,,∇,·,(cid:126)ς),+,S((cid:126)ς),,∇,(1),where,F,represents,the,spatial,component,of,a,general,hyperbolic,conservation,law,,G,denotes,the,diffusion,terms,and,S((cid:126)ς),is,the,source,term,for,the,components,of,the,state,vector.,We,provide,the,particular,system,of,equations,of,interest,in,Section,2.1.1.,In,this,report,we,present,details,of,the,implementation,and,present,initial,results,of,an,exhaust,relevant,ADR,solver,constructed,from,NEKTAR++,and,NESO-PARTICLES(NP),and,realised,as,part,of,the,NESO,framework.,The,targeted,real-world,scenario,is,the,scrape-off-layer,(SOL),which,connects,the,tokamak,divertor,(heat,sink),to,the,bulk,plasma,at,the,outer,mid-plane,edge,(heat,source).,Furthermore,the,plasma,is,coupled,to,a,system,of,neutral,particles.,The,fluid,species,exists,as,a,finite,element,representation,provided,by,the,NEKTAR++,library.,The,neutral,particle,species,is,implemented,using,the,NP,framework.,Within,NESO,there,exists,a,tight,coupling,between,finite,element,functions,within,NEKTAR++,and,particles,in,NP.,The,computationally,challenging,dynamics,are,largely,aligned,to,the,magnetic,field,,which,pred-,icates,the,use,of,a,long,,thin,computational,domain,also,aligned,with,the,magnetic,field.,In,this,domain,,the,plasma,source,is,located,in,the,middle,and,the,upper,and,lower,divertors,are,located,at,the,ends.,This,proxyapp,implements,a,system,with,an,element,of,complexity,from,many,areas,of,the,parameter,space:,advection-diffusion,in,NEKTAR++,,performance,portable,particles,via,NP,,domain,decomposition,via,MPI,,mixed,boundary,conditions,,tight,coupling,of,fluids,and,particles,,IO,,visualisation,,and,crucially,2,spatial,dimensions,and,3,velocity,dimensions.,The,final,element,of,complexity,is,critical,for,providing,functionality,to,the,tokamak,edge,modelling,community,who,wish,to,calculate,the,heat,fluxes,required,for,so-called,mean,field,equations,in,2-D,with,neutral,physics,via,particles,with,all,three,velocity,components.,The,structure,of,the,rest,of,the,report,is,as,follows.,In,Section,2,we,describe,the,technical,work,done,for,this,report,and,describe,the,equations,that,are,solved,for,the,fluid,dynamics,and,particles,,how,they,are,coupled,and,the,how,the,ionisation,rates,are,obtained.,We,also,present,an,exposition,of,the,ionisation,algorithm.,We,provide,a,description,of,the,basis,functions,in,NEKTAR++,and,their,efficient,evaluation,in,NESO,via,SYCL.,In,Section,3,,we,summarise,and,give,context,to,two,associated,grantee,reports,on,this,topic.,The,final,section,provides,a,summary,of,this,work.,2,2,Task,Work,2.1,Plasma,,particle,and,ionisation,equations,The,plasma,exhaust,is,modelled,in,2-D,in,a,long,and,thin,domain,,where,the,centre,represents,the,outer,mid-plane,edge,of,the,plasma,,ie.,broadly,where,the,hot,plasma,enters,the,scrape,off,layer,domain,,and,the,ends,represent,the,divertor,,ie.,the,plasma,from,the,centre,travels,along,the,domain,in,both,directions,towards,the,exhaust,region.,Critical,to,fusion,power,plant,design,is,the,concept,of,a,detached,plasma,whereby,a,region,of,ionisation,and,recombination,exists,between,the,plasma,and,the,divertor.,This,detached,state,is,beneficial,as,the,divertor,is,in,contact,with,merely,a,hot,gas,rather,than,a,plasma,carrying,orders,of,magnitudes,larger,heat,fluxes.,The,region,of,ionisation,and,heat,flux,near,the,ends,of,this,long,and,thin,2-D,domain,is,represented,in,this,mini-app,as,locations,of,neutral,particle,sources,(whose,presence,in,a,real,tokamak,would,arise,via,gas,puffs,,recycling,from,the,wall,and,the,processes,of,ionisation,and,recombination,between,the,neutral,gas,and,the,detached,plasma).,In,this,instance,,we,focus,on,ionisation,only,and,leave,recombination,as,an,item,of,further,work.,The,following,sections,describe,the,fluid,and,particle,dynamical,equations,and,the,origin,of,the,ionisation,rates.,2.1.1,Dynamics,equations,The,plasma,exhaust,relevant,fluid,equations,,in,the,form,of,an,advection-diffusion-reaction,(ADR),equation,,are,expressed,in,2-D,by,Eulerian,methods,,where,the,source,terms,in,the,fluid,equa-,tions,arise,from,the,ionisation,of,neutral,particles.,We,treat,these,neutral,particles,as,Lagrangian,markers.,We,describe,the,time,evolution,of,number,density,n,,velocity,(cid:126)u,and,double,energy,E,=,(g,2)nT,+,nu2,as,−,Ur,n,=,∂,∂t,∂,∂t,2)T,+,u2)),=,n(cid:126)u,=,Ur,(n(cid:126)u),+,Sn,,−∇,·,−∇,·,(n(cid:126)u,⊗,(cid:126)u,+,nT,(cid:126)I),+,Su,,(n(cid:126)u(gT,+,u2)),κd,∇,−,2T,+,SE,,−∇,·,(2),(3),(4),Ur,∂,∂t,(n((g,−,where,T,and,g,are,the,plasma,temperature,and,equation,of,state,constant,respectively,and,Si,is,the,source,term,for,variable,i,arising,from,the,ionisation,of,neutral,particles.,The,neutral,com-,putational,macroparticles,follow,straight,trajectories,from,birth,until,either,loss,due,to,complete,ionisation,or,because,they,have,travelled,past,the,ends,of,the,domain,in,the,long,direction,ie.,hit,the,divertor,target.,We,sample,initial,neutral,particle,properties,from,a,fixed,density,and,temperature,profile.,We,ad-,vect,each,neutral,particle,collisionlessly,,ie.,without,external,forces,or,inter-particle,interaction,,and,at,each,time,step,perform,an,ionisation,procedure.,This,ionisation,procedure,requires,evaluation,of,the,plasma,density,and,temperature,at,the,location,of,each,particle,to,calculate,an,ionisation,rate.,The,ionisation,rate,determines,the,proportion,of,the,weight,of,each,neutral,particle,which,will,be,deposited,onto,the,grid,as,plasma.,3,The,deposition,of,newly,ionised,plasma,takes,the,form,of,source,terms,in,the,plasma,density,,mo-,mentum,density,and,energy,density,equation.,For,conservation,of,mass,,momentum,and,energy,,the,source,terms,in,the,plasma,equations,must,be,commensurate,with,the,loss,of,neutral,density,,momentum,and,kinetic,energy,from,particles.,We,weakly,equate,particle,and,finite,element,repre-,sentations,of,a,quantity,via,a,L2,Galerkin,projection.,We,describe,the,ionisation,algorithm,below,and,describe,the,ionisation,rate,calculation,the,next,section.,for,particle,species,do,∈,position(particle),weight(particle),velocity(particle),density(x),temperature(x),rate(T,),x,←,w,←,u,←,n,←,T,←,r,←,∆w,remove,←,if,w,+,∆w,<,0,then,rnw∆t,0,←,−,∆w,remove,←,−,w,1,←,end,if,w,nS,uE,end,for,←,←,←,w,+,∆w,nS,uE,−,−,∆wu,∆wmu2/2,2.1.2,Ionisation,rates,We,use,the,following,formula,,taken,from,Ref.,[1],,as,the,ionisation,cross,section,σ,=,N,(cid:88),i=1,aiqi,(cid:17),ln,(cid:16),E,Ei,EEi,[1,−,(cid:16),E,Ei,(cid:17),−1,],−ci,bie,(5),where,E,is,the,energy,of,the,impacting,electron,,Ei,is,the,binding,energy,of,electrons,in,subshell,i,and,qi,is,the,number,of,electrons,in,the,i-th,subshell.,The,constants,ai,,bi,and,ci,are,to,be,determined,from,theory,or,experiment.,Note,that,for,hydrogen,and,Helium,like,ions,N,=,1.,In,the,following,table,we,list,the,reference,values,presented,in,[1],for,hydrogen,and,i,=,1.,Parameter,a1,b1,c1,E1,q1,Value,4x10−14cm2(eV,)2,0.6,0.56,13.6,eV,1,4,Figure,1:,Cross,section,plotted,against,electron,energy.,Figure,1,plots,this,approximation,for,the,cross,section,against,the,data,from,Freeman,and,Jones,[2].,We,use,the,cross,section,applied,to,a,Maxwellian,electron,distribution,at,temperature,T,to,generate,the,ionisation,rate,in,units,of,cm3/s,(see,Ref.,[1]),,S,=,6.7x107,N,(cid:88),i=1,aiqi,T,3,2,(cid:90),∞,Ei,T,e−x,x,dx,−,=,6.7x107,−,N,(cid:88),i=1,aiqi,T,3,2,T,Ei,(cid:16),Ei,(cid:17),Ei,T,−,−,(cid:40),T,Ei,(cid:40),bieci,Ei,T,+,ci,bieci,Ei,T,+,ci,(cid:90),∞,Ei,T,+ci,(cid:41),e−y,y,dy,(cid:41),(cid:17),(cid:16),Ei,Ei,T,−,ci,−,(6),where,T,is,in,eV.,Figure,2,compares,the,approximation,against,ADAS,provided,data.,We,assess,that,the,approximation,is,adequate,for,this,report.,Furthermore,we,use,the,closed,form,approxi-,mation,to,the,exponential,integral,by,Barry,et,al.,[3].,5,Figure,2:,Comparison,between,ADAS,data,and,approximation,plotted,against,electron,energy.,2.2,NEKTAR++,Although,NEKTAR++,includes,implementations,of,the,prerequisite,algorithms,required,by,NESO,to,project,particle,data,onto,finite,element,functions,and,evaluate,these,functions,at,particle,locations,these,methods,are,not,efficiently,implemented,for,the,particle,use,case.,The,original,implementa-,tions,in,NEKTAR++,are,adequately,implemented,for,use,in,the,setup,phase,of,a,simulation,where,small,inefficiencies,are,amortised,by,the,cost,of,the,time,stepping,loop,or,main,solve.,In,the,particle,use,case,,these,once,setup,only,functions,are,called,on,a,per-particle,basis,at,every,time,step,of,the,simulation.,This,per,particle,and,per,time,step,regime,places,these,functions,on,the,critical,path,of,the,simulation,where,small,inefficiencies,in,each,function,call,contribute,significantly,to,the,overall,runtime.,Furthermore,the,implementations,of,these,methods,are,provided,only,on,the,CPU,host,hence,particle,data,,which,is,stored,on,a,compute,device,such,as,a,GPU,,must,be,copied,to,and,from,the,host,to,call,these,functions.,To,improve,the,performance,and,portability,of,NESO,we,provide,SYCL,implementations,of,these,critical,algorithms.,2.2.1,Function,Evaluation,NEKTAR++,considers,two,representations,of,a,given,function.,The,first,representation,in,the,func-,tion,values,at,the,quadrature,points,in,the,domain.,The,second,representation,is,the,classic,finite,element,degrees,of,freedom,(DOFs),αj,such,that,u((cid:126)x),=,(cid:88),j,αjψj((cid:126)ξ((cid:126)x)),6,(7),where,u,is,the,function,represented,by,the,αj,,ψj,is,the,jth,basis,function,and,(cid:126)ξ((cid:126)x),is,the,position,in,the,reference,cell,of,the,point,(cid:126)x.,Hence,evaluating,Equation,(7),at,a,point,(cid:126)x,requires,evaluating,all,the,individual,basis,functions,at,(cid:126)x.,We,note,that,in,the,case,where,a,function,u,is,evaluated,at,all,particle,locations,with,a,mesh,cell,there,are,two,properties,that,are,computationally,beneficial,on,modern,architectures;,firstly,the,DOFs,αj,are,constant,and,identical,for,all,evaluations,in,the,cell,secondly,the,functional,form,of,the,basis,functions,is,identical,for,all,the,evaluations,within,the,cell.,The,first,property,indicates,that,evaluating,a,function,at,all,points,within,a,cell,as,a,grouped,operation,should,improve,cache,reuse,and,hence,improve,the,arithmetic,intensity,of,the,implementation.,The,second,property,indicates,that,the,basis,function,evaluations,for,different,points,in,the,cell,can,occur,in,adjacent,lanes,in,vector,computing,architectures,such,as,GPUs,and,CPU,SIMD,units.,Until,now,we,have,indexed,basis,functions,and,DOFs,with,a,single,index,j,,in,practice,there,exists,at,least,n,basis,function,indices,for,an,n,dimensional,space,(in,NEKTAR++).,In,this,report,we,consider,2-D,space,where,the,computational,mesh,is,constructed,as,a,collection,of,quadrilateral,and,triangular,elements.,Functions,defined,over,quadrilaterals,and,triangles,are,constructed,as,the,tensor,product,of,a,1-D,basis,function,in,the,first,and,second,dimensions,of,the,reference,element.,For,these,two,element,types,there,are,two,1-D,basis,functions,defined,in,NEKTAR++;,the,“modified,A”,and,“modified,B”,basis,which,we,label,as,ψa,p,(z),and,ψb,For,n,“modes”,(polynomial,degree,plus,one),the,modified,A,basis,ψa,p,q(z),respectively.,p,is,defined,as,ψa,p,(z),=,,,,1,z),2,(1,−,1,2,(1,+,z),z),1,1,2,(1,−,if,p,=,0,if,p,=,1,,,,(8),2,(1,+,z)P,(1,1),p−2,(z),otherwise,for,all,p,such,that,0,of,order,p.,The,ψb,p),q,<,(n,n,,0,≤,≤,−,},p,<,n.,In,this,report,P,(α,β),p,p,q(z),basis,is,defined,over,a,triangle,number,sized,iteration,set,(z),denotes,the,standard,(α,,β),Jacobi,polynomials,p,<,(p,,q),:,0,{,≤,and,is,defined,as,ψb,p,q(z),=,,,,ψa,q,(z),(cid:0),1,2,(1,+,z)(cid:1)p,z),(cid:0),1,1,2,(1,−,2,(1,+,z)(cid:1)p,P,(2p−1,1),q−1,if,p,=,0,if,q,=,0,(z),otherwise.,,,,(9),A,scalar,valued,function,u,of,n,modes,defined,over,a,quadrilateral,by,coefficients,α,is,evaluated,at,a,point,(cid:126)x,as,the,double,sum,u((cid:126)x),=,n−1,(cid:88),n−1,(cid:88),i=0,j=0,αjn+iψa,i,((cid:126)ξ((cid:126)x)0)ψa,j,((cid:126)ξ((cid:126)x)1).,(10),≤,O,each,with,j,<,n,},(n2),complexity,by,using,the,recursion,relation,of,the,Jacobi,We,evaluate,Equation,(10),with,j,:,polynomials,to,compute,and,store,the,values,(n),complexity.,These,basis,function,evaluations,are,then,reduced,along,0,with,the,appropriate,coefficient,α,to,construct,u((cid:126)x).,This,two,stage,approach,places,the,code,that,(n2),operation,contains,multiple,conditionals,into,the,to,consist,purely,of,a,double,loop,and,a,reduction.,We,exploit,the,local,memory,feature,of,SYCL,to,allocate,temporary,memory,for,the,basis,function,evaluations,in,each,dimension.,(n),basis,evaluation,stage,and,allows,the,j,((cid:126)ξ((cid:126)x)1),ψa,{,i,((cid:126)ξ((cid:126)x)0),ψa,{,i,<,n,i,:,0,and,O,O,O,≤,∀,∀,},7,The,evaluation,of,a,function,u,over,a,triangle,follows,the,same,structure,with,three,modifications;,the,modified,B,basis,is,applied,in,the,second,dimension,,a,correction,is,applied,for,the,second,mode,and,the,reference,coordinate,is,converted,into,a,collapsed,coordinate,prior,to,evaluation.,For,triangular,elements,modes,,and,hence,DOFs,,are,not,indexed,in,a,lexicographic,manner,but,are,arranged,such,that,mode(p,,q),=,q,+,(cid:88),,,(cid:88),,,.,1,0