有限實驗數據下工程結構概率建模及隨機響應分析

張瑞景, 戴鴻哲

(哈爾濱工業(yè)大學 土木工程學院,黑龍江 哈爾濱 150090)

在隨機工程問題中,合理考慮輸入參數的隨機性是獲得可靠結構響應的關鍵[1-2]。隨機輸入在工程應用中普遍存在,如材料特性、外部荷載、結構系統(tǒng)的初始和邊界條件等[3-4]。工程結構的隨機參數往往具有明顯的時空變異性和非高斯性,如巖土工程中的土壤參數和地下水高度、結構工程中的風荷載和地震激勵以及水文中的降水和蒸發(fā)量等。此時,結構參數應建模為非高斯隨機場[5-6];另一方面,在實際應用中,由于傳感器的存儲能力有限或增加觀測的成本過高等原因,往往只能獲得不確定的輸入參數有限的實驗數據[7-10]。因此,從有限的實驗數據中合理地建模非高斯輸入參數,并進一步進行后續(xù)隨機響應分析具有重要的工程價值。近年來,為了解決實驗數據下工程結構概率建模及隨機響應求解的問題,學者們發(fā)展了各類方法。其中,最具潛力的是基于多項式混沌(polynomial chaos, PC)展開的方法。此類方法的基本思想是首先對非高斯結構參數進行Karhunen-Loeve(KL)展開,然后再將得到的KL隨機系數進行PC展開。獲得結構參數的PC模型后,即可利用基于PC的響應求解技術進行結構隨機響應的求解[11-12]。

雖然有諸多優(yōu)勢,但基于PC展開的方法目前仍存在許多挑戰(zhàn)。此類方法的第1個挑戰(zhàn)為隨機建模。具體來說,要想獲得合理的輸入參數模型,必須準確建模多維KL隨機系數的非線性相關,并估計其聯(lián)合概率密度。此外,構造從KL隨機系數到PC變量的多維非線性變換還涉及到大量的高維積分。為此,有學者把KL隨機系數假設為相互獨立,從而大大簡化了隨機建模[13-14]。為了考慮KL隨機系數的非線性相關性,Ghanem等[15-16]利用極大熵分布、直方圖估計器來估計KL隨機系數的聯(lián)合概率密度函數,隨后再使用Rosenblatt變換構建KL隨機系數的Hermite PC表示。Soize[17]使用核密度估計器(kernel density estimation, KDE)估計KL隨機系數的聯(lián)合概率密度函數,然后通過蒙特卡羅積分求解KL隨機系數PC展開的系數。使用KDE的優(yōu)勢在于其可避免其他密度估計器在高維密度估計中的維數災難[18]。

除隨機建模外,PC類方法的另一挑戰(zhàn)則為后續(xù)隨機響應分析的求解效率。眾所周知,對于服從Wiener-Askey族分布的隨機輸入而言,使用Wiener-Askey PC可實現(xiàn)結構響應的最優(yōu)收斂[19]。然而,在實驗數據有限的情況下,KL隨機系數的分布通常在Wiener-Askey族之外,繼續(xù)采用Wiener-Askey PC會導致結構響應分析的效率降低;另一方面,即使已構造出能使響應最優(yōu)收斂的PC展開,由于PC變量存在復雜的相關性,后續(xù)的結構響應分析仍具挑戰(zhàn)[20-21]。不難看出,要想使PC類方法可應用于實際工程結構分析中,關鍵為以下2方面:首先,要能夠從實驗數據中精確地建立基于PC的結構參數模型;其次,后續(xù)響應傳播的效率應當足夠高,以適應處理高維和大規(guī)模問題的工程需求。

本文開發(fā)了一種實驗數據下的PC類結構隨機分析方法。提出了一種KDE來估計KL隨機系數的聯(lián)合分布,從而得到結構參數的隨機模型。進一步構造了KL隨機系數的任意多項式混沌(aPC)展開,并利用一種基于D-optimal的加權插值法來求解結構的隨機響應。

1 實驗數據下非高斯結構參數的隨機模型

實驗數據下的隨機結構分析第1步是合理建模非高斯結構輸入參數。基于KDE的隨機場建模方法,提出了一種可準確描述結構參數非高斯行為的隨機場模型。

1.1 有限實驗數據下非高斯結構參數的KL展開

令w(xi,θ),i=1,2,…,M為隨機場w(x,θ)在M個位置xi處的一條實驗數據,則N條獨立的實驗數據構成了隨機場w(x,θ)的實驗數據集W。其中,W={w(xi,θj)},i=1,2,…,M,j=1,2,…,N。隨機場W(θ)={w(xi,θ)}的m階截斷KL展開可表示為:

(1)

(2)

式中U為所有元素均為1的N維列向量。

KL隨機系數ξi(θ)則是由其實驗數據表征,即:

(3)

i=1,2,…,m;j=1,2,…,N

其中Wj=[w(x1,θj),w(x2,θj),…,w(xM,θj)]。

(4)

1.2 基于改進KDE的隨機模型

為了從實驗數據中準確地建模非高斯結構參數,本文提出了一種用于KL隨機系數聯(lián)合分布估計的新型KDE。該新型KDE的優(yōu)點為得到的隨機場模型準確可重構非高斯參數的最關注的2個概率特征,即二階相關性和邊緣分布。由于隨機場W(θ)的邊緣分布實質上是由一維KL隨機系數的線性組合而得到的,因此,為了匹配結構輸入參數的邊緣分布,式(4)中的帶寬s應根據一維KDE的選取準則來確定。一維KL隨機系數的概率密度為:

(5)

(6)

式中wi是帶寬s(i)的權重。wi的值隨i的增加而減小,即i越小的帶寬s(i)對ssh的貢獻越大,與特征值更大的KL隨機系數對隨機場邊緣分布貢獻更大這一現(xiàn)象契合。確定帶寬ssh后,KL隨機系數ξ(θ)聯(lián)合密度變?yōu)?

(7)

一旦KL隨機系數的聯(lián)合分布通過式(7)確定,輸入參數的非高斯模型即可由式(1)確定。

1.2.1 新型隨機場模型的統(tǒng)計性質

(8)

(9)

(10)

由式(10)可得,本文隨機模型的二階相關性為:

(11)

2 結構響應分析的任意多項式混沌法

盡管PC類的方法可以求解結構的隨機響應,由于在實驗數據下傳統(tǒng)的Wiener-Askey PC可能會導致較低的計算效率,PC類方法在實際工程結構隨機分析中的應用仍存在較大局限[22]。本文結合了結構響應分析的PC類方法,然后將KL隨機系數的聯(lián)合分布選為aPC的權函數,從而導出了結構輸入參數的aPC表示。隨后,進一步發(fā)展了基于aPC展開的高效隨機響應分析方法。

2.1 基于多項式混沌的結構隨機響應分析

在基于PC的結構隨機分析中,結構響應通常與結構隨機參數投影到相同的PC子空間中。結構參數隨機模型中KL隨機系數ξ(θ)的PC展開為:

(12)

式中:η(θ)為m維的PC變量;P=(m+p)!/(m!p!)為PC展開的截斷項數,p為m維標準正交多項式Ψj(·)的階數;αij則是PC系數。將式(12)代入式(1),即可得到非高斯輸入場的PC展開為:

(13)

(14)

式(14)也可寫為Y(η)?αTΨ(η)。

2.2 任意多項式混沌展開的構造

E[φi(ξ)φj(ξ)]

(15)

構造Ψ(ξ)的關鍵在于求解式(14)中的多維積分。一般地,該積分應使用蒙特卡羅積分進行近似:

(16)

式中Ξ(k)為KL隨機系數ξ(θ)的第k個樣本。雖然樣本Ξ(k)可通過MCMC方法產生,其中涉及的數以萬計次的聯(lián)合概率密度計算導致抽樣效率低。此外,MCMC樣本之間存在的自相關性會大大降低蒙特卡羅積分的精度。為了高效高精度地構造矩陣G,本文使用了文獻[24]中提出的KL隨機系數的高效抽樣方法。式(7)可寫為:

(17)

2.3 基于aPC展開的結構隨機響應分析

aPC展開構造完畢后,需開發(fā)求解結構隨機響應的算法。眾所周知,基于aPC的結構響應分析關鍵在于求解式(14)中的aPC系數αk。大體而言,求解系數αk的方法大體可分為侵入式和非侵入式2類,由于非侵入式方法中的插值法可直接調用第三方軟件進行結構分析[25],因此本文的研究主要圍繞插值法來開展。插值法的精度和魯棒性很大程度上取決于配置點的選擇(即實驗設計)。目前,獨立PC變量的實驗設計方案已日趨完善,而相關PC變量的實驗設計方案還鮮有報道。基于aPC的響應分析的關鍵在于抽取aPC變量的蒙特卡羅樣本。一旦獲得aPC變量的樣本,就可直接開發(fā)出aPC變量的實驗設計方案,并進行后續(xù)的響應求解。

注意到aPC變量本質上是KL隨機系數,而KL隨機系數服從混合高斯分布,故aPC變量的獨立樣本可直接產生。因此,各種獨立PC變量下的實驗設計方案即可直接推廣至aPC變量的實驗設計。為實現(xiàn)結構響應的高精度估計,本文提出了一種D-optimal加權插值方法。首先,本文將估計式(14)中aPC系數的問題列式為加權插值形式:

(18)

根據D-optimal理論,ΞED可通過求解優(yōu)化問題確定:

(19)

(20)

ΞED=ΞC(πED,:)

(21)

其中,πED=π(1:NED)。一旦確定了配置點ΞED,結構響應的aPC系數即可通過式(18)估計。

2.4 有限實驗數據下工程結構隨機分析方法

本文提出的實驗數據下的新型結構隨機分析方法總結如下:

1)對非高斯結構參數的實驗數據進行KL展開,得到均值向量,截斷后的特征值、特征向量及KL隨機系數ξ(θ)的樣本實現(xiàn);

3)產生KL變量的獨立樣本,構造aPC基函數;

4)再次產生KL變量(aPC變量)的樣本實現(xiàn),隨后利用選列主元的QR分解確定用于加權插值法的D-optimal配置點。通過式(18)估計隨機響應的aPC系數并利用式(14)求解結構的隨機響應。

步驟1)、2)旨在構建實驗數據下結構輸入參數的新型隨機模型。為高效確定后續(xù)的隨機響應,步驟3)將模型中的KL隨機系數進一步用aPC展開。最后,步驟4)開發(fā)了有限數據下基于aPC的隨機響應分析方法,從而形成了有限數據下工程結構的隨機建模和后續(xù)響應傳播的統(tǒng)一框架。注意到,由于步驟2)中的新型隨機場模型考慮了輸入參數分布與KL隨機系數分布之間的內在關系,因此本文的隨機模型可以準確地重構結構參數的邊緣分布和二階相關性,從而有效地構建了結構輸入參數的非高斯特征。此外,本文發(fā)展了結構響應求解的D-optimal加權插值法,實現(xiàn)了結構隨機響應的高效高精度求解。由于本文方法可對結構輸入參數進行合理地隨機建模,并進行高效的結構響應分機分析提供了一個有效的框架。

3 工程算例

在實際工程中,工程結構的不確定參數往往分為2大類:1)工程結構自身的機械和幾何特性。例如,由于成分和制造工藝的原因,復合材料的楊氏模量通常表現(xiàn)出空間變異性。在實際工程中,人們不可能完全地確定楊氏模量的概率信息,只能得到一組從有限數目的試件中識別出的楊氏模量的實驗數據;2)工程結構參數不確定性則源于外荷載。例如,地震動是一種具有強烈時空變異性的自然災害,在特定的場地條件下,記錄到的真實地震動時程往往十分有限。本文通過非線性結構在有限地震動數據下的隨機分析的來說明所發(fā)展方法的應用。有限的加速度時程記錄從太平洋地震工程研究中心建立的NGA強震數據庫中根據如下特定場地標準選擇出來:地震震級5≤M≤6,震源距離1 km≤D≤20 km,中硬土,土層的剪切波速Vs≥600 m/s。滿足此準則的地震動記錄共102條,將處理后的18 s的時程離散為1 801個點,步長Δt=0.01 s,如圖1所示。為保留95%以上的總能量,此算例截斷項數取為50。

圖1 根據準則選出的102條地震動時程

方差是反映隨機地震動統(tǒng)計特性的一個重要指標。為驗證本文模型隨機地震動的性能,圖2給出了地震動時程的方差。由于隨機地震動方差在整個時程中波動劇烈,對整個時程而言,地震動方差的值不在一個數量級,為便于觀察,本文僅展示了0～8 s的結果。可以看出,傳統(tǒng)KDE模型給出的邊緣分布明顯偏離了觀測值,而本文建模方法給出的地震動時程的方差與觀測值良好匹配。圖3給出了t=1.5 s、t=6.5 s、t=11.5 s和t=16.5 s處地震動時程的觀測邊緣分布與模型邊緣分布的Q-Q圖,即Quantile-Quantile圖。傳統(tǒng)KDE模型給出的邊緣分布明顯偏離了觀測到的邊緣分布,無法捕捉到地震動的非高斯特征。相比之下,本文的隨機模型給出的邊緣分布與觀測結果吻合良好,說明了本文建模方法的有效性。

圖2 0～8 s地震動的方差

圖3 4個典型時刻下地震動邊緣分布的Q-Q圖

為了驗證所發(fā)展的結構隨機響應分析方法的性能,本文進一步研究了20自由度滯回非線性結構在受隨機地震動激勵下的隨機響應。其中,結構的質量和層間剛度參數如表1所示。

表1 20自由度非線性結構的質量和層間剛度

結構的非線性行為由Bouc-Wen滯回模型表征,其控制方程為:

(22)

在本算例中,式(16)中用于構造aPC基函數的KL獨立樣本個數選取為K=105,D-optimal實驗設計中候選樣本容量取為NC=104。將所得系統(tǒng)響應與2×104次MCS給出的參考值進行比較,評估了所提響應分析方法的精度。分析結果顯示,2階aPC 展開即可對系統(tǒng)隨機響應進行高精度近似。此時,只需對式(22)中的非線性滯回系統(tǒng)進行(50+2)!/(50!2!)=1 326次確定性評估。因此,本文方法可對受隨機地震動激勵的強非線性系統(tǒng)進行高效高精度的響應分析。為更直觀說明本文方法的求解精度,部分隨機響應分析的結果展示如下。圖4和圖5分別給出了1階和2階aPC展開所得的第2層隨機位移響應的均值和方差及相應的絕對誤差,而圖6給出了第15層在t=1.5 s時的位移響應及相應的絕對誤差。可以看出,由于結構的強非線性行為,采用1階aPC展開并不能充分地近似結構響應。然而,當使用2階aPC展開時,本文方法對系統(tǒng)隨機響應的近似精度迅速上升,實現(xiàn)了系統(tǒng)響應的高精度預測。

圖4 第2層位移響應的均值及均值的絕對誤差

圖6 t=1.5 s第15層位移響應的概率密度函數及絕對誤差

4 結論

1)本文開發(fā)的隨機模型考慮了結構參數分布與KL隨機系數分布之間的內在關系,可準確地重構輸入參數的邊緣分布和二階相關性,從而有效地捕捉結構參數的非高斯特征。

2)本文構造的輸入參數任意多項式混沌表示配合所發(fā)展的D-optimal加權插值法可實現(xiàn)非線性系統(tǒng)隨機響應的高效高精度評估,從而為有限實驗數據下工程結構的隨機分析提供了一個有效的框架。

3)算例驗證了本文方法在隨機維數為50維的問題中的有效性。但是,此方法在更高維度隨機問題中的應用有待進一步研究。