TN-03_AdvancedQuantificationUncertaintiesInFusionModellingAtExascaleModelOrderReductio
======================================================================================

.. meta::
   :description: technical note

   :keywords: Advanced,Quantification,of,Uncertainties,In,Fusion,modelling,at,the,Exascale,with,model,order,Reduction,(AQUIFER),UKAEA,NEPTUNE,2057701-TN-03,Final,report,Serge,Guillas,University,College,London,14,March,2024,,revised,16,April,2024,Contents,1,Project,team,and,deliverables,2,Events,and,reporting,3,Software,outputs,4,Functional,emulation,for,temperature,profile:,reduced,order,modelling,5,Deep,Gaussian,processes,and,physics-aware,emulation,6,Gradient-based,sequential,design,for,localising,sharp,changes,7,Stochastic,emulation,by,heteroscedastic,Gaussian,processes,8,Application:,uncertainty,quantification,and,sensitivity,analysis,for,molecular,dy-,namics,simulations,with,kernelized,active,subspace,9,Polynomial,chaos,expansion,and,stochastic,collocation,(surrogate,models:,novel,simplex,stochastic,collocation),10,Data,assimilation,in,plasma,physics,models,11,FabNESO,12,Additional,Bayesian,calibration,tools,13,Software,implementation,and,HPC,deployment,Acronyms,References,2,2,3,4,8,13,18,19,22,23,39,44,44,48,50,1,1,Project,team,and,deliverables,Many,activities,are,joint,with,the,epsrc,seavea,project,due,to,synergies,and,cost,efficiencies,(e.g.,principal,investigator,(pi),time,,events,,high,performance,computing,(hpc),access,to,ARCHER2,,in-,ternational,interactions,and,dissemination),as,mentioned,in,the,aquifer,project,description.,The,team,includes:,•,Prof,Serge,Guillas,,pi,,•,Prof,Peter,Coveney,,co-investigator,,•,Dr.,Xiao,Xue,,research,fellow,,•,Dr.,Matthew,Graham,,senior,research,data,scientist,,•,Dr.,Duncan,Leggat,,research,software,engineer,,•,Dr.,Tuomas,Koskela,,principal,research,software,engineer,,•,Yiming,Yang,,PhD,student.,Below,we,present,our,deliverables,and,their,locations,in,the,report.,final,report,accelerated,Bayesian,calibration,stochastic,emulation,for,particles,da,approaches,integrated,approaches,multi-physics,physics-aware,and,deep,gp,seavea,tookit,and,da,software,implementation,and,hpc,deployment,active,subspace,single,models,emulation,of,a,pic,proxyapp,example,of,emulation,of,coupled,proxyapps,example,of,benefits,of,physics-aware,and,deep,gps,scalability,benchmark,da,filter,scalability,coupled,uq,and,da,D1.2,D2.2,D2.3,D2.4,D2.5,D2.6,D3.1,D3.2,M2.1,(M6),M2.3,(M8),M2.5,(M10),M2.6,(M11),D4.1,D4.2,D4.3,M4.1,(M16),uq,campaign,single-scale,,ARCHER2,M4.2,(M16),da,campaign,single-scale,,ARCHER2,this,report,see,sections,11.2,and,12,see,section,7,see,section,10,see,section,4.4,see,section,5,see,section,13,see,section,13.1,see,section,8,see,section,5,see,section,4.4,see,section,5,6,see,sections,11.1,and,10.5,see,section,10.5,see,section,10.5,see,section,11.1,see,section,10.5,Table,1:,Deliverables:,aquifer.,M2.3,(M8):,this,task,was,largely,not,possible,as,originally,stated,due,to,the,lack,of,stochasticity,in,the,current,versions,of,NESO-Particles,,so,the,work,in,5,on,stochastic,emulation,was,done,instead.,The,work,UKAEA,did,(in,conjunction,with,UCL,,at,the,December,2022,hackathon),on,NESO-UQ,(https://github.com/ExCALIBUR-NEPTUNE/NESO-UQ),applies,forward,UQ,to,non-stochastic,electrostatic,PIC,and,could,easily,be,extended,to,apply,the,EasySurrogate,part,of,the,toolkit,to,PIC.,2,Events,and,reporting,Our,level,of,direct,interactions,has,been,high,,as,we,meet,and,report,•,fortnightly,at,neptune,progress,meetings,on,Thursdays,,•,fortnightly,in,aquifer,technical,meetings,on,Mondays:,Ed,Threlfall,,Ander,Gray,and,James,Buchanan,represent,ukaea.,2,•,Hack,session,with,the,developers,on,2022-05-06,fully,capturing,the,new,work,in,seavea,toolkit,(FabNEPTUNE),•,ukaea,workshop,,05-06/09/2022,(Abingdon),•,ucl-ukaea,collaborative,workshop,,10/03/2023,,ucl,Our,events,have,been,joint,with,epsrc,seavea.,These,have,included:,Hackathons:,•,ICCS,seavea,Hackathon,(24,June,2022),at,Brunel,University,London,•,seavea,Hackathon,at,University,College,London,(8-9,December,2022),•,seavea,Hackathon,at,Brunel,University,London,(26-27,June,2023),•,CompBioMed,&,seavea,Hackathon,at,LRZ,,Garching,,Germany,(15,September,2023),•,seavea,Hackathon,at,ucl,(1,,4-5,December,2023),Workshops:,•,ICCS,seavea,Workshop,titled,‘verification,,validation,and,uncertainty,quantification,(vvuq),on,the,Exascale’,at,Brunel,University,London,(23,June,2022),•,SEAVEA,Applications,meeting,at,University,College,London,(12,December,2022),•,CompBioMed,conference,at,Science,Congress,Center,Munich,,Garching,,Germany,(12-14,Septem-,ber,2023),•,CIUK,ExCALIBUR,session,at,Manchester,Central,,Manchester,,UK,(7-8,December,2023),•,seavea,workshop,at,Physics,and,Applied,Mathematics,Unit,in,ISI,Kolkata,,India,on,21st,De-,cember,2022.,•,seavea,workshop,at,Centre,of,Excellence,in,climate,modelling,in,Indian,Institute,of,Delhi,,India,on,9th,May,2023,•,seavea,knowledge,exchange,at,hpc,for,ABMs,Workshop,Aberdeen,(Feb/March,2024),•,seavea,knowledge,exchange,at,Leogang,Workshop,(Jan,2024),3,Software,outputs,A,variety,of,different,software,outputs,were,produced,during,the,project,,including,development,of,new,and,existing,packages,and,integrations,with,other,software,outputs,produced,in,the,wider,NEPTUNE,project.,Table,2,provides,a,summary,of,the,software,repositories,worked,on,in,the,project,and,mentioned,in,this,report.,3,Repository,Team-RADDISH/ParticleDA,Team-RADDISH/,ParticleDA-UseCases,UCL/neso-calibration,UCL/calibr,UCL/FabNESO,UCL-CCS/FabNEPTUNE,djgroen/FabParticleDA,wedeling/MD-active-subspace,Description,Julia,package,for,distributed,particle-filter,based,data,as-,similation.,Example,integrations,of,ParticleDA,with,external,mod-,els,,specifically:,NektarDriftwave,,a,Nektar++,solver,for,a,plasma,physics,model;,anaet,,a,toy,ordinary,differential,equation,model,of,magnetohydrodynamic,instability.,Python,package,providing,a,simple,interface,for,executing,neso,solvers,and,associated,examples,calibrating,neso,mod-,els.,Parallelized,Bayesian,calibration,of,simulations,using,Gaus-,sian,process,emulation.,neso,plugin,for,FabSim3,,facilitating,execution,of,neso,sim-,ulations,on,both,local,and,remote,high,performance,comput-,ing,systems,via,a,unified,interface.,Plugin,for,FabSim3,,facilitating,execution,of,neptune,sim-,ulations,(Convection2D,and,Convection3D).,Plugin,for,FabSim3,,facilitating,execution,of,ParticleDA,fil-,tering,runs.,Python,implementations,of,Deep,active,subspace,and,kernel-,ized,active,subspace,and,the,application,in,high-dimensional,Molecular,Dynamics,(MD),simulations.,Table,2:,Summary,of,repositories,containing,packages,and,code,examples,developed,during,aquifer.,4,Functional,emulation,for,temperature,profile:,reduced,order,mod-,elling,4.1,Study,case:,isotropic,heat,transportation,The,model,for,time,evolution,of,the,temperature,field,T,is,thermal,diffusion,,which,in,a,plasma,after,Braginskii,parametrization,gives:,3,2,N,∂T,∂t,=,∇,·,(cid:16),(k,∥,b(b,∇,T,),+,k,⊥,T,(,∇,−,b[b,∇,(cid:17),T,]),where,b,indicates,the,magnetic,field,and,k,are,the,thermal,conductivities,with,respect,to,the,∥,direction,parallel,or,perpendicular,along,the,b.,The,homogeneous,Dirichlet,boundary,conditions,are,set,on,the,left,,top,,and,right,boundaries,of,the,spatial,domain.,,,k,⊥,In,this,case,,we,are,interested,in,emulating,the,relationship,between,direction,θ,of,b,and,the,temper-,ature,profile,Tx,on,the,bottom,boundary,at,steady,state,by,emulator,such,that:,femulator,:,θ,Tx,,θ,[0,,∈,→,π,2,],,Tx,∈,C[0,,1],4.2,Methods:,outer,product,emulator,The,high-resolution,physical,simulation,generates,high-dimensional,data,over,discretized,spatial,do-,mains,,such,as,wave,heights,in,a,particular,region,of,the,ocean,or,temperature,fields,over,a,heating,system.,However,,the,high,dimensionality,of,the,data,poses,significant,challenges,when,building,sta-,tistical,emulators.,Traditional,emulators,often,struggle,with,scalability,or,accuracy.,To,balance,scalability,and,accuracy,,the,outer,product,emulator,(ope),has,been,introduced.,The,4,Figure,1:,(left),Temperature,transfer,with,incident,angle,θ,(right),Temperature,profile,on,bottom,boundaries,given,θ.,The,occurrence,of,nonphysical,oscillations,,arises,due,to,the,discretization,of,the,domain,and,the,8th,order,polynomial,interpolation,in,the,numerical,solver,Firedrake,(Rathgeber,et,al.,2016).,This,issue,is,a,common,problem,in,polynomial,interpolation,(see,Runge’s,phenomenon,(Boyd,1992;,Wendland,2004)).,(a),Legendre,polynomials,(b),Haar,wavelet,functions,Figure,2:,Choices,for,regressor,functions,ope,creates,a,single,emulator,for,all,simulation,outputs,over,the,entire,domain,and,simplifies,the,representation,of,fitted,functions,using,products,of,component,functions.,In,general,,the,ope,takes,the,form:,f,(r,,s),=,v,(cid:88),j=1,βjgj(r,,s),+,ϵ(r,,s),Here,,f,(r,,s),is,the,simulator,output,with,input,r,at,output,location,s,,gj,are,the,regressor,functions,,βj,are,unknown,coefficients,,and,ϵ,is,the,residual,assumed,to,be,a,Gaussian,process,(gp),such,that,ϵ,GP,(0,,kλ(,,,,)).,∼,·,The,ope,has,two,main,characteristics,that,distinguish,it,from,conventional,multivariate,emula-,tors.,First,,the,covariance,function,of,the,residuals,is,separated,into,inputs,r,and,outputs,s,,that,is,,kλ(r,,s,,r′,,s),=,kr,ks,λ(s,,s′).,Second,,the,regressor,functions,gj,are,given,by,products,j,(r),(cid:78),gs,gj(r,,s),=,gr,In,this,case,,the,Legendre,Poly-,nomials,are,used,as,the,input,regressor,,with,gr,j,(θ),=,.,However,,due,to,the,sharp,},changes,in,the,temperature,profile,shown,in,Figure,1,,commonly,used,smooth,basis,functions,such,as,Legendre,polynomials,and,Fourier,basis,are,not,suitable.,Instead,,we,select,the,Haar,wavelets,as,the,basis,functions,for,our,output,regressors,,see,Figures,2a,and,2b.,λ(r,,r′),j,(s),where,(cid:78),is,the,outer,product,symbol.,1,,6θ2,1,,6θ,×,−,−,{,6,5,Figure,3:,Emulation,results,for,simulation,outputs,(first,row),and,pseudo,observations,with,two,different,θs.,The,red,curves,are,the,mean,emulation,result,,the,blue,dotted,red,curves,are,the,95%,Confidence,Interval,,the,blue,curves,are,the,true,simulation,result,,and,the,dotted,gray,curves,are,sampled,emulation,result.,4.3,Results,We,run,the,simulator,over,40,different,θs,which,are,evenly,chosen,to,cover,the,input,parameter,space.,For,each,run,,there,are,641,discretized,values,over,spatial,field.,We,implement,the,ope,in,two,settings:,•,Emulation,for,simulation:,In,this,case,,We,use,all,the,641,discretized,values,to,fit,the,ope,for,emulating,the,computational,model.,•,Emulation,for,pseudo-experimental,“observation”:,We,uniformly,sampled,80,values,from,the,641,discretized,values,and,added,with,random,noise,for,fitting,a,profile,by,observed,experi-,mental,data.,The,results,demonstrate,the,benefits,and,robustness,(for,pseudo,observations),of,the,reduced,order,modelling,approach,using,functional,emulation,,see,Figure,3.,However,,with,more,input,parameters,,higher,dimensional,outputs,and,coupled,models,,and,for,sharper,input-output,functional,forms,(e.g.,lower,angles),more,investigation,is,under,consideration:,possible,use,of,linked,emulations,(Ming,and,Guillas,2021),and,output,dimension,reduction,(Zhang,,Mak,,and,Dunson,2022),as,well,as,deep,gps,(Ming,,Williamson,,and,Guillas,2023).,4.4,Linked,emulation,for,coupled,proxyapps,In,this,section,,we,demonstrate,the,performances,of,the,Linked,Emulator,for,the,coupled,proxyapps.,The,linked,emulator,consists,of,two,individual,gp,emulators,,which,are,connected,in,a,feed-forward,way.,The,coupled,proxyapp,consists,of,the,an-isotropic,Heat,Transportation,as,shown,in,Section,4.1,and,take,the,boundary,temperature,profile,as,the,input,to,an,isotropic,heat,diffusion,model,with,the,central,temperature,as,the,quantity,of,interest,(qoi),,as,shown,in,Figure,4.,•,We,applied,the,principal,components,analysis,(pca),to,lower,the,dimensions,of,boundary,profile,and,build,the,first,GP,emulator,for,emulating,the,obtained,principal,components.,6,Figure,4:,Demonstration,of,emulated,coupled.,The,top,one,is,the,an-isotropic,heat,transportation,and,the,bottom,one,is,the,isotropic,heat,diffusion,model.,The,inputs,are,the,magnet,strength,and,incident,angle.,The,emulated,quantity,of,interest,(qoi),is,the,central,temperature,,as,shown,as,the,yellow,dot.,Figure,5:,The,emulation,results,and,structure,of,linked,gp,emulators.,•,Then,the,emulated,principal,components,are,used,as,the,input,for,emulating,the,qoi,by,the,second,gp,emulators.,The,emulation,outcomes,and,the,architecture,of,the,linked,emulator,are,illustrated,in,Figure,5.,Here,,the,inputs,X,consist,of,the,incident,angle,θ,and,magnet,strength,B,,with,the,first,model,generating,the,temperature,profile,along,the,bottom,domain,as,its,output.,To,reduce,the,dimensionality,of,the,temperature,profile,(e.g.,,the,resolution,of,the,discretization,mesh),from,600,to,4,,we,employed,pca,,also,known,as,proper,orthogonal,decomposition,(pod),,as,depicted,in,the,middle,of,Figure,4.,The,second,numerical,model,utilizes,the,obtained,temperature,profile,as,the,top,boundary,condition,to,simulate,the,heat,diffusion,process.,Consequently,,the,second,emulator,processes,the,reduced,temperature,profile,as,its,input,,with,the,central,temperature,denoted,as,output,Y,.,7,X,X,.,.,.,X,(1),1,GP,(2),1,GP,.,.,.,(P1),1,GP,(1),2,GP,(2),2,GP,.,.,.,(P2),2,GP,.,.,.,.,.,.,.,.,.,.,.,.,(1),L,GP,(2),L,GP,.,.,.,Y(1),Y(2),.,.,.,(PL),L,GP,Y(PL),Figure,6:,The,generic,deep,Gaussian,process,(dgp),hierarchy,5,Deep,Gaussian,processes,and,physics-aware,emulation,5.1,Deep,Gaussian,processes,A,deep,Gaussian,process,(dgp),is,a,probabilistic,model,that,builds,on,the,idea,of,a,gp,by,extending,it,to,multiple,layers.,The,model,is,comprised,of,multiple,layers,of,GP,models,,where,each,layer,is,a,gp,that,takes,the,output,of,the,previous,layer,as,its,input,(see,Figure,6).,One,of,the,key,advantages,of,dgps,over,traditional,gps,is,their,ability,to,model,complex,,highly,non-linear,functional,forms,more,effectively,,which,makes,them,well-suited,for,emulating,non-stationary,functions,,as,we,will,discuss,in,the,following,case,study.,Furthermore,,dgps,provide,more,informative,uncertainty,estimation,in,predictions,(Ming,,Williamson,,and,Guillas,2023),,making,them,valuable,in,applications,where,accurate,risk,assessment,is,critical.,5.2,Study,case:,Lorenz,model,The,Lorenz,model,is,derived,from,a,two-dimensional,model,of,Rayleigh–B´enard,convection,(rbc),,specified,by,the,system,of,ordinary,differential,equations,(odes),˙x,=,σ(y,˙y,=,x(ρ,˙z,=,xy,−,−,x),z),βz,y,−,(1),−,where,ρ,,σ,,β,are,respectively,the,Rayleigh,number,(Ra),,the,Prandtl,number,(P,r),,and,the,coupling,strength.,For,this,model,,the,Nusselt,number,(N,u),,which,quantifies,the,heat,flux,,is,given,by,N,u,=,1,+,2z∞,is,approximated,using,the,long-term,average,,calculated,over,100,000,N,u,surface,of,the,Lorenz,model,is,discontinuous,and,exists,many,integration,steps.,The,(Ra,,P,r),sub-regimes,behaving,differently,which,are,shown,in,Figure,7.,In,the,study,,we,fixed,the,coupling,strength,to,1.,The,domain,for,Rayleigh,number,and,Prandtl,number,is,[0,,100].,ρ,.,The,term,z,−,∞,5.3,Deep,Gaussian,process,results,−,In,this,section,,we,demonstrate,the,powerful,ability,of,dgps,to,emulate,non-stationary,function:,the,N,u,surface,of,the,Lorenz,model.,We,compare,dgps,with,classical,gps,,and,treed,gps,(Ra,,P,r),Gramacy,and,Lee,2008,which,combine,the,flexibility,of,decision,trees,with,the,smoothness,of,gps,to,emulate,functions,with,complex,structures.,The,sampling,scheme,is,the,Latin,hypercube,sampling,[0.1,,100.1]2,,which,is,implemented,by,the,multi-output,Gaussian,process,(mogp),with,domain,x,emulator,package1.,The,settings,for,each,model,are,outlined,below:,∈,•,gp:,The,gp,model,with,mean,function,m(x),=,0,and,Mat´ern-2.5,Kernel,K(x,,x′),=,(1,+,√5,|,x,−,γ,x′,|,x′)2,+,5(x,−,3γ2,),exp,(,x′,√5,|,x,−,γ,−,),with,length,parameter,γ,=,0.5.,|,•,dgp:,The,dgp,model,consists,of,two,gps,model,above,,with,length,parameters,γ1,=,1.0,,γ2,=,0.5,for,stabilizing,the,numerical,computation.,The,two,gps,are,connected,in,serial:,We,use,the,dgpsi,package2,for,implementing,the,dgp,,which,is,available,in,both,Python,and,R.,1https://github.com/alan-turing-institute/mogp-emulator,2https://github.com/mingdeyu/DGP,8,Figure,7:,The,(Ra,,P,r),extreme,non-stationarity.,Left:,Two-dimensional,contour,plot.,Right:,Three-dimensional,plot.,N,u,surface,of,Lorenz,Model.,There,are,many,level,sets,which,show,an,−,•,Treed,gp:,The,treed,gp,method,employs,a,decision,tree,algorithm,to,divide,a,domain,into,mul-,1,,...,,R,denote,the,R,non-overlapping,regions,partitioned,tiple,non-overlapping,subsets.,Let,r,drawn,from,the,tree,prior.,In,each,region,,a,gp,regression,is,carried,out,using,by,the,tree,a,default,gp,prior,,as,set,in,the,tgp,package3.,The,kernel,used,is,also,Mat´ern-2.5,with,length,parameter,γ,=,0.5,,which,is,the,same,as,the,other,models,compared,to,ensure,fairness.,The,tree,structure,and,each,component,gp,are,jointly,optimized,by,maximizing,the,likelihood.,∈,T,The,qualitative,results,are,presented,in,Figure,8.,As,shown,in,the,figure,,both,gp,and,dgp,start,to,capture,the,discontinuities,with,only,30,samples,,whereas,treed,gp,requires,100,samples.,This,could,be,due,to,the,requirement,that,splits,of,treed,gp,must,be,parallel,to,the,axis,,making,it,less,flexible,with,fewer,samples.,With,100,fitting,samples,,we,observe,that,dgp,outperforms,the,other,models,and,shows,the,discontinuous,boundaries,more,clearly.,In,Figure,9,,we,compute,the,root,mean,square,error,(rmse),to,measure,the,performance,improvements,with,an,increasing,number,of,training,samples.,Each,model,is,fitted,with,n,=,5,,10,,15,,30,,50,,100,,200,samples,and,evaluated,over,2500,test,[0.1,,100.1].,Overall,,dgp,has,the,lowest,rmse,for,points,evenly,spaced,over,the,domain,[0.1,,100.1],×,all,numbers,of,fitting,samples.,However,,rmse,averages,the,results,over,the,entire,domain,,which,may,obscure,the,emulation,performance,of,our,specific,region,of,interest,-,the,discontinuities.,In,the,future,,a,more,comprehensive,metric,should,be,developed,to,evaluate,emulation,results,for,non-,stationary,functions.,In,some,downstream,tasks,,such,as,uncertainty,quantification,and,parameter,calibration,,multiple,runs,of,the,emulator,may,be,required,over,the,same,evaluated,point.,Therefore,,it’s,important,to,consider,the,running,time,of,the,emulator.,We,measured,the,predicting,time,over,2500,test,points,and,demonstrated,the,fitting,and,predicting,time,of,each,model,in,Figure,10.,We,observed,that,the,fitting,and,predicting,time,of,gp,and,treed,gp,remained,constant,with,the,increase,of,the,number,of,fitting,samples.,On,the,other,hand,,the,fitting,time,of,dgp,increased,linearly,with,the,increase,of,the,number,of,fitting,samples,,whereas,the,predicting,time,increased,exponentially.,This,is,a,crucial,limitation,in,large-scale,and,high-resolution,emulation.,In,summary,,dgps,outperforms,standard,gps,and,treed,gps,in,emulating,non-stationary,functions,in,terms,of,accuracy,and,boundary,detection.,However,,the,heavy,computational,burden,may,hinder,its,applications,in,large-scale,and,high-resolution,emulations.,5.4,Physics-aware,emulation,In,this,section,,we,integrate,established,principles,of,physics,into,our,sampling,methodologies,to,assess,their,potential,for,emulating,non-stationary,functions.,We,again,consider,the,simplified,Lorenz,model,for,rbc,(see,Equation,1),,here,fixing,the,coupling,strength,to,β,=,8,3,.,3https://cran.r-project.org/web/packages/tgp/index.html,9,Figure,8:,The,mean,of,emulation,result:,gp,,dgp,treed,gp,and,true,surface,(left,to,right),with,5,,30,,100,,200,samples,(top,to,bottom),Figure,9:,The,root,mean,square,error,(rmse),of,emulation,means,with,respect,to,the,number,of,fitting,samples.,10,Figure,10:,Comparison,of,fitting,(left),and,prediction,(right),time,with,respect,to,the,number,of,fitting,samples,over,2500,test,samples.,Figure,11:,Left:,The,upper,and,lower,asymptotes,P,ru(Ra),and,P,rl(Ra).,Right:,Physics-aware,sampling,region.,−,In,the,Lorenz,model’s,(Ra,,P,r),N,u,surface,,there,are,known,asymptotes,that,dictate,the,stability,of,the,fixed,points,in,Lorenz,dynamics,(Dullin,et,al.,2007).,These,asymptotes,are,represented,by,the,discontinuous,boundaries,shown,in,Figure,7,,and,are,characterized,by,upper,and,lower,limits,P,ru(Ra),and,P,rl(Ra),(see,Figure,11),,P,ru,=,Ra,We,created,a,physics-aware,sampling,region,that,covers,both,asymptotes,,with,a,margin,of,15,units,(as,shown,in,Figure,11).,To,cover,the,non-cube,region,,we,utilized,two,latin,hypercube,designs,(lhds),(as,depicted,in,Figure,12).,In,the,following,investigation,,we,compared,the,physics-aware,sampling,approach,with,the,lhd,sampling,method,over,the,entire,domain.,In,particular,,we,employed,an,initial,sampling,of,one-third,of,the,total,number,of,samples,over,the,entire,domain,for,the,physics-aware,approach,,followed,by,the,remainder,of,the,samples,being,taken,from,the,physics-aware,region.,2(β,+,2),and,P,rl,=,β,+,1.,−,5.5,Physics-aware,emulation,result,The,results,presented,in,Figure,13,demonstrate,the,effectiveness,of,physics-aware,sampling,in,identi-,fying,boundary,contours.,The,figure,illustrates,that,gp,with,physics-aware,sampling,requires,fewer,samples,as,compared,to,gp,without,the,sampling,scheme.,However,,dgp,requires,more,samples,to,detect,the,boundary,contours,than,both,GP,with,and,without,physics-aware,sampling,but,dgp,out-,perform,others,in,terms,of,accuracy,with,enough,samples.,Furthermore,,Figure,14,shows,that,both,gp,and,dgp,with,the,physics-aware,sampling,have,lower,rmse,values,compared,to,the,gp,without,the,sampling,scheme.,This,indicates,that,the,physics-aware,sampling,scheme,significantly,enhances,the,emulation,of,non-stationarity,,especially,with,a,limited,number,of,samples,(see,rmse,of,30,sam-,ples,in,Figure,14).,Moreover,,the,rmses,are,also,smaller,for,larger,number,of,samples,compared,to,conventional,space-filling,schemes.,Although,our,setup,has,some,evident,flaws,,such,as,an,overlapping,sampling,region,by,two,lhds,11,Figure,12:,Left:,Physics-aware,sampling,scheme.,Right:,An,example,for,the,sampling,scheme,Figure,13:,Contour,level,plots,of,the,predicted,means,(red,curves),and,true,discontinuous,boundaries,(black,curves).,Columns,represent,methods,used,,namely,gp,,gp,with,physics-aware,sampling,,and,dgp,with,physics-aware,sampling,(left,to,right).,Rows,correspond,to,the,number,of,fitting,samples,used,,with,30,,60,,120,,and,150,samples,(top,to,bottom).,12,Figure,14:,The,rmse,for,gp,,gp,with,physics-aware,sampling,,and,dgp,with,physics-aware,sampling,and,a,small,uncovered,area,at,the,left-bottom,corner,(see,Figure,12,,left),,the,results,show,that,the,sampling,scheme,improves,accuracy,and,sample,efficiency.,These,findings,confirm,the,effectiveness,of,incorporating,prior,knowledge,of,physics,into,the,emulation,approaches.,Although,this,prior,knowl-,edge,is,very,strong,and,can,not,be,available,for,every,problem,,we,can,use,such,applications,under,a,careful,design.,In,conclusion,,the,results,of,this,study,demonstrate,the,potential,benefits,of,combining,machine,learning,with,physical,priors.,Moving,forward,,further,research,in,this,area,will,be,critical,for,advancing,our,understanding,in,complex,systems,and,developing,more,effective,tools,for,solving,real-world,problems.,6,Gradient-based,sequential,design,for,localising,sharp,changes,gp,emulators,typically,assume,the,computer,models,are,stationary,in,the,input,domain,,which,is,usually,a,overly,strong,assumption,in,practice,(Gramacy,2020).,The,presence,of,sharp,changes,,or,even,discontinuities,in,the,output,of,the,computer,model,can,result,in,deterioration,of,gp,models,in,terms,of,both,the,accuracy,and,efficiency,(See,Figure,15).,However,,such,sharp,changes,widely,exist,in,real,world,and,often,indicate,bifurcations,or,critical,transitions,within,the,systems,under,investigation.,These,instances,can,be,found,in,switch-like,behavior,in,genomics,(Gardner,,Cantor,,and,Collins,2000),,bifurcations,in,fluid,dynamics,(Rybin,et,al.,2015;,Dullin,et,al.,2007),,and,steep,transition,of,competing,mechanical,phenomena,in,nuclear,safety,analysis,(Lee,and,McCormick,2011).,In,this,work,,we,consider,the,localization,of,sharp,variations,with,double,purposes:,•,Physical,finding.,Identifying,these,sharp,variations,holds,importance,as,they,can,result,in,significant,implications,of,the,investigated,systems,with,better,understanding.,•,Better,emulation.,Knowing,the,locations,of,the,sharp,changes,can,aid,in,constructing,the,emulators.,Thus,,the,downstream,tasks,can,also,be,better,solved.,6.1,Gradient,of,Gaussian,process,emulator,One,key,property,of,gps,in,our,favor,is,that,any,linear,transformation,,such,as,differentiation,and,integration,,of,a,gp,remain,a,gp,(Rasmussen,,Williams,,et,al.,2006).,Gradient,is,a,linear,differential,ˆf,(x0),is,a,D-dimensional,operator,denoted,as,ˆf,(x0)]d,is,the,partial,derivative,with,respect,to,dth,dimension,of,GP,emulator.,The,dth,element,[,∇,ˆf,is,driven,by,derivatives,of,posterior,predictive,µ0,and,variance,σ2,x0d.,The,predictive,posterior,of,0,(Scheuerer,2010;,McHutchon,2015;,Rasmussen,,Williams,,et,al.,2006),as,,.,Given,a,constructed,gp,emulator,ˆf,,,the,gradient,∇,∇,∇,ˆf,(x0)),=,E(,∇,∇,µ0(x0),and,cov(,ˆf,(x0)),=,∇,∇,⊗,∇,T,σ2,0(x0),(2),From,the,Equation,2,and,predictive,posterior,of,gp,emulator,,we,can,obtain,the,gradient,of,the,predictive,posterior,explicitly,with,mean,and,variance,(Graepel,2003;,Rasmussen,,Williams,,et,al.,13,Figure,15:,Demonstration,of,prediction,error,concentration,in,regions,of,High,variations.,Upper:,functions,with,sharp,variations,(sample,paths,of,a,non-stationary,gp),are,interpolated,using,10,fitting,points.,Predictions,from,a,gp,model,with,mean,zero,and,squared,exponential,kernel,with,length,scale,γ,=,1.,Lower:,absolute,prediction,error.,2006;,Solak,et,al.,2002):,∇,µ0(x0),=,T,σ2,∇,⊗,∇,r(x0)T,R(x)−,∇,0(x0),=,σ2(cid:0)H,1y,and,1,r(x0)T,R(x)−,−,∇,r(x0)(cid:1),,∇,(3),∈,RN,is,the,set,of,emulated,targets,,[R(X)]ij,=,k(Xi,where,y,the,ith,column,of,the,design,matrix,X),,∗,r(x0),=,[,∂k(x0,x),∈,D,is,the,Hessian,matrix,with,(i,,j),element,[H(x0)]ij,=,∂2k(x0,x0),∂x0i∂x0j,∂x01,,,Xj,),+,η1,{,,,...,,∂k(x0,x),∂x0D,Xi∗=Xj∗,]T,RD,×,∇,∗,(the,Xi,},RD,and,H,=,∗,∇,⊗,and,will,be,a,indicates,T,k(x0,,x0),∈,∇,,,constant,if,k(,·,),is,shift,invariant,(Wendland,2004).,·,6.2,Extension,to,linked,Gaussian,process,emulation,In,order,to,deal,with,the,non-stationarity,,in,this,section,,we,extend,the,gradient,information,from,gp,into,linked,Gaussian,process,(lgp),by,taking,the,derivative,of,its,closed-form,expectation,(Ming,and,Guillas,2021;,Ming,,Williamson,,and,Guillas,2023).,Assume,that,the,lgp,is,constructed,with,RD,can,be,M,design,points,obtained,directly,by,taking,the,derivative,of,˜µ0(x0),such,that:,m=1.,The,gradient,of,lgp,at,new,position,x0,=,(xm,,wm,,ym)M,∈,X,⊂,D,where,A,=,R(w)−,with,[,∈,I(x0)]ij,=,∂Ii(x0),∂x0j,1y,∇,˜µ0(x0),=,∇,∇,I(x0)T,A,,(4),RM,,,the,I(x0),.,∈∈,RM,(See,(Ming,and,Guillas,2021)),,and,I(x0),∇,∈,RM,D,×,6.3,Deep,Gaussian,process,gradient,estimation,The,dgp,emulator,can,be,viewed,as,a,lgp,with,latent,internal,Input/Output,(I/O),w.,We,can,impute,N,sets,of,w,iteratively,from,the,conditional,Gaussian,distribution,,and,construct,N,copies,of,lgp,to,approximate,the,dgp,(see,Ming,,Williamson,,and,Guillas,2023,for,more,details).,Given,the,N,imputed,l,)i,,Y,=,y),,we,can,derive,the,approximate,data,denoted,with,(,predictive,mean,and,variance,of,the,dgp,gradient,through,a,mixture,of,the,gradients,from,the,N,constructed,lgp,emulators,,i,=,(X,=,x,,W,=,(wp,i=1,where,i)N,D,D,˜µ0(x0),=,∇,1,N,N,(cid:88),∇,i=1,˜µ0,i(x0),14,(5),By,exploiting,the,estimations,,several,indicators,can,be,used,to,quantify,the,local,variations,,reflecting,the,degree,of,sharp,changes,and,discontinuities.,In,this,work,,we,choose,the,gradient,norm,which,is,defined,by,Q(x0),=,=,≈,=,∥∇,(cid:113),f,(x0),∥,f,(x0)T,∇,(cid:113),∇,˜µ0(x0)T,∇,˜µ0(x0),∥,∥∇,∇,f,(x0),˜µ0(x0),Similarly,,the,estimation,of,gradient,norm,by,dgp,emulator,can,be,obtained,by,,Q(x0),≈,1,N,N,(cid:88),i=1,˜µ0,i(x0),∥∇,∥,where,each,∇,˜µ0,i(x0),is,computed,by,each,ith,imputed,lgp.,6.4,Gradient-based,acquisition,function,(6),(7),In,the,sequential,design,,a,gp-based,emulator,of,the,computer,experiment,is,firstly,constructed,based,on,the,space-filling,initial,design,with,N0,points,i=1.,The,commonly,used,space-filling,D,designs,are,optimized,lhds,and,minimax-distance,designs;,see,(Santner,et,al.,2003),for,an,overview.,Once,the,emulator,is,constructed,,the,sequential,design,begins,a,loop,over,computational,budget,N,(ie.,the,maximum,number,of,computer,model,runs).,In,n,+,1th,step,,there,are,two,steps,:,(a),Using,the,emulator,to,solve,an,acquisition,criteria,Jn,to,determine,the,next,run,of,the,computer,model,such,as,n+1.,More,precisely,,the,next,evaluation,xn+1,is,determined,such,that,(xn+1,,yn+1).,(2),Update,the,emulator,with,data,0,=,(xi,,yi)N0,n+1,=,n,D,D,D,∪,xn+1,=,arg,max,∈X,x,Jn(x),(8),{,x,c,},∈,X,|,gf,(x),where,c,In,this,paper,,we,consider,the,high-variation,region,as,the,super,level,set,of,the,gradient,norm,,such,that,B,=,.,Then,,equivalently,,we,formulate,this,problem,as,a,level,set,estimation,(lse),or,contour,localization,(cl),problem.,The,most,popular,sampling,criterion,for,lse,is,the,entropy,(Oakley,2004;,Gotovos,et,al.,c),2013;,Cole,et,al.,2023;,Booth,,Renganathan,,and,Gramacy,2023).,Denote,the,px,=,p(,as,the,probability,of,x,S,,i.e,x,located,with,in,the,high-variation,region.,The,entropy,acquisition,function,in,lse:,0,is,certain,threshold,and,gf,(x),=,f,(x),∥,f,(x),∥,≥,∥∇,∥∇,≥,R,∈,∈,≥,The,S,will,be,high,in,the,boundary,of,region,B,,∂B,(i.e.,,px,we,are,looking,for.,≈,0.5).,It,is,indeed,the,transition,region,S(x),=,px,log,(px),−,(1,−,−,px),log,(1,px),−,xn+1,=,arg,max,S(x),x,For,the,general,lse,problem,,the,px,can,be,easily,evaluated,as,the,cumulative,distribution,function,(cdf),of,the,predictive,distribution,of,gp,emulator.,In,our,sequential,learning,case,,we,hope,to,allocate,more,samples,around,the,region,B,but,also,achieves,good,explorations,of,the,input,domain.,So,we,do,not,need,to,accurately,evaluate,the,full,distribution,of,gradient,norm.,We,only,need,an,indication,for,the,degree,of,variation,and,the,uncertainty,of,the,input,domain.,We,combine,the,gradient,norm,estimation,of,dgp,and,the,predictive,variance,as,these,indications,to,form,the,new,acquisition,function,,∈X,cand,S(x),=,1,N,N,(cid:88),i=1,x,log,(pi,pi,x),−,(1,−,−,pi,x),log,(1,pi,x),−,and,pi,x,=,1,c,Φtrun(,−,˜µ0,i(x),−,∥∇,˜σ0,i(x),∥,),,15,(9),(a),The,two-dimensional,plateau,function,(b),Emulation,results,of,gradient,norm,(c),Quantitative,analysis,of,emu-,lation,results,in,function,and,gra-,dient,Figure,16:,The,demonstration,of,gradient-norm,emulation.,(a),indicates,the,two-dimensional,plateau,function.,The,black,dashed,line,indicates,the,diagonal,of,the,input,domain,,in,which,we,evaluated,the,gradient,in,(b),,red,points,are,used,to,fit,the,gp,and,dgp.,(b),gp,and,dgp,emulation,results,of,gradient,norm,(mean,+/-,2,standard,deviation,,where,negative,standard,deviation,truncated,at,0).,(c),Quantitative,analysis,of,emulation,results,over,150,points,on,diagonal.,The,mean,and,standard,deviations,of,these,results,are,obtained,over,10,random,initial,designs.,where,Φtrun,is,the,cdf,of,truncated,normal,distribution.,In,most,cases,,the,threshold,c,is,not,known.,We,take,a,substitution,of,the,explicit,threshold,level,by,an,implicit,level,as,cimp,=,(1,α),max,x,∈X,−,˜µ0(x),∥∇,,,∥,where,α,[0,,1],and,is,pre-determined.,We,show,a,demonstrative,figures,in,the,qualitative,and,quantitative,performance,of,the,proposed,acquisition,function,in,Figure,16.,The,Figure,shows,a,simple,case,of,a,two-dimensional,plateau,function,,∈,f,(x),=,2Φ,(cid:16),√2(,4,−,−,6(x1,+,x2),(cid:17),−,1,and,x,1,,1]2,[,−,∈,We,constructed,the,emulators,using,a,dataset,of,10,points,sampled,randomly,via,lhd,,as,shown,in,Figure,16,(left).,For,visualization,purposes,,we,present,the,emulation,results,for,the,gradient,norm,along,the,domain’s,diagonal,line,,indicating,the,next,design,points,selected,by,both,gp,and,dgp,methods.,The,design,point,chosen,by,dgp,more,close,to,the,high,variation,region,,aligning,with,our,objectives.,The,right,panel,provides,a,quantitative,analysis,of,the,gradient,norm,emulation’s,accuracy,,obtained,by,fitting,the,emulators,with,20,different,sets,of,lhd,points.,6.5,Study,case:,phase,transitions,in,Rayleigh-B´enard,convection,We,again,consider,the,simplified,Lorenz,model,for,rbc,(see,Equation,1),,here,fixing,the,coupling,strength,to,β,=,8,3,.,6.5.1,Settings,•,dgp,emulator:,Two-layered,dgp,with,squared,exponential,kernel,•,Design:,Initialize,with,three,lhd,samples,and,adding,designs,(one,in,each,iteration),sequential,up,to,50.,•,Evaluation,metrics:,1.,Localization:,Log,likelihood,of,Bernoulli,distribution,LLK,=,nt(cid:88),i=1,yi,log,(pxi),+,(1,yi),log,(1,pxi),−,−,16,Figure,17:,The,(Ra,,P,r),of,the,phase,transition.,More,details,could,be,found,in,Dullin,et,al.,2007.,N,u,surface,of,Lorenz,Model.,The,dashed,blue,lines,are,known,asymptotes,−,where,the,yi,=,1,2.,Approximation:,normalized,root,mean,square,error,(nrmse),and,the,pxi,=,1,Φ(,xi,−,B,∈,{,},ˆc,µQ(xi),σQ(xi),).,−,NRMSE,=,(cid:113),1,nt,(cid:80)nt,i=1(µ0(xi),maxi=1:nt,f,(xi),−,f,(xi))2,−,mini=1:nt,f,(xi),where,(xi,,f,(xi)nt,i=1,is,the,gridded,test,data,set,with,nt,=,2500.,•,Compared,methods:,1.,gp,with,proposed,sequential,design,,2.,dgp,with,‘mutual,information,for,computer,experiments’,(MICE),Beck,and,Guillas,2016,,3.,dgp,with,‘active,learning,Cohn’,(ALC),Cohn,,Ghahramani,,and,Jordan,1994.,6.5.2,Qualitative,analysis,Figure,18,demonstrates,the,sequential,localization,procedure,of,our,proposed,algorithms.,The,color,indicates,the,probability,that,the,gradient,norm,beyond,the,implicit,threshold,cimp,,which,is,the,phase,transition,region.,We,can,see,that,the,boundaries,are,clearly,captured,by,using,50,samples,,especially,the,diagonal,one.,However,,empirically,,we,found,the,algorithms,are,not,robust,enough,and,sensitive,to,the,hyper,settings,such,as,initial,design,and,cimp.,Our,next,step,will,figure,out,how,to,stabilize,the,algorithms,with,external,physical,constraint.,One,possible,way,is,to,incorporate,the,numerical,parameter,continuation,method,,since,most,of,such,phase,transitions,can,be,attributed,to,the,bifurcations,and,instabilities,of,dynamical,systems,,see,Ma,and,Wang,2005.,6.5.3,Quantitative,analysis,In,addition,to,qualitative,analysis,,we,also,made,a,quantitative,analysis,to,compare,with,other,se-,quential,design,methods,in,this,rbc,problem.,The,results,are,obtained,from,re-running,20,times,with,different,lhd,initial,designs.,From,Figure,19,,due,to,the,existences,of,sharp,variations,,our,designed,methods,outperform,than,other,sequential,design,methods,in,the,accuracy,of,both,localization,and,global,approximation.,17,Figure,18:,Localization,of,high,variation,region,B.,The,color,indicates,the,probability,that,the,gradient,norm,beyond,the,implicit,threshold,cimp.,The,four,surfaces,are,estimated,with,3,,10,,30,,50,samples,from,left,to,right.,Figure,19:,Error,analysis,of,rbc,parameter,surfaces,with,20,random,initial,designs.,Left:,Number,of,design,points,versus,negative,log-likelihood.,Right:,Number,of,design,points,versus,nrmse.,7,Stochastic,emulation,by,heteroscedastic,Gaussian,processes,7.1,Study,Case:,parameter,study,in,Tritium,breeding,ratio,Tritium,breeding,ratio,(tbr),,which,describes,the,number,of,Tritium,atoms,formed,per,incident,neutron,,is,an,essential,quantity,in,sustaining,the,fusion,reaction.,tbr,depends,on,the,large,amount,of,factors,in,a,complex,probabilistic,way.,In,this,case,study,,We,simplify,the,problems,by,only,focusing,on,the,following,two,factors,,•,Material,property:,Lithium-6,enrichment,(%),(percentage,of,Lithium-6,in,the,material,of,breeding,blanket),•,Reactor,geometry:,Blanket,radius,with,fixed,neutron,point,source.,The,neutron,transportation,can,only,be,characterized,with,respect,to,certain,stochastic,process,given,design,parameters,(i.e.,the,material,property,and,reactor,geometry).,The,stochastic,nature,requires,numerous,Monte,Carlo,samples,for,tbr,estimation,,making,it,computationally,difficult,to,explore,the,design,parameter,space,,especially,in,high,dimensional,space.,An,accurate,emulator/surrogate,model,with,reliable,uncertainty,estimation,is,desired,to,capture,its,output,variance.,7.2,Heteroscedastic,Gaussian,processes,gp,regression,provides,a,Gaussian,distribution,approximation,ˆP,(µ,,Σ),to,the,distribution,P,.,However,,the,homogeneity,of,variance,in,standard,gp,is,problematic.,When,design,parameter,θ,changes,,the,physical,system,will,change,a,lot,,as,well,the,inside,stochasticity.,Standard,gp,regression:,y,=,f,(x),+,ϵ,18,Figure,20:,Left:,Standard,gp,regression.,Right:,hgp,regression.,The,dark,dot,is,the,mean,of,Monte,Carlo,samples.,The,fitting,data,consists,of,OpenMC,mean,and,OpenMC,variance.,GP,(µf,(x),,Kf,(x,,x′)),•,f,(x),•,ϵ,∼,N,∼,(0,,σ2),Heteroscedastic,gp,regression:,y,=,f,(x),+,ϵ(x),•,f,(x),∼,GP,(µf,(x),,Kf,(x,,x′)),•,ϵ(x),(0,,r(x)),s.t,,r(x),=,exp,(cid:0)g(x)(cid:1),and,g(x),∼,N,GP,(µg(x),,Kg(x,,x′)),∼,However,,the,marginal,likelihood,and,predictive,density,is,no,longer,analytically,tractable.,There,are,massive,past,works,in,approximate,inference:,variational,inference,,Markov,chain,Monte,Carlo,(mcmc),sampling,Ming,,Williamson,,and,Guillas,2023.,In,Figure,20,,we,compared,the,performances,between,standard,and,heteroscedastic,gps.,Both,of,the,two,models,are,fitted,with,mean,and,variance,estimation,of,tbr,from,batch,of,500,particles,in,the,OpenMC,simulation,with,a,simple,onion,tokamak,configuration.,From,Figure,20,,we,can,clearly,observe,that,the,heteroskedastic,Gaussian,process,(hgp),(right,panel),captures,the,variation,of,variances,better,than,standard,gp,(left,panel).,8,Application:,uncertainty,quantification,and,sensitivity,analysis,for,molecular,dynamics,simulations,with,kernelized,active,sub-,space,The,objective,of,the,present,project,is,to,perform,a,high-dimensional,active-subspace,based,uncer-,tainty,quantification,analysis,to,assess,the,influence,of,aleatoric,,simulation,and,particularly,force-field,parameter,uncertainty,on,classical,molecular,dynamics,(md),predictions,(see,Figure,21).,The,original,f,(x),,active,subspace,approach,requires,the,explicit,estimation,of,the,derivatives,of,model,response,which,are,not,available,in,the,md,simulation.,Instead,,we,applied,gradient-based,kernel,dimension,re-,duction,(gkdr),for,the,dimension,reduction,with,a,gp,denoted,as,kernelized,active,subspace,Gaussian,process,(kas-gp),in,the,following,discussion,,as,the,low-dimensional,surrogate,and,make,a,compari-,son,with,the,recently,proposed,neural-network,based,active,subspace,called,the,deep,active,subspace,(das),method.,Both,of,the,two,approaches,are,derivative-free,and,achieve,similar,performances,in,the,dimension,reduction,,as,well,the,associated,tasks:,uncertainty,quantification,(uq),and,sensitivity,analysis,(sa).,∇,8.1,Investigated,molecular,dynamics,simulations,•,Epoxy-resin,thermosetting,polymer,materials,,predicting,mechanical,properties,,namely,E,,the,Young’s,modulus,,an,indicator,of,the,stiffness,of,the,material,,and,the,Poisson,ratio,,which,is,the,ratio,of,lateral,deformation,to,axial,deformation,when,straining,the,material,in,a,given,direction.,19,Figure,21:,Positions,of,the,atoms,in,the,binding,site,of,the,ligand-protein,complex,,of,which,the,force-field,parameters,are,most,sensitive,in,enhanced,sampling,of,molecular,dynamics,with,approximation,of,continuum,solvent,(esmacs),study.,The,ligand,is,represented,as,bond,,and,the,protein,is,shown,as,ribbon,in,white.,The,residues,at,the,binding,site,are,shown,as,ball,and,stick.,The,sp3,carbon,atoms,(coloured,orange),and,attached,hydrogen,atoms,(yellow),from,protein,have,the,most,sensitive,parameters,,along,with,hydrogen,atoms,(green),from,ligand.,The,parameters,for,the,oxygen,atoms,(purple),from,two,tyrosine,residues,are,also,important,to,the,sensitivity.,•,protein-ligand,biomolecular,systems,,where,binding,free,energies,are,computed.,Two,types,of,binding,free,energies,are,investigated:,the,absolute,binding,free,energy,(shortened,as,“binding,free,energy”,hereafter),which,is,a,quantitative,measure,of,the,strength,of,protein-ligand,binding,,and,the,relative,binding,free,energy,which,is,the,difference,of,binding,free,energies,between,two,ligands.,Reliable,prediction,of,such,properties,plays,an,important,role,in,drug,discovery,and,personalized,medicine,Within,this,study,,we,have,used,two,different,molecular,dynamics,engines,,LAMMPS,for,epoxy-,resin,polymer,materials,and,NAMD,for,protein-ligand,biomolecular,systems.,The,absolute,and,relative,binding,free,energies,are,calculated,using,the,enhanced,sampling,of,molecular,dynamics,with,approxi-,mation,of,continuum,solvent,(esmacs),and,thermodynamic,integration,with,enhanced,sampling,(ties),protocols,,respectively.,8.2,Results,The,epoxy,kas-gp,and,das,surrogates,for,E,and,the,Poisson,ratio,are,plotted,vs,the,d,=,1,dimensional,active,subspace,y1,in,Figure,22a,and,22b.,Note,that,the,one-dimensional,function,captures,the,overall,Importantly,,the,variation,of,the,md,data,f,(x),at,a,given,trend,of,the,data,well,,for,both,qois.,location,in,the,active,subspace,,i.e.,Var[f,(x),y1],,is,heavily,concentrated,around,the,prediction,with,∥,quantified,uncertainty,from,gp.,This,holds,for,both,the,training,data,and,the,test,data,(10%,of,the,data,set),,which,was,not,used,in,constructing,the,surrogates.,Hence,,while,the,original,md,model,is,a,function,of,a,103-dimensional,input,space,,,it,is,well,approximated,by,a,one-dimensional,surrogate.,The,one-dimensional,esmacs,binding-energy,surrogates,are,shown,in,Figure,22c.,A,one-dimensional,active,subspace,is,clearly,visible,,although,the,variance,of,the,training,and,test,data,around,the,surrogate,is,larger,than,for,the,epoxy,surrogates.,The,ties,one-dimensional,relative,binding,free,energy,surrogates,are,shown,in,Figure,22d.,The,one-dimensional,active,subspace,is,again,visible,,although,(like,the,esmacs,case),it,is,not,of,the,same,quality,as,for,the,epoxy,surrogates.,8.3,Global,derivative-based,sensitivity,analysis,To,assess,which,inputs,are,most,influential,,commonly-used,options,are,global,variance-based,sensi-,tivity,methods,such,that,,νi,:=,(cid:19)2,(cid:90),(cid:18),∂f,∂xi,p,(x),dx.,(10),These,indices,measure,the,(average),sensitivity,of,f,to,small,perturbations,in,the,inputs,x,,and,are,especially,suited,for,identifying,non-influential,parameters.,To,connect,(10),to,the,active,subspace,20,(a),E11,(b),Poisson,Ratio,(c),Binding,free,energy,(d),Relatively,binding,free,energy,Figure,22:,The,das,and,kas-gp,surrogates,of,all,qois,plotted,along,the,first,active,variable,y1.,method,,note,that,the,νi,are,the,diagonal,elements,of,the,C,matrix;,[ν1,,·,·,·,,,νd]T,=,diag,(C),.,(11),∇,f,will,not,be,available.,In,das,method,,the,derivative,can,be,estimated,by,the,In,most,cases,,automatic,differentiation,via,the,deep,neural,network.,To,establish,a,connection,between,Equation,10,and,the,gkdr,method,,we,can,view,the,νi,as,the,global,derivative-based,sensitivity,index,in,the,embedded,reproducing,kernel,Hilbert,space,(rkhs).,By,the,reproducing,property,of,k,,,we,can,write,X,=,x(cid:3),and,represents,the,partial,derivative,of,the,embedded,conditional,Y,the,ϕi,=,∂,,,Y,),∂xi,|,expectation.,Then,we,can,estimate,the,sensitivity,index,νi,by,gkdr,such,that,E(cid:2)k,Y,(,·,i,,,ϕx,ϕx,i,⟩H,⟨,Y,dP,(x),n,(cid:88),Diag(,ˆM,(xj))i,(12),(cid:90),˜νi,=,1,n,≈,j=1,=,Diag(,˜Mn)i,This,approach,mirrors,the,methodology,employed,in,both,active,subspace,and,das,methods,,where,the,sensitivity,index,νi,is,estimated,using,the,i-th,diagonal,element,of,the,matrix,˜Mn.,The,derivative,information,is,implicitly,provided,through,the,derivative,kernel,∂k(,for,i,=,1,,..,,d,,facilitating,a,·,∂xi,derivative-free,approach.,A,comparative,analysis,of,the,sensitivity,analysis,results,obtained,through,both,gkdr,and,das,methods,(refer,to,Figure,23),shows,a,high,degree,of,similarity.,This,empirical,evidence,substantiates,the,theoretical,link,between,gkdr,and,active,subspace,methods,,demonstrating,their,agreement,in,practical,applications.,,x),21,(a),E11,(b),Poisson,Ratio,(c),Binding,free,energy,(d),Relatively,binding,free,energy,Figure,23:,The,25,largest,global,derivative-based,sensitivity,indices,for,the,das,and,kas-gp,(stands,for,gkdr,with,low,dimensional,gp,surrogate),method,,ordered,as,ν1,ν25.,The,indices,are,normalized,by,ν1,and,ordered,according,to,the,das,ranking,,although,the,kas-gp,ranking,is,similar.,≥,·,·,·,≥,ν2,≥,9,Polynomial,chaos,expansion,and,stochastic,collocation,(surrogate,models:,novel,simplex,stochastic,collocation),As,stated,in,the,proposal,,in,D2.1,(Surrogate,models,,Model,Order,Reduction,and,verification),we,have,extended,EasyVVUQ,with,capability,to,handle,qois,which,display,discontinuities,or,high,gradient,in,the,stochastic,input,space.,This,is,important,since,surrogate,models,based,on,global,polynomials,(polynomial,chaos,expansion,(pce),/,stochastic,collocation,(sc)),will,display,nonphysical,oscillations,,see,Figure,24.,In,particular,we,implemented,the,simplex,stochastic,collocation,(ssc),method,(Edeling,,Dwight,,and,Cinnella,2016),,which,(unlike,sc),employs,the,Delaunay,triangulation,to,discretize,the,probability,space,into,simplex,elements.,Using,such,a,multi-element,technique,has,the,advantage,that,local,mesh,adaptation,can,be,performed,,such,that,only,regions,of,interest,are,adaptively,refined.,The,ssc,method,is,‘physics-aware’,(related,to,D2.6),in,the,sense,that,it,can,automatically,detect,regions,in,the,input,space,where,the,function,is,irregular.,It,achieves,this,by,enforcing,the,so-called,Figure,24:,Discontinuous,test,function,and,standard,SC,surrogate.,22,Figure,25:,Local,mesh,and,resulting,ssc,surrogate,for,the,test,function,in,24,a.,local,Extremum,Conserving,(LEC),limiter,in,all,simplex,elements.,Skipping,details,for,brevity,,this,limiter,flags,elements,through,which,a,discontinuity,runs,,increasing,the,probability,that,the,code,is,evaluated,within,that,element,at,the,next,iteration.,Simultaneously,,the,local,polynomial,order,of,the,interpolation,point,stencil,of,that,element,is,reduced,,which,avoids,nonphysical,oscillations,as,in,Figure,24,b.,The,SSC,sampling,plan,and,surrogate,are,shown,in,Figure,25,,which,can,be,replicated,by,running,the,Jupyter,notebook,tutorial,in,the,EasyVVUQ,‘tutorials’,folder,4.,Next,steps,involving,applying,the,ssc,method,to,a,Nektar++,case.,10,Data,assimilation,in,plasma,physics,models,Data,assimilation,encompasses,algorithmic,approaches,for,combining,assumptions,about,a,system,formulated,as,a,numerical,model,,with,observations,of,the,system,,to,estimate,how,the,state,of,the,system,evolves,over,time.,As,a,distinct,field,,data,assimilation,has,its,origins,in,numerical,weather,prediction,,where,combining,physical,models,with,data,is,the,central,task.,More,generally,data,assimilation,can,be,considered,a,form,of,(Bayesian),statistical,inference,in,which,we,combine,prior,information,(the,numerical,model),with,observed,data,to,infer,a,posterior,distribution,on,the,system,state,given,both,observations,and,prior,model.,Numerous,data,assimilation,approaches,have,been,proposed.,For,linear-Gaussian,models,—,that,is,models,with,linear,state,dynamics,,observations,that,linearly,depend,on,the,state,and,Gaussian,observation,and,state,noise,—,the,Kalman,filter,(Kalman,1960),offers,an,efficient,approach,for,ex-,actly,computing,the,sequence,of,distributions,on,the,system,state,at,each,observation,time,given,all,observations,up,to,that,time,(the,filtering,distributions),,by,recursively,updating,estimates,of,the,mean,and,covariance,of,the,Gaussian,filtering,distributions.,For,models,with,non-linear,state,tran-,sitions,or,observation,operators,,the,extended,Kalman,filter,can,be,used,,which,use,linearizations,of,the,state,and,observation,updates,computed,using,the,Jacobians,(matrices,of,partial,derivatives),of,the,corresponding,operators.,O,(d2,x),memory,and,As,(extended),Kalman,filters,require,updating,a,covariance,matrix,which,is,of,size,dx,dx,where,×,dx,is,the,state,dimension,,and,solving,a,linear,system,in,a,matrix,of,size,dy,dy,,where,dy,is,the,observation,dimension,,they,become,infeasible,to,apply,to,models,for,which,dx,or,dy,are,very,large,,y),floating,point,operation,costs.,Ensemble,Kalman,filters,(Evensen,1994;,with,Evensen,2006),offer,an,alternative,approach,that,form,a,Monte,Carlo,approximation,to,the,Kalman,filter,updates,using,an,ensemble,of,N,particles,,with,typically,N,(N,dx),memory,(min(dy,,N,)3,+,min(dy,,N,)2,max(dy,,N,),+,dx,min(dy,,N,)N,),floating,point,operations,costs,of,and,ensemble,Kalman,filter,updates,(Mandel,2006;,Roth,et,al.,2017),allowing,application,to,models,with,much,larger,state,dimensions,,including,use,in,operational,numerical,weather,prediction,systems.,dx,,with,the,(d3,≪,O,O,O,×,4https://github.com/UCL-CCS/EasyVVUQ/blob/dev/tutorials/simplex_stochastic_collocation_tutorial.,ipynb,23,While,they,can,be,applied,to,models,with,non-linear,state,transitions,and,observation,operators,,the,ensemble,Kalman,filter,updates,will,only,converge,to,exactly,computing,the,filtering,distributions,as,N,for,linear-Gaussian,models,,and,the,quality,of,their,estimates,can,be,poor,in,models,with,strongly,non-linear,dynamics,or,observations,,or,non-Gaussian,noise.,→,∞,Variational,methods,(Rabier,and,Liu,2003),are,another,widely,used,approach,for,data,assimilation.,Here,,the,problem,of,estimating,a,sequence,of,distributions,on,the,system,state,,is,relaxed,to,instead,computing,a,point,estimate,of,the,‘most,likely’,state,at,each,observation,time.,This,is,implemented,by,iteratively,minimizing,a,cost,function,which,combines,terms,accounting,for,our,prior,knowledge,about,the,system,(based,on,short-term,forecasts,from,previous,state,estimates,and,an,assumed,forecast,error,covariance),and,our,current,(noisy,and,partial),observations,of,the,system,state.,Two,main,variants,exist,-,3D-Var,,which,finds,a,system,state,which,minimizes,a,cost,function,including,terms,only,for,observations,of,the,initial,state,,and,4D-Var,,which,includes,terms,for,obser-,vations,of,both,the,initial,state,and,the,future,states,within,some,assimilation,window.,In,both,cases,,for,non-linear,observation,models,,computation,of,the,gradient,of,the,cost,function,require,being,able,to,compute,the,action,of,the,adjoint,of,the,Jacobian,of,the,observation,operator,(sometimes,referred,to,as,the,tangent,linear,model,),,and,in,the,4D-Var,case,,when,the,state,transition,operator,is,non-linear,we,also,require,being,able,to,compute,the,action,of,the,adjoint,of,the,Jacobian,of,the,state,transition,operator.,Variational,methods,are,widely,used,in,practice,in,operational,numerical,weather,prediction,sys-,tems,,though,as,they,produce,only,point,estimates,,they,typically,require,separate,approaches,for,uncertainty,quantification.,The,standard,‘strong,constraint’,formulations,of,variational,methods,also,do,not,account,for,error,in,the,model’s,predictions,of,how,the,state,evolves,over,time,,effectively,assum-,ing,the,model,is,perfect,over,the,assimilation,window.,Weak,constraint,4D-Var,(Evensen,,Vossepoel,,and,Leeuwen,2022),relaxes,the,perfect,model,assumption,,at,a,cost,of,increasing,the,dimension,of,the,space,being,optimized,over,,with,the,approach,jointly,optimizing,both,an,initial,state,and,values,for,model,error,‘forcing,terms’.,Particle,filters,(Gordon,,Salmond,,and,Smith,1993),offer,an,alternative,to,both,variational,methods,and,ensemble,Kalman,filter.,Similarly,to,ensemble,Kalman,filters,,particle,filters,use,an,ensemble,of,state,particles,to,represent,the,estimates,of,the,sequence,of,filtering,distributions.,Rather,than,form,a,Monte,Carlo,approximation,to,the,Kalman,filter,updates,for,linear-Gaussian,state,space,models,,the,particle,filter,instead,use,a,(importance,sampling),Monte,Carlo,approximation,to,the,recursions,updating,the,filtering,distribution,at,one,time,point,to,the,next,for,general,state,space,models,,and,so,in,the,limit,of,the,ensemble,size,N,give,consistent,estimates,of,the,true,filtering,distributions,for,arbitrary,state,space,models,with,potentially,non-linear,dynamics,and,observation,operators,or,non-Gaussian,noise,distributions.,→,∞,Particle,filters,iteratively,alternate,between,a,proposal,step,in,which,the,new,values,for,each,par-,ticle,in,the,ensemble,are,independently,generated,from,a,proposal,distribution,and,a,resampling,step,in,which,the,particles,are,resampled,according,to,a,set,of,importance,weights,,with,highly,weighted,particles,being,duplicated,and,low,weighted,particles,being,discarded.,The,proposal,step,and,compu-,tation,of,(unnormalized),importance,weights,can,be,parallelized,across,all,particles,,and,for,complex,state,space,models,will,typically,dominate,the,computational,cost,of,filtering,,making,particle,filters,ideally,suited,to,deployment,on,highly,parallel,hpc,systems,,and,so,for,performing,data,assimila-,tion,at,the,exascale.,While,normalization,of,the,importance,weights,and,the,resampling,step,do,require,synchronization,across,particles,,the,operation,cost,remains,linear,in,the,ensemble,size,N,and,only,computation,of,the,(scalar),normalized,importance,weights,requires,collective,operations,across,all,processors,,with,redistribution,of,the,particles,during,resampling,able,to,be,done,using,sparse,point-to-point,transfers.,In,comparison,,for,the,ensemble,Kalman,filter,the,operation,costs,scale,as,yN,+,dxdyN,),otherwise,,and,the,linear,algebra,oper-,O,ations,in,the,filtering,updates,require,collective,operations,on,the,full,state,vectors,(particles),when,distributing,across,multiple,processors,,resulting,in,a,much,higher,communication,overhead.,(N,3,+,dyN,2,+,dxN,2),if,dy,>,N,and,y,+,d2,(d3,O,Although,particle,filters,offer,promise,in,terms,of,their,ability,to,give,consistent,estimates,for,non-linear,and,non-Gaussian,state,space,models,,and,ease,of,scaling,on,hpc,systems,,a,key,drawback,of,particle,filters,is,that,they,can,require,an,ensemble,size,which,scales,exponentially,with,the,dimen-,24,sionality,of,the,model,(Snyder,2011;,Fearnhead,and,K¨unsch,2018),,to,avoid,a,degeneracy,sometimes,terms,ensemble,collapse,,by,which,all,members,of,the,ensemble,except,one,are,assigned,negligible,weights.,While,use,of,improved,proposal,distributions,can,help,combat,this,issue,to,a,degree,(Snyder,2011;,Slivinski,and,Snyder,2016),,algorithmic,extensions,such,as,tempering,(Frei,and,K¨unsch,2013;,Beskos,,Crisan,,and,Jasra,2014;,Bunch,and,Godsill,2016;,Beskos,,Crisan,,Jasra,,et,al.,2017),or,spatial,localization,(Rebeschini,and,Handel,2015;,Farchi,and,Bocquet,2018;,Graham,and,Thiery,2019),,are,likely,to,be,required,for,statistically,robust,use,of,particle,filters,for,data,assimilation,in,state,space,models,of,very,high,dimension,,such,as,for,the,numerical,models,of,plasma,physics,of,interest,in,the,NEPTUNE,project.,10.1,State,space,models,Particle,(and,ensemble,Kalman),filters,assume,the,system,of,interest,is,represented,as,a,state,space,model.,Specifically,we,assume,we,have,a,sequence,of,observations,y1:T,=,(yt)T,t=1,of,the,system,at,Rdy,.,The,state,of,the,T,>,0,times,with,each,observation,a,dy,dimensional,real-valued,vector,,yt,∈,system,at,each,time,index,t,1:T,is,also,assumed,to,be,representable,by,a,dx,dimensional,real-valued,Rdx,.,vector,xt,∈,The,state,dynamics,are,assumed,to,be,Markovian,,with,the,state,at,each,time,index,depending,only,on,the,previous,state,,with,initial,and,transition,distributions,specified,by,densities,p0,and,p1:t,respectively,,∈,x0,p0(,);,·,∼,xt,pt(,·,|,∼,xt,1),−,1:T,;,t,∈,and,the,observations,at,each,time,index,t,are,assumed,to,depend,only,on,the,current,state,,with,condition,distributions,specified,by,densities,g1:t,yt,∼,gt(,·,|,xt),t,1:T.,∈,Here,the,notation,z,distribution,with,density,q.,q(,·,∼,),denotes,that,a,random,vector,z,is,distributed,according,to,a,probability,A,key,point,here,is,that,the,state,transitions,are,assumed,to,be,stochastic.,This,may,seem,initially,restrictive,for,the,setting,of,interest,here,,simulation,of,plasma,physics,,where,typically,models,are,specified,as,time-dependent,systems,of,partial,differential,equations,(pdes),,with,the,state,evolution,in,such,models,inherently,deterministic.,In,any,real,setting,however,the,model,will,only,be,a,partial,description,of,the,underlying,phenomena,being,modelled,,and,if,we,do,not,account,for,this,by,specifying,we,have,some,uncertainty,in,how,the,state,evolves,over,time,compared,to,in,the,true,system,,we,will,struggle,to,get,useful,outputs,from,any,data,assimilation,algorithm.,As,a,simple,example,we,could,consider,trying,to,fit,a,simple,harmonic,motion,(shm),model,to,observations,of,a,pendulum,state.,Under,this,deterministic,model,the,pendulum,state,at,all,future,(and,past),times,would,be,entirely,determined,by,its,initial,(angular),position,and,momentum,,so,the,problem,reduces,down,to,finding,an,initial,state,such,that,the,simulated,trajectories,are,consistent,with,the,observations,,with,the,simulated,trajectories,being,sinusoids,of,constant,frequency,and,am-,plitude.,In,practice,the,real,world,observed,pendulum,state,trajectories,will,not,be,exactly,sinusoidal,–,depending,on,the,amplitude,of,the,oscillations,the,neglected,non-linear,effects,in,the,small,angle,approximation,underlying,shm,may,be,significant,,and,probably,more,importantly,,frictional,effects,not,accounted,for,in,the,shm,model,will,always,mean,the,oscillations,will,decrease,in,amplitude,over,time.,When,we,try,to,fit,the,model,to,the,observed,data,we,will,therefore,not,be,able,to,find,some,initial,state,for,the,system,which,leads,to,trajectories,consistent,with,all,the,observed,data,,with,we,likely,to,instead,only,be,able,to,find,states,which,are,consistent,with,short,contiguous,periods,of,the,observations,for,which,the,shm,assumptions,roughly,hold.,If,we,instead,fit,,for,example,,a,stochastic,differential,equation,model,to,the,observations,which,includes,terms,both,describing,a,deterministic,evolution,according,to,our,shm,model,and,a,white,noise,term,representing,unmodelled,aspects,of,the,dynamics,like,frictional,effects,(and,when,discretized,in,time,can,be,formulated,as,a,state,space,model),,then,the,additional,flexibility,in,the,model,would,allow,us,to,find,state,sequences,that,are,consistent,with,all,the,observational,data,,while,still,allowing,us,to,make,useful,inferences,about,variables,of,the,deterministic,model,which,will,provide,a,naturally,25,parsimonious,description,of,key,aspects,of,the,dynamics,(near,constant,frequency,oscillations),with,the,noise,model,then,only,needing,to,account,for,the,features,not,described,by,our,model,(decay-,ing,amplitude,,slightly,non-sinusoidal,oscillations,,slow,variations,in,frequency).,If,we,increase,the,complexity,of,the,deterministic,component,of,model,to,include,additional,aspects,of,the,system,–,for,example,including,a,damping,term,to,partially,account,for,frictional,effects,–,then,more,of,the,observed,behaviour,will,be,able,to,explained,by,the,deterministic,component,of,the,model,,with,the,stochastic,component,therefore,needing,to,account,for,less,unmodelled,aspects,,potentially,reflected,for,example,by,we,inferring,a,smaller,scale,parameter,for,the,injected,white,noise,process.,These,considerations,applies,equally,to,plasma,physics,,with,models,likely,to,differ,from,the,true,systems,of,interest,both,due,to,unmodelled,physics,and,limitations,in,the,ability,of,the,numerical,methods,used,to,simulate,models,to,resolve,behaviour,at,all,temporal,and,spatial,scales.,In,such,models,,a,key,question,is,how,to,represent,our,uncertainties,about,the,evolution,of,the,system,state,,and,importantly,ensure,that,the,perturbations,we,introduce,to,model,our,uncertainty,do,not,disrupt,the,key,underlying,physical,behaviours,we,are,seeking,to,model,,including,any,physical,constraints,on,the,system,state,such,as,boundary,conditions.,10.2,Particle,filters,The,objective,of,the,particle,filter,is,to,estimate,the,filtering,distribution,for,each,time,index,t,which,is,the,conditional,probability,distribution,of,the,state,xt,given,observations,y1:t,up,to,time,index,t,,with,the,the,filtering,distribution,at,time,index,t,denoted,πt.,As,a,slight,abuse,of,notation,we,denote,distributions,and,their,densities,(with,respect,to,the,Lebesgue,measure,unless,otherwise,noted),with,the,same,symbol,,so,that,πt,denotes,both,the,filtering,distribution,at,time,index,t,and,its,density.,An,ensemble,of,particles,(x(n),)N,n=1,represents,an,approximation,to,the,filtering,distribution,at,(cid:80)N,each,time,index,t,),,where,δx,is,a,Dirac,·,measure,at,x.,For,the,proposal,step,at,each,time,index,t,,new,values,for,the,particles,are,sampled,from,a,proposal,distribution,with,density,qt,given,the,particle,values,at,the,previous,time,index,and,current,observations,1,as,an,empirical,distribution,πt(,n=1,δ,y1:t),x(n),t,1,N,·,|,≥,≈,(,t,·,|,and,(unnormalized),importance,weights,computed,for,each,proposed,particle,value,∈,˜x(n),t,∼,qt(,x(n),t,−,1,,yt),,n,1:N,;,1),gt(yt,|,x(n),1,,yt),t,|,−,Both,of,the,sampling,of,proposals,and,computation,of,the,unnormalized,importance,weights,can,be,performed,independently,(and,so,in,parallel),for,each,particle.,pt(˜x(n),x(n),t,t,|,−,qt(˜x(n),t,˜x(n),t,t,=,w(n),,,n,1:N.,∈,),In,the,resampling,step,,the,weighted,proposed,particle,ensemble,is,resampled,from,the,weighted,empirical,distribution,x(n),t,∼,(cid:80)N,n=1,w(n),t,δ,˜x(n),t,n′=1,w(n′),(cid:80)N,t,(,),·,,,n,1:N,;,∈,to,produce,a,new,uniformly,weighted,ensemble,to,use,as,the,input,to,the,next,filtering,step.,This,step,requires,synchronization,across,the,particles.,A,particular,simple,choice,of,proposal,distribution,,is,the,na¨ıve,‘bootstrap’,proposals,which,ignore,the,observations,and,sample,the,proposals,from,the,state,transition,distributions,qt(xt,|,xt,1,,yt),=,pt(xt,−,1),,xt,−,|,1:T,;,t,∈,),,n,with,the,expression,for,the,unnormalized,importance,weights,then,simplifying,to,w(n),|,˜x(n),1:N,.,The,bootstrap,proposals,can,be,straightforwardly,implemented,any,state,space,model,t,for,which,we,can,sample,from,the,state,transition,distributions,and,evaluate,the,observation,densities,,however,the,resulting,particle,filter,can,perform,poorly,when,the,observations,are,informative,about,t,=,gt(yt,∈,26,the,state,,with,the,proposed,particles,typically,ending,up,far,from,the,regions,containing,the,most,of,the,mass,of,the,trued,filtering,distributions,,resulting,in,importance,weights,with,high,variance,,and,weight,degeneracy,with,all,but,one,normalized,importance,weight,close,to,zero.,An,alternative,to,the,bootstrap,proposals,which,can,reduce,the,tendency,for,weight,degeneracy,and,ensemble,collapse,,is,to,use,the,locally,optimal,proposal,distribution,(Doucet,,Godsill,,and,Andrieu,2000),qt(xt,|,xt,1,,yt),=,−,pt(xt,(cid:82),pt(x,xt,xt,−,1),gt(yt,|,1),gt(yt,|,−,|,|,xt),x),dx,,,1:T,;,t,∈,t,=,(cid:82),pt(x,which,minimizes,the,variance,of,the,importance,weights,,with,the,importance,weights,in,this,case,simplifying,to,w(n),1:N,which,are,independent,of,the,sampled,x),dx,,n,values,of,the,particle,proposals,˜x(1:N,),.,Sampling,from,this,locally,optimal,proposal,distribution,and,computing,the,corresponding,importance,weights,,is,intractable,in,general,as,the,integral,in,the,definition,of,the,proposal,distribution,does,not,usually,have,a,closed,form,solution.,1),gt(yt,|,x(n),t,−,t,∈,|,However,,for,an,important,subclass,of,state,space,models,with,additive,Gaussian,state,and,obser-,vation,noise,and,linear,observation,operators,,that,is,1),,Qt),,yt,∼,the,locally,optimal,proposals,have,the,tractable,form,Normal(,Ft(xt,xt,·,|,∼,−,Normal(,Htxt,,Rt),,·,|,1:T,;,t,∈,qt(xt,xt,−,|,1,,yt),=,Normal,(xt,|,mt(xt,1,,yt),,Ct),,,−,with,and,mt(xt,1,,yt),=,Ft(xt,−,1),+,QtH,T,t,(HtQtH,T,t,+,Rt)−,1(yt,−,HtFt(xt,−,1)),,−,and,corresponding,importance,weights,Ct,=,Qt,−,QtH,T,t,(HtQtH,T,t,+,Rt)−,1HtQt;,w(n),t,=,Normal(yt,|,HtFt(xt,−,1),,HtQtH,T,t,+,Rt).,10.3,ParticleDA,ParticleDA,(Giles,et,al.,2023),is,a,Julia,package,implementing,particle,filtering,algorithms,for,data,assimilation,,with,a,focus,on,allowing,the,particle,filtering,updates,to,be,computed,at,scale,on,hpc,systems,using,both,thread-based,parallelism,and,distributed,processing,using,a,message,passing,in-,terface,(mpi),implementation.,It,was,developed,initially,as,part,of,a,project,real-time,advanced,data,assimilation,for,digital,simulation,of,numerical,twins,on,HPC,(raddish),,with,an,initial,focus,on,geoscientific,applications,such,as,tsunami,modelling,and,numerical,weather,prediction.,As,well,as,its,support,for,both,thread-based,and,mpi,based,parallelism,,it,provides,implementations,of,particle,filters,using,both,bootstrap,and,locally,optimal,proposal,distributions.,As,part,of,UCL’s,work,on,AQUIFER,,we,made,various,improvements,to,ParticleDA,,both,to,extend,the,range,of,models,it,can,be,used,with,and,improve,performance:,•,The,locally,optimal,proposal,implementation,was,generalised,to,allow,use,with,a,more,general,class,of,state,space,models,(https://github.com/Team-RADDISH/ParticleDA.jl/pull/218).,•,ParticleDA’s,model,and,filter,interfaces,were,refactored,to,remove,the,previous,implicit,as-,sumption,of,models,being,finite,difference,based,spatial,discretization,of,a,time-dependent,pde,,with,the,state,vector,correspond,to,the,values,of,the,state,field(s),at,a,regular,grid,of,points,(https://github.com/Team-RADDISH/ParticleDA.jl/pull/232).,Post,refactoring,ParticleDA,now,requires,only,that,the,model’s,state,can,be,represented,as,a,flat,vector,of,fixed,dimension,,allowing,ParticleDA,to,be,applied,to,non-spatial,models,such,as,systems,of,odes,or,to,time-dependent,pde,models,using,for,example,finite,element,based,spatial,discretizations,on,potentially,arbitrarily,structured,spatial,meshes.,27,•,The,parallelization,of,the,particle,proposal,steps,across,multiple,threads,was,changed,from,a,static,scheme,in,which,groups,of,particles,were,pre-assigned,to,each,thread,,to,a,scheme,in,which,the,particles,updates,are,chunked,across,a,user,specifiable,number,of,tasks,,with,these,tasks,then,being,scheduled,across,threads,dynamically,,with,this,improving,load,balancing,and,reducing,the,chance,of,threads,sitting,idle,(https://github.com/Team-RADDISH/ParticleDA.jl/pull/,247).,10.4,ANAET,state,space,model,As,an,initial,test,case,illustrating,how,ParticleDA,can,be,used,to,perform,data,assimilation,,we,consider,a,toy,model,of,magnetohydrodynamic,instability,coupled,to,a,slow,dissipative,background,evolution,(Arter,2012).,The,resulting,axissymmetric,non-axissymmetric,extended,(anaet),model,is,described,by,the,system,of,odes,¨a,=,γra,−,−,(µ1,+,µ2b)a3,µ6a6,˙a,,−,˙b,=,ν1,ν2b2,−,−,(δ0,+,δ1b)a2,,where,a,is,a,non-axissymmetric,ideal,mode,coupled,to,an,axissymmetric,mode,b,,and,(γr,,mu1,,µ2,,µ6,,ν1,,ν2,,δ0,,δ1),are,free,parameters.,Defining,x,=,(a,,˙a,,b),,these,equations,can,be,written,as,a,first,order,ode,system,,=,,dx(τ,),dτ,−,x2(τ,),,γrx1(τ,),(µ1,+,µ2x3(τ,))x1(τ,)3,µ6x1(t)6x2(τ,),,.,−,ν1,ν2x3(τ,)2,(δ0,+,δ1b)x1(τ,)2,−,−,Here,,to,formulate,as,a,state,space,model,,we,consider,a,simple,stochastic,evolution,variant,of,this,−,model,defined,by,the,stochastic,differential,equation,(sde),system,,,dx(τ,),=,−,γrx1(τ,),(µ1,+,µ2x3(τ,))x1(τ,)3,µ6x1(τ,)6x2(τ,),,dt,+,βdw(τ,),,x2(τ,),−,−,ν1,−,ν2x3(τ,)2,(δ0,+,δ1b)x1(τ,)2,−,,with,w(t),a,Wiener,process,,with,the,magnitude,of,the,noise,introduced,into,the,state,dynamics,controlled,by,a,shared,scalar,parameter,β,>,0.,We,assume,only,observations,of,first,(a),state,component,and,independent,Gaussian,observation,noise,,that,is,an,observation,model,yt,∼,Normal(,·,|,x1(τt),,σ2),,1:T.,t,∈,Our,state,space,model,formulation,uses,a,Gaussian,approximation,to,the,state,transition,distributions,for,the,sde,system,based,on,a,splitting,numerical,scheme,which,uses,an,adaptive,ode,solver,to,solve,for,the,deterministic,component,of,the,dynamics,and,Euler–Maruyama,discretization,for,the,driving,Wiener,noise,processes.,The,code,for,the,Julia,implementation,of,the,model,,along,with,installation,instructions,and,exam-,ple,usages,with,ParticleDA,is,available,at,https://github.com/Team-RADDISH/ParticleDA-UseCases/,tree/main/ANAET.,As,an,example,,we,will,use,a,particle,filter,to,estimate,the,state,trajectories,x1:T,with,xt,=,x(τt),,τt,=,t,and,T,=,100,,given,simulated,observed,data,y1:T,generated,from,the,model.,The,model,parameters,are,fixed,to,the,values,γr,=,1,,µ1,=,0,,µ2,=,2,,µ6,=,0.001,,ν1,=,0.005,,ν2,=,0.0001,,δ0,=,0.0001,,δ1,=,0,,σ,=,0.5,and,β,=,0.02.,We,generate,a,set,of,simulated,observed,data,from,the,model,using,a,fixed,initial,state,x0,=,(1,,2.5,,0.01),when,simulating,the,observations,,but,when,filtering,use,an,initial,state,distribution,with,x0,(0,,2,,0),,I),,where,I,is,the,identity,matrix.,·,|,We,use,ParticleDA,to,compute,particle,approximations,to,the,filtering,distributions,π1:T,,,using,its,implementations,of,both,the,bootstrap,and,locally,optimal,proposals.,Figures,26,,28,and,30,show,results,for,the,locally,optimal,proposals,with,N,=,125,,N,=,250,and,N,=,500,particles,respectively,,and,Figures,27,,29,and,31,show,results,for,the,bootstrap,proposals,with,N,=,125,,N,=,250,and,N,=,500,particles,respectively.,In,all,cases,,for,the,filtering,estimates,,the,blue,line,shows,the,Normal(,−,∼,28,Figure,26:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,125,and,locally,optimal,proposals.,estimated,filtering,distribution,mean,and,the,blue,filled,region,the,estimated,mean,standard,deviation,interval.,±,three,estimated,From,Figures,30,and,31,see,that,particle,filters,using,both,the,locally,optimal,and,bootstrap,proposals,appear,to,give,reasonable,estimates,of,the,filtering,distribution,with,N,=,500,particles,,with,the,true,state,trajectories,generally,within,the,intervals,estimated,to,contain,most,of,the,mass,of,the,filtering,distributions,as,expected.,For,the,smaller,ensembles,of,N,=,250,(Figures,28,and,29),and,N,=,125,(Figures,26,and,27),,we,see,that,the,estimates,of,filters,using,both,proposals,appear,to,degrade,at,some,points,(particularly,at,earlier,times,in,the,trajectories,where,the,is,greater,uncertainty,due,to,fewer,observations),,but,that,filters,using,the,locally,optimal,proposals,outperform,filters,using,the,bootstrap,proposals.,10.5,NektarDriftwave,state,space,model,As,a,second,more,complex,test,model,,we,consider,integrating,ParticleDA,with,a,Nektar++,solver,for,a,plasma,physics,relevant,time-dependent,pde,model.,Specifically,we,define,a,state,space,model,which,extends,the,nektar-driftwave,proxyapp,,which,uses,Nektar++,to,solve,the,two-dimensional,Hasegawa–Wakatani,equations,,∂3n(s1,,s2,,τ,),+,(s1,,s2,,τ,),=,α(ϕ,∂3ζ(s1,,s2,,τ,),+,ϕ,,ζ,},n)(s1,,s2,,τ,),(s1,,s2,,τ,),=,α(ϕ,ϕ,,n,},{,−,1,+,∂2,(∂2,n)(s1,,s2,,τ,),,−,κ∂2ϕ(s1,,s2,,τ,),,2,)ϕ(s1,,s2,,τ,),=,ζ(s1,,s2,,τ,);,−,{,where,the,Poisson,bracket,operator,is,defined,as,(s1,,s2,,τ,),=,∂1a(s1,,s2,,τ,)∂2b(s1,,s2,,τ,),a,,b,},{,−,∂2a(s1,,s2,,τ,)∂1b(s1,,s2,,τ,).,The,two,dimensional,spatial,domain,is,assumed,to,be,rectangular,with,periodic,boundary,conditions.,Nektar++,is,used,to,spatially,discretize,the,system,using,a,spectral,element,method,,with,here,a,29,Figure,27:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,125,and,bootstrap,proposals.,Figure,28:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,250,and,locally,optimal,proposals.,30,Figure,29:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,250,and,bootstrap,proposals.,Figure,30:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,500,and,locally,optimal,proposals.,31,Figure,31:,True,state,trajectories,(orange),,observed,data,(green,markers),and,filtering,estimates,(blue),for,stochastic,anaet,model,using,particle,filter,with,P,=,500,and,bootstrap,proposals.,regular,quadrilateral,mesh,used.,The,time,discretization,is,performed,using,an,explicit,fourth-order,Runge–Kutta,method,,with,fixed,time,step.,To,formulate,as,a,state,space,model,,we,define,the,system,state,as,the,concatenation,of,the,coefficients,of,the,spectral,element,expansions,of,the,fields,ζ,(vorticity),and,n,(number,density),,with,the,electrostatic,potential,ϕ,being,fully,determined,by,ζ.,The,state,dimension,is,therefore,equal,to,dx,=,2m1m2(p,+,1)2,where,m1,is,the,number,of,quadrilateral,elements,along,the,s1,spatial,coordinate,,m2,is,the,number,of,quadrilateral,elements,along,the,s2,spatial,coordinate,and,p,is,the,polynomial,order,of,the,spectral,element,expansion.,To,introduce,stochasticity,in,the,state,dynamics,,after,using,the,nektar-driftwave,solver,to,simulate,the,Hasegawa–Wakatani,equations,forward,in,time,by,a,fixed,inter-observation,time,interval,∆,=,τt,1,,we,perturb,the,ζ,and,n,state,fields,by,adding,zero-mean,stationary,Gaussian,random,fields,with,Mat´ern,covariance,function,τt,−,−,C(s1,,s2,|,σ,,ℓ,,ν),=,σ2,21,ν,−,Γ(ν),(cid:33)ν,(cid:32),(cid:112)s2,1,+,s2,2,ℓ,(cid:33),(cid:32),(cid:112)s2,1,+,s2,2,ℓ,Kν,where,σ2,is,the,marginal,variance,,ℓ,a,lengthscale,parameter,and,ν,a,smoothness,parameter,,Γ,the,gamma,function,and,Kν,the,modified,Bessel,function,of,the,second,kind,of,order,ν.,By,simulating,Gaussian,random,fields,with,appropriate,choices,of,the,covariance,parameters,(σ2,,ℓ,,ν),on,the,spatial,domain,of,the,problem,and,with,same,boundary,conditions,as,the,state,field,variables,ζ,and,n,,we,can,produce,perturbations,of,the,state,fields,which,conserve,key,aspects,of,the,fields,such,as,their,smoothness,properties,and,boundary,conditions.,To,simulate,such,Gaussian,random,fields,,we,use,that,they,correspond,to,(scaled),solutions,σ√4πνℓ−,2νu,to,the,stochastic,partial,differential,equation,(spde),2,(ℓ−,∂2,1,−,−,∂2,2,),ν+1,2,u(s1,,s2),=,(s1,,s2),,W,where,spatial,statistics,by,(Lindgren,,Rue,,and,Lindstr¨om,2011).,For,ν,=,2m,is,a,Gaussian,white,noise,process,,as,first,observed,by,(Whittle,1954),and,popularized,in,N,we,additionally,1,,m,W,−,∈,32,Figure,32:,Examples,of,Gaussian,random,fields,Mat´ern,covariance,generated,using,Nektar++,by,solving,a,spde.,All,fields,are,on,a,40,40,spatial,domain,with,periodic,boundary,conditions.,From,left,to,right,the,covariance,function,lengthscale,and,smoothness,parameters,used,where,(ℓ,=,1,,ν,=,1),,(ℓ,=,1,,ν,=,3),and,(ℓ,=,2,,ν,=,3).,×,have,that,u,=,v(m),where,v(m),is,recursively,defined,by,2,(ℓ−,∂2,1,−,−,2,)v(m)(s1,,s2),=,v(m,∂2,1)(s1,,s2),,−,W,with,base,case,v(0),=,.,To,simulate,Gaussian,random,fields,with,Mat´ern,covariance,on,a,spatial,domain,in,Nektar++,,we,can,therefore,iteratively,solve,a,Helmholtz,equation,using,the,advection–,diffusion–reaction,solver,,with,forcing,(righthand,side),function,set,to,Gaussian,white,noise,for,the,base,case,,and,the,field,output,of,the,previous,iterations,for,subsequent,iterations.,We,specify,the,same,mesh,and,boundary,conditions,as,for,the,fields,used,in,the,nektar-driftwave,solver,for,the,deterministic,component,of,the,dynamics.,Figure,32,shows,examples,of,Gaussian,random,fields,generated,using,different,values,for,the,lengthscale,ℓ,and,smoothness,ν,parameters.,The,observed,data,at,each,observation,time,τt,=,t∆,is,then,assumed,to,correspond,to,noisy,observations,of,the,ζ,field,at,a,set,of,spatial,points,(s(i),1,,,s(i),2,),Normal,yt,∼,(cid:18),(cid:16),·,|,ζ(s(i),1,,,s(i),2,,,t∆),dy,i=1,,that,is,(cid:19),(cid:17)dy,,,σ2,yI,.,i=1,The,code,for,the,Julia,implementation,of,the,model,,along,with,installation,instructions,and,example,usages,with,ParticleDA,is,available,at,https://github.com/Team-RADDISH/ParticleDA-UseCases/,tree/main/NektarDriftwave.,To,test,the,scaling,of,the,ParticleDA,filter,implementation,on,a,hpc,system,,we,ran,multiple,filtering,runs,with,differing,combinations,of,thread,and,mpi,based,parallelism,on,ARCHER2,,the,UK,national,supercomputing,service.,ARCHER2,provides,a,cluster,of,5860,compute,nodes,each,with,dual,64-core,processors,,and,so,128,cores,per,node.,×,10,,[,−,20,,20],We,first,simulated,a,fixed,sequence,of,T,=,100,observations,from,the,NektarDriftwave,state,space,model,,using,a,quadrilateral,mesh,of,dimensions,m1,=,32,,m2,=,32,on,a,spatial,domain,20,,20],with,periodic,boundary,conditions.,The,polynomial,order,used,in,the,spectral,[,−,expansion,was,p,=,3,,resulting,in,an,overall,state,dimension,of,dx,=,32768.,The,ζ,field,was,observed,at,dy,=,9,points,(,10,,10),,(0,,10),,(10,,10),with,observation,noise,standard,deviation,σy,=,0.1.,The,zero-mean,independent,Gaussian,random,fields,used,to,perturb,the,ζ,and,n,state,fields,between,each,observation,time,used,a,Mat´ern,covariance,function,with,parameters,ν,=,3,,ℓ,=,1,and,σ,=,0.05,,with,a,fixed,inter-observation,interval,of,∆,=,0.5,and,a,time,integration,step,size,for,the,nektar-driftwave,solver,of,0.0005.,The,adiabacity,operator,and,background,density,gradient,lengthscale,parameters,of,the,Hasegawa–Wakatani,equations,were,set,to,respectively,α,=,2,and,κ,=,1.,The,state,was,initialized,with,independent,Gaussian,random,fields,for,ζ,and,n,with,the,same,covariance,function,and,parameters,as,the,state,noise,and,means,10,,0),,(0,,0),,(10,,0),,(,10),,(10,,10),,(0,,10),,(,−,−,−,−,−,−,E[ζ(s1,,s2,,0)],=,exp,(cid:0),1,+,s2,(s2,E[n(s1,,s2,,0)],=,exp,(cid:0),2)/4(cid:1),(s2,−,1,+,s2,(s2,2,−,1,+,s2,4)/4,,2)/4(cid:1),.,−,33,Figure,33:,Weak,scaling,efficiency,for,thread-based,parallelization,of,particle,filter,based,data,assim-,ilation,using,ParticleDA,with,NektarDriftwave,state,space,model,on,a,single,ARCHER2,node,We,first,consider,how,using,multiple,threads,to,parallelize,the,particle,updates,on,a,single,ARCHER2,node,performs,,running,a,particle,filter,with,bootstrap,proposals,and,one,particle,per,thread.,We,fil-,ter,using,the,sequence,of,T,=,100,observations,simulated,as,described,above,,using,the,same,model,parameters,for,filtering,as,used,for,generating,the,observations.,Figure,33,shows,the,weak,scaling,efficiency,for,varying,number,of,threads,,with,weak,scaling,measuring,parallel,performance,as,both,the,compute,resource,(number,of,threads,here),and,problem,size,(number,of,particles,in,ensemble),are,increased.,Here,we,specifically,defined,the,weak,scaling,efficiency,as,the,(compute,time,on,for,2,particles,on,2,threads),/,(compute,time,for,N,particles,on,N,threads).,We,can,see,that,the,weak,scaling,efficiency,shows,a,slow,but,consistent,drop,off,as,the,number,of,threads,increases,,with,an,efficiency,of,0.96,at,16,threads,and,0.67,at,128,threads.,For,reference,the,compute,time,with,2,threads,was,6590s.,As,each,compute,node,on,ARCHER2,has,128,cores,,if,the,application,was,compute,bound,we,would,expect,to,see,the,efficiency,remain,close,to,unity,up,to,128,threads.,The,drop-off,here,suggests,the,performance,bottleneck,on,a,single,node,is,not,processor,availability.,The,NektarDriftwave,state,space,model,implementation,uses,file-based,input,/,output,to,interact,with,the,Nektar++,solvers,used,to,solve,both,the,deterministic,dynamics,and,generate,the,Gaussian,random,fields,used,to,perturb,the,states,,resulting,in,multiple,writes,and,reads,of,the,state,to,and,from,hdf5,files,on,disk,for,each,particle,and,at,each,filtering,step.,As,all,threads,on,a,single,node,will,be,competing,for,the,same,bandwidth,for,reading,from,and,writing,to,disk,,we,believe,this,may,be,the,performance,bottleneck,in,this,case.,To,avoid,this,bottleneck,,it,would,be,preferable,to,share,memory,for,storing,the,states,between,the,ParticleDA,Julia,code,and,Nektar++,C++,code,to,avoid,the,round-tripping,to,and,from,disk,,however,we,did,not,have,time,to,explore,this,option.,ParticleDA,also,supports,distributing,the,particle,updates,using,mpi,,allowing,us,to,test,scaling,when,using,multiple,ARCHER2,nodes.,We,first,again,consider,weak,scaling,efficiency,,using,128,particles,and,mpi,ranks,per,node,(one,per,core).,Figure,34,show,the,weak,scaling,efficiency,for,varying,number,of,nodes,,here,defined,as,the,(compute,time,for,128,particles,on,1,node),/,(compute,time,for,128n,particles,on,n,nodes).,Unlike,the,multi-threading,case,,the,weak,scaling,efficiency,remains,close,to,unity,as,we,increase,the,number,of,nodes,,dropping,only,to,0.96,at,32,nodes,,the,maximum,tested,here,,with,again,for,reference,the,compute,time,on,a,single,node,here,8890s.,We,speculate,that,the,improved,scaling,efficiency,observed,here,compared,to,the,multi-threading,case,is,due,to,the,distribution,of,the,computation,across,multiple,nodes,increasing,both,the,number,34,128643216842Number,of,threads0.00.20.40.60.81.0Weak,scaling,efficiency,Figure,34:,Weak,scaling,efficiency,for,mpi-based,parallelization,of,particle,filter,based,data,assimilation,using,ParticleDA,with,NektarDriftwave,state,space,model,on,ARCHER2.,of,processors,available,but,also,bandwidth,for,reading,and,writing,to,disk.,We,also,note,that,the,multi-node,weak,scaling,performance,observed,here,is,significantly,better,than,previously,shown,in,(Giles,et,al.,2023),for,filtering,runs,with,another,state,space,model,(a,two-dimensional,linear,long,wave,model,of,tsunamis),on,ARCHER2,,where,the,weak,scaling,efficiency,dropped,off,sharply,when,distributing,across,multiple,nodes.,This,can,be,explained,by,the,relative,cost,of,the,proposal,and,resampling,steps,between,the,two,models,,with,the,former,(embarrassingly,parallelizable),step,dominating,the,run,time,for,the,model,here,(the,proposal,steps,constituted,96.5%,of,the,run,time,on,one,node,and,94.1%,on,32,nodes).,In,comparison,the,state,updates,in,the,tsunami,model,were,much,cheaper,to,compute,and,the,total,ensemble,sizes,used,much,larger,,resulting,in,a,much,higher,proportion,of,time,being,spent,in,the,resampling,step,,which,due,to,the,requirement,for,communication,across,ranks,becomes,more,costly,as,the,number,of,nodes,(and,ranks),is,increased.,As,the,plasma,physics,models,of,interest,to,neptune,will,be,better,reflected,by,the,computational,demands,of,the,NektarDriftwave,example,here,,we,believe,the,regime,of,interest,will,generally,be,one,were,the,cost,of,the,model,simulation,in,the,proposal,steps,dominate,,and,the,results,here,indicate,ParticleDA,should,be,able,to,achieve,high,scaling,efficiency,in,this,setting.,As,a,final,performance,benchmark,,we,also,measured,the,strong,scaling,performance,of,particle,filtering,with,ParticleDA,on,the,NektarDriftwave,model,,which,measures,the,speed-up,achieved,as,the,problem,size,(number,of,particles),is,held,constant,and,the,compute,resource,(number,of,processors),increased.,We,ran,particle,filters,(with,bootstrap,proposals),on,the,same,T,=,100,simulated,observa-,tion,sequence,as,previously,with,a,fixed,ensemble,size,of,N,=,1024,particles,,distributing,the,filtering,updates,using,mpi,across,one,or,more,nodes,with,128,ranks,per,node,(one,per,core).,Figure,35,shows,the,strong,scaling,speed-ups,(compute,time,for,N,=,1024,particles,on,n,nodes,/,compute,time,for,N,=,1024,particles,on,1,node),measured,for,runs,on,up,to,8,nodes.,We,see,that,we,achieve,close,to,ideal,of,a,linear,speed-up,as,we,increase,the,the,number,of,nodes,,suggesting,ParticleDA,is,also,able,to,achieve,good,strong,scaling,for,this,model.,The,previous,results,related,to,the,computational,performance,and,scalability,of,the,particle,filter,algorithms,implemented,in,ParticleDA.,Another,import,metric,is,how,they,perform,statistically,,in,terms,of,the,accuracy,of,their,estimation,of,the,true,filtering,distributions.,As,we,have,no,way,of,computing,the,true,filtering,distributions,for,the,state,space,model,here,,we,cannot,easily,directly,measure,the,estimation,performance.,An,indicative,proxy,for,whether,a,particle,filter,is,likely,to,35,128256512102420484096Number,of,ranks12481632Number,of,nodes0.00.20.40.60.81.0Weak,scaling,efficiency,Figure,35:,Strong,scaling,efficiency,for,mpi-based,parallelization,of,particle,filter,based,data,assim-,ilation,using,ParticleDA,with,NektarDriftwave,state,space,model,on,ARCHER2.,Dotted,line,shows,ideal,linear,speed-up.,Number,of,particles,N,mint,(cid:100)ESS(t,,N,),mediant,(cid:100)ESS(t,,N,),maxt,(cid:100)ESS(t,,N,),128,256,512,1024,2048,4096,3.17,6.16,3.99,3.44,23.5,2.54,66.5,142,251,486,954,1580,120,244,477,940,1800,3810,Table,3:,Estimated,effective,sample,size,(ess),statistics,for,particle,filter,runs,on,NektarDriftwave,model,with,bootstrap,proposals,for,varying,ensemble,size,N,.,36,1282565121024Number,of,ranks1248Number,of,nodes1248Strong,scaling,speed-up,be,giving,poor,estimates,of,the,filtering,distributions,is,the,effective,sample,size,(ess),,a,measure,of,the,number,of,independent,samples,from,the,filtering,distribution,that,would,give,a,Monte,Carlo,estimator,with,equivalent,variance,to,that,computed,using,the,importance,weighted,samples,from,the,particle,filter.,An,estimate,of,the,ess,for,the,filtering,distribution,at,time,index,t,can,be,computed,from,the,particle,importance,weights,as,(cid:100)ESS(t,,N,),=,1/,N,(cid:88),(cid:16),n=1,w(n),t,(cid:17)2,.,Table,3,shows,summary,statistics,of,the,ess,estimates,computed,from,runs,of,the,bootstrap,particle,filter,on,the,NektarDriftwave,model,(using,the,same,T,=,100,simulated,observation,sequence,and,model,set,up,as,previously),for,varying,ensembles,sizes,N,.,The,ideal,case,would,be,if,the,ess,estimates,remained,close,to,the,ensemble,size,N,.,Low,ess,estimates,are,instead,indicative,of,high,variance,in,the,importance,weights,and,a,possible,tendency,to,the,ensemble,collapse,degeneracy,discussed,previously.,Here,we,see,that,the,minimum,ess,estimates,across,the,filtering,runs,are,consistently,low,for,all,ensemble,sizes,tested,,suggesting,the,filtering,distribution,estimates,from,the,bootstrap,particle,filter,used,here,will,have,high,variance.,This,is,also,evident,when,comparing,the,particle,filter,estimates,of,the,mean,of,the,filtering,distributions,to,the,true,state,fields,used,to,generate,the,observations,used,for,filtering,,with,the,latter,corresponding,to,samples,from,the,true,filtering,distributions,in,this,simulated,data,setting.,The,true,state,fields,and,corresponding,estimates,of,the,filtering,distribution,(0,,50,,100),in,Figure,36,,with,the,particle,filter,estimates,means,are,shown,for,three,time,indices,t,here,computed,using,an,ensemble,size,of,N,=,2048.,It,can,be,seen,that,the,true,and,estimated,mean,state,fields,differ,significantly,at,later,time,indices,,and,while,the,filtering,distribution,means,would,not,be,expected,to,exactly,match,the,true,states,,the,large,differences,here,and,small,estimated,filtering,distribution,variances,(not,shown),are,indicative,of,weight,degeneracy,causing,the,ensemble,to,collapse,and,no,longer,‘track’,the,true,state,,as,also,indicated,by,the,low,minimum,ess,estimates.,The,poor,statistical,performance,here,is,not,unexpected,as,the,dimensionality,of,the,state,space,here,is,relatively,large,dx,=,32768,,and,while,the,observation,dimension,dy,=,9,is,small,the,signal,to,noise,ratio,in,the,observations,is,relatively,high.,While,we,may,be,able,to,achieve,somewhat,improved,performance,using,the,locally,optimal,proposal,,while,the,current,formulation,of,the,NektarDriftwave,state,space,model,technically,should,allow,tractable,computation,of,the,locally,optimal,proposal,dis-,tribution,,this,is,not,yet,supported,in,the,current,implementation.,More,generally,,as,noted,previously,to,achieve,robust,statistical,performance,on,complex,high-dimensional,state,space,models,,additional,approaches,such,as,tempering,and,spatial,localization,are,likely,to,be,required.,Implementing,support,for,such,approaches,in,ParticleDA,would,be,an,important,avenue,to,consider,in,any,future,work.,∈,37,(a),True,state,,t,=,0,(b),True,state,,t,=,50,(c),True,state,,t,=,100,(d),Estimated,state,,t,=,0,(e),Estimated,state,,t,=,50,(f),Estimated,state,,t,=,100,Figure,36:,True,state,fields,at,three,time,indices,t,(0,,50,,100),for,NektarDriftwave,state,space,model,∈,and,corresponding,particle,filter,estimates,of,the,mean,of,the,filtering,distributions,using,an,ensemble,with,N,=,2048,particles.,38,-0.999-0.711-0.423-0.1350.153,zeta-0.05560.208,0.471,0.734,0.997,n-0.857-0.452-0.04790.357,0.761,zeta-0.784-0.432-0.08030.271,0.623,n-2.33,-1.11,0.107,1.33,2.55,zeta-1.92,-0.993-0.06310.867,1.80,n-1.01,-0.724-0.433-0.1430.148,zeta-0.05450.212,0.479,0.746,1.01,n-0.857-0.452-0.04790.357,0.761,zeta-0.578-0.2340.111,0.455,0.800,n-2.98,-1.54,-0.1011.34,2.78,zeta-1.78,-0.7660.244,1.25,2.26,n,11,FabNESO,As,part,of,ucl’s,work,on,aquifer,,we,developed,FabNESO,(https://github.com/UCL/fabneso),a,NEPTUNE,exploratory,software,(neso),plug-in,for,FabSim3,(Groen,et,al.,2023).,FabNESO,is,designed,to,facilitate,the,execution,of,neso,simulations,on,both,local,and,remote,high,performance,computing,systems,via,a,unified,interface.,It,defines,a,series,of,tasks,to,support,different,applications,and,analyses,with,neso,simulations,,with,tasks,defined,for,running,single,and,ensembles,of,simulations,,forward,uncertainty,quantification,using,EasyVVUQ,(Richardson,et,al.,2020),and,accelerated,Bayesian,calibration,of,model,parameters,using,PyVBMC,(Huggins,et,al.,2023).,Instructions,for,installing,and,configuring,FabNESO,are,provided,on,the,repository,README,and,html,documentation,for,the,interface,for,all,tasks,and,utility,functions,is,provided,at,http://github-pages.ucl.ac.uk/,FabNESO/.,FabNESO,can,be,used,with,any,solver,executable,provided,by,neso,(or,Nektar++,more,generally),which,follows,the,standard,Nektar++,solver,command,line,interface,—,that,is,it,will,accept,as,arguments,paths,to,one,or,more,xml,files,specifying,the,conditions,of,the,equation,system,to,solve,and,geometry,and,mesh,discretization,of,the,spatial,domain,to,solve,on.,Tasks,take,a,solver,argument,which,can,be,used,to,specify,the,name,of,the,solver,executable,within,the,path,specified,by,the,(machine,specific),neso,bin,dir,configuration,option,,which,should,point,to,the,local,path,containing,the,built,neso,binaries.,Pre-defined,configurations,are,provided,for,the,Electrostatic2D3V,and,H3LAPD,solvers,based,on,respectively,the,two,stream,and,2Din3D-hw,examples,for,these,solvers,in,the,main,NESO,repository,,with,these,configuration,directories,containing,the,conditions,and,mesh,XML,files,to,use,,with,tasks,taking,a,config,argument,to,specify,the,configuration,to,use.,Users,can,add,additional,configurations,by,adding,further,sub-directories,with,the,necessary,files,to,the,FabNESO,plugin,config,files,sub-directory.,Most,tasks,also,provide,some,support,for,overriding,the,default,values,of,the,parameters,specified,in,these,configuration,files,-,either,directly,for,single,simulations,,or,as,part,of,a,scan,over,multiple,parameter,values,as,part,of,ensembles,of,runs,for,uncertainty,quantification,or,calibration.,11.1,Forward,uncertainty,quantification,-,neso,pce,ensemble,task,The,neso,pce,ensemble,task,allows,running,an,ensemble,of,neso,simulations,for,performing,a,pce,of,one,or,more,(functional),outputs,of,the,solver,,building,on,the,support,for,pce,in,EasyVVUQ,and,chaospy,(Feinberg,and,Langtangen,2015).,The,results,of,the,pce,ensemble,can,then,be,analysed,to,perform,forward,uncertainty,quantification,on,the,outputs,and,sensitivity,analyses.,For,example,the,shell,command,fabsim,archer2,neso_pce_ensemble:\,config=2Din3D-hw,solver=H3LAPD,\,nodes=4,processes=512,\,polynomial_order=3,variant=pseudo-spectral,\,HW_alpha=0:2,Te_eV=5:20,∈,[0,,2],and,Te,eV,would,schedule,an,ensemble,of,(polynomial,order,+,1)**n,parameters,=,(3,+,1)2,=,16,jobs,on,ARCHER2,running,H3LAPD,solver,instances,with,parameter,values,corresponding,to,quadrature,nodes,over,a,uniform,distribution,on,the,two-dimensional,parameter,space,spanned,by,HW,alpha,[5,,20],,with,the,specific,pce,variant,here,being,pseudo-spectral,projection.,Each,∈,ARCHER2,jobs,would,be,scheduled,to,run,on,512,mpi,ranks,(processes),across,4,nodes.,Once,this,command,is,executed,,we,can,use,fabsim,archer2,job,stat,to,check,on,the,job,statuses,and,fabsim,archer2,fetch,results,to,pull,the,job,results,back,to,our,local,machine,once,the,jobs,have,completed.,A,separate,neso,pce,analysis,task,is,provided,for,analysing,the,outputs,of,a,previous,pce,en-,semble,task,run.,For,example,if,the,output,of,the,fetch,results,task,indicated,the,results,had,been,pulled,locally,to,the,directory,/path/to/results,we,could,analyse,the,outputs,of,the,previous,neso,pce,ensemble,task,on,our,local,system,by,running,fabsim,localhost,neso_pce_analysis:config=2Din3D-hw,results_dir=/path/to/results,39,Figure,37:,Estimated,first,order,Sobol,indices,of,E(t),outputted,by,H3LAPD,solver,run,using,2Din3D-hw,example,configuration,with,respect,to,the,parameters,HW,alpha,and,Te,eV.,with,then,outputting,a,pickle,file,containing,a,serialized,version,of,the,EasyVVUQ,pce,analysis,results,object,corresponding,to,the,ensemble.,This,can,then,be,used,to,compute,various,statistics,of,the,model,outputs,,construct,a,surrogate,model,or,perform,sensitivity,analyses,of,the,outputs,with,respect,to,the,parameters.,Example,of,output,plots,produced,using,the,pce,analysis,results,object,for,a,pce,ensemble,run,executed,as,shown,above,on,ARCHER2,are,shown,in,Figures,37,and,38.,11.2,Calibrating,model,parameters,-,neso,vbmc,task,The,neso,vbmc,task,allows,accelerated,Bayesian,calibration,of,the,parameters,of,neso,solvers,using,variational,Bayesian,Monte,Carlo,(vbmc),(Acerbi,2018),,specifically,using,the,Python,implementa-,tion,in,PyVBMC,(Huggins,et,al.,2023).,The,vbmc,algorithm,is,an,approximate,Bayesian,inference,method,for,efficient,estimation,of,the,posterior,on,the,parameters,of,computationally,demanding,simu-,lator,models,given,observations,of,the,model,outputs.,Compared,to,alternative,approximate,inference,approaches,such,asmcmc,,it,allows,robust,estimation,of,the,posterior,over,moderately,dimensioned,parameter,spaces,(up,to,around,10),with,typically,far,fewer,model,evaluations,(often,in,the,low,hun-,dreds,rather,than,tens,of,thousands,typical,of,mcmc,methods),and,treating,the,model,as,a,‘black-box’,-,for,example,not,requiring,the,ability,to,compute,derivatives,with,respect,to,the,model,parameters.,These,characteristics,make,it,well,suited,the,calibration,of,neso,simulations,,which,are,expensive,to,evaluate,,have,no,gradient,information,available,and,(currently),have,a,relatively,small,number,of,parameters,we,may,wish,to,calibrate.,vbmc,combines,variational,inference,,an,efficient,approach,for,fitting,an,approximation,to,target,posterior,distribution,from,a,prescribed,parametric,family,,with,Bayesian,quadrature,,an,approach,for,estimating,the,value,of,integrals,where,the,integrand,is,expensive,to,evaluate,by,computing,a,Gaussian,process,surrogate.,It,employs,a,iterative,sequential,design,strategy,,alternating,between:,actively,sampling,new,parameter,values,to,evaluate,the,log,posterior,density,(model),at,based,on,a,current,Gaussian,process,model,for,the,log,posterior,,balancing,exploration,of,areas,of,the,parameter,space,where,we,are,uncertain,about,the,log,posterior,density,with,exploiting,of,regions,containing,high,posterior,mass;,updating,the,Gaussian,process,surrogate,for,the,log,density,based,on,the,newly,40,0100200300400500600Time,step0.00.20.40.60.81.0First,Order,Sobol,IndexHW_alphaTe_eVhigher,orders,Figure,38:,Estimated,moments,of,energy,trajectory,E(t),outputted,by,H3LAPD,solver,run,using,2Din3D-hw,example,configuration,with,a,uniform,distribution,over,the,parameters,HW,alpha,and,Te,eV.,acquired,log,posterior,density,values;,updating,the,variational,approximation,to,the,target,posterior,using,the,Gaussian,process,surrogate.,The,neso,vbmc,task,currently,supports,calibration,of,the,Electrostatic2D3V,solver,,with,a,pre-,defined,function,for,extracting,the,values,of,the,simulated,electrostatic,potential,field,along,a,line,at,the,final,simulation,time,from,the,hdf5,files,output,by,the,solver,provided.,Integration,with,additional,solvers,would,require,identification,of,suitable,outputs,to,calibrate,against,and,implementation,of,a,similar,function,to,extract,the,necessary,values,from,the,solver,outputs.,The,task,can,be,run,locally,by,specifying,the,localhost,machine,in,the,fabsim,command,or,on,a,remote,system,such,as,ARCHER2,(specified,using,archer2,machine,argument).,The,following,command,fabsim,archer2,neso_vbmc:\,config=two_stream,solver=Electrostatic2D3V,\,processes=16,reference_field_file=observations.txt,\,particle_initial_velocity=0:2:0.5:1.5,\,particle_charge_density=0:200:50:150,\,particle_number_density=0:200:50:150,would,run,PyVBMC,to,estimate,the,posterior,on,the,parameters,of,the,Electrostatic2D3V,solver,particle,initial,velocity,,particle,charge,density,and,particle,number,density,given,a,set,of,(noisy),observations,of,the,final,electrostatic,potential,field,along,a,line,recorded,in,a,data,file,observations.txt.,The,colon,delimited,lists,after,the,parameter,names,specify,the,bounds,of,the,parameters,(potentially,infinite),and,‘plausible’,bounds,specifying,a,range,over,which,it,is,believed,likely,most,of,the,posterior,mass,lies.,The,task,here,would,sequentially,schedule,jobs,running,the,solver,with,the,actively,sampled,parameter,values,on,ARCHER2,(with,the,example,command,here,specifying,to,run,each,job,on,16,mpi,ranks),,waiting,for,the,completion,of,the,job,and,computation,of,the,corresponding,log,posterior,density,value,before,scheduling,a,new,run,at,the,next,actively,sampled,parameter,set.,41,0100200300400500600Time,step0100200300400500Estdmean1%99%,Figure,39:,Simulated,noisy,observations,(orange,points),and,true,underlying,values,(blue,curve),,for,electrostatic,potential,field,along,line,at,final,timestep,of,a,simulation,using,neso,Electrostatic2D3V,solver,with,two,stream,example,configuration.,Once,PyVBMC,has,identified,certain,convergence,criteria,have,been,met,,the,final,variational,posterior,approximation,is,returned,to,allow,use,in,subsequent,analyses.,The,neso,vbmc,task,saves,a,serialization,of,this,variational,posterior,object,to,a,pickle,file.,Figure,39,shows,an,example,set,of,simulated,observations,of,the,final,electrostatic,potential,field,along,a,line,at,100,regularly,spaced,points,,generated,using,parameters,particle,initial,velocity,=,0.96,,particle,charge,density,=,99,,and,particle,number,density,=,103,,with,independent,Gaussian,observation,noise,with,standard,deviation,0.1.,Figure,40,shows,the,posterior,estimate,on,the,parameters,parameters,particle,initial,velocity,,particle,charge,density,and,particle,number,density,of,the,Electrostatic2D3V,neso,solver,with,two,stream,example,configuration,,produced,by,PyVBMC,for,the,simulated,observations,in,Fig-,ure,39.,The,parameters,bounds,were,set,as,specified,in,the,neso,vbmc,task,command,example,above,,with,a,uniform,distribution,assumed,over,the,parameter,space.,In,this,case,the,vbmc,algorithm,con-,verged,after,evaluation,of,the,model,at,130,different,set,of,parameter,values,,with,the,true,parameter,values,used,to,simulate,the,observations,contained,within,the,bulk,of,the,estimated,posterior,mass,as,expected.,With,a,relatively,limited,number,of,model,evaluations,the,vbmc,algorithm,here,has,been,able,to,capture,a,plausible,estimate,of,the,posterior,distribution.,42,40:,on,the,Estimated,posterior,three,parameters,particle,initial,velocity,,Figure,particle,charge,density,and,particle,number,density,(other,parameters,fixed),of,the,neso,Electrostatic2D3V,solver,with,the,two,stream,example,configuration,,given,the,simulated,obser-,vations,show,in,Figure,39,,estimated,using,PyVBMC,using,130,evaluations,of,the,model.,Posterior,density,estimate,is,shown,in,blue,and,true,parameter,values,used,to,generate,the,observations,in,orange.,43,12,Additional,Bayesian,calibration,tools,In,addition,to,the,FabNESO,integration,with,PyVBMC,,we,also,developed,a,basic,Python,package,nesopy,,,available,in,the,repository,https://github.com/UCL/neso-calibration,,which,provides,a,simple,interface,for,executing,neso,solvers,with,specified,parameter,values,and,extracting,the,gener-,ated,outputs,,to,allow,easy,integration,with,existing,calibration,tools,such,as,PyVBMC.,The,neso,solvers,used,can,either,be,built,natively,on,the,local,system,or,within,a,Docker,container,,with,nesopy,providing,wrappers,for,both,options.,A,notebook,illustrating,how,to,use,the,nesopy,Python,wrapper,to,perform,calibration,of,a,neso,simulation,with,PyVBMC,is,provided,in,the,neso-calibration,repository,under,examples/two,stream/pyvbmc,calibration,example.ipynb.,A,Python,implementation,of,an,alternative,calibration,approach,exploiting,Gaussian,process,em-,ulation,of,the,posterior,density,was,also,developed,during,AQUIFER,and,is,available,at,https:,//github.com/UCL/calibr,,with,the,initial,intention,of,exploring,this,as,an,alternative,to,PyVBMC,for,calibration,of,neso,simulations,,though,time,constraints,meant,that,this,work,was,not,completed.,calibr,is,a,Python,implementation,of,the,algorithm,described,in,Parallel,Gaussian,process,surrogate,Bayesian,inference,with,noisy,likelihood,evaluations,(J¨arvenp¨a¨a,et,al.,2021).,It,is,designed,to,allow,estimation,of,the,posterior,distribution,on,the,unknown,parameters,of,expensive,to,evaluate,simulator,models,given,observed,data,,using,a,batch,sequential,design,strategy,which,iterates,fitting,a,Gaussian,process,emulator,to,a,set,of,evaluations,of,the,(unnormalized),posterior,density,for,the,model,and,using,the,emulator,to,identify,a,new,batch,of,model,parameters,at,which,to,evaluate,the,posterior,density,which,minimize,a,measure,of,the,expected,uncertainty,in,the,emulation,of,the,posterior,density.,The,posterior,density,can,be,evaluated,at,the,parameter,values,in,each,batch,in,parallel,,provid-,ing,the,opportunity,for,speeding,up,calibration,runs,on,multi-core,and,multi-node,high,performance,computing,systems.,The,acquisition,functions,used,to,choose,new,parameter,values,to,evaluate,are,im-,plemented,using,the,high-performance,numerical,computing,framework,JAX,,,with,the,gradient-based,optimization,of,these,acquisition,functions,exploiting,JAX’s,support,for,automatic,differentiation.,13,Software,implementation,and,HPC,deployment,Note,than,Archer2,cycles,were,provided,to,ukaea,,and,to,two,of,the,non-UCL,neptune,grantees,from,the,seavea,allocation.,Deliverable,D3.1,and,D3.2:,D3.1:,The,new,release,28,March,2023,of,the,seavea,toolkit,featuring,advanced,surrogate,mod-,elling,(see,above,new,ssc,method).,Further,releases,will,integrate,more,of,the,methods,developed,in,aquifer.,D3.2:,Release,of,the,data,assimilation,(da),platform,was,achieved,as,FabParticleDA.,The,link,to,the,seavea,uq,platform,is,within,the,seavea,uq,platform,released,28,March,2023:,https://github.com/djgroen/FabParticleDA/tree/master.,The,test,case,works,locally.,hpc,de-,ployment,at,scale,is,part,of,Activity,4.,13.1,HPC,deployment,As,part,of,our,efforts,to,support,neptune,,we,have,been,working,on,coupling,codes,using,MUSCLE3,as,outlined,in,the,proposal.,Firstly,,we,replaced,the,mpi,implementation,with,MUSCLE3,and,verified,a,test,case,with,both,MUSCLE3,and,mpi,versions.,We,found,that,the,results,from,the,MUSCLE3,implementation,matched,those,obtained,with,the,MPI,version.,We,then,conducted,a,test,case,on,single,desktop,,which,showed,that,MUSCLE3,was,10,times,slower,than,the,mpi,version.,However,,we,expected,this,due,to,the,smaller,size,of,the,test,case,,and,communications,time,on,MUSCLE3,took,too,much,time,comparing,with,the,computational,time.,Next,,we,deployed,MUSCLE3,on,ARCHER2,using,a,larger,test,case,and,ran,both,an,mpi,job,and,a,MUSCLE3,job,on,the,same,resources.,We,found,that,they,completed,almost,at,the,same,time.,This,gives,us,confidence,that,MUSCLE3,has,the,potential,to,efficiently,couple,large,HPC,codes.,Right,now,,we’ve,been,allocated,54,000,computing,units,(cus),hours,on,ARCHER2.,One,cu,indicate,a,44,Figure,41:,Performance,scaling,plot,of,HemeLB,Mazzeo,and,Coveney,2008,deploys,on,Frontier,and,Crusher.,node,with,128,processor,cores,running,for,one,hour.,Moving,forward,,our,next,step,is,to,evaluate,the,performance,of,MUSCLE3,with,larger,test,cases,with,more,complex,scenarios,which,should,help,the,neptune,project,perform,hpc,simulations,across,different,scales.,13.1.1,HemeLB:,Performance,benchmark,on,Frontier,HemeLB-GPU5,has,undergone,extensive,testing,on,both,the,Frontier,and,Crusher,machines.,The,per-,formance,tests,conducted,on,these,machines,provide,valuable,insights,into,the,scalability,of,HemeLB.,Figure,29,illustrates,the,strong,scaling,performance,test,results,of,HemeLB-GPU,(pressure,boundary,conditions),using,different,vascular,models,on,both,Frontier,and,Crusher.,Three,vascular,models,were,tested:,the,right,leg,vascular,model,and,the,circle,of,Willis,vascular,model,at,resolutions,of,15,µm,and,6.4,µm,,respectively.,These,models,consist,of,lattice,sites,ranging,from,5.2x107,to,1.1x1010,and,ran,for,10,000,time-steps.,To,exploit,the,potential,for,exascale,computing,here,,we,have,further,tested,pipe,flow,scaling,with,3.9x1010,lattice,sites.,In,the,optimal,scaling,region,,HemeLB,demonstrated,a,single-node,performance,of,approximately,4000,million,lattice,cell,updates,per,second,(mlups).,This,performance,exhibited,linear,strong,scaling,with,the,number,of,nodes,until,degradation,sets,in,due,to,the,problem,size.,We,have,also,demonstrated,continued,strong,scaling,behaviour,on,up,to,87%,of,the,full,Frontier,system,(8192,nodes).,13.1.2,HemeLB:,Validation,of,turbulent,model,using,large-eddy,simulation,techniques,To,validate,the,channel,flow,simulation,,a,high-resolution,stenotic,channel,flow,simulation,is,estab-,lished.,As,depicted,in,Fig.,42,,the,dimensions,of,the,stenotic,channel,flow,are,configured,as,Lx,Lz,,where,Lx,=,266δ,,Ly,=,22δ,,and,Lz,=,2δ,correspond,to,the,streamwise,,spanwise,,and,vertical,di-,rections,,respectively.,Here,,δ,denotes,the,turbulent,boundary,layer,thickness,,set,at,δ,=,0.0045m.,The,lattice,resolution,is,set,to,100µm,,resulting,in,a,simulation,domain,composed,of,approximately,1.0,×,In,Fig.,43,,the,experimental,data,begins,from,y+,=,10,,whereas,the,y+,for,the,first,cell,near,the,wall,in,the,lattice,Boltzmann,method,(lbm)–large-eddy,simulation,(les),simulation,is,approximately,y+,=,1.5.,Although,there,is,a,minor,deviation,in,the,first,cell,,the,lbm–les,simulation,aligns,well,with,both,experimental,and,direct,numerical,simulation,(dns),references.,For,y+,>,30,,the,lbm–les,results,deviate,slightly,from,the,DNS,data,but,remain,in,close,agreement,with,the,experimental,results.,109,lattice,cells.,Ly,×,×,5https://github.com/UCL-CCS/HemePure-GPU,45,100101102103104Numberofnodes104105106107CumulativeMLUPSOptimalScaling75%ScalingCrusherRightLeg5.2£107latticesitesCrusherCoW15µm7.8£108latticesitesFrontierCoW6.4µm1.1£1010latticesitesFrontierExaPipe3.9£1010latticesites,Figure,42:,Illustration,of,a,stenotic,channel,flow,simulation,setup.,The,green,zone,demonstrates,the,semi-cylindrical,obstacle’s,configuration.,The,blue,zone,shows,a,snapshot,where,the,channel,flow,has,evolved,into,a,fully,developed,turbulent,velocity,profile.,The,orange,zone,represents,the,sponge,zone,,designed,to,absorb,reflective,waves,emanating,from,the,outlet.,Figure,43:,Representation,of,u+,as,a,function,of,y+,,where,the,black,dotted,line,illustrates,the,dns,reference,data,Kim,,Moin,,and,Moser,1987.,The,blue,triangle,line,depicts,the,piv,experiment,data,from,Ding,et,al.Ding,et,al.,2021.,The,red,triangle,line,represents,the,data,obtained,from,the,lbm–les,simulation.,Overall,,the,lbm–les,implementation,demonstrates,good,concordance,with,both,experimental,and,dns,simulations.,Furthermore,,we,checked,the,dimensionless,root,mean,square,(rms),for,two,velocity,components,on,streamwise,and,vertical,direction.,We,didn’t,include,the,spanwise,direction,due,to,the,comparison,with,the,particle,image,velocimetry,(piv),experiment.,u+,′rms,are,the,dimensionless,rms,velocity,components,that,are,normalized,with,the,shear,velocity,uτ,:,′rms,,v+,(cid:112),)2,.,−,⟨,u,⟩,u+,′rms,=,(u(x),uτ,Spatial,averaging,for,rms,is,performed,only,during,post-processing,due,to,the,high,resolution,of,the,simulation.,In,Figure,43,,the,red,dots,represent,the,lbm–les,simulation,results,,while,the,blue,dots,correspond,to,the,piv,experimental,reference,,and,the,black,lines,denote,the,dns,reference.,Both,experimental,and,lbm–les,results,successfully,capture,the,peak,value,of,u+,′rms,compared,with,the,dns,data.,Beyond,the,peak,,both,lbm–les,and,experimental,results,gradually,deviate,from,the,dns,data,,which,may,be,attributed,to,statistical,issues,and,the,confinement,of,the,channel,flow.,Regarding,′rms,,lbm–les,results,align,well,with,the,experimental,data,but,are,the,vertical,velocity,component,u+,slightly,lower,than,the,dns,results,for,y+,>,20.,(13),Considering,both,Figure,43,and,Figure,44,it,is,evident,that,the,lbm–les,implementation,aligns,closely,with,experimental,and,dns,results,,demonstrating,its,capability,to,capture,turbulent,statistical,quantities,accurately.,The,comparison,reveals,that,the,lbm–les,method,effectively,replicates,the,46,100101102y+05101520u+DNSresultsExperiment2023LBMresults,′,as,functions,of,y+,,with,the,black,and,black,dotted,lines,+,Figure,44:,Presentation,of,urms,representing,the,dns,reference,data,Kim,,Moin,,and,Moser,1987.,The,blue,round,and,square,lines,′,from,the,experiment,Ding,et,al.,2021,,respectively.,Additionally,,the,correspond,to,urms,red,round,and,square,lines,depict,the,lbm,simulation,outcomes,for,urms,′,,respectively.,′,and,vrms,′,and,vrms,′,and,vrms,+,+,+,+,+,key,features,of,turbulent,flows,,confirming,its,reliability,in,modeling,complex,flow,dynamics.,This,alignment,shows,the,potential,of,lbm–les,in,contributing,valuable,insights,into,the,understanding,of,turbulence,,particularly,in,the,context,of,fluid,dynamics,simulations.,47,020406080100y+0.00.51.01.52.02.5urms+′,vrms+′DNSu+′rmsDNSv+′rmsExperimentu+′rmsExperimentv+′rmsLBMu+′rmsLBMv+′rms,Acronyms,anaet,axissymmetric,non-axissymmetric,extended,aquifer,advanced,quantification,of,uncertainties,in,fusion,modelling,at,the,exascale,with,model,order,reduction,cdf,cumulative,distribution,function,cl,contour,localization,cu,computing,unit,da,data,assimilation,das,deep,active,subspace,dgp,deep,Gaussian,process,dns,direct,numerical,simulation,epsrc,Engineering,and,Physical,Sciences,Research,Council,esmacs,enhanced,sampling,of,molecular,dynamics,with,approximation,of,continuum,solvent,ess,effective,sample,size,gkdr,gradient-based,kernel,dimension,reduction,gp,Gaussian,process,hdf5,hierarchical,data,format,,version,5,hgp,heteroskedastic,Gaussian,process,hpc,high,performance,computing,html,hypertext,markup,language,kas-gp,kernelized,active,subspace,Gaussian,process,lbm,lattice,Boltzmann,method,les,large-eddy,simulation,lgp,linked,Gaussian,process,lhd,latin,hypercube,design,lse,level,set,estimation,mcmc,Markov,chain,Monte,Carlo,md,molecular,dynamics,mlups,million,lattice,cell,updates,per,second,mogp,multi-output,Gaussian,process,mpi,message,passing,interface,neptune,neutrals,and,plasma,turbulence,numerics,for,exascale,neso,NEPTUNE,exploratory,software,48,nrmse,normalized,root,mean,square,error,ode,ordinary,differential,equation,ope,outer,product,emulator,pca,principal,components,analysis,pce,polynomial,chaos,expansion,pde,partial,differential,equation,pi,principal,investigator,piv,particle,image,velocimetry,pod,proper,orthogonal,decomposition,qoi,quantity,of,interest,raddish,real-time,advanced,data,assimilation,for,digital,simulation,of,numerical,twins,on,HPC,rbc,Rayleigh–B´enard,convection,rkhs,reproducing,kernel,Hilbert,space,rms,root,mean,square,rmse,root,mean,square,error,sa,sensitivity,analysis,sc,stochastic,collocation,sde,stochastic,differential,equation,seavea,software,environment,for,actionable,and,VVUQ-evaluated,exascale,applications,shm,simple,harmonic,motion,spde,stochastic,partial,differential,equation,ssc,simplex,stochastic,collocation,tbr,Tritium,breeding,ratio,ties,thermodynamic,integration,with,enhanced,sampling,ucl,University,College,London,ukaea,United,Kingdom,Atomic,Energy,Association,uq,uncertainty,quantification,vbmc,variational,Bayesian,Monte,Carlo,vvuq,verification,,validation,and,uncertainty,quantification,xml,extensible,markup,language,49,References,Acerbi,,Luigi,(2018).,“Variational,Bayesian,Monte,Carlo”.,In:,Advances,in,Neural,Information,Pro-,cessing,Systems,31.,Arter,,Wayne,(2012).,Blue,Sky,Solutions,to,the,Magnetohydrodynamic,Trigger,Problem.,UK-Germany,National,Astronomy,Meeting,NAM2012.,url:,https://doi.org/10.13140/RG.2.2.35052.,77449.,Beck,,Joakim,and,Serge,Guillas,(2016).,“Sequential,design,with,mutual,information,for,computer,experiments,(MICE):,Emulation,of,a,tsunami,model”.,In:,SIAM/ASA,Journal,on,Uncertainty,Quantification,4.1,,pp.,739–766.,Beskos,,Alexandros,,Dan,Crisan,,and,Ajay,Jasra,(2014).,“On,the,stability,of,sequential,Monte,Carlo,methods,in,high,dimensions”.,In:,The,Annals,of,Applied,Probability,24.4,,pp.,1396–1445.,doi:,10.1214/13-AAP951.,url:,https://doi.org/10.1214/13-AAP951.,Beskos,,Alexandros,,Dan,Crisan,,Ajay,Jasra,,et,al.,(2017).,“A,stable,particle,filter,for,a,class,of,high-dimensional,state-space,models”.,In:,Advances,in,Applied,Probability,49.1,,pp.,24–48.,Booth,,Annie,S,,S,Ashwin,Renganathan,,and,Robert,B,Gramacy,(2023).,“Contour,Location,for,Reliability,in,Airfoil,Simulation,Experiments,using,Deep,Gaussian,Processes”.,In:,arXiv,preprint,arXiv:2308.04420.,Boyd,,John,P,(1992).,“Defeating,the,Runge,phenomenon,for,equispaced,polynomial,interpolation,via,Tikhonov,regularization”.,In:,Applied,mathematics,letters,5.6,,pp.,57–59.,Bunch,,Pete,and,Simon,Godsill,(2016).,“Approximations,of,the,optimal,importance,density,using,Gaussian,particle,flow,importance,sampling”.,In:,Journal,of,the,American,Statistical,Association,111.514,,pp.,748–762.,Cohn,,David,,Zoubin,Ghahramani,,and,Michael,Jordan,(1994).,“Active,learning,with,statistical,mod-,els”.,In:,Advances,in,neural,information,processing,systems,7.,Cole,,D,Austin,et,al.,(2023).,“Entropy-based,adaptive,design,for,contour,finding,and,estimating,reliability”.,In:,Journal,of,Quality,Technology,55.1,,pp.,43–60.,Ding,,Guanghui,et,al.,(2021).,“Transitional,pulsatile,flows,with,stenosis,in,a,two-dimensional,channel”.,In:,Physics,of,Fluids,33.3.,Doucet,,Arnaud,,Simon,Godsill,,and,Christophe,Andrieu,(2000).,“On,sequential,Monte,Carlo,sampling,methods,for,Bayesian,filtering”.,In:,Statistics,and,computing,10,,pp.,197–208.,Dullin,,Holger,R,et,al.,(2007).,“Extended,phase,diagram,of,the,Lorenz,model”.,In:,International,Journal,of,Bifurcation,and,Chaos,17.09,,pp.,3013–3033.,Edeling,,Wouter,Nico,,Richard,P,Dwight,,and,Pasquale,Cinnella,(2016).,“Simplex-stochastic,collo-,cation,method,with,improved,scalability”.,In:,Journal,of,Computational,Physics,310,,pp.,301–,328.,Evensen,,G.,(2006).,Data,Assimilation:,The,Ensemble,Kalman,Filter.,Springer,Berlin,Heidelberg.,isbn:,9783540383017.,url:,https://books.google.co.uk/books?id=VJ2oOecHhOYC.,Evensen,,Geir,(1994).,“Sequential,data,assimilation,with,a,nonlinear,quasi-geostrophic,model,using,Monte,Carlo,methods,to,forecast,error,statistics”.,In:,Journal,of,Geophysical,Research:,Oceans,99.C5,,pp.,10143–10162.,Evensen,,Geir,,Femke,C,Vossepoel,,and,Peter,Jan,van,Leeuwen,(2022).,“Weak,Constraint,4DVar”.,In:,Data,Assimilation,Fundamentals:,A,Unified,Formulation,of,the,State,and,Parameter,Estimation,Problem.,Springer,,pp.,49–61.,Farchi,,Alban,and,Marc,Bocquet,(2018).,“Comparison,of,local,particle,filters,and,new,implementa-,tions”.,In:,Nonlinear,Processes,in,Geophysics,25.4,,pp.,765–807.,Fearnhead,,Paul,and,Hans,R,K¨unsch,(2018).,“Particle,filters,and,data,assimilation”.,In:,Annual,Review,of,Statistics,and,Its,Application,5,,pp.,421–449.,Feinberg,,Jonathan,and,Hans,Petter,Langtangen,(2015).,“Chaospy:,An,open,source,tool,for,designing,methods,of,uncertainty,quantification”.,In:,Journal,of,Computational,Science,11,,pp.,46–57.,issn:,1877-7503.,doi:,https,:,/,/,doi,.,org,/,10,.,1016,/,j,.,jocs,.,2015,.,08,.,008.,url:,https,:,/,/,www,.,sciencedirect.com/science/article/pii/S1877750315300119.,Frei,,Marco,and,Hans,R,K¨unsch,(2013).,“Bridging,the,ensemble,Kalman,and,particle,filters”.,In:,Biometrika,100.4,,pp.,781–800.,50,Gardner,,Timothy,S,,Charles,R,Cantor,,and,James,J,Collins,(2000).,“Construction,of,a,genetic,toggle,switch,in,Escherichia,coli”.,In:,Nature,403.6767,,pp.,339–342.,Giles,,D.,et,al.,(2023).,“ParticleDA.jl,v.1.0:,A,real-time,data,assimilation,software,platform”.,In:,Geoscientific,Model,Development,Discussions,2023,,pp.,1–20.,doi:,10.5194/gmd-2023-38.,url:,https://gmd.copernicus.org/preprints/gmd-2023-38/.,Gordon,,Neil,J,,David,J,Salmond,,and,Adrian,FM,Smith,(1993).,“Novel,approach,to,nonlinear/non-,Gaussian,Bayesian,state,estimation”.,In:,IEE,proceedings,F,(radar,and,signal,processing).,Vol.,140.,2.,IET,,pp.,107–113.,Gotovos,,Alkis,et,al.,(2013).,“Active,Learning,for,Level,Set,Estimation”.,In:,Proceedings,of,the,Twenty-,Third,International,Joint,Conference,on,Artificial,Intelligence.,AAAI,Press.,Graepel,,Thore,(2003).,“Solving,noisy,linear,operator,equations,by,Gaussian,processes:,Application,to,ordinary,and,partial,differential,equations”.,In:,ICML.,Vol.,3,,pp.,234–241.,Graham,,Matthew,M,and,Alexandre,H,Thiery,(2019).,“A,scalable,optimal-transport,based,local,particle,filter”.,In:,arXiv,preprint,arXiv:1906.00507.,Gramacy,,Robert,B,(2020).,Surrogates:,Gaussian,process,modeling,,design,,and,optimization,for,the,applied,sciences.,CRC,press.,Gramacy,,Robert,B,and,Herbert,K,H,Lee,(2008).,“Bayesian,treed,Gaussian,process,models,with,an,application,to,computer,modeling”.,In:,Journal,of,the,American,Statistical,Association,103.483,,pp.,1119–1130.,Groen,,Derek,et,al.,(2023).,“FabSim3:,An,automation,toolkit,for,verified,simulations,using,high,performance,computing”.,In:,Computer,Physics,Communications,283,,p.,108596.,Huggins,,Bobby,et,al.,(2023).,“PyVBMC:,Efficient,Bayesian,inference,in,Python”.,In:,Journal,of,Open,Source,Software,8.86,,p.,5428.,doi:,10.21105/joss.05428.,url:,https://doi.org/10.21105/,joss.05428.,J¨arvenp¨a¨a,,Marko,et,al.,(2021).,“Parallel,Gaussian,Process,Surrogate,Bayesian,Inference,with,Noisy,Likelihood,Evaluations”.,In:,Bayesian,Analysis,16.1,,pp.,147–178.,doi:,10.1214/20-BA1200.,url:,https://doi.org/10.1214/20-BA1200.,Kalman,,R.,E.,(Mar.,1960).,“A,New,Approach,to,Linear,Filtering,and,Prediction,Problems”.,In:,Journal,of,Basic,Engineering,82.1,,pp.,35–45.,issn:,0021-9223.,doi:,10.1115/1.3662552.,eprint:,https,:,/,/,asmedigitalcollection,.,asme,.,org,/,fluidsengineering,/,article,-,pdf,/,82,/,1,/,35,/,5518977/35\_1.pdf.,url:,https://doi.org/10.1115/1.3662552.,Kim,,John,,Parviz,Moin,,and,Robert,Moser,(1987).,“Turbulence,statistics,in,fully,developed,channel,flow,at,low,Reynolds,number”.,In:,Journal,of,fluid,mechanics,177,,pp.,133–166.,Lee,,John,C,and,Norman,J,McCormick,(2011).,Risk,and,safety,analysis,of,nuclear,systems.,John,Wiley,&,Sons.,Lindgren,,Finn,,H˚avard,Rue,,and,Johan,Lindstr¨om,(2011).,“An,explicit,link,between,Gaussian,fields,and,Gaussian,Markov,random,fields:,the,stochastic,partial,differential,equation,approach”.,In:,Journal,of,the,Royal,Statistical,Society,Series,B:,Statistical,Methodology,73.4,,pp.,423–498.,Ma,,Tian,and,Shouhong,Wang,(2005).,Bifurcation,theory,and,applications.,Vol.,53.,World,Scientific.,Mandel,,Jan,(2006).,Efficient,implementation,of,the,ensemble,Kalman,filter.,Tech.,rep.,231.,Center,for,Computational,Mathematics,,University,of,Colorado,at,Denver,and,Health,Sciences,Center.,Mazzeo,,Marco,D,and,Peter,V,Coveney,(2008).,“HemeLB:,A,high,performance,parallel,lattice-,Boltzmann,code,for,large,scale,fluid,flow,in,complex,geometries”.,In:,Computer,Physics,Com-,munications,178.12,,pp.,894–914.,McHutchon,,Andrew,(2015).,“Nonlinear,modelling,and,control,using,Gaussian,processes”.,PhD,thesis.,University,of,Cambridge.,Ming,,Deyu,and,Serge,Guillas,(2021).,“Linked,Gaussian,process,emulation,for,systems,of,computer,models,using,Mat´ern,kernels,and,adaptive,design”.,In:,SIAM/ASA,Journal,on,Uncertainty,Quan-,tification,9.4,,pp.,1615–1642.,Ming,,Deyu,,Daniel,Williamson,,and,Serge,Guillas,(2023).,“Deep,gaussian,process,emulation,using,stochastic,imputation”.,In:,Technometrics,65.2,,pp.,150–161.,Oakley,,Jeremy,(2004).,“Estimating,percentiles,of,uncertain,computer,code,outputs”.,In:,Journal,of,the,Royal,Statistical,Society,Series,C:,Applied,Statistics,53.1,,pp.,83–93.,51,Rabier,,Florence,and,Zhiquan,Liu,(2003).,“Variational,data,assimilation:,theory,and,overview”.,In:,Proc.,ECMWF,Seminar,on,Recent,Developments,in,Data,Assimilation,for,Atmosphere,and,Ocean,,Reading,,UK,,September,8–12,,pp.,29–43.,Rasmussen,,Carl,Edward,,Christopher,KI,Williams,,et,al.,(2006).,Gaussian,processes,for,machine,learning.,Vol.,1.,Springer.,Rathgeber,,Florian,et,al.,(2016).,“Firedrake:,automating,the,finite,element,method,by,composing,abstractions”.,In:,ACM,Transactions,on,Mathematical,Software,(TOMS),43.3,,pp.,1–27.,Rebeschini,,Patrick,and,Ramon,van,Handel,(2015).,“Can,local,particle,filters,beat,the,curse,of,dimen-,sionality?”,In:,The,Annals,of,Applied,Probability,25.5,,pp.,2809–2866.,doi:,10.1214/14-AAP1061.,url:,https://doi.org/10.1214/14-AAP1061.,Richardson,,Robin,A.,et,al.,(Apr.,2020).,“EasyVVUQ:,A,Library,for,Verification,,Validation,and,Un-,certainty,Quantification,in,High,Performance,Computing”.,In:,Journal,of,Open,Research,Software.,doi:,10.5334/jors.303.,Roth,,Michael,et,al.,(2017).,“The,Ensemble,Kalman,filter:,a,signal,processing,perspective”.,In:,EURASIP,Journal,on,Advances,in,Signal,Processing,2017,,pp.,1–16.,Rybin,,Mikhail,V,et,al.,(2015).,“Phase,diagram,for,the,transition,from,photonic,crystals,to,dielectric,metamaterials”.,In:,Nature,communications,6.1,,p.,10102.,Santner,,Thomas,J,et,al.,(2003).,The,design,and,analysis,of,computer,experiments.,Vol.,1.,Springer.,Scheuerer,,Michael,(2010).,“A,comparison,of,models,and,methods,for,spatial,interpolation,in,statistics,and,numerical,analysis”.,In.,Slivinski,,Laura,and,Chris,Snyder,(2016).,“Exploring,practical,estimates,of,the,ensemble,size,necessary,for,particle,filters”.,In:,Monthly,Weather,Review,144.3,,pp.,861–875.,Snyder,,Chris,(2011).,“Particle,filters,,the,“optimal”,proposal,and,high-dimensional,systems”.,In:,Proceedings,of,the,ECMWF,Seminar,on,Data,Assimilation,for,atmosphere,and,ocean.,Citeseer,,pp.,1–10.,Solak,,Ercan,et,al.,(2002).,“Derivative,observations,in,Gaussian,process,models,of,dynamic,systems”.,In:,Advances,in,neural,information,processing,systems,15.,Wendland,,Holger,(2004).,Scattered,data,approximation.,Vol.,17.,Cambridge,university,press.,Whittle,,Peter,(1954).,“On,stationary,processes,in,the,plane”.,In:,Biometrika,,pp.,434–449.,Zhang,,Ruda,,Simon,Mak,,and,David,Dunson,(2022).,“Gaussian,Process,Subspace,Prediction,for,Model,Reduction”.,In:,SIAM,Journal,on,Scientific,Computing,44.3,,A1428–A1449.,52

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