计算力学学报

随机变量	物理意义	概率分布	均值	变异系数
m_p	主质量	对数正态	1	0.1
m_s	从质量	对数正态	0.01	0.1
k_p	主刚度	对数正态	1	0.2
k_s	次刚度	对数正态	0.01	0.2
ζ_p	主阻尼比	对数正态	0.05	0.4
ζ_s	次阻尼比	对数正态	0.02	0.5
F_s	次弹簧承载力	对数正态	20	0.1
S₀	白噪声的强度	对数正态	100	0.1

随机变量	物理意义	概率分布	均值	变异系数
m_p	主质量	对数正态	1	0.1
m_s	从质量	对数正态	0.01	0.1
k_p	主刚度	对数正态	1	0.2
k_s	次刚度	对数正态	0.01	0.2
ζ_p	主阻尼比	对数正态	0.05	0.4
ζ_s	次阻尼比	对数正态	0.02	0.5
F_s	次弹簧承载力	对数正态	20	0.1
S₀	白噪声的强度	对数正态	100	0.1

随机变量	分布类型	均值	变异系数
E_c/Pa	对数正态分布	3×10¹⁰	0.1
m₁、m₂、m₃/kg	对数正态分布	1.2×10⁵	0.1
m₄/kg	对数正态分布	1.3×10⁵	0.1
m₅/kg	对数正态分布	1.1×10⁵	0.1
m₆、m₇/kg	对数正态分布	1.0×10⁵	0.1
m₈、m₉/kg	对数正态分布	1.1×10⁵	0.1
m₁₀/kg	对数正态分布	0.5×10⁵	0.1

随机变量	分布类型	均值	变异系数
E_c/Pa	对数正态分布	3×10¹⁰	0.1
m₁、m₂、m₃/kg	对数正态分布	1.2×10⁵	0.1
m₄/kg	对数正态分布	1.3×10⁵	0.1
m₅/kg	对数正态分布	1.1×10⁵	0.1
m₆、m₇/kg	对数正态分布	1.0×10⁵	0.1
m₈、m₉/kg	对数正态分布	1.1×10⁵	0.1
m₁₀/kg	对数正态分布	0.5×10⁵	0.1

一种改进的功能函数概率密度估计的等效期望法

PDF下载

邹长星 , 翁叶耀 , 张玄一 , 赵衍刚

计算力学学报 | 研究论文 2025,42(5): 846-851

收起

计算力学学报 | 研究论文 2025, 42(5): 846-851

一种改进的功能函数概率密度估计的等效期望法

全屏

邹长星, 翁叶耀, 张玄一, 赵衍刚

作者信息

北京工业大学　城市建设学部，北京　100124

赵衍刚^*（1963-），男，博士，教授（E-mail：zhao@kanagawa-u.ac.jp）.

An improved equivalent expectation method for evaluate probability density of performance functions

Changxing ZOU, Yeyao WENG, Xuanyi ZHANG, Yangang ZHAO

Affiliations

Faculty of Architecture, Civil Engineering, Beijing University of Technology, Beijing 100124, China

出版时间: 2025-10-28 doi: 10.7511/jslx20240202002

文章导航

摘要

收起

功能函数的概率分布计算是不确定性量化及可靠度设计中的核心问题，最近提出的等效期望法是解决该问题的有效途径。本文提出了一种改进的等效期望法。针对既有等效期望法中的辅助变量取值对其概率分布有重要影响的问题，提出了辅助随机变量标准差系数的合理计算公式，得到了更为精确的辅助函数的概率分布。同时，针对既有等效期望法中由两个辅助函数的概率分布计算功能函数的概率分布时存在的计算精度问题，本文仅用一个辅助函数推导了功能函数概率分布的计算公式，提出了一种由辅助函数到功能函数的概率分布精确理论变换，使得功能函数概率分布的计算更为精确。最后通过三个算例验证了该方法的有效性和准确性。结果表明，该方法适用于高维非线性或隐式功能函数的概率分布近似计算。

关键词

结构安全 / 结构可靠度 / 功能函数 / 随机变量 / 概率密度函数

Abstract

收起

The calculation of the probability distribution of performance functions is a core issue in uncertainty quantification and reliability design, and the recently proposed equivalent expectation method (EEM) is an effective way to solve this problem. This paper proposes an improved EEM. Putting forward an empirical calculation formula for the standard deviation coefficient of auxiliary random variablesand obtaining a more accurate probability distribution of the auxiliary function. Meanwhile, aiming at the accuracy issue in calculating the probability distribution of theperformance function is proposed, the calculation formula for the PDF of the performance function is derived using only one auxiliary function. In the process of calculating the PDF, proposing an exact theoretical transformation of probability distribution from auxiliary functions to performance functions is proposed, resulting in a more accurate PDF of the performance function. Finally, the effectiveness and accuracy of the method are verified through three numerical examples. The results indicate that this method is suitable for computing the probability distribution of high-dimensional nonlinear or implicit performance functions.

Key words

structural safety / structural reliability / performance function / random variable / probability density

引用本文

邹长星, 翁叶耀, 张玄一, 赵衍刚. 一种改进的功能函数概率密度估计的等效期望法. 计算力学学报, 2025 , 42 (5) : 846 -851 . DOI: 10.7511/jslx20240202002

Changxing ZOU, Yeyao WENG, Xuanyi ZHANG, Yangang ZHAO. An improved equivalent expectation method for evaluate probability density of performance functions[J]. Chinese Journal of Computational Mechanics, 2025 , 42 (5) : 846 -851 . DOI: 10.7511/jslx20240202002

正文

收起

1　引言

收起

功能函数或模型输出的概率分布的计算是正向不确定性量化^[1]中的核心问题。如在基于可靠性的设计优化^[2]中，为了降低计算成本，需要计算功能函数的分位数，即功能函数的累积分布函数（CDF）的反函数。在敏感性分析^[3]中，使用Borgonovo指数评估敏感性需要模型输出的概率密度函数（PDF）。在实际复杂结构系统的设计中，工程师希望确定对应于不同失效概率的结构阻力^[4]。此外，失效概率仅显示功能函数小于或等于零的统计信息。而概率分布作为功能函数统计信息的完备表达为结构安全评估提供了比失效概率更有用的信息。因此，计算功能函数的概率分布非常重要。

基于概率守恒原理，文献[5，6]开发了概率密度演化方法（PDEM）来计算功能函数的PDF。在该方法中，需要使用有限差分法求解广义概率密度演化方程，计算难度和计算工作量均较大。概率变换方法（PTM）提供了实现该过程的替代方法。根据随机变量变换定理，Falsone等^[7]提出了概率变换方法（PTM），以获得线性结构响应的PDF。经过众多学者的研究，该方法扩展到解决具有随机参数和非高斯激励的问题。由于利用了线性系统的特征函数和PDF之间的关系，PTM仅适用于线性和特殊非线性系统。受PDEM和PTM的启发，Chen等^[8]提出了一种直接概率积分法（DPIM）来计算功能函数的概率分布，在这种方法中，功能函数的PDF表示为与Dirac delta函数相关的概率密度积分，并采用高斯函数对Dirac delta函数进行平滑。然而，高斯函数中平滑参数的选择对DPIM结果有显著影响^[9]，但是目前对最优平滑参数的研究还不够深入。

Cai等^[10]提出了一种计算功能函数概率分布的等效期望法（EEM）。相比于矩方法，等效期望法无需计算高阶矩，避免了计算高阶矩不准确的问题；相比于概率密度演化方法，等效期望法无需求解复杂的广义概率密度演化偏微分方程。在该方法中，当功能函数关于某一随机变量的单调反函数（MIF）可以得到时，概率分布可以精确地表示为与该变量和MIF的分布有关的期望值。而如果这样的MIF不存在或不容易获得时，则通过简单地添加辅助随机变量构造辅助函数来获得MIF，并将功能函数的概率分布表示为与辅助函数分布相关的期望，最后采用Sobol序列和点估计方法确定期望值。经研究发现，辅助随机变量标准差的选取对辅助函数概率分布的准确性有较大影响，但该方法并没有给出辅助随机变量标准差的合理取值，同时，该方法构造了两个辅助函数，并把功能函数的PDF表示为与辅助函数相关的期望，但在计算过程中采用的变换不是精准的变换，存在较大的计算误差，这限制该方法的计算精度。

因此，本文的主要目的就是在等效期望法（EEM）的基础上，通过对辅助随机变量的进一步研究来确定一个合适的辅助变量标准差，使得计算的辅助函数的概率分布更为精准；同时，仅使用一个辅助函数进行功能函数PDF计算公式的推导，而计算过程采用的变换是一个精确的变换，从而提高所提方法的计算精度。

2　等效期望估计的基本原理

收起

一般情况下，功能函数Z=G（X）的CDF可定义为

(1)

其中f_X（x）为随机变量构成的向量X={X₁，X₂，…，X_n}^T的联合概率密度函数。当随机变量X_i与其他输入随机变量X_i相关时，可以首先执行正态变换，如Rosenblatt变换^[11]、Nataf变换^[12]或高阶矩变换^[13,14]。

如果功能函数Z=G（X）中存在一个随机变量X_i，使得该随机变量对应的反函数的偏导数

（X_－_i，Z）/∂Z>0，即功能函数中存在单调反函数（MIF），则功能函数的CDF可表示为一个期望。

然而，工程实际中功能函数通常是复杂的、非线性的，甚至是隐式的，功能函数中不存在这样的一个随机变量使其对应的功能函数的反函数单调，即使存在这样的随机变量，功能函数的单调反函数的计算通常是很困难的。等效期望法中通过添加正态的辅助随机变量Θ构造了两个辅助函数，构造的辅助函数为

(2)

(3)

而辅助函数的MIF是比较容易获得的，Θ对应的反函数分别为Θ=H－G（X）和Θ=G（X）－K，其单调性可以轻松确定：

(4)

(5)

根据式（4，5），辅助函数H和K的PDF和CDF可确定为

(6)

(7)

(8)

(9)

其中F_Θ（•）和f_Θ（•）分别为Θ的CDF和PDF，假设Θ是一个正态分布的随机变量，平均值μ_Θ=0，标准偏差σ_Θ=σ_z，其中σ_z是功能函数的标准差。以上等式采用Sobol序列来计算期望值，通过计算发现，简单选取σ_Θ=σ_z计算得到的结果并不精准。

根据式（2，3），功能函数Z=（H+K）/2，功能函数Z的PDF和CDF可表示为一个由辅助函数H和K的联合PDF构造的积分，推导得到的最终结果为

(10)

其中u_i和ω_i分别是Gauss-Hermite积分的横坐标和权重，c（•，•）是一个由H和K的CDF拟合的二元copula密度函数，k_i由变换F_K（k_i）=Φ（u_i）计算得到，Φ为标准正态分布函数。需要特别说明的是，该变换需要使用式（8）进行计算，但解该方程存在一定的计算困难；另一种方法是根据矩原理，利用辅助函数K的前几阶矩的信息构造一个立方正态分布后获得k_i，但利用矩方法时高阶矩的计算并不是精确的，且存在一定的范围限制。同时，不管采用何种方法，F_K（k_i）的计算本身就不是准确的，所以，运用变换F_K（k_i）=Φ（u_i）计算得到的k_i的计算结果存在比较大的误差。

3　改进的等效期望法

收起

3.1　辅助随机变量标准差系数的确定

利用Sobol序列进行点估计，可以确定H的PDF和CDF为

(11)

(12)

式中N是采用Sobol序列点集的个数，x_i为由Sobol序列产生的样本点，G（x_i）为样本点代入功能函数后计算所得的功能函数值，如果涉及有限元分析则为结构响应。

在等效期望估计方法的研究中，为简单起见，辅助随机变量Θ~N（0，σ_Z），然而经过对辅助函数CDF尾部的计算发现，辅助随机变量标准差的选取对计算辅助函数概率分布的准确性有较大影响。根据式（6，7），辅助函数H的CDF和PDF的计算精度与功能函数值的分布范围有关，为提高计算精度，将辅助随机变量Θ的标准偏差修正为σ_Θ=kσ_Z，选取不同的标准差来计算CDF，并观察其尾部拟合效果，发现标准差系数k与功能函数的取值范围及标准差有关，通过大量试算，得到标准差系数k的经验计算公式。

(13)

式中maxG（x_i）为样本点代入功能函数后计算所得的最大的功能函数值，minG（x_i）为最小的功能函数值。

3.2　功能函数概率分布计算公式的推导

根据式（2），功能函数Z可以表示为

(14)

值得注意的是，在式（14）中，辅助函数H与辅助变量Θ是存在相关性的，根据功能函数Z的CDF的定义，可得

(15)

此时，功能函数Z的CDF实际上为一个变上限积分，对该变上限函数进行求导即可得到功能函数Z的PDF。

(16)

根据Sklar定理，H和Θ的联合PDF可以表示为

(17)

其中c（•，•）是一个由H和Θ的CDF拟合的二元copula密度函数。在这项研究中，高度稳健的TLL方法用于计算二元copula密度函数。TLL方法的详细算法参见文献[15]。功能函数Z的PDF可重新表示为

(18)

注意f_Θ（Θ）是一个均值为0，标准差为kσ_Z的正态分布的PDF，可以轻松地将其转换为一个权重函数exp（－u²/2），使其满足Gauss-Hermite型积分公式的使用条件，本研究采用七点Gauss-Hermite积分计算功能函数Z的PDF：

(19)

其中u_i和ω_i分别是Gauss-Hermite积分的横坐标和权重；F_Θ（Θ_i）=Φ（u_i），式中Φ为标准正态概率分布函数，F_Θ（•）也为一正态分布函数，该变换是一个精确的变换，对于每一个u_i，都有一精确的Θ_i与之对应。相比于等效期望法中的变换F_K（k_i）=Φ（u_i），该变换克服了其计算困难且不精确的缺点，使得该方法能够得到一个更为精确的结果。确定PDF后，对PDF进行积分即可得到功能函数的CDF：

(20)

4　数值算例与研究

收起

4.1　算例1　强非线性功能函数

为了讨论提出方法在计算强非线性功能函数概率分布时的适用性，本算例研究了一个受白噪声激励的两自由度系统，如图1所示。

副弹簧失效可靠性的功能函数可表示为

(21)

其中γ=m_s/m_p为质量比，ω_a=（ω_p+ω_s）/2和ζ_a=（ζ_p+ζ_s）/2分别为平均频率和平均阻尼比，θ_sr=（ω_p－ω_s）/ω_a是一个调谐参数。

根据式（13）计算所得标准差系数k=1.9，图2所示为辅助随机变量的标准差σ_Θ=1.9σ_Z时辅助函数H的PDF和对数坐标CDF。为了进行比较，还使用10⁷样本的MCS获得了辅助函数的概率分布。可以看出，采用本文方法得到的辅助函数的PDF和CDF与MCS结果非常吻合。

在得到比较精确的辅助函数的概率分布后，即可计算功能函数的概率分布。同时，功能函数的概率分布也可以通过使用和不使用功能函数的MIF的等效期望方法得到。图3所示为采用1024个Sobol点获得的功能函数的PDF和对数坐标CDF。为了进行比较，还使用了1024个Sobol点的EEM以及10⁷样本的MCS获得了功能函数的概率分布。可以看出，（1）采用本文方法得到的功能函数的PDF和CDF与MCS结果基本一致。（2）与EEM相比，本文提出的方法能得到更精确的尾部CDF，具有更高的计算精度和计算效率。（3）对于更高维度的非线性问题，由于该方法只需要得到功能函数的响应值，而与功能函数的维度无关。即使功能函数的维度继续增长，Sobol序列也能够生成对应的样本点。

4.2　算例2　涉及动态可靠性的隐式功能函数

本算例对非平稳地震动作用下的10自由度集中质量非线性框架结构进行研究，如图4所示。结构的几何尺寸为，高度h₁=4 m，h₂=h₃=…=h₁₀=3 m，立柱截面=500 mm×500 mm，梁截面的抗弯刚度为E_bI→∞。采用瑞利阻尼，即C=ζ₁M+ζ₂K，其中ζ₁=0.01 Hz，ζ₂=0.005 s，M和K分别为集总质量矩阵和初始刚度矩阵。此外，恢复力采用Bouc-Wen模型，

(22)

其中k为初始刚度，u为层间位移，z为滞回位移，满足

(23)

其中参数设置为α=0.01，A=1.2，β=1.4，γ=0.2，r=1。

在本算例中，施加在目标结构上的非平稳地面运动加速度由以下谱密度函数建模：

(24)

其中S（ω）是平稳情况下的双侧功率谱密度函数，并且使用Clough-Penzien谱^[16]，g（t）为调制函数：

(25)

其中t₁和t₂分别为地震激励的静止部分的开始时刻和结束时刻，T为地震激励的持续时间，c为衰减系数。在该算例中，参数取为t₁=0.8 s，t₂=7 s，T=17 s和c=0.35。

采用四阶龙格-库塔法对结构进行每时程地震动加速度的确定性动力响应分析。第一层位移的最大绝对值认为是评估该结构可靠性的极值。然后，对该建筑的可靠性进行研究。假设结构破坏定义为结构位移超过允许极限的极值。因此，功能函数为

(26)

计算所得标准差系数k=1.65，图5所示为辅助随机变量的标准差σ_Θ=1.65σ_Z时辅助函数H的PDF和对数坐标CDF。为了进行比较，还使用10⁶样本的MCS获得了辅助函数的概率分布。可以看出，采用本文方法得到的辅助函数的PDF和CDF与MCS结果基本一致。

根据表2中随机变量的概率分布信息，可以使用Sobol序列产生样本点，对结构进行有限元计算获得结构的动力响应。图6所示为采用1024个Sobol点获得的功能函数的PDF和对数坐标CDF。为了进行比较，还使用了2048个Sobol点的EEM以及10⁶样本的MCS获得了功能函数的概率分布。可以看出，（1）采用本文方法得到的功能函数的PDF和CDF与MCS结果基本一致。（2）与EEM相比，本文提出的方法能得到更精确的尾部CDF，具有更高的计算精度。

5　结论

收起

本文在等效期望法（EEM）的基础上，针对既有等效期望法中的辅助变量取值影响其概率分布计算结果准确性的问题，提出了标准差系数k的经验计算公式，从而确定了辅助随机变量的标准差，使得计算的辅助函数的概率分布与MCS结果非常吻合，相比于EEM中确定的辅助变量更为精准；同时，仅利用一个辅助函数进行了功能函数概率分布计算公式的推导，计算过程中利用辅助随机变量进行了精确的变换，克服了等效期望法中利用构造的辅助函数进行变换的计算缺点，提高了计算精度和计算效率，得到如下结论。

（1）与EEM相比，该方法只需对功能函数进行较少的计算次数即可获得与MCS基本一致的结果，该方法在计算功能函数概率分布时有较好的精度和效率；

（2）与EEM相比，通过提出的标准差系数经验计算公式，得到了更为恰当的辅助随机变量标准差，同时能够得到比较精确的辅助函数的PDF和CDF。

基金

收起

国家自然科学基金(52278135)

参考文献

收起

文献

收起

参考文献引证文献

排序方式：

[1]

Wang

Z Q

, Broccardo

. A novel active learning-based Gaussian process metamodelling strategy for estimating the full probability distribution in forward UQ analysis[J]. Structural Safety, 2020, 84: 101937.

[2]

Dubourg

, Sudret

, Bourinet

J M

. Reliability-based design optimization using Kriging surrogates and subset simulation[J]. Structural and Multidisciplinary Optimization, 2011, 44(5): 673-690.

[3]

Saltelli

. Sensitivity analysis for importance assessment[J]. Risk Analysis, 2002, 22(3): 579-590.

[4]

Zhang

X D

, Low

Y M

, Koh

C G

. Maximum entropy distribution with fractional moments for reliability analysis[J]. Structural Safety, 2020, 83: 101904.

[5]

, Chen

J B

. Probability density evolution method for dynamic response analysis of structures with uncertain parameters[J]. Computational Mechanics, 2004, 34(5): 400-409.

[6]

. Probability density evolution method: Background, significance and recent developments[J]. Probabilistic Engineering Mechanics, 2016, 44: 111-117.

[7]

Falsone

, Settineri

. Explicit solutions for the response probability density function of linear systems subjected to random static loads[J]. Probabilistic Engineering Mechanics, 2013, 33: 86-94.

[8]

Chen

G H

, Yang

D X

. Direct probability integral method for stochastic response analysis of static and dynamic structural systems[J]. Computer Methods in Applied Mechanics and Engineering, 2019, 357: 112612.

[9]

X L

, Chen

G H

, Cui

H C

, et al. Direct probability integral method for static and dynamic reliability analysis of structures with complicated performance functions[J]. Computer Methods in Applied Mechanics and Engineering, 2021, 374: 113583.

[10]

Cai

C H

, Zhao

Y G

, Lu

Z H

, et al. An equivalent expectation evaluation method for approximating the probability distribution of performance functions[J]. Structural Safety, 2022, 95: 102180.

[11]

Hohenbichler

, Rackwitz

. Non-normal dependent vectors in structural safety[J]. Journal of the Engineering Mechanics Division, 1981, 107(6): 1227-1238.

[12]

Liu

P L

, Der Kiureghian

. Multivariate distribution models with prescribed marginals and covariances[J]. Probabilistic Engineering Mechanics, 1986, 1(2): 105-112.

[13]

Z H

, Cai

C H

, Zhao

Y G

. Structural reliability analysis including correlated random variables based on third-moment transformation[J]. Journal of Structural Engineering, 2017, 143(8): 04017067.

[14]

Z H

, Cai

C H

, Zhao

Y G

, et al. Normalization of correlated random variables in structural reliability analysis using fourth-moment transformation[J]. Structural Safety, 2020, 82: 101888.

[15]

Nagler

. Kdecopula: An R package for the kernel estimation of bivariate copula densities[J]. Journal of Statistical Software, 2018, 84(7): 1-22.

[16]

赵岩, 张亚辉, 林家浩. 基于虚拟激励法的大跨度桥梁非平稳随机振动分析闭合解[J]. 计算力学学报, 2024, 41(1): 202-208.

(ZHAO

Yan

, ZHANG

Yahui

, LIN

Jiahao

. Closed solution for non-stationary random vibration analysis of long-span bridge based on pseudo excitation method[J]. Chinese Journal of Computational Mechanics, 2024, 41(1): 202-208. (in Chinese))

2025年第42卷第5期

PDF下载

引用本文

BibTeX

文章信息

doi: 10.7511/jslx20240202002

接收时间：2024-02-02
首发时间：2026-03-24
出版时间：2025-10-28

补充材料

相关文章

文章信息

作者

出版历史

收稿日期：2024-02-02
修回日期：2024-04-21

基金

国家自然科学基金(52278135)

作者信息

北京工业大学　城市建设学部，北京　100124

参考文献

分享链接

https://castjournals.cast.org.cn/joweb/jslxxb/CN/10.7511/jslx20240202002

分享至

全文二维码

扫描看全文

引用本文

BibTeX

本文的引用情况

2种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT