基于改進統(tǒng)計矩點估計法和最大熵原理的結構整體可靠度分析

王 濤，李正良,2，范文亮,2

(1. 重慶大學土木工程學院，重慶 400045；2. 山地城鎮(zhèn)建設與新技術教育部重點實驗室(重慶大學)，重慶 400045)

由于材料、荷載不可避免的隨機性，結構可靠度分析成為結構設計的基礎和重要保障。經(jīng)過長期的發(fā)展，結構構件可靠度方法已形成了數(shù)值近似法(如一次二階矩法[1])、隨機模擬法(如Monte Carlo 法[2])、代理模型法(如響應面法[3]) 和數(shù)值積分法(如矩方法[4]) 四大類方法，其基本理論已經(jīng)相對成熟且能夠應用于工程實際[1]。然而，結構整體可靠度分析仍在發(fā)展階段。大體上，結構整體可靠度方法可分為三類，即Monte Carlo 法[5]、失效模式識別法[5]以及等價描述方法[6-8]。Monte Carlo 法由于其計算昂貴難以應用于大型復雜結構系統(tǒng)，常作為一種校準法；而經(jīng)典的失效模式識別法則存在組合爆炸和相關失效問題[6]；對于復雜工程結構進行整體可靠度分析，Monte Carlo 法計算成本難以接受，失效模式識別法應用困難，相比而言，由于等價描述方法將結構整體可靠度用一個等價的功能函數(shù)描述，有效地避免了組合爆炸與相關失效問題，形成了一種較為實用可行的思路。

功能函數(shù)經(jīng)過等價描述后，其函數(shù)形式變得較為復雜，由于基于驗算點的方法處理函數(shù)形式復雜的功能函數(shù)會存在多驗算點問題，采用此類方法進行整體可靠度分析可能致使其計算精度不理想[9]。將等價描述后的功能函數(shù)與不依賴驗算點的方法結合進行結構整體可靠度分析是一種行之可效的思路。矩方法是獲得結構可靠度的一種有效途徑，既不依賴于驗算點，也無需復雜的積分計算和求導運算，因而得到廣泛應用[10-15]。其主要思路是先獲取功能函數(shù)的統(tǒng)計矩，再而基于統(tǒng)計矩信息求解功能函數(shù)的概率密度函數(shù)，最后通過對功能函數(shù)的概率密度函數(shù)積分獲取失效概率或可靠指標。

功能函數(shù)的統(tǒng)計矩估計通常通過數(shù)值積分方法實現(xiàn)，其主要包括降維近似統(tǒng)計矩點估計法[10-13]、稀疏網(wǎng)格配點法[14]以及容積積分法[15]。其中，研究者對降維近似統(tǒng)計矩點估計法進行了大量的研究，通過引入不同的降維近似模型，分別發(fā)展了單變量、雙變量和三變量降維近似統(tǒng)計矩點估計法[10-13]。對于結構整體可靠度問題，等價功能函數(shù)蘊含了多個失效模式的相關性，變量間的交互影響一般較為復雜，采用單變量與雙變量降維近似統(tǒng)計矩點估計法往往不能達到滿意的精度[13]。研究發(fā)現(xiàn)，三變量降維近似統(tǒng)計矩點估計法能夠較為準確地考慮結構整體可靠度問題中各隨機變量的交互影響，適用于結構整體可靠度的評估[13]。但三變量降維近似統(tǒng)計矩點估計法所需的結構分析次數(shù)較多，如何在保證精度的前提下改善其計算效率有待進一步提高。

一般而言，在獲取功能函數(shù)的統(tǒng)計矩后，采用Pearson 分布[16]、立方正態(tài)[4]、鞍點近似[13]、移位廣義對數(shù)分布[17]以及最大熵原理[18]等方法可獲得功能函數(shù)的概率密度函數(shù)。其中，最大熵原理具有統(tǒng)一形式的最大熵概率密度函數(shù)，能夠較為精確地擬合單峰與多峰概率密度函數(shù)，被認為是擬合概率密度函數(shù)最具有無偏估計的一類方法[18]，因而學者們將該方法廣泛應用于工程結構可靠度分析[18-22]。

為此，本文在結構整體可靠度問題分類的基礎上給出其對應等價功能函數(shù)的統(tǒng)一描述，進而提出了將有效維度兩步分析法與共軛無跡變換法結合的改進統(tǒng)計矩點估計法，并結合最大熵原理發(fā)展了一類兼顧精度與效率的結構整體可靠度分析方法。

1 結構整體可靠度問題描述

根據(jù)多失效模式產(chǎn)生的來源不同，可將結構整體可靠度問題分為兩類：結構經(jīng)典整體可靠度問題和結構一般整體可靠度問題[23]。結構經(jīng)典整體可靠度問題主要研究理想彈塑性結構倒塌或形成可變機構的概率，即結構經(jīng)典整體可靠度問題的失效模式由結構倒塌這一物理機制進行識別[7]；而結構一般整體可靠度問題的失效模式可以由失效準則直接確定。例如考察某多層結構(層數(shù)n≥2)的可靠度問題，任意某層的位移大于某限值時即認為結構失效，則該結構體系具有n個失效判別準則，相應地對應著n個失效模式[23]。因此，結構一般整體可靠度的多失效模式源于可靠度問題本身。

1.1 結構經(jīng)典整體可靠度

當理想彈塑性結構承受完全相關的隨機荷載作用時，隨機向量Θ={ΘL,ΘS}T，其中ΘS和ΘL分別表示結構的隨機參數(shù)和荷載的隨機參數(shù)。定義ΘL,1為參考荷載，則ΘL可表示為：

其中，r(ΘL)為荷載向量。根據(jù)非線性發(fā)展過程[24]，彈塑性結構整體可靠度可等效為結構的極限承載力Fmax大于所施加荷載ΘL,1的概率。因此，結構經(jīng)典整體可靠度對應的功能函數(shù)可等價描述為：

對于非完全相關荷載下的結構經(jīng)典整體可靠度分析，詳見文獻[25]。

1.2 結構一般整體可靠度

當失效模式已知或已識別出時，此類結構整體可靠度問題稱為結構一般整體可靠度。對于結構一般整體可靠度問題，其涉及的各失效模式及其邏輯關系均可根據(jù)預定的失效準則直接確定。

根據(jù)等價極值事件[8]，串聯(lián)系統(tǒng)、并聯(lián)系統(tǒng)及混聯(lián)系統(tǒng)這三類邏輯關系的結構一般整體可靠度問題的等價功能函數(shù)可表示為：

式中：Gi(·)為串聯(lián)系統(tǒng)或并聯(lián)系統(tǒng)中第i個單元的功能函數(shù)；Gij(·)為第i個串聯(lián)子系統(tǒng)中第j個單元的功能函數(shù)。

綜上，結構經(jīng)典整體可靠度與結構一般整體可靠度對應的等價功能函數(shù)可統(tǒng)一表示為：

2 基于有效維度兩步分析法和共軛無跡變換法的改進統(tǒng)計矩點估計法

2.1 變量的獨立標準正態(tài)變換

對于隨機向量Θ={Θ1, Θ2,···,ΘN}，通過引入Nataf 變換[26]，可將其轉換到獨立標準正態(tài)空間，得到新的隨機向量U={U1,U2,···,UN}。于是，兩類結構整體可靠度問題對應的功能函數(shù)可改寫為：

式中，N-1(·)表示Nataf 變換的逆變換。

2.2 傳統(tǒng)的三變量降維近似統(tǒng)計矩點估計法

根據(jù)文獻[11]，Z的原點矩MZ,k可近似計算如下：

雖然式(7)將高維積分轉換為多個低維積分使得結構分析次數(shù)得到縮減，但三維和二維積分所需的結構分析次數(shù)仍較多。據(jù)此，本文提出有效維度兩步分析法，并結合共軛無跡變換法進一步提高式(7)的計算效率。

2.3 有效維度兩步分析法

其中，?為衡量影響程度的閾值，本文中取?=5%[15]。

采用有效維度分析法進行交互項評估僅需2N+2 次結構分析，具有較高的效率。然而，該方法存在交互項誤判問題[21]。鑒于此，本文首先通過有效維度分析法進行交互項預判，在此基礎上引入一次交叉項判別準則[12]進行二次判定。定義Us和Ut的交互項判定的示性函數(shù)Qh,2(Us,Ut)為[12]：

本文在有效維度分析法的基礎上引入二次判定，故將此方法稱為有效維度兩步分析法。

2.4 改進統(tǒng)計矩點估計法

相比于傳統(tǒng)的高斯求積公式，研究者近年來在非線性濾波領域發(fā)展的共軛無跡變換法具備更高的計算效率[27-29]。本文引入9 階代數(shù)精度的共軛無跡變換法進行分量函數(shù)積分計算，其基本原理與過程如下。

其次，通過Isserlis 定理[27]，可建立如下矩約束方程：

式中：Uj為隨機變量；為求積節(jié)點的坐標；ν為隨機變量數(shù)目；κj為非負的整數(shù)且1≤κ1+κ2+···+κν≤9；n為除零點外求積節(jié)點的數(shù)量。求解式(19)，可獲取縮放變量ri和權系數(shù)wi。特別的，當ν=2 或3 時，縮放變量ri和權系數(shù)wi如表1所示。

表1 二維及三維系統(tǒng)共軛無跡變換方法對應的ri 和wiTable 1 ri and wi in conjugate unscented transformation method for two-dimensional and three-dimensional systems

據(jù) 此，令ui={us,i,ut,i}或ui={us,i,ut,i,uv,i}，ui與αi分別為二維或三維共軛無跡變換方法的節(jié)點與權系數(shù)，則式(7)中的二維及三維積分可由下式確定：

其中，N1=4N(N-1)(N-2)/3。

對于式(7)中的一維積分，直接采用高斯求積公式進行計算：

式中：uq,i和ωi分別為Gauss-Hermite 積分的求積節(jié)點和權系數(shù)；本文取m=5。

對于式(7)中的二維積分，采用下式計算：

對于式(7)中三維積分，若Qh(Us,Ut,Uv)=1，則：

若Qh(Us,Ut,Uv)=0，則存在如下三種情況：

1) 當任意一對變量的交叉項不存在，如Qh,2(Us,Ut)=0 時，經(jīng)推導，式(7)中三維積分可簡化為：

2) 當任意兩對變量的交叉項不存在，如Qh,2(Us,Ut)=0 且Qh,2(Us,Uv)=0 時，經(jīng)推導，式(7)中三維積分可簡化為：

3) 當三對變量的交叉項均不存在，經(jīng)推導，式(7)中三維積分可簡化為：

3 基于最大熵原理的結構整體可靠度分析

3.1 最大熵原理

結構響應統(tǒng)計矩獲取后，通過最大熵原理進行結構整體可靠度評估[18-19]。以功能函數(shù)Z的熵取最大值為目標函數(shù)，以功能函數(shù)Z的前4 階原點矩MZ,k(k=0,1,···,4)為約束條件，建立優(yōu)化模型為：

對熵函數(shù)H引入Lagrange 乘子，構造Lagrange函數(shù)，得到：

通過求解非線性方程可獲取Lagrange 系數(shù)λ0, λ1, ···, λ4，從而估計響應函數(shù)Z的密度函數(shù)：

則結構的失效概率與可靠指標為：

3.2 算法流程

基于本文方法的結構整體可靠度分析基本步驟如下：

1) 通過式(5)將結構整體可靠度問題的功能函數(shù)統(tǒng)一表示。

2) 通過式(6)將任意分布轉換到標準獨立正態(tài)空間。

3) 根據(jù)式(7)計算功能函數(shù)的前四階原點矩：

a) 根據(jù)式(21)和式(22)分別計算單變量函數(shù)和雙變量函數(shù)的積分值；

b) 根據(jù)Qh(Us,Ut,Uv)的結果通過式(23)～式(26)計算三變量函數(shù)的積分值；

c) 將單變量函數(shù)、雙變量函數(shù)及三變量函數(shù)的積分值代入式(7)中，計算功能函數(shù)的前四階矩。

4) 基于求得的功能函數(shù)前四階矩，按照式(27) ～式(30)計算結構的失效概率和可靠指標。

4 算例分析

本節(jié)通過結構一般整體可靠度與結構經(jīng)典整體可靠度兩個算例驗證本文方法的計算效率、精度以及適用性，分別將本文方法分別與Monte Carlo 模擬方法(MCS)、寬界限法(WBM)[30]、窄界限法(NBM)[31]、5 點Gauss-Hermite 積分的雙變量降維近似統(tǒng)計矩點估計法(BDRM)[11]及5 點Gauss-Hermite 積分的三變量降維近似統(tǒng)計矩點估計法(TDRM)[11]計算的統(tǒng)計矩及可靠度進行了對比。本文將MCS 計算結果視為標準解，其他各方法的相對誤差為：

式中：Value為對應方法的計算結果；MCV為Monte Carlo 方法計算的結果。

4.1 算例1. 結構一般整體可靠度問題

考察一個如圖1 所示的六層鋼筋混凝土框架結構，框架的跨度l=7.5 m，層高h=3 m，梁與柱的截面尺寸分別為300 mm×400 mm 和500 mm×500 mm，水平荷載F可表示為：

圖1 六層鋼筋混凝土框架結構Fig. 1 Six story reinforced concrete frame structure

考慮梁柱的彈性模量Eb和Ec，豎向荷載P1和P2以及水平荷載F0為隨機變量，各隨機變量的統(tǒng)計特性如表2 所示。考慮結構體系的最大層間位移角超過某限值為失效狀態(tài)，該可靠度問題可以描述為：

表2 算例1 中隨機變量的統(tǒng)計特征Table 2 Statistical characteristics of random variables for Example 1

式中：Xj(j=1,2,···,6)代表第j層與第j-1 層的相對位移； φB為層間位移角限值，本文中 φB取值為1/50。

該問題可歸類為結構一般整體可靠度問題，根據(jù)式(6)，其單一等價功能函數(shù)可以表示為：

通過有效維度兩步分析法易知，Qh(U1,U3,U4) =Qh(U2,U3,U4) =Qh(U3,U4,U5) =Qh,2(U3,U4) =0，進而，可計算Z的前四階統(tǒng)計矩和可靠指標，其計算結果如表3 和表4 所示，表中?表示該方法與Monte Carlo 方法的相對誤差。

表3 算例1 中極限狀態(tài)函數(shù)前四階矩結果Table 3 Results of first four order moments of limit state function for Example 1

表4 算例1 中各類方法的可靠指標與有限元分析次數(shù)Table 4 Reliability index and number of finite element analysis of each method for Example 1

由表3 可以發(fā)現(xiàn)，TDRM 和本文方法估計的前四階矩均具有較高的精度，其最大相對誤差低于2%；雖然BDRM 亦能精確地估計前兩階矩，然而其估計的三階矩具有6.32%的相對誤差。就可靠指標的計算精度與效率而言，由表4 易知：采用經(jīng)典的WBM 僅能給出較大的可靠指標界限值，其上、下界限可靠指標均具有較大的相對誤差；雖然NBM 計算的可靠指標的界限值較窄，但其仍難以定量評估結構整體可靠度；通過BDRM估計的可靠指標相對誤差達到12.5%，亦難以精確評估結構整體可靠度；幸運的是，本文方法和TDRM 能夠精確地估計可靠指標，其最大誤差低于1%，但本文方法所需的有限元分析次數(shù)為650 次，相較于TDRM 所需的1528 次有限元分析，本文方法計算效率提高2.35 倍左右。綜合計算精度與效率，推薦采用本文方法進行結構一般整體可靠度的評估。

4.2 算例2. 結構經(jīng)典整體可靠度問題

考察如圖2 所示的理想彈塑性單層單跨鋼框架結構的整體可靠度問題，其中，材料的彈性模量E為2.0×105MPa，各構件截面積均為0.09 m2。考慮節(jié)點荷載F0以及各構件截面抗彎承載力M1、M2和M3為相互獨立的隨機變量，其統(tǒng)計特性如表5 所示。當鋼框架結構變成可變機構時認為結構失效。

圖2 單層單跨的理想彈塑性鋼框架結構 /mFig. 2 One-story and one-bay perfectly elastoplastic steel frame structure

表5 算例2 中隨機變量的統(tǒng)計特征Table 5 Statistical characteristics of random variables for Example 2

根據(jù)問題描述可知，該問題為結構經(jīng)典整體可靠度問題。根據(jù)式(6)，該問題的單一等價功能函數(shù)可表示為：

通過有效維度兩步分析法可得，Qh(U1,U2,U3) =Qh(U1,U3,U4) =Qh,2(U1,U3) = 0，進而，采用各類方法計算出Z的前四階統(tǒng)計矩和可靠指標如表6 和表7 所示。

表6 算例2 中極限狀態(tài)函數(shù)前四階矩結果Table 6 Results of first four order moments of limit state function for Example 2

表7 算例2 中各類方法的可靠指標與有限元分析次數(shù)Table 7 Reliability index and number of finite element analysis of each method for Example 2

通過表6 和表7 可以看出，采用BDRM 估計的三階矩和可靠指標的相對誤差均較大，其三階矩的相對誤差高達38.97%；采用經(jīng)典的WBM 和NBM 僅能給出可靠指標的界限值難以精確量化，且需進行結構失效模式識別步驟；而采用TDRM及本文方法估計的統(tǒng)計矩和可靠指標均能夠得到較為精確的結果，其相對誤差均小于2%，但本文方法的計算效率較TDRM 提升2.5 倍多，效果顯著。

5 結論

結構整體可靠度的定量評估是可靠度領域的基本問題之一。本文結合有效維度兩步分析法和共軛無跡變換法，發(fā)展了高效的改進統(tǒng)計矩點估計法；進而，結合最大熵原理和結構整體可靠度等價描述法，提出了適用于結構經(jīng)典和一般整體可靠度的分析方法。結果表明：

(1) 本文方法對兩類結構整體可靠度問題均表現(xiàn)出較高的精度水平。

(2) 相比于傳統(tǒng)的三變量降維近似統(tǒng)計矩點估計法，本文方法能在保證精度的前提下具有更高的效率；

(3) 雙變量降維近似統(tǒng)計矩點估計法對于兩類整體可靠度問題均具有較高的相對誤差，難以對結構整體可靠度問題進行精確評估。