A time-dependent uncertainty design and optimization strategy using an enhanced Kriging surrogate model

Abstract

Uncertainties exist widely in the time-dependent performance degeneration processes of engineering systems in practice. Generally, in order to simplify the calculation, the random processes of uncertainty information are usually treated as time-independent or monotonic processes. The corresponding uncertainty analysis approaches are time independent. However, in this situation, the failure probabilities of performance can only be considered at the end of structure lifetime. To deal with the above challenge, an enhanced outcrossing rate method using the Kriging interpolation method is proposed in this study. The proposed method can utilize the correlation information between two design variables to predict the stress level and describe the time-dependent degeneration process. An uncertainty-based design and optimization problem of machine tool spindle is utilized to illustrate the effectiveness of the proposed strategy.

Introduction

In practical engineering, designers need to consider and deal with the uncertainty information which varies in space and time.^1–10 This time-dependent information includes stochastic loads, component degeneration, and fatigue process. ¹¹ Consequently, the uncertainty evaluations for such engineering problems are usually time dependent.^12–15 The prediction and evaluation of time-dependent uncertainty are important for designers to make maintenance schedules of complex structural components. In the past decades, a lot of research has been developed and proposed in the time-dependent field to enhance safety and reliability.^16–19 However, with the increasing complexity of engineering, how to perform the precise evaluation of time-dependent reliability is still a challenge.

Recently, approximate reliability approaches have attracted more and more attention when complex structures are involved. ²⁰ Being one of these approximate reliability approaches, the crossing rate method is commonly utilized to approximate the failure rate of important performance, instead of the expensive simulation method. However, in some special cases, it is not easy or even hard to obtain the analytical solutions for estimating the outcrossing rate. To solve this problem, the response surface methods (RSMs) are introduced in time-dependent design and optimization problems.^21–27 In RSM, to match a reasonable digital model, the sample size should be at least three times that of the model coefficients. Consequently, with the increase of the number of input design variables, RSM is conservative and less efficient for highly nonlinear design problems.

As one of the popular RSMs, the Kriging method has been introduced in engineering design widely. At first, this method was utilized in deterministic optimization problems. Uncertainties are not taken into consideration during this period. ²⁸ Cai et al. ²⁹ developed a multi-point sampling method based on Kriging to solve the efficiency problem in global optimization problem. Zhang et al. ³⁰ proposed a framework based on Kriging model to solve the design problem considering both design variables and model parameters. Zhang et al. ³¹ introduced the Kriging model to enhance the calculation efficiency and solve a multi-objective multidisciplinary design optimization problem. Li et al. ³² introduced a local Kriging approximation method in most probable point (MPP)-based uncertainty design optimization. In the above methods, there is a presumption that the correlation function is exponentially or normally distributed. Furthermore, the parameters in the above distributions are estimated by maximum likelihood estimator. However, these assumptions are not suitable for time-dependent reliability evaluations. The accuracy of design and optimization under uncertainty generally relies on the accuracy of the digital model which captures the performance variations.^33–42 To solve the above problems, this study focuses on the combination of Kriging method and Bayes theory to calculate the stress spectrum and strength degeneration in full life cycle. A sequential sampling approach is also introduced to update the digital model during the optimization process.

This study is organized as follows. In section “Time-dependent stress calculation method,” a time-dependent stress calculation method is presented in detail. In section “The enhanced Kriging method,” information on enhanced Kriging method is provided. In section “The time-independent uncertainty design and optimization,” the time-independent uncertainty design and optimization is presented. In section “Example,” a structure design problem of machine tool spindle is presented to show the effectiveness of the given strategy. Section “Conclusion” provides the conclusions.

Time-dependent stress calculation method

In general, the formula of uncertainty-based design and optimization can be given as

$$\begin{array}{c}\min :f(\mathbf{d},\mathbf{x},\mathbf{p})\hfill \\ \mathrm{s}.\mathrm{t}.\mathrm{Pr}[g_i(\mathbf{d},\mathbf{x},\mathbf{p})\le 0]\ge \left[R_i^{\mathrm{T}}\right]\hfill \\ {\mathbf{d}}^L\le \mathbf{d}\le {\mathbf{d}}^U,{\mathbf{\mu }}_{\mathbf{x}}^{\mathbf{L}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{x}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{x}}^{\mathbf{U}}\hfill \\ \mathbf{i}\mathbf{=}1\mathbf{,}2\mathbf{,}\dots \mathbf{,}\mathbf{n}\hfill \end{array}$$

(1)

where $\mathbf{x}$ is the vector of random design variables, $\mathbf{d}$ is a vector of deterministic design variables, $f(•)$ is the objective function of the design optimization problem, ${\mathbf{\mu }}_{\mathbf{x}}$ is a vector of the mean values of $\mathbf{x}$ , $\mathbf{p}$ is a vector of random design parameters, $g_i(•)$ is the limit state function of the $i\mathrm{th}$ reliability constraints, and $[R_i^{\mathrm{T}}]$ is the acceptable reliability value.

In practical engineering, random information in space and time exists in many design problems, such as time-dependent loads and component degeneration. Consequently, the time-variant uncertainties should be considered and evaluated in these situations. Generally, the formula of time-dependent uncertainty-based design optimization model is as follows

$$\begin{array}{c}\min :f(\mathbf{d},\mathbf{x},\mathbf{p},\mathbf{y}(t))\hfill \\ \mathrm{s}.\mathrm{t}.⪻[g_i(\mathbf{d},\mathbf{x},\mathbf{p},\mathbf{y}(t)|t\in [0,T])\le 0]\ge \left[R_i^{\mathrm{T}}\right]\hfill \\ {\mathbf{d}}^L\le \mathbf{d}\le {\mathbf{d}}^U,{\mathbf{\mu }}_{\mathbf{x}}^{\mathbf{L}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{x}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{x}}^{\mathbf{U}}\hfill \\ \mathbf{i}\mathbf{=}1\mathbf{,}2\mathbf{,}\dots \mathbf{,}\mathbf{n}\hfill \end{array}$$

(2)

where $\mathbf{y}(t)$ is the vector of time-independent random design variables and $\mathrm{Pr}[g_i(\mathbf{d},\mathbf{x},\mathbf{p},\mathbf{y}(t)|t\in [0,T])\ge 0]$ is the probability of failure at time $t$ in the period $[0,T]$ .

Theoretically, the value of reliability can be evaluated using the following equation

$$R=1-\exp (-\underset{0}{\overset{t}{\int }}\lambda (t)\mathrm{d}t)$$

(3)

where $\lambda (t)$ is the failure rate at $t$ . The conditional probability can be utilized to denote this failure rate as

$$\lambda (t)=\underset{\mathrm{\Delta }t\to 0}{lim}\frac{p\{g_i(\mathbf{d},\mathbf{x},\mathbf{p},\mathbf{y}(t+\mathrm{\Delta }t))>0|g_i(\mathbf{d},\mathbf{x},\mathbf{p},\mathbf{y}(t))\le 0\}}{\mathrm{\Delta }t}$$

(4)

Here, we introduce the crossing rate $\psi (t)$ to approximate the above failure rate. However, it should be noted that only a few analytical solutions can exist in some situations to estimate the value of $\psi (t)$

$$\lambda (t)\approx \psi (t)$$

(5)

Generally, it is not easy to obtain the analytical solution for the degeneration process or stress spectrum under multi-loadings. In this situation, simulation-based methods are effective to calculate reliability, which is shown in equation (6)

$$P_f=\underset{G}{\int }\cdots \int I[g(\mathbf{x})\le 0]f(\mathbf{x})d\mathbf{x}$$

(6)

where $G$ denotes the failure domain.

In this study, to keep the trade-offs between computational efficiency and accuracy, some uncertainty information is assumed to be available to predict the time-dependent reliability. During the optimization procedure, the updating of parameters is performed using the sequential sampling method and the Bayes theory. The stress level can be evaluated with sequential samplings, which is shown in equation (7)

$$P({\sigma }_s|\mathrm{data})=P(\mathrm{data}|{\sigma }_s)\pi ({\sigma }_s)$$

(7)

where $P({\sigma }_s|\mathrm{data})$ is the posterior distribution of sequential information and $\pi ({\sigma }_s)$ is the prior distribution of stress.

The enhanced Kriging method

The global approximation Kriging model can be re-denoted as

$$y_{\mathrm{pre}}(\mathbf{x})=\sum _{i=1}^{n_s}{y}_i(\mathbf{x})b_i+\epsilon (\mathbf{x})$$

(8)

where $y_{\mathrm{pre}}(•)$ is the undetermined prediction value, $n_s$ is the number of basic functions in the model, $\mathbf{b}=[b_0,b_1,\dots ,b_n{]}^{\mathrm{T}}$ is the vector of coefficients with the parameters of the above digital model, and $\epsilon (\mathbf{x})$ is a normal stationary process with mean $\mu =0$ and covariance $\sigma =1$

$$\mathrm{Cov}({\mathbf{x}}_i,{\mathbf{x}}_j)={\sigma }_{\mathrm{kriging}}^{2}R({\mathbf{x}}_i,{\mathbf{x}}_j),i,j=1,\dots ,n$$

(9)

where n is the number of sample points, $R({\mathbf{x}}_i,{\mathbf{x}}_j)$ is a correlation function of points ${\mathbf{x}}_i$ and ${\mathbf{x}}_j$ , and ${\sigma }_{\mathrm{kriging}}^2$ is the process variance

$$\mathbf{R}=\left[\begin{array}{ccc}R({\mathbf{x}}_1,{\mathbf{x}}_2)\hfill & \cdots \hfill & R({\mathbf{x}}_1,{\mathbf{x}}_n)\hfill \\ \multicolumn{1}{c}{⋮}& \ddots \hfill & ⋮\hfill \\ \multicolumn{1}{c}{R({\mathbf{x}}_1,{\mathbf{x}}_n)}& \cdots \hfill & R({\mathbf{x}}_n,{\mathbf{x}}_n)\hfill \end{array}\right]$$

(10)

The correlation function $R({\mathbf{x}}_i,{\mathbf{x}}_j)$ can be determined by the designers. The prediction variance should be minimized to guarantee the accuracy of prediction

$$\begin{array}{c}\multicolumn{1}{c}{\min :\mathrm{Var}[y_{\mathrm{pre}}({\mathbf{x}}_0)]=\min :2\sum _{i=1}^n{b}_{i}R({\mathbf{x}}_i,{\mathbf{x}}_0)}\\ \multicolumn{1}{c}{-\sum _{i=1}^n\sum _{j=1}^n{b}_i{b}_{j}R({\mathbf{x}}_i,{\mathbf{x}}_0)}\end{array}$$

(11)

The sum of coefficients should be equal to one to ensure that the estimation is unbiased, as shown in equation (12)

$$\sum _{i=1}^n{b}_i=1$$

(12)

Here, we take partial derivatives with respect to $b_i,i=1,2,\dots ,n$ , to determine the coefficients. The values of these equations are set to be zero. The coefficients can be evaluated as

$$\mathbf{b}={\left[\begin{array}{ccc}R({\mathbf{x}}_1,{\mathbf{x}}_2)\hfill & \cdots \hfill & R({\mathbf{x}}_1,{\mathbf{x}}_n)\hfill \\ \multicolumn{1}{c}{⋮}& \ddots \hfill & ⋮\hfill \\ \multicolumn{1}{c}{R({\mathbf{x}}_1,{\mathbf{x}}_n)}& \cdots \hfill & R({\mathbf{x}}_n,{\mathbf{x}}_n)\hfill \end{array}\right]}^{-1}\left[\begin{array}{c}R({\mathbf{x}}_1,{\mathbf{x}}_0)\hfill \\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{R({\mathbf{x}}_n,{\mathbf{x}}_0)}\end{array}\right]$$

(13)

In the original Kriging method, the gradient-free optimization strategy is utilized to search the optimal correlation parameter $\mathbf{\theta }$ for the Kriging method. The maximization process for $\mathbf{\theta }$ is shown as

$$max:F(\mathbf{\theta }\mathbf{)}\mathbf{=}\frac{1}{2}\mathbf{ln}\mathbf{(}\mathbf{\left|}\mathbf{R}\mathbf{\right|}\mathbf{)}\mathbf{+}\frac{\mathbf{n}}{2}\mathbf{ln}\mathbf{(}{\mathbf{\sigma }}_{\mathrm{kriging}}^2\mathbf{)}$$

(14)

In this study, the uncertainty of distribution parameters and statistical model selection for correlation function are quantified. In the case where various correlation models are possible candidates, the mean value of correlation matrix with incorporating both model and parameter uncertainties is

$$E(\mathbf{R})=\sum _{i=1}^{m}P(M_i)\underset{{\theta }_i}{\int }P(\mathbf{R}|{\theta }_i,M_i)f({\theta }_i|M_i)\mathrm{d}{\theta }_i$$

(15)

where $M_i$ is the $i\mathrm{th}$ correlation model, $P(M_i)$ is the model uncertainty, and ${\theta }_i$ represent the distribution parameters within $M_i$ . Then, the posterior probabilistic distribution can be calculated by Bayes theory using the sequential observation data.

The time-independent uncertainty design and optimization

Both the loads on the structure and the performance of structural components may vary with time. Consequently, the extreme stress level and strength degeneration during the life cycle should be taken into consideration to satisfy the acceptable reliability. The Bayes theory and Kriging method are combined and utilized to evaluate the value of time-dependent reliability.

In this study, the uncertainty evaluation approaches assume that the reliabilities of different stress levels are either fully correlated $(C(s_1,s_2)=1)$ or uncorrelated $(C(s_1,s_2)=0)$ random variables. For stress prediction, the correlation between different stress levels is considered with Kriging interpolation approach.

Here, the prior distribution of regression coefficients ${\mathbf{b}}_0$ is provided in the Kriging model. Then the stress levels at different time points $t_1,t_2,\dots ,t_k$ are analyzed. As shown in equations (16)–(19), the correlation of gradients can be introduced to evaluate the correlation matrix in the case when the samples are insufficient

$$R({\sigma }_{t_l},{\sigma }_{t_m})=M(\theta ,|{\mathbf{F}}_l^k-{\mathbf{F}}_m^k|)$$

(16)

$$R=[\begin{array}{cccc}R({\sigma }_{t_1},{\sigma }_{t_1})& R({\sigma }_{t_1},{\sigma }_{t_l})R({\sigma }_{t_1},\partial {\sigma }_{t_1})& \cdots & R(\partial {\sigma }_{t_1},\partial {\sigma }_{t_l})\\ ⋮& ⋮& ⋮& ⋮\\ R({\sigma }_{t_1},{\sigma }_{t_1})& R({\sigma }_{t_1},{\sigma }_{t_1})& \cdots & R({\sigma }_{t_1},{\sigma }_{t_1})\end{array}]$$

(17)

$${\sigma }_{t_k}=f(F_1(t_k),F_2(t_k),\dots ,F_n(t_k))$$

(18)

$${\sigma }_{t_{l+1}}=\sum _{k=1}^l{\sigma }_{t_k}{b}_k+\epsilon (\mathbf{F})$$

(19)

Then, the time-dependent uncertainty design and optimization model can be given as

$$\begin{array}{c}\min :f(\mathbf{x},\mathbf{d},\mathbf{p})\hfill \\ \mathrm{s}.\mathrm{t}.P[{\sigma }_i(\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F}(t_k),t_k|k=1,2,\dots ,l,\dots ,m)\hfill \\ \multicolumn{1}{c}{-{\sigma }_{max,i}\le 0]\ge [R_i^{\mathrm{T}}]}\\ {\sigma }_i(\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F}(t),t_l)=\sum _{k=1}^{l-1}{b}_k{\sigma }_i(t_k)+\epsilon (\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F})\hfill \\ {\sigma }_i(\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F}(t),t_{l+m})=\sum _{k=1}^{m+l-1}{b}_k{\sigma }_i(t_k)+\epsilon (\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F})\hfill \\ i=1,2,\dots ,n\hfill \end{array}$$

(20)

Generally, there are two methods that can be utilized to solve the probabilistic constraints in equation (20): the MPP-based reliability method as in equation (21) and the advanced mean value method as in equation (22)

$$\begin{array}{c}\min :\sum _{k=1}^{l-1}{b}_{i,k}{\sigma }_{i,t_k}+{\epsilon }_i(\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F})-{\sigma }_{max,i}\hfill \\ \mathrm{s}.\mathrm{t}.\Vert \mathbf{u}\Vert =\left[{\beta }_i^{\mathrm{target}}\right]\hfill \end{array}$$

(21)

$$P(g_i{\sigma }_i(\mathbf{x},\mathbf{d},\mathbf{p},\mathbf{F}(t),t_l)-{\sigma }_{max,i}\le 0)\ge R_i^{\mathrm{T}}$$

(22)

In practical engineering, the assumption that the distribution of the limit state function is known is not always available. In this situation, the MPP can give necessary information to choose sampling. In this study, the MPP-based reliability analysis approaches are introduced and utilized to calculate reliability.

Example

In this section, a time-dependent uncertainty design and optimization problem of shift in computer control machine tools is introduced. The VC++ software is utilized here as the implementation environment. With the given failure data of shift, two challenges should be taken into consideration in the design stage of the shift: the acceptable stiffness requirement and the weight of shift. Consequently, we treat the inner diameter of the three parts, the external diameter of the first part, and the length of extended end segment as the input design information in this problem. The formulation of this uncertainty-based design optimization problem is given in equation (23). The objective in this case is to minimize the weight of shift. The details of the uncertainty information are given in Table 1

$$\begin{array}{c}\min :f=\pi ({\mu }_1^2-{\mu }_2^2){\mu }_4+\pi ({({\mu }_1+4\times {10}^{-3})}^2-{\mu }_2^2)\hfill \\ \multicolumn{1}{c}{\times 11\times {10}^{-3}}\\ +\pi ({({\mu }_1+3\times {10}^{-3})}^2-{\mu }_2^2)(200\times {10}^{-3}-{\mu }_4-{\mu }_5)\hfill \\ +\pi ({({\mu }_1+8\times {10}^{-3})}^2-{\mu }_3^2)({\mu }_5-16\times {10}^{-3})\hfill \\ \mathrm{s}.\mathrm{t}P\{g_1=y-y^{\mathrm{accept}}\le 0\}\ge 0.99\hfill \\ 30\le {\mu }_1\le 40;20\le {\mu }_2\le 30;15\le {\mu }_3\le 25;\hfill \\ 50\le {\mu }_4\le 60;68\le {\mu }_5\le 78\hfill \end{array}$$

(23)

Table 1.

Design details of input random uncertainty information.

Variables		Mean value	Deviation	Distribution
Design variables	$x_1$	${\mu }_{x_1}(\mathrm{m})$	$6\times {10}^{-5}\mathrm{m}$	Normal
	$x_2$	${\mu }_{x_2}(\mathrm{m})$	$6\times {10}^{-5}\mathrm{m}$	Normal
	$x_3$	${\mu }_{x_3}(\mathrm{m})$	$6\times {10}^{-5}\mathrm{m}$	Normal
	$x_4$	${\mu }_{x_4}(\mathrm{m})$	$6\times {10}^{-5}\mathrm{m}$	Normal
	$x_5$	${\mu }_{x_5}(\mathrm{m})$	$6\times {10}^{-5}\mathrm{m}$	Normal
Design parameters	$F_t$	6758 N	13514 N	Interval
	$F_r$	2456 N	4487 N	Interval
	$F_n$	2325 N	4214 N	Interval
	$y^{\mathrm{accept}}$	$1.6\times {10}^{-3}\mathrm{m}$	$1.6\times {10}^{-4}\mathrm{m}$	Normal

Because of the uncertainty in practical engineering, the stiffness is a stochastic process as shown in Figure 1. The finite element strategy is used to simulate the variation of stiffness. Here, we assume that the correlation function follows the normal distribution. The design solutions from the Kriging interpolation approach are shown in Figure 2. In this case, the optimal design solutions are x = [x₁, x₂, x₃, x₄, x₅] = [0.0336, 0.031, 0.0207, 0.058, 0.084]. It should be noted that the accuracy of the Kriging interpolation method can be accepted.

Figure 1.

The stochastic loads of shift in this case study.

Figure 2.

The mean value of deflection using the Kriging interpolation method.

Conclusion

In this study, an enhanced time-dependent uncertainty evaluation approach is developed using the Bayes theory and Kriging interpolation method. Based on the proposed method, the time-variant uncertainty design and optimization problem for complex structures can be solved effectively. Furthermore, engineering simulations during the design process for structural design are performed to search the optimal design solutions. For time-variant uncertainty in practical engineering, the multilevel extreme stress is taken into consideration in this study to satisfy the reliability requirements during the entire life period of the structure. The correlation between loads and design variables is utilized by Kriging approximation approach to predict the strength degeneration process and load spectrum. Then, the value of reliability can be obtained. In the case where the observation information is insufficient, the sensitivity of reliability to the input variables and loads is introduced and analyzed by Bayes theory. The purpose of the above process is to update the stochastic process for the correlation and stress function. The updating procedure with multi-loading is performed during the design and optimization procedure. In the case study, using Kriging interpolation method, the correlations between different stress levels are assumed to follow either normal or interval distribution. Also, the parameter uncertainties in the correlation function are taken into consideration. It should be noted that, in this study, only aleatory uncertainties are taken into consideration. In practical engineering, both aleatory and epistemic uncertainties exist widely. In the future work, we will focus on this challenge to solve the reliability-based design and optimization with mixed uncertainties.