Inscrição na biblioteca: Guest
International Journal for Uncertainty Quantification

Publicou 6 edições por ano

ISSN Imprimir: 2152-5080

ISSN On-line: 2152-5099

The Impact Factor measures the average number of citations received in a particular year by papers published in the journal during the two preceding years. 2017 Journal Citation Reports (Clarivate Analytics, 2018) IF: 1.7 To calculate the five year Impact Factor, citations are counted in 2017 to the previous five years and divided by the source items published in the previous five years. 2017 Journal Citation Reports (Clarivate Analytics, 2018) 5-Year IF: 1.9 The Immediacy Index is the average number of times an article is cited in the year it is published. The journal Immediacy Index indicates how quickly articles in a journal are cited. Immediacy Index: 0.5 The Eigenfactor score, developed by Jevin West and Carl Bergstrom at the University of Washington, is a rating of the total importance of a scientific journal. Journals are rated according to the number of incoming citations, with citations from highly ranked journals weighted to make a larger contribution to the eigenfactor than those from poorly ranked journals. Eigenfactor: 0.0007 The Journal Citation Indicator (JCI) is a single measurement of the field-normalized citation impact of journals in the Web of Science Core Collection across disciplines. The key words here are that the metric is normalized and cross-disciplinary. JCI: 0.5 SJR: 0.584 SNIP: 0.676 CiteScore™:: 3 H-Index: 25

Indexed in

FRAMEWORK FOR CONVERGENCE AND VALIDATION OF STOCHASTIC UNCERTAINTY QUANTIFICATION AND RELATIONSHIP TO DETERMINISTIC VERIFICATION AND VALIDATION

S. Maysam Mousaviraad

IIHR-Hydroscience and Engineering, The University of Iowa, Iowa City, Iowa 52242, USA

Wei He

Visiting scholar from NAOCE, Shanghai Jiao Tong University, Shanghai, China

Matteo Diez

Visiting scholar from CNR-INSEAN, Via di Vallerano 139, 00128 Rome, Italy

Frederick Stern

IIHR-Hydroscience and Engineering, The University of Iowa, Iowa City, Iowa 52242, USA

Abstract

A framework is described for convergence and validation of nonintrusive uncertainty quantification (UQ) methods; the relationship between deterministic verification and validation (V&V) and stochastic UQ is studied, and an example is provided for a unit problem. Convergence procedures are developed for Monte Carlo (MC) without and with metamodels, showing that in addition to the usual user-defined acceptable confidence intervals, convergence studies with systematic refinement ratio are required. A UQ validation procedure is developed using the benchmark UQ results and defining the comparison error and its uncertainty to evaluate validation. A stochastic influence factor is defined to evaluate the effects of input variability on the performance expectation and four possibilities are identified in making design decisions. The unit problem studies a two-dimensional airfoil with variable Re and normal distribution using high-fidelity Reynolds-averaged Navier–Stokes (RANS) simulations. Deterministic V&V studies achieve monotonic grid convergence and validation at the validation uncertainty interval of 2.2% D, averaged between lift and drag, with an average error of 0.25% D. For MC with Latin hypercube sampling the converged results are obtained with 400 computational fluid dynamics (CFD) simulations and are used as validation benchmark in the absence of experimental UQ. The stochastic influence factor is small such that the output expected value is not distinguishable from the deterministic solution. The output uncertainty is one order of magnitude smaller for lift than drag, implying that lift is only weakly dependent on Re. Several metamodels are used with MC, reducing the number of CFD simulations to a minimum of 4. The results are converged and validated at the average intervals of 0.1% for expected value (EV) and 10.7% for standard deviation (SD). The Gauss quadrature and the polynomial chaos (PC) method are validated using 8 and 7 CFD simulations, respectively, at the average intervals of 0.08% for EV and 7.5% for SD. The error values are smallest for the metamodels, followed by the PC method and then the Gauss quadratures.

KEYWORDS: deterministic verification and validation, nonintrusive stochastic uncertainty quantification, convergence criteria, UQ validation, RANS CFD simulation


1. INTRODUCTION

Uncertainty analysis is essential to guarantee and evaluate the fidelity of simulation-based design (SBD). Herein the focus is physic-based SBD using high-fidelity computational fluid dynamics (CFD). Until recently, uncertainty analysis focused on deterministic verification and validation (V&V) for deterministic numerical and modeling uncertainties. Recent focus has been stochastic uncertainty quantification (UQ) as an essential first step for robust and reliability-based design optimization. These methods focus on stochastic-based uncertainties, assessing the stochastic properties of simulation outputs including expected value (EV), standard deviation (SD), and cumulative distribution function (CDF) or probability density function (PDF) with respect to a set of uncertain input variables with known distribution functions (aleatory). For this purpose, the uncertain input variables for CFD are categorized as geometrical, operational, and numerical/modeling uncertainties [1]. Most studies are carried out for geometrical and operational uncertainties, while some consider numerical/modeling uncertainties, e.g., [2, 3]. Some studies (e.g., [4]) have considered both deterministic V&V and stochastic UQ. For complete reviews of the applications of UQ methods (intrusive and nonintrusive) in fluids refer to [5, 6]

The Mach number (Ma) variability is studied for NACA0012 airfoil at 3° angle of attack (AoA) [1]. The input distribution is truncated Gaussian distribution (99% interval) with the input mean and standard deviations of μMa = 0.3 (Reynolds number Re = 3M) and σMa = 5% μMa, respectively. The Reynolds-averaged Navier–Stokes (RANS) solver FINE/Hexa with Spalart–Allmaras turbulence model is used on a 50K grid. A probabilistic collocation method is used with five deterministic CFD solutions. The results show SD = 10% EV for both drag coefficient (CD) and lift coefficient (CL). The pressure and skin friction distributions on the airfoil surface are plotted along with 99% confidence level uncertainty bars. Transonic RAE 2822 airfoil at AoA = 2.79° with variable Ma, μMa = 0.734, σMa = 0.68% μMa (Re = 6.5M), and uniform distribution is studied [7]. The in-house RANS solver Joe with the Spalart–Allmaras turbulence model is used on a 112K grid. A stochastic collocation UQ method based on quadrature and polynomial interpolation is used with a maximum of 10 deterministic CFD simulations. The percent differences of EV versus deterministic solution at mean Ma are 0.5 for CD and 0.3 for CL. SD values for CD and CL are 4.06 and 0.7% EV, respectively. NACA0012 airfoil with variable AoA, μAoA = 3°, σAoA = 1.7% μAoA, and beta distribution is studied at high Ma = 0.73 (Re = 6M) [8]. The RANS solver elsA is used for deterministic simulations. The Monte Carlo (MC) method with N = 100 CFD solutions is used, and then improved by using a metamodel generating up to NM = 1M points. The polynomial chaos (PC) method is also used for 2, 3, 4, and 5 deterministic points. EV and SD values for CL with their 90%, 95%, and 99% confidence intervals from the metamodel are shown versus NM and compared against MC results with N = 100. The difference between MC with N = 100 and the metamodel with NM = 1M is 0.08% for EV and 0.96% for SD. The SD value for NM = 1M is 1.04% EV for CL. PC results are closer to NM = 1M results than MC with N = 100 for both EV and SD.

Multivariable studies are carried out in [1], and further in [9, 10] for the same problem. Truncated Gaussian distributions (99% interval) are assumed for all variables. The first multivariable study considers two geometrical input uncertainties at 3° AoA: relative thickness and maximum camber with 0.425% and 0.4472% coefficients of variation, respectively. The second study considers four joint UQ variables: Ma, AoA, relative thickness, and maximum camber on a 80K grid. The same conditions as the first study are considered for geometrical uncertainties, with μMa = 0.3, σMa = 10% μMa, = 5°, and σAoA = 10% μAoA. The UQ method uses 35 deterministic solutions based on the four-dimensional Halton sequence for probabilistic radial function approach. The EV for CL/CD is 1% lower than the deterministic value. For the surface pressure distribution, EV values with one SD upper and lower uncertainty bars are shown. The single-variable study in [7] is also extended to include a multivariable UQ for Ma (μMa = 0.734, σMa = 0.68% μMa), AoA (μAoA = 2.79°, σAoA = 3.58% μAoA), and thickness-to-chord ratio (μt/c = 0.1211, σt/c = 4.12% μt/c) with uniform distributions. The mean and variations are stated to be selected “based on expert opinions of realistic variations in practical operating conditions” [7]. A maximum of 93 = 729 deterministic CFD using 9 points for each variable is used. The PDFs are obtained and are nonuniform. The 95% uncertainty bars for surface pressure coefficient are shown, with asymmetric upper and lower bounds. The relative importance of each variable on EV and SD is also evaluated. Another study considered the same problem with normal distribution and two variables, Ma and AoA, using the RANS solver TAU with the k−ω turbulence model [11]. Additionally, a second strategy is considered where Ma and AoA variability are obtained based on two turbulent velocity vectors with normal distributions. Both MC and sparse Gauss–Hermite grid methods are used. MC results are obtained for up to N = 1500 (almost 39 points for each variable) for the first strategy and N = 17,000 (almost 130 points for each variable) for the second. SD values for CD and CL, respectively, are 14.6 and 2.0% EV for the first strategy and 15.9 and 1.4% EV for the second strategy. The sparse grid method does not converge for up to 647 deterministic simulations for the second strategy with error up to 190% against MC with N = 17,000 for SD of CL.

The AGARD 445.6 transonic wing is studied solving the inviscid steady Euler equations with the finite volume code CFL3D [12]. A C-H topology grid with 527K points is used. Two joint input variables are considered: Ma with uniform distribution between 0.8 and 1.1 (mean Ma = 0.95, coefficient of variation 9.1%) and AoA with uniform distribution between −2.0° and 2.0° (mean AoA = 0.0, standard deviation = 1.155). MC-LHS (Latin hypercube sampling) results are obtained with 1000 deterministic CFD solutions and the EV and SD of CD and CL are calculated with 95% confidence intervals using the Bootstrap method. PC results are then obtained for up to 5th degree (42 deterministic CFD solutions) with a point-collocation method and Hammersley sampling using twice the minimum number of points required and the least-squares method. The 5th degree PC results are within the 95% confidence interval of MC-LHS results. The authors later extended their work to study uncertainty on the first aeroelastic mode of the wing due to the same input variables and conditions [13]. A linear structural mode is considered for the modal aeroelastic analysis built in to the CFD code, and the uncertainty is evaluated for the damping factor used to reconstruct the transient aeroelastic response. PC results are obtained for up to 8th degree (90 deterministic simulations) and are within 95% confidence intervals of MC-LHS results with 1000 deterministic simulations. The same transonic aeroelastic wing with random free stream velocity is studied for unsteadiness extending the unsteady adaptive stochastic finite elements approach to resolve the effects of randomness [14]. The inviscid Euler equations were solved using an unstructured hexahedral mesh. The free stream velocity was assumed to follow a beta distribution. The results showed that although the mean free stream velocity was a safety margin of 5% below the deterministic flutter velocity, a 3.5% coefficient of variation still resulted in a nonzero flutter probability of 6.19%.

The convergence and validation of UQ results are not addressed quantitatively in the previous studies. In [8] the confidence intervals of EV and SD versus the number of MC samples are plotted and it is shown that the confidence intervals become negligible for very large samples. It is also shown that the results for MC with metamodels fall in the confidence intervals of MC with CFD points/items, without defining the error values or discussing the validation. For the MC method, quantitative convergence is defined in the literature as achieving an “acceptable” confidence interval for EV with a user-defined confidence level [15]. The random errors/uncertainties are estimated for certain metamodels/surrogate models using the training points [16]. For the PC method, seven different metrics are used to evaluate the accuracy of PC approximations against the exact solutions for the designed examples with analytical functions [17]. The accuracy of PC improved in some metrics as additional terms were retained, but did not exhibit this behavior in all metrics.

The relationship between CFD deterministic V&V and stochastic UQ is rarely studied. An error quantification procedure is developed in which both deterministic and stochastic errors are identified and evaluated, with sensitivity analysis to rank the contribution of each [16]. Three sources of error are recognized as numerical, model form, and uncertainty quantification errors, with the numerical error consisting of input, discretization, and surrogate model prediction errors. The deterministic errors are corrected by adding its value to the prediction, and the random errors are included in the model output through sampling. The corrected model outputs together with observed data are then used to estimate model form error through MC sampling, where MC error is also considered. The sensitivity analysis methods are provided for both local and global sensitivity measures.

The objective of the present research is to address issues of convergence and validation of UQ results and the relationship between CFD deterministic V&V and stochastic UQ and present an example for a unit problem. A framework is described that can be used for the convergence and validation of nonintrusive UQ methods. In addition to the usual statistical criteria of acceptable user-defined confidence intervals, convergence studies with systematic refinement ratio are included, similar to those in the deterministic verification [18]. Convergence procedures are provided for MC methods with and without metamodels focusing on EV, SD, CDF, their uncertainties, and the significance level of the chi-square test if the type of the output distribution is identifiable. For other UQ methods, the procedure simply studies the convergence of EV and SD versus the number of deterministic simulations. The UQ validation uses benchmark UQ results and defines the comparison error and its uncertainty to evaluate validation similar to the deterministic studies [19]. The procedure can be applied both to validate the UQ results from different methods against a benchmark result at the simulation or experimental level, and to validate the simulation-based UQ solutions against experimental UQ data. The relationship between the deterministic V&V and the stochastic UQ is studied in terms of a stochastic influence factor and its uncertainty. Four possibilities are identified in making design decisions. The unit problem is high-fidelity RANS simulation of the NACA0012 airfoil at 3° AoA and variable Re with normal input distribution.

2. DETERMINISTIC V&V METHODS

Deterministic verification and validation of CFD simulations is conducted at the individual user, code, model, grid type, etc. level for specified applications with available benchmark experimental validation data and uncertainties. It is assumed that all input variables/models are deterministic. The error estimates in deterministic V&V are based on fixed values; thus, the uncertainty estimates at the 95% level of confidence require similar reasoning as used for single realizations and zero-order-replication level systematic uncertainty estimation in experimental fluid dynamics (EFD). This is the smallest uncertainty interval that could be achieved with the CFD code for the specified application. For multiple realizations and Nth-order replication level uncertainty analysis in EFD, the systematic uncertainty is root sum squared (RSS) with the random uncertainty and multiplied by the coverage factor to obtain the total uncertainty. Note however that recent verification procedures have provided statistical evidence for the 95% level of confidence [20], and validation procedures take into consideration the numerical and experimental uncertainties.

Herein, the overall V&V approach in [18] along with detailed procedures for verification and validation in [20, 19], respectively, are followed. The methodology and resulting definitions are based on equations derived for deterministic simulation numerical and modeling errors and uncertainties, which provide the overall mathematical framework. Deterministic simulation modeling and numerical errors are assumed additive such that simulation uncertainties are root sum squared. For solutions in the asymptotic range, the estimated numerical error and its estimated error are used to obtain a corrected solution (numerical benchmark) and its uncertainty. Verification procedures identify the most important numerical error sources (such as iterative, grid size, and time-step errors) and provide error and uncertainty estimates. Validation methodology and procedures use benchmark experimental data and properly take into account both deterministic numerical and experimental uncertainties in estimating deterministic modeling errors and validation uncertainty, including the option of using corrected solutions. Simulation-based design requires verification of all simulations, but hopefully validation only for the final or unusual designs.

Verification procedures estimate numerical uncertainties (USN) based on iterative (UI), grid (UG), and time-step (UT) uncertainties:

(1)

UI is evaluated using the procedures described in [21]. Iterative convergence requires UI to be at least one order of magnitude smaller than UG and UT. Grid and time-step convergence studies are carried out for three solutions (S) with systematic refinement ratio

(2)

where 3, 2, and 1 represent the coarse, medium, and fine grids, respectively. A sufficiently large refinement ratio should be chosen since for small values (i.e., very close to 1.0) solution changes are small and sensitivity may be difficult to identify compared to iterative errors. Solution changes ε and the convergence ratio R are defined by

(3)

The convergence is evaluated as

(4)

When monotonic convergence is achieved, the Richardson extrapolation (RE) is used to estimate the order of accuracy pRE and error estimate δRE:

(5)

The factor of safety method [20] is used for estimations of UG and UT . For the solutions to be close to the asymptotic range, pRE needs to be close to the theoretical order of accuracy. However, the criteria for attainment of the asymptotic range are lacking. A possible criterion is that monotonic convergence should be established based on the convergence ratio R for S1 (toward Sc, the corrected solution), P (toward 1), UG, and UT for multiple (at least three) grid/time-step triplets with the same refinement ration r along with UIUG. In some cases it may be that damped oscillatory convergence is deemed acceptable for S1, P, UG, and/or UT, although this would require many grid triplets [22].

For oscillatory convergence, uncertainties are estimated simply based on oscillation maxima and minima, provided that the mean value is not drifting [18]. For nonsmooth or nonmonotonic convergence, the least-squares [23] or response-surface [24] methods may be used, which require solutions for more than three grids. However, there are some issues in using these methods as discussed in [20].

A validation procedure defines the comparison error (E) and the validation uncertainty (UV) using experimental benchmark data (D) and its uncertainty (UD). If UV bounds E, the combination of all the errors in D and S is smaller than UV and validation is achieved at the UV interval, which is a more stringent requirement than |E| ≤ UE:

(6)

If UV ≪ |E|, the sign and magnitude of E ≈ δSM can be used to make modeling improvements, where δSM is simulation modeling error.

3. STOCHASTIC UQ METHODS

3.1 Introduction

Stochastic UQ methods using CFD simulations are also conducted at the individual user, code, model, grid type, etc. level for specified applications, but in this case for stochastic input variables with known distributions (aleatoric) and typically without available benchmark EFD validation data and uncertainties. The focus is typically on output variable distribution parameters, namely EV, SD, and possibly CDF, along in some cases with their confidence intervals. Typical input variables are for geometrical and operational uncertainties, although some models may have stochastic numerical/modeling input variables [1].

The set of independent input stochastic variables is denoted by the vector ξ, with probability density function P(ξ). The stochastic properties of interest for the output function J(ξ) are [25]

(7)
(8)
(9)

where

(10)

The PDF is obtained from CDF as

(11)

In general, output functions evaluated using CFD (S) have deterministic (systematic) modeling and numerical errors (δS). Output functions evaluated using EFD (D) have both systematic and random error components (δD). For CFD solutions in the asymptotic range a corrected solution SC can be used along with its error estimate δSC [21]. The corrected values are used in the procedures and applications to the mechanical stress UQ study in [16], but the corrected solution error estimate is not included. Thus in general we have

(12)

where J and δJ represent the corrected or uncorrected CFD or the uncorrected EFD value and its error, respectively. Substituting Eq. (12) in Eqs. (7) and (8) and assuming that J and δJ are uncorrelated:

(13)
(14)

The UQ methods in Section 3 are described generally for any deterministic function of ξ, simply denoted as J, and may be applied straightforwardly when the δJ are systematic (epistemic) errors. For random (aleatoric) errors/uncertainties, the distribution parameters of δJ may be used at each ξ to generate samples based on the distribution function, and the methods of Section 3 may be then applied (see, e.g., [26]). Equations (13) and (14) indicate that a separate UQ procedure could be pursued for the error terms than the procedure for the output function itself. In that case, different number of items at different ξ values may be used for δJ than J, depending on the UQ methods chosen.

CFD deterministic numerical and modeling errors/uncertainties are systematic in nature and may or may not depend on the UQ independent variable ξ. For cases such as the present unit problem where they only weakly depend on ξ, their stochastic effects may be neglected. For other cases, Section 2 methods may be used to estimate δJ (ξ) values. EFD deterministic systematic and random errors also may or may not depend on ξ. For the random errors, however, the stochastic effects may not be neglected even if uncertainties are independent of ξ. Again the usual EFD uncertainty assessment procedures could be used to estimate δJ (ξ).

For the unit problem with RANS simulations, various nonintrusive UQ methods are employed as described in the following. Nonintrusive approaches are most attractive for high-fidelity computations since they require little modifications to the existing simulation code, and allow the deterministic simulations to run in parallel.

3.2 Monte Carlo with Latin Hypercube Sampling

The MC-LHS methodology follows [25, 27]. In MC analysis, a sampling procedure is used to construct a set of independent realizations of ξ from P(ξ). The collocation of the corresponding unique solutions is called the sample solution set:

(15)

The convergence rate of MC method for a randomly selected sample set is small, 1/√N. Latin hypercube sampling improves convergence by sampling from equiprobable bins in the parameter range

(16)

where Φ−1 is the inverse of the CDF function called the quantile function, and μ and σ are mean and standard deviation of the input variable, respectively.

The integrals in Eqs. (7)–(9) for an unbiased sampling, including random and LHS, can be approximated by

(17)
(18)
(19)

The y values in the current study are determined using k = √N number of bins in the sample range. After the cumulative histogram is calculated, PDF is obtained by numerical differentiation using a simple three-point finite difference formula.

The output distribution obtained by Eq. (19) needs to be identified since it is generally not the same as the input distribution. The chi-square (χ2) test may be used to identify the type of output distribution from well-known distributions such as normal, beta, gamma, etc. The test procedure is based on multinomial distribution, i.e., the number of deterministic solutions in the ith bin, denoted by xi, is the product of the probability of the ith group from the well-known distribution, pi, and the total number of simulations, N [28]. The quantity χ2 defined as

(20)

is shown to have the chi-square distribution with (k − 1) degrees of freedom. The hypothesized distribution is accepted if the χ2 value from Eq. (20) is smaller than or equal to the critical value, χ2k−1critical). The most usual value to assign αcritical is 0.05, that is, the test is conducted with 5% risk [28]. It is recommended that Npi should be at least as large as 5. For the current study, k = √N number of bins is considered initially, and then several groups with smaller probabilities are combined to satisfy the requirement if needed. If the test is passed, the actual significance level α is determined from the chi-square distribution with (k − 1) degrees of freedom using the χ2 value from Eq. (20).

After the type of output distribution is identified, one can estimate the 95% confidence interval of J, herein called output stochastic uncertainty UJ

(21)

where CL and CU are the lower and upper coverage factors and SD is the standard deviation of the output function obtained from Eq. (18). UJ is the uncertainty propagated to the output from the input stochastic uncertainty Uξ, which can be obtained similarly using Eq. (21) with input standard deviation σ. For normal distribution and large N, for example, the symmetric 95% confidence intervals are obtained with CU = CL = 2.0. If the type of output distribution cannot be identified, the coverage factors may be obtained using the Chebyshev inequality [29].

3.3 MC-LHS with Metamodels

A limited number of deterministic CFD solutions, M, are used to build the metamodels. Herein the M points are chosen from the deterministic solutions used for MC-LHS and therefore follow the LHS sample set as per Eq. (16). The sample solution set in Eq. (15) is then predicted by the metamodels. The MC-LHS procedure is used to obtain the UQ results following Eqs. (16)–(19).

Three global metamodels are used for this single-variable unit problem which show good results; therefore no local metamodels are employed. Two response surface functions, polynomials and power law, and a polynomial kernelbased learning model, least-squares support vector machines (LS-SVM), are used. The polynomial metamodels are used without least-squares fitting, while the other two metamodels are based on least-squares methods.

Polynomial metamodels with degree d are used with the following form:

(22)

For each polynomial degree, M = d + 1 deterministic points are used to obtain the coefficients in Eq. (22).

The power law response surface metamodel has the following form:

(23)

where a and k are constant, and o(xk) is an asymptotically small function of xk. A least-squares fitting is used to find a and k using M > 2 samples.

LS-SVM methodology follows [30], which finds the solution for support vector machines (SVMs) by solving a set of linear equations. SVMs basically map design variables x to a higher dimensional space by defining kernel function. The LS-SVM model for function estimation using a training set {xk, yk} is

(24)

where K(x, xk) is kernel function, herein a polynomial kernel of K(x, xk) = (xTkx + 1)d with d a tuning parameter, and αk and b solutions from the following system:

(25)

where y = [y1,…, yM], 1v = [1,…, 1], α = [α1,…, αM], Ωkl = K (xk, xl) : k, l = 1,…, M, and γ is an additional tuning parameter.

Note that the selection of metamodels depends on the specific problem as a large number of choices are available. For complex problems with a large number of uncertain variables, the construction of certain metamodels may be more expensive than the computational cost of a direct MC-LHS approach.

3.4 Quadrature Formulas

For estimation of the integrals in Eqs. (7) and (8), numerical integration techniques can be used based on a number of collocation points, I. Gauss quadrature, trapezoidal, and Simpson′s integration methods are used. For the current untruncated problem, good results are not expected with trapezoidal and Simpson′s, but these methods may be useful for finite functions such as truncated input distributions.

Gauss quadratures provide high-order accuracy for the integrations, with the only drawback requiring a node distribution that depends on I. This means that the deterministic solutions obtained for the MC-LHS method cannot be used for Gauss quadrature evaluations, and that a different set of CFD simulations is required for each I. The Gauss quadrature estimation in the case of normal distribution is

(26)

where ξi and wi are the quadrature nodes and the corresponding weights of the Gauss–Hermite polynomials, which can be found in [27] for example.

3.5 Polynomial Chaos Method

The PC methodology used herein follows [31]. One-dimensional PC expansion of order P is

(27)

where uk are the deterministic PC coefficients and Ψk(ξ) are the PC basis functions corresponding to the kth mode. For the basis functions, Hermite polynomials which are orthogonal in the random space are used. To find the uk values with a given PC order P, P + 1 deterministic simulations are used to form a linear system of equations based on Eq. (27). The P + 1 points are chosen from the deterministic points used for MC-LHS and therefore follow the LHS rule as per Eq. (16). After the polynomial coefficients are determined, the evaluations of statistical properties are carried out on JP instead of J as

(28)

4. FRAMEWORK FOR CONVERGENCE AND VALIDATION OF UQ METHODS

4.1 Convergence Studies for MC Results

For EV and SD estimates from MC method, their respective uncertainties UEV and USD are estimated at the 95% level of confidence. The confidence intervals are evaluated based on the methodology and procedures summarized in [32] assuming normal distribution. The procedure for EV assumes that the true value at the ε confidence level is bounded by EV − UEV and EV + UEV such that

(29)

UEV is evaluated using the Student′s t-distribution with (N − 1) degree of freedom:

(30)

The confidence level in the current framework is considered to be ε = 95%. For large N, t-distribution converges to normal distribution such that UEV = (2/ √N)SD.

For SD, the confidence interval is evaluated using the chi-square distribution such that for the 95% confidence level we have

(31)

where χ2(N−1);0.025 and χ2(N−1);0.975 are the 2.5% and 97.5% percentiles, respectively, of the χ2 distribution with (N−1) degree of freedom. Note that since the χ2 distribution is nonsymmetric, the upper and lower limits of the interval have different values:

(32)

The upper bound is larger for small N and converges to the lower bound value for large N. For the current study, the upper bound value is considered as the uncertainty, assuming large N:

(33)

Note that the procedures in Eqs. (30) and (33) are obtained assuming normal distribution. However, in general, the following conditions are possible:

  1. Normal distribution with large or small N: the procedures may be applied regardless of N.
  2. Non-normal distribution with large N: the procedures can be still used according to the central limit theorem. If the skewness of the distribution is small, at least N = 20 should be used for UEV and N = 50 for USD [32]. For the purpose of the current framework the central limit theorem is expected to almost always hold true since the χ2 test is not valid for small samples and the convergence criteria described later usually require a large number of samples.
  3. Non-normal distribution with small N: for this unlikely case, the Chebyshev inequality can be used which provides the confidence interval of expected value only [29].

The uncertainty in the empirical CDF is evaluated following [16]. The MC estimate of the empirical CDF from Eq. (19) follows a binomial distribution for any specific value y. The binomial distribution approaches the normal distribution for large N, e.g., N > 20. Now the 95% confidence bounds are evaluated similar to Eq. (21), which for normal distribution gives

(34)

Note that the CDF uncertainty is a function of y, and it goes to zero at both ends of the CDF curve. The maximum value max(UCDF) for CDF(y) = 0.5 is

(35)

In addition to the usual statistical criteria of acceptable user-defined confidence intervals, convergence studies are required similar to deterministic studies described in Eq. (3). The convergence studies are deemed necessary since UEV% EV1 could be very small before the solutions are converged as shown later for the current unit problem, and USD% SD1 and max(UCDF) are only functions of N and not sensitive to the solution changes.

For a single triplet convergence study, the MC samples for large, medium, and small N values are denoted as N1, N2, and N3, respectively, with a systematic refinement ratio:

(36)

The r value should be sufficiently large since for small values very close to 1.0 solution changes are small and sensitivity may be difficult to identify. According to the experience from deterministic grid studies, r = 2 or 21/2 is usually reasonable. The corresponding EV, SD, and α values estimated from the MC analyses are used to determine the convergence ratios. For EV, for example, we will have

(37)

The convergence criteria are described based on the four possible convergence conditions in Eq. (4) for a single triplet study as a minimum requirement, although multiple triplet study is desirable to establish asymptotic range. Divergence (monotonic or oscillatory) is not acceptable, i.e., when the absolute value of the convergence ratio R is larger than 1.0. Monotonic convergence (0 < R < 1) is preferred as a precondition to approaching asymptotic range, but if ε21 is very small (usually a couple of magnitudes smaller than the value itself), oscillatory convergence (−1 < R < 0) may have to be accepted. This provides the MC convergence criteria as

(38a)
(38b)

For the user-defined values, the choices are associated with the degree of risk considered acceptable for each problem. For the current unit problem for example, UEV% EV1 smaller than 1% and USD% SD1 and max (UCDF) × 100 values smaller than 10% are specified.

Note that if some of the criteria in Eq. (38) are not applicable, for example when the output distribution does not follow any standard distribution and α is not available, the convergence is still evaluated based on the other criteria.

4.2 Convergence for MC with Metamodels

Since metamodels are used to predict the output functions using a limited number of deterministic points M, the goodness of predictions needs to be evaluated. A few metrics are used in the literature to study the fitting accuracy [33]. Herein the unbiased average absolute error is used:

(39)

where l is the cross-validation set excluding the training points. If the total number of deterministic simulation points are T, M of which used to train the metamodel, the number of points for the cross-validation set is l = TM. The training points are excluded in the cross validation to avoid a favorable bias. Acceptably small |E|, for example 0.01 (less than 1% prediction error) for the current unit problem, is required for all metamodels. Also, convergence of |E| versus M needs to be studied following Eqs. (36) and (37). After a good fit is obtained, convergence of MC results versus N is required following the criteria in Eq. (38). Finally, convergence studies of EV, SD, and α estimates versus M are carried out following Eqs. (36) and (37).

In summary, the procedure for convergence studies of MC with metamodels requires a minimum of three training sets based on M1, M2, and M3 number of points, each with a minimum of three predicted MC points N1, N2, and N3. However, note that again multiple triplet studies are desirable to establish asymptotic range. The minimum convergence criteria for MC with metamodels include

(40a)
(40b)
(40c)
(40d)

Note that similar to Eq. (38), if some criteria are not applicable the convergence is still evaluated based on the remaining criteria. For the current unit problem, for example, each polynomial degree (d) is trained with a constant number of points (M = d + 1) and the criteria for convergence versus M, i.e., Eqs. (40b) and (40d), are not applicable. For some problems, evaluation of a fitting error may not be feasible due to the computational costs. The random uncertainty may be estimated for some metamodels (see, e.g., [16]) and used to evaluate the goodness of predictions.

4.3 Validation of UQ Results

Validation is evaluated for EV and SD using benchmark results EVB and SDB. The comparison error values are defined as

(41)

Similar to Eq. (6) for deterministic V&V studies, the uncertainties in the comparison error values in Eq. (41) are obtained and used as validation uncertainties:

(42)

Note that for small N, the upper and lower bounds of USD are different, as per Eq. (32). For convergence studies, the upper bound is considered since it is larger. For validation studies, however, the smaller bound may be used to impose more strict criteria if N is small.

Validation is achieved at the “noise level” UV imposed by both benchmark and solution uncertainties, if the absolute value of E is less than UV:

(43)

Ideally we would have the converged experimental UQ values as a benchmark to validate CFD UQ results including MC. However, since experimental UQ is typically not available, the converged MC results may be used as a benchmark to evaluate the other UQ methods. For MC with metamodels, UEV and USD are estimated along with UBEV and UBSD to obtain the validation uncertainties. For quadrature formulas and PC, confidence intervals are not defined and validation is carried out only based on benchmark uncertainties. This results in smaller UV and therefore more stringent validation criteria, as per Eq. (43).

5. RELATIONSHIP BETWEEN DETERMINISTIC V&V AND STOCHASTIC UQ

Deterministic V&V and stochastic UQ using CFD simulations are different since the former is based on fixed input and output values with focus on simulation S and its comparison error E and numerical USN and validation UV uncertainties, whereas the latter is based on stochastic input and output values with the usual focus on EV, SD, and CDF. Additionally, for stochastic UQ Section 4 describes a framework for evaluating convergence and validation of UQ methods.

The relationship between deterministic V&V and stochastic UQ can be assessed in terms of a stochastic influence factor λ and its uncertainty Uλ:

(44)
(45)

The simulation uncertainty US can be estimated by USN or UV depending on whether deterministic modeling uncertainties are considered or not. Note that UEV and US are independently evaluated using different uncorrelated procedures such that even if the same CFD code or experimental equipment is used for both of their evaluations there is no cancellation of errors.

In making design decisions there are several possibilities of interest:

  1. |λ| < Uλ
    1. Uλ < user-defined value: Uλ acceptable and negligible stochastic influence on the performance expectation
    2. Uλ > user-defined value: Uλ needs to be reduced in order to insure negligible stochastic influence on the performance expectation
  2. |λ| > Uλ
    1. Uλ < user-defined value: Uλ acceptable and UQ needed, i.e., EV is a better estimate of performance than S
    2. Uλ > user-defined value: UQ needed but Uλ needs to be reduced

Reductions in Uλ must consider the relative importance of UEV versus US and additionally confidence issues associated with levels of convergence and validation as discussed in Section 4. The assessment will most likely be done using simulations; since UQ experiments are largely not available. As discussed in Section 3.1, if the error values δJ* are included then λ, EV, and UEV would include additional terms resulting from Eqs. (12)–(14). In particular, UEV would include a contribution from the stochastic numerical and modeling errors/uncertainties.

6. UNIT PROBLEM

The unit UQ problem studies the effects of stochastic Re variations at AoA = 3° on CD and CL output variables for NACA0012. The mean Re value is μRe = 1.7 M, according to the tripped experimental data, with a standard deviation of σRe = 10% μRe and nontruncated normal distribution.

Tripped [34] and untripped [35] wind tunnel experimental studies are available. The results show limited transition effects at AoA = 3° while variations are large and nonlinear at high AoA (especially around critical AoA). Small AoA = 3° is selected to avoid the critical conditions and reduce deterministic computational costs accordance with previous UQ studies. The tripped data are applicable for the current study since no transition model is used in the CFD code and the small Ma = 0.3 ensures negligible compressible effects.

7. CFD METHODS

Incompressible RANS computations are performed in the relative inertial coordinates with the CFDShip-Iowa code [36] using a blended k−ε/k−ω shear stress transport (SST) turbulence model for fully turbulent simulations according to the tripped data. The grid used for the UQ simulations has 357K grid points with γ+ = 0.35 at the mean Re. Numerical methods incorporate finite difference discretization, with a second-order upwind scheme for convection terms and second-order centered for viscous terms. Temporal terms are discretized using a second-order backward Euler scheme. Incompressibility is enforced by strong pressure/velocity coupling, using a projection algorithm. Each deterministic CFD run takes 20 min clock time and 10 CPU hours on the HPCMP IBM Power 6 machine.

To accomplish the large number of high-fidelity simulations required by UQ methods, a Python program is developed to automate and parallelize the CFD simulations. The program creates directory, name-list file, boundary condition file, and job submission file and then submits the job for each simulation by reading Re values from an input file. A separate Python program is developed for postprocessing the CFD results and is run after all deterministic simulations are completed. This program modifies, compiles, and runs a Fortran code for each simulation to calculate the converged force components.

8. RESULTS

8.1 Deterministic V&V

Deterministic V&V is performed at the mean Re = 1.7M. Grid convergence studies are conducted for five grids (three grid triplets) ranging from 22.2K to 5.67M grid points with a grid refinement ratio of r = 2.0. An O-type grid topology is used and the grid extends to 15 times the chord length in each direction, as shown in Fig. 1. The γ+ values for all grids are reported in Table 1. To achieve two-dimensional results, five grid points are specified in the z direction extending one chord length. Zero-gradient boundary conditions are used for lateral boundaries and far-field boundary conditions are specified for the far-field boundaries.

FIG. 1: O-type grid design for NACA0012 shown for grid #3.


TABLE 1: Summary of grids used for the deterministic V&V study at AoA = 3°

Grid Dimensions y+ Computation time (min) on IBM P6 32 processors CD CL E% D
x × y × z Size D = 0.0128 D = 0.303 CD CL
1 1610 × 704 × 5 5.67M 0.1 750 0.0128 0.304 0.20 0.30
2 805 × 352 × 5 1.42M 0.18 110 0.0129 0.305 0.61 0.63
3 405 × 176 × 5 356.4K 0.35 20 0.0130 0.308 1.63 1.59
4 201 × 88 × 5 88.4K 0.7 4 0.0132 0.312 3.09 2.81
5 101 × 44 × 5 22.2K 1.4 2 0.0151 0.350 18.51 15.51

Iterative convergence studies are conducted with the flow traveling 50 chord lengths. The UI values for all triplets are at least one order of magnitude smaller than ε12. The UI / UG ratios are small (between 0.007 and 0.054) such that the iterative uncertainties are negligible compared with the grid uncertainties and are not considered in the validation studies.

Grid verification studies show that for both CD and CL monotonic convergence is achieved for all grid triplets, as per Table 2. The solutions are closer to the asymptotic range (P = 1) and grid uncertainties are smaller (between 0.3 and 1.6% S1) compared to the authors′ experience for most of the previous ship hydrodynamics problems. The multiple grid-triplet results are studied considering RS1, RP, and RUG. For S1, monotonic convergence is achieved for both CD and CL, i.e., the grid triplet (1, 2, 3). The P values monotonically converge for CD, but are monotonically divergent for CL. The UG values monotonically diverge for both CD and CL. Therefore the asymptotic range is not yet established and finer grids are required, as can also be seen in Fig. 2. The grid convergence studies are also carried out for the minimum and the maximum Re in the N = 400 LHS range, using the grid triplet (3, 4, 5). The results (not included) show that R, P, and UG versus Re are almost constant (less than 1% variation). The reason is perhaps that for this thin boundary layer and lift dominated problem with no transition effects, the flow physics does not change significantly in the current range of Re. This verifies that neglecting δJ* (ξ) in Eq. (12) is a reasonable assumption for this unit problem.

TABLE 2: Deterministic grid verification and validation results

Parameter Grids Refinement ratio UI / ε21 R P RP UG% S1 RUG UI / UG UD% D UV% D E% D
CD 3, 4, 5 2 0.096 0.08 1.80 0.94 1.47 6.58 0.007 2.28 2.71 1.63
2, 3, 4 2 0.164 0.17 1.27 1.50 0.013 2.28 2.73 0.61
1, 2, 3 2 0.203 0.34 0.77 0.38 0.022 2.28 2.36 0.20
CL 3, 4, 5 2 0.103 0.10 1.69 1.14 1.60 1.86 0.008 1.99 2.55 1.59
2, 3, 4 2 0.483 0.18 1.26 1.18 0.039 1.99 2.32 0.63
1, 2, 3 2 0.505 0.34 0.77 0.31 0.054 1.99 2.02 0.30

FIG. 2: Grid convergence results at Re = 1.7M and AoA = 3°.

The validation results show that validation is achieved for all grid triplets. The comparison error values averaged between CD and CL are 1.61% D for (3, 4, 5), 0.62% D for (2, 3, 4), and 0.25% D for (1, 2, 3). The average UV intervals are 2.63% D for (3, 4, 5), 2.52% D for (2, 3, 4), and 2.19% D for (1, 2, 3). The average UD / UG ratio is 1.4 for the grid triplet (3, 4, 5) and increases to 6.2 for (1, 2, 3), such that reducing UV requires a reduction in UD. Nonetheless, the UD values are relatively small (average 2.1%) compared to the ship hydrodynamics experiments in the recent G2010 workshop [37]. The velocity and pressure contours obtained with grid #3 are shown in Fig. 3.

FIG. 3: Pressure (left) and velocity (right) contours at Re = 1.7M and AoA = 3°.

8.2 MC-LHS Results

The convergence is studied for MC-LHS results by systematically increasing the number of deterministic simulations N with r = 2 until the convergence criteria in Eq. (38) are satisfied. The results are shown in Table 3 for N ≥ 50. The usual criteria of reasonably small uncertainties are studied first as per Eq. (38a). As mentioned earlier, the USD% SD and max(UCDF)% values are only functions of N and do not depend on the solutions, with USD% SD always being greater than max(UCDF)%. The USD% SD values are 25.5% for N = 50 and 7.53% for N = 400 for both CD and CL. The values for max(UCDF)% are 14% for N = 50 and 5% for N = 400. The UEV% EV values are reasonably small for all MC-LHS studies with N ≥ 50. The average UEV% EV values are 0.2% for N = 50 and 0.08% for N = 400. The convergence studies are conducted based on Eq. (38b) and the R values are reported in Table 3. The R value for each N corresponds to the triplet for which N is the largest number of deterministic points. Note that χ2 tests could not be performed for small N < 50. For EV and SD, monotonic convergence is achieved for N ≥ 200 deterministic CFD simulations, while χ2 test results monotonically converge for N = 400. Since monotonic convergence is achieved, the ε12 values do not need to be examined.

TABLE 3: Convergence studies for the MC-LHS results

Output function N Values Uncertainties Convergence Studies
EV SD α
EV SD% EV UEV% EV USD% SD max(UCDF)% R E% EVB R E% SDB α R E% αB
CD 50 0.013154 1.35 0.38 25.49 14.14 0.49 0.011 0.51 5.92 96.55 0.99
100 0.013155 1.40 0.28 16.53 10.0 −0.97 0.005 −0.67 2.89 96.87 0.67
200 0.013156 1.42 0.20 11.05 7.07 0.60 0.002 0.59 1.09 97.35 1.50 0.17
400 0.013156 1.44 0.14 7.53 5.0 0.59 0.000 0.60 0.00 97.52 0.35 0.00
CL 50 0.307914 0.17 0.05 25.49 14.14 0.47 −0.0014 0.51 5.83 94.73 0.90
100 0.307912 0.18 0.04 16.53 10.0 −1.02 −0.0007 −0.65 2.84 94.65 0.98
200 0.307911 0.18 0.03 11.05 7.07 0.56 −0.0002 0.59 1.06 95.23 −7.25 0.38
400 0.30791 0.18 0.02 7.53 5.0 0.61 0.0000 0.60 0.00 95.59 0.62 0.00

The convergence trends for EV, SD, and α versus N are plotted in Fig. 4 for CD and CL. Although monotonic convergence is achieved, the solutions are still changing versus N for SD and more significantly for α. Therefore the asymptotic values are not reached yet and larger N is required. The average difference for N = 200 versus 400 is 0.001α for EV, 1.07% for SD, and 0.28% for α.

FIG. 4: UQ convergence trends for MC-LHS results.

For output distributions, the histogram analysis results and the Npi values for each bin of the χ2 tests are shown in Fig. 5. Both CD and CL distributions are found closest to normal, although not exactly normal since the average confidence level is α = 96%. The normal input and output distributions for this single-variable problem imply that the output functions are close to linear, as will be shown later with the metamodels. Studies of the other force components are carried out (not included) and show that Cfx, Cpx, and Cpy are close to the normal, while Cfy is close to the gamma distribution.

FIG. 5: Histogram and Npi values (normal distribution) for MC-LHS χ2 test with N = 400.

The relationship between deterministic V&V and MC-LHS stochastic UQ results is studied in Table 4. Since the UQ simulations are performed with grid #3, the grid triplet (3, 4, 5) is selected for the relationship studies. The λ values are 0.09% S for CD and −0.01% S for CL, implying that the Re variability causes an increase in CD and a decrease in CL. For both CD and CL, the λ values are smaller than both UEV and USN and therefore |λ| < Uλ, i.e., condition 1a or 1b. The average Uλ values are 1.5% S for US = USN and 2.6% S for US = UV, which are reasonably small. This implies condition 1a, i.e., negligible stochastic influence on the performance expectation. USN/UEV ratios are 10.5 for CD and 80 for CL, such that reduction in USN is required to reduce Uλ when US = USN. Similarly, UV needs to be reduced for reduction in Uλ when US = UV, i.e., both UD and USN need to be reduced. Comparing CD and CL, the λ and UEV values are one order of magnitude smaller for CL, while the USN, UV, and Uλ values are of the same order. The smaller UEV for CL is correlated to its smaller UJ, which implies that CL is less dependent on Re variations than CD. This is consistent with findings of [7] for Ma variability. The opposite is expected for variable AoA.

TABLE 4: Relationship between deterministic V&V and stochastic UQ with converged MC-LHS

λ % S UEV % S USN % S UD % S UV % S Uλ = (√ U2EV + U2SN) % S Uλ = (√ U2EV + U2V) % S UJ % EV
CD 0.09 0.14 1.47 2.28 2.71 1.48 2.72 2.88
CL −0.01 0.02 1.60 1.99 2.55 1.60 2.55 0.37
Avg. 0.05 0.08 1.53 2.13 2.63 1.54 2.63 1.63

Figure 6 shows the deterministic V&V and stochastic UQ results. Deterministic values include the N = 400 CFD solutions used for MC-LHS with grid #3 (S3), CFD solution at the mean Re with grid #1 (S1), and EFD (D) value at the mean Re, with the deterministic uncertainty bounds USN−3, USN−1, and UV. The figure shows that the UV bars (for S3) bound the D points such that validation is achieved for both CD and CL. The S1 solutions are closer to the D values, with an average error of 0.25% D between CD and CL, than S3 solutions with an average error of 1.61% D. The uncertainties are also smaller for S1 with an average UG1 = 0.3% S1 than S3 with an average UG3 = 1.5% S3. The stochastic UQ results include EV, its uncertainty UEV, and output uncertainty UJ. The UEV values are very small compared to the US values such that Uλ is almost the same as USN or UV (average 0.05% difference) and are not shown in Fig. 6. The EV and S3 values are very close with average λ = 0.05% S. The Uλ bounds are obviously much larger, with average 1.5% S for US = USN and 2.6% S for US = UV. The stochastic influence on the performance expectation is therefore insignificant for this unit problem. It is expected as per the V&V studies at the minimum and maximum Re that if grid #1 was used for the UQ studies, the EV values would still be close to the deterministic values such that the same condition 1a is obtained. The S3 values versus Re show that the output function is steeper for CD than for CL. The UJ bound is 2.9% S for CD, almost the same order as UV, while for CL it is only 0.4% S and is 86% smaller than UV. Therefore for CL, the output uncertainties due to the input Re variations are less significant than the deterministic numerical and modeling uncertainties.

FIG. 6: Deterministic V&V and stochastic UQ results: MC-LHS for CD and CL of NACA0012 at 3° AoA and variable Re with μRe = 1.7M and σRe = 10% μRe.

8.3 MC-LHS with Metamodels Results

The results for polynomial metamodels are shown in Table 5. Convergence studies are conducted based on Eq. (40) excluding the convergence versus M, i.e., Eqs. (40b) and (40d), since each polynomial degree (d) is trained with a constant number of points (M = d + 1). The |E| values are smallest for 3rd to 6th degree polynomials with average |E| = 0.001, while for 1st degree (linear) for example the average is |E| = 0.34. The MC-LHS analyses require N = 400 prediction points with all polynomials to satisfy the convergence criteria, as expected. The resulting UEV and USD values are comparable to the benchmark values UBEV and UBSD. The validation studies show that the 1st degree polynomial is not validated for EV, while all other polynomials are validated at the average intervals of UVEV = 0.2% EVB for CD, UVEV = 0.025% EVB for CL, and UVSD = 10.6% SDB for both CD and CL. The error values are larger for SD than EV, with an average EEV = 0.04% EVB and ESD = 0.9% SDB. The UQ errors are almost correlated with the fitting errors, i.e., smallest for 3rd to 6th degree polynomials.

TABLE 5: Validation studies for polynomial metamodels

Output function Polynomial degree M |E| N EV SD
UEV % EV UVEV % EVB E % EVB USD % SD UVSD % SDB E % EVB
CD 1 2 0.6418 400 0.1458 0.2050 −0.63878 7.53 10.65 −1.89
2 3 0.0305 400 0.1473 0.2060 −0.00291 7.53 10.65 −2.24
3 4 0.0024 400 0.1441 0.2037 −0.00217 7.53 10.65 −0.03
4 5 0.0004 400 0.1441 0.2037 0.000001 7.53 10.65 −0.02
5 6 0.001 400 0.1444 0.2039 0.00011 7.53 10.65 −0.22
6 7 0.0026 400 0.1438 0.2036 0.00158 7.53 10.65 0.14
7 8 0.0025 400 0.1436 0.2034 0.00026 7.53 10.65 0.32
8 9 0.0102 400 0.1422 0.2024 −0.00016 7.53 10.65 1.31
CL 1 2 0.0769 400 0.01859 0.02612 0.07653 7.53 10.65 −1.20
2 3 0.0028 400 0.01863 0.02615 −0.00019 7.53 10.65 −1.52
3 4 0.0004 400 0.01834 0.02595 −0.00031 7.53 10.65 0.05
4 5 0.0003 400 0.01832 0.02593 0.00005 7.53 10.65 0.16
5 6 0.0005 400 0.01832 0.02593 −0.00015 7.53 10.65 0.14
6 7 0.0009 400 0.01831 0.02592 0.00074 7.53 10.65 0.22
7 8 0.0038 400 0.01908 0.02647 −0.00003 7.53 10.65 −3.95
8 9 0.0007 400 0.01826 0.02589 −0.00058 7.53 10.65 0.50

Results for the power law and the LS-SVM metamodels are shown in Table 6. For the |E| values, monotonic convergence versus M requires at least M = 16 training points for the power law and M = 28 for the LS-SVM. The average fitting errors with the same M = 28 points are 0.005 for the power law and 0.009 for the LS-SVM. The MC-LHS analyses require N = 400 predictions to achieve convergence, similar to the polynomial metamodels. The EV values monotonically converge with M = 16, except for CL with the LS-SVM metamodel, which requires at least M = 28. The SD values require at least M = 28 to achieve monotonic convergence, except for CL with the power law metamodel, which converges with M = 16. The converged results are validated at the average intervals of UVEV = 0.2% EVB for CD, UVEV = 0.02% EVB for CL, and UVSD = 10.6% SDB for both CD and CL. The average error for EV is 0.01% EVB, while the average error for SD is 0.3% SDB. Comparing the two metamodels, the EV values are better predicted with the LS-SVM, while the SD results are closer to the benchmark with the power law metamodel.

TABLE 6: Convergence and validation studies for power law and LS-SVM metamodels

Metamodel Output function M Fitting N EV SD
|E| R R UEV % EV UVEV % EVB E % EVB R USD % SD UVSD % SDB E % SDB
Power law CD 3 0.027 400 0.14404 0.20370 −0.027 7.53 10.65 −0.03
7 0.012 400 0.14409 0.20374 −0.010 7.53 10.65 −0.05
16 0.01 0.19 400 0.25 0.14408 0.20373 −0.006 −0.58 7.53 10.65 −0.04
28 0.008 0.6 400 0.74 0.14405 0.20373 −0.003 0.92 7.53 10.65 −0.02
CL 3 0.007 400 0.01832 0.02593 0.0071 7.53 10.65 0.177
7 0.003 400 0.01834 0.02594 0.0025 7.53 10.65 0.080
16 0.003 0.14 400 0.20 0.01835 0.02595 0.0016 0.96 7.53 10.65 −0.013
28 0.002 0.79 400 0.97 0.01835 0.02595 0.0007 0.02 7.53 10.65 −0.015
LS-SVM CD 3 0.03 400 0.1473 0.2060 −0.0028 7.53 10.65 −2.24
7 0.026 400 0.1467 0.2056 −0.0012 7.53 10.65 −1.85
16 0.02 1.63 400 0.50 0.1458 0.2049 −0.0004 1.70 7.53 10.65 −1.20
28 0.016 0.53 400 0.06 0.1452 0.2045 −0.0003 0.60 7.53 10.65 −0.81
CL 3 0.0028 400 0.01863 0.02615 −0.00021 7.53 10.65 −1.52
7 0.0024 400 0.01858 0.02611 −0.00018 7.53 10.65 −1.24
16 0.002 1.25 400 4.69 0.01851 0.02606 −0.00003 1.36 7.53 10.65 −0.85
28 0.0017 0.61 400 0.17 0.01846 0.02603 −0.00006 0.72 7.53 10.65 −0.58

Overall the smallest fitting error with metamodels is obtained for the 4th degree polynomials with the following equations:

(46)

where x = (Re − 1.7 × 106)/(1.7 × 106), CD = ( + 10.1)/800, and CL = ( + 87.4) × 0.0035. The higher order terms in Eq. (46) are very small compared with the first orders, such that the functions are close to linear, as can also be seen in Fig. 6. This explains the results of the χ2 tests showing that output distributions are close to the normal.

8.4 Quadrature Results

The validation studies for the Gauss quadratures are shown in Table 7. Note that LHS points cannot be used and additional CFD simulations are carried out for each I. Validation is achieved at the same intervals of UVEV = 0.14% EVB for CD, UVEV = 0.02% EVB for CL, and UVSD = 7.53% SDB for both CD and CL, for all the quadrature rules. The error values are of the same order for all I values, with an average EEV = 0.001% EVB and ESD = 1.1% SDB.

TABLE 7: Validation studies for Gauss quadratures

Output function I EV SD
E % EVB UVEV % EVB E % SDB UVSD % SDB
CD 2 −0.00128 0.14 −0.03 7.53
3 −0.00230 0.14 −1.31 7.53
4 −0.00224 0.14 −1.34 7.53
5 −0.00210 0.14 −1.32 7.53
6 −0.00233 0.14 −1.31 7.53
7 −0.00232 0.14 −1.34 7.53
8 −0.00224 0.14 −1.30 7.53
CL 2 0.00007 0.02 −0.32 7.53
3 0.00025 0.02 −1.31 7.53
4 0.00028 0.02 −1.29 7.53
5 0.00028 0.02 −1.21 7.53
6 0.00032 0.02 −1.23 7.53
7 0.00027 0.02 −1.28 7.53
8 0.00039 0.02 −1.30 7.53

The Gauss quadrature integrations are also carried out using the LS-SVM metamodel with M = 28 to predict the integral points (not included). The average error values are 7.06% larger for the SD of CD and 0.54% larger for the SD of CL than without the metamodel; while for EV the differences are not significant.

The trapezoidal integrations are carried out (not included) with the systematic refinement ratio r = 2 for I. The results show that monotonic convergence is achieved with I = 200. The error values, however, are relatively large especially for the SD of CL with E = 77% SDB. The errors remain large for I = 400, where validation is achieved only for the SD of CD with 7.5% interval.

For the standard Simpson′s integration, LHS points cannot be used since equally spaced abscissas are required. To avoid additional CFD simulations, the LS-SVM metamodel with M = 28 is used to predict the required points. The convergence studies with r = 2 show the same trend as for the trapezoidal integrations, i.e., monotonic convergence is achieved (with I = 33), but the error values are relatively large especially for SD of CL (with E = 90% SDB). The errors remain large for I values up to 513, where validation is achieved only for the SD of CD with 7.5% interval. This is perhaps due to the finite integration range neglecting the tails of the infinite functions.

8.5 PC Results

The PC validation results for P = 2−10 are presented in Table 8. Validation is achieved for all PC orders except for EV with P = 10 and SD with P = 9 and 10. Similar to the Gauss quadratures, the validation intervals are the same for all PC orders, i.e., UVEV = 0.14% EVB for CD, UVEV = 0.02% EVB for CL, and UVSD = 7.53% SDB for both CD and CL. The error values for both CD and CL are smallest for P = 6, with the average EEV = 0.001% EVB and ESD = 1.2% SDB. This shows that higher order PC may not always achieve more accurate results, consistent with the polynomial metamodels in Table 5.

TABLE 8: Validation studies for polynomial chaos

Output function P No. det. CFD EV SD
E % EVB UVEV % EVB .5 E % SDB UVSD % SDB
CD 2 3 −0.0053 0.14 −3.48 7.53
3 4 −0.0045 0.14 −1.36 7.53
4 5 −0.0024 0.14 −1.35 7.53
5 6 −0.0023 0.14 −1.51 7.53
6 7 −0.0011 0.14 −1.23 7.53
7 8 −0.0022 0.14 −1.35 7.53
8 9 −0.0025 0.14 −1.99 7.53
9 10 −0.0018 0.14 −302.02 7.53
10 11 0.3528 0.14 −2355.32 7.53
CL 2 3 0.00010 0.02 −2.75 7.53
3 4 −0.00003 0.02 −1.24 7.53
4 5 0.00033 0.02 −1.13 7.53
5 6 0.00014 0.02 −1.15 7.53
6 7 0.00084 0.02 −1.09 7.53
7 8 0.00025 0.02 −7.70 7.53
8 9 −0.00012 0.02 −1.72 7.53
9 10 0.00015 0.02 −292.51 7.53
10 11 0.25881 0.02 −13706.16 7.53

8.6 Comparison of UQ Methods

Table 9 summarizes the results of the different UQ methods for the current unit problem. For the polynomial metamodels, the same degree as PC (6th degree) is selected for comparison even though the 4th degree polynomial is slightly better. For the Gauss quadratures, the results for I = 8 are selected, since I = 2 seems unreasonably small. For the Simpson′s integrations, the results for I = 33 are selected after which the errors remain almost constant.

TABLE 9: Validation studies for polynomial chaos

UQ method No. det. CFD EV SD
|E| % EVB |E| % SDB
CD CL Avg. CD CL Avg.
MC-LHS N = 400 400 0.0 0.0 0.0 0.0 0.0 0.0
MC-LHS with polynomial metamodel M = 7, N = 400 7 0.0016 0.0007 0.0012 0.00.14 0.22 0.18
MC-LHS with power law metamodel M = 28, N = 400 28 0.0028 0.0007 0.0018 0.02 0.02 0.02
MC-LHS with LS-SVM metamodel M = 28, N = 400 28 0.0003 0.0001 0.0002 0.81 0.58 0.70
Average for metamodels 0.0016 0.0005 0.0010 0.32 0.27 0.30
Polynomial chaos method P = 6 7 0.0011 0.0008 0.0010 1.23 1.09 1.16
Gauss quadratures I = 8 8 0.0022 0.0004 0.0013 1.3 1.3 1.30
Trapezoidal integration I = 400 400 0.4787 0.4768 0.4778 4.02 177.46 90.74
Simpson′s integration with LS-SVM metamodel I = 33, M = 28 28 0.5032 0.5015 0.5024 5.35 190.25 97.80

For the MC-LHS results with different metamodels, the same N = 400 is required for MC convergence. The convergence studies versus M are conducted for the power law and the LS-SVM metamodels, for both of which at least M = 28 training points are required.

Validation is achieved for all the methods except for the trapezoidal and Simpson′s integrations, which are therefore excluded in the comparison discussions hereafter. The error values are relatively smaller for the metamodels, but none of them can be considered superior. For EV the average error values are 0.00019% for the LS-SVM and 0.001% EVB for the other two metamodels. For SD the average error values are 0.02% for the power law, 0.18% for the polynomial, and 0.7% SDB for the LS-SVM. The average error values between the metamodels are considered for the following comparisons. For EV the average error values are 0.001% for the metamodels and the PC method and 0.0013% EVB for the Gauss quadratures. For SD the average error values are 0.3% for the metamodels, 1.2% for the PC method, and 1.3% SDB for the Gauss quadratures.

Overall, the error values for the unit problem are smallest for the metamodels, followed by the PC method and then the Gauss quadratures. The numbers of CFD simulations are 7 for the polynomial metamodel and the PC method, 8 for the Gauss quadratures, and 28 for the power law and the LS-SVM metamodels.

9. CONCLUSIONS AND FUTURE WORK

A framework is described for convergence and validation of nonintrusive UQ methods, the relationship between deterministic V&V and stochastic UQ is studied, and an example is provided for a unit problem. Convergence procedures are developed for MC without and with metamodels, showing that in addition to the usual user-defined acceptable confidence intervals, convergence studies with systematic refinement ratio is required. A UQ validation procedure is developed using benchmark UQ results and defining the comparison error and its uncertainty to evaluate validation similar to deterministic studies. In the absence of experimental UQ, the converged simulation-based MC results may be used as the validation benchmark.

Deterministic V&V is appropriate for estimating iterative and grid numerical errors/uncertainties and deterministic modeling errors/uncertainties due, e.g., to turbulence models. Stochastic UQ methods are appropriate for estimating output uncertainties due to geometrical and operational stochastic input variables. The variations of the deterministic errors with respect to the input variables contribute to the stochastic UQ results and are considered if significant. A stochastic influence factor is defined to evaluate the effects of the stochastic input variability on the performance expectation, and four possibilities in making design decisions are identified.

The unit problem is high-fidelity RANS simulations of NACA0012 airfoil at AoA = 3° and variable Re with normal distribution, mean Re = 1.7M, and 10% variation coefficient. Deterministic V&V is conducted for the mean Re. Monotonic convergence is achieved for three grid triplets with the grids ranging from 22.2K to 5.67M grid points. The numerical uncertainties for this unit problem are relatively small compared with the previous ship hydrodynamics problems. However, the multiple grid-triplet studies show that the asymptotic range is not yet established and requires finer grids. Validation is achieved at an average validation uncertainty interval of 2.2% D, averaged between lift and drag coefficients, with an average error value of 0.25% D. The grid convergence studies are also conducted for the minimum and maximum Re and it is verified that the variations of the deterministic V&V results are negligible. For the MC-LHS results it is shown that the UEV values could be reasonably small before the solutions are convergent. The MC convergence criteria are satisfied with N = 400 deterministic CFD simulations. The distributions for CD and CL output functions are closest to the normal with an average 96% significance level of the chi-square test, implying weakly nonlinear outputs. The stochastic influence factors are smaller than the Uλ values of 1.5% S for US = USN and 2.6% S for US = UV. The condition 1a is determined, i.e., stochastic effects are negligible on the performance expectation. The UJ values are relatively small for this unit problem, but in general can be large. The UJ value is 2.88% S for CD, while it is only 0.37% S for CL, implying that CL is only weakly dependent on Re. For the MC-LHS with metamodels, at least N = 400 for all metamodels and M = 28 for the power law and the LS-SVM metamodels are required to satisfy the convergence criteria. Validation is achieved for the metamodels, the Gauss quadratures, and the PC method against the benchmark solution, i.e., the converged MC results. Future work should consider finer grids to establish the asymptotic range and obtain the numerical benchmark solution for the deterministic studies. Similarly, the asymptotic solution for the MC-LHS method should be obtained as a numerical benchmark for UQ validations.

Note that some of the methods used for the unit problem may not be feasible for the industrial applications. The sources of complexity could be large number of input variables, nonmonotonic and/or nonsmooth output functions, significant δ J* variations versus ξ, or large metamodel prediction errors which need to be considered in the UQ evaluations. For some problems the computational resources may not allow a direct MC, or obtaining converged MC results may not be possible. The construction of certain metamodels may be more expensive than a direct MC and therefore unreasonable.

Currently work is underway to extend the framework and the relationship studies for an industrial problem with both single and multivariable UQ input variables. A Delft catamaran free to sinkage and trim in calm water with speed, geometry, and speed and geometry input variables is studied using high-fidelity unsteady Reynolds averaged Navier- Stokes (URANS) and lower fidelity potential flow. The output functions are more challenging than the unit problem such that more robust metamodels are included. The two-variable study requires efficient sampling methods with minimum correlation between the variables. The preliminary results show that the current framework is applicable and the convergence and validation criteria are mostly met. The stochastic influence factor is found significant, unlike the unit problem, with combined 2a and 2b conditions. For some of the outputs the χ2 tests do not pass for the well-known distributions and therefore α is not available. For the geometric variations the outputs are noisy and convergence is not always easy to achieve. The dependence of the V&V results to the input variables are found significant such that the simulation numerical errors should be considered in the UQ studies. The next step will include added-resistance/powering in regular and irregular waves for the Delft catamaran. Also it will be considered to carry out UQ experiments for validations, if possible. Finally, the UQ results will be used as the precursory step to robust and reliability-based design optimization applications.

ACKNOWLEDGMENTS

This research is supported by the Office of Naval Research through grant no. N00014-01-1-0073 under the administration of Dr. Patrick Purtell, and grant no. N00014-08-1-0491 under the administration of Dr. Ki-Han Kim and is carried out in collaboration with AVT-191, Application of Sensitivity Analysis and Uncertainty Quantification to Military Vehicle Design. Computations were performed at the NAVY DoD Supercomputing Resource Center.

REFERENCES

1. Hirsch, C., Non-deterministic simulation for CFD-based design methodologies, Proc. of NATO RTO-AVT-147 Symp. on Computational Uncertainty in Military Vehicle Designi, Paper No. 29, Athens, Greece, Dec. 3–6, 2007.

2. Evans, T. P., Tattersall, P., and Doherty, J. J., Identification and quantification of uncertainty sources in aircraft related CFD computations–An industrial perspective, Proc. of NATO RTO-AVT-147 Symp. on Computational Uncertainty in Military Vehicle Design, Paper No. 6, Athens, Greece, Dec. 3–6, 2007.

3. Couaillier, V. and Delbos, M., Accurate boundary conditions with classical high order CFD methods, Proc. of NATO RTOAVT- 147 Symp. on Computational Uncertainty in Military Vehicle Design, Paper No. 29, Athens, Greece, Dec. 3–6, 2007.

4. Peter, J., Lazareff, M., and Couaillier, V., Verification, validation and error estimation in CFD for compressible flows, Int. J. Eng. Syst. Model. Simul., 2(1/2):75–86, 2010.

5. Najm, H., Uncertainty quantification and polynomial chaos techniques in computational fluid dynamics, Annu. Rev. Fluid Mech., 41:35–52, 2009.

6. Knio, O. and Le Maître, O., Uncertainty propagation in CFD using polynomial chaos decomposition, Fluid Dyn. Res., 38:616–640, 2006.

7. Witteveen, J. A. S., Doostan, A., Chantrasmi, T., Pecnik, R., and Iaccarino, G., Comparison of stochastic collocation methods for uncertainty quantification of the transonic RAE 2822 airfoil, Proc. of Workshop on Quantification of CFD Uncertainties, Brussels, Belgium, Oct. 29–30, 2009.

8. Peter, J. and Bompard, M., Uncertainty quantification for robust design using polynomial chaos and metamodels, presented in Workshop SOUMO, Sophia-Antipolis, France, Oct. 5–7, 2009.

9. Loeven, G. J. A., Witteveen, J. A. S., and Bijl, H., Probabilistic collocation: An efficient non-intrusive approach for arbitrarily distributed parametric uncertainties, Proc. of 45th AIAA Aerospace Sciences Meeting and Exhibit, Reno, NV, Jan. 8–11, 2007.

10. Loeven, G. J. A., Witteveen, J. A. S., and Bijl, H., A probabilistic radial basis function approach for uncertainty quantification, Proc. of NATO RTO-AVT-147 Symp. on Computational Uncertainty in Military Vehicle Design, Paper No. 29, Athens, Greece, Dec. 3–6, 2007.

11. Litvinenko, A., Matthies, H. G., Eisfeld, B., and Kroll, N., Quantification of uncertainties in numerical aerodynamics, Proc. of NODESIM-CFD Workshop on Quantification of CFD Uncertainties, Brussels, Belgium, Oct. 29–30, 2009.

12. Hosder, S., Walters, R. W., and Balch, M., Efficient sampling for non-intrusive polynomial chaos applications with multiple input uncertain variables, Proc. of 48th Structures, Structural Dynamics, & Materials Conf., Paper No. 2007-1939, Honolulu, HI, Apr. 23–26, 2007.

13. Hosder, S.,Walters, R.W., and Balch, M., Efficient uncertainty quantification applied to the aeroelastic analysis of a transonic wing, Proc. of 46th AIAA Aerospace Sciences Meeting and Exhibit, Paper No. 2008-729, Reno, NV, Jan. 7–10, 2008.

14. Witteveen, J. A. S. and Bijl, H., Effect of randomness on multi-frequency aeroelastic responses resolved by unsteady adaptive stochastic finite elements, J. Comput. Phys., 228:7025–7045, 2009.

15. Ata, M. Y., A convergence criterion for the Monte Carlo estimates, Simul. Model. Pract. Theory, 15(3):237–246, 2007.

16. Liang, B. and Mahadevan, S., Error and uncertainty quantification and sensitivity analysis in mechanics computational models, Int. J. Uncertainty Quantification, 1(2):147–161, 2011.

17. Field, Jr., R. V. and Grigoriu, M., On the accuracy of the polynomial chaos approximation, J. Probab. Eng. Mech., 19:65–80, 2004.

18. Stern, F., Wilson, R., and Shao, J., Quantitative approach to V&V of CFD simulations and certification of CFD codes, Int. J. Numer. Methods Fluids, 50:1335–1355, 2006.

19. Coleman, H. W. and Stern, F., Uncertainties and CFD code validation, J. Fluids Eng., 119(4):795–803, 1997.

20. Xing, T. and Stern, F., Factors of safety for Richardson extrapolation, J. Fluids Eng., 132(6):061403-1–061403-13, 2010.

21. Stern, F., Wilson, R. V., Coleman, H. W., and Paterson, E. G., Comprehensive approach to verification and validation of CFD simulations—Part 1: Methodology and Procedures, ASME J. Fluids Eng., 123(4):793–802, 2001.

22. Coleman, H., Stern, F., Di Mascio, A., and Campana, E., Technical note: The problem with oscillatory behavior in grid convergence studies, ASME J. Fluids Eng., 123(2):438–439, 2001.

23. Eça, L. and Hoekstra, M., An evaluation of verification procedures for CFD applications, Proc. of 24th Symp. on Naval Hydrodynamics, Fukuoka, Japan, July 8–13, 2002.

24. Logan, R. W. and Nitta, C. K., Comparing 10 methods for solution verification and linking to model validation, J. Aerosp. Comput. Inf. Commun., 3(7):354–373, 2006.

25. Helton, J. C. and Davis, F. J., Latin hypercube sampling and the propagation of uncertainty in analyses of complex systems, Reliab. Eng. Syst. Saf., 81:23–69, 2003.

26. Du, X. and Chen, W., Methodology for managing the effect of uncertainty in simulation-based design, AIAA J., 38(8):1471–1478, 2000.

27. Le MaÎtre, O. P. and Knio, O. M., Spectral Methods for Uncertainty Quantification: With Applications to Computational Fluid Dynamics, Springer, New York, 2010.

28. Ochi, M. K., Applied Probability and Stochastic Processes in Engineering and Physical Sciences, Wiley, New York, 1990.

29. Bendat, J. S. and Piersol, A. G., Measurement and Analysis of Random Data, John Wiley and Sons, Inc., New York, 1966.

30. Suykens, J. A. K., Gestel, T. V., Brabanter, J. D., Moor, B. D., and Vandewalle, J., Least Squares Support Vector Machines, World Scientific, Singapore, 2002.

31. Hosder, S., Walters, R. W., and Perez, R., A non-intrusive polynomial chaos method for uncertainty propagation in CFD simulations, Proc. of 44th AIAA Aerospace Sciences Meeting and Exhibit, Paper No. 2006-891, Reno, NV, Jan. 9–12, 2006.

32. Kreyszig, E., Advanced Engineering Mathematics, Wiley, New York, 1993.

33. Jin, R., Chen, W., and Simpson, T. W., Comparative studies of metamodelling techniques under multiple modeling criteria, Struct. Multidiscip. Optim., 23:1–13, 2001.

34. Gregory, N. and Wilby, P. G., NPL 9615 and NACA 0012–A comparison of aerodynamic data, Tech. Report ARC/CP–1261, National Physical Laboratory, London, UK, 1973.

35. Sheldahl, R. E. and Klimas, P. C., Aerodynamic characteristics of seven symmetrical airfoil sections through 180-degree angle of attack for use in aerodynamic analysis of vertical axis wind turbines, Tech. Report SAND80-2114, Sandia National Laboratories, Albuquerque, NM, 1981.

36. Carrica, P. M., Huang, J., Noack, R., Kaushik, D., Smith, B., and Stern, F., Large-scale DES computations of the forward speed diffraction and pitch and heave problems for a surface combatant, Comput. Fluids, 39(7):1095–1111, 2010.

37. Simonsen, C. and Stern, F., CFD simulation of KCS sailing in regular head waves, Proc. of Gothenburg 2010: A Workshop on CFD in Ship Hydrodynamics, Chalmers University of Technology, Gothenburg, Sweden, Dec. 8–10, 2010.

Portal Digital Begell Biblioteca digital da Begell eBooks Diários Referências e Anais Coleções de pesquisa Políticas de preços e assinaturas Begell House Contato Language English 中文 Русский Português German French Spain