結構可靠度指標的最優化方法及MATLAB實現

李 艷，劉肖凡，黃永超

(武漢工業學院土木工程與建筑學院，湖北武漢430023)

在結構可靠度計算過程中，需要建立大量的概率統計模型以及優化相關的算法。這一系列算法的實現可以通過一些基礎的計算編程語言如Fortran和C語言來完成。然而要想通過這一方法來得到較為準確的計算結果，則必需要熟練掌握這兩門編程語言，這一點對廣大的非計算機專業人員來說具有較大的難度。除此之外，目前市面上無可靠度計算的專用軟件，因此要解決這一計算難題必需求助于功能強、效率高、綜合性較強的計算軟件。美國Math Works公司推出的Matlab軟件語法簡單、用戶界面友善、矩陣運算功能強大，其已經成為數值分析計算、數值仿真、信號處理等高級課程的基本運算工具［1］。在可靠度指標優化計算過程中需要建立大量與概率、統計和最優化方法相關的數值計算方法在Matlab環境中均可很容易實現。

1 從一次二階矩到最優化方法

可靠度指標的常用計算方法有一次二階矩法、蒙特卡羅法和響應面法［2］。本文從獨立正態分布變量在極限狀態方程為線性時可靠度指標β的幾何涵義出發，運用最優化原理，提出了不求導數的最優化分析方法。在標準正態坐標系中可靠度指標的幾何圖解如圖1所示。

圖1(a)為二維坐標系中標準正態坐標系中可靠度指標β的幾何圖解。假設R和S是相互獨立并且服從正態分布的變量，其極限狀態方程如下：

Z=R-S=0.

圖1 在標準正態坐標系中可靠度指標的幾何圖解

其中，μR、σR、μS、σS分別為R和S的均值和標準差。由式(1)對這一簡單的原坐標系進行標準化，可知(直線距離)，其中點是標準正態坐標系的原點，(μR、μS)是在原坐標系中的變量均值點，設計驗算點為P*并滿足極限狀態方程［3］：

在極限狀態方程具有多個正態變量的條件下，變量如下：

同理對其進行標準化的處理如下：

標準化處理結果為：

標準正態坐標系下，可靠度指標等于原點到極限狀態面的最小距離：

圖1(b)所示為三個變量的情況。

可靠度指標的最優化問題無論在什么樣的前提下均可轉換成原點到極限狀態曲面最短距離的問題。且這一過程可以方便又快捷的在Matlab中實現。

2 優化方法求可靠度指標數學模型的建立

首先，設n個影響結構的可靠度服從任意分布的獨立隨機變量，并假設該極限狀態函數為：

利用笛卡爾隨機空間內一次二階矩理論常用處理辦法，先通過拉克維茨—菲斯萊法(簡稱R—F法)或根據變量特性確定的其他方法將非正態變量正態化。前者的基本假設是正態化處理后的當量正態分布滿足如下條件：在設計驗算點的概率分布函數值以及概率密度函數值分別等于原隨機變量在處所對應的概率分布函數值和概率密度函數值［4］。如圖2所示。

圖2

通過上述對當量的正態化處理即可得到相應的等效分布的均值和標準差：

由于驗算點的未知性，只能通過對可靠度指標的優化處理才能確定。先設β為極限狀態曲面上點P(X1…Xn)的函數，再求得β在該區域的最小值，從而驗算點和可靠度指標的值均可得出［5］。于是，求式(6)的最小值可轉化為求下式的最小值：

因此可以得出如下優化數學模型：

從式(5)中可得出通過拉克維茨—菲斯萊法正態化后的均值和標準差對于不同的分布有不同的簡化方式，具體如下。

(1)正態分布

(2)對數正態分布

由于極限狀態函數中的一個變量可以由其他變量來表示：

于是約束優化問題式(7)可轉化為無約束優化問題：

程序框圖如圖(3)所示。

圖3 程序框圖

3 優化數學模型的Matlab實現

3.1 約束優化法

該問題可靠度指標的求解可歸結為如下約束優化數學模型［6］：

MATLAB程序如下。

global Mu Sigma

Mu=［102.520000000.008］，Sigma=［0.405000000.0015］

X0=［4.5］

A=［］，b=［］，Aeq=［］，beq=［］，lb=［］，ub=［］

［X，fval，exitflag ，output］=fmincon(@bata2，A，b，X0，Aeq，beq，lb，ub，@st)

bata=sqrt(fval

Pf=cdf(‘norm’，-bata，0，1)

Function CC=bata2(X)

global Mu Sigma

CC=((X(1)-Mu(1))/Sigma(1))^2+((X(2)-Mu(2))/Sigma(2))^2+((X(3)-Mu(3))/Sigma(3))^2+((X(4)-Mu(4))/Sigma(4)^2 function［c，ceq］=st(X)

c=［］

ceq=X(2)/360-0.0069*X(1)*X(2)^4/X(3)*X(4

程序中要說明的是：(1)將global函數中的變量Mu，Sigma設為全局變量;(2)初始迭代點為均值;(3)必須將目標函數(β2)寫入bata2.m子函數［7］。計算結果見表1。

3.2 無約束優化法

無約束優化法求解可靠度指標需建立如下約束優化數學模型：

取消約束條件子函數如下：

global Mu Sigma

Mu=［102.520000000.008］，Sigma=［0.405000000.0015］

X0=［4.5］

［X，fval，exitflag，output］ =fminsearch(@bata2w，X0)

bata=sqrt(fval);

Pf=cdf('norm'，-bata，0，1);

X1=X3*X4/(2.484X2.^3);

Function CC=bata2w(X);

global Mu Sigma

X1=X3*X4/(2.484X2.^3)

CC=((X(1)-Mu(1))/Sigma(1))^2+((X(2)-Mu(2))/Sigma(2))^2+((X(3)-Mu(3))/Sigma(3))^2+((X(4)-Mu(4))/Sigma(4)^2;

計算結果如表1。

表1 計算結果

由表1可知，約束優化法和無約束優化法計算結構可靠指標，兩者計算結果精度接近。

4 結論

最優化計算方法不用求任何函數的導數即可直接求得結構可靠度指標，該方法易于編制計算程序，計算結果精度高，通用性好，比較適用于巖土工程可靠度問題的分析［8］。而采用Matlab語言進行結構可靠度的數值計算，可充分發揮利用矩陣運算功能及各種工具箱的快速和便捷，大大提高編程和計算的效率。這將對結構可靠度理論的研究起到積極推進作用。

［1］ 張亮，趙娜.用MATLAB實現JC法計算結構可靠度程序［J］.電腦知識與科技，2009(10).

［2］ 貢金鑫.工程結構可靠度計算方法［M］.大連：大連理工大學出版社，2003(9).

［3］ 冷伍明.基礎工程可靠度分析與設計理論［M］.長沙：中南大學出版社，2000(4).

［4］ 趙國藩，金偉良，貢金鑫.結構可靠度理論［M］.北京：中國建筑工業出版社，2000.

［5］ 李志華，張光海，康海貴.基于Matlab優化工具箱的工程結構可靠度計算［J］.四川建筑科學研究，2005(6).

［6］ 石博強，趙金.MATLAB數學計算與工程分析范例教程［M］.北京：中國鐵道出版社，2005(5).

［7］ 徐華，朝澤剛，劉勇.用MATLAB實現結構可靠度計算［J］.科協論壇，2009.

［8］ 陸舟.可靠度指標優化計算法在巖土工程可靠度研究中的應用［J］.巖土工程界，2004(1).