# Uncertainty Quantification and Sensitivity Analysis of Transonic Aerodynamics with Geometric Uncertainty

Authors: Xiaojing Wu; Weiwei Zhang; Shufang Song
Journal: International Journal of Aerospace Engineering
 (2017)

---

## Abstract

Airfoil geometric uncertainty can generate aerodynamic characteristics fluctuations. Uncertainty quantification is applied to compute its impact on the aerodynamic characteristics. In addition, the contribution of each uncertainty variable to aerodynamic characteristics should be computed by the uncertainty sensitivity analysis. In the paper, Sobol’s analysis is used for uncertainty sensitivity analysis and a nonintrusive polynomial chaos method is used for uncertainty quantification and Sobol’s analysis. It is difficult to describe geometric uncertainty because it needs a lot of input parameters. In order to alleviate the contradiction between the variable dimension and computational cost, a principal component analysis is introduced to describe geometric uncertainty of airfoil. Through this technique, the number of input uncertainty variables can be reduced and typical global deformation modes can be obtained. By uncertainty quantification, we can learn that the flow characteristics of shock wave and boundary layer separation are sensitive to the geometric uncertainty in transonic region, which is the main reason that transonic drag is sensitive to the geometric uncertainty. The sensitivity analysis shows that the model can be simplified by eliminating unimportant geometric modes. Moreover, which are the most important geometric modes to transonic aerodynamics can be learnt. This is very helpful for airfoil design.


---

## Body

## 1. Introduction

There inherently exist a vast of uncertainties in the practical aircraft design and application, thus causing fluctuations of aircraft performance. Therefore, it is important to take these uncertainties into account at the beginning of aircraft design. With the development of computer technology, computational fluid dynamics (CFD) technology has been widely used to solve problems in aerodynamic mechanics. The traditional CFD simulation is deterministic. However, a variety of uncertainties are inevitable in CFD simulation with the increasing complexity of the fluid problem in reality, leading to the mismatch between the CFD simulation results and the actual results [1, 2].The sources and classifications of uncertainty in CFD were described in [3, 4]. Several uncertainty quantification (UQ) methods have been used in CFD, including Monte Carlo simulation (MCS), sensitivity analysis (SA), moment methods, and polynomial chaos expansions (PCE) in [5]. MCS is a statistical method, which needs a large number of samples to accurately analyze the uncertainty quantification. SA and moment methods are suitable to solve the problem of small parameter uncertainty or linear model. Recently, PCE has been widely applied to UQ of fluid problems. PCE methods can be divided into intrusive and nonintrusive methods according to the coupling ways with CFD solvers. In general, an intrusive approach calculates the unknown polynomial coefficients by projecting the resulting equations into basis functions for different modes, and it requires modifying the CFD codes, which may be difficult and time-consuming for complex problems such as N-S simulation. Mathelin et al. [6] used the intrusive PCE to research the quasi one-dimensional duct flow uncertainty propagation problems. Xiu and Karniadakis [7] extended the original polynomial methods to the general polynomial chaos (gPC) and applied it to the incompressible channel flow and the flow around a cylinder. To overcome the shortcomings of intrusive polynomial chaos, nonintrusive polynomial chaos (NIPC) has been developed. CFD is regarded as a black box without changing the CFD program code in nonintrusive methods. There are two different methods to build a PCE model: projection method and regression method. The projection method is based on the orthogonality of PCE to evaluate the coefficients by a multidimensional integral. There are two approaches to evaluate this integral: random sampling method and Gaussian quadrature-based method. Random sampling method uses MCS to evaluate the unknown coefficients, but its convergence rate is low. Deterministic sampling method uses the quadrature to evaluate the unknown coefficients. The quadrature-based methods are more efficient than random sampling methods in low-dimensional problems. However, they are inefficient for relatively high-dimensional problems because of the exponential rising in quadrature points with the increasing dimensions. The other nonintrusive collocation method evaluates PCE using regression. Hosder et al. [8] developed this regression method to NIPC for stochastic computational fluid dynamics. UQ based on PCE and its applications in fluid mechanics were comprehensively reviewed in [9, 10].In addition to UQ, global sensitivity analysis (GSA) plays an important role in quantifying the relative importance of each uncertainty source. By means of this technique, uncertainties can be systemically studied to measure their effects on the system outputs so as to screen out the uncertainties with negligible contributions to simplify the problem. A comprehensive review of uncertainty sensitivity analysis can be found in [11]. The variance-based method, also called Sobol’s analysis [12, 13], is one of the most popular practices in many disciplines. Sobol’s analysis can measure the relative importance of one input variable to the model’s output by computing the partial variance. There are many methods to compute Sobol’s indices, such as Fourier Amplitude Sensitivity Test (FAST) and sampling-based method [14, 15]. Sudret [16] introduced PCE to build surrogate models that allow one to compute Sobol’s indices analytically as postprocessing of PCE coefficients.In stochastic aerodynamic analysis, most research on UQ studies the flight condition parameter uncertainty based on NIPC methods. Loeven et al. [17] conducted a subsonic aerodynamic analysis around a NACA0012 airfoil with an uncertain free stream velocity using a commercial flow solver. From their study, an uncertain free stream velocity leads to the highest variation in pressure on the upper surface near the leading edge. Simon et al. [18] focused on the transonic stochastic response of two-dimensional airfoil to parameter uncertainty (Ma and α) using gPC. Two kinds of nonlinearities are critical to transonic aerodynamics in their study: the leeward shock movement characteristics and boundary layer separation on the aft part of airfoil downstream the shock. Chassaing and Lucor [19] conducted a stochastic investigation of flows about NACA0012 airfoil at transonic speeds. Liu et al. [20] conducted a stochastic fluid analysis on a 3D wind blade case considering the wind speed as an uncertain parameter. From their study, when the flow separation appears, the separation vortex region corresponds to the maximum variation area, and the maximum variation extends to the trailing edge even to the whole suction side. Moreover, Sobol’s analysis was also introduced to stochastic aerodynamic analysis. Simon et al. [18] applied Sobol’s analysis to the stochastic transonic aerodynamic analysis considering the uncertainty of Mach number and angle of attack. Hosder and Bettis [21] conducted an uncertainty and sensitivity analysis for reentry flows with inherent and model-form uncertainties. However, rare research can be found on stochastic aerodynamics analysis considering geometric uncertainty. The geometric uncertainty on aerodynamic surfaces resulting from manufacturing errors has significant effect on the aerodynamic performance. Due to the high cost of precise surface manufacturing techniques, it is often impractical to remove the impact of these geometric variations by improving the manufacturing tolerance. In other words, the geometric uncertainty resulting from manufacturing errors is unavoidable. Therefore, it is important to conduct an UQ for aerodynamics evaluation considering this geometric uncertainty. Moreover, a GSA for aerodynamics evaluation considering geometric uncertainty is also necessary, because it can obtain important geometric variation parameter to aerodynamics by sensitivity analysis. It is difficult to describe geometric variation because it needs a lot of input parameters. The high dimensions of input uncertainty can lead to the high computational cost to build PCE model, which makes it more difficult to conduct a stochastic aerodynamic analysis considering geometric uncertainty.In this paper, we conduct UQ and GSA of transonic aerodynamics of airfoil RAE2822 considering geometric uncertainty. The principle component analysis (PCA) technique combining with an airfoil parameterization method is used to describe the airfoil shape deformation and obtain the typical deformation modes. Then, a stochastic PCE model is built for UQ and GSA. The paper is structured as follows. Section2 introduces UQ and the variance-based GSA (Sobol’s analysis). In Section 3, two methods, MCS and PCE, are introduced to UQ and GSA in detail. In Section 4, UQ and GSA of aerodynamics are conducted. Section 5 outlines several useful conclusions.


## 2. Uncertainty Quantification and Global Sensitivity Analysis

The aim of UQ is to compute the influence of input parameters on the model’s output quantitatively. The input uncertainty may be expressed in a number of ways, for example, interval bounds or probability density functions. In the paper, the uncertainty is quantified by the probabilistic analysis and the robustness (mean and variance) of outputs is focused on. The mean and variance are shown by the following equation:(1)EY=Efx=∫-∞+∞fxpxxdxVY=Efx-μx2=∫-∞+∞fx-μf2pxxdx.Although the influence of the input uncertainty on output variables can be identified by UQ, it cannot obtain the contribution of each uncertainty. In this case, A GSA is applied to study how the output of a model can be apportioned to different sources of input uncertainty. Since most of the sensitivity analyses are based on the partial derivate, this approach is sometimes called local sensitivity. However, the fatal limitation of a derivative-based approach is that when the model is uncertain is unwarranted. In other words, derivatives are only informative at the base point where they are computed and do not provide an exploration of the rest of the uncertainty space. The alternative uncertainty sensitivity analysis is GSA. The main goal of GSA is to demonstrate the relative importance of each input uncertainty among the overall uncertainty in the output. It is often useful to rank the contribution of each uncertainty to the overall uncertainty in the output. A detailed review of a variety of GSA methods can be found in [15]. In the current study, a variance-based global sensitivity method is used to accomplish this task. The method is based on Sobol’s variance decomposition, by which f-function can be uniquely decomposed into 2n functional terms of increasing dimensions:(2)Y=f0+∑ifi+∑i∑j>ifij+⋯+f12⋯k,where f0=E(Y), fi=V(E(Y∣Xi))-g0, and fij=V(E(Y∣Xi,Xj))-fi-fj-f0. Adding variance to both sides of (2), we can obtain the following equation:(3)VY=∑iVi+∑i∑j>iVij+⋯+V12⋯k,where Vi=V[E(Y∣Xi)],Vij=V[E(Y∣Xi,Xj)]-Vi-Vj,…. The first-order partial variance Vi can be regarded as the average reduction of model output variance resulting from fixed Xi; that is, Vi measures the individual contribution of Xi to the total variance V(Y). The larger Vi is, the more reduction of output variance can be obtained by reducing the uncertainty of Xi. The second-order partial variance Vij quantifies the interaction effect between Xi and Xj. Similar interpretations can be given to the higher-order partial variances. Dividing both sides of (3) by VY, we obtain (4)Si1,…,is=Vi1,…,isV∑iSi+∑i∑jSij+∑i∑j>i∑l>jSijl+⋯+S123⋯k=1.S i 1 , … , i s are Sobol’s indices, which measure the relative contribution to the total variance. Si are the main sensitivity indices. And Sij,Sijl,…,S123⋯k measure the interaction effects.(5)Si=VEY∣XiVY.Another commonly used measurement is the total partial varianceVTi, which is defined as the summation of all terms in (6) with subscripts including i. In other words, VTi incorporates both the individual effect of Xi and its interaction effects with all the other input variables. (6)VTi=∑℘iVi1,…,is,℘i=i1,…,is:∃k,1≤k≤s,ik=i(7)STi=VTiVY=VY-VEY∣X~iVY=EVY∣X~iVY.UsuallySi is used to select important variables, while STi is more suitable to screen noninfluential variables.


## 3. The Nonintrusive Polynomial Chaos

There are many uncertainty analysis methods available for UQ [22, 23]. In the paper, The NIPC method is used for the current study. The PCE is a stochastic method based on the spectral representation of uncertainty. According to the spectral representation, the random function can be decomposed into deterministic and stochastic components. For example, a random variable (X) can be represented by the following equation:(8)Xx,ξ=∑j=0∞αjxΨjξ≈∑j=0PαjxΨjξ,where αj(x) is the deterministic component and Ψj(ξ) is the random basis function corresponding to the jth mode. From (8), the random variable X is the function of deterministic independent variable vector x and the n-dimensional standard random variable vector ξ=(ξ1,ξ2,…,ξn). The polynomial chaos expansion given by (8) contains an infinite number of terms. In a practical computational context, the terms of PCE can be truncated by both order p and dimension n. The number of terms is finite, which is given by the following equation:(9)Nt=P+1=n+p!n!p!,where n is the number of random dimensions and p is the order of polynomial chaos. When the input uncertainty obeys Gauss distribution, the basis function is the multidimensional Hermite polynomial.There are two ways to solve the coefficients of (8): intrusive and nonintrusive methods. The intrusive method computes the unknown polynomial coefficients by projecting the resulting equations into basis functions for different modes. It requires the modification of the deterministic code, which may be difficult, expensive, and time-consuming for complex computational problems such as complex N-S equations. Alternatively, the nonintrusive method treats the CFD as a black box without changing the program code when propagating uncertainty. Now, we pay close attention to how to solve the coefficients in nonintrusive methods. There are two methods: projection method and regression method.The main idea of projection method is to use Galerkin projection to solve coefficients. Equation (8) can be transformed to (10) by inner product: (10)Xx,ξ,Ψkξ=∑j=0PαjxΨjξ,Ψkξ.· represents the inner product, which can be expressed by the following formula:(11)fξgξ=∫Rfξgξwξdξ.For Gaussian random variable, the basis functions are Hermite orthogonal polynomials. Because of orthogonality, (10) can be transformed to the following equation:(12)Xx,ξ,Ψkξ=αkxΨk2ξ.And then the following can be derived:(13)αkx=Xx,ξ,ΨkξΨk2ξ=1Ψk2ξ∫Rα∗x,ξΨkξωξdξ.The sampling approaches can be divided into random sampling method and deterministic sampling method. Random sampling method uses MCS, LHS, and other methods to compute the projection integrals. However, its convergence rate is low. Deterministic sampling method uses the quadrature for the numerical evaluation of the unknown coefficients. Usingn-dimensional Gauss-Hermite quadrature, with q points in each dimension, we can compute the unknown coefficients by the following equation:(14)αk=1Ψk2∑i1=1q⋯∑in=1quxi1,…,xinΨkxi1,…,xin∏k=1nωik.In (14), one-dimensional integral is expanded into a high-dimensional form by tensor product, the calculation times of which require p+1n points for pth-order polynomial chaos. For low-dimensional problems, the sampling efficiency of deterministic sampling method has been greatly improved compared with random sampling method. However, the calculation times grow exponentially with the increasing dimensions.The regression method is another approach to compute the coefficients using a linear equation system based on a selected set of points. A linear system of equations can be obtained by(15)α=ΨTΨ-1·ΨT·α∗,where Ψ is given by Ψij=Ψj(ξi),i=1,…,N;j=0,…,P-1. The accuracy of regression method depends on the selection of sample points in design space. According to [16], the optimal design is given by the roots of Hermite polynomial, in which the optimal sample points are p+1n. The sample point of NIPC with projection method is p+1n, which grows exponentially with the increasing number of input dimensions. Compared with projection based method, the sample points can be reduced by regression based method. The sample points are 2(P+1), which gives a better approximation at each polynomial degree in [8]. To study how the number of samples influences the analysis results, the oversampling ratio np is defined as(16)np=number  of  samplesP+1.After computing the unknown coefficients, the approximate PCE model is built, and then the model is used to compute the first two statistical moments analytically as follows:(17)EXξ=∫RXξWξdξ=α0σPC2=Var∑j=0P-1αj2Ψjξ=∑j=1P-1αj2EΨj2ξ.In order to compute Sobol’s indices, Sobol’s decomposition based on PCE should be derived. The detailed derivation process is in [16]. The important equations are shown:(18)Ψjx≡Ψj:Ψjx=∏i=1nHαixi,where Hk(x) is the kth Hermite polynomial. Let us define the set α, where the indices {i1,…,is} are nonzero:(19)℘i1,…,is=α:αk>0∀k=1,…,n,k∈i1,…,isαk=0∀k=1,…,n,k∉i1,…,isfPCx=f0+∑i=1n∑α∈℘ifαΨαxi+∑1≤i1<i2≤n∑α∈℘i,jfαΨαxi1,xi2+⋯+∑1≤i1<⋯<is≤n∑α∈℘i,jfαΨαxi1,xi2+⋯+∑α∈℘1,2,…,nfαΨαx1,…,xnVi1,…,is=∑α∈℘i1,…,isfα2Ψα2.Sobol’s analysis of PCE is shown by (19). Combining (4)–(7), Sobol’s indices are obtained.


## 4. Uncertainty Quantification and Sensitivity Analysis for Transonic Aerodynamics

### 4.1. The Description of Geometric Uncertainties

In engineering, the geometric uncertainty on aerodynamic surfaces resulting from manufacturing errors or wearing is unavoidable, which has significant influence on the aerodynamic performance. Therefore, UQ and GSA considering geometric uncertainty to aerodynamics are conducted in this section. To carry out this work, the first step is to describe the geometric uncertainty in the computing environment.In [24, 25], the main geometric variation modes are obtained by principle component analysis (PCA) with a large amount of geometric statistical manufacturing error data of airfoil. The main variation modes based on the statistical manufacturing error data are used to achieve the geometric variation. The description of the geometric uncertainty of airfoil is shown in the following equation: (20)g=gn+g-+∑i=1nsσizivi,where gn is the nominal geometry; g- is the average geometric variation; vi is the geometric mode shape; ns is the number of mode shapes used to represent the variation in geometry. The geometric mode can be computed by PCA based on manufacturing samples. σi is the ith singular value of the measurement snapshot matrix, which represents the geometric variability attributable to the ith mode. zi is a random parameter which obeys the standard normal distribution; thus, the product σizi is the stochastic contribution of ith mode.It is difficult to describe this geometric variation in the computing environment. In [26], a Gaussian random process simulation is used to obtain the geometric data. Then PCA is used to obtain main geometric modes. In addition, a variety of parametric methods have been employed to describe the geometric variation in aerodynamic design so far, likePARSEC-11 geometry parameterization, class-shape function transformation (CST) method, free-form deformation (FFD) method, Hicks-Henne bump functions, and so forth [27–31]. By changing the parameter of these parametric methods, the geometric variation is realized. Generally, these parametrization methods need a lot of parameters to represent the airfoil shape. To reduce the dimensions of the variables, the PCA technology combined with airfoil parameterization is used to describe the geometric variation used in the paper. Firstly, a parametric method is used to generate a set of sample data. Then the PCA is used to obtain the main deformation modes based on the generated sample data. In this paper, we use a parameterized representation of the airfoil RAE2822 by CST method with 24 parameters. The data of the measurement points on the airfoil surface are obtained by random perturbation of CST parameters. In this way, the PCA based sample data is conducted. Figure 1 shows the eigenvalue with the number of modes, which represents the geometric variation attributable to each mode. The eigenvalue of one mode is smaller; the proportion of the mode to geometric variation is smaller. It can be observed that the proportion of the first 12 modes is 99.32%, which means that they can be used to express the geometric variation well. Through the PCA, the number of input parameters is reduced from 24 to 12. The first 12 modes obtained by PCA are shown in Figure 2. The 12 modes include 6 modes of the upper surface (mode 2, mode 3, mode 5, mode 7, mode 9, and mode 10) and 6 modes of the lower surface (mode 1, mode 4, mode 6, mode 8, mode 11, and mode 12). It can be seen that the lower-order modes in these modes are global geometric deformation modes and present some typical geometric deformation. Specifically, mode 1 and mode 2 are the scale modes in the thickness direction; mode 3 and mode 4 are translation modes of the maximum thickness in the axial direction; and mode 5 and mode 6 are the extrusion modes of the upper surface. So, the number of parameters used to describe the geometric variation is reduced to a large extent and the typical deformation modes are obtained through the PCA technique.Figure 1
The eigenvalue with the number of modes.Figure 2
The first 12 deformation modes of RAE2822 airfoil.


### 4.2. Deterministic Flow Computations (CFD Simulations)

The deterministic steady-state flow solutions are computed by the Reynolds-Averaged Navier-Stokes (RANS) equations. The two-dimensional Navier-Stokes equations can be written in the integral form as(21)∂∂t∬ΩQdΩ+∫∂ΩFQndΓ=∫∂ΩGQndΓ,where Ω is the control volume; ∂Ω is the boundary of control volume; n=(nx,ny)T is the outer normal vector of the control volume boundary; Q is the vector of conservative variables; F(Q) is the inviscid flux; and G(Q) is the viscous flux.For the cell-centered finite volume method, the computational domain is divided into nonoverlapping control volumes that completely cover the domain. The interface variables are derived from the average values of the grid cells in order to calculate the flux of control volumes. Through spatial discretization, the equations of the integral form are transformed to the linear ordinary differential equations, and the flow variables can be obtained by the time marching method. We adopt AUSM scheme to evaluate the numerical flux which can also be evaluated by the upwind scheme [32]. The geometry of airfoil is changing in the process of UQ, and radial basis functions (RBF) are used for mesh deformation [33]. Figure 3 shows the computing grids of RAE2822 and the reliability verification of the CFD program.Figure 3
Computing grids of RAE2822 and the reliability verification of the CFD program.
(a)
Global grid (b)
Local grid (c)
The comparison ofCp


### 4.3. Uncertainty Quantification of Aerodynamics

UQ of aerodynamics is conducted in this section. The uncertain variables arezi(i=1,2,…,12) in (20), which obey Gaussian distribution and are independent from each other. The flow state is selected in a transonic region (Ma = 0.73, α=2.5°, and Re=3.0×106); NIPC and MCS are used to quantify uncertainty. MCS is used to aim at validating the accuracy and efficiency of NIPC. The computational cost of MCS is selected as 5000 samples. Figure 4 shows the range of geometric uncertainty in the paper. It can be seen that the range of fluctuation reaches 0.006 near the maximum thickness, which is 10% of the maximum thickness. In this paper, the fluctuation range of the geometry should be greater than the fluctuation caused by the manufacturing error.Figure 4
The range of geometric uncertainty.


#### 4.3.1. The Difficulty of NIPC for the High-Dimensional Problem

For a high-dimensional problem, the number of PCE terms is very high. Now, we study the number of terms withp increasing for the 12-dimensional problem. From Table 1, it can be seen that the number of terms increased rapidly with p increasing for the 12-dimensional problem.Table 1
The number of PCE terms withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 p = 6 N t 13 91 455 1820 6188 18564It is mentioned in (14) that the projection based NIPC method needs p+1n samples for UQ, which grows exponentially with the increasing number of input dimensions. Table 2 shows the number of samples of projection based NIPC method with p for the 12-dimensional problem. It can be seen that the computational cost is too high and unacceptable.Table 2
The number of samples of the projection based NIPC method withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 N s 4096 531441 16777216 244140625 2176782336The regression based NIPC method adopts the linear regression method to calculate the coefficients of PCE. To study how the number of samples influences the analysis results, the oversampling rationp is defined as (16).Compared with the projection based NIPC method, the computational cost of the regression based NIPC method is much lower. Therefore, the regression based NIPC method is used for the stochastic aerodynamic analysis considering geometric uncertainty.


#### 4.3.2. The Verification of NIPC for the UQ considering Geometric Uncertainty

The verification of the regression based NIPC method for the UQ considering geometric uncertainty is conducted.Firstly, we study the convergence of NIPC by increasing the polynomial chaos orderp. Figures 5–7 show the standard deviation of aerodynamics with the order of NIPC and the comparison with that of MCS. It can be seen that when p=2 or 3, the results of NIPC can coincide with those of MCS. However, there is a little difference between NIPC and MCS. It can also be observed that the results of NIPC are wrong when p=4. We think that the number of terms of PCE is 1820 for p=4 and n=12, which makes it difficult to compute the coefficients correctly by linear regression method. This leads to the wrong results for high-order accuracy.Figure 5
Standard deviation distribution ofCp with the order of NIPC.Figure 6
Standard deviation distribution ofCf with the order of NIPC.Figure 7
Standard deviation distribution of the aerodynamic coefficients with the order of NIPC.Moreover, the convergence of NIPC by increasing the oversampling rationp is also studied. Figures 8–10 show the standard deviation of aerodynamics with np of NIPC and the comparison with that of MCS. When np=4, the results of NIPC can coincide well with those of MCS. Through comparison, it can be observed that the results of NIPC with p=2 and np=4 are even better than those of NIPC with p=3 and np=2.Figure 8
Standard deviation distribution ofCp with oversampling ratio np.Figure 9
Standard deviation distribution ofCf with oversampling ratio np.Figure 10
Standard deviation distribution of the aerodynamic coefficients oversampling rationp.Lastly, the computational costs are compared with MCS. Table3 shows the number of samples to UQ by NIPC and MCS. It can be seen that the number of samples is 364 with p=2 and np=4, which is much lower than that of MCS.Table 3
The number of samples to UQ by NIPC and MCS.
NIPC MCS p = 2, np=2 p = 3, np=2 p = 4, np=2 p = 2, np=4 p = 2, np=6 — N 182 910 3640 364 546 5000


#### 4.3.3. The Results of UQ for Transonic Aerodynamics

Figure11 shows the error bars of Cp and Cf distribution by UQ. It can be observed that the aerodynamic loads are sensitive to the geometric variation of the upper surface. In transonic region, the location of shock fluctuates with the change of the geometric shapes, which can result in intensive fluctuations of the flow variables. We call the region in which the standard variations become much larger suddenly shock disturbance region. It can also be seen that the skin-friction behavior displays discrepancies from the Cp distribution. The variation of Cf is similar to that of Cp before and in the shock disturbance region. However, the difference occurs in the downstream region of the shock disturbance region; the standard deviation keeps larger magnitudes below the shock disturbance region and then gradually decreases. This indicates that not only the shock wave but also the boundary layer separation is sensitive to geometric uncertainty in the transonic region.Figure 11
The error bars ofCp and Cf distribution by UQ.


### 4.4. Global Sensitivity Analysis of Aerodynamics

The important deformation modes can be found out through GSA, which can provide a useful reference for subsequent design work and manufacturing process. In the paper, the NIPC method (p=2; np=4) is used for Sobol’s analysis. The aerodynamic load distribution and aerodynamic coefficients are included for Sobol’s analysis. Figures 12 and 13 show the results of Sobol’s analysis of the pressure and friction coefficient distribution. The partial standard deviation, total standard deviation, and standard deviation caused by the coupling effect for each mode are given in Figures 12 and 13. It can be seen that the deformation modes of the upper surface are much more important than those of lower surface for transonic aerodynamic loads (Cp and Cf). Specifically, the first 3 modes are important for aerodynamic load distribution. For Cp distribution, the main variation of modes 2, 3, and 5 is in shock disturbance region. For Cf distribution, the main variation of modes 2 and 5 is in both shock disturbance region and boundary layer separation region, while the main variation of mode 3 is only in shock disturbance region. In addition, it can also be observed that the coupling effect among geometric deformation modes is weak in Figures 12 and 13. In order to visually observe the effects of modes, the pressure coefficient and friction coefficient of error bars are given in Figures 14 and 15. From Section 4.1, we can observe that mode 1 and mode 2 are the scale modes in the thickness direction; mode 3 and mode 4 are translation modes of the maximum thickness in the axial direction; and mode 5 is the extrusion mode of the upper surface. We can draw the following conclusions: the scale mode and extrusion mode are important for the flow characteristic of shock wave and boundary layer separation, but the translation mode is only important for flow characteristic of shock wave.Figure 12
The results ofCp distribution by Sobol’s analysis (σi is the partial standard variation with uncertain variable i; σTi is the total partial standard variation with uncertain variable i; σcouplingi is the nonlinear coupling effect among uncertain variable i and other uncertainty variables).Figure 13
The results ofCf distribution by Sobol’s analysis (σi is the partial standard variation with uncertain variable i; σTi is the total partial standard variation with uncertain variable i; σcouplingi is the nonlinear coupling effect among uncertain variable i and other uncertainty variables).Figure 14
The error bars ofCp distribution by Sobol’s analysis.Figure 15
The error bars ofCf distribution by Sobol’s analysis.The results of GSA for aerodynamic coefficients including lift coefficient (Cl), drag coefficient (Cd), moment coefficient (Cm), and lift-to-drag ratio (L/D) are shown in Figure 16. For Cl, mode 1 and mode 2 are important, which indicate that the scale modes are important for Cl. For Cd and L/D, modes 2, 3, and 5 are important; the first 3 deformation modes of the upper surface have bigger contribution. For Cm, the scale mode of upper surface is most important. In addition, it can also be seen that the coupling effect among deformation modes is very weak.Figure 16
The results of Sobol’s analysis of aerodynamics ((a) Sobol’s indices forCl; (b) Sobol’s indices for Cd; (c) Sobol’s indices for Cm; (d) Sobol’s indices for L/D).
(a)
Sobol’s indices forCl (b)
Sobol’s indices forCd (c)
Sobol’s indices forCm (d)
Sobol’s indices forL/DFrom the results of GSA, it can be observed that modes 1, 2, 3, and 5 are important for aerodynamics. Therefore, the UQ of aerodynamics is conducted considering the 4 modes. Moreover, because the coupling effects among the modes are small from GSA, we adopt regression method with parametersp=2 and np=4 to build the PCE model without the interaction terms for UQ. The comparison of the UQ results is shown in Table 4. It can be observed that the UQ results are relatively small different from each other; and the number of sample points is reduced from 364 to 32. That is to say, the unimportant modes are eliminated and PCE model can be modified by GSA, which can reduce the computational cost in UQ and robust design.Table 4
The comparison of the UQ results for aerodynamic coefficients between the selected 4 modes and the full 12 modes.
n = 4 (NIPC) n = 12 (NIPC) C l Mean 0.7234 0.7232 StD 0.0096 0.0097 C d Mean 0.0163 0.0163 StD 0.0012 0.0013 C m Mean 0.2718 0.2720 StD 0.0053 0.0057 L / D Mean 44.679 44.627 StD 3.1699 3.3608 N  64 364


## 4.1. The Description of Geometric Uncertainties

In engineering, the geometric uncertainty on aerodynamic surfaces resulting from manufacturing errors or wearing is unavoidable, which has significant influence on the aerodynamic performance. Therefore, UQ and GSA considering geometric uncertainty to aerodynamics are conducted in this section. To carry out this work, the first step is to describe the geometric uncertainty in the computing environment.In [24, 25], the main geometric variation modes are obtained by principle component analysis (PCA) with a large amount of geometric statistical manufacturing error data of airfoil. The main variation modes based on the statistical manufacturing error data are used to achieve the geometric variation. The description of the geometric uncertainty of airfoil is shown in the following equation: (20)g=gn+g-+∑i=1nsσizivi,where gn is the nominal geometry; g- is the average geometric variation; vi is the geometric mode shape; ns is the number of mode shapes used to represent the variation in geometry. The geometric mode can be computed by PCA based on manufacturing samples. σi is the ith singular value of the measurement snapshot matrix, which represents the geometric variability attributable to the ith mode. zi is a random parameter which obeys the standard normal distribution; thus, the product σizi is the stochastic contribution of ith mode.It is difficult to describe this geometric variation in the computing environment. In [26], a Gaussian random process simulation is used to obtain the geometric data. Then PCA is used to obtain main geometric modes. In addition, a variety of parametric methods have been employed to describe the geometric variation in aerodynamic design so far, likePARSEC-11 geometry parameterization, class-shape function transformation (CST) method, free-form deformation (FFD) method, Hicks-Henne bump functions, and so forth [27–31]. By changing the parameter of these parametric methods, the geometric variation is realized. Generally, these parametrization methods need a lot of parameters to represent the airfoil shape. To reduce the dimensions of the variables, the PCA technology combined with airfoil parameterization is used to describe the geometric variation used in the paper. Firstly, a parametric method is used to generate a set of sample data. Then the PCA is used to obtain the main deformation modes based on the generated sample data. In this paper, we use a parameterized representation of the airfoil RAE2822 by CST method with 24 parameters. The data of the measurement points on the airfoil surface are obtained by random perturbation of CST parameters. In this way, the PCA based sample data is conducted. Figure 1 shows the eigenvalue with the number of modes, which represents the geometric variation attributable to each mode. The eigenvalue of one mode is smaller; the proportion of the mode to geometric variation is smaller. It can be observed that the proportion of the first 12 modes is 99.32%, which means that they can be used to express the geometric variation well. Through the PCA, the number of input parameters is reduced from 24 to 12. The first 12 modes obtained by PCA are shown in Figure 2. The 12 modes include 6 modes of the upper surface (mode 2, mode 3, mode 5, mode 7, mode 9, and mode 10) and 6 modes of the lower surface (mode 1, mode 4, mode 6, mode 8, mode 11, and mode 12). It can be seen that the lower-order modes in these modes are global geometric deformation modes and present some typical geometric deformation. Specifically, mode 1 and mode 2 are the scale modes in the thickness direction; mode 3 and mode 4 are translation modes of the maximum thickness in the axial direction; and mode 5 and mode 6 are the extrusion modes of the upper surface. So, the number of parameters used to describe the geometric variation is reduced to a large extent and the typical deformation modes are obtained through the PCA technique.Figure 1
The eigenvalue with the number of modes.Figure 2
The first 12 deformation modes of RAE2822 airfoil.


## 4.2. Deterministic Flow Computations (CFD Simulations)

The deterministic steady-state flow solutions are computed by the Reynolds-Averaged Navier-Stokes (RANS) equations. The two-dimensional Navier-Stokes equations can be written in the integral form as(21)∂∂t∬ΩQdΩ+∫∂ΩFQndΓ=∫∂ΩGQndΓ,where Ω is the control volume; ∂Ω is the boundary of control volume; n=(nx,ny)T is the outer normal vector of the control volume boundary; Q is the vector of conservative variables; F(Q) is the inviscid flux; and G(Q) is the viscous flux.For the cell-centered finite volume method, the computational domain is divided into nonoverlapping control volumes that completely cover the domain. The interface variables are derived from the average values of the grid cells in order to calculate the flux of control volumes. Through spatial discretization, the equations of the integral form are transformed to the linear ordinary differential equations, and the flow variables can be obtained by the time marching method. We adopt AUSM scheme to evaluate the numerical flux which can also be evaluated by the upwind scheme [32]. The geometry of airfoil is changing in the process of UQ, and radial basis functions (RBF) are used for mesh deformation [33]. Figure 3 shows the computing grids of RAE2822 and the reliability verification of the CFD program.Figure 3
Computing grids of RAE2822 and the reliability verification of the CFD program.
(a)
Global grid (b)
Local grid (c)
The comparison ofCp


## 4.3. Uncertainty Quantification of Aerodynamics

UQ of aerodynamics is conducted in this section. The uncertain variables arezi(i=1,2,…,12) in (20), which obey Gaussian distribution and are independent from each other. The flow state is selected in a transonic region (Ma = 0.73, α=2.5°, and Re=3.0×106); NIPC and MCS are used to quantify uncertainty. MCS is used to aim at validating the accuracy and efficiency of NIPC. The computational cost of MCS is selected as 5000 samples. Figure 4 shows the range of geometric uncertainty in the paper. It can be seen that the range of fluctuation reaches 0.006 near the maximum thickness, which is 10% of the maximum thickness. In this paper, the fluctuation range of the geometry should be greater than the fluctuation caused by the manufacturing error.Figure 4
The range of geometric uncertainty.


### 4.3.1. The Difficulty of NIPC for the High-Dimensional Problem

For a high-dimensional problem, the number of PCE terms is very high. Now, we study the number of terms withp increasing for the 12-dimensional problem. From Table 1, it can be seen that the number of terms increased rapidly with p increasing for the 12-dimensional problem.Table 1
The number of PCE terms withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 p = 6 N t 13 91 455 1820 6188 18564It is mentioned in (14) that the projection based NIPC method needs p+1n samples for UQ, which grows exponentially with the increasing number of input dimensions. Table 2 shows the number of samples of projection based NIPC method with p for the 12-dimensional problem. It can be seen that the computational cost is too high and unacceptable.Table 2
The number of samples of the projection based NIPC method withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 N s 4096 531441 16777216 244140625 2176782336The regression based NIPC method adopts the linear regression method to calculate the coefficients of PCE. To study how the number of samples influences the analysis results, the oversampling rationp is defined as (16).Compared with the projection based NIPC method, the computational cost of the regression based NIPC method is much lower. Therefore, the regression based NIPC method is used for the stochastic aerodynamic analysis considering geometric uncertainty.


### 4.3.2. The Verification of NIPC for the UQ considering Geometric Uncertainty

The verification of the regression based NIPC method for the UQ considering geometric uncertainty is conducted.Firstly, we study the convergence of NIPC by increasing the polynomial chaos orderp. Figures 5–7 show the standard deviation of aerodynamics with the order of NIPC and the comparison with that of MCS. It can be seen that when p=2 or 3, the results of NIPC can coincide with those of MCS. However, there is a little difference between NIPC and MCS. It can also be observed that the results of NIPC are wrong when p=4. We think that the number of terms of PCE is 1820 for p=4 and n=12, which makes it difficult to compute the coefficients correctly by linear regression method. This leads to the wrong results for high-order accuracy.Figure 5
Standard deviation distribution ofCp with the order of NIPC.Figure 6
Standard deviation distribution ofCf with the order of NIPC.Figure 7
Standard deviation distribution of the aerodynamic coefficients with the order of NIPC.Moreover, the convergence of NIPC by increasing the oversampling rationp is also studied. Figures 8–10 show the standard deviation of aerodynamics with np of NIPC and the comparison with that of MCS. When np=4, the results of NIPC can coincide well with those of MCS. Through comparison, it can be observed that the results of NIPC with p=2 and np=4 are even better than those of NIPC with p=3 and np=2.Figure 8
Standard deviation distribution ofCp with oversampling ratio np.Figure 9
Standard deviation distribution ofCf with oversampling ratio np.Figure 10
Standard deviation distribution of the aerodynamic coefficients oversampling rationp.Lastly, the computational costs are compared with MCS. Table3 shows the number of samples to UQ by NIPC and MCS. It can be seen that the number of samples is 364 with p=2 and np=4, which is much lower than that of MCS.Table 3
The number of samples to UQ by NIPC and MCS.
NIPC MCS p = 2, np=2 p = 3, np=2 p = 4, np=2 p = 2, np=4 p = 2, np=6 — N 182 910 3640 364 546 5000


### 4.3.3. The Results of UQ for Transonic Aerodynamics

Figure11 shows the error bars of Cp and Cf distribution by UQ. It can be observed that the aerodynamic loads are sensitive to the geometric variation of the upper surface. In transonic region, the location of shock fluctuates with the change of the geometric shapes, which can result in intensive fluctuations of the flow variables. We call the region in which the standard variations become much larger suddenly shock disturbance region. It can also be seen that the skin-friction behavior displays discrepancies from the Cp distribution. The variation of Cf is similar to that of Cp before and in the shock disturbance region. However, the difference occurs in the downstream region of the shock disturbance region; the standard deviation keeps larger magnitudes below the shock disturbance region and then gradually decreases. This indicates that not only the shock wave but also the boundary layer separation is sensitive to geometric uncertainty in the transonic region.Figure 11
The error bars ofCp and Cf distribution by UQ.


## 4.3.1. The Difficulty of NIPC for the High-Dimensional Problem

For a high-dimensional problem, the number of PCE terms is very high. Now, we study the number of terms withp increasing for the 12-dimensional problem. From Table 1, it can be seen that the number of terms increased rapidly with p increasing for the 12-dimensional problem.Table 1
The number of PCE terms withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 p = 6 N t 13 91 455 1820 6188 18564It is mentioned in (14) that the projection based NIPC method needs p+1n samples for UQ, which grows exponentially with the increasing number of input dimensions. Table 2 shows the number of samples of projection based NIPC method with p for the 12-dimensional problem. It can be seen that the computational cost is too high and unacceptable.Table 2
The number of samples of the projection based NIPC method withp increasing (n=12).
p = 1 p = 2 p = 3 p = 4 p = 5 N s 4096 531441 16777216 244140625 2176782336The regression based NIPC method adopts the linear regression method to calculate the coefficients of PCE. To study how the number of samples influences the analysis results, the oversampling rationp is defined as (16).Compared with the projection based NIPC method, the computational cost of the regression based NIPC method is much lower. Therefore, the regression based NIPC method is used for the stochastic aerodynamic analysis considering geometric uncertainty.


## 4.3.2. The Verification of NIPC for the UQ considering Geometric Uncertainty

The verification of the regression based NIPC method for the UQ considering geometric uncertainty is conducted.Firstly, we study the convergence of NIPC by increasing the polynomial chaos orderp. Figures 5–7 show the standard deviation of aerodynamics with the order of NIPC and the comparison with that of MCS. It can be seen that when p=2 or 3, the results of NIPC can coincide with those of MCS. However, there is a little difference between NIPC and MCS. It can also be observed that the results of NIPC are wrong when p=4. We think that the number of terms of PCE is 1820 for p=4 and n=12, which makes it difficult to compute the coefficients correctly by linear regression method. This leads to the wrong results for high-order accuracy.Figure 5
Standard deviation distribution ofCp with the order of NIPC.Figure 6
Standard deviation distribution ofCf with the order of NIPC.Figure 7
Standard deviation distribution of the aerodynamic coefficients with the order of NIPC.Moreover, the convergence of NIPC by increasing the oversampling rationp is also studied. Figures 8–10 show the standard deviation of aerodynamics with np of NIPC and the comparison with that of MCS. When np=4, the results of NIPC can coincide well with those of MCS. Through comparison, it can be observed that the results of NIPC with p=2 and np=4 are even better than those of NIPC with p=3 and np=2.Figure 8
Standard deviation distribution ofCp with oversampling ratio np.Figure 9
Standard deviation distribution ofCf with oversampling ratio np.Figure 10
Standard deviation distribution of the aerodynamic coefficients oversampling rationp.Lastly, the computational costs are compared with MCS. Table3 shows the number of samples to UQ by NIPC and MCS. It can be seen that the number of samples is 364 with p=2 and np=4, which is much lower than that of MCS.Table 3
The number of samples to UQ by NIPC and MCS.
NIPC MCS p = 2, np=2 p = 3, np=2 p = 4, np=2 p = 2, np=4 p = 2, np=6 — N 182 910 3640 364 546 5000


## 4.3.3. The Results of UQ for Transonic Aerodynamics

Figure11 shows the error bars of Cp and Cf distribution by UQ. It can be observed that the aerodynamic loads are sensitive to the geometric variation of the upper surface. In transonic region, the location of shock fluctuates with the change of the geometric shapes, which can result in intensive fluctuations of the flow variables. We call the region in which the standard variations become much larger suddenly shock disturbance region. It can also be seen that the skin-friction behavior displays discrepancies from the Cp distribution. The variation of Cf is similar to that of Cp before and in the shock disturbance region. However, the difference occurs in the downstream region of the shock disturbance region; the standard deviation keeps larger magnitudes below the shock disturbance region and then gradually decreases. This indicates that not only the shock wave but also the boundary layer separation is sensitive to geometric uncertainty in the transonic region.Figure 11
The error bars ofCp and Cf distribution by UQ.


## 4.4. Global Sensitivity Analysis of Aerodynamics

The important deformation modes can be found out through GSA, which can provide a useful reference for subsequent design work and manufacturing process. In the paper, the NIPC method (p=2; np=4) is used for Sobol’s analysis. The aerodynamic load distribution and aerodynamic coefficients are included for Sobol’s analysis. Figures 12 and 13 show the results of Sobol’s analysis of the pressure and friction coefficient distribution. The partial standard deviation, total standard deviation, and standard deviation caused by the coupling effect for each mode are given in Figures 12 and 13. It can be seen that the deformation modes of the upper surface are much more important than those of lower surface for transonic aerodynamic loads (Cp and Cf). Specifically, the first 3 modes are important for aerodynamic load distribution. For Cp distribution, the main variation of modes 2, 3, and 5 is in shock disturbance region. For Cf distribution, the main variation of modes 2 and 5 is in both shock disturbance region and boundary layer separation region, while the main variation of mode 3 is only in shock disturbance region. In addition, it can also be observed that the coupling effect among geometric deformation modes is weak in Figures 12 and 13. In order to visually observe the effects of modes, the pressure coefficient and friction coefficient of error bars are given in Figures 14 and 15. From Section 4.1, we can observe that mode 1 and mode 2 are the scale modes in the thickness direction; mode 3 and mode 4 are translation modes of the maximum thickness in the axial direction; and mode 5 is the extrusion mode of the upper surface. We can draw the following conclusions: the scale mode and extrusion mode are important for the flow characteristic of shock wave and boundary layer separation, but the translation mode is only important for flow characteristic of shock wave.Figure 12
The results ofCp distribution by Sobol’s analysis (σi is the partial standard variation with uncertain variable i; σTi is the total partial standard variation with uncertain variable i; σcouplingi is the nonlinear coupling effect among uncertain variable i and other uncertainty variables).Figure 13
The results ofCf distribution by Sobol’s analysis (σi is the partial standard variation with uncertain variable i; σTi is the total partial standard variation with uncertain variable i; σcouplingi is the nonlinear coupling effect among uncertain variable i and other uncertainty variables).Figure 14
The error bars ofCp distribution by Sobol’s analysis.Figure 15
The error bars ofCf distribution by Sobol’s analysis.The results of GSA for aerodynamic coefficients including lift coefficient (Cl), drag coefficient (Cd), moment coefficient (Cm), and lift-to-drag ratio (L/D) are shown in Figure 16. For Cl, mode 1 and mode 2 are important, which indicate that the scale modes are important for Cl. For Cd and L/D, modes 2, 3, and 5 are important; the first 3 deformation modes of the upper surface have bigger contribution. For Cm, the scale mode of upper surface is most important. In addition, it can also be seen that the coupling effect among deformation modes is very weak.Figure 16
The results of Sobol’s analysis of aerodynamics ((a) Sobol’s indices forCl; (b) Sobol’s indices for Cd; (c) Sobol’s indices for Cm; (d) Sobol’s indices for L/D).
(a)
Sobol’s indices forCl (b)
Sobol’s indices forCd (c)
Sobol’s indices forCm (d)
Sobol’s indices forL/DFrom the results of GSA, it can be observed that modes 1, 2, 3, and 5 are important for aerodynamics. Therefore, the UQ of aerodynamics is conducted considering the 4 modes. Moreover, because the coupling effects among the modes are small from GSA, we adopt regression method with parametersp=2 and np=4 to build the PCE model without the interaction terms for UQ. The comparison of the UQ results is shown in Table 4. It can be observed that the UQ results are relatively small different from each other; and the number of sample points is reduced from 364 to 32. That is to say, the unimportant modes are eliminated and PCE model can be modified by GSA, which can reduce the computational cost in UQ and robust design.Table 4
The comparison of the UQ results for aerodynamic coefficients between the selected 4 modes and the full 12 modes.
n = 4 (NIPC) n = 12 (NIPC) C l Mean 0.7234 0.7232 StD 0.0096 0.0097 C d Mean 0.0163 0.0163 StD 0.0012 0.0013 C m Mean 0.2718 0.2720 StD 0.0053 0.0057 L / D Mean 44.679 44.627 StD 3.1699 3.3608 N  64 364


## 5. Conclusion

In the paper, the regression based NIPC method was used for UQ and GSA for transonic aerodynamics with geometric uncertainty. It is difficult to describe the variation of geometric shapes, which always requires a large number of parameters. The high dimensions of input parameters can result in the high computational cost of NIPC. To alleviate the contradiction between the number of input dimensions and the computational cost, the PCA technique combined with CST method is used to describe the geometric uncertainty, which reduces the number of input variables. Moreover, typical geometric deformation modes can be obtained by the method.From the results of UQ, the fluctuation of aerodynamic characteristics is mainly centered in the shock disturbance region and boundary layer interference region in transonic region. The characteristics of shock wave are sensitive to geometric uncertainty, and since the transonic drag mainly includes wave drag, the transonic drag is sensitive to geometric uncertainty. From the results of GSA of aerodynamic loads, it is shown that the first three deformation modes of the upper surface have the greatest impact on the aerodynamic characteristics. Specifically, the scale mode and extrusion mode are important for the flow characteristics of shock wave and boundary layer separation, but the translation mode is only important for the flow characteristics of boundary layer separation. Moreover, from the results of GSA of aerodynamic coefficients, we can observe that the scale modes (modes 1 and 2) are important for lift characteristics, while the contribution of the transition modes (modes 3 and 4) is negligible. The first 3 modes of the upper airfoil surface are important to drag and lift-to-drag characteristics. By GSA, we can eliminate unimportant geometric modes and simplify model for UQ and learn which typical modes are important for transonic aerodynamics, which is helpful for transonic airfoil design.



---

*Source: 8107190-2017-02-26.xml*