Abstract
In the present work, supersonic flows over an axisymmetric base and a 24-deg compression ramp are investigated using the generalized - (GEKO) model implemented in the commercial software, ANSYS FLUENT. GEKO is a two-equation model based on the - formulation, and some specified model coefficients can be tuned depending on the flow characteristics. Uncertainty quantification (UQ) analysis is incorporated to quantify the uncertainty of the model coefficients and to calibrate the coefficients. The Latin hypercube sampling (LHS) method is used for sampling independent input parameters as a uniform distribution. A metamodel is constructed based on general polynomial chaos expansion (gPCE) using ordinary least squares (OLS). The influential coefficient closure is obtained by using Sobol indices. The affine invariant ensemble algorithm (AIES) is selected to characterize the posterior distribution via Markov chain Monte Carlo sampling. Calibrated model coefficients are extracted from posterior distributions obtained through Bayesian inference, which is based on the point-collocation nonintrusive polynomial chaos (NIPC) method. Results obtained through calibrated model coefficients by Bayesian inference show superior prediction with available experimental measurements than those from original model ones.
1. Introduction
Recently, uncertainty quantification (UQ) has been applied to computational fluid dynamics (CFD) by improvement of computing performance and reduction of computational cost. UQ is used to estimate not only the mean value and standard deviation but also the probability distributions of QoIs (Quantities of Interest) by considering the probability distributions of input variables such as boundary conditions or applied model coefficients. The Monte Carlo (MC) technique has been the general method used to compute output variables by sampling the random inputs with a specific distribution. Its basic concept is relatively simple but too expensive to achieve reasonable convergence, which means that application to CFD is difficult. The generalized polynomial chaos (gPC) method that was recently proposed by Xiu [1] requires much fewer samples than the MC technique so that studies using UQ can be performed efficiently.
Uncertainties in CFD can be divided into aleatory uncertainty and epistemic uncertainty. Aleatory uncertainty relates to the random components contained in a real system; from the perspective of CFD, examples would include uncertainties in geometry and boundary conditions such as velocity and pressure. On the other hand, epistemic uncertainty occurs due to a lack of knowledge regarding physical phenomena and modeling hypotheses. Thus, uncertainties caused by the process of model formulation and estimation of model coefficients, such as those used in modeling turbulence or multiphase flow, are typical examples of epistemic uncertainty.
The methodology of the UQ can be classified into two approaches, the forward problem and the inverse problem, depending on how the random inputs are treated [2]. In the case of a forward problem, the distribution of random inputs is given with a fixed mean and interval, and then one needs to compute the stochastic outputs with mean and deviation values based on the selected UQ methodology. However, in the inverse problem, the Bayesian inference is adopted to estimate the distribution of random inputs using observation data such as experimental results or synthetic data. Through the process of model calibration, the distribution of the random inputs is updated, and then QoIs are computed statistically; this process is called model prediction.
Schaefer et al. [3] applied point-collocation nonintrusive polynomial chaos (NIPC) to various turbulence models to study the effect of the random inputs, including aleatory and epistemic uncertainties, on the aerodynamic performance coefficients. They assumed model coefficients of the Spalart-Allmaras (SA) model and - family models as a uniform distribution and obtained the statistical distributions of pressure and aerodynamic coefficients over an RAE 2822 airfoil and an axisymmetric bump. Additionally, influential coefficients among the model coefficients of the corresponding models are investigated through Sobol indices.
For high-speed flows, Huan et al. [4] conducted a global sensitivity analysis to identify effective input parameters, which are inflow, fuel inflow, and wall boundary conditions, and turbulence model parameters in the CFD of a scramjet system by using UQ methods. Additionally, the estimated model errors for models of different fidelity were investigated. Burt and Josyula [5] considered aleatory and epistemic uncertainties in sensitivity analysis/UQ calculations. Global sensitivity analysis and UQ are integrated with a direct simulation, which is a Monte Carlo gas flow simulation code for a hypersonic double-cone flow.
Recently, the generalized - (GEKO) model, which is a two-equation model based on the - model formulation, was proposed by Menter et al. [6]. The GEKO model itself includes model coefficients that can be controlled by users depending on the flow type, which means that the coefficients are able to be tuned in a variety of flows. Figat and Piatkowska applied the GEKO model to investigate the interaction between the propeller and the fuselage of an aircraft [7]. Brezgin and Aronson adjusted the value of for the sensitivity regarding adverse pressure gradients and spreading rates, which were overpredicted by conventional models (standard -, SST -) [8]. Additionally, Li et al. [9] compared the GEKO model with other turbulence models to investigate the influence of wall conduction on the heat transfer of supercritical n-decane in active regenerative cooling channels.
In the present work, supersonic flows over an axisymmetric base and 24-deg compression ramp are investigated. In the simulation of a 24-deg compression ramp, predicting shock wave/boundary layer interaction (SWBLI) accurately is an important factor in the design of high-speed flight vehicles. On the exterior flow over the aircraft, SWBLI can cause loss of control, some peaks on surface thermal loading. In internal flows, it can enhance distortion and pressure losses and can even cause catastrophic events leading to engine unstart. Rizzetta [10] investigated supersonic flow over a 24-deg compression ramp by using explicit algebraic Reynolds stress models and - models, but the model that gave the best overall results was the Speziale-Abid [11] realizable - model. Gerolymos et al. [12] also used Reynolds stress models and concluded that the WNF-LSS RSM model gave the best prediction of skin friction at reattachment and in the relaxation region. In terms of experimental studies, Settles et al. [13] tested a 24-deg compression ramp and mentioned that two-dimensionality can be considered because the three-dimensional perturbations of the flow are minor. However, aerodynamic fences are required to achieve two-dimensionality. Horstman et al., Dolling and Murphy, and Settles and Dodson also experimentally studied compression ramps [14–16].
In the case of base flow, the accurate prediction of drag has been a challenge problem in CFD. A supersonic body experiences major drag from skin friction drag, wave drag, and pressure drag. If the drag is able to be controlled, stability and control of vehicles can be enhanced. As for an experiment with supersonic base flow, Herrin and Dutton [17] experimented on a supersonic axisymmetric base flow at . Many researchers have tried to analyze the supersonic base flow by using CFD. Forsythe et al. [18] used Reynolds-averaged Navier-Stokes (RANS) and detached eddy simulation (DES) with a compressibility correction to analyze supersonic base flow. Both models using the compressibility correction showed improved pressure distribution compared to the models that did not use the correction. In general, it is known that large eddy simulation and hybrid RANS/LES simulation are able to predict the pressure distribution and recirculation zone better than traditional RANS, but when the computational time and cost are considered, RANS is still the most widely used method in engineering.
In the present work, UQ with Bayesian inference is applied to the GEKO turbulence model to estimate the optimized model coefficients and QoIs in supersonic flows over an axisymmetric base and 24-deg compression ramp. The present work demonstrates Bayesian uncertainty quantification of GEKO turbulence model coefficients for the first time, to the best of our knowledge. The point-collocation NIPC technique [1, 19], which is one of the typical nonintrusive methods in UQ, is adopted for the UQ framework. Three model coefficients are considered the random input variables. As for QoIs, the pressure coefficient along the base radius () and reattachment point (RP) are set for the axisymmetric base case, whereas the separation point (SP) and RP are set for the 24-deg compression ramp case. The process of UQ, including input variable sampling, surrogate model construction, global sensitivity analysis, and Bayesian inference, is performed using the UQLab framework [20].
2. Numerical Method
2.1. Numerical Method
The numerical simulation is performed using the commercial software ANSYS FLUENT v.19.1 [21] to solve compressible Navier-Stokes equations and the energy equation. The advection upstream splitting method (AUSM+) scheme is selected because it provides exact resolution of contact and shock discontinuities, preserves positivity of scalar quantities, and gives less dissipation than the Roe scheme does. As for spatial discretization, the third-order monotone upstream-centered (MUSCL) scheme is incorporated to improve spatial accuracy by reducing numerical diffusion. The Courant-Friedrichs-Lewy (CFL) number is set to gradually increase from 0.25 to 5 depending on the convergence circumstance. To compensate for the defect predicting low pressure on the base, a compressibility correction lowering the turbulent eddy viscosity production is applied to the computational models [18, 22–24]. On the other hand, for the 24-deg compression ramp, a compressibility correction is not used due to degraded prediction in near-wall regions [10, 12, 25–28]. Results obtained with the compressibility correction are denoted by “-CC” in the manuscript. As mentioned in the previous section, the GEKO model is adopted in the present work. Other models including baseline (BSL) - model, shear stress transport (SST) model, and - realizable model in FLUENT [21] are also used for comparison.
2.2. Turbulence Model
The characteristic of the GEKO model is that the model coefficients can be controlled to tune the model to various flow scenarios. In this present work, three model coefficients are considered for compressible flows.
The GEKO model formulation is as follows:where
The coefficients of the GEKO model are implemented through the functions , which can be controlled to achieve different goals in a variety of parts of the computational domain. Although the details of formulations , including values of , , and , are not released, the range of the coefficient values and characteristics of each coefficient are given below:
(1). Main parameter for adjusting separation prediction for boundary layers. Affects all flows. Increasing reduces eddy viscosity, leading to greater sensitivity to adverse pressure gradients for boundary layers and to lower spreading rates for free shear flows.
(2). Affects mostly the inner part of wall boundary layers (no impact on free shear flows). Increasing leads to higher wall shear stress and wall heat transfer rates in nonequilibrium flows.
(3). Affects only free shear flows (boundary layer shielded due to function , which is discussed later). Increasing increases spreading rates of free shear flows. For each value of , an optimal value of exists, which maintains optimal free shear flow. This value is given by the correlation , which is the default.
(4). Is active in a submodel of (no impact when is equal to 0). Affects mostly jet flows. Increasing while is active decreases the spreading rate for jets. Allows adjustments to the spreading rate of jet flows while maintaining the spreading rate of the mixing layer.
(5). Nonlinear stress-strain term to account for secondary flows in corners (e.g., wing-body junctions).
(6). An existing model for curvature correction, which can be combined with the GEKO model.
All coefficients can be accessed globally or locally by using user-defined functions (UDFs), which allows a global or zonal model optimization.
The coefficients and are designed for free shear flows, whereas and have an effect on boundary layers. To avoid any influence of and on boundary layers, a blending function is incorporated, which deactivates and in the boundary layer. The blending function is given by
This function activates following the shear flow parameters:
There are two important aspects to . Firstly, inside boundary layers, the function is equal to 1, which means that the function becomes equal to 0. Secondly, the parameter is a subparameter of . As mentioned above, it only affects the simulation in cases where is not equal to 0. The Min and Max values of are suggestions from [8]. There might be situations where values lower than Min (0.5) or higher than Max (1.0) would be appropriate for specific flows. To avoid negative effects on free mixing layers due to changes in , however, the use of (Equation (3)) is recommended. has an effect on the secondary flows at the corner in the 3D domain, and affects the cylinder geometry as a cyclone for the curvature correction. Hence, in the present work, model coefficients (, , and ) are varied within the ranges given in Table 1.
3. Uncertainty Quantification
3.1. Point-Collocation Nonintrusive Polynomial Chaos (NIPC)
In this current work, point-collocation NIPC, which was proposed by Hosder et al. [29], is applied. The point-collocation NIPC technique requires the minimum number of random input variables calculated by Equation (6) consisting of the polynomial order (), the number of random input variables (), and the oversampling ratio (). The minimum amount of random data is extracted by using the Latin hypercube sampling (LHS) method. Extracted random input variables are used for Equation (7).
Here, are the polynomial chaos coefficients and is an element of an orthogonal family, and represents the output corresponding to a given random input. Note that and in the above equation correspond to given values, and is obtained by applying the Gram-Schmidt technique to the following matrix:
The gPC technique utilizes the properties of a base function orthogonal. The basis functions of the gPC used according to the distribution of random input variables are shown in Table 2. In the present work, the Legendre basis function is used due to the uniform distribution.
The least square minimization method is adopted, which has the advantage that an arbitrary number of points can be used to calculate the coefficients, as long as they are a representative sample of the random input vector. Least square minimization is a method that minimizes the truncation error. Truncated PCE and a residual are expressed as follows:where is the order of PCE, is the truncation error, is a vector containing the coefficients, and is the matrix that assembles the values of all the orthonormal polynomials in . The least square minimization can be expressed as follows:
Equation (10) is solved by ordinary least squares (OLS). The ordinary least square solution iswhere
is the model response corresponding to a random input vector.
3.2. Leave-One-Out Cross-Validation Error
The method of the leave-one-out cross-validation error is adopted to evaluate the error of the constructed PCE. The is designed to overcome overfitting limitations by using cross-validation. When the least square minimization is used for calculating the coefficients, the formulation to calculate iswhere is the th component of the vector, which is given below:
In Equation (14), is the experimental matrix. is the model response and is the PCE result. Here, is the mean of the model response.
3.3. Bayesian Inference
Bayesian inference is a process of updating the probability of based on the frequency at which datum is observed after assuming a prior distribution of the random variable . Unlike the frequentist methodology, which estimates conditional probabilities based on a large number of samples, Bayesian inference estimates the posterior probability of a random variable with a small number of samples. The basic formula for Bayesian inference is shown below. Here, is the prior distribution of the input parameters. This could be set up heuristically or could be set up through a predetermined tendency. Similarly, is the likelihood. The likelihood measures the suitable statistical model with respect to given data that consider observations and prior distribution. The factor is used to normalize the right-hand side. As a result, the posterior distribution is finally calculated.
In the case of the independent observations given, the likelihood function is modeled as follows:
If the independent observations are used, the likelihood function is expressed as a form of multiplication (Equation (16)). In Equation (18), is an observation, is a computational forward model, and is a set of input parameters. The discrepancy term represents the effects of the measurement error, which is obtained from experimental references [4, 25–29].
4. Computational Setup
4.1. Geometry
The geometry of the base was selected from Forsythe et al. [18]. The cylinder length was set to , where is the base radius (31.75 mm). This cylinder length is determined to match the experimental momentum thickness [18]. The outlet is located downstream from the base, while the far-field boundary is at from the axis of symmetry.
In the case of the 24-deg compression ramp, the geometry was adopted from Gerolymos et al. [12]. A no-slip wall is placed on the top and bottom. The upper wall is located 0.2 m away from the bottom wall of the downstream section. The inlet is 0.2 m from the beginning of the 24-deg ramp, which is 0.15 m long. Various frames were used in the compression ramp experimental setups [13, 30] for the profile measurements. The details of frames and specific locations where experimental data were measured are thoroughly explained in [12, 13, 30].
4.2. Grid
A structured grid is generated by using Pointwise [31]. To resolve the turbulent boundary layer near the wall and satisfy the wall unit within 1.0, which is averaged in the computational domain, the first cell height is set to (m). For the grid sensitivity test, a total of five types of grid were generated, which is discussed in detail in the result section. The total number of cells was 37,370 ( upstream and downstream from the base), while 180,000 cells () were generated for the compression ramp case. The grids are shown in Figures 1 and 2.
(a)
(b)
(a)
(b)
4.3. Boundary Conditions
For axisymmetric base flow, the experimental conditions of Herrin and Dutton [17] were used in the present computation. supersonic flow comes from the inlet with (a unit Reynolds number of ). Adiabatic wall conditions were adopted at the wall.
In the case of the compression ramp, profiles obtained from the 2D channel simulation were specified for all dependent variables (e.g., , , and ) at the upstream inlet boundary of the computational domain. Supersonic inflow condition was computed at the inlet boundary (a unit Reynolds number of ). The wall temperature was fixed to , in accordance with measurements [13, 14]. Extrapolation was applied at the downstream outlet boundary, and test conditions are shown in Table 3.
5. Results
5.1. Grid Test
The grid sensitivity was investigated using the grid convergence index (GCI) using five types of grid (tiny, coarse, medium, fine, and extra fine), which are characterized in Tables 4 and 5. There are several GCI methods, but simplified least square GCI (SLS-GCI), which was developed by Eça and Hoekstra [32], is used. The results obtained by SLS-GCI can be explained using three components with respect to QoIs:
(1). A measure of how far the computed value is from the value of the asymptotic numerical value.
(2). The local error of estimation meaning the difference between the numerical and estimated results.
(3) Extrapolation Values. An estimate of the value of each QoI at zero grid spacing.
To evaluate the grid test, RP and (Equation (19)) are used for the axisymmetric base case. On the 24-deg compression ramp, RP and SP are considered. These are QoIs, which are used not only in the grid test but also for UQ. Table 6 shows and values. The values in the brackets show , and both are very close to zero, which means the computation is within the asymptotic range and the grids are reasonably generated. Table 7 and Figure 3 show extrapolation values with respect to the QoIs. The values inside the brackets in Table 7 mean the error (%) between each extrapolation value and the result of the medium mesh. Combination has the biggest spacing ratio, but the difference from the result of the medium mesh is less than 5% (Max 4.5%). Hence, the medium mesh is selected for the computational domain. Figure 4 shows the results of all grids.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
5.2. Deterministic Solver Validation
5.2.1. Axisymmetric Base
In the case of the base flow, the velocity profile was obtained 1 mm upstream of the base corner. Figure 5 shows the velocity profiles, which are compared with the profile obtained from the experimental data. Whether compressibility correction is used or not, the profiles of all models appear in a similar shape. Overall, all models predicted the boundary layer profile reasonably well under the current mesh.
Pressure distributions along the base surface are shown in Figure 6. The base pressure distributions obtained without the compressibility correction are quite a bit lower than the experimental data [17]. On the other hand, the models including the compressibility correction show higher base pressure distributions but have greater variations along the base surface. Qualitative flat pressure distribution along the base surface is not predicted at all by RANS models [18, 22–24, 33, 34]. The reverse flow near the wake axis stagnates at the center of the base surface where the relatively high pressure appears (see Figures 6 and 7).
One of the QoIs for axisymmetric base flow is the shear layer RP. Figure 7 shows the axial velocity along the centerline of the axis. It is clearly shown that results from the models without compressibility correction match the experimental data well. If a smaller recirculation region is predicted, as in the models without compressibility correction, it causes the flow to turn sharply behind the base, leading to a more enhanced expansion wave and a reduction in pressure. A small separation region, therefore, causes larger pressure drag than a large separated region. Models with compressibility correction overpredicted not only the RP but also the peak reverse velocity. When the reverse flow is overpredicted, the flow is more accelerated outward along the base radius, which can cause a pressure reduction. The failure of RANS models to predict a flat pressure distribution along the base surface due to the large turbulent eddy viscosity is a well-known issue [34].
The contours of the Mach number and turbulent viscosity ratio are shown in Figure 8. The upper and lower parts of the figure correspond to GEKO-CC and GEKO, respectively. It is clear that the GEKO-CC model predicts RP further from the base than the GEKO model. As mentioned in the introduction, compressibility correction reduces turbulent eddy viscosity, which can lead to higher pressure levels, but it also increases the reattachment length. Thus, the size and length of the recirculation region behind the base are influenced by overpredicting or underpredicting the turbulent eddy viscosity. This factor can also be controlled by the coefficient of GEKO model, and it is indirectly indicated through turbulent viscosity ratio contours (Figure 8) how big an influence could have. The relationship between overprediction and/or underprediction of the RP and compressibility correction has also been investigated by many researchers [18, 22–24, 33, 34].
(a)
(b)
5.2.2. 24-deg Compression Ramp
For the 24-deg compression ramp case, the profile including all dependent variables is computed at the inlet. The boundary layer thickness (), displacement thickness (), and momentum thickness () of the results obtained from the models (BSL, SST, and realizable models) are shown in Table 8. The profiles are extracted at a position where the profiles were measured in the experiment [30]. The SST model predicts a relatively lower value of than the other models do.
Figure 9 shows nondimensional pressure distributions along the bottom surface of the computational domain. To compare computational results, experimental data obtained from [13, 15, 16] were used. (The experimental data was found in [16], but the exact source was not given [12].) The result of the GEKO model is quite similar to the result obtained from the SST model, as indicated [6]. Both models overpredicted the length of the reversed flow region, which means that a separation shock occurs too early. The models underpredicted even lower pressure compared to the experimental data after the corner. The BSL and realizable models performed the best in the pressure distribution. The initial pressure rise of the results obtained from both models matches the measured rise. This means that the separation shock predicted from numerical flow fields occurred in the same location as observed in the experimental flow fields.
Figure 10 shows the skin friction on the bottom surface of the compression ramp domain. As mentioned before, it is seen that the GEKO and SST models both overpredicted the separation length. Further, the GEKO and SST models predicted quite low skin friction in the relaxation region and more delayed RP than measured. In contrast with the GEKO and SST models, again, the BSL and realizable models show better agreement when it comes to the separation shock, SP, and skin friction distribution after the RP. The BSL and realizable models accurately predicted the point that skin friction started decreasing, which means the separation shock location is predicted quite well. Looking at it more closely, the BSL model predicted the separation shock location slightly more accurately than the realizable model did, but the realizable model predicted the SP a little bit closer to the experimental data. As for the RP, all models failed to predict the RP perfectly. Having said that, however, the BSL and realizable models predicted the RP much more closely than the GEKO and SST models did, which means that the BSL and realizable models predicted much shorter separation lengths.
The Mach number contours and static pressure are shown in Figures 11 and 12, respectively. In the Mach number contours, all models predict separation shock and reattachment shock, but the GEKO and SST models predict much larger separation bubbles than the BSL and realizable models do, as seen above in Figures 9 and 10. Additionally, the GEKO and SST models underpredict pressure compared to the realizable and BSL models in the relaxation region.
5.3. Uncertainty Quantification
5.3.1. Sensitivity Analysis
The random input variables in the present work are , , and . The LHS method was used for sampling these input variables, assuming a uniform distribution. The selection of the interval of random inputs is critical in the UQ process and has many effects on the distribution of QoIs and the posterior one of RVs. In the present work, the interval of each model coefficient was set with the interval of coefficients of the GEKO model proposed by Menter et al. [6], which is shown in Table 1. The total number of samples can be set if the order of polynomial chaos and oversampling rate are set (Equation (6)). Once the total samples are computed in a deterministic solver and calculated, a metamodel (surrogate model) can be constructed by using the results obtained from the deterministic solver. One of the several standards validating the robustness of the constructed metamodel is .
Even though the same number of samples was used, of the axisymmetric base case is much higher than that of the 24-deg compression ramp case, which means that the axisymmetric base case is less robust than the 24-deg compression ramp case. To investigate what components affect , we carried out a sensitivity analysis. As mentioned above, the number of samples depends on the order of polynomial chaos () and oversampling rate . The left side of Table 9 presents the number of samples required for each combination of and . The right side of Table 9 shows the corresponding to each QoI. decreases only when the increases, which means increasing the could be a better choice if the data is limited or computational cost is expensive. However, an appropriate order for the polynomial chaos () is also crucial when it comes to Bayesian inference. Schaefer et al. [3] set and mentioned that most of the turbulence uncertainty can be quantified by using . However, in the present work considering Bayesian inference, the oversampling rate was not sufficient to construct the surrogate model. In this study, is set to 3, and is selected considering QoI distribution. Figure 13 shows the distributions of and RP corresponding to the polynomial order, . All distributions look quite similar to each other. In order to decide the specific order for the surrogate model and Bayesian inference, the 5th order is selected for the base flow. In the case of ramp flow, the 3rd order is chosen due to the limited data.
5.3.2. Global Sensitivity Analysis
A total of 168 samples for base flow and 60 samples for ramp flow are used to construct the metamodel and to analyze which model closure coefficient is dominant in each flow. Table 10 shows Sobol indices of closure coefficients for QoIs of each test case. The largest contributors to uncertainty in each flow case are typed in italic (closure coefficients with Sobol indices of less than were not considered significant for reduced dimensionality analyses [3]). It is apparent that is the most dominant closure coefficient. As Menter et al.[6] noted, is the main parameter for adjusting separation prediction and has effects on the eddy viscosity for boundary layers and the spreading rate for free shear flows. The next dominant closure coefficient is , which affects wall shear stress. Sobol indices of in the 24-deg compression ramp case are higher than those in the axisymmetric base case. This could be presumed to be affected by the fact that the SWBLI phenomena occurred on the surface of the compression ramp. has not made as much contribution to either case as the other closure coefficients have.
5.3.3. Bayesian Inference
Based on the metamodel constructed using OLS, a likelihood function was calculated using the observation data. In the present work, experimental data of QoIs such as the pressure coefficient, RP in the base flow, and SP and RP in the ramp flow were applied for the observation data. Experimental data that is assumed to have a Gaussian distribution including an error of the experimental measurement is considered for the observation data. For MCMC sampling, the affine invariant ensemble algorithm (AIES) is used to overcome the issue that MCMC algorithms result in poor convergence when there is a strong correlation shown in the posterior distribution between parameters. It performed 1500 steps and 300 parallel chains. The first half of the sample points generated by all chains are removed as burn-in, and postprocessing is performed after burn-in. It took about 30 minutes for MCMC sampling in each case. Figures 14 and 15 show prior and posterior distributions of QoIs for the axisymmetric base, while Figures 16 and 17 show the distributions for the 24-deg compression ramp. The common finding from these four figures is that the prior distribution, which was assumed as a uniform distribution, produces a posterior distribution that does not show a uniform distribution any longer. This means that the results composed with the input variables assumed as a uniform distribution are not uniform distributions, and also, each input variable presents a specific distribution by the observation data, which were used to calculate the likelihood function. However, each model coefficient shows a different posterior distribution and different correlated values (red lines) depending on the selection of QoI, even though in the same flow, axisymmetric base flow, or compression ramp one.
The correlated coefficients were extracted from the posterior distributions in Table 11, which shows the correlated coefficients according to QoIs in both the base and ramp flows. The correlated coefficients of are predicted as 2.290 for and 0.723 for RP in the same axisymmetric base flow. When the minimum and maximum values of 0.7 and 2.5 are considered, two correlated values are located near the minimum or maximum values for different QoIs, and RP. This means that if a better prediction of is required, should be set to 2.290, and if the RP is predicted well, should be set to 0.723. Further, , the second dominant parameter, is correlated differently and with a different sign. However, in the compression ramp flow, the values of are predicted similarly: 0.864 for SP and 0.713 for RP. When the default value of (1.75) is considered, the optimized values for SP and RP in the ramp flow are located near the minimum value (0.7).
Based on the correlated coefficients by Bayesian inference, the deterministic simulations with the GEKO model are conducted and the results are shown in Figures 18 and 19. The lines labeled UQ- and UQ-RP in the legend of Figure 18 present the results with the coefficients calibrated in terms of and RP, respectively. Similarly, UQ-SP and UQ-RP in Figure 19 adopted the correlated coefficients for SP and RP, respectively. As shown in Figure 18, the pressure distribution of UQ- shows the best agreement with the experimental data, with improvement over the results of the original GEKO-CC model and the UR-RP model since three model coefficients, , , and , are calibrated using observation data for Bayesian inference. However, the RP is predicted best by UQ-RP, whose coefficients are calibrated with the RP, even though the lowest peak of the reverse velocity is overpredicted in UQ-RP relative to the others. As for the 24-deg compression ramp case in Figure 19, UQ-SP, which is calibrated with the SP data, shows better agreement in the pressure distribution and skin friction than the original GEKO model. UQ-RP predicts a delayed SP but more accurate RP relative to the other results.
(a)
(b)
(a)
(b)
Tables 12 and 13 show the values of QoIs and the error between computational results and experimental data for each case. The values in the brackets indicate the percent error, and UQ- shows a 4% error for , which is an improvement of 2.5% over GEKO-CC. UQ-RP has a 7.5% error with respect to RP, whereas GEKO-CC shows a 25% error. In the compression ramp case, UQ-SP shows better agreement than GEKO-CC and UQ-RP, with a 3% error in the prediction of the SP. However, the error in RP is still high. UQ-RP is able to predict both SP and RP within 40%.
6. Conclusion
In the present study, the GEKO model was investigated to quantify the uncertainty of its coefficients for supersonic flows over an axisymmetric base and a 24-deg compression ramp. For the axisymmetric base flow, compressibility correction aided the models to predict the pressure along the base surface more accurately than the models without any correction. However, using the compressibility correction increased the radial variation of the pressure due to the increased centerline velocity, which also increased the distance to the shear layer RP.
As for the 24-deg compression ramp case, in which the compressibility correction was not applied, the GEKO model failed to predict the separation shock point and separation length. It even underpredicted pressure and skin friction in the relaxation region.
To overcome these issues and quantify the uncertainties, PC-NIPC and Bayesian inference were applied to the coefficients of the GEKO model, which can be tuned for specific flows. It was found that the proper value of the oversampling rate in Bayesian inference with PC-NIPC might be 3 than 2 in the forward problem with PC-NIPC. The LHS method was used for sampling random input variables with a uniform distribution. Based on the surrogate model constructed by OLS, the dominant coefficient was able to be figured through Sobol indices. was the most influential coefficient in both test cases. Observation data made of Gaussian distribution considering experimental measurement errors was used to calculate the likelihood function. The MCMC was calculated by using AIES for better convergence in cases where a strong correlation between parameters is expected. It was found that the calibrated coefficients were computed differently depending on the flow type and QoIs, even in the same flow type. After extracting calibrated coefficients from the posterior distribution, they were computed in the deterministic solver and the results clearly showed what QoIs they were calibrated for. The results with calibrated coefficients show better agreement than GEKO-CC for both test cases.
The results obtained from the GEKO model using calibrated coefficients based on Bayesian inference show superior prediction for two types of compressible flow than the results from the default coefficients of the GEKO model.
Data Availability
Data are available on request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was supported by the National Research Foundation of the Korean Government (NRF-2019R1F1A1059564).