Fréchet-derivative-based global sensitivity analysis of the physical random function model of ground motions

Physical random function model of ground motions

The StoModel studied in this research is based on the source-path-site mechanisms. Specifically, this model consists of two physical models: the amplitude spectrum model $$ A_R(\boldsymbol{\Theta}, \omega) $$ and the phase spectrum model $$ \Phi_R(\boldsymbol{\Theta}, \omega) $$. Then, at a point on the surface of a specific local engineering site with the epicentral distance of $$ R $$, the acceleration of ground motion can be generated by^[23]

where the amplitude spectrum model $$ A_R(\boldsymbol{\Theta}, \omega) $$ is given by

where $$ K=10^{-5} $$ s/km is the attenuation parameter, and the phase spectrum model reads

In the StoModel, $$ A_0 $$ is the amplitude parameter of the source, $$ \tau $$ is the Brune source parameter describing the decay process of the fault rupture, $$ \zeta_g $$ is the equivalent damping ratio of the site, and $$ \omega_g $$ is the equivalent predominate circular frequency of the site. On the consideration of physical interpretations of the parameters $$ \{A_0, \tau, \zeta_g, \omega_g\} $$, it is appropriate to denote $$ \boldsymbol{\Theta} = [A_0, \tau, \zeta_g, \omega_g] $$ as the basic random source that characterizes the uncertainty of ground motions. It should be emphasized that the remaining parameters $$ \{a, b, c, d, R\} $$, which have a crucial impact on the phase spectrum model, may also possess randomness. However, to clearly illustrate the proposed method in this paper, these parameters are assumed to be deterministic by setting $$ a=1.02 $$, $$ b=403 $$, $$ c=1.89 $$, $$ d=130 $$, and $$ R = 20 $$ km, as referenced from Wang & Li (2012)^[24]. The reader interested in the randomness of $$ \{a, b, c, d, R\} $$ is directed to Ding et al. (2018, 2022)^[25,26] for further information.

According to the site classification recommended in the Chinese code for seismic design of buildings (GB 50011-2010)^[44], the marginal probability density functions (PDFs) of $$ \{A_0, \tau, \zeta_g, \omega_g\} $$ are estimated^[32] using a database of 4438 seismic ground motions from the Pacific Earthquake Engineering Research Center (PEER). The assumption of Independence is made by $$ \{A_0, \tau, \zeta_g, \omega_g\} $$, and the probabilistic information of $$ \{A_0, \tau, \zeta_g, \omega_g\} $$ is summarized in Table 1. The PDFs of $$ \{A_0, \tau, \zeta_g, \omega_g\} $$ for Site I, Site II, Site III, and Site IV are shown in Figure 1.

Figure 1. PDFs of parameters of physical random function model of ground motions according to the site classification in the Chinese design code (GB 50011-2010)^[44].

It is evident that the distribution parameters of basic random sources vary greatly for different site classes. This variation can be attributed to both the physical characteristics of different site classes and the statistical uncertainty originating from the earthquake ground motions. In other words, the distribution parameters derived from the analysis of $$ 4438 $$ seismic ground motions^[32] may possess epistemic uncertainty, as demonstrated in the study by Li & Liu (2015)^[33] where different distribution parameters are estimated from the records of 2008 Wenchuan earthquakes. To this regard, it is valuable to study how these uncertainties may affect the uncertainty of structural responses, such as maximum inter-story drift angle, top displacement, etc.

In Table 1, the PDFs of Lognormal distribution and Gamma distribution are given by

Besides, for the sake of simplicity, the distribution parameters in Table 1 are numbered in order as follows: $$ \{a_{1}, b_{1}\} $$ for $$ A_0 = \Theta_1 $$, $$ \{a_{2}, b_{2}\} $$ for $$ \tau = \Theta_2 $$, $$ \{a_{3}, b_{3}\} $$ for $$ \zeta_g = \Theta_3 $$, and $$ \{a_{4}, b_{4}\} $$ for $$ \omega_g = \Theta_4 $$, for Site I to Site IV.

Uncertainty propagation via the probability density evolution method

In this section, We introduce the basic theory and numerical algorithm of the PDEM^[45], which is adopted to estimate the PDF of the QoI.

Without loss of generality, Let us consider a MDOF structure with the equation of motion given by:

where $$ \bf{M} $$ and $$ \bf{C} $$ are mass and damping matrices of $$ n \times n $$, respectively, $$ \bf{f} $$ is a linear or nonlinear force vector of $$ n \times 1 $$. The displacement, velocity, and acceleration of $$ n\times 1 $$ are denoted by $$ \boldsymbol{X} $$, $$ \dot{\boldsymbol{X}} $$, and $$ \ddot{\boldsymbol{X}} $$, respectively. $$ \boldsymbol{\xi} $$ is the loading influence matrix of $$ n \times 1 $$. In this paper, the uncertainty denoted by $$ \boldsymbol{\Theta} $$ from the earthquake ground motions is taken into account, while the structure is considered to be deterministic.

Let $$ X(t) $$ be the QoI that is a function of structural responses, i.e., $$ X(t) = g(\boldsymbol{X}, \dot{\boldsymbol{X}}, \ddot{\boldsymbol{X}}) $$ where $$ g(\cdot) $$ is a linear or nonlinear mapping. For instance, The maximum inter-story drift or the top displacement may be of interest in seismic reliability assessment for building structures. In this paper, the QoI is defined as $$ X = \max_{t=0}^T\{\vert X_{\rm top}(t)\vert\} $$, where $$ X_{\rm top}(t) $$ is the time history of the structural top displacement. For most well-posed engineering systems, $$ X(t) $$ is a unique function that exists as a function of $$ \boldsymbol{\Theta} $$, which can be expressed by

where $$ h = \partial H / \partial t $$ stands for the generalized velocity. Then, based on the principle of probability preservation^[46], the joint PDF of $$ (X, \boldsymbol{\Theta}) $$ is governed by^[2]

which is referred to as the generalized density evolution equation (GDEE). Finally, the PDF of QoI can be calculated by integrating $$ \boldsymbol{\theta} $$ after solving Equation (7), i.e.,

where $$ \Omega_{\boldsymbol{\Theta}} $$ is the sample space of $$ \boldsymbol{\Theta} $$.

In general, the numerical algorithm of the PDEM consists of the following four steps:

Step 1.1.  Generation of representative points. Denote $$ \Omega_{\boldsymbol{\Theta}} $$ be partitioned into $$ N $$ disjoint subdomains satisfying $$ \cup_{q=1}^N \Omega_q = \Omega_{\boldsymbol{\Theta}} $$ and $$ \Omega_p \cap \Omega_q = \emptyset $$ for $$ p\ne q $$ and $$ p, q=1, 2, \cdots, N $$. For each subdomain $$ \Omega_q $$, select one representative point $$ \boldsymbol{\theta}_q \in \Omega_q $$ and calculate its assigned probability $$ P_q $$ defined by

The way to partition the sample space can be referred to Chen et al. (2009)^[47]. Besides, to minimize the point discrepancy of representative points, the GF-discrepancy minimization strategy^[48] is adopted in this work.

Step 1.2.  For each $$ \boldsymbol{\Theta} = \boldsymbol{\theta}_q $$, generate the stochastic ground motion $$ a_{R}(\boldsymbol{\theta}_q, t) $$ by Equations (1) to (3). Then, solve Equation (5) and Equation (6) to obtain the generalized velocity $$ h(\boldsymbol{\theta}_q, t) $$.

Step 1.3.  For each $$ \boldsymbol{\Theta} = \boldsymbol{\theta}_q $$, solve Equation (7) in a discretized version, i.e.,

with the initial condition $$ p_{X\boldsymbol{\Theta}}^{(q)}(x, \boldsymbol{\theta}_q, t_0) = \delta_{\rm D}[x - x_0] P_q $$ where $$ \delta_{\rm D}[\cdot] $$ is the Dirac delta function. Equation (10) is a typical partial differential equation that can be numerically solved via finite difference methods^[45].

Step 1.4.  Assemble the results in Step 1.3, i.e., $$ p_X(x, t) = \sum_{q = 1}^N p_{X\boldsymbol{\Theta}}^{(q)}(x, \boldsymbol{\theta}_q, t) $$.

Uncertainty propagation via the change of probability measure

The aforementioned PDEM is available only if the input PDF is precisely determined. In other words, when the input PDF denoted by $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ is already known, the corresponding output PDF denoted by $$ p_{X}^{(1)}(x, t) $$ can be accurately estimated by the PDEM. What if the PDF of $$ \boldsymbol{\Theta} $$ is changed from $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ to $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$? To obtain the PDF of $$ X $$ denoted by $$ p_X^{(2)}(x, t) $$ in terms of $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$, one may complement the probability density evolution analysis again, which undoubtedly requires another loop of deterministic analyses. To avoid additional model evaluations, the COM^[4] is briefly introduced in this section. The combination of PDEM and COM is essential for a quick Fre-GSA in Section 2.4.

The backbone of the COM is based on the Radon-Nikodým theorem, which ensures that

where $$ \circ $$ means an operator on a function and $$ \mathcal{T} $$ is the Radon-Nikodým derivative defined by

where $$ P_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$ and $$ P_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ are the cumulative distribution functions (CDFs) in terms of $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$ and $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$, respectively. Note that Equation (11) means one can obtain $$ p_X^{(2)}(x, t) $$ directly by using the Radon-Nikodým derivative $$ \mathcal{T} $$ rather than adding additional deterministic analyses.

For some simple stochastic systems, the analytical formula of Radon-Nikodým derivative can be found in Chen & Wan (2019)^[4]. Nevertheless, it is always impossible to obtain an exact expression of Radon-Nikodým derivative for complex and nonlinear stochastic systems, but the COM can be numerically accomplished with the aid of the PDEM.

The numerical algorithm of the PDEM-COM is summarized as follows:

Step 2.1.  Complete one round of probability density evolution analysis via PDEM introduced in Section 2.2. Store the point set $$ \mathcal{M}^{(1)} = \{\boldsymbol{\theta}_q^{(1)}, P_q^{(1)}\}_{q = 1}^N $$ and the corresponding generalized velocity $$ h(\boldsymbol{\theta}_q^{(1)}, t) $$, where $$ \boldsymbol{\theta}_q^{(1)} $$ is the $$ q $$-th representative point with respect to the PDF $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$.

Step 2.2.  Considering the input PDF is changed from $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ to $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$, recalculate the assigned probability by

where $$ \Omega_q^{(1)} $$'s are the Voronoi cells determined by $$ \boldsymbol{\theta}_q $$'s. This generates a new point set $$ \mathcal{M}^{(2)} = \{\boldsymbol{\theta}_q^{(1)}, P_q^{(2)}\}_{q = 1}^N $$. By doing so, the location of the $$ q $$-th representative point is unchanged (still $$ \boldsymbol{\theta}_q^{(1)} $$), which means the generalized velocity $$ h(\boldsymbol{\theta}_q^{(1)}, t) $$ can be reused.

Step 2.3.  Solve the GDEE in Equation (10) with a new initial condition $$ p_{X\boldsymbol{\Theta}}^{(q)}(x, \boldsymbol{\theta}_q^{(1)}, t_0) = \delta_{\rm D}[x - x_0] P_q^{(2)} $$ and then assemble the results.

It should be emphasized that the accuracy of the PDEM-COM depends on whether the support of $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ mostly covers that of $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$. When the supports of $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ and $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$ are largely coincident, the PDEM-COM provides relatively accurate results. Nonetheless, if the supports of $$ p_{\boldsymbol{\Theta}}^{(1)}(\boldsymbol{\theta}) $$ and $$ p_{\boldsymbol{\Theta}}^{(2)}(\boldsymbol{\theta}) $$ are exclusive, for example, $$ \Omega_{\Theta}^{(1)} = (-1, 0) $$ while $$ \Omega_{\Theta}^{(2)} = (0, 1) $$, the accuracy of the PDEM-COM may immediately collapse because there is no universally perfect method. To address this issue, the accuracy of the PDEM-COM can be remarkably improved by an augmenting method. For more details, please refer to Wan et al. (2023)^[49].

Fréchet-derivative-based global sensitivity analysis

The Fre-GSA provides a quantitative approach to identify the most influential variables of input for stochastic systems. In this analysis, a Fre-GSI is calculated, and its parametric form is defined by^[43]

and the $$ j $$-th Fre-GSI reads

where $$ \xi_j $$ is the $$ j $$-th distribution parameter of the $$ i $$-th random variable $$ \Theta_i $$, and the norm term in the denominator is defined as $$ \Vert p \Vert = \frac{1}{2} \int \vert p \vert $$. It should be emphasized that Equation (15) holds on the assumption that $$ \xi_j $$'s are independent. Specifically, in the studied StoModel, $$ \Theta_i $$ is associated with $$ \boldsymbol{\xi} = \left\{a_{i}, b_{i}\right\} $$ where $$ i=1, 2, 3, 4 $$ for Site I, Site II, Site III, and Site IV.

For the $$ j $$-th Fre-GSI, the corresponding IM is given by^[50]

which theoretically satisfies that $$ 0\le \mathcal{S}_j \le 1 $$.

The $$ j $$-th Fre-GSI can be calculated via the PDEM-COM in Section 2.3, which mainly consists of three steps:

Step 3.1.  Firstly, calculate $$ p_X(x, t;\boldsymbol{\xi}) $$ with respect to $$ p_{\boldsymbol{\Theta}}(\boldsymbol{\theta} ; \boldsymbol{\xi}) $$ via the PDEM, in which the point set $$ \mathcal{M} = \{\boldsymbol{\theta}_q, P_q\}_{q = 1}^N $$ and the corresponding generalized velocity $$ h(\boldsymbol{\theta}_q, t) $$ are stored.

Step 3.2.  Calculate $$ p_X(x, t;\boldsymbol{\xi}+\boldsymbol{e}_j\Delta{\xi_j}) $$ and $$ p_X(x, t;\boldsymbol{\xi}-\boldsymbol{e}_j{\xi_j}) $$ in terms of $$ p_{\boldsymbol{\Theta}}(\boldsymbol{\theta} ; \boldsymbol{\xi} + \boldsymbol{e}_j\Delta{\xi_j}) $$ and $$ p_{\boldsymbol{\Theta}}(\boldsymbol{\theta} ; \boldsymbol{\xi} - \boldsymbol{e}_j\Delta{\xi_j}) $$, respectively, where $$ \boldsymbol{e}_j $$ is a selection vector whose elements are zeros except its $$ j $$-th location equal to one, and $$ \Delta{\xi} $$ is a small perturbation, e.g., $$ \Delta{\xi} = {0.01} $$ for $$ \xi_j = {0} $$ and $$ \Delta{\xi_j} = 0.01{\xi_j} $$ for $$ {\xi_j} \ne {0} $$. Note that this step can be speedily accomplished by adopting the PDEM-COM, as mentioned in Section 2.3.

Step 3.3. Approximate the $$ j $$-th Fre-GSI via a central difference scheme, i.e.,

where the norm term in the denominator can be numerically or analytically computed^[4].

Acknowledgments

The financial support provided by the National Natural Science Foundation of China (NSFC Grant No. 52208206) and the Fundamental Research Funds for the Central Universities (Grant Nos. G2022KY05103) is highly appreciated.

Authors' contributions

Conceptualization, Investigation, Methodology, Formal analysis, Software, Writing - original draft, and Writing - review & editing: Wan Z

Resources, Funding acquisition, and Data curation: Tao W, Ding Y, Xin L

Availability of data and materials

Some or all of the data, models, or code generated or used during the study are available from the author upon request.

Financial support and sponsorship

The National Natural Science Foundation of China (NSFC Grant No. 52208206);

The Fundamental Research Funds for the Central Universities (Grant Nos. G2022KY05103).

Conflicts of interest

All authors declared that there are no conflicts of interest.

Ethical approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Copyright

© The Author(s) 2023.