利用多项式混沌展开的结构可靠性分析
张明1, 王恩志1, 刘耀儒1, 齐文彪2, 王德辉3    
1. 清华大学 水沙科学与水利水电工程国家重点实验室, 北京 100084;
2. 吉林省水利水电勘测设计研究院, 长春 130021;
3. 福州大学 土木工程学院, 福州 350116
摘要:实际工程设计优化、设计空间搜索、灵敏度分析、可靠性分析等问题,若单次模拟比较费时,直接用原模型进行数千、甚至数百万次模拟是不可能完成的任务。多项式混沌展开方法是解决这类问题的有效方法,其方法表达和程序实现是应用中关注的重点问题。该文介绍了多项式混沌展开方法的数学理论,并将之用于结构可靠性分析。首先,将结构可靠性分析的功能响应函数以多项式混沌展开表示,其中统一采用Hermite多项式。给出Hermite多项式的一种适合计算机程序生成的通项形式,实现多项式混沌展开的计算程序的通用化,以及多项式次数的自适应选择。其次,利用具有显式功能函数的结构可靠性分析算例,考察所构建的代理模型的正确性和适用性。结果表明,多项式混沌展开的次数越高,模型的精度就越高,具有良好的收敛性,同时表明仅就考察代理模型而言,利用显式功能函数是最简便的方式。
关键词优化设计    结构可靠性    多项式混沌展开    Hermite多项式    功能函数    
Reliability analysis of structures using polynomial chaos expansions
ZHANG Ming1, WANG Enzhi1, LIU Yaoru1, QI Wenbiao2, WANG Dehui3    
1. State Key Laboratory of Hydroscience and Engineering, Tsinghua University, Beijing 100084, China;
2. Jilin Province Water Resource and Hydropower Consultative Company, Changchun 130021, China;
3. College of Civil Engineering, Fuzhou University, Fuzhou 350116, China
Abstract: Practical engineering issues such as design optimization, design space exploration, sensitivity analyses, and reliability analyses need many simulations. If a single simulation is very time-consuming, engineers cannot perform the thousands or even millions of simulations needed for such analyses. The polynomial chaos expansion (PCE) method is an effective method that allows analyses of complex problems. This paper introduces the mathematical theory of the PCE method and presents a structural reliability analysis example. The performance response function for the structural reliability analysis is expressed as a PCE using Hermite polynomials. A general form of the Hermite polynomial, which is suitable for use in a computer program, is used to generalize the PCE analysis program and the adaptive selection of the polynomial order. Then, the accuracy and applicability of the surrogate model are verified using structural reliability analysis examples with explicit performance functions. The results show that the model has an excellent convergence rate with higher order PCE giving higher accuracy. The examples also show that the direct use of explicit performance functions is the easiest way to investigate PCE surrogate models.
Key words: optimal design    structural reliability    polynomial chaos expansion (PCE)    Hermite polynomial    performance function    

在实际工程优化设计、灵敏性分析等问题中,经常遇到这样的场合:一方面需要进行复杂的高精度仿真分析,如结构的有限元分析、计算流体动力学分析;另一方面这种费时的计算分析还要重复进行多次,如可靠性分析、易损性分析[1-3]。这就给大型复杂的仿真分析带来了挑战,利用所建立的仿真模拟模型即使可以完成数次分析过程,但是大量地进行重复却因耗时过长成为不可能。例如用最简单的直接Monte Carlo方法来计算失效概率,如果用原模型模拟一次的机时为半个小时,符合工程实际需要的失效概率最低为10-4量级,则所需随机模拟次数至少10万次,整个计算需要2 000多天!

为了降低计算量,可以对计算方法进行改进,如利用方差缩减技术的各种改进Monte Carlo方法(如重要性抽样)[4],但是其效果也都十分有限。为应对上述挑战,从根本上克服始终利用原模型的局限,产生了基于代理模型的优化设计,用代理模型取代这些费时的仿真模型,在代理模型的基础上进行优化[5-14]。其基本思想是对原性能函数构建代理模型,代理模型不但能产生已有的性能响应,还能预测其他参数下的性能响应,代理模型足够简单,可极大降低计算量,保证精度的同时提高效率。

代理模型是原始真实模型的虚拟和近似,其复杂性、逼近程度跟原模型的复杂程度没有必然的联系,无需对原模型做任何改动,对原模型是无扰动、非侵入的[1, 13]。采用代理模型方法时,数值模拟分析与概率分析分开进行,不需要像传统的随机有限元方法[15]那样修改有限元源代码,数值分析程序没有改变,这种便利性有利于代理模型的应用和与商业数值分析软件的结合。

这样的代理模型已有不少,如利用响应面方法[4]、kriging方法[3]、人工神经网络(artificial neural network, ANN)方法[4]等建立模型,但这些方法在响应面试验设计、kriging方法的相关函数、ANN的深度等方面均比较灵活,没有统一的做法。近年焕发生机的多项式混沌(展开)方法[9-14],是当系统参数具有随机不确定性时,用来确定动力系统中不确定性演化的有效方法。通过对随机变量输入到响应函数输出的映射过程建立代理模型,利用正交多项式的优异性能,使得多项式混沌代理模型方法极大地降低了其他诸多不确定性分析方法的意义,成为代理模型的代表,也成为不确定性研究的发展趋势。

目前多项式混沌代理模型的研究仍然很少,在结构分析领域也几乎是空白。李典庆等[13-14]利用简单的多项式混沌展开建立非侵入式随机分析方法,随机分析无需修改有限元代码,能够有效与商业有限元软件结合,在岩土工程领域得到了应用。但就其多项式混沌展开而言,多项式都是手动输入计算机程序中,无法或难以对任意次多项式进行程序化,且多项式混沌与有限元方法混在一起也不够清晰。

多项式混沌展开的理论基础较为深奥,加之多元正交多项式本身的形式很复杂,这既限制了其推广应用,也在计算机程序实现上造成了困难。因此,给出便于程序实现的多项式的通项形式,实现任意次多项式混沌的程序化,以及多项式混沌次数的自动选择,最大程度地降低理解上的难度,减少出错几率,并利用显式性能函数以最简单的方式检验模型的正确性,这些都是在多项式混沌展开理论和应用中需要解决的问题,也是一些基础性的问题。此外,以工程中的实际算例,考察多项式混沌代理模型的适用性也是必要的。

1 多项式混沌展开

多项式混沌(polynomial chaos, PC),或称Wiener混沌(Wiener chaos),是指采用多项式基构成随机空间,来描述和传播(合成)随机变量的不确定性。“混沌”只是简单地指不确定性,与奇怪吸引子等无关,即一个动力系统因输入参数的不确定性,其输出所具有的不确定性[16]。“多项式”是指用多项式展开来传播不确定性,常采用正交多项式进行展开和作为空间的正交基。有时更明确地将多项式混沌称作多项式混沌展开(polynomial chaos expansion, PCE)。

设系统的输入为Ξ=[Ξ1; Ξ2; …; Ξn]为独立随机变量(有时又称为芽(germ)),输出为gΞ(Ξ),用多项式混沌表示为[16]

$ g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi})=\sum\limits_{i=0}^{\infty} c_{i} P_{i}(\boldsymbol{\varXi}). $ (1)

其中:Pi(Ξ)为i次正交多项式,常称为i次同质混沌(homogeneous chaos);ci为待定系数。多项式混沌展开在均方意义下是收敛的[7-8]

正交多项式Pi(Ξ) (i=0, 1, …)构成Hilbert函数空间的完备正交基,其内积为

$ \begin{aligned} \left\langle P_{i}(\boldsymbol{\varXi}), P_{j}(\boldsymbol{\varXi})\right\rangle &=\int_{\mathbb{R}^{n}} P_{i}(\xi) P_{j}(\xi) W(\xi) \mathrm{d} \xi=\\ &\left\langle P_{i}^{2}(\boldsymbol{\varXi})\right\rangle \delta_{i j} \end{aligned} $ (2)

其中,W(ξ)为权函数,δij为Kronecker delta。

对于随机空间,W(ξ)为独立随机变量Ξ的概率密度函数$ f_{\boldsymbol{\varXi}}(\xi)=f_{\varXi_{1}}\left(\xi_{1}\right) f_{\varXi_{2}}\left(\xi_{2}\right) \cdots f_{\varXi_{n}}\left(\xi_{n}\right)$,此时式(2)可用数学期望函数E表示为

$ E\left[P_{i}(\boldsymbol{\varXi}) P_{j}(\boldsymbol{\varXi})\right]=E\left[P_{i}^{2}(\boldsymbol{\varXi})\right] \delta_{i j}. $ (3)

可注意到P0(Ξ)=1。由正交性式(3)可知,对于i≥1,Pi(Ξ)的均值E[(Pi(Ξ)]=0,于是E[Pi2(Ξ)] 即Pi(Ξ)的方差。因此,内积是半正定的,它给出了具有有限方差的随机变量空间的一个度量。

由式(3)知,要保证正交性,正交多项式须与不同分布的随机变量Ξ相对应。例如,与均匀分布变量对应的是Legendre多项式,Jacobi多项式相应于beta分布变量,Laguerre多项式相应于gamma分布变量[1, 7-8]

如果Ξ=[Ξ1; Ξ2; …; Ξn]为独立标准正态随机变量,Pi(Ξ)为Hermite多项式,此时式(1)可写成[1, 7-8]

$ \begin{gathered} g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi})=a_{0} \varGamma_{0}+\sum\limits_{i_{1}=1}^{n} a_{i_{1}} \varGamma_{1}\left(\boldsymbol{\varXi}_{i_{1}}\right)+ \\ \sum\limits_{i_{1}=1}^{n} \sum\limits_{i_{2}=1}^{i_{1}} a_{i_{1} i_{2}} \varGamma_{2}\left(\boldsymbol{\varXi}_{i_{1}}, \boldsymbol{\varXi}_{i_{2}}\right)+ \\ \sum\limits_{i_{1}=1}^{n} \sum\limits_{i_{2}=1}^{i_{1}} \sum\limits_{i_{3}=1}^{i_{2}} a_{i_{1} i_{2} i_{3}} \varGamma_{3}\left(\boldsymbol{\varXi}_{i_{1}}, \boldsymbol{\varXi}_{i_{2}}, \boldsymbol{\varXi}_{i_{3}}\right)+\cdots . \end{gathered} $ (4)

其中,Γq(Ξi1, Ξi2, …, Ξiq)为自由度为qq次Hermite多项式,a0, ai1, ai1i2, ai1i2i3, …为待定系数。

多项式混沌通常采用与独立标准正态随机变量相对应的Hermite多项式。如果多项式混沌采用与非正态随机变量相应的正交多项式,则称为广义多项式混沌(generalized polynomial chaos, GPC)[17]

多项式混沌展开中各正交多项式基函数是相互独立的标准随机变量的函数,输入变量的非正态性和相关性会影响代理模型的构建。若应用广义多项式混沌于实际的不确定分析,理论上不需要随机变量的正态化处理,但是广义多项式混沌由于使用了不同的正交多项式,加大了问题解决办法一般化的难度。因此,在本文的多项式混沌展开中只使用Hermite多项式,而这只需要在之前解决非正态随机变量的相关性并作等效正态化处理即可。

式(4)可以用系数bi和正交多项式Ψi(Ξ) (i=0, 1, …)重写成以下简洁的紧凑形式:

$ g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi})=\sum\limits_{i=0}^{\infty} b_{i} \varPsi_{i}(\boldsymbol{\varXi}). $ (5)

实际应用时,可采用有限项的多项式混沌展开表达式。取至p次多项式,对于式(1)则∞变成p,而式(4)和式(5)分别为

$ \begin{gathered} g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi}) \approx a_{0} \varGamma_{0}+\sum_{i_{1}=1}^{n} a_{i_{1}} \varGamma_{1}\left(\boldsymbol{\varXi}_{i_{1}}\right)+ \\ \sum_{i_{1}=1}^{n} \sum_{i_{2}=1}^{i_{1}} a_{i_{1} i_{2}} \varGamma_{2}\left(\boldsymbol{\varXi}_{i_{1}}, \boldsymbol{\varXi}_{i_{2}}\right)+\cdots+ \\ \sum_{i_{1}=1}^{n} \sum_{i_{2}=1}^{i_{1}} \cdots \sum_{i_{p}=1}^{i_{p-1}} a_{i_{1} i_{2} \cdots i_{p}} \varGamma_{p}\left(\boldsymbol{\varXi}_{i_{1}}, \boldsymbol{\varXi}_{i_{2}}, \cdots, \boldsymbol{\varXi}_{i_{p}}\right), \end{gathered} $ (6)
$ g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi}) \approx \sum\limits_{i=0}^{m-1} b_{i} \varPsi_{i}(\boldsymbol{\varXi})=\sum\limits_{\alpha \in \mathbb{N}^{n}} a_{\boldsymbol{\alpha}} \varGamma_{\boldsymbol{\alpha}}(\boldsymbol{\varXi}). $ (7)

其中,bjΨj分别对应于式(4)中的ai1i2ijΓj(Ξi1, Ξi2, …, Ξij);α为多元指标。

对于多元指标α,可通过分级字典序列(graded lexicographic order)来确定[8]。在每一级分级字典序列中,多元指标α按字典序列进行排列。在分级字典序列中,第0级对应|α|=0,只有1个多元指标,在多项式混沌展开中表现为常数项;而第$i(i \in \mathbb{N}^{n} $i≥1)级对应|α|=i,所包含多元指标的个数mi

$ m_{i}=\frac{1}{i !} \prod\limits_{j=0}^{i-1}(n+j) . $ (8)

对于n元随机变量、截断至p次的多项式混沌展开,式(6)和式(7)中多项式的个数均为

$ m=1+\sum\limits_{i=0}^{p} m_{i}=\left(\begin{array}{c} n+p \\ p \end{array}\right)=\frac{(n+p) !}{n ! p !}. $ (9)

n元变量x=[x1; x2; …; xn]的q次Hermite多项式由下式确定[1, 18]

$ \begin{gathered} \varGamma_{q}\left(x_{i_{1}}, x_{i_{2}}, \cdots, x_{i_{q}}\right)=(-1)^{q} \exp \left(\frac{1}{2} \boldsymbol{x}^{\mathrm{T}} \boldsymbol{x}\right) \cdot \\ \frac{\partial^{q}}{\partial x_{i_{1}} \partial x_{i_{2}} \cdots \partial x_{i_{q}}} \exp \left(-\frac{1}{2} \boldsymbol{x}^{\mathrm{T}} \boldsymbol{x}\right) . \end{gathered} $ (10)

特别地,一元概率学者的Hermite多项式(probabilists' Hermite polynomials) Hei(x)与其次数i的奇偶性一致,据式(10)还可写成以下递归公式:

$ \left\{\begin{array}{l} H e_{0}(x)=1 , \\ H e_{1}(x)=x , \\ H e_{i}(x)=x H e_{i-1}(x)-(q-1) H e_{i-2}(x), \quad i>1 . \end{array}\right. $ (11)

多元q次Hermite多项式Γq(xi1, xi2, …, xiq)为某些单变量Hermite多项式Hei(xj) (i=1, 2, …; j=1, 2, …, n)的乘积。本文将其写成以下形式,以便在计算机程序上实现:

$ {\varGamma _q}(\overbrace {\underbrace {{x_{{i_1}}}, \cdots , {x_{{i_1}}}}_{{n_1}{\rm{个 }}}, \underbrace {{x_{{i_2}}}, \cdots , {x_{{i_2}}}}_{{n_2}个}, \cdots , \underbrace {{x_{{i_m}}}, \cdots , {x_{{i_m}}}}_{{n_m}个}}^{q{\rm{ 个}}}) = \\ H e_{n_{1}}\left(x_{i_{1}}\right) H e_{n_{2}}\left(x_{i_{2}}\right) \cdots H e_{n_{m}}\left(x_{i_{m}}\right). $ (12)

结合式(11),多元q次Hermite多项式能便捷地利用式(12)得到,也易知其关于某个变量的奇偶性。

2 结构可靠性分析的多项式混沌代理模型 2.1 多项式混沌展开代理模型

结构可靠性分析依赖于结构的性能响应,需要已知结构的功能函数。式(6)或(7)提供了一种可能性,可以设计合适的p次多项式混沌展开gΞ(Ξ)来建立结构功能函数Z=gX(X)的一个代理(surrogate),来代替原功能函数gX(X)进行可靠度分析。由此建立的代理模型(surrogate model)又称元模型(metamodel)、仿真器(emulator)。

X=T(Ξ),其中T为转换函数,此时有

$ Z=g_{\boldsymbol{X}}(\boldsymbol{X})=g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi})=\sum\limits_{i=0}^{m} b_{i} \varPsi_{i}(\boldsymbol{\varXi}). $ (13)

其中,$g_{\boldsymbol{\varXi}}=g_{\boldsymbol{X}} \circ \boldsymbol{T} $

在建立某功能函数的多项式混沌代理时,向量XΞ的维数不必相同。例如,若功能函数Z~χ2(k),则Z可以表示成k个独立标准正态变量Ξi的平方和,即$Z=\sum\limits_{i=1}^{k} \boldsymbol{\varXi}_{i}^{2} $,这也是其严格的多项式混沌形式,其中k不一定等于X的维数n。一般地,一个有意义的有限项多项式混沌代理模型的Ξ的维数应等于或大于X的维数,#Ξ≥#X

biΨi(Ξ)或ai1i2iqΓq(Ξi1, Ξi2, …, Ξiq)称作模态(mode),其系数biai1i2iq就是模态的强度。如果结构的响应信息即功能函数Z是已知的,则式(13)中的待定系数[1]

$b_{i}=\frac{E\left[Z \varPsi_{i}(\boldsymbol{\varXi})\right]}{E\left[\varPsi_{i}^{2}(\boldsymbol{\varXi})\right]}. $ (14)

可以注意到b0=μZ,其中μZ为功能函数Z的均值。

在利用式(14)时,还可注意到

$E\left(\boldsymbol{\varXi}_{i}^{k}\right)= \begin{cases}0, & k \text { 为奇数; } \\ \frac{k !}{2^{k / 2}(k / 2) !}=(k-1) ! !, & k \text { 为偶数. }\end{cases} $ (15)

以上确定待定系数的方法称为Galerkin投影方法。确定多项式混沌中的待定系数,也常采用随机配点的方法。选取独立标准正态空间中的配点ξi(i=1, 2, …, N),将其映射为原始空间中的配点xi=T(ξi),计算功能函数的值gX(xi),由式(7)和式(13),可得确定多项式混沌中的系数的线性方程组如下:

$ \boldsymbol{A} b=\boldsymbol{g}_{X}. $ (16)

其中:$ g_{\boldsymbol{X}}=\left[g_{\boldsymbol{X}}\left(\boldsymbol{x}_{1}\right) ; g_{\boldsymbol{X}}\left(\boldsymbol{x}_{2}\right) ; \cdots ; g_{\boldsymbol{X}}\left(\boldsymbol{x}_{N}\right)\right], \boldsymbol{b}=$ $\left[b_{0} ; b_{1} ; \cdots ; b_{m}\right]$ 为待定系数, A为相应的系数矩阵。即

$\boldsymbol{A}=\left[\begin{array}{ccccc} 1 & \varPsi_{1}\left(\xi_{1}\right) & \varPsi_{2}\left(\xi_{1}\right) & \cdots & \varPsi_{m}\left(\xi_{1}\right) \\ 1 & \varPsi_{1}\left(\xi_{2}\right) & \varPsi_{2}\left(\xi_{2}\right) & \cdots & \varPsi_{m}\left(\xi_{2}\right) \\ \vdots & \vdots & \vdots & & \vdots \\ 1 & \varPsi_{1}\left(\xi_{\mathrm{N}}\right) & \varPsi_{2}\left(\xi_{\mathrm{N}}\right) & \cdots & \varPsi_{m}\left(\xi_{\mathrm{N}}\right) \end{array}\right]. $ (17)

配点数N须多于待定系数的数目m,式(16)可用基于最小二乘回归分析的方法求解。

可以注意到,随机配点方法中并不要求显式功能函数一定存在,所需要的只是在配点处功能函数的值。通过有限次试验或试算,得到结构的多组输入与响应,只要能够确定多项式混沌展开中的系数,即使结构输入与响应的映射是隐式的,分析也可进行下去。

基于得到的代理函数gΞ(Ξ),当利用一次可靠度方法时,会用到梯度$ \nabla g_{\boldsymbol{\varXi}}(\boldsymbol{\varXi})$,为此需要计算Hermite多项式$ \varGamma_{q}\left(x_{i_{1}}, x_{i_{2}}, \cdots, x_{i_{q}}\right)$的导数。考虑到式(12),可以借助以下关系式来完成:

$ \mathrm{He}_{i}^{\prime}(x)=i H e_{i-1}(x). $ (18)

至此,已设法将功能函数表示为多项式混沌展开,建立了代理模型。在构建完成的代理模型上可以得到原始模型的所有概率特征,不确定性分析可以很方便地进行,所有的分析步骤都是程式化的。

2.2 多项式混沌展开的程序实现

多项式混沌展开方法首先要构建多项式混沌展开模型,即选择与输入随机变量相对应的正交多项式基函数,将输出随机变量表示为正交多项式基函数的线性组合,然后确定多项式混沌展开系数,即通过计算得到各正交多项式基函数的对应权重。对该方法计算机程序化时的几个问题做出说明和评述如下。

1) 概率学者的Hermite多项式计算。MATLAB的内置函数hermite H,可计算任意次一元物理学者的Hermite多项式(physicists' Hermite polynomials) Hq(x),由此可转换成概率学者的Hermite多项式Heq(x),但当处理大维数的向量时计算速度很慢。本文根据需要,利用式(11)编制计算任意次数的一元概率学者的Hermite多项式的函数hermite He。在本文方法中,该步也是计算多元多项式Γq的基础。

2) 多元多项式混沌的多重指标序列的确定。对于有n个随机变量、截断次数为p次的多项式混沌,利用分级字典序列来生成全部m个多元序列(参见式(8)和(9)),按照多元指标序列的顺序,利用式(12)依次生成(计算)所对应的Hermite多项式。此步也是多项式混沌展开方法实施中最为关键的一步。靠手动将多项式逐个输入进计算机程序的办法是不可取的,这极大地增加了出错的几率,当随机变量个数多、截断次数高时,也是不可行的。

3) 多项式混沌展开次数p的自动选取。这相对容易。可以先选定某一p值,根据前后两次计算得到结果(如可靠指标β)是否满足收敛准则来确定,如果前后两次计算结果的相对误差不满足要求,则可使p增加1次,继续重复相同的计算过程,直到该误差满足要求为止。代理模型是关于模型的模型,原模型的差异性决定了代理模型恐难有统一的收敛性标准,这里确定多项式混沌次数的做法是一种简便、有效的实用方法。

3 算例分析 3.1 简支梁问题[1]

图 1所示简支梁受集中力P1P2作用。P1P2均服从正态分布,均值分别为μP1=10 kN和μP2=15 kN,标准差分别为σP1=2 kN和σP2=3 kN。梁的允许最大弯矩Ma服从均值为μMa=50 kN·m、标准差为σMa=5 kN·m的正态分布。求P1P2相互独立及相关系数为ρP1P2=0.25两种情况下梁的失效概率。

图 1 集中力作用的简支梁

梁的最大弯矩为Mm=(P1+2P2)l/9,其中梁长l=9 m。因此梁的极限状态函数为

$\begin{gathered} Z=g\left(M_{\mathrm{a}}, P_{1}, P_{2}\right)=M_{\mathrm{a}}-M_{\mathrm{m}}= \\ M_{\mathrm{a}}-\frac{1}{9}\left(P_{1}+2 P_{2}\right) l \end{gathered} $

根据结构可靠性理论[4],梁的可靠指标β存在解析解,即

$ \beta=\frac{\mu_{M_{\mathrm{a}}}-\left(\mu_{P_{1}}+2 \mu_{P_{2}}\right) l / 9}{\sqrt{\sigma_{M_{\mathrm{a}}}^{2}+\left(\sigma_{P_{1}}^{2}+4 \rho_{P_{1} P_{2}} \sigma_{P_{1}} \sigma_{P_{2}}+4 \sigma_{P_{2}}^{2}\right) l^{2} / 81}} . $

考虑到输入随机变量均服从正态分布,且原功能函数次数最高为1次,故选择一次多项式混沌展开构建代理模型。令相关正态随机变量X=[Ma; P1; P2],通过简单的线性变换[4]可将其变成独立标准正态变量Ξ=[Ξ1; Ξ2; Ξ3],因此原功能函数准确的多项式混沌代理函数为gΞ(Ξ)=b0+b1Ξ1+b2Ξ2+b3Ξ3。据此,可靠指标$ \beta=b_{0} / \sqrt{b_{1}^{2}+b_{2}^{2}+b_{3}^{2}}$

在计算待定系数bi时,分别使用了Galerkin投影和随机配点的方法。相应地,可靠指标β的计算结果,当ρP1P2=0时,分别为1.240 4和1.240 3,解析解为1.240 3;当ρP1P2=0.25时,分别为1.186 0和1.186 8,解析解为1.186 8。

3.2 重力坝问题[19-20]

砌石重力坝如图 2所示。砌石容重γ服从正态分布,均值21.8 kN/m3,变异系数0.03;混凝土容许抗压强度fc服从对数正态分布,均值2 570 kN/m2,变异系数0.278 6;作用水头H服从正态分布,均值58.1 m,变异系数0.017;摩擦系数f服从正态分布,均值0.72,变异系数0.131。

图 2 重力坝问题示意图(单位:m)

考虑坝基平面抗滑稳定的极限状态,得到相应的功能函数为$ Z=f \sum W-P$,其中$ \sum W$为坝自重和铅直方向水压力之和,P为坝水平方向水压力之和。代入相关数据得到

$Z_{1}=1643.5 \gamma f+58.7 H f-1147 f-5 H^{2} . $

功能函数Z1对应的失效概率pf1可通过以下在失效区域上的积分得到,即

$ {{p}_{\text{f}1}}=\iiint\limits_{\gamma \le \left( 1147+5{{H}^{2}}/f-58.7H \right)/1643.5}{{}}{{f}_{\gamma }}(\gamma ){{f}_{H}}(h){{f}_{f}}(f)\text{d}\gamma \text{d}h~\text{d}f= $
$ \begin{aligned} &\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} F_{\gamma}\left[\left(1147+5 H^{2} / f-58.7 H\right) / 1643.5\right] \cdot\\ &f_{H}(h) f_{f}(f) \mathrm{d} h \mathrm{~d} f . \end{aligned} $

其中,fγ(γ)、fH(h)和ff(f)分别表示γHf的概率密度函数,Fγ(γ)表示γ的累积分布函数。

完成此积分,得坝基平面抗滑稳定的失效概率的精确解为pf1=2.043 1×10-3

将独立正态随机变量[γ; H; f]标准化为独立标准正态随机变量[Ξ1; Ξ2; Ξ3],多项式混沌展开次数p分别取1, 2, …, 6,基于相应的p次多项式混沌展开代理模型,利用Monte Carlo方法计算坝的失效概率,得到失效概率的估计值$\hat{p}_{\mathrm{f} 1} $依次为2.215 7×10-3, 2.042 8×10-3, 2.044 1×10-3, 2.042 9×10-3, 2.046 1×10-3, 2.042 1×10-3。可见,二次多项式混沌展开即可得到满意的计算结果,多项式混沌代理模型收敛速度很快。

考虑坝踵容许抗压强度的极限状态,得到相应的功能函数$ Z=f_{\mathrm{c}}-\left(\frac{\sum W}{T}+\frac{6 \sum M}{T^{2}}\right)$,其中$\sum M $为各力对于坝底面形心的力矩之和,T为坝底的宽度。代入相关参数得到

$ Z_{2}=f_{\mathrm{c}}-8.474 \gamma+0.1856 \mathrm{H}-0.0036 \mathrm{H}^{3}-34.04 . $

功能函数Z2对应的失效概率pf2可通过以下在失效区域上的积分得到,即

${{p}_{\text{f}2}}=\iiint\limits_{{{f}_{\text{c}}}\le 8.474\gamma -0.1856H+0.036{{H}^{3}}+34.04}{{}}\left( {{f}_{\text{c}}} \right)\cdot $
$ \begin{gathered} f_{\gamma}(\gamma) f_{H}(h) \mathrm{d} f_{\mathrm{c}} \mathrm{d} \gamma \mathrm{d} h= \\ \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} F_{f_{\mathrm{c}}}\left(8.47 \gamma-0.1856 H+0.036 H^{3}+34.04\right) \\ f_{\gamma}(\gamma) f_{H}(h) \mathrm{d} \gamma \mathrm{d} h \end{gathered} $

其中,ffc(fc)和Ffc(fc)分别表示fc的概率密度函数和累积分布函数。

完成此积分,得坝踵容许抗压强度极限状态的失效概率的精确解为pf2=3.063 4×10-5

将正态随机变量γH分别化为标准正态随机变量Ξ1Ξ2,将对数正态变量fc正态标准化[4]成为Ξ3,基于多项式混沌展开代理功能函数,自适应地选择多项式混沌展开次数,利用Monte Carlo方法计算坝的失效概率,当p=6时可满足精度要求,得到失效概率的估计值为$ \hat{p}_{\mathrm{f} 2}$ =3.069 0×10-5。加大p无疑使构建的代理模型更准确,增加失效概率估计值收敛的稳定性,但必然带来更大的运算量,使模拟时间大幅增加,这时就面临着精度和效率平衡。

4 结论

本文在介绍和分析多项式混沌展开的基础上,给出多元Hermite多项式通项的正确形式,以及多项式混沌用于不确定性分析的一般格式,考察在结构可靠性分析中的正确性与适用性。主要结论如下:

1) Hermite多项式表示的多项式混沌方法,适用于所有结构可靠性分析问题,此时要求式中的变量Ξ为独立正态随机变量。在应用时须先将实际的非正态相关随机变量X进行相应的独立标准正态化处理。本文的做法便于程序的规范化。

2) 多项式混沌展开计算程序的通用化、一般化,有赖于多元Hermite多项式通式的正确表达。有了式(12)即可方便地利用一元Hermite多项式得到多元Hermite多项式,使得选择任意次数的多元多项式混沌成为可能。

3) 多项式混沌展开作为一种收敛快且成本低的不确定性量化方法,在保证计算精度的前提下,应尽量利用自适应次数选择,采用最低的多项式混沌次数。

4) 多项式混沌展开模型适用于任何存在输入—响应映射的情形,但是考察其正确性,利用简单的存在显式功能函数的实例进行考察是最适宜的。

参考文献
[1]
CHOI S K, GRANDHI R V, CANFIELD R A. Reliability-based structural design[M]. London: Springer, 2007.
[2]
侯少康, 刘耀儒. 双护盾TBM掘进数值仿真及护盾卡机控制因素影响分析[J]. 清华大学学报(自然科学版), 2021, 61(8): 809-817.
HOU S K, LIU Y R. Numerical simulations of double-shield TBM tunneling for analyzing shield jamming control factors[J]. Journal of Tsinghua University (Science and Technology), 2021, 61(8): 809-817. (in Chinese)
[3]
熊芬芬, 杨树兴, 刘宇, 等. 工程概率不确定性分析方法[M]. 北京: 科学出版社, 2015.
XIONG F F, YANG S X, LIU Y, et al. Probabilistic uncertainty analysis methods for engineering problems[M]. Beijing: Science Press, 2015. (in Chinese)
[4]
张明, 金峰. 结构可靠度计算[M]. 北京: 科学出版社, 2015.
ZHANG M, JIN F. Structural reliability computations[M]. Beijing: Science Press, 2015. (in Chinese)
[5]
赵威. 结构可靠度分析代理模型方法研究[D]. 哈尔滨: 哈尔滨工业大学, 2012.
ZHAO W. Study on surrogate model methods for structural reliability analysis[D]. Harbin: Harbin Institute of Technology, 2012. (in Chinese)
[6]
BAROTH J, SCHOEFS F, BREYSSE D. Construction reliability: Safety, variability and sustainability[M]. Hoboken: John Wiley & Sons, 2013.
[7]
PETTERSSON M P, IACCARINO G, NORDSTRÖM J. Polynomial chaos methods for hyperbolic partial differential equations[M]. Cham: Springer, 2015.
[8]
PELLEGRINI G. Polynomial chaos expansion with applications to PDEs[D]. Verona: University of Verona, 2013.
[9]
XU J, WANG D. Structural reliability analysis based on polynomial chaos, Voronoi cells and dimension reduction technique[J]. Reliability Engineering & System Safety, 2019, 185: 329-340.
[10]
YU H, GILLOT F, ICHCHOU M. A polynomial chaos expansion based reliability method for linear random structures[J]. Advances in Structural Engineering, 2012, 15(12): 2097-2111. DOI:10.1260/1369-4332.15.12.2097
[11]
HAWCHAR L, EL SOUEIDY C P, SCHOEFS F. Time-variant reliability analysis using polynomial chaos expansion[C]//12th International Conference on Applications of Statistics and Probability in Civil Engineering. Vancouver, Canada: ICASP, 2015.
[12]
CLARICH A, MARCHI M, RIGONI E, et al. Reliability-based design optimization applying polynomial chaos expansion: Theory and applications[C]//10th World Congress on Structural and Multidisciplinary Optimization. Orlando, USA, 2013.
[13]
李典庆, 蒋水华. 边坡可靠度非侵入式随机分析方法[M]. 北京: 科学出版社, 2016.
LI D Q, JIANG S H. Non-intrusive stochastic analysis method of slope reliability[M]. Beijing: Science Press, 2016. (in Chinese)
[14]
蒋水华, 冯晓波, 李典庆, 等. 边坡可靠度分析的非侵入式随机有限元法[J]. 岩土力学, 2013, 34(8): 2347-2354.
JIANG S H, FENG X B, LI D Q, et al. Reliability analysis of slope using non-intrusive stochastic finite element method[J]. Rock and Soil Mechanics, 2013, 34(8): 2347-2354. (in Chinese)
[15]
HALDAR A, MAHADEVAN S. Reliability assessment using stochastic finite element analysis[M]. New York: John Wiley & Sons, 2000.
[16]
O'HAGAN A. Polynomial chaos: A tutorial and critique from a statistician's perspective[EB/OL]. (2013-03-24). http://www.tonyohagan.co.uk/academic/pdf/Polynomial- chaos. pdf.
[17]
XIU D B. Generalized (Wiener-Askey) polynomial chaos[D]. Providence: Brown University, 2004.
[18]
WüNSCHE A. Hermite and Laguerre 2D polynomials[J]. Journal of Computational and Applied Mathematics, 2001, 133(1-2): 665-678. DOI:10.1016/S0377-0427(00)00681-6
[19]
蒋水华, 张曼, 李典庆. 基于Hermite正交多项式逼近法的重力坝可靠度分析[J]. 武汉大学学报(工学版), 2011, 44(2): 170-174.
JIANG S H, ZHANG M, LI D Q. Reliability analysis of gravity dam using Hermite orthogonal polynomials approximation method[J]. Engineering Journal of Wuhan University, 2011, 44(2): 170-174. (in Chinese)
[20]
吕泰仁, 吴世伟. 用几何法求构件的可靠指标[J]. 河海大学学报(自然科学版), 1988, 16(5): 86-93.
LYU T R, WU S W. Geometric method for solving reliability index of structures[J]. Journal of Hohai University (Natural Sciences), 1988, 16(5): 86-93. (in Chinese)