POD-Based Model-Order Reduction for Discontinuous Parameters

: Reduced-order models (ROMs) based on proper orthogonal decomposition (POD) are widely used in industry. Due to the rigid requirements on the input data, these methods struggle with discontinuous parameters, e.g., optional rear spoiler on a car. In order to also include these types of parameters, a new method is presented that splits the full-order model (FOM) domain with its discontinuous parameters into multiple ROM subdomains. The resulting subdomains then again comply with the ROM requirements, and the established and proven ROM methods can be applied. The steps involved in computing a ROM based on the proposed method, by setting up the subdomains, mapping the FOM data into the domains, as well as computing the ROMs on the domains, are shown in detail in this paper. The method is employed on two use cases. The academic one-dimensional use case focuses on how the steps involved are employed and analyzes the introduced errors. The second use case’s FOM is based on the DrivAer body with an optional rear spoiler computed using computational ﬂuid dynamics (CFD) and demonstrates the usage in an industrial environment.


Introduction
Reduced-order models (ROMs) nowadays play an important role in industry and significantly improve how full-order model (FOM) results are used with respect to optimization problems, merging experimental and numerical data as well as within an interactive design process [1][2][3]. By computing FOM solution snapshots for various parameter combinations in the so-called offline phase, these results can be evaluated and exploited in near real time during the online phase. In this paper, the order reduction is performed using proper orthogonal decomposition (POD) [4] coupled with a radial basis function (RBF) [5] interpolation method (POD+I) for predictions of unknown parameter combinations. The main FOM input in mind is a computational fluid dynamic (CFD) solution.
Utilizing POD+I, especially for industrial use cases, is often not straight forward, since many aspects of the POD computation needs to be taken into account. Due to the POD requirement of conform meshes across the input snapshots, these methods struggle with: • Discontinuous parameters. • Complex geometry transformations. • Nonlinear effects across the parameter domain often times not being modeled accurately. In this paper, a new approach in dealing with the aforementioned shortcomings is introduced. The proposed method is referred to as Overset POD (oPOD) since it shares some similarities with the overset mesh CFD approach in which multiple overlapping meshes are used within the CFD context.
In the oPOD approach, the ROM domain is separated into multiple subdomains. For each subdomain a POD+I model is computed. The individual predictions for an unknown parameter combination are then merged together for the complete ROM domain.
In this section the oPOD method, its terminology, as well as the steps involved in 63 creating and using an oPOD are described. First, the key problems oPOD addresses are 64 explained in detail and a brief overview over the method is given followed by an in-depth 65 explanation in the following sub sections and a summary. In order to better describe 66 and illustrate potential use cases and the oPOD method itself, a schematic CFD case of a 67 parameterized airfoil with and without a flap is utilized as displayed in figure 1.  A hard requirement for using POD is a congruent ordering of values across all snap-73 shots. For snapshots resulting from CFD simulation the ordering is determined by the 74 CFD index space which is unique for each generated CFD grid. Hence, utilizing snapshots 75 for different geometries is problematic. Utilizing mesh morphing (see [7]) keeps the CFD 76 element order and deforms the mesh with its geometries to a new configuration. Snapshots 77 resulting from these transformations share the same index space and can be used in a POD 78 context. This strategy works well for a wide range of deformation problems, but has its 79

Problem Definition
A hard requirement for using POD is a congruent ordering of values across all snapshots. For snapshots resulting from CFD simulacvtion, the ordering is determined by the CFD index space which is unique for each generated CFD grid. Hence, utilizing snapshots for different geometries is problematic. Utilizing mesh morphing (see [7]) keeps the CFD element order and deforms the mesh with its geometries to a new configuration. Snapshots resulting from these transformations share the same index space and can be used in a POD context. This strategy works well for a wide range of deformation problems but has its short comings. Large or nonlinear deformations can lead to mesh quality degradation and therefore negatively affect the CFD calculations. The same problem applies to deformations close to a wall or stationary object. This is schematically depicted in Figure 2 for different parameters p f . The cells between the stationary airfoil and the morphed flap are significantly stretched and a good resolution of the flow characteristics in this region is impaired.
Version July 6, 2022 submitted to Journal Not Specified 3 of 19 short comings. Large or non linear deformations can lead to mesh quality degradations 80 and therefor negatively effect the CFD calculations. The same problem applies to defor-81 mations close to a wall or stationary object. This is schematically depicted in figure 2 for  While the problems with mesh deformations are mostly gradual with increasing 86 magnitude of the deformation, it is not possible to model discontinuous parameters with 87 POD. Figure 1 shows the use case with and without a flap. For both cases a new grid needs 88 to be generated and the resulting snapshots do not share the same CFD index space. In general, oPOD uses a set of CFD snapshots, that are all based on a consistent set of 91 parameters. Dependent on the use case, the CFD domain is split into multiple oPOD sub 92 domains (step 1). For each sub domain a new mesh is generated, which can, depending on 93 the oPOD application, be significantly coarser than a CFD mesh. The CFD snapshots are 94 then mapped into the created oPOD sub domains (step 2) and a sub domain specific local 95 POD+I reduced order model is created (step 3). In order to perform an oPOD prediction 96 (step 4) the relevant POD+I models in the corresponding sub domains are queried and the 97 individual predictions are combined to again fill the complete oPOD domain. The criteria on how the CFD snapshots are split into the oPOD domains are use case 100 specific and one of the most critical steps in the oPOD process. During this step it is 101 important to keep in mind how an oPOD prediction (step 4) is performed and how the 102 multiple oPOD domains interact with one another.

103
The primary reason for creating an oPOD sub domain is to handle snapshots that are 104 based on different CFD index spaces. As discussed in section 2.1, the main reasons this 105 occurs are a CFD input deck where no mesh morphing was applied and each snapshot is 106 based on a unique CFD index space.

107
The sub domain creation is demonstrated below for the discontinuous parameter p d 108 for the CFD flow example above.

109
Discontinuous parameter 110 In general, a separate oPOD domain decouples certain regions from one another. It 111 is therefore possible to embed a discontinuous CFD snapshot parameter in a single or 112 multiple overlapping domains. While the problems with mesh deformations are mostly gradual with increasing magnitude of the deformation, it is not possible to model discontinuous parameters with POD. Figure 1 shows the use case with and without a flap. For both cases, a new grid needs to be generated and the resulting snapshots do not share the same CFD index space.

General oPOD Overview
In general, oPOD uses a set of CFD snapshots that are all based on a consistent set of parameters. Dependent on the use case, the CFD domain is split into multiple oPOD subdomains (step 1). For each subdomain, a new mesh is generated, which can, depending on the oPOD application, be significantly coarser than a CFD mesh. The CFD snapshots are then mapped onto the created oPOD subdomains (step 2) and a subdomain-specific local POD+I reduced-order model is created (step 3). In order to perform an oPOD prediction (step 4), the relevant POD+I models in the corresponding subdomains are queried and the individual predictions are combined to again fill the complete oPOD domain.

oPOD Domains (Step 1)
The criteria on how the CFD snapshots are split into the oPOD domains are use-casespecific and one of the most critical steps in the oPOD process. During this step, it is important to keep in mind how an oPOD prediction (step 4) is performed and how the multiple oPOD domains interact with one another.
The primary reason for creating an oPOD subdomain is to handle snapshots that are based on different CFD index spaces. As discussed in Section 2.1, the main reasons this occurs are a CFD input deck where no mesh morphing was applied and each snapshot is based on a unique CFD index space.
The subdomain creation is demonstrated below for the discontinuous parameter p d for the CFD flow example above.

Discontinuous Parameter
In general, a separate oPOD domain decouples certain regions from one another. It is therefore possible to embed a discontinuous CFD snapshot parameter in a single or multiple overlapping domains. This is visualized in Figure 3 in which the discontinuous CFD snapshot parameter p d is used to model a single airfoil or an airfoil with a trailing flap. The oPOD domains created based on the CFD domain is a background domain that only contains the airfoil for p d = w/o FLAP and for p d = w/ FLAP a separate foreground domain . In Figure 3, the background domain almost shares the same bounds as the CFD domain and only models the flow around as well as the wake of airfoil. The foreground domain models the region around the flap and the trailing wake of it.
). An interpolation region can be used on the boundary of the foreground 123 domain to smooth out the predictions of the background and foreground models (more on 124 this is discussed in step 4).

125
In this example, it is quite obvious that the setup of the domains is very dependent on 126 the use case. Since, depending on p d , the wake of the flow would have significantly different 127 features, the foreground domains need to extend as far into the background domain as 128 necessary in order to accurately model these features. Also, any other parameterization or 129 a larger value range for the parameter α would lead to a different setup of the foreground 130 domains. After setting up the oPOD domains resulting in new meshes for each domain in step 133 1, the CFD snapshots need to be mapped onto these meshes. If any geometric deformation 134 on the CFD snapshot was applied, it needs to be applied to the oPOD domain as well by 135 using, among others, radial basis functions [8,9] or the linear elasticity equations [10]. The 136 resulting elements of the deformed oPOD meshes may end up outside of the CFD domain. 137 By using the method of signed distances [11,12] these elements can be identified and are 138 left out of the mapping process.

139
With the prepared and aligned oPOD domains to the corresponding CFD snapshots, 140 the modeled variables can be mapped from the CFD snapshots to the oPOD domains. 141 Possible mesh interpolation methods are a consistent interpolation scheme such as the 142 advancing front algorithm [13], an R-tree spatial search algorithm [14] or a conservative 143 interpolation method by local Galerkin projection [15]. 144 Another important feature of the mapping process is that not all CFD snapshots 145 need to be mapped into all oPOD domains. Visiting the example in figure 3 again shows, 146 that each domain only receives snapshots for the corresponding value of p d . Hence, this 147 parameter is not part of the domain specific POD+I model since it is constant for all domain 148 specific snapshots. Also, the parameter p f is irrelevant for the background domain.  A prediction for a sought parameter p * = α * , p * f , p * d now would query the model of the background domain and depending on the value of p * d the foreground model (if p * d = w/ FLAP). An interpolation region can be used on the boundary of the foreground domain to smooth out the predictions of the background and foreground models (more on this is discussed in step 4).
In this example, it is quite obvious that the setup of the domains is very dependent on the use case. Since, depending on p d , the wake of the flow would have significantly different features, the foreground domains need to extend as far into the background domain as necessary in order to accurately model these features. Moreover, any other parameterization or a larger value range for the parameter α would lead to a different setup of the foreground domains.

Mapping Snapshots from CFD to oPOD Domains (Step 2)
After setting up the oPOD domains resulting in new meshes for each domain in step 1, the CFD snapshots need to be mapped onto these meshes. If any geometric deformation on the CFD snapshot was applied, it needs to be applied to the oPOD domain as well by using, among others, radial basis functions [8,9] or the linear elasticity equations [10]. The resulting elements of the deformed oPOD meshes may end up outside of the CFD domain. By using the method of signed distances [11,12] these elements can be identified and are left out of the mapping process.
With the prepared and aligned oPOD domains to the corresponding CFD snapshots, the modeled variables can be mapped from the CFD snapshots to the oPOD domains. Possible mesh interpolation methods are a consistent interpolation scheme such as the advancing front algorithm [13], an R-tree spatial search algorithm [14] or a conservative interpolation method by local Galerkin projection [15].
Another important feature of the mapping process is that not all CFD snapshots need to be mapped onto all oPOD domains. Visiting the example in Figure 3 again shows that each domain only receives snapshots for the corresponding value of p d . Hence, this parameter is not part of the domain-specific POD+I model since it is constant for all domain-specific snapshots. Moreover, the parameter p f is irrelevant for the background domain.

POD+I Creation in oPOD Domains (Step 3)
Once the oPOD domain-specific snapshots are prepared in steps 1 and 2, the POD+I model of each domain can be computed. In the case that all points/cells/elements of the oPOD subdomain are part of the CFD domain, a POD+I model based on the mapped snapshots can be computed. Otherwise, possibly due to transformations, it is possible that the domain for a number of snapshots lie outside of the corresponding CFD domain. In this case, not all variables of the CFD snapshot could be mapped on the corresponding oPOD snapshot domain. In order to also use these snapshots for the local POD+I models, a snapshot reconstruction algorithm needs to be employed. This reconstruction is discussed in detail in Section 3.3.
The Algorithm 1 summarizes steps 2 and 3. It shows how, based on the CFD snapshots W i CFD and a domain-specific reference grid G, a local oPOD domain POD+I model can be created. By looping over all available CFD snapshots, the relevant domain-specific snapshots are identified and the corresponding parameter vector is transformed to the local domain-specific parameters (lines 3 to 5). If a geometry transformation is applied, the reference grid G is transformed, resulting in a transformed grid G x i . Otherwise, the reference grid G is used (line 6). After mapping the CFD snapshot on the transformed mesh G x i , the resulting snapshots with the corresponding parameters are stored in separate sets (lines 7 to 9). After the loop finishes, the local POD+I model is computed based on the snapshots and parameters present in the stored sets (line 10).

Algorithm 1 oPOD domain-specific algorithm for local POD+I model creation.
Require: if W i CFD is relevant for local POD+I then 5:

oPOD Prediction (Step 4)
An oPOD prediction for a sought parameter p * is performed by querying the local POD+I models in only the relevant oPOD domains. Since the parameters used by the local POD+I models differ, it is necessary to filter or possibly transform the oPOD parameters to the parameters used by the local models. Depending on the construction of the oPOD domains, the predictions of various oPOD domains can be smoothed out. This is especially necessary if the domains are not large enough to model domain-specific phenomena and a more visualization-driven application is focused. Here, the signed distances from Section 2.4 can be used to smooth out the predictions of the multiple domains.
With the signed distances d i , computed at each element of the foreground domain with respect to the edge of the domain or a geometric object, and the interpolation distance r, the interpolation region ( Figure 3) is defined. By mapping the results of the background domain onto the foreground domain (similar to step 2), the predictions W * f and W * b of the fore and background model, respectively, are now present in each element of the foreground domain. The interpolated predictionŴ * is now computed witĥ Fluids 2022, 7, 242 6 of 18

Summary
An alternative visualization of Algorithm 1 for the creation of the POD+I models in the separate oPOD subdomains is displayed in Figure 4 as a flow chart.
Transformation and mapping of oPOD sub domain reference grid to a specific CFD snapshot.  Based on a set of already created oPOD subdomains, Figure 4a displays how, for a specific subdomain G j , the CFD snapshots and parameter are mapped, transformed and the resulting snapshots W x and parameters p x are used to compute the domain-specific (POD+I) j model. The details of the "map and transform CFD snapshots" block are shown in Figure 4b. In a loop over all CFD snapshots, it is checked if a specific CFD snapshot W i CFD with the parameter p i is relevant for the current oPOD subdomain G j . If this is the case, then the parameter and reference grid are transformed and the CFD solution is mapped onto the transformed grid G x . The resulting snapshot W x and parameter p x are stored and used to compute the domain-specific (POD+I) j model in Figure 4a.
All resulting subdomain models (POD+I) j are stored and used again for an oPOD prediction of a sought-after parameter combination as outlined in Section 2.6.

Theory
The aim to use or enhance predictions of reduced-order models is a major focus in the field and ". . . has seen tremendous development in the past two decades, especially in the broad domain of computational mechanics" [16]. An overview of the latest snapshot-based methods for parameterized partial differential equations are available in [1,16].
Related approaches to oPOD can be found in the literature. Shifted POD [17] ". . . extends the POD by introducing time-dependent shifts of the snapshot matrix" [17]. With this shift, it is possible to improve the POD predictions considerably by reducing the nonlinear effects across the snapshots. While it is not the key focus of this paper, an oPOD subdomain can be transformed and therefore "shifted" analogous the shifted POD, making use of a similar improvement as a consequence. The Domain-Decomposition POD (DD-POD) [7] utilizes POD in the context of CFD in order to model turbulent flows around parameterized geometries. Here, the ROM subdomain is coupled to CFD through an overlapping region using a Schwarz-type method [18], resulting in a nonlocal CFD boundary condition. This is comparable to oPOD in the sense that oPOD subdomains are also coupled to the surrounding domain.
The remaining section reviews the methods involved in computing a POD+I model with a special focus on snapshot reconstruction, which is needed in oPOD in the case the an oPOD subdomain transformation leads to elements outside the FOM domain.

Radial Basis Functions/Kriging
Radial basis function (RBF) or Kriging models are suitable to model complex functions. A comprehensive guide can be found in [5,19], and a short overview over the core method and features is given here.
With samples X = [x 1 , x 2 , . . . , x n ] T ∈ R n×m and responses y = [y 1 , y 2 , . . . , y n ] T ∈ R n , the model approximationf of an unknown sample x * can be approximate bŷ with n, m ∈ N, the basis function ψ and the weights w. The weights are determined by solving the linear system of equations Ψw = y, with Ψ ∈ R n×n being the Gram matrix which is defined as Popular choices for the RBF model basis function ψ are thin plate spine (TPS) with ψ(r) = r 2 ln r and cubic with ψ(r) = r 3 , with r being the euclidean distance of the two sample points r(x 1 , x 2 ) = x 1 − x 2 . A characteristic feature of the model basis functions is that their response decreases or increases monotonically with respect to the distance of the two samples [5].
The Kriging basis function is defined with ψ(x 1 , and introduces a varying exponent p, which typically is set to either 1 or 2, and a hyperparameter vector θ, which can be determined by a maximum likelihood estimation [19] (Section 2.4).

Proper Orthogonal Decomposition
In the field of fluid dynamics, three established reduced-order model-based proper orthogonal decomposition (POD) approaches exist by either using an interpolation coupled method (POD+I) [6,20,21], a CFD flux residual minimization scheme [22,23] or a Galerkin projection-based framework [24]. A comprehensive introduction can be found in [4,16,24]. A brief review of POD is given below.
For a set of independent parameters p = [p 1 , . . . , p d ] ∈ R d and CFD flow solution snapshots W i = W p i ∈ R n s , i = 1, . . . , m; the snapshot matrix is defined by By solving the m × m dimensional snapshot correlation eigenvalue problem the normalized eigenvectors V k ∈ R m can be found. By ordering the eigenvalues and corresponding eigenvectors based on λ 1 ≥ λ 2 ≥ · · · ≥ λ m , the m POD modes form an orthonormal basis, which is given by The coefficients a i k are determined with Equations (2) and 3 by Initial snapshots W i can be expressed with the newly found basis and an approximate snapshot W * at an untried parameter combination p * can be computed by a coefficient vector a * analogous to Equation (4) with The coefficient vector a * can be obtained by interpolating its components with a RBF or Krigin model at the parameter combination p * based on the data set p i , a i m i=1 . The POD coefficient matrix A is then defined with A = a 1 , . . . , a m T ∈ R m×m and the POD modes with U = U 1 , . . . , U m ∈ R n s ×m . A mode reduction of m r is applied to the POD by stripping the last m r rows and columns from the POD coefficient and modes matrix, respectively.

POD+I with Snapshot Reconstruction
Based on a snapshot reconstruction method introduced in [6,25], it is possible to iteratively guess and fill possibly missing elements of the snapshots. For a fixed number of iterations, a POD based on the current snapshot ensemble is computed, and the missing elements are repaired by using the resulting POD modes of the current iteration. In [6], a constant number of modes p is proposed. Setting p too low results in only low-level features being set in the missing elements, while setting p too high does not average out the to-be-detected features sufficiently.
In this paper, an improved method is introduced by determining p separately in each iteration of the snapshot reconstruction with the R 2 -Score of the corresponding interpolation model of the POD coefficients. The R 2 -Score is a measure for the accuracy of a surrogate model. An R 2 -Score of 1 denotes a perfect accuracy in the prediction, while an R 2 -Score of 0 denotes a model that constantly predicts the expected value of the responses regardless of the model input. Ref. [26] suggests to reduce the modes of a POD based on the R 2 -Score of the POD coefficients surrogate models by reducing the modes of the POD to the point that all coefficient surrogates have an R 2 -Score > 0. By using the R 2 -Score measure in the snapshot reconstruction procedure, it is possible to include the snapshot features depending on the current state of the snapshot ensemble during this iterative process.
The Algorithm 2 enlists in detail the snapshot reconstruction procedure with the initial setup in lines 1 and 2 where the initial snapshot matrix 1 Y is constructed by setting the missing values to the mean of the snapshot matrixW. In k max iterations, the reduced POD modes and coefficients are obtained with the R 2 -Score reduction and an R 2 -Score cut-off set to 0 (lines 4-5). Lines 7-9 apply the mask matrix to the modes and a best guess for the POD coefficient vector a is computed and applied to the corresponding column of the snapshot matrix of the next iteration (line 10). At last, the final POD modes and coefficients are computed after the loop finishes based on the final snapshot matrix k max Y resulting from the last iteration.

Require:
• n snapshots W i for parameter p i ∈ P with i ∈ [1, . . . , n] and P = [p 1 , . . . , Create POD coefficients and modes for k Y

5:
A, U ← R 2 -SCOREREDUCTION(A , U , P) Reduce coefficients and modes based on R 2 -Score reduction 6: Apply mask on j th mode 8: Apply mask on j th snapshot 9: a ← u T uu T −1 f Solve for coefficient vector 10: Create POD coefficients and modes for k max Y

Results
In this paper, oPOD is applied on two test cases. Both test cases have a discrete parameter p d that allows to optionally add or remove a feature to the overall process, splitting the domain into a separate background and a foreground POD+I. The first test case is a two-parameter one-dimensional process of an analytic function. It is well-suited to demonstrate the investigated oPOD method and analyze the effects each step has on the model prediction as well as the introduced errors. The second test case is the DrivAer body [27] with four parameters computed in a CFD environment.

Analytic Test Case
The analytic test case is split into two subsections. Both sections model the same process y in the domain x ∈ [0, 1]. The first section analyzes the performance of oPOD when all snapshot data is available to the model, while in the second section, only a subset of the information for the foreground model is available. Here, the reconstruction of the snapshot matrix, as depicted in Section 3.3, is applied. The to-be-modeled process consists of three distinct parameters. With a discrete parameter p d , and two continuous parameters p 1 and p 2 , the analytic test case is described by   As displayed in figure 5b the parameter p 1 introduces a scaling factor to the Gaussian 300 function, while the parameter p 2 shifts a local bell curve which can be seen in figure 301 5a. While the background model spans the complete domain, the foreground domain is 302 limited to x f = [p 2 − 0.2, p 2 + 0.2] ⊂ x. Also, the foreground model is dynamically shifted 303 according to parameter p 2 . 304 Figure 8 shows how a prediction with oPOD compares with a standard POD+I model. 305 The POD+I model was trained solely on the snapshots with p d = ON for the complete 306 domain x. As it can be seen in figure 8a the strong non linearities in the snapshots are 307 difficult to model for the POD+I model. It specifically struggles in the region of the peak 308 of the bell curve. oPOD is better capable of modeling this region, since the foreground 309 domain x f is shifted along with the curve and therefore differences in the snapshots are 310 smaller and easier to model (analogous to [17]). Figure 8b shows how the interpolation 311 in oPOD for the background POD b +I and the foreground POD f +I model is performed. 312 With an interpolation distance set to r = 0.075 from the domain boundary the output 313 of both models is linearly smoothed out. At the domain boundary for x = 0.15 oPOD 314 uses the values of the background model POD b +I while at x = 0.15 + r the values of the 315 foreground model POD f +I are used. In between the values are linearly interpolated with 316 respect to the distance from the domain boundary. This interpolation scheme results in a 317 good approximation to the to be model process.

318
The L 2 prediction errors across the parameter space for p d = ON are depicted in figure 319 6. The error e is given for parameter p * 1 and p * 2 by the model prediction oPOD i x i , p * 1 , p * 2 320 compared to the process true value Y i = y x i , p * 1 , p * 2 for varying x i ∈ x with e p * 1 , p * 2 = 321 oPOD − Y . Figure 6a shows the L 2 error for the standard POD+I model while figure 6b 322 shows the L 2 errors for the oPOD predictions. It is important to note the different color 323 ranges in the color bars of the two figures. As expected, the errors close to the snapshot 324  As displayed in Figure 5b, the parameter p 1 introduces a scaling factor to the Gaussian function, while the parameter p 2 shifts a local bell curve which can be seen in Figure 5a. While the background model spans the complete domain, the foreground domain is limited to x f = [p 2 − 0.2, p 2 + 0.2] ⊂ x. Moreover, the foreground model is dynamically shifted according to parameter p 2 . Figure 8 shows how a prediction with oPOD compares with a standard POD+I model. The POD+I model was trained solely on the snapshots with p d = ON for the complete domain x. As it can be seen in Figure 8a, the strong nonlinearities in the snapshots are difficult to model for the POD+I model. It specifically struggles in the region of the peak of the bell curve. oPOD is better capable of modeling this region, since the foreground domain x f is shifted along with the curve, and therefore differences in the snapshots are smaller and easier to model (analogous to [17]). Figure 8b shows how the interpolation in oPOD for the background POD b +I and the foreground POD f +I model is performed. With an interpolation distance set to r = 0.075 from the domain boundary, the output of both models is linearly smoothed out. At the domain boundary for x = 0.15, oPOD uses the values of the background model POD b +I, while at x = 0.15 + r the values of the foreground model POD f +I are used. In between, the values are linearly interpolated with respect to the distance from the domain boundary. This interpolation scheme results in a good approximation to the to-be-modeled process. Version July 6, 2022 submitted to Journal Not Specified  The L 2 error comparison of oPOD with and without snapshot reconstruction in figure 342 9 shows that the error increases due to the reconstruction of the snapshots. Compared to the 343 errors for oPOD without snapshot reconstruction the largest errors are no longer in between 344 the snapshot samples, but on the border of the parameter domain at p * = [100, 0.3]. This is 345 due to the fact that in those regions errors due to the reconstruction of the snapshots are 346 introduced. 347 Figure 10 displays how the error introduced in the snapshot reconstruction increases with increasing l co and the number of missing values in the snapshots data. The plotted cut-off error e co (l co ) is the maximum oPOD prediction error for a specific cut-off length l co e co (l co ) = max with e l co being the error evaluation for an oPOD with the corresponding l co . For this 348 academic use case the error for l co < 0.3, in which up to 50% of data is reconstructed, re-349 mains very small. For larger l co the error increases up to the point where the reconstruction 350 The L 2 prediction errors across the parameter space for p d = ON are depicted in Figure 6. The error e is given for parameter p * 1 and p * 2 by the model prediction oPOD i x i , p * 1 , p * 2 compared with the process true value Y i = y x i , p * 1 , p * 2 for varying Figure 6a shows the L 2 error for the standard POD+I model while Figure 6b shows the L 2 errors for the oPOD predictions. It is important to note the different color ranges in the color bars of the two figures. As expected, the errors close to the snapshot sample locations in the parameter space are small and increase with greater distance, but the errors for oPOD compared with POD+I are significantly smaller. sample locations in the parameter space are small and increase with greater distance, but 325 the errors for oPOD compared to POD+I are significantly smaller. In contrast to the previous section only the foreground model is now computed on 328 the bases of incomplete model data. Setting the cut-off length l co = 0.3, all data for x < 0.3 329 and x > 0.7 is stripped from the snapshot matrix and reconstructed using algorithm 2. 330 As depicted in figure 7 with the areas the loss in data is for some snapshots in the 331 foreground model significant and reduces the available data to half compared with the full 332 data.

Incomplete Model Data
In contrast to the previous section, only the foreground model is now computed on the bases of incomplete model data. Setting the cut-off length l co = 0.3, all data for x < 0.3 and x > 0.7 is stripped from the snapshot matrix and reconstructed using Algorithm 2. As depicted in Figure 7 with the areas, the loss in data is for some snapshots in the foreground model significant and reduces the available data to half compared with the full data. In contrast to the previous section only the foreground model is now computed on 328 the bases of incomplete model data. Setting the cut-off length l co = 0.3, all data for x < 0.3 329 and x > 0.7 is stripped from the snapshot matrix and reconstructed using algorithm 2. 330 As depicted in figure 7 with the areas the loss in data is for some snapshots in the 331 foreground model significant and reduces the available data to half compared with the full 332 data.  The error introduced with the snapshot reconstruction for the chosen l co is negligible 334 as the oPOD prediction with missing data almost matches the oPOD without missing data 335 (see figure 8).  The error introduced with the snapshot reconstruction for the chosen l co is negligible as the oPOD prediction with missing data almost matches the oPOD without missing data (see Figure 8). The L 2 error comparison of oPOD with and without snapshot reconstruction in Figure 9 shows that the error increases due to the reconstruction of the snapshots. Compared with the errors for oPOD without snapshot reconstruction, the largest errors are no longer in between the snapshot samples but on the border of the parameter domain at p * = [100, 0.3]. This is due to the fact that in those regions errors due to the reconstruction of the snapshots are introduced.  The L 2 error comparison of oPOD with and without snapshot reconstruction in figure 342 9 shows that the error increases due to the reconstruction of the snapshots. Compared to the 343 errors for oPOD without snapshot reconstruction the largest errors are no longer in between 344 the snapshot samples, but on the border of the parameter domain at p * = [100, 0.3]. This is 345 due to the fact that in those regions errors due to the reconstruction of the snapshots are 346 introduced. 347 Figure 10 displays how the error introduced in the snapshot reconstruction increases with increasing l co and the number of missing values in the snapshots data. The plotted cut-off error e co (l co ) is the maximum oPOD prediction error for a specific cut-off length l co e co (l co ) = max p * ,p * e l co (p * 1 , p * 2 ),  with e l co being the error evaluation for an oPOD with the corresponding l co . For this academic use case, the error for l co < 0.3, in which up to 50% of data is reconstructed, remains very small. For larger l co , the error increases up to the point where the reconstruction of snapshot data is infeasible, since the portion of data that needs to be reconstructed is too large (75% for l co = 0.4). of snapshot data is infeasible, since the portion of data that needs to be reconstructed is too 351 large (75% for l co = 0.4).
352 Figure 10. L 2 prediction errors in the parameter space P for oPOD foreground POD+I with increasing snapshot reconstruction size for p d = ON.

DrivAer Test Case
The DrivAer body [27] is used to demonstrate the oPOD capabilities for industrial CFD simulations. All snapshots are computed using OpenFOAM ® [28] with the Reynoldsaveraged Navier-Stokes (RANS) equations, a k-ω SST turbulence model and a free stream velocity of 30 m s . Each simulation consists of approximately 50 million cells with the geometries being deformed using mesh morphing based on radial basis functions. The simulation is based on five parameters. The four continuous parameters (see Figure 11    The fifth parameter p d ∈ [w/ SPOILER, w/o SPOILER] is discontinuous and controls 362 whether or not the spoiler is part of the CFD domain resulting in a hyper cube parameter 363  The fifth parameter p d ∈ [w/ SPOILER, w/o SPOILER] is discontinuous and controls whether or not the spoiler is part of the CFD domain resulting in a hypercube parameter domain P = p x × p z × p θ × p di f f × p d . The to-be-modeled variables are the windowed average velocity field U and pressure coefficient c p of the last 2000 iterations of the RANS simulations. Figure 12a shows the oPOD domain decomposition used for the test case. Within the CFD domain (not displayed in its correct size for visualization purposes), the background domain encapsulates the complete DrivAer body as well as its wake and is responsible for modeling the front of the body, and in case of p d = w/o SPOILER, also the wake. For p d = w/ SPOILER, the foreground domain models the wake and the flow around at the rear of the body. It contains all deformation-affected geometry parts and expands upstream as far as necessary in order to model the upstream effects of the spoiler. Moreover, it extends to the back of the ROM domain due to the dominant effect the spoiler has on the wake of the flow field.   The fifth parameter p d ∈ [w/ SPOILER, w/o SPOILER] is discontinuous and controls 362 whether or not the spoiler is part of the CFD domain resulting in a hyper cube parameter 363 domain P = p x × p z × p θ × p di f f × p d . The to be modeled variables are the windowed 364 average velocity field U and pressure coefficient c p of the last 2000 iterations of the RANS 365 simulations.
366 Figure 12a shows the oPOD domain decomposition used for the test case. Within the 367 CFD domain (not displayed in its correct size for visualization purposes) the background 368 domain encapsulates the complete DrivAer body as well as its wake and is responsible 369 for modeling the front of the body and in case of p d = w/o SPOILER also the wake. For 370 p d = w/ SPOILER the foreground domain models the wake and the flow around at the 371 rear of the body. It contains all deformation affected geometry parts and expands upstream 372 as far as necessary in order to model the upstream effects of the spoiler. Also it extends to 373 the back of the ROM domain due to the dominant effect the spoiler has on the wake of the 374 flow field.  For the results presented below, the foreground domain is modeled using 74 snapshots which are placed in the parameter space P f = p x × p z × p θ × p di f f × (p d = w/ SPOILER) using the quasi-random number Halton sequence [29]. The background domain is modeled using five equidistantly distributed samples in the parameter space P b = p di f f × (p d = w/o SPOILER). While the background model only depends on p di f f as its parameter, the foreground model is built using all continuous parameters.
In Figure 12b, the midslice of the applied signed distance field ratio d r (step 4, Section 2.6) is shown. Setting the interpolation region radius to r = 0.15, all domain boundaries except the rear boundary are interpolated according to Equation (1). The signed distance field is not affected by the DrivAer body since the foreground model should be used directly in the wake.

oPOD Prediction Evaluation
For the predictions of U and c p below a single prediction point p * = p * x , p * z , p * θ , p * di f f = (0.01, 0.06, −13.36, 0.43), which is not part of the oPOD sampling, is chosen. Depicted in Figures 13 and 14, the high-fidelity CFD solutions as well as the oPOD predictions of the magnitude of the velocity field and c p are shown. Compared with the high-fidelity CFD solutions, the differences of the oPOD predictions with and without spoiler are negligible. The wake of the body is modeled correctly for both states of the discontinuous parameter p d , and the transition of the background and foreground models at the edge of the foreground model is smooth without any discontinuities. Especially in the spoiler region, the model predictions are in very good agreement with the CFD solutions. For the results presented below, the foreground domain is modeled using 74 snapshots 376 which are placed in the parameter space using the quasi-random number Halton sequence [28]. The background domain is mod-378 eled using five equidistantly distributed samples in the parameter space P b = p di f f × 379 (p d = w/o SPOILER). While the background model only depends on p di f f as its parame-380 ter the foreground model is built using all continuous parameters.

381
In figure 12b, the mid slice of the applied signed distance field ratio d r (step 4, section 382 2.6) is shown. Setting the interpolation region radius to r = 0.15 all domain boundaries 383 except the rear boundary are interpolated according to equation 1. The signed distance field 384 is not affected by the DrivAer body since the foreground model should be used directly in 385 the wake.  For the results presented below, the foreground domain is modeled using 74 snapshots 376 which are placed in the parameter space using the quasi-random number Halton sequence [28]. The background domain is mod-378 eled using five equidistantly distributed samples in the parameter space P b = p di f f × 379 (p d = w/o SPOILER). While the background model only depends on p di f f as its parame-380 ter the foreground model is built using all continuous parameters.

381
In figure 12b, the mid slice of the applied signed distance field ratio d r (step 4, section 382 2.6) is shown. Setting the interpolation region radius to r = 0.15 all domain boundaries 383 except the rear boundary are interpolated according to equation 1. The signed distance field 384 is not affected by the DrivAer body since the foreground model should be used directly in 385 the wake.  For the results presented below, the foreground domain is modeled using 74 snapshots 376 which are placed in the parameter space P f = p x × p z × p θ × p di f f × (p d = w/ SPOILER) 377 using the quasi-random number Halton sequence [28]. The background domain is mod-378 eled using five equidistantly distributed samples in the parameter space P b = p di f f × 379 (p d = w/o SPOILER). While the background model only depends on p di f f as its parame-380 ter the foreground model is built using all continuous parameters.

381
In figure 12b, the mid slice of the applied signed distance field ratio d r (step 4, section 382 2.6) is shown. Setting the interpolation region radius to r = 0.15 all domain boundaries 383 except the rear boundary are interpolated according to equation 1. The signed distance field 384 is not affected by the DrivAer body since the foreground model should be used directly in 385 the wake.  The plots in section 4.2.1 show oPOD predictions with an enabled interpolation region 398 between the background and foreground oPOD sub domains. As mentioned above, the 399 prediction of a significantly large foreground sub domain model, which encapsulates all 400 variances introduced by the to be modeled parameters, does not need to be interpolated 401

oPOD Analysis
The plots in Section 4.2.1 show oPOD predictions with an enabled interpolation region between the background and foreground oPOD subdomains. As mentioned above, the prediction of a significantly large foreground subdomain model, which encapsulates all variances introduced by the to-be-modeled parameters, does not need to be interpolated with the background model. However, since it is beneficial to limit the subdomain to the region of interest due to memory and performance reasons as well as possible complex geometry transformations, an interpolation region can be required.
The effect of the interpolation region is displayed in Figure 15. Marked by the white circles are the areas in which the interpolation effect of the predicted values at the border of the foreground and background domain are most pronounced.  with the background model. But since it is beneficial to limit the sub domain to the region 40 of interest due to memory and performance reasons as well as possible complex geometry 40 transformations, an interpolation region can be required.

40
The effect of the interpolation region is displayed in figure 15. Marked by the white 40 circles are the areas, in which the interpolation effect of the predicted values at the border 40 of the foreground and background domain are most pronounced.

40
The final analysis of oPOD is a comparison against an alternative approach on model-40 ing the discontinuous spoiler parameter p d . In this approach, a domain with a fixed mesh 40 on which the solutions of the snapshots are mapped, is used. As a consequence, the cells in 41 the spoiler area occasionally are on the upper or lower side of the spoiler, depending on the 41 snapshots parameter, or do feature a flow without spoiler at all. This can be problematic for 41 the POD+I model, since high non linearities of the modeled variables in the corresponding 41 cells are introduced due to the fact that these cells are subject to high fluctuations compared 41 to the parameters. The final analysis of oPOD is a comparison against an alternative approach on modeling the discontinuous spoiler parameter p d . In this approach, a domain with a fixed mesh on which the solutions of the snapshots are mapped is used. As a consequence, the cells in the spoiler area occasionally are on the upper or lower side of the spoiler, depending on the snapshots parameter, or do feature a flow without spoiler at all. This can be problematic for the POD+I model since high nonlinearities of the modeled variables in the corresponding cells are introduced due to the fact that these cells are subject to high fluctuations compared with the parameters. Figure 16 compares the alternative method for modeling discontinuous parameters with oPOD. Displayed in Figure 16a is the oPOD prediction for a solution with spoiler. It can be seen that the stagnation point in the front as well as the separation point at the trailing edge of the spoiler are accurately modeled. Moreover, the recirculation area beneath the spoiler is closely attached to the spoiler geometry. For the alternative approach, the model struggles with predicting these specific features, as the nonlinearities in the cells are too severe. trailing edge of the spoiler are accurately modeled. Also the recirculation area beneath the 419 spoiler is closely attached to the spoiler geometry. For the alternative approach the model 420 struggles with predicting these specific features, as the non linearities in the cells are too 421 severe.

423
A novel method for handling discontinuous parameters for POD+I ROMs is presented. 424 The steps involved are discussed and the method is demonstrated for two use cases. The 425 industrial DrivAer use case shows significant improvements compared to an alternative 426 way of modeling discontinuous parameters. It is possible to perform complex geometry 427 trailing edge of the spoiler are accurately modeled. Also the recirculation area beneath t spoiler is closely attached to the spoiler geometry. For the alternative approach the mod struggles with predicting these specific features, as the non linearities in the cells are t severe.

Conclusion
A novel method for handling discontinuous parameters for POD+I ROMs is present

Conclusions
A novel method for handling discontinuous parameters for POD+I ROMs is presented. The steps involved are discussed and the method is demonstrated for two use cases. The industrial DrivAer use case shows significant improvements compared with an alternative