孫建勛 潘 雷 谷良賢

西北工業(yè)大學(xué)，西安，710072

0 引言

敏感性分析是多學(xué)科設(shè)計(jì)優(yōu)化的重要工具，可分為一維函數(shù)和多維函數(shù)的敏感性分析。目前，常用的一維敏感性分析方法有差商法、復(fù)步長(zhǎng)法、逆波蘭表示法、符號(hào)分析法和自動(dòng)微分法等5種[1]，而多維敏感性分析則有敏感性矩陣法、滯后耦合伴隨法等[2]。這些方法在進(jìn)行函數(shù)敏感性計(jì)算時(shí)都需要進(jìn)行大量函數(shù)值計(jì)算，從而使得敏感性分析的代價(jià)極其巨大。研究表明，在一次優(yōu)化循環(huán)中，有70%以上的計(jì)算費(fèi)用用于函數(shù)敏感性分析[3]。因而，降低敏感性分析費(fèi)用，尋找快速有效的敏感性分析方法一直是敏感性分析的一個(gè)重要研究方向[4]。

響應(yīng)面法是一種常用的近似技術(shù)，它通過(guò)對(duì)函數(shù)的一系列已知點(diǎn)進(jìn)行擬合從而推測(cè)函數(shù)的變化趨勢(shì)[5]。本文結(jié)合響應(yīng)面法特點(diǎn)，提出一種基于響應(yīng)面法的函數(shù)敏感性分析方法。該方法計(jì)算簡(jiǎn)單，適于快速敏感性分析。在優(yōu)化過(guò)程中引入該方法能夠減少目標(biāo)函數(shù)計(jì)算次數(shù)，加快優(yōu)化收斂過(guò)程，減少優(yōu)化費(fèi)用。

1 多項(xiàng)式響應(yīng)面與函數(shù)敏感性的關(guān)系

基于多項(xiàng)式的響應(yīng)面法通常采用高階多項(xiàng)式，其中最常用的是二階響應(yīng)面，也稱二次響應(yīng)面。用η表示二次響應(yīng)面模型的響應(yīng)量，則具有n個(gè)變量的二次響應(yīng)面函數(shù)為

當(dāng)給定一組已知數(shù)據(jù)時(shí)，則可利用最小二乘法對(duì)式（1）進(jìn)行求解。

將待擬合函數(shù)y＝f（X）在X0點(diǎn)展開，得

其中，▽f（X0）為函數(shù)f在X0點(diǎn)的梯度，即敏感性；Hf（X0）為函數(shù)f在X0點(diǎn)的Hessian矩陣；o（‖X‖3）表示余項(xiàng)為三階無(wú)窮小量。對(duì)比式（1）與式（2）可知，多項(xiàng)式響應(yīng)面（response surface method，RSM）近似的主要思想仍為Taylor展開的思想，響應(yīng)面的系數(shù)從某種程度上近似反映了原函數(shù)的一階梯度與二階梯度。

實(shí)際問(wèn)題中，函數(shù)往往以隱函數(shù)的形式表示，因而有必要研究隱函數(shù)情況下的敏感性與響應(yīng)面系數(shù)的關(guān)系?？紤]最一般的情況，狀態(tài)變量u為該學(xué)科設(shè)計(jì)變量Xd、耦合設(shè)計(jì)變量Xs以及耦合狀態(tài)變量Ys的函數(shù)，即u＝f（Xd，Xs，Ys），對(duì)應(yīng)的響應(yīng)面可描述為

式中，β為待定系數(shù)。

對(duì)式（3）求全微分：

則在任一點(diǎn)（Xd，Xs，Ys）處，根據(jù)全微分關(guān)系du＝▽fdx可知

式（5）即為響應(yīng)面法計(jì)算任一點(diǎn)梯度的近似計(jì)算公式。

需要注意的是，在進(jìn)行系統(tǒng)全局敏感性計(jì)算時(shí)，耦合狀態(tài)變量Ys的微分dYTs需通過(guò)各個(gè)學(xué)科的響應(yīng)面進(jìn)行求解。分別按上述方法對(duì)每個(gè)子系統(tǒng)進(jìn)行處理，并假定一個(gè)系統(tǒng)只含有一個(gè)狀態(tài)變量，則n個(gè)子系統(tǒng)中共含有n個(gè)dYs，構(gòu)成一個(gè)n元線性方程組。解此方程組則可求得各學(xué)科的敏感性。

顯然，這種基于響應(yīng)面法的敏感性分析方法的優(yōu)點(diǎn)在于計(jì)算簡(jiǎn)單，只需要根據(jù)響應(yīng)面的系數(shù)求解若干個(gè)由式（5）組成的線性方程組，避免了大量調(diào)用原函數(shù)，尤其是對(duì)于存在狀態(tài)變量耦合的系統(tǒng)，可以避免狀態(tài)變量的迭代求解過(guò)程。因而該方法簡(jiǎn)單快捷，可用于快速敏感性分析。

2 敏感性近似能力分析

采用一個(gè)三元顯函數(shù)算例（式（6））來(lái)測(cè)試基于RSM的函數(shù)敏感性計(jì)算能力：

圖1為部分參數(shù)固定情況下的函數(shù)圖像，由圖1可知，該函數(shù)的敏感性變化范圍比較大，變化情況多樣，適于敏感性計(jì)算性能測(cè)試。

圖1 x2、x3 固定時(shí)測(cè)試函數(shù)y（x1）圖像

采用二階分塊響應(yīng)面近似方法[6]，首先檢驗(yàn)該RSM曲面擬合精度。在設(shè)計(jì)空間中按照均勻試驗(yàn)設(shè)計(jì)方法選取45個(gè)試驗(yàn)點(diǎn)，其中x1維5個(gè)測(cè)試點(diǎn)，其他兩維都是3個(gè)測(cè)試點(diǎn)。圖2為測(cè)試結(jié)果，其中，大部分計(jì)算結(jié)果都在解析解附近，近似結(jié)果與解析結(jié)果的最大誤差不超過(guò)5%。因而，采用分塊RSM有效地?cái)M合了原設(shè)計(jì)空間，可以用于梯度近似計(jì)算。

圖2 RSM近似性能圖像

在整個(gè)設(shè)計(jì)空間內(nèi)按照均勻試驗(yàn)設(shè)計(jì)的方法選取112個(gè)測(cè)試點(diǎn)，其中x1維7個(gè)測(cè)試點(diǎn)，其他兩維都是4個(gè)測(cè)試點(diǎn)。圖3給出了RSM在測(cè)試點(diǎn)處計(jì)算得到的近似敏感性值和解析值的比值。表1給出了圖3部分測(cè)試點(diǎn)的解析結(jié)果、近似結(jié)果和計(jì)算誤差。

表1 圖3中的典型點(diǎn)敏感性結(jié)果

圖3 RSM梯度近似結(jié)果和解析值比值

圖3表明，RSM方法近似梯度計(jì)算結(jié)果整體性能較好，但存在部分點(diǎn)近似結(jié)果誤差較大。結(jié)合表1可以看出，當(dāng)解析梯度數(shù)值較大時(shí)，RSM方法近似梯度計(jì)算得到的結(jié)果和解析結(jié)果誤差不大，如點(diǎn)35、點(diǎn)89處的結(jié)果。而當(dāng)解析結(jié)果接近0時(shí)，即當(dāng)原函數(shù)接近其極點(diǎn)時(shí)，RSM方法計(jì)算結(jié)果和原函數(shù)的解析結(jié)果會(huì)有較大的相對(duì)誤差，甚至出現(xiàn)符號(hào)反向的情況，但此時(shí)絕對(duì)誤差卻很小，如點(diǎn)67處的x1維梯度值，解析結(jié)果為0.01，近似結(jié)果為－0.009，兩者相對(duì)誤差高達(dá)190%，而絕對(duì)誤差小于0.02，不及該點(diǎn)其他維梯度絕對(duì)誤差的1%。同樣，點(diǎn)73處x1維梯度值的解析結(jié)果和近似結(jié)果有著較大的相對(duì)誤差，但兩者的絕對(duì)誤差則很小，不足該點(diǎn)處其他維梯度絕對(duì)誤差的0.5%。與其他維敏感性相比，該維的計(jì)算誤差可以忽略。

圖4給出了在x2、x3維固定時(shí)，RSM方法沿x1維對(duì)3個(gè)自變量的梯度計(jì)算結(jié)果與解析結(jié)果的比值變化圖。x2、x3維的近似梯度計(jì)算結(jié)果比較滿意，除x3維個(gè)別點(diǎn)誤差在15%附近外，其余結(jié)果誤差都在10%以內(nèi)。對(duì)該結(jié)果進(jìn)行分析可得出與圖3同樣的結(jié)論。

圖4 x1維近似梯度和解析梯度的比值

3 基于響應(yīng)面近似的函數(shù)優(yōu)化

敏感性分析的一個(gè)重要用途是用于優(yōu)化設(shè)計(jì)，下面用一個(gè)實(shí)例來(lái)檢驗(yàn)基于響應(yīng)面近似的敏感性快速分析策略的性能。該算例用于計(jì)算在兩沖量作用下共面橢圓軌道轉(zhuǎn)移問(wèn)題，如圖5所示。選取優(yōu)化變量為轉(zhuǎn)移軌道的近地點(diǎn)幅角ω、半長(zhǎng)軸a和偏心率e，優(yōu)化目標(biāo)為最小化軌道轉(zhuǎn)移的能量需求，模型表述為

式中，下標(biāo)0表示初始軌道參數(shù)；下標(biāo)f表示終端軌道參數(shù)。

圖5 兩沖量共面橢圓軌道轉(zhuǎn)移

在上述條件下每次軌道轉(zhuǎn)移需要的沖量大小計(jì)算過(guò)程見(jiàn)文獻(xiàn)[7]。

對(duì)問(wèn)題歸一化處理，設(shè)初始軌道與目標(biāo)軌道的半長(zhǎng)軸長(zhǎng)度均為0.5、近地點(diǎn)幅角分別為0、－π／6，偏心率均為0.1。選用序列二次規(guī)劃法進(jìn)行優(yōu)化計(jì)算，并分別采用有限差分法和基于RSM的敏感性快速計(jì)算法獲得優(yōu)化過(guò)程需要的梯度信息，當(dāng)函數(shù)敏感性小于0.1時(shí)采用有限差分法獲得敏感性。優(yōu)化得到的最優(yōu)轉(zhuǎn)移軌道參數(shù)為（a，w，e）＝ （5.438，0.499，0.001），此時(shí)需要的轉(zhuǎn)移能量大小為0.176。

圖6所示為兩種方法的收斂過(guò)程，該圖表明，兩種方法中目標(biāo)函數(shù)的收斂過(guò)程很接近，這表明基于RSM的敏感性分析法可有效用于優(yōu)化過(guò)程。

圖6 優(yōu)化收斂過(guò)程

由圖6b可知，基于RSM的敏感性近似優(yōu)化過(guò)程的目標(biāo)函數(shù)計(jì)算次數(shù)要少于基于有限差分優(yōu)化過(guò)程的目標(biāo)函數(shù)計(jì)算次數(shù)。主要原因是，響應(yīng)面法在迭代開始大量計(jì)算目標(biāo)函數(shù)，生成響應(yīng)面。在迭代后期采用快速敏感性分析法，不需要計(jì)算目標(biāo)函數(shù)。有限差分法則自始至終都需要進(jìn)行目標(biāo)函數(shù)計(jì)算?；赗SM近似的優(yōu)化過(guò)程總共循環(huán)了39次，目標(biāo)函數(shù)計(jì)算次數(shù)為237；基于有限差分的優(yōu)化過(guò)程總共循環(huán)了44次，目標(biāo)函數(shù)計(jì)算次數(shù)為244?？偟膩?lái)看，采用RSM敏感性近似計(jì)算的優(yōu)化過(guò)程目標(biāo)函數(shù)計(jì)算次數(shù)要少。

4 結(jié)論

（1）RSM在大部分情況下能夠有效近似目標(biāo)曲面的敏感性。但在目標(biāo)曲面極點(diǎn)附近，即目標(biāo)曲面的敏感性較小時(shí)，RSM敏感性近似結(jié)果較差。

（2）基于RSM的敏感性計(jì)算方法能夠有效用于優(yōu)化過(guò)程，減少優(yōu)化過(guò)程的目標(biāo)函數(shù)計(jì)算次數(shù)，加快收斂過(guò)程，并得到理想結(jié)果。

（3）基于RSM方法的函數(shù)敏感性分析計(jì)算簡(jiǎn)單，在非極點(diǎn)附近計(jì)算精度較高，是一個(gè)良好的快速敏感性分析工具。

[1]蔡毅，邢巖，胡丹.敏感性分析綜述[J].北京師范大學(xué)學(xué)報(bào)（自然科學(xué)版），2008，44（1）：9-16.

[2]鐘毅芳，陳柏鴻，王周宏.多學(xué)科綜合優(yōu)化設(shè)計(jì)原理與方法[M].武漢：華中科技大學(xué)出版社，2006.

[3]Edmond P，Pelletier D，Etienne S.A Sensitivity E－quation Method for Compressible Subsonic Laminar Airfoil Flows[C]／／23rd AIAA Applied Aerodynamics Conference.Ontario，Canada，2005：AIAA 2005－4601.

[4]Martins J，Kroo I，Alonso J.An Automated Method for Sensitivity Analysis Using Complex Variables[C]／／38th Aerospace Sciences Meeting and Exhibit.Reno，NV，2000：AIAA－2000－0689.

[5]余大勝，康海貴.響應(yīng)面法計(jì)算結(jié)構(gòu)可靠度的回顧[J].工業(yè)建筑，2006，36：238－242.

[6]潘雷，谷良賢.分塊響應(yīng)面法研究[J].計(jì)算機(jī)工程與應(yīng)用，2009，45（19）：37－40.

[7]楊嘉墀.航天器軌道動(dòng)力學(xué)與控制[M].北京：宇航出版社，2001.