There are two main avenues to design space ex- ploration. In the ﬁrst approach, a simulation is run, analyzed, the problem modiﬁed, and the simulation run again. In the second approach, an ensemble simulation is performed and the battery of results is leveraged to construct a surrogate model for a given quantity of interest (QoI). The ﬁrst ap- proach allows a practitioner to methodically move through the design space and analyze a solution ﬁeld. A disadvan- tage of this technique is that each new simulation requires time consuming setup. The second approach provides the practitioner with a global view of the problem, but requires a priori design space limits and the QoI speciﬁcation. In this work we introduce an immersive simulation software frame- work that enables practitioners to maintain the ﬂexibility of the ﬁrst approach, while eliminating the burden of setting up new simulations. Immersive simulation can also be used to inform the second approach, establishing limits and clarify- ing QoI selection prior to the launch of an ensemble sim- ulation. We demonstrate live, reconﬁgurable visualization of on-going simulations coupled with live, reconﬁgurable problem deﬁnition that guides users in determining problem parameters. Ultimately, an immersive simulation framework enables more efﬁcient design space exploration that reduces the gap between simulations, data analysis and insight ex- traction. Keywords HPCScientiﬁc visualization Immersive simulationin situ simulationComputational steering Felix NewberryCorey Wetterer-Nelson John A EvansAlireza DoostanKenneth E Jansen Ann and H.J. Smead Aerospace Engineering Sciences Department, University of Colorado, Boulder, CO 80303, USA, E-mail: kenneth.jansen@colorado.edu1 Introduction The solution of systems of partial differential equations (P- DEs) provides scientists, engineers, and other practitioners with insight into a wide variety of physical phenomena. High performance computing (HPC) is the driving force behind increasingly accurate simulations that are capable of rapidly solving multi-physics, three-dimensional systems, or higher order approximations. The traditional simulation workﬂow for PDE solvers consists of sequential steps: meshing, deﬁ- nition of boundary conditions, solver and visualization. This structure leads to numerous inefﬁciencies. Problem redeﬁ- nition and the role of visualization as a post-processing step delays engineering/scientiﬁc discovery and design. In HPC applications, the practitioner’s pursuit of insight is hampered not only by time devoted to each new problem redeﬁnition, but also by additional queue time. These traditional work- ﬂow inefﬁciencies inhibit design space exploration and sub- sequent ensemble simulation. One method to address workﬂow inefﬁciencies is with notiﬁcation and monitoring systems that query simulation data during runtime [39, 37, 44]. While this output does pro- vide practitioners with basic solver behavior, such as con- vergence, it does not permit live reconﬁgurable simulation. Co-processing, the concept of performing in situ visualiza- tion and analysis while a simulation is ongoing, allows users to access simulation data immediately and mitigates the bur- den of writing large simulation data to disk [22, 46, 17, 29]. A summary of in situ methods, infrastructures and applica- tions to HPC is presented in the state-of-the-art (STAR) re- port [7]. The authors observe that while numerous infras- tructures are available, only a limited number of produc- tion quality frameworks have emerged. For instance, there are co-processing implementations within VisIt [11] using ADIOS [26], Lib-Sim [31], GLEAN [47], and in ParaView Catalyst (Catalyst ) [1, 4]. A concluding remark in [7] is 2 Felix Newberry et al. that given the entry cost for new in situ frameworks is rela- tively high, a move towards a generic data interface would be of great value. More recently, attempts have been made to address this concern. The introduction of the in situ in- terface SENSEI [6] and infrastructure ALPINE [27] provides an interface between a given solver and a variety of visu- alization and analysis routines. While SENSEI is designed solely to support Catalyst ,Lib-Sim andADIOS infrastruc- tures, ALPINE plays both the supporting role and provides its own execution model and visualization. In [5], a perfor- mance analysis of SENSEI shows the infrastructure is highly ﬂexible and has low overhead. A powerful extension of co-processing is computational steering, where users can take advantage of real time simula- tion feedback to rapidly explore the design space and extract insight. Computational steering has been previously demon- strated in the ﬂow solver PHASTA with Catalyst [50, 54], though steering was accomplished via live edits of a solver parameter ﬁle, and geometric deformations were not imple- mented. More recently, Catalyst has been used to model turbidity currents [9], where steering is introduced through application speciﬁc LibMesh -sed-imentation code on solver parameters such as time step and tolerances. Another ge- ological application utilized in situ visualization in earth- quake modeling but did not use computational steering [33]. In [51], computational steering is enabled via a user inter- face, however this software is closed source. These exam- ples show the current state of computational steering either requires additional software, the in situ editing of simulation solver parameter inputs or are closed source. Ensemble simulation, informed via prior design space exploration, is a useful means of further insight extraction. Through visualizing and analyzing solution output across a range of different realizations of input parameters, prac- titioners can better understand how some output quantity of interest (QoI) is a function of input variables. The ﬁeld of uncertainty quantiﬁcation (UQ) addresses the challenge of mapping uncertain inputs to output QoIs [19, 28, 52]. Large ensembles allow the collection of useful statistics, cre- ation of surrogate models for a given QoI, sensitivity anal- ysis, and much more. Further, the evaluation of ensembles has traditionally been performed in serial, with parallel im- plementations that are predominantly ad hoc . While there has been some commercial software development in this area [2], in this work we employ the open source Python library libEnsemble , developed at Argonne National Lab- oratory as part of the Department of Energy Exascale Com- puting Project to coordinate the concurrent evaluation of dy- namic ensembles of calculations on massively parallel re- sources [25]. In [10], libEnsemble is coupled with mul- tiobjective optimization to achieve better resource utiliza- tion for HPCs. Having performed an ensemble simulation, the practitioner may access one of several open source UQPython libraries such as Uncertainpy [45], chaospy [16], andUQ-PyL [48]. Additional Python libraries for speciﬁc UQ ﬁelds such as active subspaces are also in use [12]. The priority research directions for in situ data management iden- tiﬁed by the authors in [36] include ensemble analysis, UQ, and surrogate models. We note that there is a distinction be- tween random input variables addressed by UQ and design variables examined in design space exploration. In this work we apply UQ methodologies to design variables by treating them as uniform random variables. 1.1 Contributions of this work While numerous studies have demonstrated co-processing, to the authors’ knowledge, the valuable link between co- processing capabilities, extraction of relevant data from en- semble generation, and subsequent application of UQ tech- niques such as sensitivity analysis, has not been explored. In this document, we introduce an immersive simulation soft- ware framework that enables more rapid design space ex- ploration with computational steering and subsequently in- forms ensemble generation for more efﬁcient insight extrac- tion. Our work builds directly upon previous implementations ofCatalyst in the ﬂow solver PHASTA [54, 4]. Where ear- lier demonstrations of co-processing functioned with Cata- lyst hooked directly into the PHASTA solver, our work dif- fers in that we implement SENSEI to interface between the solver and Catalyst . The main contribution of this work is to introduce and demonstrate an immersive simulation software framework. We describe the relationship between a PDE solver and immersive simulation software compo- nents that enable computational steering and modiﬁcation of physical, geometric and solver parameters. We provide soft- ware linkages that others can leverage or modify as needed. Equipped with an immersive simulation framework, a prac- titioner may interact with their simulation in situ via the ParaView graphical user interface (GUI) . We complement these immersive simulation tools with a demonstration of libEnsemble for straightforward ensem- ble simulation and chaospy to implement UQ methods to analyze the ensemble results. We show the ability of im- mersive simulation and ensemble tools to complement and enhance one another. The remainder of this paper is organized as follows. Sec- tion 2 describes the simulation workﬂow, ﬁrst addressing the traditional approach in Section 2.1 followed by an immer- sive simulation workﬂow in Section 2.2. Next, Section 3 de- tails the software elements that make up an immersive sim- ulation software framework. Section 4 provides a numerical demonstration of immersive simulation tools on a 2D ag- gressive subsonic diffuser with the ﬂuid solver PHASTA . Sec- tion 5 describes ensemble simulation software tools before Software Tools to Enable Immersive Simulation 3 Section 6 applies the ensemble tools, coupled with an im- mersive simulation framework, to a global sensitivity analy- sis. Finally, Section 7 summarizes this study’s conclusions. 2 The simulation workﬂow We now describe the traditional simulation workﬂow, high- lighting aspects that inhibit scientists’ and engineers’ abil- ity to gain insight from their simulations. We then present an immersive simulation workﬂow that alleviates standard workﬂow inefﬁciencies. 2.1 Traditional simulation workﬂow The traditional simulation workﬂow, whether a lone simula- tion (Figure 1 (a)), or an ensemble simulation (Figure 1 (b)), consists of a series of sequential steps. Considering Figure 1 (a), we refer to components of the workﬂow that occur be- fore and after the solver step as belonging to the pre- and post-process, respectively. First, in the pre-process we de- velop an idea into a solvable problem through deﬁning an appropriate geometry, mesh, and boundary conditions. Next, we calculate the solution to our problem through our chosen PDE solver. Finally, in the post-process section we perform visualization and data analysis to determine the utility of the solution data. If the present results are satisfactory, then we have obtained the insight we seek and the simulation work- ﬂow is complete. Alternatively, if the results are unsatisfac- tory and we want to redeﬁne the problem to explore the de- sign space, then we return to the pre-process region of the workﬂow. Difﬁculties with the traditional simulation workﬂow a- rise due to the computational expense of individual work- ﬂow components coupled with the persistent need for prob- lem redeﬁnition. Solution data may persuade a practitioner to return to the pre-process steps either to ask a new question of the simulation, or to improve the method of answer to an existing question. As an example, let us consider the visual- ization of a ﬂuid ﬂow. Inspection of high gradient regions in the ﬂow, such as boundary and shear layers, may reveal that further mesh reﬁnement is needed to accurately resolve ﬂow behavior. Alternatively, the practitioner may decide geomet- ric adjustments would be a productive line of inquiry and be motivated to update the geometry and associated mesh. The role of problem redeﬁnition as a post-process step in the tra- ditional workﬂow delays engineering and scientiﬁc insight. A vital concern for HPC applications that will grow with the arrival of exascale computing is the bottleneck of IO. The PDE introduced in [15], in which an unsteady ﬂow prob- lem is addressed, was scaled up to 3 :1 million processes and 92 billion elements in [38]. To save the complete so- lution of a simulation this size would generate O(1)ter-abyte of data per second. Although the code can determine a write frequency and produce numerous ﬂow statistics, anal- ysis of unsteady ﬂow structures is not practical with the tra- ditional simulation workﬂow. To undertake post-processing steps that involve writing, and re-reading of data at this rate is unreasonable. Additionally, we ﬁnd conventional meth- ods impractical even if we consider specialized libraries that achieve close to machine IO bandwidth limit (240 Gb=s, or more than 4 seconds to write 1 second of simulation data) [18, 30]. Both the aforementioned traditional workﬂow concerns, i.e., the delay of insight caused by problem redeﬁnition and impracticality of visualization in HPC applications, are ex- acerbated when we consider ensemble simulation; see Fig- ure 1 (b). Catering to an array of pre-process settings and tai- loring post-process steps to each ensemble member is cum- bersome for the practitioner. Additionally, the traditional work- ﬂow delays the initial setup of an ensemble due to the prob- lem redeﬁnition necessary to ﬁgure out design space limits and adjust a given QoI. 2.2 Immersive simulation workﬂow An immersive simulation workﬂow intends to mitigate some of the standard workﬂow inefﬁciencies through live, recon- ﬁgurable simulation. In essence, a new co-process stage in- tegrates steps from both the pre-process, such as problem deﬁnition, and post-process, such as visualization and data- analysis, to be contemporary with the ongoing/live advance- ment of the ﬂow solver; see Figure 2. Co-processing methods aim to vastly reduce the quantity of data storage through in situ visualization that the user is able to interact with and calibrate to their speciﬁc needs [22, 46]. In co-processing, the visualization pipeline is shifted from being a post-processing step within the PDE solver workﬂow to be concurrent with the simulation. The progres- sion of HPC from present petascale to future exascale com- puting, coupled with corresponding increase in simulation complexity, such as multi-scale and multi-physics, further enhance the desire for co-processing as opposed to tradi- tional post-processing capabilities. Immersive simulation couples co-processing that allows live visualization with computational steering that allows live problem redeﬁnition. When immersed, a practitioner is capable not only of visualizing their problem solution in situ, but also of intuitively redeﬁning the problem deﬁnition for rapid exploration of the design space. The practitioner may choose to steer solver parameters, boundary conditions, and even the problem geometry depending on what aspects of the problem space engage them. In the HPC setting, an immersive simulation practitioner can evade the queue time associated with each new problem redeﬁnition. In this man- ner, the behavior of different combinations of input parame- 4 Felix Newberry et al. Idea Geometry Meshing Boundary Conditions Solver Visualization Data Analysis Problem Redeﬁnition? InsightNoYesPre-processPost-process (a) Traditional simulation workﬂow. Idea Geometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-processGeometry Meshing Boundary Conditions Solver Visualization Problem Redeﬁnition? YesPre-processPost-process Data Analysis Problem Redeﬁnition? InsightYes NoPost-process (b) Traditional ensemble simulation workﬂow. Fig. 1: Traditional singular (a) and ensemble (b) simulation workﬂows. Idea Parametric Geometry Parametric Meshing Parametric Boundary Conditions SolverLive, Reconﬁgurable VisualizationLive, Reconﬁgurable Problem Redeﬁnition Data Analysis Problem Redeﬁnition? InsightNoYesPre-processPost-processCo-process Fig. 2: Immersive simulation workﬂow. Software Tools to Enable Immersive Simulation 5 ters can be quickly understood and inform appropriate limits for subsequent or concurrent ensemble simulation. 3 Software tools for immersive simulation In this section we describe the core components of an im- mersive simulation software framework; see Figure 3. This framework enables multiple forms of immersive simulation that work with a traditional PDE solver to perform both live, reconﬁgurable visualization and live, reconﬁgurable prob- lem deﬁnition. The colored boxes indicate software com- ponents while the grey boxes show functionality. Thin ar- rows out of the in situ infrastructure Catalyst , or inter- face SENSEI , represent data extraction/compression, while thin arrows into these components relay analyst input via theParaView GUI. Large data streams, shown by thick ar- rows, correspond to the full solution state (problem mesh and solution). The analyst obtains insight, shown by a wide arrowhead, through engagement in an immersive simula- tion through the ParaView GUI. A crucial aspect of im- mersive simulation is that the commands related by the ana- lyst to reconﬁgure the simulation are conveyed with small data streams. With this setup, an analyst can operate the live visualization and problem deﬁnition for swift simula- tion feedback and insight extraction. We now provide some context on the integration of core immersive simulation soft- ware components SENSEI ,Catalyst , and ParaView with the PDE solver. While we use the ﬂuid solver PHASTA , many aspects of the links between software and all of the resulting advantages are general. We note that SENSEI ,Catalyst , andParaView are all open source products of Kitware. 3.1 SENSEI in situ interface SENSEI is a generic in situ interface best summarized by the underlying premise of “write once, use everywhere” [6]. The interface supports a number of in situ infrastructures such as VisIt /Libsim ,Catalyst , and ADIOS . Simulation output is mapped to the VTK data model via a data adaptor and an analysis adaptor, while an in situ bridge links the respective adpators and prompts in situ analysis. An in depth descrip- tion of SENSEI is provided in [6] as well as a brief view from the simulation code, and the perspective of different choices ofin situ infrastructure. In the SENSEI repository mini-apps are provided as a guide for different PDE solver applica- tions. The mini-app topics span SENSEI Python bindings, C++solvers for dynamic periodic oscillators, the Mandel- brot set and vortex simulation, and the PHASTA mini-app that we build upon for this work. In general the linking software to enable immersive sim- ulation appears in three places. First, the SENSEI mini-appbridge must be updated to match any parameters the practi- tioner seeks to pass from their PDE solver. As shown in [6], the bridge is custom to the implementation code and requires initialization, analysis, and ﬁnalization steps. In this work, we focus on how this bridge communicates with the PDE solver from the solver’s view, and also on the practical im- plementation in terms of necessary scripts at runtime. From the perspective of the PHASTA solver the abbreviated sim- ulation code to highlight the SENSEI footprint is presented in Figure 4. Within the simulation script we make calls to subroutines that initialize, co-process and ﬁnalize the inter- action of the interface with the solver. These subroutines are deﬁned in a second PDE solver script sensei interface.f shown in Figure 5. Figure 5 essen- tially reﬁnes the main simulation routine calls to SENSEI by providing additional context. The routines called here are deﬁned in the SENSEI PHASTA bridge. For instance, the ini- tialization subroutine calls any in situ analysis adaptor, such asCatalyst , before communicating mesh details vital to live visualization and a set of ﬁelds determined by the user. For the full context of these scripts we refer the inter- ested reader to our GitHub repository www.github.com/ SimNautilus/nautilusFlow that documents the links be- tween the PDE solver PHASTA andin situ interface SENSEI . In our implementation of an immersive simulation frame- work we initially setup the main solution variables, such as mesh coordinates, velocity and pressure ﬁelds, for live visu- alization. Next, we improve the framework with the intro- duction of steering parameters or more complex ﬁelds that are passed to SENSEI . Having setup the PDE solver to communicate ﬁelds to SENSEI , we next address the implementation at runtime. As noted, several in situ infrastructures are supported. The prac- titioner can setup the desired analysis type with ease in the equivalent of senseiPhasta.xml ; see Figure 6. This script is called via sensei adaptors init() during the SENSEI initial- ization steps in Figure 5. While the SENSEI bridge is created during initialization, the analysis it facilities, for instance via Catalyst , is reconﬁgurable during a simulation run. 3.2 Catalyst in situ infrastructure Catalyst , formally the ParaView co-processing library [15], is an in situ visualization library with an adaptable API that is built on VTK[40] and ParaView [1]. Through building the library with VTK,Catalyst can access a plethora of use- ful algorithms such as IO writers, visualization ﬁlters, and graphics rendering. The visualization pipelines can be de- scribed in C ++or Python, and the ParaView GUI facilitates the generation of Catalyst co-processing Python scripts both directly, and via a Python trace and shell. Additionally, Catalyst is designed to connect visualization pipelines re- motely through server-client architecture. While the addi- 6 Felix Newberry et al. Analyst ParaView Graphical User Interface Catalyst In Situ InfrastructureSENSEI In Situ InterfaceLive, Reconﬁgurable Problem DeﬁnitionLive, Reconﬁgurable Visualization PDE SolverLarge Data Streams Small Data Streams Large Insight Streams Fig. 3: Immersive simulation software framework. s u b r o u t i n e i t r d r v ( . . . ) c a l l p h s e n s e i i n i t ( . . . ) ! advance s i m u l a t i o n c a l l p hs e n s e i c o p r o c e s s ( . . . ) ! advance s i m u l a t i o n c a l l p h s e n s e i f i n i ( ) ! c l e a n u p end Fig. 4: Pseudocode of the SENSEI footprint in PDE solver simulation script. s u b r o u t i n e p hs e n s e i c o p r o c e s s ( . . . ) ! read i n f i e l d s from s o l v e r c a l l s e n s e i c o p r o c e s s ( . . . ) end s u b r o u t i n e p h s e n s e i i n i t ( . . . ) c a l l s e n s e i a d a p t o r s i n i t ( ) c a l l c r e a t e p o i n t s a n d a l l o c a t e c e l l s ( . . . ) ! i t e r a t e PDE s o l v e r p a r t i t i o n s do ib l o c k =1 , n b l o c k s ! c o u n t v e r t i c e s per element , ! e l e m e n t s i n block , ! and pass t o SENSEI c a l l i n s e r t b l o c k o f c e l l s ( . . . ) enddo c a l l a d d f i e l d s ( . . . ) end s u b r o u t i n e p h s e n s e i f i n i ( ) c a l l s e n s e i a d a p t o r s f i n i ( ) end Fig. 5: Pseudocode of the PDE solver link between simula- tion script in Figure 4 and SENSEI . tion of the Catalyst adapter to solvers must to some degree inhibit solver performance, this inﬂuence has been shown to be negligible when compared to the savings in ﬁle I/O and data management [15]. In our numerical example we implement Catalyst with the ﬂuid solver PHASTA via the<s e n s e i > <! A s i m p l e C a t a l y s t A n a l y s i s > <a n a l y s i s t y p e =” c a t a l y s t ” p i p e l i n e =” s l i c e ” a r r a y =” p r e s s u r e ” a s s o c i a t i o n =” p o i n t ” imagef i l e n a m e =” p r e s s u r e %t s . png ” imagewidth =” 1920 ” imageh e i g h t =” 1080 ” s l i c enormal =” 0 , 0 , 1 ” e n a b l e d =” 0 ” / > <! C a t a l y s t cop r o c e s s s c r i p t > <a n a l y s i s t y p e =” c a t a l y s t ” p i p e l i n e =” p y t h o n s c r i p t ” f i l e n a m e =” / p a t h t of i l e / m y c a t a l y s t p i p e l i n e . py ” e n a b l e d =” 1 ” / > <s e n s e i > Fig. 6: Pseudocode of the runtime XML description that conﬁgures in situ analysis routines [6]. SENSEI in situ interface, speciﬁcally through analysis adap- tor in Figure 6. A powerful function of Catalyst that has been alluded to is the setup of ParaView pipelines. The entire pipeline script can be generated intuitively from a ParaView GUI session that writes out a Catalyst co-processing script as shown in Figure 7. The co-processing script can be tailored to speciﬁc problems and is loaded automatically when live visualization is launched to facilitate the user’s analysis. The practitioner can update the pipeline live through the Catalyst connection in the ParaView GUI. Having described the implementation of the interface SENSEI and infrastructure Catalyst we are equipped with an immersive simulation framework of Figure 3. With these immersive simulation tools in hand a practitioner can lever- age interactive visualization and problem redeﬁnition for ef- ﬁcient design space exploration. Software Tools to Enable Immersive Simulation 7 # CoProcessor d e f i n i t i o n  def C r e a t e C o P r o c e s s o r ( ) : def C r e a t e P i p e l i n e ( . . . ) : c l a s s P i p e l i n e : # ParaView p i p e l i n e # Sequence o f f i l t e r s return P i p e l i n e ( ) # P r o c e s s i n g method  def DoCoProcessing ( . . . ) : # Update c o p r o c e s s o r w i t h newly # g e n e r a t e d s i m u l a t i o n data # I f t h e p i p e l i n e hasn ’ t been s e t u p # yet , t h i s w i l l s e t i t up . c o p r o c e s s o r . U p d a t e P r o d u c e r s ( . . . ) # W r i t e o u t p u t data , # i f a p p r o p r i a t e . c o p r o c e s s o r . WriteData ( . . . ) ; # W r i t e image c aptur e , # i f a p p r o p r i a t e . c o p r o c e s s o r . WriteImages ( . . . ) # L i v e V i s u a l i z a t i o n , i f e n a b l e d . c o p r o c e s s o r . D o L i v e V i s u a l i z a t i o n ( . . . ) Fig. 7: Pseudocode of the Catalyst co-processing script that enables in situ visualization and loads a chosen ParaView pipeline [4]. 4 Numerical demonstration: immersive simulation tools We demonstrate the utility of immersive simulation tools on a 2D aggressive subsonic diffuser simulated in PHASTA [50]. In this section we ﬁrst describe the diffuser model equipped with parametric geometry deformation and then illustrate the application of immersive simulation tools described in Section 3. Code examples of the implementation of immer- sive simulation tools in PHASTA can be found in our GitHub repository www.github.com/SimNautilus/nautilusFlow . 4.1 2D aggressive subsonic diffuser Aircraft frequently require air intake systems that reroute air from the free stream velocity and slow it to speeds appropri- ate for the engine [3]. The intake arrives at the compres- sor on the aerodynamic interface plane (AIP). It is desirable that the AIP air exhibits low swirl, low distortion, and high pressure recovery, with minimal increase of the intake drag, lest engine performance is inhibited [32]. Therefore, under- standing of intake ﬂow is imperative to effective simulation of engine performance. Typical intake geometries consist of a duct to direct the ﬂow connected to a diffuser that, via an increase in cross- sectional area, trades upstream kinetic energy for a rise in downstream static pressure. This increase in area producesan adverse pressure gradient that encourages separated ﬂow. Furthermore, the ﬂow dynamics may be complicated by duct curvature that introduces spanwise pressure gradients. The diffuser model we study is motivated by a 3D tran- sonic diffuser with two signiﬁcant simpliﬁcations made in order to reduce computational expense: the model is a 2D slice of the 3D geometry midplane and the ﬂow is treated as incompressible as opposed to compressible. Shifting from 3D to 2D greatly reduces the number of elements in the simulation mesh, and treating the ﬂow as incompressible permits a much larger time step while maintaining a con- verging simulation. We emphasize the severity of the model simpliﬁcations on the diffuser; the purpose of this simulation is to demonstrate an immersive simulation software frame- work and to capture quantitatively the inﬂuence of upper and lower blowers on core ﬂow behavior. The diffuser geometry, depicted in Figure 8, is consid- ered compact with a length to diameter ratio of L=D1:03 and a high expansion ratio of ER=2:2. Active ﬂow control is implemented through two tangential blowers, denoted as the upper blower (UB) and lower blower (LB). The blow- ers are designed to inhibit or delay ﬂow separation on their respective surfaces. We prescribe a trapezoid waveform in- ﬂow condition to both the UB and LB, where we set the wave’s mean, amplitude, and period. The time the trapezoid wave spends rising, falling, and at the maximum value of the wave is set to 1 :5e3 s, 1:5e3 % and 4e3s, respectively. We have ﬁxed the period to 8e 3sfor the study and ensured a temporal resolution of 80 time steps for each blower period with a time step of dt=1e4s. The incompressible ﬂow through the diffuser is modeled with an inﬂow Mach number M¥0:052 and bulk Reynolds number Reb4e4. This Reynolds number is similar to the 2D plane diffuser experiments [34, 8]. The computational mesh made up of 241k tetrahedral elements is shown in Fig- ure 8. The mesh is only one element deep (into the page as displayed in Figure 8), and periodic boundary conditions are employed in the depth direction to model 2D ﬂow. Every boundary layer has ﬁrst point off the wall is set to 2e 6 m and a stretching ratio of 1 :25 to maintain y+<1 at the wall boundaries. The mesh is further reﬁned at the entrance of the blowers and downstream to resolve the recirculating region. In this study, we examine the QoI measuring the recir- culating region length deﬁned as `r=rZ ¯ux<0dx d wall; (1) where ¯ uxdenotes the time-averaged stream-wise velocity anddwallthe distance from the wall. We integrate with re- spect to dwallas opposed to dyto penalize recirculating re- gions that extend further into the core body of the ﬂow. In the case where separation occurs in both the lower and up- 8 Felix Newberry et al. (a) (b) Fig. 8: 2D aggressive subsonic diffuser geometry (a) entire geometry, (b) close up of the blower control region showing the UB and LB. per walls, we add the respective integrals together prior to taking the square root to ﬁnd the length. 4.1.1 Parametric geometry deformation In order to explore the effect of UB location along the curved upper extent of the diffuser, a custom parametric mesh mod- iﬁcation tool was developed which slides the UB tangen- tially along the curved surface, moving the mesh vertices smoothly, allowing for geometric design exploration in situ without the need for re-meshing or multiple hand-generated CAD conﬁgurations. This procedure relies on Sederberg’s freeform deformation algorithm [41] with a geometry-ﬁtted collection of 493 background non-uniform rational B-Spline (NURBS) patches deformed in unison with a custom geom- etry management algorithm. The geometry modiﬁcation is driven by a single parameter jvarying continuously from 0 to 1 representing locations far down the curve of the diffuser and past midway up the diffuser as depicted in Figure 9 (b). These extents were chosen to limit mesh distortion under ge- ometry modiﬁcation. Figure 9 (a) illustrates the collection of background NURBS patches already deformed to the j=1 extent. This geometry deformation algorithm was custom de- signed for the diffuser problem, employing low level Open- CASCADE commands. The general procedure here is extensi- ble to other 2D geometric considerations. Further, this pro- cedure is a proof of concept for how geometric variables can be included in immersive simulation design space ex- ploration procedures. Future improvement to this procedure may incorporate higher level abstractions for parametric ge- ometry modiﬁcation such as the programmable CAD system Engineering Sketchpad [21]. More general tools for geomet- ric modiﬁcation to HPC applications are presented in [49].4.2 Computational steering We employ Catalyst ,SENSEI , and the ParaView GUI to implement live visualization and computational steering on the 2D diffuser. The parameters for which we enable steer- ing are the mean and amplitude of the UB and LB velocities, the UB position parameterized between j=0 and j=1, as shown in Figure 9 (b), and the weight term w, useful for as- sessing time averaged quantities. The QoI we are interested in, the recirculating region deﬁned in Section 4.1, is mea- sured as a time average to account for the transient boundary conditions on the UB and LB. We update the time averaged stream-wise velocity with ¯ux=wux+(1w)¯ut0x; where uxis the present time step’s stream-wise velocity and ¯ut0xdenotes the previous time step’s averaged stream-wise velocity. Increasing wweighs recent steps more heavily and “shortens” the time-average interval, while decreasing wre- duces the weight of new time steps and “lengthens” the time- average interval. In general, we set w=0:01. In this manner, we can calibrate the memory of the time-averaged velocity ﬁeld. To equip PHASTA with the steering parameters the cho- sen ﬁelds must be passed from the PHASTA solver to SENSEI , and the SENSEI PHASTA adaptor updated for these ﬁelds as described in Section 3. Once implemented, the steering pa- rameters can be interacted with during a Catalyst live vi- sualization; see Figure 10. The user can adjust a given steer- ing parameter, for instance the UB position, and immedi- ately see how this alteration impacts the ﬂow behavior. In Figure 10 we see a complex Catalyst VTK pipeline that determines regions in the diffuser ﬂow with negative time averaged stream-wise velocity and calculates the associated recirculation length. Note, the practitioner’s interaction with theParaView GUI allows them to create, visualize, and ex- plore different VTKpipelines and easily examine other QoIs. In this manner, the practitioner can directly observe the ef- fect of variations in key parameters on their QoIs, gain in- Software Tools to Enable Immersive Simulation 9 (a) (b) Fig. 9: Parametric 2D aggressive subsonic diffuser geometry: (a) spline track used to guide the deformation of the background NURBS surfaces (b) close up of the blower control region showing blower position set to j=0 and j=1. sight more rapidly, and do so without saving the visualized ﬂow data to disk. We leverage computational steering to clarify QoI se- lection and setup design space limits for subsequent ensem- ble simulation. Our recirculation length QoI, deﬁned in Sec- tion 4.1, initially accounted for the upper recirculation area alone and integrated over dyrather than dwall. Computa- tional steering allows us to observe the lower separation area and trial integration with respect to distance from the wall. Live adjustments of the blower control steering parameters suggest that low values of UB velocity mean, i.e., <1:43% total diffuser mass ﬂow rate, produce negligible reduction in the upper recirculation area. High values, >2:85%, result in fully attached ﬂow, indicate additional velocity would be a poor use of ﬂow mass. Similar testing of the LB mean produces an interval of [0:00;2:34]. The LB mean is use- ful at small percentages too because the diffuser geometry means that lower separation occurs only when the core ﬂow is pulled up by large UB action. This tug of war effect can be seen in Figure 10. Finally, we decide to examine the interval [0;1]for blower position j. While a higher blower position improves this QoI, our study here is 2D and incompressible, hence we stay conservative in our design space constraints. 5 Software tools for ensemble sensitivity analysis In Section 4.2, we employed immersive simulation to ex- plore the design space of the 2D diffuser problem and to determine a QoI. Our next objective is to leverage this infor- mation to spawn an ensemble simulation. We ﬁrst, in Section 5.1, provide background on global sensitivity analysis. Next, in Sections 5.2 and 5.3 we present software useful for efﬁ- cient ensemble generation and the application of UQ tech- niques, respectively. 5.1 Sensitivity analysis Parameter sensitivity is central to achieving an optimal de- sign and performing UQ. Through visualization of a solu-tion subject to deviations from nominal parameter values in boundary conditions, material properties, or geometric pa- rameters, a user can navigate the design space with greater efﬁcacy. In variance-based sensitivity analysis, the variance of a given QoI is decomposed into fractions that are at- tributed to speciﬁc input parameters, or interactions of those parameters. Monte Carlo simulation is a popular candidate for obtaining sensitivity metrics but can become infeasible for computationally demanding models. A suitable alterna- tive is to use a UQ technique to approximate the QoI with an expansion in multivariate orthogonal polynomials, known as the polynomial chaos (PC) expansion [19, 53]. We next de- scribe PC expansions followed by their application to global sensitivity analysis. We note that our PC expansion is de- scribed with random input variables general to the UQ ﬁeld. Our implementation of the PC expansion treats design vari- ables as uniform random variables. 5.1.1 Polynomial chaos expansions We consider the scalar QoI, u(XXX), assumed to have ﬁnite variance, as a function of a d-dimensional vector of random inputs XXX:= (X1;:::;Xd)with joint probability density func- tion f(xxx). The PC expansion approximates the QoI as u(XXX) =¥ å j=1cjyj(XXX); where yj(XXX)is a multivariate orthogonal polynomial eval- uated at the random inputs and weighted by deterministic coefﬁcients cj. The polynomials yj(XXX)are chosen to be or- thogonal with respect to the probability measure f(xxx). For instance, if XXXfollows a jointly uniform or Gaussian distri- bution, then yj(XXX)are multivariate Legendre or Hermite polynomials, respectively [53]. We assume yj(XXX)are nor- malized such that E[y2 j(XXX)] = 1, where E[]represents the mathematical expectation operator. The expansion is trun- cated to a ﬁnite number of terms as u(XXX)P å j=1cjyj(XXX): (2) 10 Felix Newberry et al. Fig. 10: ParaView GUI displaying the Catalyst live visualization and SENSEI computational steering. On the left, the pipeline browser shows the visualization pipeline and steering parameters. An expansion with total order pand dimension dhasP= (p+d)! p!d!terms. As P!¥, for a sufﬁciently smooth u(XXX), the PC expansion converges in the mean-square sense to u. The PC coefﬁcients can be used to determine QoI statistics, construct a surrogate model, or perform sensitivity analy- sis. For instance, the QoI mean and variance can be com- puted from the PC expansion coefﬁcients as mPC=c1and s2 PC=åP j=2c2 j, respectively. To identify PC coefﬁcients c= (c1;:::; cP)Tin (2) we generate realizations of XXXdenoted byxxx. We seek to solve for cin the linear system from N independent samples of xxxas uYYYc; (3) where YYY(i;j):=yj(xxxi)andu:= (u(xxx1);:::; u(xxxN))T: (4) If the resulting regression problem is over-determined, with N>P, we employ least squares approximation [23, 20], whereas if it is under-determined, with N<P, we apply compressed sensing [14, 35, 24, 13]. Compressed sensing methods are advantageous because they reduce the number of sample evaluations to determine PC coefﬁcients. 5.1.2 Global sensitivity analysis via PC expansions The ﬁrst order Sobol’ indices, SXk, measure the effect of varying one parameter alone, Xk, averaged over the vari- ations in other parameters. Total indices, ST Xk, refer to the variation in the QoI accounted for by one parameter and any interactions it has with other input parameters. In this work,we apply global sensitivity analysis analytically as a post- processing of PC expansion coefﬁcients [42]. The ﬁrst order Sobol’ indices are given by SXk=å i2Ikc2 i s2 PC; (5) where, Ikis a set of indices for which yj(XXX),j=1;:::; P, is a function of only Xk. Similarly, the total Sobol’ indices are given by ST Xk=å i2IT kc2 i s2 PC; (6) where, IT kis a set of indices for which yj(XXX),j=1;:::; P, is a function of Xkand any other inputs. 5.2 Ensemble generation via libEnsemble In this work, we employ libEnsemble , an open source Py- thon library developed at Argonne National Laboratory as part of the Department of Energy Exascale Computing Project to coordinate the concurrent evaluation of dynamic ensem- bles of calculations on massively parallel resources [25]. libEnsemble aims to achieve extreme scaling, resilience, task monitoring and resource recovery, portability and ﬂexi- bility, and exploit persistent data ﬂow. The ensemble is coor- dinated via a manager/worker scheme that operates on sev- eral communication protocols. Software Tools to Enable Immersive Simulation 11 libEnsemble is implemented through four components. The generator function gen fproduces values for simula- tions, the simulator function sim fperforms a simulation for given values from gen f, the allocation function determines whether gen forsim fshould be called, and the calling script deﬁnes the parameters for the previous three functions and executes libEnsemble . The libEnsemble user manual [25] provides more detail on each of these components, and var- ious tutorials. In Figure 11 we show an abbreviated version of the calling script. Note, in our application we have no generator function because sim freads inputs directly from a ﬁle. In the calling script we specify the number of work- ers, and the number of simulations. For example, an HPC job may use 50 nodes in total, but require 2 nodes per sim- ulation. In this case, the number of workers would be 25. In Figure 11 we ﬁrst setup the PHASTA executable. Next we create an ID number for each simulation before deﬁn- ing parameters for the simulation and generation functions. An exit criteria is speciﬁed to terminate the ensemble when the maximum number of simulations is reached. Finally, the ensemble is launched. As simulations complete and work- ers become available new simulations are launched. Other ensemble manager tools are also available, for instance, UQ Pipeline [43]. 5.3 Sensitivity analysis with chaospy We use the Python library chaospy to perform UQ PC based sensitivity analysis [16]. More generally, chaospy includes a suite of UQ methods such as low-discrepancy sampling, quadrature creation, polynomial manipulations, among oth- ers. A collection of Jupyter notebooks lets users explore the variety of chaospy applications. We present pseudocode of the chaospy enabled sensi- tivity analysis in Figure 12. The main idea is to construct a PC expansion and use the coefﬁcients to calculate Sobol indices following Sections 5.1.1 and 5.1.2, respectively. We note that in this pseudocode we do not include the validation error. In practice, we recommend reserving 20% of sam- ples to compute the PC expansion validation error, and use this to calibrate the polynomial order p. 6 Numerical demonstration: ensemble simulation tools We have now, courtesy of an immersive simulation frame- work, established design space limits and conﬁdence in our choice of QoI. We are also equipped with software tools to perform ensemble simulation and analyze the results. In this section, we demonstrate the utility of Catalyst during the collection of ensemble data, and carry out a global sensitiv- ity analysis on the 2D agressive subsonic diffuser.import numpy as np # Import l i b E n s e m b l e c o n t e n t from l i b e n s e m b l e . l i b E import l i b E from l i b e n s e m b l e . e x e c u t o r s . m p i e x e c u t o r import MPIExecutor from l i b e n s e m b l e . t o o l s import p a r s e a r g s # Import p h a s t a s i m f and a l l o c f from p h a s t a s i m import p h a s t a e v a l u a t e as s i m f from a l l o c a t i o n f u n c import p h a t a a l l o c as a l l o c f nworkers = 25 nsim = 50 sim app = p a t h t of i l e / p h a s t a . exe e x c t r = MPIExecutor ( c e n t r a l m o d e = True ) e x c t r . r e g i s t e r c a l c ( f u l l p a t h =sim app , c a l c t y p e = ’ sim ’ ) # S e t up H0used as run ID H0 = np . z e r o s ( n sim , d t y p e = [ ( ’ x ’ , f l o a t ) , ( ’ s i m i d ’ , i n t) ] ) H0 [ ’ s i m i d ’ ] = range ( nsim ) l i b E s p e c s =f’ nworkers ’ : nworkers , ’comms ’ : ’ l o c a l ’ g s i m s p e c s =f’ s i m f ’ : sim f , ’ i n ’ : [ ’ x ’ ] , ’ o u t ’ : [ ( ’ y ’ , f l o a t ) , ( ’ c a l c s t a t u s ’ , s t r) ]g g e n s p e c s =fg a l l o c s p e c s =f’ a l l o c f ’ : a l l o c f , ’ o u t ’ : [ ( ’ x ’ , f l o a t ) ]g e x i t c r i t e r i a =f’ sim max ’ : l e n( H0)g H, p e r s i s i n f o , f l a g = l i b E ( s i m s p e c s , g e n s p e c s , e x i t c r i t e r i a , a l l o c s p e c s = a l l o c s p e c s , l i b E s p e c s = l i b E s p e c s , H0=H0 ) Fig. 11: Pseudocode of the libEnsemble calling script. 6.1 Catalyst pipeline QoI extraction TheCatalyst library is successfully interfaced with PHA- STAviaSENSEI to enable real-time QoI extraction. First, the PHASTA adaptor ﬁle is edited to include additional ﬁelds that facilitate QoI computation. This adaptor ﬁle determines the ﬁelds that are passed from the PHASTA solver to ParaView during a simulation. Next, PHASTA is built with SENSEI en- abled and the updated adaptor ﬁle. Finally, a ParaView pipe- line is saved to a Python co-processing script that is read dur- ing the PHASTA run at a speciﬁed frequency. The execution of the Catalyst co-processing step results in some com- putational overhead which we assess via a strong scalability study. 12 Felix Newberry et al. # S e t u p d i s t r i b u t i o n w i t h s t o c h a s t i c # d i m e n s i o n d v a r i a b l e s d i s t r i b u t i o n = cp . I i d ( cp . Uniform (1 ,1) , d ) # C o n s t r u c t PC e x p a n s i o n p o l y n o m i a l o r d e r = 4 p o l y n o m i a l e x p a n s i o n = cp . o r t h t t r ( p o l y n o m i a l o r d e r , d i s t r i b u t i o n ) # E v a l u a t e e x p a n s i o n a t samples t o # c o n s t r u c t p s i p s i = p o l y n o m i a l e x p a n s i o n (samples ) . T # S o l v e f o r pc c o e f f i c i e n t s v i a l e a s t # s q u a r e s or compressed s e n s i n g sigma = 0 . 0 0 1 # s p g l 1 s o l v e r t o l e r a n c e ch a t , r e s i d , grad , i n f o = spg . spg bpdn ( p s i , e v a l u a t i o n s , sigma ) # C a l c u l a t e Sobol i n d i c e s from PC e x p a n s i o n m o d e l a p p r o x i m a t i o n = np . sum ( ch a t s p gp o l y n o m i a l e x p a n s i o n ) # F i r s t o r d e r s o b o l i n d i c e s1 = cp . Sens m ( m o d e l a p p r o x i m a t i o n , d i s t r i b u t i o n ) # T o t a l s o b o l i n d i c e st = cp . S e n s t ( m o d e l a p p r o x i m a t i o n , d i s t r i b u t i o n ) si n t e r a c t i o n s = s ts1 Fig. 12: Pseudocode of the chaospy sensitivity analysis. Samples is an array of samples fxxxigN i=1, and evaluations are the corresponding simulated QoI ufrom equation (4). We performed a strong scalability study on the Argonne Leadership Computing Facility’s Theta computational re- source with a2572K element mesh of the diffuser, par- titioned to run on 6, 12, 24 and 48 cores, with 1 MPI process per core. Theta is a 11 :69-petaﬂops Intel-Cray supercom- puter. Each simulation computes 25 time steps, and saves output data at 2 step intervals. In Figure 13 we plot the scal- ability study timing results for no write, conventional write andCatalyst write cases. The results are computed as the average of three repetitions. No write denotes a PHASTA run that neither saves data to restarts or via the Catalyst pipe- line, conventional write indicates that a restart ﬁle contain- ing the entire ﬂow solution is written every 2 time steps and Catalyst write employs the Catalyst QoI pipeline to ef- ﬁciently extract data every 2 time steps. In Figure 13 (a) we observe that the performance of the varying write methods and the optimal scaling curve in terms of the wall-clock time are indistinguishable. In Figure 13 (b)we plot the strong scalability factor, deﬁned as s f=tre fncre f tnc; (7) where tandncare a given execution time and number of cores, and tre fandncre fare the execution time and num- ber of cores for the 6 core, no write, reference case. Optimal scaling would be a horizontal line with a scaling factor of 1. We observe in Figure 13 (a) that for all cases the paral- lel overhead increases with number of cores. Importantly, while Catalyst write is marginally less efﬁcient than not writing at all, its performance is comparable to the conven- tional write method. In the traditional workﬂow approach, where the QoI cal- culation is performed as a post-processing step, the entire solution ﬁeld must be saved at a given time step to a restart ﬁle of approximately 205 Mb. In the Catalyst write ap- proach, the ﬁle containing the QoI is just 4 Kb, indicating a storage reduction of approximately ﬁfty thousand times. These results demonstrate that the signiﬁcant memory ad- vantages of writing output data via the Catalyst pipeline come at minimal computational overhead. 6.2 Ensemble generation We have enabled the straightforward evaluation of ensem- bles of PHASTA simulations through the implementation of thelibEnsemble Python library described in Section 5.2 [25]. libEnsemble facilitates two levels of parallelism. In the ﬁrst, an individual PHASTA simulation may employ nu- merous processors, in this demonstration we partition a 241 k element mesh and conduct each simulation with 24 pro- cessors. In the second, libEnsemble manages individual PHASTA evaluations and can launch new simulations imme- diately after previous simulations are ﬁnalized and resource becomes available. The user selects or provides functions to generate simulation input, carry out a simulation, and man- age the ensemble of simulations. This setup allows the user to tailor libEnsemble to their speciﬁc solver requirements. For instance, to obtain an ensemble of PHASTA simulations we create a custom simulation script that for a given input realization creates a run directory with the necessary geom- etry and restart ﬁles, reads the inputs from ﬁle, writes the inputs to the solver.inp script and executes PHASTA . We have developed and documented a set of PHASTA speciﬁc libEnsemble scripts and that enable the parallel evaluation of ensemble members, allowing much greater ﬂex- ibility in job submission. libEnsemble with PHASTA has been implemented on local resources, and additionally on the Argonne Leadership Computing Facility’s Cooley clus- ter. Cooley has 126 compute nodes, each with 12 cores and 1 NVIDIA Tesla K80 GPU. Software Tools to Enable Immersive Simulation 13 (a) The wall-clock time. (b) Strong scaling factor (7). Fig. 13: Strong scalability study comparing the performance of no write, where no output data is saved, conventional write, where output data is written to restart ﬁles and Catalyst write where PHASTA is equipped with the Catalyst pipeline for efﬁcient QoI extraction. In each case we plot the average of three repetitions of a 25 time step simulation with output data saved at every second step. 13 (b) presents the corresponding strong scalability factor 6.3 Global sensitivity analysis We apply global sensitivity analysis via PC expansions, de- tailed in Section 5.1.2, to the 2D aggressive subsonic dif- fuser. Our objective is to understand the inﬂuence of the UB and LB velocities and the UB position on the length of the recirculation area (1). While the UB and LB provide us with ﬂow control, the mass ﬂow rate diverted to these blowers in- creases intake drag of the entire diffuser. Therefore, we want to determine a blower conﬁguration that reduces the recircu- lation area while expending minimal diverted ﬂow. To explore the design space we treat the UB and LB velocities and amplitudes and the UB position as uniform random variables. The variables and their respective inter- vals are summarized in Table 1. The limits were determined through computational steering via in situ exploration of the parameter space. Figure 14 shows the PC based Sobol in- dices constructed from 495 samples using polynomial order p=4. The PC validation error was 10 :1%. In Figure 14 we see that the most inﬂuential parameter is the UB velocity mean, while the LB velocity amplitude has minimal effect on the QoI variance. The UB position jalso plays a signiﬁ- cant role. Further, it is apparent that the recirculation length is sensitive to the interactions between all the parameters. This matches our expectation given the tug of war dynamic between the upper and lower blower’s over the core ﬂow. We plot in Figure 15 the entire ensemble in terms of mass cost and recirculation length, colored by the sampled input parameters. We exclude AL(see Table 1) from Fig- ure 15 due to its small sensitivity found in Figure 14. Mass cost is computed simply as the sum of the UB and LB meanParameter Symbol Distribution UB velocity mean mU U[1:43;2:85]% UB velocity amplitude AU U[0:00;2:13]% LB velocity mean mL U[0:00;2:34]% LB velocity amplitude AL U[0:00;1:17]% UB Position j U[0;1] Table 1: Upper blower (UB) and lower blower (LB) input parameters. Velocity metrics are given as a percentage of the total diffuser mass ﬂow rate. velocities. Consequently, it is no surprise that in Figures 15 (a) and (c) we ﬁnd that increasing the mean mass ﬂow rates is correlated with higher mass cost. However, this increase in mass is correlated with our aim of reducing the recirculation length. The best performing samples are those that produce the smallest recirculation length while expending the least mass. In general, these samples correspond to higher blower positions as illustrated in Figure 15 (d). While Figure 15 (d) implies that moving the blower position even further down- stream may achieve better results, we refrain from doing this given the severe simpliﬁcations of this study, namely incom- pressible and 2D ﬂow. Future work could transition exami- nation of the blower position into 3D and compressible ﬂow. We select a sample, identiﬁed by a red circle in Figure 15, as ﬂow conﬁguration that yields low recirculation length and minimizes mass cost. This conﬁguration reduces the re- circulation length to 0, at a mass cost of approximately 4 :5% the core mass ﬂow. In Figure 16 we compare our chosen ﬂow conﬁguration to the diffuser no active ﬂow control. The 14 Felix Newberry et al. Fig. 14: Global sensitivity analysis of the 2D aggressive sub- sonic diffuser. This study was carried out in Python with chaospy and an order p=4 PC expansion from 495 sam- ples with validation error error . The sum of the interactions and ﬁrst order is the total order. time averaged stream-wise velocity is plotted. We observe in Figure 16 (a) that the recirculation area, outlined in black, is very large. In Figure 16 (b) the core ﬂow has been bent upwards by the action of the UB and the recirculation area is reduced to 0. 7 Conclusions In this work we introduce an immersive simulation software framework that enables practitioners to evade many of the typical workﬂow inefﬁciencies in solving systems of PDEs. We describe the in situ interface SENSEI and infrastructure Catalyst and provide software linkages to show their inter- action with the PDE solver. In addition we demonstrate the use of complementary ensemble simulation software tools libEnsemble , for intuitive evaluation of simulation ensem- bles, and chaospy , for straightforward implementation of UQ techniques. We demonstrate the efﬁcacy of this immersive simula- tion framework on a 2D diffuser problem in the PHASTA . With SENSEI andCatalyst we implement computational steering on the diffuser problem, including parametric ge- ometry deformation, that allows for rapid design space ex- ploration and eliminates the need of writing data to disk. We show that using a Catalyst pipeline to save output QoI requires minimal computational overhead, scales well, and reduces memory demands by approximately ﬁfty thousand times. A global sensitivity analysis integrates all the immer- sive simulation components together. Computational steer-ing allows the practitioner to test input parameters’ inﬂu- ence, guide the selection of appropriate parameter limits and determine a useful QoI. The ensemble evaluation with lib- Ensemble is fast, intuitive, and takes full advantage of mas- sively parallel resources. Moreover, the saving of output QoI with pre-constructed Catalyst pipeline drastically reduces the total memory footprint. Finally, we use chaospy to cal- culate Sobol’ indices and determine the optimal diffuser con- ﬁguration. Our future research includes a plethora of advancements to immersive simulation capabilities. While the application of UQ techniques as a post-processing analysis is useful, the ability to establish sample limits during computational steer- ing and launch ensembles interactively from within Para- View GUI would lead to yet more rapid insight extraction. In this work we provide a proof of concept for parametric geometry deformation that is designed speciﬁcally for this problem’s geometry. More general geometric deformation tools are an interesting line of inquiry. Acknowledgments The work was supported by National Science Foundation grant 1740330 (SSE) and Department of Energy grant DE- SC0021411 (FASTMath SciDAC). The SENSEI enhancements developed by Kitware were supported by the Ofﬁce of Science, Ofﬁce of Advanced Sci- entiﬁc Computing Research, of the U.S. Department of En- ergy under a subcontract from Lawrence Berkely National Laboratory No. 7566608, through the grant “Scalable Anal- ysis Methods and In Situ Infrastructure for Extreme Scale Knowledge Discovery.”