全局灵敏性分析方法的研究_谷良贤

全局灵敏性分析方法的研究
谷良贤,仇理宽,龚春林
(西北工业大学航天学院,陕西西安
710072)
摘 要:全局灵敏性分析是解决多学科优化中复杂耦合系统的重要方法,对全局灵敏性分析技术进行了论述,基于一阶模型,建立了二阶以及高阶导数模型,发展了全局敏感性计算方法。

编制了基于梯度的优化算法,并与全局灵敏度方程(GSE )求得的导数信息融合,进行了实例应用。

与多种优化方法进行了比较,证明了设计方法的可行性、优越性。

关键词:多学科优化;全局灵敏度分析;耦合系统
中图分类号:V 421.1 文献标识码:A 文章编号:1671 654X (2011)01 0001 04
Study on G lobal Sensitivity AnalysisM ethod
GU L iang x i a n ,Q I U L i kuan ,GONG Chun lin
(S chool of Aeronautics ,N orthw estern Poly technical University ,X i an 710072,Ch ina)
Abst ract :Syste m Sensiti v ity A na l y sis(SSA)is an i m portant m ethod to so l u te co m p lex ,interna ll y coupled syste m s inM ultidisciplinary design opti m ization .The study presents an overvie w of the SSA,deve loped the
solution of Syste m Sensiti v ity and constructed the Second O r der and H i g her O rder Sensiti v ity Derivatives M odel based on the F irst O r derM ode.l The research w or ks out several opti m ization m et h ods based on de rivatives .The application of the co m binati o n o f the GSE m ethod and severa l opti m izati o n m ethods are g i v en .The desi g n m et h od is feasi b le and advantageous by co m paring w ith o ther opti m izati o n m ethods .K ey w ords :mu lti d isc i p li n ar y desi g n opti m izati o n ;syste m sensiti v ity analysis ;coupled syste m s
引言
在复杂系统的多学科设计优化(M u lti d iscipli n ary Design Opti m izati o n ,MDO )中,由于系统各学科/子系统间耦合关系复杂,相互作用,导致了系统分析的计算困难和子系统间信息交换复杂等难题,如何解决学科
耦合问题成为MDO 亟待解决的难点之一。

现任A I A A MDO 技术委员会主席、美籍波兰裔数学家Sobiesk i 提出全局灵敏度方程GSE (G lobal Sensi ti v ity Equation )用于耦合系统的灵敏度分析
[1]
,GSE
利用每个子系统输出对其输入的偏导数构造整个系统的全局敏感度方程,求解该方程得到系统状态变量对设计变量的全导数,即可采用常规优化方法求解整个系统的优化设计问题。

GSE 将子系统与整个系统的灵敏度分析联系起来,与传统差分方法相比避免了全局差分,不需要对X 每个扰动求解全局导数,大大提高了计算效率和准确性,既体现了大系统中各子系统的
相互耦合作用,又实现了并行处理。

在此方法基础上,许多研究人员应用到工程设计中,取得了丰硕成果:H a j e la 等利用全局敏感度方程探索了通用飞机气动/性能/结构/控制一体化设计的方法[2]
;北航的方卫国、
郦正能等人基于全局灵敏度方程与约束变尺度法结合用于4座公务机方案优化的设计
[3];北航的姜欢等人
利用全局灵敏度方程求解多层次弹道优化设计[4]
;西工大的姜亚楠等人采用全局灵敏度方程计算准则函数及约束函数对设计变量的灵敏度,得到飞行器稳健设计最优解
[5]
;国防科大的王振国、颜力等人对适用于飞
行器MDO 的灵敏度分析方法进行了系统研究[6]。

但
是这些研究都是基于GSE 方法求得的一阶导数(灵敏度)信息加以利用,二阶及其高阶导数信息还少有应用,而许多算法需要二阶甚至高阶导数信息。

为此,本文建立了一阶、二阶以及高阶全局灵敏度计算模型,并结合基于梯度的优化算法进行了验证,证明了所建立一阶、二阶全局灵敏度方程的合理性,正确性。

对于
收稿日期:2010 09 01
基金项目:国家自然科学基金资助项目(20092012082900)
作者简介:谷良贤(1957-),女,河南沁阳人,教授,博士生导师,研究方向为飞行器总体设计、多学科优化设计。

第41卷 第1期航空计算技术
Vo.l 41N o .12011年1月
A eronauti ca l Co m puti ng T echn i que
Jan .2011
GSE方法的广泛推广具有一定意义。

1 灵敏度计算
对于耦合的多学科优化问题,其数学描述为:
y1=a(X,y2)
y2=b(X,y1)
(1)其中学科(子系统)1、2间存在相互影响,必须通过迭代求解才能解耦。

求解y1,y2的一阶全导数:
d y1 d X =
y1
y2
d y2
d X
+
y1
X
d y2 d X =
y2
y1
d y1
d X
+
y2
X
(2)
整理(2)得:
1- y1 y2
- y2
y11d y1
d X
d y2
d X
=
y1
X
y2
X
(3)
方程(3)即为全局灵敏性方程(GSE),可以推广到多维形式。

设X=[x1,x2, ,x k]为设计变量,Y=[Y1,Y2, ,Y m],Y i表示第i(1 i m)个子系统的输出变量,则全局灵敏性方程为:
I - Y1
Y2
- Y1
Y m
- Y2
Y1I - Y2 Y m
- Y m Y1- Y m
Y2 I
d Y1
d X
d Y2
d X
d Y m
d X
=
Y1
X
Y2
X
Y m
X
(4)
式(4)可简写为AB=C
其中I为单位矩阵,C称为局部灵敏性导数(1ocal
sensiti v ity derivatives,LSD),它表示不考虑其他影响,
各子系统的输出响应对输入变量的偏导数;A系数矩
阵称为全局灵敏性矩阵,它代表了每个子系统对其他
子系统输出的灵敏性;B称为全局灵敏性导数,表示考
虑子系统间的耦合关系,各子系统输出响应对任意输
入变量的导数信息。

全局灵敏性导数B可以通过A-1和C求得,即:
d Y1
d X
d Y2
d X
d Y m
d X
=
I
- Y1
Y2
- Y1
Y m
- Y2
Y1I
- Y2
Y m
- Y m
Y1
- Y m
Y2
I
-1 Y
1
X
Y2
X
Y m
X
(5)
对式(1)推导其二阶导数:
y 1=
y1
y2y 2+(
y1
y2) y
2+(
y1
X
) =
y1
y2y 2+
2y1
y2 y2y
2y 2+
2y1
y2
X
y 2+
2y1
X y2y 2+
2y1
X X
y 2=
y2
y1y 1+(
y2
y1)
y 1+(
y2
X) =
y2
y1
y 2+
2y2
y1 y1y 1y 1+
2y2
y1 X y 1+
2y2
X y1y 2+
2y2
X X
(6)
整理式(6)得:
1
- y1
y2
-
y2
y11
y 1
y 2
=
2y1
y2 y2
2y2
y1 y10
y 1y 1
y 2y 2
+
2y1
y2 X+
2y1
X y2
2y2
y1 X+
2y2
X y10
y 1
y 2
+
2y1
X X
2y2
X X
(7)上式同样可以推广到多维形式,其二阶全局灵敏性方程为:
I
- Y1
Y2
- Y1
Y m
- Y2
Y1I
- Y2
Y m
- Y m
Y1
- Y m
Y2
I
Y 1
Y 2
Y m
=
2
航空计算技术 第41卷第1期
m
i=1,j=1 2
Y 1
Y i Y j
Y i Y
j m
i=1,j=1 2
Y 2
Y i Y j Y i Y j
m
i=1,j=1 2
Y m
Y i Y j
Y i Y j +0 2Y 1 Y 2 X + 2
Y 1
X Y 2
2
Y 1 Y m X +
2
Y 1
X Y m
2
Y 2 Y 1 X + 2
Y 2
X Y 1
0 2
Y 2 Y m X + 2
Y 2
X Y m
2
Y m Y 1 X + 2Y m
X Y 1
2
Y m Y 2 X + 2
Y m
X Y 2
Y 1Y 2 Y m
+
2
Y 1
X
X 2
Y 2 X X 2Y m X X
计算二阶导数所需要的一阶导数信息可以通过式(5)得到。

高阶导数信息可以类似推导得出:1 Y 1 Y 2
- Y 1
Y m - Y 2
Y 1I
- Y 2
Y m - Y m Y 1
- Y m Y 2
I
Y
(n )
1Y (n )2
Y
(n )
m
=
2
Y 1
Y 2 X 2
Y 1 Y m X 2
Y 2 Y 1 X 0 2Y 2 Y m X 2
Y m
Y 1 X
2
2
Y m Y 2
X
Y (n-1)
1Y (n-1)
2
Y
(n-1)m
m i=1,j=1 2
Y 1
Y i Y j Y (n-1)
i Y
j m
i=1,j=1 2
Y 2 Y i Y j Y (n-1)
i Y j
m
i=1,j=1 2
Y m
Y i Y j Y (n-1)
i Y j + (n)
Y 1
X X (n)
Y 2 X X (n )Y m X X
+L L 依据不同的函数有不同的表达式。

2 灵敏度计算在优化设计中的应用
图1 基于G SE 优化方法流程图
基于GSE 优化方法的流程如图1。

设函数f (X ),在X (k)
点泰勒展开有:
f (X ) f (X
(k )
)+f 1(X
(k)
) X +
12f 2(X (k )) X 2+ +1n!
f n (X (k)) X n
取二阶项,以牛顿法为例,令f (X (k)
)=0,构造迭
代式如下:
X
(k +1)
=X
(k )
f 2(X
(k)
)
-1f 1
(X
(k)
)
上式中目标函数在X (k)
处的一阶导数f 1
(X
(k)
)、
二阶导数f 2
(X
(k)
)信息可通过前文推导的GSE 方法
一阶、二阶模型求得。

设优化的目标函数是Y i (1
k ),在灵敏度方向寻优,得到此灵敏度方向最优点然后在此最优点X 0用以上方法计算出新的全局灵敏度,并在新的灵敏度方向进行寻优;如此反复直到灵敏度趋近0时停止,即可求出系统全局最优解。

3 算例
下面用两个算例对本文建立的方法进行验证。

算例1为一简单线性耦合问题,首先用GSE 方法与差分方法分别计算灵敏度,验证一阶、二阶模型的正确性,然后结合牛顿法进行优化求解,并与单利用GSE 一阶导数信息的变尺度算法优化的结果进行了比较。

算例2是一非线性耦合问题,通过用GSE 方法计算一阶、二阶导数信息,并结合优化算法进行优化,证明了其在多学科优化中的实用性。

算例1:M i n i m ize
F 1=X 3
1-4X 2
2-X 1X 2+Y 1+Y 2+20Sub ject to
Y 1=X 2
1-4X 2+Y 2+20
Y 2=-10X 1+X 2
2-Y 1+40
其中X =[X 1,X 2]代表输入变量,Y =[Y 1,Y 2]代表各子系统的输出变量。

取初始点(0,0),用GSE 方法得到初始点处的一阶灵敏度为
-10-4,二阶灵敏度为2-1
-12
;用差分方法得到一阶灵敏度-9.99
-3.99
,二阶灵敏度
,说明了一阶、二阶导
3 2011年1月谷良贤等:全局灵敏性分析方法的研究
数计算模型的合理性、正确性。

用GSE 方法求出约束函数和目标函数的一阶、二阶导数信息,融合到牛顿法中,分别与牛顿法利用差分方法求得的导数信息得到的结果和只用GSE 一阶信息的变尺度法得到的优化结果比较如表1。

表1 三种基于梯度优化算法的计算结果
方法最小值点最小值迭代次数耗时/s 差分牛顿法(8.00478,6.00318)8.000031.07G SE 牛顿法(8.6)820.61GSE 变尺度法(7.9951,5.9604)8.0014
10
0.82
由表中数据可以看出,GSE 牛顿法与差分牛顿法相比,耗时少,精度高;与只利用GSE 一阶导数信息的变尺度法相比,精度更高,迭代次数更少。

当运用到耦合严重系统复杂的飞行器设计优化时,GSE 方法耗时
少精度高的优点将更加突出。

算例2:
M i n f =X 2
2+X 3+Y 1+e -Y 2
学科1分析
Y 1=X 2
1+X 2+X 3-0.2Y 2学科2分析Y 2=
Y 1+X 1
+X 2
Subject to g 1=
Y 1
8
-1 0g 2=10-Y 2
10 0
-10 X 1 100 X 2 100 X 3 10
其中X =[X 1,X 2,X 3]为设计变量,Y =[Y 1,Y 2]为状态变量。

此算例为一非线性耦合问题,且含有不等式约束条件。

为此,采用罚函数法将约束问题转化成无约束问题,构造新的目标函数为:
F (X ,r )=f (X )+r
m i=1
1
g i (X )
其中:r m
i=1
1/g i (X )为障碍项,r 为障碍因子,f (X )为
原目标函数,g i (X )为不等式约束。

取初始点为[3,0.1,0.1],障碍因子初始取0.1,
缩小率 取0.1。

利用GSE 方法求得原目标函数f (X )在初始点处的一阶导数为:
5.6509
0.77141.9654二阶导数为:2.01140.0160
0.0128
0.01602.00500.00220.0128
0.0022
0.0021
现目标函数F (X )在初始点处障碍因子r =0.1时
的一阶导数为:
103
-2.9250
-0.4133-0.5133
二阶导数为:
106 2.37390.32760.40950.32760.04540.0565
0.40950.0565
0.0708
融合到牛顿法中求得优化结果与理论结果比较见
表2。

图2 迭代次数与目标函数的变化
表2 基于G SE 的牛顿法优化结果与理论值比较方法最优值点f
理论值
(3.0283,0.0012,0)
8.0029
基于GSE 的牛顿法(3.0279,0.0016,0.0017)8.0046
由表2计算结果可以看出,GSE 融合牛顿法与理论结果相比基本相同,证明了利用GSE 方法求解一阶、二阶导数信息然后融合基于梯度的优化方法在多学科优化中的实用价值。

4 结论
对于复杂耦合系统的多学科优化问题,全局灵敏度分析方法具有很高的实用价值。

由于各学科间的高
(下转第13页)
4 航空计算技术 第41卷第1期
ing,1995.
[3] L ieu B i nh K,M oa rre,f R ashad,Jovanov ic,M i ha ilo R.P reven
ti ng transiti on to t urbulence usi ng strea mw ise trave li ng w aves:d irect nu m er ica l si m u l a ti ons[Z].Am erican Control Con ference(ACC),J une2010.
[4] Spa lart P R.T rends i n turbu lence treat m ents[R].A I AA-00
-2306,2000.
[5] S iH a i qing,W ang T ongguang.A Fu lly I mp licit So l ve r of3-
D Euler Equations on M ultiblock Curv ilinear G r i ds[Z].
M odern Physics L etters B.2005,19(N os.28&29):1483-
1486,38(5):595-599.
[6] 司海青,王同光.基于全隐式无分裂算法求解三维N S
方程[J].计算力学学报,2009,26(2):252-257.
[7] 常书平,王永生,庞之洋.用基于模型的方法数值模拟圆
柱绕流[J].舰船科学技术,2009,31(2).
(上接第4页)
度耦合关系,传统的方法不能取得很好的效果,全局灵敏度分析方法避免了差分方法带来的效率低、准确度低的弊端,而且系统越复杂、耦合越严重,GSE方法的优势越明显。

本文基于一阶模型,建立了二阶以及高阶模型,通过算例验证了一阶、二阶模型的正确性、优越性。

并通过一非线性耦合算例验证了利用GSE方法计算一阶、二阶导数信息结合基于梯度的优化算法在多学科优化中的合理性、实用性。

参考文献:
[1] Sob i eszczansk.i Sob ies k i J Sensiti v ity o f co m plex,interna l y
coupled system s[J].A I AA Journa,l1990,28(1).
[2] H aje l a P.A ppli cation of G loba l Sensiti v ity Equa ti ons i n M u l
tidisc i p li nary A ircraft Synthesis[J].Journa l o f A ircraft,
1992,27(12):1002-1010.
[3] 方卫国,郦正能.基于全局敏度方程方法的飞机方案设
计[J].航空学报,1998,19(3).
[4] 姜欢,周浩,陈万春,等.基于全局灵敏度方程的多层次弹
道优化设计[J].北京航空航天大学学报,2009,35(12).
[5] 姜亚楠,王和平.基于灵敏度分析的飞行器稳健设计优化
方法[J].航空计算技术,2009,39(5).
[6] 颜力,陈小前,王振国.飞行器多学科设计优化中的灵敏
度分析方法研究[J].航空计算技术,2005,35(1).
(上接第8页)
2)通过对哈尔滨终端区的容量评估实例,验证了本文模型和方法的正确性,评估结果得到了管制员的认可;
3)本文通过对实例结果的分析,提出了本方法的应用范围。

参考文献:
[1] 国际民用航空组织.9426-AN/924-1984,空中交通服务
计划手册[S].蒙特利尔:国际民用航空组织,1999.
[2] 陈亚青,孙宏.进近管制员工作进程分类及工作负荷研究
[J].中国安全科学学报,2006(2):65-68.
[3] R e i d G B,N ygren T E.T he subjecti ve workload assess ment
techn i que:A sca li ng pro cedure fo r m easur i ng m enta l wo rk
l oad[M].N orth-H o lland:E lsev i er Sc i ence Pub i shers,
1988:185-218.
[4] H art S G,Stave land L E.D ev elopment of NA S A-TLX
(T as k Load Index):R esults of emp irical and theo re ti ca l re
search.In Human M enta lW ork l oad[M].N o rt h-H o lland:
Elsev ier Sc i ence Pub ishe rs,1988:139-184.
[5] 万莉莉,胡明华.管制员工作负荷及扇区容量评估问题研
究[J].交通运输工程与信息学报,2006,4(2):70-75. [6] 刘欢.区域动态容量评估技术研究[D].南京:南京航空
航天大学民航学院,2009.
[7] Berts m as D,Patterson S S.T he a ir traffic fl ow m anage m ent
prob le m w ith enroute capac iti es[J].O pera ti ons R esearch,
1994,42(2):249-261.
13
2011年1月陈 都等:脱体涡模拟方法研究及其应用。