Multidisciplinary systems with transient behavior under time-varying inputs and coupling variables pose significant computational challenges in reliability analysis. Surrogate models of individual disciplinary analyses could be used to mitigate the computational effort; however, the accuracy of the surrogate models is of concern, since the errors introduced by the surrogate models accumulate at each time-step of the simulation. This paper develops a framework for adaptive surrogate-based multidisciplinary analysis (MDA) of reliability over time (A-SMART). The proposed framework consists of three modules, namely, initialization, uncertainty propagation, and three-level global sensitivity analysis (GSA). The first two modules check the quality of the surrogate models and determine when and where we should refine the surrogate models from the reliability analysis perspective. Approaches are proposed to estimate the potential error of the failure probability estimate and to determine the locations of new training points. The three-level GSA method identifies the individual surrogate model for refinement. The combination of the three modules facilitates adaptive and efficient allocation of computational resources, and enables high accuracy in the reliability analysis result. The proposed framework is illustrated with two numerical examples.

## Introduction

Multidisciplinary analysis (MDA) has been intensively studied during the past two decades [1,2]. Current MDA methods can be roughly classified into three groups, namely, the field elimination method, the monolithic method, and the partitioned method [3,4]. Based on these MDA methods, approaches have been developed for the reliability analysis of multidisciplinary systems by considering the uncertainty sources [57]. While most of the current multidisciplinary reliability analysis (MDRA) methods have only considered aleatory uncertainty, a couple of studies have also included model uncertainty [2,8,9]. The above MDRA methods have only considered time-independent reliability analysis, whereas time-dependent behavior and uncertainty sources are of interest in many practical multidisciplinary systems [10].

In the past decades, reliability analysis methods have been investigated for engineering systems with time-dependent uncertainty sources [1114]. These methods, however, are developed for single disciplinary systems and cannot be directly applied to multidisciplinary systems, due to the complicated couplings between different disciplinary simulations. Motivated by solving this problem, Zhu et al. [15] recently developed a time-dependent reliability analysis method for a multidisciplinary system using the first-order reliability method (FORM). Since FORM could have large errors for problems with nonlinear limit-state functions, the method presented in Ref. [15] can only be applied to multidisciplinary systems where FORM is accurate [15].

This work aims to develop a general time-dependent MDRA method by exploring the surrogate modeling approach, where expensive physics simulation models are replaced by surrogate models. Monte Carlo simulation (MCS) is then performed on the surrogate models to estimate the failure probability. To guarantee the accuracy of reliability analysis, two types of strategies have been pursued for surrogate model refinement, namely, uncertainty minimization and region of interest (RoI). In the former, the surrogate model is refined by minimizing either the bias or variance in the prediction [16]. In the latter, the surrogate model is refined mainly in the RoI (i.e., the limit state), which is important for reliability analysis [17]. The limit-state surrogate modeling methods have been developed for both time-independent [1820] and time-dependent reliability analysis [11,21], and have shown excellent performance. However, current limit-state surrogate methods focus only on single disciplinary systems and are inapplicable to multidisciplinary systems.

In this paper, we are interested in pursuing the RoI strategy to improve the efficiency and accuracy of the surrogate model in time-dependent MDRA. Several challenges as detailed in Sec. 2.4 are encountered. The main challenge comes from the fact that the surrogate models built for individual disciplinary systems are not directly connected with the quantities of interest (QoIs) with respect to reliability analysis, which brings challenges in determining the RoI. Some questions that need to be answered include when to refine the surrogate models, which disciplinary surrogate models to refine, and where to refine the individual disciplinary surrogate models. This work solves these challenges through the development of a systematic adaptive surrogate modeling framework. In this framework, approaches are proposed to determine whether the surrogate models need to be refined and at what time instant to refine. Along with that, a three-level global sensitivity analysis (GSA) method is developed to adaptively allocate the computational resources to refine the individual disciplinary surrogate models. The developed framework is generic and applicable to both time-dependent and time-independent MDRA. Two numerical examples are used to demonstrate the effectiveness of the proposed framework.

The remainder of the paper is organized as follows: Section 2 provides background concepts on MDRA under time-dependent uncertainty. Section 3 describes the proposed adaptive surrogate modeling framework. Section 4 illustrates the application of the proposed methodology with two numerical examples. Concluding remarks are provided in Sec. 5.

## Background

In this section, we first introduce motivation problems of the proposed research. We then briefly review the differences between time-dependent and time-independent multidisciplinary analyses, and the general procedure of RoI-based surrogate modeling methods for reliability analysis. Based on that, we summarize the challenges in the adaptive surrogate modeling for time-dependent MDRA.

### Motivation Problems.

Coupled multidisciplinary analysis is very common in the analysis of many engineering systems. The coupling of individual disciplinary analysis is represented as the exchanging physical quantities between different computational models. The required computational effort is very large when sophisticated simulation models are used in the analysis, even for deterministic MDA [22]. One example is the aero-elastic analysis of an aircraft wing, which involves coupled finite element analysis (FEA) and computational fluid dynamics simulation as shown in Fig. 1 [4]. Since time-dependent MDRA requires repeated MDAs under different realizations of uncertainty sources (i.e., random variables and stochastic processes), it is computationally unaffordable when the original simulation models are used.

Another motivation problem is the reliability analysis of a panel structure on a conceptual hypersonic aircraft vehicle [23]. As the vehicle is subjected to hypersonic flow, an attached oblique shock is created at the forebody leading edge [24]. The resulting aerodynamic pressure causes elastic deformation of the panel, which feeds back to alter the aerodynamic pressure on the panel. The panel is also subjected to aerothermal effects (aerodynamic heating and heat transfer). The aerothermal effects are coupled with the aero-elastic effects [24]. A detailed description of this panel structure is available in Refs. [23] and [25]. Figure 2 shows the four disciplinary analysis models of the panel (aerodynamics, aerodynamic heating, heat transfer, and structural deformation) [23].

In the four-disciplinary analysis, finite element models are used in each individual disciplinary analysis to predict the response of coupling variables and response variables. A partitioned approach is employed to solve the four-disciplinary analysis over time since the partitioned method enables us to use relatively independent solution procedures for each discipline [23]. The partition sequence in the analysis is: (1) aerodynamic pressure, (2) aerodynamic heating, (3) thermal, and (4) structural [23]. Since the simulation models in each disciplinary analysis are computationally expensive, directly performing time-dependent MDRA is computationally prohibitive. New strategies are required to reduce the required computational effort in time-dependent MDRA without sacrificing the accuracy of reliability analysis.

The nodal responses of FEA or computational fluid dynamics analyses are usually used as coupling variables and response variables in MDA and MDRA. As a result, the coupling and response variables may have both spatial and temporal variability. In this paper, we only focus on coupled systems with temporal variability alone. The developed method in this paper could be integrated with the method presented in Ref. [25] to handle spatial and temporal variability simultaneously. In what follows, we first provide a generalized model of multidisciplinary analysis under time-dependent uncertainty. Based on the generalized model, we present the proposed new method.

### Generalized Multidisciplinary Analysis Under Time-Dependent Uncertainty.

Figure 3 shows a generalized multidisciplinary system with three individual disciplines [26]. The notations of different variables are also given in the figure. In Secs. 24, we use $X$ to represent input random variables and $Y(t)$ to represent input stochastic processes. For the ith discipline, the inputs include random variables ($Xi, Xs$), stochastic processes ($Yi(t), Ys(t)$), and time ($t$). The subscript s refers to shared inputs across multiple disciplines. There are also coupling variables ($Lij(t)$) between the ith and jth disciplines. The response variables ($Zi(t)$) are functions of the coupling variables and input variables.

In time-independent MDRA, MDA is a deterministic analysis for a given realization of uncertainty sources, and will converge after several iterations. This situation also holds for time-dependent MDRA if the time-dependent inputs only affect the system response variables but not the coupling variables [15]. However, when the time-dependent inputs affect the coupling variables, the coupling variables will change over time; thus we do not have distributions of converged values of these coupling variables and will not be able to decouple the disciplinary analyses. As a result, the computational effort becomes prohibitive in time-dependent MDRA, due to the need to run the fully coupled MDA over the entire duration of interest for each realization of the inputs. In this situation, the partitioned method of MDA needs to be employed. As mentioned previously, a partitioned approach to MDA will enable the use of separate solution procedures for each disciplinary analysis and also handle different time intervals in coupled multidisciplinary analysis, i.e., the same time-step or identical discretization for individual disciplines is not required. In addition, note that the partition sequence may affect the result of MDA; the proposed MDRA method is based on the assumption of a predetermined partition sequence from deterministic MDA.

As shown in Fig. 4, at any time instant $ti$, the MDA given in Fig. 3 can be reformulated as Fig. 4. Based on the formulation given in Fig. 4, the MDA becomes a multilevel analysis problem at a given time instant and the analysis is coupled together between time instants.

According to Fig. 4, the models of coupling variables and response variables are given by
$Lkj(ti)=gLkj(Xs,Xk,Ys(ti),Yk(ti),Lj>k,k(ti−1),Lj
(1)
$Zk(ti)=gZk(Xs, Xk, Ys(ti), Yk(ti), Lk•(ti),Lj>k,k(ti−1),Lj
(2)

where $gLkj(⋅)$ is a vector of models of coupling variables $Lkj(ti)$, $gZk(⋅)$ is a vector of models of response variables $Zk(ti)$, $Lk•(ti)$ is a vector of all coupling variables from the kth discipline, $Lj>k,k(ti−1)$ is a vector of coupling variables from the jth discipline to the kth discipline with j > k, and $Lj are coupling variables from the jth discipline to the kth discipline with j < k.

Next, we briefly review the procedure of adaptive surrogate-based reliability analysis. After that, we discuss the challenges in adaptive surrogate modeling for time-dependent MDRA.

### Adaptive Surrogate Modeling for Reliability Analysis.

In adaptive surrogate-based reliability analysis, the original computer simulation models are substituted with surrogate models. Training points are added adaptively in the RoI (i.e., around the limit state) to guarantee the efficiency and accuracy of reliability analysis. As discussed in Sec. 1, adaptive surrogate modeling approaches were first developed for time-independent reliability analysis and have been extended to time-dependent reliability analysis of single disciplinary systems [11]. The general procedure of adaptive surrogate modeling for time-dependent reliability analysis can be summarized as follows:

• Step 1: Build the initial surrogate model and generate random realizations of random variables and stochastic processes.

• Step 2: Perform predictions using the surrogate model and samples generated in step 1.

• Step 3: Check convergence of the failure probability estimate. If the convergence criterion is satisfied, obtain the failure probability estimate. Otherwise, go to the next step.

• Step 4: Identify a new training point in the region of interest, update the surrogate model by adding the new training point into current training points, and then go to step 2.

In step 4, the identification of new training points in the RoI can be achieved based on different learning functions (e.g., the expected feasibility function defined in Ref. [17] or the U function defined in Ref. [27]) or performing global sensitivity analysis [28]. Even if the detailed implementation procedures may vary with methods, most current adaptive surrogate modeling approaches can be generalized into the above summarized procedure.

### Adaptive Surrogate Modeling for MDRA.

Depending on the system configuration defined in the reliability block diagram (RBD) included in Fig. 3, the time-dependent reliability of multidisciplinary system can be expressed as a function of the QoIs, $Zk(t)$, $k=1, 2,… , Nd$, where $Nd$ is the number of disciplines and $Zk(t)$ are the response variables given in Eq. (2). For systems with response variables defined in series and parallel system configurations in the RBD, the system reliability is expressed as
$Rseries(t0, te)=Pr{∩k=1NdZk(ti)≥0, ∀ti∈[t0, te]}$
(3)
$Rparallel(t0, te)=Pr{∪k=1NdZk(ti)≥0, ∀ti∈[t0, te]}$
(4)

where $Zk(ti)≥0$ is the safe state of the kth disciplinary response variables, “$∀$” means “for all,” $[t0, te]$ is the time interval of interest, and $Rseries(t0, te)$ and $Rparallel(t0, te)$ are the reliability values of series and parallel system, respectively.

For systems with other configurations, the time-dependent reliability needs to be defined according to the reliability block diagram. The objective of this paper is to develop a framework to efficiently and accurately estimate the time-dependent multidisciplinary system reliability by pursuing the adaptive surrogate modeling strategy. However, the procedure summarized in Sec. 2.3 cannot be directly applied to MDRA due to the difference between multidisciplinary system and single disciplinary system. A major difference is that the QoIs in a single disciplinary system are only affected by the uncertainty sources of that discipline, whereas the QoIs in a multidisciplinary system are affected by not only the uncertainty sources in their own discipline, but also by those in the other disciplines. The uncertainty sources are present in each discipline and the uncertainty accumulates as the analysis proceeds from one model to the other and over time. This is also why current available surrogate modeling method for single disciplinary system reliability analysis [29] cannot be applied to solve the MDRA problem. To efficiently and accurately estimate the time-dependent failure probability, a question that needs to be answered is how to effectively allocate the computational resources to train the individual disciplinary surrogates. To answer this question, the following challenges need to be addressed:

• Whether: After we replace the original simulation models with surrogate models, the first question that needs to be answered is whether the surrogate models are good enough to be used to estimate the failure probability. This question is answered by quantifying the uncertainty in the failure probability estimate and assessing whether the uncertainty is larger than a threshold.

• Where: At what realizations of the inputs and uncertainty sources to refine the surrogate models? If the surrogate models of the individual disciplinary analyses are not directly connected with the failure probability estimate, currently available learning functions cannot be used to determine where to refine the surrogate models.

• When: As shown in Fig. 4, the multidisciplinary system simulation under time-dependent uncertainty is performed iteratively over time. The uncertainty due to various sources will accumulate over time. At what time instant to refine the surrogate models needs to be addressed?

• Which: After we identified when and at what realization of the uncertainty sources to refine the surrogate models, the next challenge is to determine which disciplinary surrogates to refine and which surrogate model within the discipline to refine (since each disciplinary analysis is substituted with multiple surrogate models, one for each output of the disciplinary analysis).

In Sec. 3, we develop a framework called adaptive surrogate-based multidisciplinary analysis of reliability over time (A-SMART) to address the above questions.

## Proposed Method

In this section, we first provide an overview of the proposed A-SMART framework. Following that, we explain the framework in detail.

### Overview of the Proposed Framework.

In the proposed A-SMART framework, for the purpose of generalization, we assume that all the simulation models are expensive. We substitute individual disciplinary simulation models with surrogates. In other words, we build surrogate models for all $Lkj(t)$ and $Zk(t)$ and determine how to allocate computational resources to refine these surrogate models. Figure 5 shows the flowchart of the main steps of the proposed framework.

The proposed framework (as given in Fig. 5) consists of three main modules, namely, initialization, uncertainty propagation, and three-level GSA. The initialization and uncertainty propagation modules inherit the basic idea of the single loop Kriging surrogate modeling approach presented by Hu and Mahadevan [11]. In these two modules, initial surrogate models are built for the response and coupling variables. Then the surrogate model uncertainty sources are propagated through the time-dependent MDA to quantify the effects of surrogate model uncertainty on the time-dependent failure probability estimate. Based on the results of uncertainty propagation, we can quantify the uncertainty in the failure probability estimate from the current set of surrogate models and thus determine whether the surrogate models need to be refined, where to refine, and when to refine. In the three-level GSA module, we overcome the challenge that the surrogate model uncertainty sources are not directly connected to the time-dependent failure probability by performing GSA at three levels. The first two modules address the challenges of whether, where, and when summarized in Sec. 2.4. The last module solves challenge of which. A brief summary of the three modules corresponding to the procedure given in Sec. 2.3 is given as follows:

• Module 1—Initialization: In this module, random realizations of $Xk, Xs$ and $Yk(t), Ys(t)$, $k=1, 2, …, Nd$ are generated. Initial surrogate models are built for $Lkj(t)$ and $Zk(t)$. This is similar to step 1 of the procedure summarized in Sec. 2.3.

• Module 2—Uncertainty propagation: This module includes steps 2–4 shown in Fig. 5. The basic idea is similar to steps 2–3 discussed in Sec. 2.3. Due to the multilevel simulation at each time instant and the coupling between different time instants, the implementation details are quite different from that of single disciplinary analysis even if the basic idea is similar.

• Module 3—Three-level GSA: This module covers steps 5–8 of the flowchart presented in Fig. 5. In this module, we first identify the response variable that makes the highest contribution to the uncertainty in the failure probability estimate (this is level 1 GSA) since the response variables are directly connected to the failure probability. After the response variable is identified, we identify which surrogate model (i.e., response surrogate model or coupling variable surrogate model) to refine to reduce the uncertainty in the response variable through level 2 and 3 GSAs.

In Secs. 3.23.4, we explain the three modules in details. For the sake of illustration, we use a multidisciplinary system with two disciplines to explain the three modules.

### Module 1: Initial Surrogate Modeling.

According to the description of multidisciplinary system given in Sec. 2.2, at a given time instant $ti$, a two disciplinary system is simulated as shown in Fig. 6(a). For each discipline, there are models of coupling variables and response variables. For example, as depicted in Fig. 6(b), for the kth discipline, there are models of $Lk·(t)$ and $Zk(t)$, $∀k=1, …, Nd$.

In module 1, we build initial surrogate models for $Lk·(t)$ and $Zk(t)$, $∀k=1, …, Nd$. In this paper, the surrogate models are built using the Kriging surrogate method. For the sake of illustration, in Sec. 3.2, we use $XLk=[Xk, Xs, Yk(ti),Ys(ti), Ljk,k(ti−1)]$ to represent all the inputs of $Lk·(t)$.

When we are building the surrogate models, we do not know the domains of the coupling variables. To generate training points, we need to get a rough guess of the domains by performing several deterministic MDAs over the time duration of interest with fixed random realizations of random variables and stochastic processes. Once we have a rough guess of the domains, we generate training points of $XLk$, $Lk·$, and t. We then evaluate the models of individual disciplines separately at these training points. Based on the training points, we build the initial surrogate models for $Lk·(t)$ and $Zk(t)$, $∀k=1, …, Nd$. Note that the surrogate models are built for one time-step MDA as indicated in Fig. 4. Here, the generated training points are therefore for partitioned MDA at one time-step. In addition to the generated training points, the MDA iterations used to get the initial guess of coupling variables ranges can also be added into the pool of initial training points for surrogate modeling. Even if the initial training points may not fully cover the ranges of $XLk$, $Lk·$, and t, more training points will be added in important regions later in the proposed method through adaptive sampling approach to improve the accuracy of reliability analysis.

From the surrogate modeling, we have the initial surrogate models as
$L̂kj(t)=ĝLk·(XLk, t)$
(5)
$Ẑk(ti)=ĝZk(XLk,Lk•(t), t)$
(6)

where $ĝLk·$ and $ĝZk$ represent the surrogate model approximations of $gLk·$ and $gZk$, respectively.

These surrogate models are then used to substitute the original simulation models in time-dependent MDA. Along with the surrogate models, in the initialization module, we also generate random samples for $Xk, Xs, Yk(ti), Ys(ti)$, $∀k=1, …,Nd; i=1, 2, …,Nt$, where $Nt$ is the number of time instants in the time interval $[t0, te]$. We denote the generated random samples as $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, $k=1, …,Nd; i=1, 2, …, Nt; q=1,2, …, NMCS$, where $xk(q)$ and $xs(q)$ represent the qth sample of $Xk$ and $Xs$, $yk(q)(ti)$ and $ys(q)(ti)$ are the qth trajectory of $Yk(t)$ and $Ys(t)$ at time instant $ti$, and $NMCS$ is the number of MCS samples.

For a given realization $XLk=xLk$ and $Lk•(t)=lk•(t)$, we have surrogate model predictions as
$L̂kj(t)=[L̂kj(1)(t), L̂kj(2)(t), …, L̂kj(NkjL)(t)], where L̂kj(p)(t)∼N(μLkj(p)(t), σLkj(p)2(t))$
(7)
in which $NkjL$ is the number of coupling variables from the kth discipline to the jth discipline, $L̂kj(p)(t)$ is the prediction of the pth coupling variable from the kth discipline to the jth discipline, $μLkj(p)(t)$ and $σLkj(p)2(t)$ are the mean and variances of the prediction of $L̂kj(p)(t)$ at time instant $t$, and
$Ẑk(t)=[Ẑk(1)(t), Ẑk(2)(t), …, Ẑk(NkZ)(t)], where Ẑk(p)(t)∼N(μZk(p)(t), σZk(p)2(t))$
(8)

where $NkZ$ is the number of response variables from the kth discipline, and $μZk(p)(t)$ and $σZk(p)2(t)$ are the mean and variances of the prediction of $Ẑk(p)(t)$ at time instant $t$. The mean and variance are obtained from Kriging prediction.

The above equations imply that after substituting the original simulation models with surrogate models, we introduce surrogate model uncertainty into MDA. Next, we will discuss the propagation of the surrogate model uncertainty and investigate how to determine the quality of the surrogate models based on the uncertainty propagation.

### Module 2: Uncertainty Propagation in Surrogate Model-Based MDA

#### Uncertainty Analysis of Surrogate Model-Based MDA.

In this work, to investigate the effect of surrogate model uncertainty on the failure probability estimate, we first quantify the uncertainty in the response variables due to the uncertainty of the coupling and response surrogates for any given realization of the input variables.

Taking the multidisciplinary system shown in Fig. 6(a) as an example, for the qth sample of the input variables (i.e., $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, $k=1, 2; i=1, 2, …, Nt;$) generated in Sec. 3.2, the MDA under uncertainty is a multimodel uncertainty quantification (UQ) problem as shown in Fig. 7. In the network, square nodes represent deterministic nodes and elliptical nodes represent random nodes. At any time instant $ti$, a multimodel UQ is performed and the uncertainty quantified at $ti$ will be propagated to $ti+1$ due to the coupling between time instants. This process continues until the maximum number of iterations is reached.

Even if the output of each random node follows a normal distribution as given in Eqs. (7) and (8) for a fixed input, it is hard to get analytical expressions for the forward UQ over time. In this work, the sampling-based method is used to quantify the uncertainty of response variables and coupling variables over time. At the first time instant $t1$, $L21(t0)$ is a deterministic node with initial condition as its value. The deterministic value is then passed to node $L̂12(t1)$ as shown in Fig. 7. $L̂12(t1)$ is a vector of random variables given by $L̂12(t1)=[L̂12(1)(t1), L̂12(2)(t1), …, L̂12(N12L)(t1)], where L̂12(p)(t1)∼N(μL12(p)(t1), σL12(p)2(t1))$ due to the surrogate model uncertainty. $L̂21(t1)$ and $Ẑ1(t1)$ are vectors of random variables conditioned on $L̂12(t1)$. In order to obtain the distributions of $L̂21(t1)$ and $Ẑ1(t1)$, we first generate $NMCS2$ samples of $L̂12(t1)$. Defining the samples of $L̂12(t1)$ as $l12(w)(t1)$, $w=1, 2, …, NMCS2$, the samples of $L̂21(t1)$ and $Ẑ1(t1)$ are then generated based on the following distributions:
${L̂21(p)(t1)|l12(w)(t1)}∼N(μL21(p)(t1)|l12(w)(t1), σL21(p)2(t1)|l12(w)(t1)), ∀p=1, 2, …, N21L$
(9)
${Ẑ1(p)(t1)|l12(w)(t1)}∼N(μZ1(p)(t1)|l12(w)(t1), σZ1(p)2(t1)|l12(w)(t1)), ∀p=1, 2, …, N1Z$
(10)

where $μL21(p)(t1)|l12(w)(t1)$ and $σL21(p)2(t1)|l12(w)(t1)$ are the conditional mean and variance predicted from the Kriging surrogate model by using $[x2(q), xs(q), y2(q)(t1), ys(q)(t1), l12(w)(t1)]$ as inputs.

Similarly, we obtain the random samples of $Ẑ2(t1)$. The samples of $L̂21(t1)$ are then passed to the next time instant, and the above process is repeated until the end of the duration of interest is reached. Note that the above forward UQ is based on a given realization of uncertain inputs $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, $k=1,2; i=1, 2, …, Nt$. From the uncertainty analysis, we have the random samples of response variables as
$zk(p, w)(ti), k=1, …, Nd; p=1, …, NkZ; w=1, 2, …, NMCS2; i=1, 2, …, Nt$
(11)

where $zk(p, w)(ti)$ is the wth sample of the pth response variable of the kth discipline at $ti$.

For the system given in Fig. 7, $Nd=2$. Next, we study how to determine whether the surrogate models need to be refined or not based on the random samples given in Eq. (11).

#### Error Analysis of the Multidisciplinary System Failure Probability Prediction.

As indicated in Eq. (10), for given $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$ and time instant $ti$, $Ẑk(p)(ti), ∀p=1, 2, …, NkZ$ is a random variable. Based on the sign of $μZk(p)(ti)$, we can predict the safety state of discipline $k$ at time instant $ti$. However, due to the uncertainty in $Ẑk(p)(ti)$, we may make an error on the sign of the actual $Zk(p)(ti)$ (where $Zk(p)(ti)$ is a deterministic value from the original simulation), if the sign of the mean prediction is used as the sign of actual response. To quantify the probability of making an error on the sign of the actual $Zk(p)(ti)$, the following probability is defined based on the concept of U function presented in Ref. [27]:
$Perror(p, q, k, i)={Pr{Ẑk(p)(ti)<0}, if μZk(p)(ti)>01−Pr{Ẑk(p)(ti)<0}, otherwise$
(12)

where $Perror(p, q, k, i)$ is the probability of making an error on the sign of the pth response of the kth discipline at time instant i for given $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, and $Pr{⋅}$ stands for probability.

Using the samples generated in Sec. 3.3.1, Eq. (12) is estimated as
$Perror(p, q, k, i)≈{∑w=1NMCS2I(zk(p, w)(ti))/NMCS2, if μZk(p)(ti)>01−∑w=1NMCS2I(zk(p, w)(ti))/NMCS2, otherwise$
(13)
in which $I(zk(p, w)(ti))=1$, if $zk(p, w)(ti)<0$ and $I(zk(p, w)(ti))=0$, otherwise.
Equation (13) gives the probability of making an error at a given time instant. In time-dependent reliability analysis, the safety state for given $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, $i=1, 2, …, Nt$, is predicted using the mean prediction of $Zk(p)(ti)$ as follows:
$Ĩk,p(q)=I(mini=1, 2, …, Nt{μZk(p)(ti)})$
(14)

where $Ĩk,p(q)$ is the safety state indicator of the pth response of the kth discipline for the qth realization of random variables and stochastic processes.

If uncertainty in the mean prediction is considered, the safety state is computed using the surrogate model predictions as follows:
$Ik,p(q)=I(mini=1, 2, …, Nt{Ẑk(p)(ti)})$
(15)
In the above equation, $Ik,p(q)$ is a random variable due to the uncertainty in the surrogate model prediction $Ẑk(p)(ti)$. Due to the surrogate modeling uncertainty, we may make an error on the sign of the true safety state when we use $Ĩk,p(q)$ to substitute $Ik,p(q)$. The probability of making an error is given by
$Perrort(p, q, k)={Pr{Ik,p(q)=1}=Pr{∪i=1NtẐk(p)(ti)<0}, if Ĩk,p(q)=0Pr{Ik,p(q)=0}=1−Pr{∪i=1NtẐk(p)(ti)<0}, otherwise$
(16)
In the above equation, $Ẑk(p)(ti)$, $i=1, 2, …, Nt$, are correlated random variables due to the property of Kriging surrogate model [28]. Solving the probability given in Eq. (16) is not straightforward since it is also a time-dependent reliability analysis. Even if methods have been developed to estimate the probability in the form of Eq. (16) [30], we directly use the samples generated in Sec. 3.3.1 to approximate this probability to reduce the complexity of the proposed method. Based on Eq. (11), we approximate Eq. (16) as
$Perrort(p, q, k)≈{1NMCS2∑w=1NMCS2I(mini=1, 2, …, Nt{zk(p, w)(ti)}), if Ĩk,p(q)=01−1NMCS2∑w=1NMCS2I(mini=1, 2, …, Nt{zk(p, w)(ti)}), otherwise$
(17)

In addition, if there exists such a $ti$ that $I(μZk(p)(ti))=1$ and $Perror(p, q, k, i)<0.05$, we can directly assign $Perrort(p, q, k)=0$. This is due to the fact that once the response is classified as failed at a time instant and the probability of making an error on the classification is low, the corresponding trajectory can be classified as failed.

The above analysis is for only the qth response of the kth discipline. For a multidisciplinary system with $Nd$ disciplines, the safety state of the system needs to be predicted based on the reliability block diagram and the safety state indicators of individual disciplinary responses. For example, if the reliability block diagram is a series system, the system safety state computed using the mean prediction is given by
$ĨS(q)={1, ∑k=1Nd∑p=1NkZĨk,p(q)>00, otherwise$
(18)
Similar to the safety state of individual disciplinary response (i.e., Eq. (15)), the system safety state $IS(q)$ is a random variable when surrogate model predictions are used in Eq. (18). The corresponding probability of making an error on the sign of the system state for the series system is computed by
$Perrorsys(q)≈{Pr{IS(q)=1}=1-∏k=1Nd∏p=1NkZ[1−1NMCS2∑w=1NMCS2I(mini=1, 2, …, Nt{zk(p, w)(ti)})], if ĨS(q)=0Pr{IS(q)=0}=∏k=1Nd∏p=1NkZ[1−1NMCS2∑w=1NMCS2I(mini=1, 2, …, Nt{zk(p, w)(ti)})], otherwise$
(19)
For other types of reliability block diagram, $Perrorsys(q)$ needs to be computed based on the system configuration and Eq. (16). The above discussed procedure (i.e., Secs. 3.3.1 and 3.3.2) is performed for all $xk(q), xs(q), yk(q)(ti), ys(q)(ti)$, $q=1, 2, …, NMCS$. From that, we obtain $Perror(p, q, k, i)$, $∀q=1, 2, …, NMCS$; $k=1, …, Nd; p=1, 2, …, NkZ; i=1, 2, …, Nt$, $IS(q)$, and $Perrorsys(q)$ using Eqs. (13), (18), and (19), respectively. Based on the values of $Perrorsys(q)$ and $IS(q)$, we then estimate the potential relative error of the time-dependent failure probability estimate by implementing the convergence criterion developed in Ref. [11] as follows:
$εrmax=maxNf2*∈[0, N2]{|Nf2−Nf2*|Nf1+Nf2*×100%}$
(20)
where $Nf1=∑q=1NMCSIN1(q)$ and $Nf2=∑q=1NMCSIN2(q)$, and
$IN1(q)={1, if ĨS(q)=1 and Perrorsys(q)<0.05 0, otherwise, ∀q=1, 2, …, NMCS$
(21)
$IN2(q)={1, if ĨS(q)=1 and Perrorsys(q)≥0.050, otherwise, ∀q=1, 2, …, NMCS$
(22)
Based on Eq. (20), we determine whether the surrogate models need to be refined or not. In this paper, if $εrmax>5%$, we refine the surrogate models. Otherwise, we get the time-dependent failure probability estimate as
$pf(t0, te)≈∑q=1NMCSĨS(q)/NMCS$
(23)
If the surrogate models need to be refined, the input setting at which the surrogate models need to be refined is identified by
$q*=argmaxq=1, 2, …, NMCS{Perrorsys(q)}$
(24)
The basic principle of Eq. (24) is to identify the realization of input random variables and stochastic processes that give the highest probability of making an error on the safety state prediction of the system. After $xk(q*), xs(q*), yk(q*)(ti), ys(q*)(ti)$, $i=1, 2, …, Nt$ are identified using Eq. (24), the time instant on the trajectory that the surrogate models need to be refined is obtained by
$ik*=argmaxi=1, 2, …, Nt{Perror(p, q*, k, i)}$
(25)

The above equation indicates that if the kth disciplinary response needs to be refined, it should be refined at the input setting $xk(q*), xs(q*), yk(q*)(tik*), ys(q*)(tik*), tik*$. From Eqs. (20)(25), the questions of whether, where, and when (summarized in Sec. 2.4) have been answered. Next, we investigate which surrogate to refine through the development of a three-level GSA scheme.

### Module 3: Three-Level GSA.

As shown in Fig. 5, there are three levels of analysis: coupling variables→response variables→system safety indicator. To minimize the maximum probability of making an error $Perrorsys(q*)$ by reducing the surrogate model uncertainty, a three-level GSA scheme is developed.

#### Level 1 GSA: Which Discipline to Refine?.

Global sensitivity analysis ranks the contribution of each input variable (A) on the variance Var(B) of a quantity of interest (B). Sobol' indices can be used to quantify the uncertainty contributions, with two types of indices: first-order indices and total indices. The first-order index measures the contribution of variables without considering its interactions with the other variables and is given by [28]
$SiI=VarAi(EA∼i(B|Ai))/Var(B)$
(26)

where $Ai$ is the ith input variable, $A∼i$ is the vector of variables excluding variable $Ai$, $Var(B)$ is the variance of quantity of interest (B), and $EA∼i(B|Ai)$ is the expectation by freezing $Ai$.

The total effects index of $Ai$ considers its interactions with other variables and is given by
$SiT=1−VarA∼i(EAi(B|A∼i))/Var(B)$
(27)
In level 1 GSA, our quantity of interest is the system safety indicator $IS(q*)$ and the corresponding inputs are the safety indicators $Ik,p(q*), k=1, 2, …, Nd; p=1, 2, …, NkZ$ of individual disciplinary response variables. Due to the surrogate model uncertainty, $IS(q*)$ and $Ik,p(q*), k=1, 2, …, Nd; p=1, 2, …, NkZ;$ are random variables. Using the first-order indices given in Eq. (26), the first-order index of the safety indicator of the pth response of the kth discipline is computed by
$Sk,pI=[VarIk,p(EIk,p∼(IS(q*)|Ik,p(q*)))]/Var(IS(q*))$
(28)
Since both $Ik,p(q*)$ and $IS(q*)$ are binary random variables, Eq. (28) is rewritten as
$Sk,pI=VarIk,p(EIk,p∼(IS(q*)|Ik,p(q*)))Pr{IS(q*)=1}−Pr{IS(q*)=1}2$
(29)
For $EIk,p∼(IS(q*)|Ik,p(q*))$, we have
$EIk,p∼(IS(q*)|Ik,p(q*)=0)=Pr{IS(q*)=1|Ik,p(q*)=0}$
(30)
$EIk,p∼(IS(q*)|Ik,p(q*)=1)=Pr{IS(q*)=1|Ik,p(q*)=1}$
(31)
Based on Eqs. (30) and (31), we have
$VarIk,p(EIk,p∼(IS(q*)|Ik,p(q*)))=Pr{IS(q*)=1|Ik,p(q*)=0}2Pr{Ik,p(q*)=0}+Pr{IS(q*)=1|Ik,p(q*)=1}2Pr{Ik,p(q*)=1}−Pr{IS(q*)=1}2$
(32)
Combining Eqs. (29) and (32), we have
$Sk,pI=1Pr{IS(q*)=1}−Pr{IS(q*)=1}2{Pr{IS(q*)=1|Ik,p(q*)=0}2Pr{Ik,p(q*)=0}+Pr{IS(q*)=1|Ik,p(q*)=1}2Pr{Ik,p(q*)=1}−Pr{IS(q*)=1}2}$
(33)

where $Pr{IS(q*)=1}$ is computed using Eq. (19), $Pr{Ik,p(q*)=1}$ is computed using Eq. (16), and $Pr{IS(q*)=1|Ik,p(q*)=0}$ is computed by combining the RBD with Eq. (16).

With the first-order sensitivity indices, the discipline response that needs to be refined is obtained as
$[k*, p*]=argmaxk=1, 2, …, Nd; p=1, 2, …, NkZ{Sk,pI}$
(34)
At the same time, according to Eq. (25), the time instant that the $p*$th response of the $k*$ discipline needs to be refined is $tik*$. After the response that needs to be refined is identified, we go to the level 2 GSA to determine which surrogate model to refine since there are response surrogate and several coupling surrogates as indicated in Fig. 6(b).

#### Level 2 GSA: Response Surrogate or Coupling Surrogate?.

At the $ik*$th time instant, as shown in Fig. 6(b), the surrogate model of the $p*$th response of the $k*$ discipline is given by
$Ẑk*(p*)(tik*)=ĝZk*(p*)(xk*(q*), xs(q*), yk*(q*)(tik*), ys(q*)(tik*), Ljk*,k*(tik*−1), Lk·(tik*), tik*)$
(35)
In the above equation, $q*$ is obtained using Eq. (24), the random input variables are $Lj, $Lj>k*,k*(tik*−1)$, and $Lk·(tik*)$ due to the surrogate model uncertainty of the coupling variables, and the QoI is $Ẑk*(p*)(tik*)$. The distributions of $Lj, $Lj>k*,k*(tik*−1)$, and $Lk·(tik*)$ are determined by the random samples generated in Sec. 3.3.1. For a given realization of $l(tik*)=[Ljk*,k*(tik*−1)=lj>k*,k*(tik*−1), Lk·(tik*)=lk·(tik*)]$, as discussed in Eq. (8), $Ẑk*(p*)(tik*)∼N(μZk*(tik*)|l(tik*), σZk*(tik*)|l(tik*))$ is a random variable due to the surrogate model uncertainty of $ĝZk*(⋅)$. We are interested in analyzing the contributions of the coupling surrogate model uncertainty and response surrogate model uncertainty to the uncertainty of $Ẑk*(p*)(tik*)$. To separate these two sources of uncertainty, we introduce the cumulative distribution function (CDF) value $UZk*(p*)$ of $Ẑk*(p*)(tik*)$ as an input variable into Eq. (35). We then have
$Ẑk*(p*)(tik*)=μZk*(p*)(tik*)|L¯(tik*)+UZk*(p*)σZk*(p*)(tik*)|L¯(tik*)$
(36)
in which $L¯(tik*)=[Ljk*,k*(tik*−1), Lk·(tik*)]$, and $μZk*(p*)(tik*)|L¯(tik*)$ and $σZk*(p*)(tik*)|L¯(tik*)$ are the conditional mean and standard deviation predicted from the response surrogate model $ĝZk*(⋅)$.
Based on Eq. (36), we compute the first-order sensitivity index of coupling variables $L¯(tik*)$ by grouping them together and use Eq. (26) as follows:
$SL¯I=VarL¯(EL¯∼(Zk*(p*)(tik*)|L¯(tik*)))Var(Zk*(p*)(tik*))=VarL¯(μZk*(p*)(tik*)|L¯(tik*))Var(Zk*(p*)(tik*))$
(37)
In the above equations, $VarL¯(μZk*(p*)(tik*)|L¯(tik*))$ and $Var(Zk*(p*)(tik*))$ are computed directly using the samples generated in Sec. 3.3.1 (i.e., Eqs. (10) and (11)). If $SL¯I, where $Se1$ is a threshold for the first-order sensitivity index ($Se1=0.5$ for the numerical examples given in Sec. 4), it means that the uncertainty in $Zk*(p*)(tik*)$ mainly comes from the surrogate model uncertainty of $ĝZk*(⋅)$. In that case, we will refine $ĝZk*(⋅)$. The new training point of $ĝZk*(⋅)$ is identified using the weighted mean square error criterion [31] as below:
$w*=argmaxw=1, 2, …, NMCS{f(l¯(w)(tik*))[σZk*(p*)2(tik*)|l¯(w)(tik*)]}$
(38)

where $l¯(w)(tik*)=[ljk*,k*(w)(tik*−1), lk·(w)(tik*)]$ is the wth sample of $L¯(tik*)$ generated in Sec. 3.3.1, $f(l¯(w)(tik*))$ is the joint probability density function of $l¯(w)(tik*)$ computed using kernel smoothing density function, and $σZk*(p*)2(tik*)|l¯(w)(tik*)$ is the conditional variance.

From Eq. (38), we have the new training point of $ĝZk*(⋅)$ as $[xk*(q*), xs(q*), yk*(q*)(tik*), ys(q*)(tik*), ljk*,k*(w*)(tik*−1), lk·(w*)(tik*), tik*]$. Based on the new training point, surrogate model $ĝZk*(⋅)$ is updated. After that, the level 2 GSA is performed again until we have $SL¯I≥Se1$.

If we get $SL¯I≥Se1$ in the first round level 2 GSA, it implies that the uncertainty in $Zk*(p*)(tik*)$ mainly comes from the coupling surrogate models. In that case, we will go to the level 3 GSA to decide, which coupling surrogate to refine.

#### Level 3 GSA: Which Coupling Surrogate?.

Since the coupling variables are replaced with surrogate models in MDA, the surrogate models of the coupling variables are coupled together over time. During the analysis, the surrogate model prediction uncertainty of one coupling surrogate model (for a given realization of input random variables) will be propagated to the uncertainty in the prediction of another coupling surrogate model. As a result, the uncertainty in the output of any coupling surrogate model comes from prediction uncertainty of all coupling surrogate models. This makes GSA to decide which coupling surrogate model to refine very difficult. To separate the contribution of each coupling surrogate model uncertainty, we introduce a vector of auxiliary variables $Ujk(ti)=[Ujk(1)(ti), Ujk(2)(ti), …, Ujk(NjkL)(ti)]$ (i.e., CDF values) for each coupling variables $Ljk(ti)$ at each time instant $ti$. After introducing the auxiliary variables, the network given in Fig. 7 at $ti$ changes to the network shown in Fig. 8. In this network, a triangular node represents a functional node (i.e., the variable indicated by this node is a deterministic function of its parent nodes). For example, the functional node $l12(ti)$ is defined as
$l12(ti)=[l12(1)(ti), l12(2)(ti), …, l12(N12L)(ti)], where l12(p)(ti)=μl12(p)(ti)|l12(ti−1)+U12(p)(ti)σl12(p)(ti)|l12(ti−1), ∀p=1, 2, …, N12L;$
(39)
Based on this transformation, the $p*$th response of the $k*$th disciplinary at time instant $tik*$ becomes a function of the auxiliary variables as follows:
$Zk*(p*)(tik*)=gZk*(p*)new(c, Ujk(ti),UZk*(ti), ∀j, k=1, 2, …, Nd; i=1, …, ik*)$
(40)

where $Ujk(ti)$ are the auxiliary variables corresponding to the coupling variables $Ljk(ti)$ from the jth discipline to the $k$th discipline at $ti$, $gZk*(p*)new(⋅)$ is a new function defined through a network defined similar to Fig. 8, $UZk*(p*)(ti)$ is the auxiliary variable of the $p*$th response of the $k*$th discipline at $ti$, and $c$ includes all other deterministic variables. For the sake of illustration, we do not list all of the involved deterministic variables, which include the realizations of random variables, stochastic processes, and parameters of the Kriging models.

For example, considering Fig. 8, we have $Z2(1)(ti)$ as
$Z2(1)(ti)=gZ2new(c, U21(tj), U12(tj), UZ2(tj), ∀j=1, …, i)$
(41)
In Eq. (40), each auxiliary variable represents the surrogate model uncertainty of a coupling variable at a time instant. These auxiliary variables (i.e., CDF values) are statistically independent from each other. This enables us to separate the contributions of different coupling surrogates to the uncertainty of QoI, which is the $p*$th response of the $k*$th discipline identified in level 1 GSA. Based on Eq. (40), we analyze the contribution of each coupling surrogate to the uncertainty of $Zk*(p*)(tik*)$ by grouping the associated auxiliary variables together. For instance, the first-order sensitivity index of the surrogate model prediction of $Ljk(p)$ is computed by
$SUjk(p)I=VarUjk(p)(EUjk(p)∼(Zk*(p*)(tik*)|Ujk(p)))/Var(Zk*(p*)(tik*))$
(42)

where $Ujk(p)=[Ujk(p)(t1), Ujk(p)(t2), …, Ujk(p)(tik*)]$, $Ujk(p)∼$ represent the other auxiliary variables, $EUjk(p)∼(Zk*(p*)(tik*)|Ujk(p))$ is computed by performing uncertainty propagation using the network given in Fig. 8, and $Var(Zk*(p*)(tik*))$ is computed in Eq. (37).

With $SUjk(p)I$, the coupling surrogate that needs to be refinedis identified by
$[p, j, k]=argmaxj,k∈[1, 2, …, Nd]; p∈[1, …, NjkL]{SUjk(p)I}$
(43)

After the coupling surrogate model is identified, the new training point for the surrogate model is identified using the weighted mean square error criterion given in Eq. (38) and the samples of the coupling variables at time instant $tik*$. After the coupling surrogate models are updated, the distributions of the coupling variables are updated. The updated samples of the coupling variables are then used to perform level 2 GSA again. If $SL¯I, we go to the next step. Otherwise, level 3 GSA is performed again to further refine the coupling surrogate.

Until now, we have discussed all the implementation procedure of the proposed A-SMART framework. Next, we will use numerical examples to demonstrate the effectiveness of the proposed method.

## Numerical Examples

In this section, a mathematical example and an engineering example are used to illustrate the proposed method. In each example, the proposed method is compared with two other methods, namely, variance-minimization method and MCS. For MCS, 5 × 105 realizations are generated for random variables and stochastic processes and the Karhunen–Loève (KL) expansion method is used to generate samples for stochastic processes [11]. In the variance-minimization method, the surrogate models are refined by minimizing the maximum variance over iterations. In the Kriging surrogate, the zeroth order trend function is used.

### A Mathematical Example.

A multidisciplinary system with three disciplines as shown in Fig. 9 is employed as our first example. In this example, $X1=∅$, $X2=∅$, $X3=∅$, $Xs=[X1, X2]$, $Y1(t)=∅$, $Y2(t)=∅$, $Y3(t)=∅$, and $Ys(t)=[Y(t)]$. $X1$ and $X2$ are standard normal random variables and the initial conditions of $L21(t0)$, $L31(t0)$, and $L32(t0)$ are random variables follow Lognormal distributions $L21(t0)∼LN(11, 0.1)$, $L31(t0)∼LN(22, 0.1)$, and $L32(t0)∼LN(0.5, 0.02)$, respectively. $Y(t)$ is a Gaussian stochastic process with distribution given by $Y(t)∼N(10, 22)$ and correlation given by
$ρY(t1, t2)=exp (−((t2−t1)/0.2)2)$
(44)
The time-dependent failure probability of the system is defined as
$pf(t0, te)=Pr{(Z1(t1)<0∩Z2(t2)<0)∪Z3(t3)<0, ∃t1, t2, t3∈[t0, te]}$
(45)

where $t0=0$ and $te=1.5$. The time interval is discretized into 20 time instants in MDA.

We first perform time-dependent system reliability analysis using the proposed method. In the initial surrogate modeling, 35 training points are generated for each model (i.e., $L12$, $L21$, $L13$, $L31$, $L32$, $L23$, $Z1$, $Z2$, $Z3$). Table 1 gives the results comparison of the proposed method and MCS, which include the estimated $pf(t0, te)$, the number of function evaluations (NOF), and the relative estimation error ($ε%$). $ε%$ is computed by
$ε%=[|pf(t0, te)−pfMCS(t0, te)|/pfMCS(t0, te)]×100%$
(46)

where $pfMCS(t0, te)$ is the results estimated from MCS.

It shows that the proposed method is able to estimate the time-dependent system failure probability efficiently and accurately. For the 651 NOF of the proposed method, the NOFs are allocated as $NL12=72$, $NL13=124$, $NL21=37$, $NL23=35$, $NL31=35$, $NL32=180$, $NZ1=69$, $NZ2=58$, and $NZ3=41$, where $NLij$ means the NOF allocated to $Lij$. It implies that more computational resources are allocated to models $L13$, $L32$, and $Z1$ than the other models. Figure 10 shows the comparison of the convergence history of $pf(t0, te)$ from the proposed method and the variance minimization-based method with respect to NOF. It indicates that the proposed method converges much faster than the variance minimization-based method.

### A Compound Cylinder.

Our second example is a compound cylinder as shown in Fig. 11, which is modified from Ref. [1]. This compound cylinder is actually a two-component system with a single discipline. However, each component of the system has different analysis models and the outputs of one model will act as inputs of the other model. This is mathematically the same as the multidisciplinary system. This example can be extended to have multiple disciplines by including bending and modal analyses of the cylinder. In that case, FEA analysis needs to be employed. In this paper, for the sake of illustration, the model presented in Ref. [1] is used. Figure 12 shows the MDA models and variables of the MDA models. In this example, $X1=[a, bin, S1]$, $X2=[bout, c, S2, δ]$, $Xs=[E, ρ]$, $Y1(t)=∅$, $Y2(t)=∅$, and $Ys(t)=[p(t)]$. The QoI are the safety sates ($Z1(t)=[Z1(1)(t), Z1(2)(t)]$, $Z2(t)=[Z2(1)(t), Z2(2)(t)]$) of the cylinders with respect to the tangential stresses ($σ11(t), σ12(t)$) at the internal radius and external radius of the inner cylinder and their counterparts ($σ21(t), σ22(t)$) of the outer cylinder. The expressions of the stresses are given in Eqs. (47)(50)
$σ11(t)=−2L21(t)(bin−k1t)2(bin−k1t)2−a2+(c2+a2)(c2−a2)p(t)$
(47)
$σ12(t)=−L21(t)[(bin−k1t)2+a2](bin−k1t)2−a2+a2(c2+(bin−k1t)2)(c2−a2)p(t)$
(48)
$σ21(t)=L21(t)(bout+k2t)2+c2c2−(bout+k2t)2+a2[(bout+k2t)2+c2](c2−a2)(bout+k2t)2p(t)$
(49)
$σ22(t)=2(bout+k2t)2L21(t)c2−(bout+k2t)2+2a2p(t)c2−a2$
(50)

where $k1=5×10−3 in/yr$ and $k2=1.5×10−3 in/yr$, representing the effect of corrosion on the thickness of cylinders over time.

The random variables $a$ (in), $bin$ (in), $bout$ (in), $c$ (in), $δ$ (in), $ρ$ (psi), $S1$ (psi), $S2$ (psi), $E$ (psi), and $L210$ (psi) are assumed to follow Gaussian distributions with mean values of 7, 11, 11, 13, 4 × 10−3, 0.3, 1 × 107, 1 × 107, 3 × 107, and 2 × 103, respectively, and standard deviations of 0.1, 0.1, 0.1, 0.1, 1 × 10−4, 0.01, 1 × 106, 1 × 106, 3 × 106, and 10, respectively. $p(t)$ (psi) is a Gaussian stochastic process with mean 4.2 × 106, standard deviation 3 × 105, and correlation given by
$ρp(t1, t2)=exp (−((t2−t1)/ζp)2)$
(51)

where $ζp=0.5 years$ is the correlation length of pressure.

The time-dependent failure probability of the system is defined as
$pf(t0, te)=Pr{Z1(1)(t1)<0∪Z1(2)(t2)<0∪Z2(1)(t3)<0∪Z2(2)(t4)<0, ∃t1, t2, t3, t4∈[t0, te]}$
(52)

where $t0=0 years$ and $te=5 years$.

We perform time-dependent system reliability analysis for the compound cylinder. In the initial surrogate modeling, 13 training points are generated for each model. Table 2 gives the results comparison between the proposed method and MCS. Table 3 lists the allocation of NOFs to different simulation models. Following that, Fig. 13 shows the comparison of convergence history of different methods. Similar conclusions are obtained as that from the first example.

### Discussion.

The above two numerical examples demonstrate that the proposed method can effectively allocate the computational resources to efficiently and accurately perform time-dependent MDRA. In the proposed method, only temporal variability is considered in MDRA.

As discussed in our motivation problems (i.e., Sec. 2.1), the coupling variables and response variables may have both spatial and temporal variabilities. When the spatial variability is also considered, the coupling and response variables will be high-dimensional, which brings significant challenge to the surrogate modeling. In that situation, dimension reduction methods need to be employed to reduce the dimension and thus make the surrogate modeling possible. For instance, the singular value decomposition or proper orthogonal decomposition methods can be adopted first to map the high-dimensional coupling and response variables into low-dimensional variables in latent space. Surrogate models are then built in the low-dimensional latent space. Such a method has been presented in Ref. [25].

Based on the dimension reduction methods, the proposed adaptive surrogate modeling framework can be employed to reduce the required computational effort for surrogate modeling to achieve an accurate prediction of the system failure probability. The integration of the proposed adaptive surrogate modeling method with the dimension reduction method, however, is not straightforward. The convergence criterion and the three-level GSA methods need to be integrated with dimension reduction methods, and need to be further extended to account for the spatial variability. Integration of the proposed adaptive surrogate modeling method with the method presented in Ref. [25] will be studied in our future work.

## Conclusion

This paper developed a new A-SMART framework for time-dependent multidisciplinary system reliability analysis. Surrogate models are built separately for each coupling variable and response of each disciplinary analysis. Computational resources are then adaptively allocated to the surrogate models to improve the efficiency and accuracy of reliability analysis, through the integration of three modules: initialization, uncertainty propagation, and three-level GSA. A mathematical example and an engineering application example demonstrated that the proposed method can effectively allocate the computational resources to efficiently and accurately perform time-dependent MDRA.

In the proposed framework, sampling-based approaches are employed to perform uncertainty propagation and global sensitivity analysis. For time-dependent MDRA over a very long time period, the sampling-based approaches might be computationally expensive even if it is much cheaper than using the original simulation models. In addition, the correlation between the predictions of individual samples as discussed in Ref. [28] is not considered in the developed framework. Accounting for the correlation effect will further improve the accuracy and efficiency of the proposed method.

## Funding Data

• Air Force Office of Scientific Research (Grant No. FA9550-15-1-0018).

## References

1.
Du
,
X.
,
Guo
,
J.
, and
Beeram
,
H.
,
2008
, “
Sequential Optimization and Reliability Assessment for Multidisciplinary Systems Design
,”
Struct. Multidiscip. Optim.
,
35
(
2
), pp.
117
130
.
2.
Liang
,
C.
,
,
S.
, and
Sankararaman
,
S.
,
2015
, “
Stochastic Multidisciplinary Analysis Under Epistemic Uncertainty
,”
ASME J. Mech. Des.
,
137
(
2
), p.
021404
.
3.
Felippa
,
C. A.
,
Park
,
K.
, and
Farhat
,
C.
,
2001
, “
Partitioned Analysis of Coupled Mechanical Systems
,”
Comput. Methods Appl. Mech. Eng.
,
190
(
24
), pp.
3247
3270
.
4.
Rangavajhala
,
S.
,
Sura
,
V. S.
,
Hombal
,
V. K.
, and
,
S.
,
2011
, “
Discretization Error Estimation in Multidisciplinary Simulations
,”
AIAA J.
,
49
(
12
), pp.
2673
2683
.
5.
Du
,
X.
, and
Chen
,
W.
,
2005
, “
Collaborative Reliability Analysis Under the Framework of Multidisciplinary Systems Design
,”
Optim. Eng.
,
6
(
1
), pp.
63
84
.
6.
,
S.
, and
Smith
,
N.
,
2006
, “
Efficient First-Order Reliability Analysis of Multidisciplinary Systems
,”
Int. J. Reliab. Saf.
,
1
(
1–2
), pp.
137
154
.
7.
Sankararaman
,
S.
, and
,
S.
,
2012
, “
Likelihood-Based Approach to Multidisciplinary Analysis Under Uncertainty
,”
ASME J. Mech. Des.
,
134
(
3
), p.
031008
.
8.
Jiang
,
Z.
,
Li
,
W.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2015
, “
A Spatial-Random-Process Based Multidisciplinary System Uncertainty Propagation Approach With Model Uncertainty
,”
ASME J. Mech. Des.
,
137
(
10
), p.
101402
.
9.
Jiang
,
Z.
,
Chen
,
S.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2016
, “
Reduction of Epistemic Model Uncertainty in Simulation-Based Multidisciplinary Design
,”
ASME J. Mech. Des.
,
138
(
8
), p.
081403
.
10.
DeCarlo
,
E. C.
,
Smarslok
,
B. P.
, and
,
S.
,
2016
, “
Segmented Bayesian Calibration of Multidisciplinary Models
,”
AIAA J.
,
54
(
12
), pp.
3727
3741
.
11.
Hu
,
Z.
, and
,
S.
,
2016
, “
A Single-Loop Kriging Surrogate Modeling for Time-Dependent Reliability Analysis
,”
ASME J. Mech. Des.
,
138
(
6
), p.
061406
.
12.
Wang
,
Z.
, and
Wang
,
P.
,
2012
, “
A Nested Extreme Response Surface Approach for Time-Dependent Reliability-Based Design Optimization
,”
ASME J. Mech. Des.
,
134
(
12
), p.
121007
.
13.
Singh
,
A.
,
Mourelatos
,
Z. P.
, and
Li
,
J.
,
2010
, “
Design for Lifecycle Cost Using Time-Dependent Reliability
,”
ASME J. Mech. Des.
,
132
(
9
), p.
091008
.
14.
Singh
,
A.
,
Mourelatos
,
Z.
, and
Nikolaidis
,
E.
,
2011
, “
Time-Dependent Reliability of Random Dynamic Systems Using Time-Series Modeling and Importance Sampling
,”
SAE Int. J. Mater. Manuf.
,
4
(
1
), pp.
929
946
.
15.
Zhu
,
Z.
,
Hu
,
Z.
, and
Du
,
X.
,
2015
, “Reliability Analysis for Multidisciplinary Systems Involving Stationary Stochastic Processes,”
ASME
Paper No. DETC2015-46168.
16.
Hombal
,
V.
, and
,
S.
,
2011
, “
Bias Minimization in Gaussian Process Surrogate Modeling for Uncertainty Quantification
,”
Int. J. Uncertainty Quantif.
,
1
(
4
), pp.
321
349
.
17.
Bichon
,
B. J.
,
Eldred
,
M. S.
,
Swiler
,
L. P.
,
,
S.
, and
McFarland
,
J. M.
,
2008
, “
Efficient Global Reliability Analysis for Nonlinear Implicit Performance Functions
,”
AIAA J.
,
46
(
10
), pp.
2459
2468
.
18.
Zhu
,
Z.
, and
Du
,
X.
,
2016
, “
Reliability Analysis With Monte Carlo Simulation and Dependent Kriging Predictions
,”
ASME J. Mech. Des.
,
138
(
12
), p.
121403
.
19.
Dubourg
,
V.
,
Sudret
,
B.
, and
Deheeger
,
F.
,
2013
, “
Metamodel-Based Importance Sampling for Structural Reliability Analysis
,”
Probab. Eng. Mech.
,
33
, pp.
47
57
.
20.
Wang
,
Z.
, and
Wang
,
P.
,
2014
, “
A Maximum Confidence Enhancement Based Sequential Sampling Scheme for Simulation-Based Design
,”
ASME J. Mech. Des.
,
136
(
2
), p.
021006
.
21.
Wang
,
Z.
, and
Chen
,
W.
,
2016
, “
Time-Variant Reliability Assessment Through Equivalent Stochastic Process Transformation
,”
Reliab. Eng. Syst. Saf.
,
152
, pp.
166
175
.
22.
Liang
,
C.
, and
,
S.
,
2016
, “
Stochastic Multidisciplinary Analysis With High-Dimensional Coupling
,”
AIAA J.
,
54
(
4
), pp.
1209
1219
.
23.
Culler
,
A. J.
, and
McNamara
,
J. J.
,
2011
, “
Impact of Fluid-Thermal-Structural Coupling on Response Prediction of Hypersonic Skin Panels
,”
AIAA J.
,
49
(
11
), pp.
2393
2406
.
24.
Smarslok
,
B.
,
Culler
,
A.
, and
,
S.
,
2012
, “Error Quantification and Confidence Assessment of Aerothermal Model Predictions for Hypersonic Aircraft,”
AIAA
Paper No. 2012-1817.
25.
Hu
,
Z.
, and
,
S.
,
2017
, “
A Surrogate Modeling Approach for Reliability Analysis of a Multidisciplinary System With Spatio-Temporal Output
,”
Struct. Multidiscip. Optim.
,
56
(
3
), pp.
553
569
.
26.
Du
,
X.
, and
Chen
,
W.
,
2002
, “
Efficient Uncertainty Analysis Methods for Multidisciplinary Robust Design
,”
AIAA J.
,
40
(
3
), pp.
545
552
.
27.
Echard
,
B.
,
Gayton
,
N.
, and
Lemaire
,
M.
,
2011
, “
AK-MCS: An Active Learning Reliability Method Combining Kriging and Monte Carlo Simulation
,”
Struct. Saf.
,
33
(
2
), pp.
145
154
.
28.
Hu
,
Z.
, and
,
S.
,
2016
, “
Global Sensitivity Analysis-Enhanced Surrogate (GSAS) Modeling for Reliability Analysis
,”
Struct. Multidiscip. Optim.
,
53
(
3
), pp.
501
521
.
29.
Wang
,
Z.
, and
Wang
,
P.
,
2015
, “
An Integrated Performance Measure Approach for System Reliability Analysis
,”
ASME J. Mech. Des.
,
137
(
2
), p.
021406
.
30.
Hu
,
Z.
, and
Du
,
X.
,
2013
, “
Time-Dependent Reliability Analysis With Joint Upcrossing Rates
,”
Struct. Multidiscip. Optim.
,
48
(
5
), pp.
893
907
.
31.
Arendt
,
P. D.
,
Apley
,
D. W.
, and
Chen
,
W.
,
2013
, “
Objective-Oriented Sequential Sampling for Simulation Based Robust Design Considering Multiple Sources of Uncertainty
,”
ASME J. Mech. Des.
,
135
(
5
), p.
051005
.