The above procedure will continue until the convergent criterion is satisfied. Figure 10 shows the case that the sample points in the previous optimization iteration are reused in the two successive iterations. As is seen that two points are located in the overlapping area of two successive sampling regions, thus are reused in the next iteration for regression to identify the significant polynomials/calculate the PCE coefficients. In this way, the computational cost can be further reduced.

Figure 10. Illustration of the reuse of sample points.

4.1. The trust region scenario

40 Uncertainty Quantification and Model Calibration

of the last design point μk

obtained new design point μk

q

j ¼ abs � Yðμk xÞ �

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Xd i¼1 ðμk xi Þ 2

functions with respect to the iteration number k.

where jμk

at μk

xj 2 ¼

x, i.e., jYk

The trust region method is a traditional approach that has been widely used in nonlinear numerical optimization [28]. The basic idea of the trust region method is that in the trust region of the current iteration design point, the second-order Taylor expansion is used to approximate the original objective function. If the accuracy of the current second-order Taylor expansion is satisfied, the size of the trust region is increased to speed up the convergence, and if not it is reduced to improve the accuracy of approximation. To reduce the computational cost of design optimization, the idea of the trust region technique has been extended and applied to reliability-based wing design optimization [29], multifidelity wing aero-structural optimization [30], and multifidelity surrogate-based wing optimization [31], which has been widely believed as an efficient strategy in design optimization. For example, when the trust region technique is applied to meta-model-based design optimization, during optimization, the sample points are sequentially generated in the trust region and the radius of the trust region is

dynamically adjusted based on the accuracy of the meta-model in the local region.

The scenario of trust region is extended here to reduce the computational cost of sDD-PCEbased robust optimization. The basic idea is that the radius of a trust region is determined by the distance between two successive design points and the variation of the corresponding objective function values. If the updated design point μkþ1 x lies in the current trust region, it is considered that the insignificant terms of its PCE model remain unchanged compared to those

Meanwhile, the sample points lying in the overlapping area of two adjacent sampling regions are reused for the PCE coefficient regression for the updated design point to further save the computational cost. Generally, for a practical engineering optimization problem, there is only one performance function that is computationally expensive. Therefore, only one PCE model is required to be constructed and the UP for the rest of the functions can be conveniently implemented by MCS. In this study, it is assumed that the PCE model is only constructed for

the objective function and the general steps of the proposed method is as below.

the Latin Hypercube sample points are generated around μ0

r1 ¼ max ζ1

and |Yk

� � �μk x � � � 2

Step 0: Set the iteration number as k = 1 and the initial staring design point μ0

optimization with sDD-PCE without trust region and obtain a new design variable μk

Step 1: After the kth optimization iteration, define/update the trust region at the current

, ζ2 , r2 ¼ max ζ1

n o n o

x as a rectangle with each length as

X, i.e., the inner loop is eliminated at the updated design point.

� � � Yk � � �, ζ2

; ζ1 and ζ2 are user-defined parameters, which can be constants or

x, do robust

x to calculate the PCE coefficients.

x, where

ð24Þ

4.2. Robust design using sparse data-driven PCE and trust region

#### 4.3. Comparative studies

The first example is the Ackley Function:

$$f(\mathbf{X}) = -20 \exp\left(-0.2\sqrt{\frac{1}{d}\sum\_{j=1}^{d} \mathbf{X}\_j^2}\right) - \exp\left(\frac{1}{d}\sum\_{j=1}^{d} \cos\left(2\pi\mathbf{X}\_j\right)\right) + 22.71282,\\ d = 10\qquad(25)$$

The robust design optimization of this example is:

$$\begin{aligned} \min F &= \mu\_f + k \sigma\_f \\ -10 \le \mu\_{\chi\_j} &\le 10, j = 1, 2, \dots, 10 \end{aligned} \tag{26}$$

All the design variables are considered to follow uniform distribution with variation of �0.2 around their mean values and k in Eq. (26) is set as k = 20. In this study, the fmincon function in Matlab is used for optimization, and ζ1 and ζ2 in Eq. (45) are set as ζ1 = 0.5 and ζ2 = 0.5. Meanwhile, the obtained optimal design variables of sDD-PCE-based robust design with and without trust region scenario as well as the deterministic design without considering any uncertainties are respectively substituted into Eq. (26) through MCS (with 1e6 runs) to calculate the mean μf and standard deviation σf of the objective function.

The results are shown in Table 31, from which it is found that compared to the robust optimization without the trust region scenario (denoted by without), the obtained performance results (μf, σf, and F) of the robust optimization with the trust region scenario (denoted by with) are comparable. However, the number of function calls (denoted as Funcall) is clearly reduced. The decrease in computational cost is attributed to the application of trust region scenario and the reusing of sample points. Meanwhile, the optimal designs of the two robust designs are both less sensitive to uncertainties (smaller σf) compared to the results of deterministic design (denoted by DD). These results demonstrate the effective and advantage of the proposed method.

The second example is the robust design optimization of an automobile torque arm, shown in Figure 11.

In this problem, the four geometrical parameters (a, d1, d2, and l) are considered as design variables, and the yielding strength Sy, Young's modulus E, and the applied force Q are deterministic parameters.

$$\begin{aligned} \min f(a, d\_1, d\_2, l) &= \frac{\pi a d\_2^2}{4} + 2\left(l - \frac{d\_1}{2} - \frac{d\_2}{2}\right) a^2 \\\\ \text{s.t. } g\_1(a, d\_2, l) &= \frac{Q(2l - d\_2) d\_2}{4IS\_y} - 1.0 \le 0 \\\\ g\_2(a, d\_1, d\_2, l) &= 1.0 - \frac{\pi^2 E a^4}{3(2l - d\_1 - d\_2)^2} \frac{d\_2 - d\_1}{Ql} \le 0 \\\\ 5 \le a \le 15, 45 \le d\_1 \le 55, 55 \le d\_2 \le 65, 110 \le l \le 210 \end{aligned} \tag{27}$$

Table 31. Results of the Ackley Function.

4.3. Comparative studies

Figure 11.

deterministic parameters.

The first example is the Ackley Function:

42 Uncertainty Quantification and Model Calibration

0 @

The robust design optimization of this example is:

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 d X d

1

A � exp

min F ¼ μf þ kσf

�10 ≤ μXj

1 d X d

0 @

All the design variables are considered to follow uniform distribution with variation of �0.2 around their mean values and k in Eq. (26) is set as k = 20. In this study, the fmincon function in Matlab is used for optimization, and ζ1 and ζ2 in Eq. (45) are set as ζ1 = 0.5 and ζ2 = 0.5. Meanwhile, the obtained optimal design variables of sDD-PCE-based robust design with and without trust region scenario as well as the deterministic design without considering any uncertainties are respectively substituted into Eq. (26) through MCS (with 1e6 runs) to calcu-

The results are shown in Table 31, from which it is found that compared to the robust optimization without the trust region scenario (denoted by without), the obtained performance results (μf, σf, and F) of the robust optimization with the trust region scenario (denoted by with) are comparable. However, the number of function calls (denoted as Funcall) is clearly reduced. The decrease in computational cost is attributed to the application of trust region scenario and the reusing of sample points. Meanwhile, the optimal designs of the two robust designs are both less sensitive to uncertainties (smaller σf) compared to the results of deterministic design (denoted by

The second example is the robust design optimization of an automobile torque arm, shown in

In this problem, the four geometrical parameters (a, d1, d2, and l) are considered as design variables, and the yielding strength Sy, Young's modulus E, and the applied force Q are

> 2 4 þ 2 l � d1

4ISy

5 ≤ a ≤ 15, 45 ≤ d1 ≤ 55, 55 ≤ d2 ≤ 65, 110 ≤ l ≤ 210

g2ða, d1, d2, lÞ ¼ 1:0 � π2Ea4

2 � d2 2

a2

d2 � d1 Ql ≤ 0 ð27Þ

2

� �

� 1:0 ≤ 0

3ð2l � d1 � d2Þ

DD). These results demonstrate the effective and advantage of the proposed method.

j¼1

cos ð2πXjÞ

1

≤ 10, j ¼ 1, 2, …, 10 ð26Þ

A þ 22:71282, d ¼ 10 ð25Þ

j¼1 X2 j

late the mean μf and standard deviation σf of the objective function.

min fða, d1, d2, lÞ ¼ πad2

s:t: g1ða, d2, lÞ ¼ Qð2l � d2Þd2

vuut

fðX޼�20 exp �0:2

Figure 11. Automobile torque arm.

where the objective function f represents the volume of the arm, the first constraint g1 denotes the yielding failure at section A-A, the second constraint g2 denotes the buckling failure at the two connecting rods, and I ¼ a2ðd2 � aÞ 2 =2 þ a4=6.

The distribution parameters of the four design variables and design parameters are shown in Table 32.

Table 32. Distribution parameters for design variables and design parameters.

The corresponding robust design optimization model is formulated as

$$\begin{aligned} \min F &= \omega\_1 \frac{\mu\_f}{\mu\_f^\*} + \omega\_2 \frac{\sigma\_f}{\sigma\_f^\*}, \omega\_1 = 0.5, \omega\_2 = 0.5\\ \text{s.t.} \quad G\_1(a, d\_2, l) &= \mu\_{\varrho\_1} + k \sigma\_{\varrho\_1} \le 0\\ &\quad G\_2(a, d\_1, d\_2, l) = \mu\_{\varrho\_2} + k \sigma\_{\varrho\_2} \le 0\\ &\quad 5 \le \mu\_a \le 15, 45 \le \mu\_{d\_1} \le 55, 55 \le \mu\_{d\_2} \le 65, 110 \le \mu\_l \le 210 \end{aligned} \tag{28}$$

As has been mentioned above, the PCE model is only constructed for the objective function and the results are shown in Table 33. It is noticed that the robust optimization designs with and without the trust region scenario yields comparable results, while the function calls (objective function calls) required by design with trust region is evidently smaller. The deterministic design cannot even obtain a feasible optimal solution with both constraint violated (>0), since it does not consider uncertainties during design. These results further demonstrate the effectiveness and advantage of the proposed method.

#### 4.4. Summary

The employment of the trust region in sDD-PCE-based robust optimization can evidently reduce the computational cost. However, the determination of the trust region in this chapter is still very subjective and a more rigorous method should be explored. In this section as well as Section 3, the scenarios of sparse PCE and trust region are only employed to DD-PCE to save the computational cost. However, the methods proposed here are also applicable to other PCE approaches, such as gPCE and GS-PCE.

In this chapter, the latest advances in PCE theory and approach for probabilistic UP are comprehensively presented in detail. However, it does not limit the application of PCE to nonprobabilistic UP to address epistemic uncertainties. Sudret and Schöbi have proposed a two-level metamodeling approach using nonintrusive sparse PCE to surrogate the exact computational model to facilitate the uncertainty quantification analysis, in which the input variables are modeled by probability-boxes (p-boxes), accounting for both aleatory and epistemic uncertainty [32]. The Fuzzy uncertainty propagation in composites has been implemented using Gram-Schmidt polynomial chaos expansion, in which the parameter uncertainties are represented by fuzzy membership functions [5]. A general framework has been proposed for a dynamical uncertain system to deal with both aleatory and epistemic uncertainty using PCE, where the uncertain parameters are described through random variables and/or fuzzy variables [33]. The mix UP approach is proposed, in which the inner loop PDFs are calculated using the PCE, and outer loop bounds can be computed with optimization-based interval

Table 33. Results of automobile torque arm.

estimation [34]. PCE has also been applied for solving Bayesian inverse problem as "surrogate posterior." However, it has been indicated that the accuracy cannot always be ensured by PCE since a sufficiently accurate PCE for this problem requires a high order, making PCE impractical compared to directly sampling the posterior [35].

#### Author details

The corresponding robust design optimization model is formulated as

, ω1 ¼ 0:5, ω2 ¼ 0:5

ð28Þ

þ kσg1 ≤ 0

As has been mentioned above, the PCE model is only constructed for the objective function and the results are shown in Table 33. It is noticed that the robust optimization designs with and without the trust region scenario yields comparable results, while the function calls (objective function calls) required by design with trust region is evidently smaller. The deterministic design cannot even obtain a feasible optimal solution with both constraint violated (>0), since it does not consider uncertainties during design. These results further demonstrate

The employment of the trust region in sDD-PCE-based robust optimization can evidently reduce the computational cost. However, the determination of the trust region in this chapter is still very subjective and a more rigorous method should be explored. In this section as well as Section 3, the scenarios of sparse PCE and trust region are only employed to DD-PCE to save the computational cost. However, the methods proposed here are also applicable to other

In this chapter, the latest advances in PCE theory and approach for probabilistic UP are comprehensively presented in detail. However, it does not limit the application of PCE to nonprobabilistic UP to address epistemic uncertainties. Sudret and Schöbi have proposed a two-level metamodeling approach using nonintrusive sparse PCE to surrogate the exact computational model to facilitate the uncertainty quantification analysis, in which the input variables are modeled by probability-boxes (p-boxes), accounting for both aleatory and epistemic uncertainty [32]. The Fuzzy uncertainty propagation in composites has been implemented using Gram-Schmidt polynomial chaos expansion, in which the parameter uncertainties are represented by fuzzy membership functions [5]. A general framework has been proposed for a dynamical uncertain system to deal with both aleatory and epistemic uncertainty using PCE, where the uncertain parameters are described through random variables and/or fuzzy variables [33]. The mix UP approach is proposed, in which the inner loop PDFs are calculated using the PCE, and outer loop bounds can be computed with optimization-based interval

X μf σf F G1 G2 Funcall

DD [8.13, 55.00, 55.00, 110.00] 2.6616e4 1.2171e3 1 0.1848 5.1509e4 82 with [8.53, 54.10, 58.67, 111.03] 3.1027e4 1.3355e3 1.1315 �0.0123 �1.1833e5 658 without [8.57, 52.68, 57.50, 110.00] 3.0332e4 1.3093e3 1.1077 �4.0000e�4 �1.2913e2 1283

þ kσg2 ≤ 0 5 ≤ μa ≤ 15, 45 ≤ μd1 ≤ 55, 55 ≤ μd2 ≤ 65, 110 ≤ μl ≤ 210

μf μ� f þ ω2 σf σ� f

G2ða, d1, d2, lÞ ¼ μg2

s:t: G1ða, d2, lÞ ¼ μg1

the effectiveness and advantage of the proposed method.

PCE approaches, such as gPCE and GS-PCE.

4.4. Summary

μ�

Table 33. Results of automobile torque arm.

min F ¼ ω1

44 Uncertainty Quantification and Model Calibration

Shuxing Yang\*, Fenfen Xiong and Fenggang Wang

\*Address all correspondence to<EMAIL_ADDRESS>
School of Aerospace Engineering, Beijing Institute of Technology, Beijing, China

#### References

[25] Hu C, Youn BD. Adaptive-sparse polynomial chaos expansion for reliability analysis and design of complex engineering systems. Structural & Multidisciplinary Optimization. 2011; 43(3):419–442

[10] Coelho R, Bouillard P. Multi-objective reliability-based optimization with stochastic

[11] Xiu D, Karniadakis GE. The wiener-askey polynomial chaos for stochastic differential

[12] Wiener N. The homogeneous chaos. American Journal of Mathematics. 1938;60(1):897–936 [13] Fan et al. Parameter uncertainty and temporal dynamics of sensitivity for hydrologic models: A hybrid sequential data assimilation and probabilistic collocation method.

[14] Guerine A, Hami AE, Walha L, et al. A polynomial chaos method for the analysis of the dynamic behavior of uncertain gear friction system. European Journal of Mechanics - A/

[15] Meecham WC, Siegel A. Wiener-Hermite expansion in model turbulence at large Reynolds numbers. Physics of Fluids (1958-1988). 1964;7(8):1178–1190. DOI: 10.1063/1.1711359

[16] Witteveen JAS, Bijl H. Modeling arbitrary uncertainties using Gram-Schmidt polynomial chaos. 44th AIAA Aerospace Sciences Meeting and Exhibit; 9–12 January 2006; Reno,

[17] Wan X, Karniadakis GE. Multi-element generalized polynomial chaos for arbitrary prob-

[18] Oladyshkin S, Nowak W. Data-driven uncertainty quantification using the arbitrary polynomial chaos expansion. Reliability Engineering & System Safety. 2012;106(4):179–190 [19] Hosder S, Walters RW, Balch M. Efficient sampling for non-intrusive polynomial chaos applications with multiple uncertain input variables. 48th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics and Materials Conference; 23–26 April 2007; Honolulu,

[20] Eldred MS. Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design. 50th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference; 4–7 May 2009; Palm Springs, California [21] Abramowitz M, Stegun I, Mcquarrie D A. Handbook of Mathematical Functions. Dover

[22] Xiong F, Greene S, Chen W, Xiong Y, Yang S A new sparse grid based method for uncertainty propagation. Structural & Multidisciplinary Optimization. 2009;41(3):335–349 [23] Efron B, Hastie T, Johnstone I, Tibshirani R. Least angle regression. Mathematics. 2004;32

[24] Tatang MA, Pan W, Prinn RG, McRae GJ. An efficient method for parametric uncertainty analysis of numerical geophysical models. Journal of Geophysics Research. 1997;102

ability measures. SIAM Journal on Scientific Computing. 2006;28(3):901–928

metamodels. Evolutionary Computation. 2011;19(4):525–560

Environmental Modelling & Software. 2016;86:30–49

Solids. 2016;59:76-84

46 Uncertainty Quantification and Model Calibration

Nevada

Hawall

(2):407–499

(D18):21925–21932

Publications, New York,1964

equations. SIAM Journal on Scientific Computing. 2002;24(2):619–644

### **State‐of‐the‐Art Nonprobabilistic Finite Element Analyses**

Wang Lei, Qiu Zhiping and Zheng Yuning

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/intechopen.68154

#### Abstract

The finite element analysis of a mechanical system is conventionally performed in the context of deterministic inputs. However, uncertainties associated with material properties, geometric dimensions, subjective experiences, boundary conditions, and external loads are ubiquitous in engineering applications. The most popular techniques to handle these uncertain parameters are the probabilistic methods, in which uncertainties are modeled as random variables or stochastic processes based on a large amount of statistical information on each uncertain parameter. Nevertheless, subjective results could be obtained if insufficient information unavailable and nonprobabilistic methods can be alternatively employed, which has led to elegant procedures for the nonprobabilistic finite element analysis. In this chapter, each nonprobabilistic finite element analysis method can be decomposed as two individual parts, i.e., the core algorithm and preprocessing procedure. In this context, four types of algorithms and two typical preprocessing procedures as well as their effectiveness were described in detail, based on which novel hybrid algorithms can be conceived for the specific problems and the future work in this research field can be fostered.

Keywords: interval finite element method, fuzzy finite element method, arithmetic approach, perturbation approach, sampling approach, optimization approach, subinterval technique, surrogate model

#### 1. Introduction

The traditional finite element analysis (FEA) was performed in the context of deterministic parameters. However, uncertainties associated with material properties, geometric dimensions, and external loads are always unavoidable in engineering. The ability to include uncertainties is of great value for a design engineer. In the last decade, criticism has arisen regarding the general application of the probabilistic concept. Especially when the statistical information

© 2017 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

on uncertainties is limited [1], the subjective probabilistic analysis result may be obtained by the probabilistic method [2, 3], which proves to be of little value and does not justify the high computational cost [3–5]. Consequently, nonprobabilistic concepts have been introduced.

In this context, interval and fuzzy approaches are gaining more and more momentum for the uncertainty analysis and optimization of numerical models in their descriptions. In the interval approach, uncertainties are considered to be contained within a predefined range and only the lower and upper bounds are necessary for each uncertain parameter. The fuzzy approach further extends this methodology by the α-level technique, where α stands for the extent that a specific value is member of the range of possible input values. From this viewpoint, a fuzzy analysis requires the consecutive solution for a number of interval analysis based on the α-level technique [6]. For this reason, current researches on nonprobabilistic uncertainty propagation mainly focus on the solution and implementation of the interval analysis. In the past decades, the interval and fuzzy concepts in FEA have been studied extensively and some typical solution schemes for the interval FEA (IFEA) and fuzzy FEA (FFEA) were developed. This chapter is to give an overview of state-of-the-art numerical implementations of IFEA and FFEA in applied mechanics.

FFEA aims to obtain a fuzzy description of an FEA result, starting from fuzzy descriptions of all uncertainties. The α-level technique subdivides the membership function range into a number of discrete α-levels. The α-cuts of the input quantities are defined as xiα ¼ {xi ∈ Xi, μ~xi ðxiÞ ≥ α} where μ~x ðxÞ is the membership function. This means that an α-cut is the interval resulting from intersecting the membership function at μ~xi ðxiÞ ¼ α. The α-level interval describes the grade of membership to the fuzzy set for each element in the domain and enables the representation of a value that is only to a certain degree member of the set. However, the confidence interval defined in statistics is the range of likely values for a population parameter, such as the population mean. The selection of a confidence level for an interval determines the probability that confidence interval produced will contain the true parameter value. The intersection with the membership function of the input uncertainties on each α-level results in an interval and IFEA is formulated, resulting in an interval for the output on the considered α-level. The fuzzy solution is finally assembled from the resulting intervals on each sublevel. The IFEA is based on the interval description of uncertainties and its goal is to capture the range of specific output quantities of interest that corresponds to a given interval description of input uncertainties. For the sake of simplicity, the static analysis of a mechanical system is adopted in this chapter to explain current IFEA schemes. The FEA equation can be expressed in a general form as follows:

$$\mathbf{K(p)}\mathbf{U(p)} = \mathbf{F(p)}\tag{1}$$

where K and F stand for the stiffness matrix and load vector, respectively; U represents the static response vector; and p is the input parameter vector of the mechanical system. In the IFEA, p is quantified as an interval input vector pI and shown in Figure 1.

where pc i is the nominal value, Δpi is the interval radius. Then, the IFEA equation is accordingly rewritten as follows:

$$\mathbf{K(p^l)}\mathbf{U(p^l)} = \mathbf{F(p^l)}\tag{2}$$

Figure 1. The diagram of interval variable p.

on uncertainties is limited [1], the subjective probabilistic analysis result may be obtained by the probabilistic method [2, 3], which proves to be of little value and does not justify the high computational cost [3–5]. Consequently, nonprobabilistic concepts have been introduced.

In this context, interval and fuzzy approaches are gaining more and more momentum for the uncertainty analysis and optimization of numerical models in their descriptions. In the interval approach, uncertainties are considered to be contained within a predefined range and only the lower and upper bounds are necessary for each uncertain parameter. The fuzzy approach further extends this methodology by the α-level technique, where α stands for the extent that a specific value is member of the range of possible input values. From this viewpoint, a fuzzy analysis requires the consecutive solution for a number of interval analysis based on the α-level technique [6]. For this reason, current researches on nonprobabilistic uncertainty propagation mainly focus on the solution and implementation of the interval analysis. In the past decades, the interval and fuzzy concepts in FEA have been studied extensively and some typical solution schemes for the interval FEA (IFEA) and fuzzy FEA (FFEA) were developed. This chapter is to give an overview

of state-of-the-art numerical implementations of IFEA and FFEA in applied mechanics.

discrete α-levels. The α-cuts of the input quantities are defined as xiα ¼ {xi ∈ Xi, μ~xi

IFEA schemes. The FEA equation can be expressed in a general form as follows:

IFEA, p is quantified as an interval input vector pI and shown in Figure 1.

KðpI

ÞUðpI

intersecting the membership function at μ~xi

50 Uncertainty Quantification and Model Calibration

where pc

rewritten as follows:

FFEA aims to obtain a fuzzy description of an FEA result, starting from fuzzy descriptions of all uncertainties. The α-level technique subdivides the membership function range into a number of

μ~x ðxÞ is the membership function. This means that an α-cut is the interval resulting from

membership to the fuzzy set for each element in the domain and enables the representation of a value that is only to a certain degree member of the set. However, the confidence interval defined in statistics is the range of likely values for a population parameter, such as the population mean. The selection of a confidence level for an interval determines the probability that confidence interval produced will contain the true parameter value. The intersection with the membership function of the input uncertainties on each α-level results in an interval and IFEA is formulated, resulting in an interval for the output on the considered α-level. The fuzzy solution is finally assembled from the resulting intervals on each sublevel. The IFEA is based on the interval description of uncertainties and its goal is to capture the range of specific output quantities of interest that corresponds to a given interval description of input uncertainties. For the sake of simplicity, the static analysis of a mechanical system is adopted in this chapter to explain current

where K and F stand for the stiffness matrix and load vector, respectively; U represents the static response vector; and p is the input parameter vector of the mechanical system. In the

i is the nominal value, Δpi is the interval radius. Then, the IFEA equation is accordingly

Þ¼ FðpI

ðxiÞ ≥ α} where

ðxiÞ ¼ α. The α-level interval describes the grade of

KðpÞUðpÞ¼ FðpÞ ð1Þ

Þ ð2Þ

where the superscript "I" hereinafter represents an interval input. The exact solution set of this interval equation can be expressed as:

$$\mathbf{U} = \left\{ \overline{\mathbf{U}} | \mathbf{K}(\overline{\mathbf{p}}) \overline{\mathbf{U}} = \mathbf{F}(\overline{\mathbf{p}}), \forall \overline{\mathbf{p}} \in \mathbf{p}^{\mathrm{I}} \right\} \tag{3}$$

It is noted that interdependencies among entries of the response vector are introduced due to sharing the common input vector and a nonconvex polyhedron is always defined [7], which makes it extremely difficult to obtain the exact solution [5]. However, only individual ranges of some components in the response vector are of interest for real-life problems. Therefore, by neglecting the aforementioned interdependencies, the smallest hypercube approximation denoted as UI around the exact solution set is an alternative object for current IFEA. The kth component of UI is expressed as follows:

$$\mathcal{U}\_k^\mathrm{I} = \left[ \mathcal{U}\_{k'}^\mathrm{L} \, \mathcal{U}\_k^\mathrm{U} \right] = \left[ \min\_{\mathbf{p} \in \mathbf{p}^\mathrm{I}} \mathcal{U}\_k(\mathbf{p}), \max\_{\mathbf{p} \in \mathbf{p}^\mathrm{I}} \mathcal{U}\_k(\mathbf{p}) \right], \; k = 1, 2, \ldots, N \tag{4}$$

where superscripts "L" and "U" represent the lower and upper bounds of an interval variable, respectively; N is the total number of response components of interest. Accordingly, the smallest hypercube solution of IFEA equation is expressed as:

$$\mathbf{U}^{\mathrm{I}} = \left[ \mathbf{U}\_{1\prime}^{\mathrm{I}} \mathbf{U}\_{2\prime}^{\mathrm{I}} \dots \mathbf{U}\_{N}^{\mathrm{I}} \right]^{\mathrm{T}} \tag{5}$$

where "T" is a transposition operator.

#### 2. Core algorithms

From published literatures, four types of algorithms for IFEA have been well established. Most of the current schemes are formulated based on these core algorithms.

#### 2.1. Arithmetic approach

The key point of arithmetic approach is to translate the complete deterministic numerical FE procedure to an interval procedure using the arithmetic operations. Each substep of the interval algorithm calculates the range of the intermediate subfunction instead of the deterministic result. Based on this principle, the interval bounds of the output can be obtained. The original solution procedure for IFEA is the interval arithmetic approaches [7–10], in which all basic deterministic algebraic operations are replaced by their interval arithmetic counterparts.

The major advantage of the arithmetic approach is its simplicity. However, the major drawback of this method is its repeated vulnerability to conservatism. It is shown that these methods suffer considerably from the overestimation effect, also referred to as the dependency problem, and for the real-life problems, the resulting overestimation may render the final result totally useless [5]. A simple example is shown as follows. Consider the function

$$f(\mathbf{x}) = \mathbf{x}^2 - \mathbf{x} + \mathbf{1} \tag{6}$$

applied on the interval x ¼ ½0, 1�. Applying arithmetic approach, both terms are assumed independently. This leads to the interval solution fðxÞ¼½0, 2�. However, the exact range of the function equals fðxÞ ¼ 3 4 , 1 � �. That is to say, an arithmetic interval operation introduces conservatism in its result if neglecting the correlation that exists between the operands. Besides, the integration of interval arithmetic approaches with software for FEA is also a challenge in real applications.

#### 2.2. Perturbation approach

The perturbation approach has been widely applied in structural response analyses and other applications. Compared to arithmetic approaches, perturbation methods are more popular due to its simplicity and efficiency in IFEA and can be available in the original, improved, and modified versions.

#### 2.2.1. Original version

The first-order Taylor expansions of the interval stiffness matrix and load vector at the nominal (mid-) values of interval parameters were firstly obtained as:

$$\begin{aligned} \mathbf{K(p^l)} &= \mathbf{K(p^c)} + \sum\_{i=1}^n \frac{\partial \mathbf{K(p^c)}}{\partial p\_i} \Delta p\_i^l = \mathbf{K^c} + \Delta \mathbf{K^l} \\ \mathbf{F(p^l)} &= \mathbf{F(p^c)} + \sum\_{i=1}^n \frac{\partial \mathbf{F(p^c)}}{\partial p\_i} \Delta p\_i^l = \mathbf{F^c} + \Delta \mathbf{F^l} \end{aligned} \tag{7}$$

where pc is the nominal (mid-) value of the interval input vector and ΔpI i ¼ ½�Δpi ,Δpi � is the interval radius of the ith interval parameter, i.e.,

State‐of‐the‐Art Nonprobabilistic Finite Element Analyses http://dx.doi.org/10.5772/intechopen.68154 53

$$\begin{aligned} \mathbf{p}^c = (\mathbf{p}^\mathsf{U} + \mathbf{p}^\mathsf{L})/2 &= \left[p\_1^c, p\_2^c, \dots, p\_n^c\right]^\mathsf{T} \\ \Delta \mathbf{p} = (\mathbf{p}^\mathsf{U} - \mathbf{p}^\mathsf{L})/2 &= \left[\Delta p\_1, \Delta p\_2, \dots, \Delta p\_n\right]^\mathsf{T} \end{aligned} \tag{8}$$

And the interval radiuses of the stiffness matrix and load vector in Eq. (7) are expressed as follows, respectively.

$$\begin{aligned} \Delta \mathbf{K}^{\mathrm{I}} &= \sum\_{i=1}^{n} \frac{\partial \mathbf{K}(\mathbf{p}^{c})}{\partial p\_{i}} \Delta p\_{i}^{\mathrm{I}} = [-\Delta \mathbf{K}, \Delta \mathbf{K}] \\ \Delta \mathbf{F}^{\mathrm{I}} &= \sum\_{i=1}^{n} \frac{\partial \mathbf{F}(\mathbf{p}^{c})}{\partial p\_{i}} \Delta p\_{i}^{\mathrm{I}} = [-\Delta \mathbf{F}, \Delta \mathbf{F}] \end{aligned} \tag{9}$$

The FEA model for the perturbed system can be rewritten as follows:

$$(\mathbf{K}^{\mathbf{c}} + \Delta \mathbf{K}^{\mathbf{l}})(\mathbf{U}^{\mathbf{c}} + \Delta \mathbf{U}^{\mathbf{l}}) = \mathbf{F}^{\mathbf{c}} + \Delta \mathbf{F}^{\mathbf{l}} \tag{10}$$

By expanding Eq. (10) and neglecting the second-order perturbed term, the following equations can be obtained.

$$\begin{aligned} \mathbf{U}^{\mathbf{c}} &= \left(\mathbf{K}^{\mathbf{c}}\right)^{-1} \mathbf{F}^{\mathbf{c}}\\ \mathbf{K}^{\mathbf{c}} \Delta \mathbf{U}^{\mathbf{l}} &= \Delta \mathbf{F}^{\mathbf{l}} - \Delta \mathbf{K}^{\mathbf{l}} (\mathbf{K}^{\mathbf{c}})^{-1} \mathbf{F}^{\mathbf{c}} \end{aligned} \tag{11}$$

Substituting Eq. (10) into Eq. (11) yields the interval radius of the response vector as:

$$\Delta \mathbf{U}^{\mathrm{I}} = (\mathbf{K}^{\mathrm{c}})^{-1} \sum\_{i=1}^{n} \frac{\partial \mathbf{F}(\mathbf{p}^{\mathrm{c}})}{\partial p\_{i}} \Delta p\_{i}^{\mathrm{I}} - (\mathbf{K}^{\mathrm{c}})^{-1} \sum\_{i=1}^{n} \frac{\partial \mathbf{K}(\mathbf{p}^{\mathrm{c}})}{\partial p\_{i}} \Delta p\_{i}^{\mathrm{I}} (\mathbf{K}^{\mathrm{c}})^{-1} \mathbf{F}^{\mathrm{c}} \tag{12}$$

And the radius vector of the response vector is estimated in the original interval perturbation method [11] as follows:

$$
\Delta \mathbf{U} = \sum\_{i=1}^{n} \left( \left| \left( \mathbf{K}^{\mathbf{c}} \right)^{-1} \right| \left| \frac{\partial \mathbf{F}(\mathbf{p}^{\mathbf{c}})}{\partial p\_{i}} \right| + \left| \left( \mathbf{K}^{\mathbf{c}} \right)^{-1} \right| \left| \frac{\partial \mathbf{K}(\mathbf{p}^{\mathbf{c}})}{\partial p\_{i}} \right| \left| \left( \mathbf{K}^{\mathbf{c}} \right)^{-1} \right| \left| \mathbf{F}^{\mathbf{c}} \right| \right) \Delta p\_{i} \tag{13}
$$

The smallest hypercube solution can thus be determined as:

$$\mathbf{U}^{\mathrm{I}} = [\mathbf{U}^{\mathrm{c}} - \Delta \mathbf{U}, \mathbf{U}^{\mathrm{c}} + \Delta \mathbf{U}] \tag{14}$$

The major drawback of this method is that a significant overestimation is introduced by the original interval perturbation method, indicating that it applies to the interval analysis of problems with "small" interval parameters.

#### 2.2.2. Improved version

2.1. Arithmetic approach

52 Uncertainty Quantification and Model Calibration

counterparts.

the function

applications.

function equals fðxÞ ¼ 3

2.2. Perturbation approach

and modified versions.

2.2.1. Original version

The key point of arithmetic approach is to translate the complete deterministic numerical FE procedure to an interval procedure using the arithmetic operations. Each substep of the interval algorithm calculates the range of the intermediate subfunction instead of the deterministic result. Based on this principle, the interval bounds of the output can be obtained. The original solution procedure for IFEA is the interval arithmetic approaches [7–10], in which all basic deterministic algebraic operations are replaced by their interval arithmetic

The major advantage of the arithmetic approach is its simplicity. However, the major drawback of this method is its repeated vulnerability to conservatism. It is shown that these methods suffer considerably from the overestimation effect, also referred to as the dependency problem, and for the real-life problems, the resulting overestimation may render the final result totally useless [5]. A simple example is shown as follows. Consider

applied on the interval x ¼ ½0, 1�. Applying arithmetic approach, both terms are assumed independently. This leads to the interval solution fðxÞ¼½0, 2�. However, the exact range of the

vatism in its result if neglecting the correlation that exists between the operands. Besides, the integration of interval arithmetic approaches with software for FEA is also a challenge in real

The perturbation approach has been widely applied in structural response analyses and other applications. Compared to arithmetic approaches, perturbation methods are more popular due to its simplicity and efficiency in IFEA and can be available in the original, improved,

The first-order Taylor expansions of the interval stiffness matrix and load vector at the nominal

∂KðpcÞ ∂pi

∂Fðpc Þ ∂pi

ΔpI

ΔpI

i ¼ Kc þ ΔKI

ð7Þ

� is the

i ¼ ½�Δpi

,Δpi

i ¼ Fc þ ΔFI

Þ þXn i¼1

i¼1

(mid-) values of interval parameters were firstly obtained as:

Þ ¼ Kðpc

Þ ¼ FðpcÞ þXn

where pc is the nominal (mid-) value of the interval input vector and ΔpI

KðpI

FðpI

interval radius of the ith interval parameter, i.e.,

fðxÞ ¼ x2 � x þ 1 ð6Þ

4 , 1 � �. That is to say, an arithmetic interval operation introduces conser-

The most typical improved interval perturbation method was proposed in Ref. [12], in which the radius vector of the response vector was calculated as follows:

$$
\Delta \mathbf{U} = \sum\_{i=1}^{n} \left| (\mathbf{K^c})^{-1} \frac{\partial \mathbf{F}(\mathbf{p^c})}{\partial p\_i} - (\mathbf{K^c})^{-1} \frac{\partial \mathbf{K}(\mathbf{p^c})}{\partial p\_i} (\mathbf{K^c})^{-1} \mathbf{F}^c \right| \Delta p\_i \tag{15}
$$

Accordingly, the smallest hypercube solution of IFEA can also be determined by Eq. (14). Although with better accuracy compared to the original one, an interval translation effect, i.e., the translation of the resulting interval w.r.t. the accurate one, is always introduced by the improved interval perturbation method.

#### 2.2.3. Modified versions

Compared with the original version of the perturbation approach where only first-order terms are considered, the main aspect of the following two modified interval perturbation methods [13, 14] is that the interval bounds are calculated by retaining part of higher order terms in Neumann series. Therefore, the modified methods can obtain more accurate response bounds. The key expressions are summarized as follows:

$$\begin{aligned} \mathbf{U}^{\mathbf{c}} &= (\mathbf{K}^{\mathbf{c}})^{-1} \left[ \mathbf{I} + \sum\_{i=1}^{n} \mathbf{E}\_{i}^{\mathbf{c}} \right] \mathbf{F}^{\mathbf{c}} \\ \Delta \mathbf{U}^{\mathbf{I}} &= \sum\_{k=1}^{n} \left\{ (\mathbf{K}^{\mathbf{c}})^{-1} \left[ \mathbf{I} + \sum\_{i=1}^{n} \mathbf{E}\_{i}^{\mathbf{c}} \right] \frac{\partial \mathbf{F}(\mathbf{p}^{\mathbf{c}})}{\partial p\_{k}} \right\} \Delta p\_{k}^{\mathbf{I}} + \sum\_{i=1}^{n} (\mathbf{K}^{\mathbf{c}})^{-1} \Delta \mathbf{E}\_{i}^{\mathbf{I}} \mathbf{F}^{\mathbf{c}} \end{aligned} \tag{16}$$

where

$$\begin{aligned} \mathbf{E}\_i^c &= \left( (\mathbf{I} + \Delta p\_i \mathbf{K}\_i)^{-1} + (\mathbf{I} - \Delta p\_i \mathbf{K}\_i)^{-1} - 2\mathbf{I} \right) / 2 \\ \Delta \mathbf{E}\_i &= \left| (\mathbf{I} + \Delta p\_i \mathbf{K}\_i)^{-1} - (\mathbf{I} - \Delta p\_i \mathbf{K}\_i)^{-1} \right| / 2 \end{aligned} \tag{17}$$

and

$$\mathbf{K}\_{i} = \frac{\partial \mathbf{K}(\mathbf{p}^{\epsilon})}{\partial p\_{i}} (\mathbf{K}^{\epsilon})^{-1} \tag{18}$$

Different estimations of the radius vector of the response vector were, respectively, obtained as follows:

$$\Delta \mathbf{U} = \left| \sum\_{k=1}^{n} \left\{ (\mathbf{K}^{\mathbf{c}})^{-1} \left[ \mathbf{I} + \sum\_{i=1}^{n} \mathbf{E}\_{i}^{\mathbf{c}} \right] \frac{\partial \mathbf{F}(\mathbf{p}^{\mathbf{c}})}{\partial p\_{k}} \right\} \Delta p\_{k} + \sum\_{i=1}^{n} (\mathbf{K}^{\mathbf{c}})^{-1} \Delta \mathbf{E}\_{i} \mathbf{F}^{\mathbf{c}} \right| \tag{19}$$

$$\Delta \mathbf{U} = \sum\_{k=1}^{n} \left| (\mathbf{K}^{\mathbf{c}})^{-1} \left[ \mathbf{I} + \sum\_{i=1}^{n} \mathbf{E}\_{i}^{\mathbf{c}} \right] \frac{\partial \mathbf{F}(p^{\mathbf{c}})}{\partial p\_{k}} \right| \Delta p\_{k} + \sum\_{i=1}^{n} \left| (\mathbf{K}^{\mathbf{c}})^{-1} \Delta \mathbf{E}\_{i} \mathbf{F}^{\mathbf{c}} \right| \tag{20}$$

It should be pointed out that significant unpredicted estimation is always introduced by Eqs. (19) and (20). A more reasonable estimation of the radius vector of the response vector is simultaneously determined herein as follows:

$$
\Delta \mathbf{U} = \sum\_{k=1}^{n} \left| (\mathbf{K^c})^{-1} \left[ \mathbf{I} + \sum\_{i=1}^{n} \mathbf{E\_i^c} \right] \frac{\partial \mathbf{F(p^c)}}{\partial p\_k} \right| \Delta p\_k + \sum\_{i=1}^{n} \left| (\mathbf{K^c})^{-1} \right| \Delta \mathbf{E\_i} \left| \mathbf{F^c} \right| \tag{21}
$$

And a slight conservatism is alternatively resulted in by Eq. (21). The smallest hypercube solution for the IFEA is finally determined as Eq. (14). It is worth mentioning that the spectral radius of ðKc Þ �1 ΔK increases with the increase in ΔKI . ðKc þ ΔKÞ �1 can be expanded with a Neumann series if and only if kðKc Þ �1 ΔKk is less than 1 based on the criteria of convergence for a Neumann series. Therefore, these methods applies to the interval analysis of nonlinear problems with "small" interval parameters and the accuracy for those with "large" interval inputs can be improved by the subinterval technique in Section 3.1. Furthermore, the integration of all interval perturbation methods with current FEA software for the system simulation remains a great challenge.

#### 2.3. Sampling approach

#### 2.3.1. Vertex method

ΔU ¼ Xn

improved interval perturbation method.

54 Uncertainty Quantification and Model Calibration

The key expressions are summarized as follows:

Uc ¼ ðKc

Xn k¼1

> Ec i ¼ �

ΔUI ¼

ΔU ¼ Xn

ΔU ¼ Xn k¼1

� � � � �

k¼1

ðKc Þ �1 Iþ

ðKc Þ �1 Iþ

� � � � �

Þ �1 Iþ

> ðKc Þ �1 Iþ

Xn i¼1 Ec i

ðI þ Δpi

ΔEi ¼ ðI þ Δpi

� � � � � Fc

Xn i¼1 Ec i

KiÞ

KiÞ

Ki ¼ ∂Kðpc

Xn i¼1 Ec i

Xn i¼1 Ec i

" #

" #

( )

" #

( )

" #

2.2.3. Modified versions

where

and

follows:

i¼1

ðKc Þ

� � � � �

�1 ∂Fðpc Þ ∂pi

� ðKc Þ

Accordingly, the smallest hypercube solution of IFEA can also be determined by Eq. (14). Although with better accuracy compared to the original one, an interval translation effect, i.e., the translation of the resulting interval w.r.t. the accurate one, is always introduced by the

Compared with the original version of the perturbation approach where only first-order terms are considered, the main aspect of the following two modified interval perturbation methods [13, 14] is that the interval bounds are calculated by retaining part of higher order terms in Neumann series. Therefore, the modified methods can obtain more accurate response bounds.

> ∂Fðpc Þ ∂pk

�1 þ ðI � Δpi

�1 � ðI � Δpi

Þ ∂pi

> ∂Fðpc Þ ∂pk

∂Fðpc Þ ∂pk

� � � � �

Different estimations of the radius vector of the response vector were, respectively, obtained as

ðKc Þ

ΔpI k þXn i¼1 ðKc Þ �1 ΔEI iFc

KiÞ

KiÞ �1

Δpk þXn

Δpk þXn

i¼1 ðKc Þ �1 ΔEiFc

� � � � �

ðKc Þ �1 ΔEiFc

i¼1

�1 � 2I � =2

> � � � � � =2

�1 ð18Þ

� � � � �

� � � � �

�1 ∂Kðpc Þ ∂pi

ðKc Þ �1 Fc

� � � � �

Δpi ð15Þ

ð16Þ

ð17Þ

ð19Þ

ð20Þ

The vertex method was originally developed in Ref. [15], which can be viewed as a sampling technique with vertices being input samples of the FEA model. This method has been popular for the implementation of IFEA [16–21] due to its main aspect of simple formulation and the black-box property. If the behavior of the target response w.r.t. uncertain parameters can be guaranteed to be monotonic, the vertex method firstly proposed in Ref. [15] yields the exact solution. It should be pointed out that the concept of monotonicity in this section means monotonic along all principal directions where only one parameter is changing at a time. However, it is very hard—if not impossible —to prove the property of monotonicity in a general way, e.g., in the application of structural dynamics [22]. The number of FEA runs necessary for the vertex method is given as:

$$N = \mathcal{Z}^n\tag{22}$$

where n is the number of interval parameters. It is noted that the computational cost for the vertex method exponentially increases w.r.t. the number of interval parameters, which results in a dimensionality curse.

#### 2.3.2. Transformation method

To promote the accuracy of the vertex method for nonmonotonic problems, transformation methods for the epistemic uncertainty propagation were developed. Its original version was firstly proposed in literature [23]. This method is based on the α-level strategy and on each α-level the interval problem is defined. The interval solution strategy then consists of a dedicated sampling strategy in the space spanned by α-cut of fuzzy parameters. This method is available in a general, a reduced, and an extended form, with the most appropriate form to be selected depending on the type of model to be evaluated [23, 24]. If the behavior of the target response w.r.t. uncertain parameters can be guaranteed to be monotonic, the reduced transformation method yields the exact solution. If it shows nonmonotonic behavior, instead, the extended transformation method can be applied, in which more observation points were added in a well-directed way to the search domain after rating the monotonicity of the response w.r.t. different uncertain parameters on the basis of a classification criterion [24].

The computational cost of the transformation method is governed by the number of FEA runs N to be performed. In the case of the general transformation method, this number is given as:

$$N = \sum\_{k=1}^{m+1} k^n \tag{23}$$

where m is the number of discrete α-levels and n is the number of fuzzy parameters. It is noted that the number of FEA runs grows exponentially w.r.t. the number of fuzzy inputs, which makes the general transformation method computational tedious for high-dimensional problems. The main aspect of the transformation method, its characteristic property of reducing fuzzy arithmetic to multiple crisp-number operations entails that this method can be implemented without major problems into an existing software environment for system simulation. Expensive rewriting of the program codes is not required [25]. Some of the most recent applications can be found in Refs. [25–32]. Besides, a program named as FAMOUS (fuzzy arithmetical modeling of uncertain systems) has been developed [25], which provides an interface to commercial software environments. Primarily developed in Matlab environment, FAMOUS actually works as a standalone application on both Windows and Linux platforms.

#### 2.4. Optimization approach

In essence, calculating the solution set expressed in Eq. (3) is equivalent to performing a global optimization, aimed at the minimization and maximization of the components of the deterministic analysis results {U}. The lower and upper bounds of the output of a classical FEA model are determined by the optimization approach through a search algorithm within the domain spanned by the interval parameters. If the global minimum and maximum of the analysis result are found by the search algorithm, it returns the smallest hypercube solution around the exact one. The optimization is performed independently on each element of the response vector. Furthermore, as the behavior of the target response w.r.t. uncertain parameters is rather unpredictable, the computational cost of the optimization approach in general is strongly problem-dependent. It is noted that the optimization approach is immune to the excessive conservatism for the interval arithmetic approaches because the optimization strategy approaches the smallest hypercube solution from its inside, which means that it does not guarantee conservatism until the actual bounds are captured. Additionally, the smooth behavior of the target response w.r.t. uncertain parameters facilitates the search for the global extrema over the space spanned by uncertain parameters. The directional search-based algorithm [16, 33, 34], linear programming [35], and genetic algorithm [36] were utilized to formulate the procedure of IFEA or FFEA. More applications can be found in [37–39]. It is worth mentioning that the optimization approach and Monte Carlo simulation can be adopted to verify the accuracy of other schemes for IFEA and FFEA.

#### 3. Preprocessing procedures

Except for the aforementioned core algorithms for IFEA/FFEA, two types of preprocessing procedures are always adopted to improve either the accuracy or efficiency.

#### 3.1. Subinterval technique

selected depending on the type of model to be evaluated [23, 24]. If the behavior of the target response w.r.t. uncertain parameters can be guaranteed to be monotonic, the reduced transformation method yields the exact solution. If it shows nonmonotonic behavior, instead, the extended transformation method can be applied, in which more observation points were added in a well-directed way to the search domain after rating the monotonicity of the response w.r.t. different uncertain parameters on the basis of a classification criterion [24].

The computational cost of the transformation method is governed by the number of FEA runs N to be performed. In the case of the general transformation method, this number is given as:

> N ¼ Xmþ1 k¼1 k

where m is the number of discrete α-levels and n is the number of fuzzy parameters. It is noted that the number of FEA runs grows exponentially w.r.t. the number of fuzzy inputs, which makes the general transformation method computational tedious for high-dimensional problems. The main aspect of the transformation method, its characteristic property of reducing fuzzy arithmetic to multiple crisp-number operations entails that this method can be implemented without major problems into an existing software environment for system simulation. Expensive rewriting of the program codes is not required [25]. Some of the most recent applications can be found in Refs. [25–32]. Besides, a program named as FAMOUS (fuzzy arithmetical modeling of uncertain systems) has been developed [25], which provides an interface to commercial software environments. Primarily developed in Matlab environment, FAMOUS actually works as a standalone application on both Windows and Linux platforms.

In essence, calculating the solution set expressed in Eq. (3) is equivalent to performing a global optimization, aimed at the minimization and maximization of the components of the deterministic analysis results {U}. The lower and upper bounds of the output of a classical FEA model are determined by the optimization approach through a search algorithm within the domain spanned by the interval parameters. If the global minimum and maximum of the analysis result are found by the search algorithm, it returns the smallest hypercube solution around the exact one. The optimization is performed independently on each element of the response vector. Furthermore, as the behavior of the target response w.r.t. uncertain parameters is rather unpredictable, the computational cost of the optimization approach in general is strongly problem-dependent. It is noted that the optimization approach is immune to the excessive conservatism for the interval arithmetic approaches because the optimization strategy approaches the smallest hypercube solution from its inside, which means that it does not guarantee conservatism until the actual bounds are captured. Additionally, the smooth behavior of the target response w.r.t. uncertain parameters facilitates the search for the global extrema over the space spanned by uncertain parameters. The directional search-based algorithm [16, 33, 34], linear programming [35], and genetic algorithm [36] were utilized to formulate the procedure of IFEA or FFEA. More applications can be found in [37–39]. It is worth

2.4. Optimization approach

56 Uncertainty Quantification and Model Calibration

n ð23Þ

For the accuracy improvement, the subinterval technique w.r.t. interval inputs is developed [11, 40] and can be integrated with the interval arithmetic and perturbation approaches. The main aspect of the subinterval technique is the ability to relax requirements of "small" or "narrow" interval inputs for nonlinear responses. However, there remain two challenges as follows:

#### 3.2. Surrogate model

To enhance the efficiency of IFEA and FFEA, many surrogate models of the real numerical model are always adopted when dealing with engineering design problems often involving large-scale FEA models. The main aspect of the surrogate model is to avoid the large amount of computational time. Apart from the conventional surrogate models always used in the optimization procedure of IFEA and FFEA, e.g. response surface models [41, 42], Kriging models [43–45], radial basis function models [46–48] and sparse grid meta-models [49–51], those for the sampling and optimization approaches including the high dimensional model representation (HDMR) and the component mode synthesis (CMS) are gaining momentum in recent years. CMS was originally introduced in Ref. [52], in which a Ritz-type transformation to each individual component of a structure was adopted. The deformation of each component is approximated using a limited number of component modes. For each of these vectors, only a single degree of freedom (DOF) was retained in the reduced component model, yielding a large reduction in DOF for each component and the entire structure. Thus, the computational cost for the FEA is drastically reduced. From this viewpoint, CMS can also be seen as a special surrogate model of the expensive numerical FEA for the improvement in the computational efficiency. The repeated FEAs required in the context of IFEA can benefit from this computational time reduction obtained by CMS.

#### 4. Hybrid algorithms

Numerous schemes for IFEA and FFEA have been developed based on the core algorithms and preprocessing procedure, which can be classified into the following three cases.

#### 4.1. Subinterval-based hybrid algorithms

Divide the large interval parameter pI i ði ¼ 1, 2, …, nÞ into Ni subintervals and its rith subinterval can be expressed as follow:

$$(p\_i^{\mathcal{I}})\_{r\_i} = \left[ p\_i^{\mathcal{L}} + \frac{\mathcal{Z}(r\_i - 1)\Delta p\_i}{N\_i}, p\_i^{\mathcal{L}} + \frac{2r\_i \Delta p\_i}{N\_i} \right], \ r\_i = 1, 2, \dots, N\_i \tag{24}$$

The number of subintervals for each interval parameter may be different. Nsub combinations can be produced by taking a subinterval out of each interval parameter.

$$N\_{\text{sub}} = \prod\_{i=1}^{n} N\_i \tag{25}$$

For each subinterval combination, the IFEA model can be rewritten as:

$$\mathbf{K(p\_{r\_1r\_2\ldots r\_n}^1)} \mathbf{U(p\_{r\_1r\_2\ldots r\_n}^1)} = \mathbf{F(p\_{r\_1r\_2\ldots r\_n}^1)}, \ r\_i = 1, 2, \ldots, N\_{i\nu}; \ i = 1, 2, \ldots, n \tag{26}$$

where pI r1r2…rn stands for a subinterval combination and is composed of the r1th subinterval of the first interval parameter, the r2th subinterval of the second one and up to the rnth subinterval of the nth one. In a conclusion, Eq. (26) stands for Nsub subinterval IFEA equations. For each subinterval IFEA equation, the response vector can be obtained by using core algorithms in Section 2, e.g., interval arithmetic approaches, perturbation approaches, and vertex method. For two adjacent subinterval vector pI r1…rr…rn and pI r1…rrþ1…rn , the following formulae hold true, i.e.,

$$\mathbf{K}(\mathbf{p}^{\mathrm{l}}\_{r\_1...r\_r...r\_n}) \cap \mathbf{K}(\mathbf{p}^{\mathrm{l}}\_{r\_1...r\_r+1...r\_n}) = \mathbf{K}(p^{\mathrm{l}}\_{r\_1}...\prime p^{\mathrm{U}}\_{r\_r} = p^{\mathrm{L}}\_{r\_r}...\prime p^{\mathrm{l}}\_{n}) \tag{27}$$

$$\mathbf{F}(\mathbf{p}^{\mathrm{I}}\_{r\_1...r\_r...r\_n}) \cap \mathbf{F}(\mathbf{p}^{\mathrm{I}}\_{r\_1...r\_r+1...r\_n}) = \mathbf{F}(p^{\mathrm{I}}\_{r\_1}......p^{\mathrm{U}}\_{r\_r} = p^{\mathrm{L}}\_{r\_r}......p^{\mathrm{I}}\_{n}) \tag{28}$$

where pU rr and pL rr are the upper bound of pI rr and lower bound of pL rrþ1, respectively. Thus, the intersection of UðpI r1…rr…rn Þ and UðpI r1…rrþ1…rn Þ does not equal to an empty set, i.e.,

$$\mathbf{U}(\mathbf{p}^{\mathrm{I}}\_{r\_{1}\ldots r\_{r}\ldots r\_{n}}) \cap \mathbf{U}(\mathbf{p}^{\mathrm{I}}\_{r\_{1}\ldots r\_{t}+1\ldots r\_{n}}) \neq \mathfrak{Q} \tag{29}$$

It is shown from Eq. (29) that the interval response vectors for each subinterval combination are simply connected. Therefore, the interval response vector can be obtained as follows by using the interval union operation.

$$\begin{split} \mathbf{U}(\mathbf{p}^{\mathrm{I}}) &= \bigcup\_{r\_{1}=1,2,\ldots,N\_{i}} \mathbf{U}(\mathbf{p}^{\mathrm{I}}\_{r\_{1}r\_{2}\ldots r\_{i}\ldots r\_{n}}) = \\ & \left[ \min\_{r\_{i}=1,2,\ldots,N\_{i}} \left( \mathbf{U}(\mathbf{p}^{\mathrm{I}}\_{r\_{1}r\_{2}\ldots r\_{i}\ldots r\_{n}}) \right) \max\_{r\_{i}=1,2,\ldots,N\_{i}} \left( \mathbf{U}(\mathbf{p}^{\mathrm{I}}\_{r\_{1}r\_{2}\ldots r\_{i}\ldots r\_{n}}) \right) \right] \end{split} \tag{30}$$

The above subinterval method is shown in Figure 2 with 50 subintervals when considering one uncertain parameter x.

The interval arithmetic approach, subinterval technique and Taylor series expansion were integrated [40]. More applications can be found in [13, 53, 54].

#### 4.2. Surrogate model-based hybrid algorithms

efficiency. The repeated FEAs required in the context of IFEA can benefit from this computa-

Numerous schemes for IFEA and FFEA have been developed based on the core algorithms

, pL i þ

The number of subintervals for each interval parameter may be different. Nsub combinations

Nsub ¼ Yn

r1r2…rn

r1…rr…rn and pI

r1…rrþ1…rn

r1…rrþ1…rn

r1…rrþ1…rn

Þ ∩ UðpI

r1…rr…rn

the first interval parameter, the r2th subinterval of the second one and up to the rnth subinterval of the nth one. In a conclusion, Eq. (26) stands for Nsub subinterval IFEA equations. For each subinterval IFEA equation, the response vector can be obtained by using core algorithms in Section 2, e.g., interval arithmetic approaches, perturbation approaches, and vertex method.

i¼1

r1r2…rn stands for a subinterval combination and is composed of the r1th subinterval of

Þ ¼ KðpI r1 , …, pU

Þ ¼ FðpI r1 , …, pU

r1…rrþ1…rn

rr and lower bound of pL

r1…rrþ1…rn

rr ¼ pL rr ,…, pI

rr ¼ pL rr , …, pI

Þ does not equal to an empty set, i.e.,

2riΔpi Ni

ði ¼ 1, 2, …, nÞ into Ni subintervals and its rith subinter-

, ri ¼ 1, 2, …, Ni ð24Þ

Ni ð25Þ

Þ, ri ¼ 1, 2, …, Ni; i ¼ 1, 2, …, n ð26Þ

, the following formulae hold true,

nÞ ð27Þ

nÞ ð28Þ

rrþ1, respectively. Thus, the

Þ 6¼ ∅ ð29Þ

and preprocessing procedure, which can be classified into the following three cases.

i

2ðri � 1ÞΔpi Ni

can be produced by taking a subinterval out of each interval parameter.

For each subinterval combination, the IFEA model can be rewritten as:

Þ ¼FðpI

r1r2…rn

Þ ∩ KðpI

Þ ∩ FðpI

ÞUðpI

� �

i þ

tional time reduction obtained by CMS.

58 Uncertainty Quantification and Model Calibration

4.1. Subinterval-based hybrid algorithms

ðpI i Þri ¼ pL

KðpI r1r2…rn

For two adjacent subinterval vector pI

KðpI

FðpI

r1…rr…rn

r1…rr…rn

r1…rr…rn

rr are the upper bound of pI

Þ and UðpI

UðpI

where pI

i.e.,

where pU

rr and pL

intersection of UðpI

Divide the large interval parameter pI

val can be expressed as follow:

4. Hybrid algorithms

Taylor series expansion was integrated with the interval arithmetic approach in [40] and a method named as Taylor expansion with extrema management was proposed by integrating the higher order Taylor series expansion and the optimization approach [55] to detect possible nonmonotonic influences.

The transformation method was integrated with HDMR in Ref. [25]. And a component mode transformation method was developed [56] by combing the CMS with the transformation method to provide a significant reduction of the computational cost for large mechanical problems with uncertain parameters. Besides, a hybrid method was proposed for the interval frequency response analysis by integrating the optimization and interval arithmetic approach in [57], which was further integrated with CMS in Ref. [22]. An acceptable computational cost and a limited amount of conservatism in the analysis result were achieved by these hybrid algorithms.

Figure 2. The diagram of subinterval method.

#### 4.3. Hybrid core algorithms

The aforementioned core algorithms can be combined together to achieve a better tradeoff between the accuracy and efficiency, e.g., frameworks [22, 57–60] formulated by the global optimization methods and interval arithmetic ones.

To improve the computational efficiency, any core algorithm in Section 2 can be integrated with reanalysis method [61], which is fundamentally an intrusive FEA. It is noted that the major computational cost of IFEA consists of repeated solutions of the deterministic FEA systems while the main goal of the re-analysis method is to accelerate this conventional FEA solution process. It is shown that the application of the re-analysis method in the context of IFEA can reduce the computational cost by one order of magnitude compared to those based on the conventional FEA strategy [5].

#### 5. Conclusions

This chapter presents the state-of-the-art and recent advances in nonprobabilistic finite element analyses. The main advantages and shortcomings of each nonprobabilistic finite element analysis method are discussed.

The arithmetic approach is the most straightforward strategy for nonprobabilistic finite element analyses. However, this chapter further shows that the interval arithmetic implementation of the finite element procedure is conservative. Therefore, the development of an adequate methodology for solving the uncertain parameter dependency problem is still the main challenge in the domain of arithmetic approach. The perturbation approach has been widely used in structural response analyses and other applications due to its simplicity and efficiency. The accuracy of the original perturbation methods can be improved by retaining part of higher order terms in Neumann series or Taylor series as shown in the improved and modified versions. The sampling approach like vertex method yields the exact solution under the condition that the behavior of the target response w.r.t. uncertain parameters can be guaranteed to be monotonic and has been popular for the implementation of IFEA due to its main aspect of simple formulation and the black-box property. However, when tackling the nonmonotonic problems, the extended transformation methods should be applied by adding more observation points in a well-directed way. The optimization approach is more and more acknowledged as standard procedure in an interval finite element context except for the high computational cost.

Moreover, in this context, two typical preprocessing procedures, e.g., subinterval technique and surrogate model to improve either the accuracy or efficiency are described in detail. Additionally, novel hybrid algorithms, including subinterval-based hybrid algorithms, surrogate model-based hybrid algorithms and hybrid core algorithms can be conceived by combining the aforementioned core algorithms and preprocessing procedures to achieve a better tradeoff between the accuracy and efficiency for the specific problems and the future work in this research field can be fostered.

#### Author details

4.3. Hybrid core algorithms

60 Uncertainty Quantification and Model Calibration

on the conventional FEA strategy [5].

analysis method are discussed.

5. Conclusions

computational cost.

this research field can be fostered.

optimization methods and interval arithmetic ones.

The aforementioned core algorithms can be combined together to achieve a better tradeoff between the accuracy and efficiency, e.g., frameworks [22, 57–60] formulated by the global

To improve the computational efficiency, any core algorithm in Section 2 can be integrated with reanalysis method [61], which is fundamentally an intrusive FEA. It is noted that the major computational cost of IFEA consists of repeated solutions of the deterministic FEA systems while the main goal of the re-analysis method is to accelerate this conventional FEA solution process. It is shown that the application of the re-analysis method in the context of IFEA can reduce the computational cost by one order of magnitude compared to those based

This chapter presents the state-of-the-art and recent advances in nonprobabilistic finite element analyses. The main advantages and shortcomings of each nonprobabilistic finite element

The arithmetic approach is the most straightforward strategy for nonprobabilistic finite element analyses. However, this chapter further shows that the interval arithmetic implementation of the finite element procedure is conservative. Therefore, the development of an adequate methodology for solving the uncertain parameter dependency problem is still the main challenge in the domain of arithmetic approach. The perturbation approach has been widely used in structural response analyses and other applications due to its simplicity and efficiency. The accuracy of the original perturbation methods can be improved by retaining part of higher order terms in Neumann series or Taylor series as shown in the improved and modified versions. The sampling approach like vertex method yields the exact solution under the condition that the behavior of the target response w.r.t. uncertain parameters can be guaranteed to be monotonic and has been popular for the implementation of IFEA due to its main aspect of simple formulation and the black-box property. However, when tackling the nonmonotonic problems, the extended transformation methods should be applied by adding more observation points in a well-directed way. The optimization approach is more and more acknowledged as standard procedure in an interval finite element context except for the high

Moreover, in this context, two typical preprocessing procedures, e.g., subinterval technique and surrogate model to improve either the accuracy or efficiency are described in detail. Additionally, novel hybrid algorithms, including subinterval-based hybrid algorithms, surrogate model-based hybrid algorithms and hybrid core algorithms can be conceived by combining the aforementioned core algorithms and preprocessing procedures to achieve a better tradeoff between the accuracy and efficiency for the specific problems and the future work in Wang Lei, Qiu Zhiping\* and Zheng Yuning

\*Address all correspondence to<EMAIL_ADDRESS>
Institute of Solid Mechanics, School of Aeronautic Science and Engineering, Beihang University, Beijing, China

#### References

[24] Hanss M. The extended transformation method for the simulation and analysis of fuzzyparameterized models. International Journal of Uncertainty, Fuzziness and Knowledge-Based Systems. 2003;11(06):711–727. doi:doi:10.1142/S0218488503002491

[11] Qiu Z, Elishakoff I. Antioptimization of structures with large uncertain-but-non-random parameters via interval analysis. Computer Methods in Applied Mechanics and Engineering. 1998;152(3):361–372. doi:http://dx.doi.org/10.1016/S0045-7825(96)01211-X [12] McWilliam S. Anti-optimisation of uncertain structures using interval analysis. Computers & Structures. 2001;79(4):421–430. doi:http://dx.doi.org/10.1016/S0045-7949(00)00143-7 [13] Xia B, Yu D, Liu J. Interval and subinterval perturbation methods for a structural-acoustic system with interval parameters. Journal of Fluids and Structures. 2013;38:146–163. doi:

[14] Wang C, Qiu Z, Wang X, Wu D. Interval finite element analysis and reliability-based optimization of coupled structural-acoustic system with uncertain parameters. Finite Elements in Analysis and Design. 2014;91:108–114. doi:http://dx.doi.org/10.1016/j.finel.2014.07.014 [15] Dong W, Shah HC. Vertex method for computing functions of fuzzy variables. Fuzzy Sets and Systems. 1987;24(1):65–78. doi:http://dx.doi.org/10.1016/0165-0114(87)90114-X [16] Rao S, Sawyer JP. Fuzzy finite element approach for analysis of imprecisely defined

[17] Chen L, Rao SS. Fuzzy finite-element approach for the vibration analysis of impreciselydefined systems. Finite Elements in Analysis and Design. 1997;27(1):69–83. doi:http://dx.

[18] Akpan UO, Koko TS, Orisamolu IR, Gallant BK. Practical fuzzy finite element analysis of structures. Finite Elements in Analysis and Design. 2001;38(2):93–111. doi:http://dx.doi.

[19] Qiu Z, Wang X, Chen J. Exact bounds for the static response set of structures with uncertain-but-bounded parameters. International Journal of Solids and Structures.

[20] Qiu Z, Xia Y, Yang J. The static displacement and the stress analysis of structures with bounded uncertainties using the vertex solution theorem. Computer Methods in Applied Mechanics and Engineering. 2007;196(49–52):4965–4984. doi:10.1016/j.cma.2007.06.022

[21] Xu M, Qiu Z, Wang X. Uncertainty propagation in SEA for structural–acoustic coupled systems with non-deterministic parameters. Journal of Sound and Vibration. 2014;333

[22] De Gersem H, Moens D, Desmet W, Vandepitte D. Interval and fuzzy dynamic analysis of finite element models with superelements. Computers & Structures. 2007;85(5–6):304–

[23] Hanss M. The transformation method for the simulation and analysis of systems with uncertain parameters. Fuzzy Sets and Systems. 2002;130(3):277–289. doi:http://dx.doi.

2006;43(21):6574–6593. doi:http://dx.doi.org/10.1016/j.ijsolstr.2006.01.012

(17):3949–3965. doi:http://dx.doi.org/10.1016/j.jsv.2014.03.003

319. doi:http://dx.doi.org/10.1016/j.compstruc.2006.10.011

systems. AIAA Journal. 1995;33(12):2364–2370. doi:10.2514/3.12910

http://dx.doi.org/10.1016/j.jfluidstructs.2012.12.003

62 Uncertainty Quantification and Model Calibration

doi.org/10.1016/S0168-874X(97)00005-X

org/10.1016/S0168-874X(01)00052-X

org/10.1016/S0165-0114(02)00045-3

[50] Klimke A, Willner K, Wohlmuth BI. Uncertainty modeling using fuzzy arithmetic based on sparse grids: Applications to dynamic systems. International Journal of Uncertainty, Fuzziness and Knowledge-Based Systems. 2004;12(6):745–759. doi:10.1142/S0218488504003181

[36] Möller B, Graf W, Beer M. Fuzzy structural analysis using α-level optimization. Compu-

[37] Moens D, Vandepitte D. Fuzzy finite element method for frequency response function analysis of uncertain structures. AIAA Journal. 2002;40(1):126–136. doi:10.2514/2.1621 [38] Farkas L, Moens D, Vandepitte D, Desmet W. Application of fuzzy numerical techniques for product performance analysis in the conceptual and preliminary design stage. Com-

puters & Structures. 2008;86(10):1061–1079. doi:10.1016/j.compstruc.2007.07.012

[39] Farkas L, Moens D, Vandepitte D, Desmet W. Fuzzy finite element analysis based on reanalysis technique. Structural Safety. 2010;32(6):442–448. doi:10.1016/j.strusafe.2010.04.004

[40] Zhou YT, Jiang C, Han X. Interval and subinterval analysis methods of the structural analysis and their error estimations. International Journal of Computational Methods.

[41] de Boor C, Ron A. On multivariate polynomial interpolation. Constructive Approxima-

[42] Myers RH, Montgomery DC, Anderson-Cook CM. Response Surface Methodology: Process and Product Optimization Using Designed Experiments. Hoboken, NJ: John Wiley &

[43] Jones DR, Schonlau M, Welch WJ. Efficient global optimization of expensive black-box functions. Journal of Global Optimization. 1998;13(4):455–492. doi:10.1023/a:1008306431147

[44] Martin JD, Simpson TW. Use of Kriging models to approximate deterministic computer

[45] Kleijnen JPC. Kriging metamodeling in simulation: A review. European Journal of Operational Research. 2009;192(3):707–716. doi:http://dx.doi.org/10.1016/j.ejor.2007.10.013 [46] Park J, Sandberg IW. Universal approximation using radial-basis-function networks.

[47] Chen S, Cowan CFN, Grant PM. Orthogonal least squares learning algorithm for radial basis function networks. IEEE Transactions on Neural Networks. 1991;2(2):302–309.

[48] T Sev, Shin YC. Radial basis function neural network for approximation and estimation of nonlinear stochastic dynamic systems. IEEE Transactions on Neural Networks. 1994;5

[49] Klimke A, Nunes RF, Wohlmuth BI. Fuzzy arithmetic based on dimension-adaptive sparse grids: A case study of a large-scale finite element model under uncertain parameters. International Journal of Uncertainty, Fuzziness and Knowledge-Based Systems.

models. AIAA Journal. 2005;43(4):853–863. doi:10.2514/1.8650

Neural Computation. 1991;3(2):246–257. doi:10.1162/neco.19<IP_ADDRESS>

tational Mechanics. 2000;26(6):547–565. doi:10.1007/s004660000204

2006;3(2):229–244. doi:10.1142/S0219876206000771

tion. 1990;6(3):287–302. doi:10.1007/bf01890412

Sons; 2011.

64 Uncertainty Quantification and Model Calibration

doi:10.1109/72.80341

(4):594–603. doi:10.1109/72.298229

2006;14(5):561–577. doi:10.1142/S0218488506004199

### **Epistemic Uncertainty Quantification of Seismic Damage Assessment**

Hesheng Tang, Dawei Li and Songtao Xue

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/intechopen.68153

#### Abstract

The damage-based structural seismic performance evaluations are widely used in seismic design and risk evaluation of civil facilities. Due to the large uncertainties rooted in this procedure, the application of damage quantification results is still a challenge for researchers and engineers. Uncertainties in damage assessment procedure are important consideration in performance evaluation and design of structures against earthquakes. Due to lack of knowledge or incomplete, inaccurate, unclear information in the modeling, simulation, and design, there are limitations in using only one framework (probability theory) to quantify uncertainty in a system because of the impreciseness of data or knowledge. In this work, a methodology based on the evidence theory is presented for quantifying the epistemic uncertainty of damage assessment procedure. The proposed methodology is applied to seismic damage assessment procedure while considering various sources of uncertainty emanating from experimental force-displacement data of reinforced concrete column. In order to alleviate the computational difficulties in the evidence theory-based uncertainty quantification analysis (UQ), a differential evolutionbased computational strategy for efficient calculation of the propagated belief structure in a system with evidence theory is presented here. Finally, a seismic damage assessment example is investigated to demonstrate the effectiveness of the proposed method.

Keywords: damage model, epistemic uncertainty, uncertainty quantification, evidence theory, differential evolution algorithm

#### 1. Introduction

With widespreading of the concept and applications of performance-based earthquake engineering (PBEE) and performance-based seismic design (PBSD), the effective measures for assessing the performance state of structural components or entire structure have been deeply

investigated in seismic engineering. In consistence with the different performance assessment criteria, the evaluation and measurements of damage states for structural components are divided into three main branches (e.g., displacement-based approach, energy-based measure and the combination of both). Due to the simplicity and convenience of observation and description for structural damage states, the displacement-based approach and corresponding damage index (e.g., inelastic displacement, maximum inter story drift ratio, and ductility demand, etc.) have been widely documented in building seismic evaluation and retrofit of existing building guidelines [1]. Notwithstanding the prevalent application of displacement method in damage assessment, the defect of lacking the influence of low cyclic fatigue of structural components is obvious. The hysteretic energy dissipation is considered as a more reasonable indicator for seismic structural damage, because it is a cumulative parameter involved cyclic-plastic deformations in a structure during earthquakes [2]. Despite the effectiveness of hysteretic energy, experimental observations demonstrate that the expression of energy would be significantly affected by the exceedance plastic deformation [3]. And the cumulative laboratory experimental data on structural members and structures indicate the fact that the structure is damaged by a combination of the excessive deformation and hysteretic energy. Park–Ang damage model [4], which takes into account the effects of both the first exceedance failure and cumulative damage failure in low-cycle-fatigue for a structural component during seismic load, is served as a baseline for many researches. Due to intrinsic simplicity as well as calibrations against a significant amount of observed seismic damages, the Park-Ang model and its modified version have been extensively implemented in seismic performance evaluation of structures [5–7].

Although the applicability and practicability of using the Park-Ang model and its modified versions have been supported by many researchers [8, 9], it should be noted that the Park-Ang-damage-index-based performance evaluation is still a challenging task due to the large uncertainties associated with the damage model parameters [10]. With the influence of these uncertainties [11, 12], the evaluation results of structural damage state are always represented with the empirical interval value (e.g., the minor damage state is represented by 0.25<D<0.4 or 0.11<D<0.4, etc. [13]). Some of these uncertainties stem from factors that are inherently random (or aleatory) in engineering or scientific analysis (e.g., material properties such as Young's modulus of steel; compression strength of concrete). Others arise from a lack of knowledge, ignorance, or modeling (e.g., simplification of mathematical model of buildings for structural analysis purposes). The large uncertainties associated with the Park-Ang damage model are derived from limited experimental data and approximate modeling (lack of knowledge) [2, 4, 5, 10]. Considering the importance of damage model in assessment of damage state for a structure, the epistemic uncertainty shall be taken into account in seismic damage state assessment with great care. Hence, it is significant to present a comprehensive uncertainty analysis methodology to quantify the epistemic uncertainty and obtain more reliable results.

The traditional probability theory, based on the sufficient statistical information, is used to model the objective uncertainty (random), which is inherent in physical variability of materials and environment. Unfortunately, the limited number of experimental data set cannot support the strong assumption of probability theory, and the process of collecting data is always costly and time consuming. These shortcomings lead the assessment result of damage state of structures are not aleatory but epistemic. In the past decades, several alternative approaches have been developed to deal with epistemic uncertainty. Some of the potential uncertainty theories are the theory of fuzzy set [14], possibility theory [15], the theory of interval analysis [16], imprecise probability theory [17], and evidence theory [18, 19]. Among these promising uncertainty representation models, evidence theory with the ability of handling aleatory and epistemic uncertainty is used for UQ, risk assessment, and reliability analysis.

With two complementary measures of uncertainty such as belief and plausibility, using evidence theory to UQ is flexible and effective. In comparison with the calculation of single probability density function (PDF) in probability theory, the computationally intensive problem involves computing the bound values over all possible discontinuous sets which is a main shackle of wide application for evidence theory. In order to break the computational barriers in the evidence theory-based UQ, the differential-evolution-based interval optimization is employed to enhance the computational efficiency as described by the authors [20].

#### 2. Sources of uncertainty in seismic damage assessment

investigated in seismic engineering. In consistence with the different performance assessment criteria, the evaluation and measurements of damage states for structural components are divided into three main branches (e.g., displacement-based approach, energy-based measure and the combination of both). Due to the simplicity and convenience of observation and description for structural damage states, the displacement-based approach and corresponding damage index (e.g., inelastic displacement, maximum inter story drift ratio, and ductility demand, etc.) have been widely documented in building seismic evaluation and retrofit of existing building guidelines [1]. Notwithstanding the prevalent application of displacement method in damage assessment, the defect of lacking the influence of low cyclic fatigue of structural components is obvious. The hysteretic energy dissipation is considered as a more reasonable indicator for seismic structural damage, because it is a cumulative parameter involved cyclic-plastic deformations in a structure during earthquakes [2]. Despite the effectiveness of hysteretic energy, experimental observations demonstrate that the expression of energy would be significantly affected by the exceedance plastic deformation [3]. And the cumulative laboratory experimental data on structural members and structures indicate the fact that the structure is damaged by a combination of the excessive deformation and hysteretic energy. Park–Ang damage model [4], which takes into account the effects of both the first exceedance failure and cumulative damage failure in low-cycle-fatigue for a structural component during seismic load, is served as a baseline for many researches. Due to intrinsic simplicity as well as calibrations against a significant amount of observed seismic damages, the Park-Ang model and its modified version have been extensively implemented in seismic performance evaluation of struc-

Although the applicability and practicability of using the Park-Ang model and its modified versions have been supported by many researchers [8, 9], it should be noted that the Park-Ang-damage-index-based performance evaluation is still a challenging task due to the large uncertainties associated with the damage model parameters [10]. With the influence of these uncertainties [11, 12], the evaluation results of structural damage state are always represented with the empirical interval value (e.g., the minor damage state is represented by 0.25<D<0.4 or 0.11<D<0.4, etc. [13]). Some of these uncertainties stem from factors that are inherently random (or aleatory) in engineering or scientific analysis (e.g., material properties such as Young's modulus of steel; compression strength of concrete). Others arise from a lack of knowledge, ignorance, or modeling (e.g., simplification of mathematical model of buildings for structural analysis purposes). The large uncertainties associated with the Park-Ang damage model are derived from limited experimental data and approximate modeling (lack of knowledge) [2, 4, 5, 10]. Considering the importance of damage model in assessment of damage state for a structure, the epistemic uncertainty shall be taken into account in seismic damage state assessment with great care. Hence, it is significant to present a comprehensive uncertainty analysis methodology to quantify the

The traditional probability theory, based on the sufficient statistical information, is used to model the objective uncertainty (random), which is inherent in physical variability of materials and environment. Unfortunately, the limited number of experimental data set cannot support the strong assumption of probability theory, and the process of collecting data is

epistemic uncertainty and obtain more reliable results.

tures [5–7].

68 Uncertainty Quantification and Model Calibration

To effectively describe the damage state of structural components or entire structure, the original Park-Ang damage model and modified model were developed. The original Park-Ang damage model was presented here to access the uncertainty influence of the evaluation on the damage state of column components. There are various methods to estimate constants in Park-Ang damage model in different studies. In addition to diverse combination measures, the empirical estimation value and calibration value dispersed in a large range. Using the classification method proposed by Oberkampf and Helton [21], the aleatory and epistemic uncertainties involved in Park-Ang damage model are listed as:

In consideration of these aleatory and epistemic uncertainties in Park-Ang damage model, the quantification influence of uncertainties is indispensable. To achieve this goal, a series of empirical expressions are summarized. Then, the Structural Performance Database of Pacific Earthquake Engineering Research Center (PEER) is used to construct the uncertain sources of parameters of damage assessment models. Using these calibration results of column set, the parameter uncertainties are represented by the fluctuation of ratio of empirical values and calibration values.

#### 2.1. Park-Ang model and empirical expression of its constants

The Park-Ang damage model [4] combines the first exceedance failure and cumulative damage failure with a linear expression as:

$$D = \delta\_{\mathbf{m}} / \delta\_{\mathbf{u}} + \beta \left[ dE / F\_{\mathbf{y}} \delta\_{\mathbf{u}} \right] \tag{1}$$

where δm is the maximum deformation under earthquake, δu is the ultimate deformation under monotonic load, ð dE is the cumulative energy under earthquake, β is the energy coefficient, and Fy is the yield strength. In order to simplify the analysis procedure, the value of Fy, δu, and β are always assumed as the constants and have nothing to do with the loads pattern. Following above assumption, the value of damage index D for per-load stage can be computed by only using the current value of δm and ð dE. Furthermore, the damage evolution of structures and components can be described and this evolution index is supported to estimate the true damage stage of structure and components.

In the last two decades of the twentieth century, a set of experimental results were conducted and some illuminate-, empirical-, or mechanical-based expression of Fy, δu, and β were successively generated. Park et al. [4] computed the value of β as given in Eq. (2):

$$\beta = (-0.447 + 0.073l/d + 0.24n\_0 + 0.314\rho\_t) \times 0.7^{\rho\_w} \tag{2}$$

where l and d denote the length span and effective height of cross section, n0 is the axial load ratio, ρt is the longitude tension steel ratio (%), and ρw is the confinement ratio (%). Kunnath et al. [5] used 260 beams and columns data to fit the value of β as given in Eq. (3):

$$\beta = [0.37n\_0 + 0.36(k\_\mathrm{p} - 0.2)^2]0.9^{v\_\mathrm{w}}\tag{3}$$

where kp ¼ ρt f y=0:85f c is normalized steel ratio and ρw is confinement ratio. Similarly, δu can be determined with statistical approach or fundamental method using the mechanics of concrete and steel. Using the typical statistical measure, Park [6] evaluated the ultimate displacement as:

$$
\delta\_{\rm u} = 0.52(l/d)^{0.93} \rho^{-0.27} \rho\_{\rm av}^{0.48} n\_0^{-0.48} f\_{\rm c}^{-0.15} \times \delta\_{\rm y} \tag{4}
$$

where ρ is normalized steel ratio and δy is the yield displacement of components that can be computed with [4] and other factors are same as above. Compared to above statistical calculation model, EU 8 [22] and Fardis and Biskinis [23] presented two different models with the mechanics of concrete and steel:

$$\delta\_{\rm u} = \frac{1}{\mathcal{V}\_{\rm el}} 0.016(0.3)^{\rm n\_0} \left[ \frac{\max(0.01, \omega')}{\max(0.01, \omega)} f\_{\rm c} \right]^{0.225} \left[ \min \left( 9, \frac{l}{h} \right) \right]^{0.35} 25 \left( ^{a\rho\_m \frac{l\_{\rm v}}{f\_{\rm c}}} \right) \times l \tag{5}$$

$$\delta\_{\rm u} = \alpha\_{\rm st} (1 - 0.4a\_{\rm yce}) (1 + 0.5\alpha\_{\rm sl}) (0.3)^{n\_l} \left[ \frac{\max(0.01, \alpha')}{\max(0.01, \alpha)} f\_c \right]^{0.175} \left( \frac{l}{h} \right)^{0.4} 25 \left( ^{a\rho\_{\rm s}^{f\_{\rm yw}}} \right) \times l \tag{6}$$

where γel is coefficient of primary and secondary elements, ω0 and ω are mechanical steel ratio of compression and tension reinforcement, respectively, h is cross-section height, α is confinement effective factor, ρsx is confinement steel ratio, fyw is yield strength of stirrup, and αst, αcyc, and αsl are coefficients for type of steel, loading, and anchorage slip. For the yield strength of concrete components, the expression is given by Panagiotakos and Fardis [24]:

$$F\_{\mathbf{Y}} = \frac{bd^3}{l} \phi\_{\mathbf{y}} \left\{ E\_{\mathbf{c}} \frac{k\_{\mathbf{y}}^2}{2} \left( 0.5(1+\delta') - \frac{k\_{\mathbf{y}}}{3} \right) + \frac{E\_{\mathbf{s}}}{2} \left[ (1-k\_{\mathbf{y}})\rho + (k\_{\mathbf{y}}-\delta')\rho' + \frac{\rho\_{\mathbf{v}}}{6}(1-\delta') \right] \right\} (1-\delta') \quad (7)$$

Conventionally, the damage index D can be obtained by using above expressions to obtain the nominal value of Park-Ang constants. Owing to limited statistical data and incomplete knowledge of mathematical model to predict these constants, the large convergence is reported as in [4–6, 23, 24]. Furthermore, these uncertainties will influence the quantification result of Park-Ang damage index. In order to verify the impact of damage quantification result derived from uncertainty of Park-Ang model constants, we present the structural performance database of PEER [25] to calibrate these constants and determine the uncertainty fluctuation range of each constant.

#### 2.2. Comparison between the calibration results and empirical results

2.1. Park-Ang model and empirical expression of its constants

ð

computed by only using the current value of δm and

estimate the true damage stage of structure and components.

sively generated. Park et al. [4] computed the value of β as given in Eq. (2):

failure with a linear expression as:

70 Uncertainty Quantification and Model Calibration

under monotonic load,

The Park-Ang damage model [4] combines the first exceedance failure and cumulative damage

where δm is the maximum deformation under earthquake, δu is the ultimate deformation

coefficient, and Fy is the yield strength. In order to simplify the analysis procedure, the value of Fy, δu, and β are always assumed as the constants and have nothing to do with the loads pattern. Following above assumption, the value of damage index D for per-load stage can be

of structures and components can be described and this evolution index is supported to

In the last two decades of the twentieth century, a set of experimental results were conducted and some illuminate-, empirical-, or mechanical-based expression of Fy, δu, and β were succes-

where l and d denote the length span and effective height of cross section, n0 is the axial load ratio, ρt is the longitude tension steel ratio (%), and ρw is the confinement ratio (%). Kunnath

where kp ¼ ρt f y=0:85f c is normalized steel ratio and ρw is confinement ratio. Similarly, δu can be determined with statistical approach or fundamental method using the mechanics of concrete and steel. Using the typical statistical measure, Park [6] evaluated the ultimate displacement as:

0:93ρ�0:27ρ0:48

where ρ is normalized steel ratio and δy is the yield displacement of components that can be computed with [4] and other factors are same as above. Compared to above statistical calculation model, EU 8 [22] and Fardis and Biskinis [23] presented two different models with the

Þ

� �0:225

f c

n0 maxð0:01,ω0

maxð0:01,ωÞ

2

ω n�0:48 0 f �0:15

min 9,

Þ

f c � �0:175 l

l h � � � � 0:35

β ¼ ð�0:447 þ 0:073l=d þ 0:24n0 þ 0:314ρt

et al. [5] used 260 beams and columns data to fit the value of β as given in Eq. (3):

δu ¼ 0:52ðl=dÞ

n0 maxð0:01,ω0

maxð0:01,ωÞ

mechanics of concrete and steel:

δu ¼ 1 γel

0:016ð0:3Þ

δu ¼ αstð1 � 0:4αcycÞð1 þ 0:5αslÞð0:3Þ

β ¼ ½0:37n0 þ 0:36ðkp � 0:2Þ

ð

ð

dE is the cumulative energy under earthquake, β is the energy

dE=Fyδu ð1Þ

dE. Furthermore, the damage evolution

Þ · 0:7ρω ð2Þ

�0:9pw ð3Þ

c · δy ð4Þ

f yw f c � �

> 25 αρs f yw f c � �

· l ð5Þ

· l ð6Þ

25 αρsx

h � �0:4

D ¼ δm=δu þ β

In this work, the calibration set is selected from the structural performance database of PEER and the selection criteria are such as (1) the cross section of column is rectangle; (2) the column is loaded cyclically until failure and the corresponding failure model is dominated by flexure; (3) the longitude bars in column should not be spliced and the column should experience more than two hysteretic cycles. In conformity with these criteria, 185 specimens are selected. Using these column load-displacement data, the performance points on the backbone curve of column under cyclic load are calibrated.

Similar to the most studies [23], the ultimate deformation under monotonic load δu is defined as a distinct reduction on the negative stiffness slope of backbone curve and 80% of maximum strength which is always assumed as Fu. Unfortunately, the missing monotonic load experiments oblige us to employ the statistical relationship of ultimate displacement under cyclic load and monotonic load to characterize the ultimate displacement. Herein, the failure displacement under typical load histories is assumed as 60% of their ultimate deformation capacity, which is firstly observed by Panagiotakos and Fardis [24]. For yield force, we defined that the value is 75% of the maximum force. Following above definitions, the energy coefficient β is computed with the assumption that damage index D is 1 at the ultimate state. In light with above definitions, the performance point is marked on the backbone of columns as depicted in Figure 1.

As shown in Figure 1, the column backbone curves are divided into two categories: one with obvious ultimate state point (the 80% maximum force) like in Figure 1a, the other with the largest displacement in backbone curve (e.g., Figure 1b). In order to yield the uncertainty distribution of empirical model, the attention is concentrated on the first category. Using the selected force-displacement data, the comparison of empirical model results and calibration results is given in Figure 2.

Figure 1. Performance point of backbone curve with obvious ultimate state point (a) and with the largest displacement point (b).

Figure 2. Comparison of predicted results and experimental results of β, δu and Fy.

As shown in Figure 2, the predicted and experimental values are scattered in a wide range, and this means the researchers should carefully handle the uncertainty derived from the empirical model in the process of evaluating damage state with Park-Ang model. Employing the parameter ε to represent the variability of predicted model deviation, the experimental value Vexp can be expressed as Vexp = Vpre · ε. Taking into account the major fluctuation range of ε and the number of experimental samples, the ε which is located in the interval [1/3, 3] is selected and the range of data points which located less than 1/3 or more than 3 are discarded. In the light of above rules, the uncertainty source of β, δu, and Fy consists of 83, 111, and 173 specimen, respectively.

Along with classical concept, probability theory plays a key role in the UQ of physical model, and the distribution type is determined by the hypothesis test and related parameter are calibrated by enough experimental data. However, the limited data of experimental set and large variation restricted the ability of probability theory. As a generalized UQ measure, evidence theory is compatible with both aleatory and epistemic uncertainties. So, the evidence theory is adopted in this work to handle the epistemic uncertainty rooted in parameters of Park-Ang damage model.

#### 3. Evidence theory and differential evolution-based UQ for seismic damage assessment

#### 3.1. Basic of the evidence theory

As shown in Figure 2, the predicted and experimental values are scattered in a wide range, and this means the researchers should carefully handle the uncertainty derived from the empirical model in the process of evaluating damage state with Park-Ang model. Employing the parameter ε to represent the variability of predicted model deviation, the experimental value Vexp can be expressed as Vexp = Vpre · ε. Taking into account the major fluctuation range of ε and the number of experimental samples, the ε which is located in the interval [1/3, 3] is selected and the range of data points which located less than 1/3 or more than 3 are discarded. In the light of above rules, the uncertainty source of β, δu, and Fy consists of 83, 111, and 173

Figure 2. Comparison of predicted results and experimental results of β, δu and Fy.

Figure 1. Performance point of backbone curve with obvious ultimate state point (a) and with the largest displacement

Along with classical concept, probability theory plays a key role in the UQ of physical model, and the distribution type is determined by the hypothesis test and related parameter are calibrated by enough experimental data. However, the limited data of experimental set and large variation restricted the ability of probability theory. As a generalized UQ measure, evidence theory is compatible with both aleatory and epistemic uncertainties. So, the evidence theory is adopted in this work to handle the epistemic uncertainty rooted in parameters of

specimen, respectively.

point (b).

72 Uncertainty Quantification and Model Calibration

Park-Ang damage model.

Evidence theory is a theoretical framework for reasoning with partial and unreliable information. It was proposed by Dempster [18] and further improved by Shafer [19]. Compared to the classical uncertain model theory, it offers the possibility to explicitly represent doubt and conflict. As the most basic concept of the evidence theory, the fame of discernment Ω is defined as a set of mutually exclusive elementary propositions. Due to limited information, the propositions can be scattered, nested, or partially overlapped. Thus, the mutually exclusive elementary propositions construct the power set F = 2Ω. Given the measureable sample space (Ω, F), the basic belief assignment (BBA) on F, m is a mapping F! [0, 1] that satisfies the following axioms:

$$m(A) \ge 0 \quad m(\mathfrak{Q}) = 0 \quad \sum m(A) = 1 \quad \text{for each } A \subseteq \Omega \,. \tag{8}$$

An element A∈F for which m(A) > 0 is named a focal element. Corresponding to the scattered, nested, or partially overlapped propositions in F, it seems more reasonable to make use only of this available information to produce two uncertain measures, the Belief (Bel) and the Plausibility (Pl) functions (Figure 3).

Similar to the additive rule in probability, belief and plausibility measures of proposition A can be calculated from following formula:

$$\text{Bel}(A) = \sum\_{B \subseteq A} m(B) \quad \text{for all } B \subseteq 2^{\Omega} \tag{9}$$

$$Pl(A) = \sum\_{B \cap A \neq \mathcal{O}} m(B) \quad \text{for all } B \subsetneq \mathcal{2}^{\Omega} \tag{10}$$

where A represents different elements in F. In terms of two complementary sets A and Ã, the sum of belief and plausible function is not required to be one. But the weaker rule Pl(A)+Bel(Ã)=1 is satisfied, and this expression is completely different from probability distribution function p in probability theory, that is, p(A)+p(Ã)=1. As the most remarkable distinction from probability theory, evidence theory allows evidence stemming from different sources and employs the rules of combination to aggregate [26]. One of the most important combination rules is Dempster's rule which has following formulation. Given two independent BBA m(B1) and m(B2), the Dempster's rule can be expressed as:

Figure 3. Belief function (Bel) and Plausibility function (Pl) of proposition A.

$$m(A) = \frac{\sum\_{B\_1 \cap B\_2 = A} m(B\_1) m(B\_2)}{(1 - K)} \quad \text{for all } A \ne \bigotimes \tag{11}$$

where K ¼ X B1 ∩ B2¼∅ mðB1ÞmðB2Þ can be viewed as contradiction or conflict among the informa-

tion given by the independent knowledge sources.

#### 3.2. Evidence theory-based UQ of seismic damage assessment using differential evolution

#### 3.2.1. Evidence-based uncertainty representation

For the purpose of UQ, the first step is the uncertainty representation of parameters using evidence theory, in which separate belief structures for each uncertain parameter should be constructed. In this work, we adopt a general methodology as described previously by Salehghaffari et al. [27] to obtain necessary information from available data and express the uncertain variables in the mathematical framework of evidence theory.

According to Salehghaffari et al. [27], two principle steps are involved in this methodology: (1) representation of uncertain parameters in several intervals through drawing bar charts by using all available data or directly from expert opinions and (2) identification of three relationships between all adjacent intervals and construction of the associated BBA structure. To further illustrate this, assuming that D1 and D2 represent the number of data points within two adjacent intervals I1 and I2, respectively, and D1 > D2, three relationships of two adjacent intervals can be identified as agreement ðD2=D1 ≥ 0:8Þ, conflict ð0:5 ≤ D2=D1 < 0:8Þ, and ignorance ðD2=D1 < 0:5Þ (see Figure 4), the corresponding belief structure and BBA value for these three relationships are calculated by Eqs. (12)–(14), respectively.

$$m(\{I\} = \{I\_1, I\_2\}) = (D\_1 + D\_2) / D\_T \tag{12}$$

$$m(\{I\_1\}) = D\_1 / D\_T, \qquad m(\{I\_2\}) = D\_2 / D\_T \tag{13}$$

Figure 4. Three relationships of uncertain intervals.

$$m(\{I\_1\}) = D\_1 / D\_{T\prime} \qquad m(\{I\_1, I\_2\}) = D\_2 / D\_{T\prime} \tag{14}$$

where DT denotes the total number of data points, following this approach, a reasonable BBA structure of uncertain parameter is constructed based on available data and knowledge, a more detailed illustration of uncertainty representation in intervals with assigned BBA value is referred in Salehghaffari et al. [27].

Employing this strategy, the uncertainty of Park-Ang model parameters can be properly represented with the evidence theory. In Figure 5, we use εA(β), εB(β), and ε(Fy) to denote the variability of the predicted models in Refs. [4, 5] for energy constant β and the one in Ref. [24] for yielding force Fy of columns. The εC(δu), εD(δu), and εE(δu) in Figure 6 represent the fluctuation of the empirical model for ultimate displacement under monotonic loading in Refs. [6, 22, 23].

Figure 5. Evidential uncertainty description of ε(β) and ε(Fy).

mðAÞ ¼

tion given by the independent knowledge sources.

3.2.1. Evidence-based uncertainty representation

Figure 4. Three relationships of uncertain intervals.

where K ¼

X

74 Uncertainty Quantification and Model Calibration

B1 ∩ B2¼∅

X

uncertain variables in the mathematical framework of evidence theory.

three relationships are calculated by Eqs. (12)–(14), respectively.

B1 ∩ B2¼AmðB1ÞmðB2Þ

3.2. Evidence theory-based UQ of seismic damage assessment using differential evolution

For the purpose of UQ, the first step is the uncertainty representation of parameters using evidence theory, in which separate belief structures for each uncertain parameter should be constructed. In this work, we adopt a general methodology as described previously by Salehghaffari et al. [27] to obtain necessary information from available data and express the

According to Salehghaffari et al. [27], two principle steps are involved in this methodology: (1) representation of uncertain parameters in several intervals through drawing bar charts by using all available data or directly from expert opinions and (2) identification of three relationships between all adjacent intervals and construction of the associated BBA structure. To further illustrate this, assuming that D1 and D2 represent the number of data points within two adjacent intervals I1 and I2, respectively, and D1 > D2, three relationships of two adjacent intervals can be identified as agreement ðD2=D1 ≥ 0:8Þ, conflict ð0:5 ≤ D2=D1 < 0:8Þ, and ignorance ðD2=D1 < 0:5Þ (see Figure 4), the corresponding belief structure and BBA value for these

ð1 � KÞ for all A 6¼ ∅ ð11Þ

mð{I} ¼ {I1,I2}Þ¼ðD1 þ D2Þ=DT ð12Þ

mð{I1}Þ ¼ D1=DT; mð{I2}Þ ¼ D2=DT ð13Þ

mðB1ÞmðB2Þ can be viewed as contradiction or conflict among the informa-

Figure 6. Evidential uncertainty description of ε(δu).

#### 3.2.2. Uncertainty propagation using differential evolution

In evidence theory community, uncertainty variable is usually expressed to be a series of focal element intervals based on limited information and the joint frame of discernment is composed of the Cartesian products of uncertain intervals, then, the BBA value of each element of joint frame of discernment is also the Cartesian product of BBA value assigned on the corresponding interval. Given two independent uncertain parameters u1 ∈ U1 and u2 ∈ U2, and corresponding focal element C1 and C2, the joint BBA structure of this problem is defined as:

$$\mathbf{C} = \mathbf{C}\_1 \otimes \mathbf{C}\_2 \qquad \forall \boldsymbol{\mu} \in \mathbf{U} \; \forall \boldsymbol{\mu}\_1 \in \mathbf{U}\_1 \; \forall \boldsymbol{\mu}\_2 \in \mathbf{U}\_2 \tag{15}$$

$$m(\mathbb{C}) = m(\mathbb{C}\_1)m(\mathbb{C}\_2) \tag{16}$$

where the symbol ⊗ denotes the Cartesian products. Using Eqs. (15) and (16), the joint uncertainty input of system can be seemed as the multidimensional hypercube. Therefore, uncertainty propagation is a progress of finding the maximum and minimum of the system response value in each hypercube interval (proposition of the joint belief structure). To propagate the represented uncertainties of Park-Ang damage model constants, the damage index D is considered as system response.

Considering epistemic uncertainty of the system, the belief and plausibility functions of the response are obtained on the basis of the combined BBAs of the input parameters from different information sources using the evidence combination rules. For the prediction response process D = f(Y), whose input parameter vector Y = (Y1,…,Yn) has n variables with epistemic uncertainty, the joint proposition C of elementary proposition is constructed for the Park-Ang damage index prediction system model as:

$$\mathbf{C} = \{ \mathbf{c}\_{k} = [\mathbf{x}\_{1i\_{1}}, \mathbf{x}\_{2i\_{2}}, \dots, \mathbf{x}\_{ni\_{n}}] : \mathbf{x}\_{1i\_{1}} \in \mathbf{X}\_{1}, \mathbf{x}\_{2i\_{2}} \in \mathbf{X}\_{2}, \dots, \mathbf{x}\_{ni\_{n}} \in \mathbf{X}\_{n} \} \tag{17}$$

where X1, X2, …, Xn denote the intervals sets (frame of discernment) of the n variables Y1, Y2, …, Yn, and the relevant numbers of the intervals are I1, I2,…, In. x1i1 , x2i2 , ⋯, xnin denote the subintervals, 0 ≤ jij ≤ Ij (j = 1,2,…,n); ck denotes the n-dimensional joint proposition set constructed by several subintervals, and there are I1, I2,…,In joint proposition sets ck in C. The BBA of the joint proposition set C is defined as:

$$m\_{\mathbf{c}}(\mathbf{c}\_{k}) = m\_{1}(\mathbf{x}\_{1i\_{1}})m\_{2}(\mathbf{x}\_{2i\_{2}})\cdots m\_{n}(\mathbf{x}\_{ni\_{n}}) \tag{18}$$

Thus every element of the Cartesian set C is required to be checked in the evaluation of the belief and plausibility functions by finding the system response bounds. That is to say the minimum and maximum responses of each joint set are needed to calculate:

$$[D\_{\min}, D\_{\max}] = [\min[f(\mathbf{c}\_k)], \max[f(\mathbf{c}\_k)]] \tag{19}$$

As uncertain variable is represented by many discontinuous set instead of smooth and continuous explicit function, time consuming is inevitable in UQ with evidence theory. There are two main approaches to find the bounds of the system response: sampling and optimization. The accuracy of sampling approach is highly dependent on the number of samples and the number of hypercubes, and the process is costly. On the contrary, optimization methods have the potential to dramatically reduce the computational work. To alleviate this computational burden, based on authors' previous work [11], the differential evolution (DE) [28] optimization approach is used to calculate the response bounds of each hypercube and compute the composite BBA of each hypercube, propagation of the represented uncertainty through Park-Ang damage model (Eq. (1)). The characteristics of derivative-free and capability of handling discrete belief structure make DE method to be a good choice for such an interval bound task.

3.2.2. Uncertainty propagation using differential evolution

76 Uncertainty Quantification and Model Calibration

is considered as system response.

Park-Ang damage index prediction system model as:

BBA of the joint proposition set C is defined as:

element C1 and C2, the joint BBA structure of this problem is defined as:

In evidence theory community, uncertainty variable is usually expressed to be a series of focal element intervals based on limited information and the joint frame of discernment is composed of the Cartesian products of uncertain intervals, then, the BBA value of each element of joint frame of discernment is also the Cartesian product of BBA value assigned on the corresponding interval. Given two independent uncertain parameters u1 ∈ U1 and u2 ∈ U2, and corresponding focal

where the symbol ⊗ denotes the Cartesian products. Using Eqs. (15) and (16), the joint uncertainty input of system can be seemed as the multidimensional hypercube. Therefore, uncertainty propagation is a progress of finding the maximum and minimum of the system response value in each hypercube interval (proposition of the joint belief structure). To propagate the represented uncertainties of Park-Ang damage model constants, the damage index D

Considering epistemic uncertainty of the system, the belief and plausibility functions of the response are obtained on the basis of the combined BBAs of the input parameters from different information sources using the evidence combination rules. For the prediction response process D = f(Y), whose input parameter vector Y = (Y1,…,Yn) has n variables with epistemic uncertainty, the joint proposition C of elementary proposition is constructed for the

where X1, X2, …, Xn denote the intervals sets (frame of discernment) of the n variables Y1, Y2, …, Yn, and the relevant numbers of the intervals are I1, I2,…, In. x1i1 , x2i2 , ⋯, xnin denote the subintervals, 0 ≤ jij ≤ Ij (j = 1,2,…,n); ck denotes the n-dimensional joint proposition set constructed by several subintervals, and there are I1, I2,…,In joint proposition sets ck in C. The

Thus every element of the Cartesian set C is required to be checked in the evaluation of the belief and plausibility functions by finding the system response bounds. That is to say the

As uncertain variable is represented by many discontinuous set instead of smooth and continuous explicit function, time consuming is inevitable in UQ with evidence theory. There are two main approaches to find the bounds of the system response: sampling and optimization. The accuracy of sampling approach is highly dependent on the number of samples and the number of hypercubes, and the process is costly. On the contrary, optimization methods have the

minimum and maximum responses of each joint set are needed to calculate:

C ¼ {ck ¼ ½x1i1 , x2i2 , ⋯, xnin � : x1i1 ∈ X1, x2i2 ∈ X2, ⋯, xnin ∈ Xn} ð17Þ

mcðckÞ ¼ m1ðx1i1 Þm2ðx2i2 Þ⋯mnðxnin Þ ð18Þ

½Dmin, Dmax�¼½min½fðckÞ�, max½fðck�� ð19Þ

C¼C1 ⊗ C2, ∀u∈ U ∀u1 ∈ U1 ∀u2 ∈ U2 ð15Þ

mðCÞ ¼ mðC1ÞmðC2Þ ð16Þ

DE is arguably one of the most powerful stochastic real-parameter optimization algorithms for solving complex and computational optimization problems in current use. As a novel evolutionary computation technique, differential evolution resembles the structure of an evolutionary algorithm (EA). However, unlike traditional EAs, the DE-variants perturb the current generation population members with the scaled differences of randomly selected and distinct population members. The characteristics together with other factors of DE make it a fast and robust algorithm and as an alternative to EA. Since late 1990s, DE started to find several significant applications to the optimization problems arising from diverse domains of science and engineering. In a recently published article, Das and Suganthan [29] provided a comprehensive survey of the DE algorithm and its basic concepts, different structures and variants for solving various optimization problems, as well as applications of DE variants to practical optimization problems.

In the context of DE, the individual trial solutions (which constitute a population) are called parameter vectors or genomes. Let S ∈ Rn be the search space of the problem. Then, the ndimensional vector can be represented by xi = (xi1, xi2,…, xin) T ∈ S, i = 1, 2, …, NP, and DE algorithm utilizes NP as a population for each iteration, called a generation of the algorithm. For the damage index assessment response process, its parameter vector is generated by the uncertainty variables (β, Fy, and δu) in ranges according to their respective belief structures. DE operates through the same computational steps as employed by a standard EA, including crossover, mutation, crossover, and selection operators, but differs from traditional EAs, DE employs difference of the parameter vectors to explore the objective function landscape. As above brief description, the pseudocode of DE is presented in Figure 7 and with a detailed survey of the DE family of algorithms can be found in Ref. [29].

Take the pseudocode of DE in mind, the illustration of DE-based computational strategy for finding the propagated belief structure by the example as shown in Figure 8 (only one uncertain parameter is considered).

The procedure of uncertainty propagation using the DE strategy is as follows:

$$\begin{aligned} \text{3.11 (Moundan Map): } & \text{the minimal vector is given as} \\ & \begin{aligned} \text{3.12 (+)} \quad & \text{X}\_{i}^{G+1} = \text{X}\_{i}^{G} + F\_{1} \left(\text{X}\_{\text{best}}^{G} - \text{X}\_{i}^{G}\right) + F \left(\text{X}\_{\text{rt}}^{G} - \text{X}\_{\text{rt}}^{G}\right) \\ & \text{3.2 (Cross step: for each mutually vector, the trial vector is generated using} \\ & \text{FOR } j=1 \text{ to } ID \\ & \text{IF rand}(j) \leq CR \text{ or } j\text{-randn}(i) \text{ THEN} \quad \mathfrak{U}\_{g}^{G+1} = \text{V}\_{g}^{G+1} \\ & \text{END } \text{IF} \\ & \text{END } \text{IF} \\ & \text{END } \text{IF} \\ & \text{3.3 } \text{Sdelection step: Evaluate each trial vector } \mathfrak{U}\_{i}^{(G+1)} \\ & \text{IF } j \left(\mathfrak{U}\_{i}^{G+1}\right) < f \left(\mathfrak{X}\_{i}^{G+1}\right) \quad \text{THEN } \mathcal{X}\_{i}^{G+1} = \mathfrak{U}\_{i}^{G+1} \\ & \text{OTHERWISE } \mathcal{X}\_{i}^{G+1} = \mathfrak{X}\_{i}^{G} \\ & \text{IF } j \left(\mathfrak{X}\_{i}^{G+1}\right) < f \left(\mathfrak{X}\_{\text{best}}^{G}\right) \text{ THEN } \mathcal{X}\_{\text{best}}^{G+1} = \mathfrak{X}\_{i}^{G+1} \\ & \text{OTHERWISE } \mathcal{X}\_{\text{best}}^{G+1} = \mathfrak{X}\_{\text{best}}^{G} \end{aligned} $$

Figure 7. Pseudocode of DE.

Figure 8. Uncertainty propagation of belief structure of system by DE.

Once the BBA structure of the Park-Ang damage index response is constructed, observed evidence on simulation responses is used in the determination of target propositions to estimate uncertainty measures, i.e., cumulative belief function (CBF) and cumulative plausibility function (CPF).

#### 3.2.3. Uncertainty measurement for seismic damage assessment

In evidence theory framework, the plausibility function Pl and belief function Bel are used to denote the uncertainty measurement. Employing the construction rule proposed by Sentz et al. [30], the CBF and CPF of Park-Ang damage index D less than the threshold value are formulated as follows:

$$Pl(D\_{\text{thr}}) = \sum\_{u\_D \cap \: l\_D \neq \mathcal{Q}} m(u) \; \mathcal{U}\_D = \{u\_D \le D\_{\text{thr}}\} \tag{20}$$

$$Bel(D\_{\text{thr}}) = \sum\_{\substack{\mu\_D \subset \downarrow \mathcal{U}\_D}} m(\mu) \quad \mathcal{U}\_D = \{\mu\_D \le D\_{\text{thr}}\} \tag{21}$$

Where uD ∩ UD 6¼ ∅ means that the joint focal element u can be entirely or partially within the threshold domain uD ≤ Dthre and uD ⊆ UD means that the joint focal element uD can be entirely within the threshold domain uD ≤ Dthre. Summarized above subparts, the separate stages of UQ framework of evidence theory using differential evolution optimization is shown in Figure 9.

Figure 9. Procedure of UQ of Park-Ang damage model.

Once the BBA structure of the Park-Ang damage index response is constructed, observed evidence on simulation responses is used in the determination of target propositions to

Figure 8. Uncertainty propagation of belief structure of system by DE.

Figure 7. Pseudocode of DE.

78 Uncertainty Quantification and Model Calibration

#### 4. Case study

In order to investigate the effectiveness and feasibility of the proposed UQ measures, the column "zahn86u7" [31] is selected to compute the Park-Ang damage index in its load step. The backbone curve and load history are shown in Figure 10.

As shown in Figure 10a, the ultimate cyclic displacement is calibrated by using the average value of 80% maximum force point on the force capacity reduction slope of positive and negative direction. The effective path in Figure 10b denotes the load path from initial state to ultimate state and the load path is the global displacement history. Using the properties of column, listed in the webpage of PEER, the nominal value of constants in Park-Ang damage model β, δu and Fy can be estimated by the empirical expressions from Eq. (2) to Eq. (7), respectively. In consistent with Section 3.2, the uncertainty distribution of model constants can be depicted as the nominal value multiply the factor ε. Taking the computed results into the evidence representation process, the BBA structures of β, δu, and Fy with different models are listed in Tables 1 and 2.

Taking above uncertain information into the differential evolution-based uncertainty propagation framework, the evidential UQ results for each load step as shown in Figure 11.

To validate the generality of evidence theory, the variability of Park-Ang model parameters is also represented by probability theory. The goodness of fit test is applied to test the distribution type and determine the related distribution parameters. The uncertainty distribution information of model B for β model C for δu and Fy is presented in Table 3.

From Table 3, the values of ε(β) and ε(δu) do not refuse the normal and lognormal distribution. We use two strategies to construct the probability input of variables. In first strategy, the lognormal distribution is applied to fit all the uncertainty inputs and the cumulative distribution

Figure 10. Backbone curve (a) and load path (b) of columns of column test.

Table 1. The BBA structure for multisource of β and Fy.

4. Case study

80 Uncertainty Quantification and Model Calibration

are listed in Tables 1 and 2.

In order to investigate the effectiveness and feasibility of the proposed UQ measures, the column "zahn86u7" [31] is selected to compute the Park-Ang damage index in its load step.

As shown in Figure 10a, the ultimate cyclic displacement is calibrated by using the average value of 80% maximum force point on the force capacity reduction slope of positive and negative direction. The effective path in Figure 10b denotes the load path from initial state to ultimate state and the load path is the global displacement history. Using the properties of column, listed in the webpage of PEER, the nominal value of constants in Park-Ang damage model β, δu and Fy can be estimated by the empirical expressions from Eq. (2) to Eq. (7), respectively. In consistent with Section 3.2, the uncertainty distribution of model constants can be depicted as the nominal value multiply the factor ε. Taking the computed results into the evidence representation process, the BBA structures of β, δu, and Fy with different models