概率功能度量求解的共軛梯度步長調節(jié)法

易 平， 謝東赤

（大連理工大學 建設工程學部，大連116024）

1 引 言

工程結構中的材料屬性、幾何尺寸和外部荷載通常都具有隨機性，基于概率參數(shù)的優(yōu)化設計不可避免，也更符合工程實際。概率結構優(yōu)化設計PSDO（Probabilistic Structural Design Optimization）［1－5］通常表述為滿足概率約束條件下目標函數(shù)的最小化。直接的PSDO本質上是一個兩層次迭代問題，外部為確定性優(yōu)化，內部為可靠性分析。內層可靠度分析可采用功能度量法PMA（Performance Measure Approach）［4］和 可 靠 指 標 法RIA（Reliability Index Approach）［5］，其中功能度量法更加穩(wěn)定高效。

功能度量法通過在標準正態(tài)空間中搜索目標可靠指標球面上的最小功能目標點MPTP（Minimum Performance Target Point）獲得其概率功能度量。改進均值法 AMV（Advanced Mean Value）［6］常用于確定 MPTP，但對于一些非線性程度較高的功能函數(shù)，AMV迭代可能會陷入周期振蕩或混沌解。因此，相繼提出多種改善其收斂性的算法。Youn等［7－9］利用連續(xù)三次迭代點處的梯度信息和功能函數(shù)凹凸性，提出了混合均值法HMV（Hybrid Mean Value method），然而該算法對某些問題的求解仍因不收斂而失效［10］；Yang等［11］從混沌動力學的角度出發(fā)，利用穩(wěn)定轉換法對AMV迭代格式實施收斂控制；Meng等［12］改進混沌控制方法CC（Chaos Control），提出了混合混沌控制策略MCC（Modified Chaos Control），兩種混沌控制方法雖然提高了效率，但在實際應用中，尤其是對于高維問題，很難選取合適的迭代因子及對合矩陣。

共軛梯度法是求解非線性無約束優(yōu)化問題的常用方法，穩(wěn)定高效的特點使其在優(yōu)化領域得到了廣泛的發(fā)展。基于 RMIL（Rivaie Mustafa Ismail and Leong）共軛梯度方法［13］的思想并結合自適應步長調節(jié)策略［10，14，15］，本文提出了共軛梯度步長調節(jié)法 CGS（Conjugate gradient step length adjustment method）。多個數(shù)值算例表明，無論是凸函數(shù)還是凹函數(shù)，與其他求解方法相比，CGS方法更加高效且穩(wěn)健。

2 概率約束評估的功能度量法

2.1 基于概率功能度量的概率結構優(yōu)化設計

概率結構優(yōu)化設計數(shù)學模型通常為［5］

式中d為設計變量向量，x為隨機變量向量，gj（d，x）為第j號功能函數(shù)，gj（d，x）≤0代表失效，失效概率定義為

式中 Fgj（＊，d）為功能函數(shù)gj（d，x）的累積分布函數(shù)，f（x）為x的聯(lián)合概率密度函數(shù)，Pj，t為給定的許可失效概率，與目標可靠指標βj，t的關系為Pj，t≈Φ（－βj，t）。

式（1）的概率約束可采用概率功能度量的方式表達，其形式為

式中 Gpj（d）為第j號功能函數(shù)gj（d，x）的概率功能度量。采用概率功能度量表達概率約束，即為功能度量法（PMA）。

2.2 概率功能度量求解

概率功能度量求解需要在標準正態(tài)空間中進行，若初始隨機向量不是標準正態(tài)分布向量，則需借助Rosenblatt／Nataf等變換［5］，即u＝T（x），相應功能函數(shù)的變換為g（d，x）＝g（d，T－1（u））＝G（d，u）。

概率功能度量的求解本質上為標準正態(tài)空間中的一個優(yōu)化問題，其數(shù)學優(yōu)化模型可表示為

式中u＊為在目標可靠指標球面上功能函數(shù)值最小的點，稱為最小功能目標點（MPTP）。從而概率功能度量Gp為MPTP處的功能函數(shù)值，也就是目標可靠指標球面上功能函數(shù)的最小值。

改進均值法（AMV）［6］常用于確定 MPTP，其迭代格式為

式中uk代表第k次迭代點，n（uk）為點uk處的最速下降方向。AMV的迭代過程如圖1所示。

從式（6）可以看出，AMV方法中通過uk處的最速下降方向來實現(xiàn)MPTP的更新。對于非線性程度較低的功能函數(shù)以及凸函數(shù)，AMV方法具有良好的收斂性能；然而，對于非線性程度較高的功能函數(shù)，AMV迭代格式可能出現(xiàn)周期振蕩或混沌等不收斂現(xiàn)象［10］。針對AMV方法求解存在的弊端，HMV，CC 和 MCC 等方法相繼提出［8，11，12］，雖然改善了收斂性能，但其操作復雜，難以滿足實際應用的需求。

3 共軛梯度步長調節(jié)法

20世紀50年代初，為求解線性方程組提出了最早的共軛梯度法［16］。隨后，F(xiàn)letcher等［17］將共軛梯度法應用于求解非線性無約束優(yōu)化問題，稱之為Fletcher-Reeves（FR）共軛梯度法；Keshtegar［14，15］將共軛搜索方向用于可靠指標的求解，相比最速下降搜索方向，采用共軛梯度方向更具有魯棒性。基于共軛搜索方向的思想，并結合步長調節(jié)策略，本文提出了共軛梯度步長調節(jié)法（CGS）以求解概率功能度量。

3.1 基本的迭代格式

CGS方法的迭代過程如圖1所示。可以看出，uk為第k次迭代點，Sk為功能函數(shù)在點uk處的共軛梯度方向。在點uk處沿著Sk方向移動距離λk（λk＞0），可得到，

所以CGS迭代格式為

CGS的共軛梯度方向Sk為

式中 共軛梯度系數(shù)θk用于控制前一次共軛搜索方向Sk－1對當前搜索方向Sk的影響，如圖2所示，通常不同共軛梯度算法將定義不同的θk。若θk的值較大，會使新的共軛搜索方向太接近前一次搜索方向，影響收斂性［18］。因此，必須設定θk的合理取值以保證算法的收斂性和高效性。本文基于（Rivaie Mustafa Ismail and Leong，RMIL［13］）確定共軛梯度系數(shù)θk，其中θkRMIL定義為

圖1 AMV和CGS方法的迭代Fig.1 Iterative procedure of AMV and CGS

圖2 CGS方法迭代參數(shù)Fig.2 Iterative parameters of CGS

從圖2可以看出，負梯度方向（－Ug（uk））和共軛搜索方向Sk之間的夾角為φk，φk的值主要由共軛梯度系數(shù)θk確定。為使迭代逐漸收斂，每次迭代時應滿足0≤φk＜π／2［14］。

經(jīng)算例分析（見3.4節(jié)中算例1）發(fā)現(xiàn)，當共軛梯度系數(shù)θk采用時，迭代過程中出現(xiàn)θk大幅波動的情況。當θk＜0時，新的共軛搜索方向Sk與前一次搜索方向Sk－1幾乎反向，會延緩迭代進程，降低收斂效率；當θk太大時，則可能導致φk≥π／2，影響算法的收斂性。因此，為了避免θk大幅波動影響收斂效率，可以假定0≤θk≤‖Ug（uk）‖／‖Sk－1‖。結合，本文共軛梯度系數(shù)θk的取值為

3.4 節(jié)算例1表明，采用如上共軛梯度系數(shù)時，θk的取值相對穩(wěn)定，保證了收斂穩(wěn)定性，提高了效率。

3.2 確定步長的方法

除共軛搜索方向，步長調節(jié)是本文所提方法的另一重要措施。步長λ的選擇通常與功能函數(shù)的非線性程度有關。然而，在實際問題求解中，步長λ的取值與函數(shù)非線性程度的量化關系無法確立，因此通常先假定一個較大的步長λ，使迭代點快速向最優(yōu)點附近靠近，隨后逐漸調小λ以精確逼近最優(yōu)點。Keshtegar［15］在可靠指標求解中假定初始步長也即最大步長為

式中 常數(shù)M滿足10≤M≤100。Keshtegar提出的步長準則看似解決了初始步長選擇的難題，但實際求解過程中發(fā)現(xiàn)，根據(jù)Keshtegar步長準則給出的初始步長經(jīng)常因為太小造成收斂緩慢的局面。因此，通過限定初始步長下限值并結合Keshtegar步長準則，最終得到初始步長取值為

同時，取λ1＝λ0。而在后續(xù)迭代過程中，若某一步出現(xiàn)了 ‖uk＋1－uk‖＞‖uk－uk－1‖，將步長λk調整為λk＋1＝ckλk，ck稱為步長調節(jié)因子［14］，ck根據(jù)每次迭代過程中梯度和共軛梯度的情況進行取值，取值如下。

綜上所述，當k＞0時，λk的取值如下：

3.3 迭代過程

CGS方法迭代過程如下。

（1）將原始隨機變量x轉化為標準正態(tài)變量u，設定常數(shù)M，迭代停止準則ε＝0.001。

（2）令k＝0，設定初始迭代點u0，算例中沒有明確說明則取u0為坐標原點，根據(jù)式（14）計算初始步長λ0，根據(jù)式（9，10）計算u1，取λ1＝λ0，k＝k＋1，式（14）的常數(shù)M本文均取35。

（3）采用有限差分法計算功能函數(shù)在點uk處的梯度Ug（uk）。

（4）根據(jù)式（12）確定θk，根據(jù)式（10）確定Sk，結合式（9）計算uk＋1，根據(jù)式（15，16）調整步長。

（5）若‖uk＋1－uk‖＜ε，停止迭代，uk＋1即為u＊，g（uk＋1）為所求Gp；否則，令k＝k＋1，返回步驟（3）。

3.4 算 例

通過多個算例驗證CGS方法求解概率功能度量時的有效性和穩(wěn)健性，并與其他常見方法如AMV，HMV，CC和MCC進行比較，考證CGS方法的計算效率。參與對比的這些方法均在Matlab R2014a軟件中實現(xiàn)，均采用同一收斂準則‖uk＋1－uk‖＜0.001，靈敏度分析均采用有限差分法。

算例1 考慮如下凹功能函數(shù)［12］，

圖3 功能函數(shù)示意圖（算例1）Fig.3 Performance function diagram （example 1）

標準正態(tài)空間中，功能函數(shù)和最小功能目標點u＊如圖3所示，可以看出，在u＊附近的功能函數(shù)為凹函數(shù)。表1～表11中，Iters為迭代次數(shù)，NFE為功能函數(shù)計算次數(shù)。表1對比了分別采用和本文所提出的兩種共軛梯度系數(shù)的計算結果，初始迭代點分別取u0＝（0.0，0.0）和u0＝（3.0，3.0）?？梢钥闯?，采用兩種共軛搜索方向雖然均能得到收斂解，但采用時，效率明顯高于。表2給出了初始迭代點為u0＝（0.0，0.0）時，兩種不同共軛梯度系數(shù)的迭代歷程?？梢钥闯?，整個迭代過程中，共軛梯度系數(shù)的取值發(fā)生大幅波動，因而共軛搜索方向的波動也較大，影響計算效率；而的取值相對平穩(wěn)，避免了連續(xù)迭代過程中共軛搜索方向大幅波動的情況；表1結果也驗證了采用的高效性和穩(wěn)健性。

表3給出了初始迭代點為u0＝（0.0，0.0）時，CGS，AMV，HMV和CC等常見方法的計算結果。文獻［12］給出MCC方法計算結果Gp＝－21.672，但經(jīng)計算發(fā)現(xiàn)該答案錯誤?？梢钥闯?，當采用CGS方法時，計算效率很高，僅經(jīng)過12次迭代和37次函數(shù)調用就得到了收斂解，概率功能度量值為Gp＝－76.035。因此，針對非線性較高的凹函數(shù)，CGS方法具有良好的收斂性。

算例2 考慮凹功能函數(shù)［12］：

表4給出了各方法的計算結果，除AMV方法外，其余方法都能收斂，且與文獻［12］給出的求解結果Gp＝－2.2293一致?？梢钥闯?，本文提出的CGS方法不僅能夠得到收斂解，其計算效率也明顯高于其他方法，僅需11次迭代和34次功能函數(shù)調用，再次驗證了CGS方法求解凹函數(shù)的高效性。

CGS方法的迭代歷程列入表5?？梢钥闯觯ㄟ^步長準則自動選取的初始步長λ＝51.356，經(jīng)兩次調節(jié)減小到3.340。由于CGS采用了自適應步長調節(jié)策略，因此無需估計功能函數(shù)非線性程度，也無需假定λ的合適初值，通過步長準則會自動選取初值λ，并在迭代過程中不斷調整，直至最終收斂。

算例3 考慮凸功能函數(shù)［7］：G（X）＝－exp（X1－7）－X2＋10

Xi～N（6.0，0.82），i＝1，2，βt＝3.0

文獻［7］給出了求解結果，最小功能目標點u＊＝（2.8963，0.7813），概率功能度量Gp＝－0.358。從表6可以看出，各種方法的求解結果與文獻所得結果均非常吻合。當采用AMV方法時，經(jīng)過8次迭代和25次功能函數(shù)調用達到收斂，說明對于凸功能函數(shù)，AMV方法具有較好的收斂性能。同時，相對于其他方法，CGS方法的計算效率也較高，經(jīng)過11次迭代和34次功能函數(shù)調用得到收斂解，即CGS方法對于凸函數(shù)的計算效率與AMV基本相當。

表1 采用不同共軛梯度系數(shù)的計算結果Tab.1 Results of different conjugate gradient coefficients

表2 不同共軛梯度系數(shù)的迭代歷程Tab.2 Iterative history of different conjugate gradient coefficients

表3 算例1中各方法的計算結果Tab.3 Results of various methods in example 1

表4 算例2中各方法的計算結果Tab.4 Results of various methods in example 2

算例4 考慮如下功能函數(shù)［9］：

該功能函數(shù)隨機變量較多，有7個隨機變量。表7給出了不同方法的計算結果，除了AMV方法以外，其他方法均能收斂，概率功能度量Gp＝0.075，與文獻［9］的結果一致。CGS方法經(jīng)過13次迭代和105次函數(shù)調用，計算效率相對較高。

算例5 考慮如下高度非線性功能函數(shù)［11］：

表8給出了各常見方法的計算結果。采用初始迭代點u0＝（0.0，0.0）時，只有CGS方法收斂，且所得結果與文獻［11］結果吻合。若采用初始點u0＝（－2.0，－2.0），CGS方法僅經(jīng)過15次迭代就收斂到最小功能目標點，與CC方法相比，CGS方法更加高效。由此可看出，CGS方法對初始點的依賴較小，且能顯著改善高度非線性功能函數(shù)概率功能度量求解的收斂性能。

算例6 如圖4所示，三跨十二層框架結構［10］由同種材料制成，楊氏模量為20GPa，P為水平荷載。梁柱單元截面面積分為5組，圖中數(shù)字編號表示不同單元截面面積Ai（i＝1，2，…，5）的組別，截面慣性矩為Ii＝αiA2i（i＝1，2，…，5），系數(shù)αi的值列入表9。Ai（i＝1，2，…，5）與P是相互獨立的隨機變量，其統(tǒng)計信息列入表9。當結構最大位移超過0.096m時認為結構失效，目標可靠指標為βt＝3.0。

文獻［10］給出概率功能度量為Gp＝－0.0134，最小功能目標點u＊＝（－0.6923，－0.2559，－1.0827，－1.1523，0.0002，2.4403）。從表10可以看出，僅AMV方法和CGS方法能夠收斂，且計算結果與文獻結果相近，其中CGS方法調用函數(shù)次數(shù)較少，效率較高。

表5 算例2中CGS方法的迭代歷程Tab.5 Iterative history of CGS in example 2

表6 算例3中各方法的計算結果Tab.6 Results of various methods in example 3

表7 算例4中各方法的計算結果Tab.7 Results of various methods in example 4

表8 算例5中各方法的計算結果Tab.8 Results of various methods in example 5

算例7 圖5所示十桿桁架，材料楊氏模量為1.0E7psi，施加兩個垂直向下荷載P1＝P2＝1.0E5lb，各桿單元橫截面面積是相互獨立的正態(tài)分布隨機變量，均值為10in2，變異系數(shù)為0.05。結構功能函數(shù)定義為桁架的最大豎向位移不超過4in，目標可靠指標βt＝3.0。

圖4 十二層框架結構Fig.4 Twelve-story frame structure

圖5 十桿桁架Fig.5 Ten-bar truss

表9 算例6中隨機變量信息Tab.9 Information of random variables in example 6

表10 算例6中各方法的計算結果Tab.10 Results of various methods in example 6

表11 算例7中各方法的計算結果Tab.11 Results of various methods in example 7

從表11可以看出，僅AMV方法和CGS方法能夠收斂，CGS方法經(jīng)18次迭代收斂，其效率遠遠高于AMV方法，再次驗證了CGS方法的高效性。

4 結 論

針對AMV方法求解高度非線性功能函數(shù)的概率功能度量難以收斂的情況，本文提出了共軛梯度步長調節(jié)法（CGS），引入新的共軛梯度搜索方向和自適應步長調節(jié)策略，以求解概率功能度量。新的共軛搜索方向保證了收斂穩(wěn)定性；自適應步長調節(jié)策略無需了解功能函數(shù)凹凸性及非線性程度等先驗信息，無需尋找步長λ的合適取值。通過限定步長準則自動選取初始步長λ，并不斷調節(jié)λ，直至最終收斂。多個數(shù)值算例表明，相比于AMV和HMV等方法，本文提出的CGS方法不僅顯著地改善了高度非線性功能函數(shù)以及凹功能函數(shù)概率功能度量求解的收斂性能，同時對于凸功能函數(shù)也保持了較高的計算效率。