多输出性能的模拟是飞机结构系统中常见的问题。一方面我们关心的输出性能通常有多个，如研究复合材料飞机结构往往要同时关注最大位移和临界强度等多个输出性能^[1]；另一方面，结构系统的输出性能模型可能随时间和空间发生变化^[2]，如结构系统的动态性能等。传统的重要性测度分析对每个输出进行产生多个重要性测度结果，很难从整体上对工程设计给出实用的指导。针对结构系统多输出情况下的重要性测度分析，Gamboa等人基于协方差分解定义了一种新的广义重要性测度指标^[3]。而Campbell等人基于输出分解方法提出在一个合适的函数基底上展开原模型再进行重要性测度分析的方法^[4]。这2个指标都由传统的基于方差的Sobol指标^[5]发展而来，而Ruan等人^[6]提出了一种改进的方差重要性测度，该改进的指标包含输入变量对输出响应量方差影响的更全面信息。本文将在文献[6]基础上利用Gamboa等人提出的协方差分解方法，定义了结构系统多输出情况下重要性测度新指标。Gamboa等人给出了其所提指标求解的Monte Carlo算法，但其计算代价十分巨大，不适合大多数工程实际问题。Garcia-Cabrejo等人^[7]总结了Gamboa和Campbell所提的2种指标，并把多项式混沌展开(PCE)方法引入到多输出情况下重要性测度的求解，相对于文献[3]中的Monte Carlo数字模拟法和文献[4]中的主成分分析方法，PCE极大地减少了模型运行次数。对于单个输出的重要性测度分析，文献[8]提到Rosenblueth建议使用一种降维的代理模型，从而极大地降低了运算量。Zhang等人在此基础上提出乘法降维模型(M-DRM),相对于PCE方法进一步减少了模型调用次数。本文将M-DRM方法扩展到结构系统多输出情况下重要性测度新指标的求解。

Sobol等人^[5]利用响应量方差的ANOVA(analysis of variance)分解给出了基于方差的重要性测度指标。根据ANOVA分解，功能函数Y=g(x)可以被分解为一系列子项之和,即

点击浏览公式

式中，x=｛x₁,x₂,…,x_n｝为n维输入变量,g₀表示Y=g(x)的均值,g_i(x_i)表示由于x_i变化引起的Y=g(x)的变化,而g_i1i2(x_i1,x_i2)给出了x_i1,x_i2相互作用对输出的影响。在基本输入变量相互独立的情况下,响应量的方差可以由各方差分解项之和表示,即

点击浏览公式

其中各方差分解项为

点击浏览公式

公式(2)表明输出响应量的方差可以分解为各个输入变量以及它们相互作用的方差项,基于方差的重要性测度指标S_i1,…,is则定义为方差分解项与响应量方差之比,即

点击浏览公式

公式(4)中对应于单个输入变量的重要性测度指标为

点击浏览公式

文献[10]中给出了各方差项的求解公式为

点击浏览公式

式中,E(·)为期望算子,V(·)为方差算子。

基于方差的重要性测度只适用于单个输出情况。针对多输出情况,Gamboa等人^[3]基于协方差分解方法定义了一种新的广义重要性测度指标。首先考虑Y₁=f(X₁,…,X_n)和Y₂=g(X₁,…,X_n)2个单模式输出。定义u和v是输入变量(X₁,…,X_n)的2个互补子集,根据ANOVA分解,2个输出响应可以分解为

点击浏览公式

进一步推导如下

点击浏览公式

点击浏览公式

公式(9)的分解表明Y₁和Y₂的协方差可以分解为:子集u中输入变量变化引起的协方差,子集v中输入变量变化引起的协方差和子集u与子集v相互作用的影响项。进一步扩展到m个输出模式,则有

点击浏览公式

为得到标量形式的指标,Gamboa在文献[3]中建议对公式左乘一个M矩阵再求取矩阵的迹,并说明M取为单位矩阵,即

点击浏览公式

进而得出多输出模式下单个输入变量的重要性测度指标

点击浏览公式

以条件矩的形式表示(12)式指标,可得

点击浏览公式

由于不同输出的响应量量纲的差异,公式(13)对于方差的直接加和运算会抹除某些输出模式的方差效应,造成集成的指标失去部分输出模式的有效信息,从而使得所定义的广义重要性测度指标失去了对工程设计的指导意义。为此,我们将使用无量纲化的输出模型进行多输出情况下的重要性测度分析,即

点击浏览公式

这样,Gamboa定义的广义重要性测度指标转化为

点击浏览公式

公式(15)是基于传统的Sobol指标定义的,而Ruan等人在文献[6]中提出可以变换公式(6)中V_i的求解公式可以得到,

点击浏览公式

从(16)式的推导可以看出:当输入变量X_i取不同实现值时,V_i实际上定义为V(Y)与V(Y|X_i)差的平均值,反映了其对响应量无条件方差的平均影响。但这种期望运算可能使V(Y)-V(Y|X_i)在X_i的随机取值范围内的正负值相互抵消而掩盖了输入变量X_i的V_i对于响应量方差的较大影响。为避免上述指标求解过程中存在的正负抵消问题,可以引入平方运算,即:。它保证了在期望运算后X_i的每个实现值对响应量方差的影响被完整的保留下来而不会由于正负值相互抵消而掩盖某些影响。对改进的指标计算式进一步推到如下：

点击浏览公式

由(17)式可以看出,改进的指标S_i^m 实际上是对原指标的一种修正,即在S_i的基础上考虑了X_i不同实现值时对响应量方差影响的变异性V[V(Y|X_i)],更加完整地反映X_i取不同实现值对响应量所造成的影响。在结构系统多输出情况下,通过X_i取不同实现值对各个输出情况下响应量方差所造成的影响的加和来权衡X_i对系统的重要性程度,此时更加完整地反映X_i对每个输出情况下响应量方差所造成的影响是十分必要且有意义的。

结合上述观点,我们可以定义结构系统多输出情况下重要性测度新指标:

点击浏览公式

用于结构系统多输出情况下重要性测度指标的求解方法很多,各有优缺点。如文献[3]中的Monte Carlo数字模拟法,文献[6]中的SDP方法,文献[7]中的PCE方法等等。其中Monte Carlo数字模拟法最为直观,且易于编程实现。Monte Carlo数字模拟法的收敛解可作为精确解,但其需要采用双层抽样计算,效率低,在本文中将Monte Carlo法的求解结果解作为参考对照解。由于乘法降维模型的高效性,本文将在此模型的基础上建立新指标的求解方法,具体过程如下。

点击浏览公式

式中,h₀表示所有基本输入变量取均值时功能函数的响应值,h(x_i,μ_-i)表示除输入变量x_i外其他基本输入变量都取均值时功能函数的响应值。文献[9]中给出了乘法降维模型应用到传统单模式Sobol′重要性测度的详细过程。本文将这种替代模型应用到新指标的求解中。为后续计算推导表示方便,引入如下计算式,

点击浏览公式

式中,h(X_i)=h(x_i,μ_-i)表示公式(19)降维替代模型中各输入变量一维计算模型,而ρ、θ、ω分别表示各输入变量一维计算模型的一阶矩、二阶矩和四阶矩。ρ、θ、ω的计算可采用Gauss型求积公式^[11],同时针对不同的输入分布采用不同类型的Gauss型求积公式。如Gauss-Chebyshev求积公式可以用于均匀分布,Gauss-Hermite求积公式可以用于正态分布,Gauss-Laguerre求积公式可以用于指数分布等。(n+1)个节点的Gauss型求积公式的代数精度可以达到(2n+1),本文的2个算例均采用5节点的Gauss型求积公式,在减少调用模型次数同时保证求解精度。

为简明推导过程,先令Y=h(x₁,x₂,…,x_n)表示任意一个单模式输出响应,则替代模型化为h(x)=带入新指标(18)各项推导如下,

响应量方差的求解:

点击浏览公式

条件方差的求解:

点击浏览公式

这样将E(Y|X_i)看成X_i的单一变量函数,求解其响应量的方差,

点击浏览公式

V[V(Y|X_i)]的求解:先求解V(Y|X_i),

点击浏览公式

这样将V(Y|X_i)看成X_i的单一变量函数,求解其响应量的方差,

点击浏览公式

结构系统多输出情况下重要性测度新指标的求解:

点击浏览公式

从公式(26)可以看出,使用乘法降维替代模型求解新指标只需要m×[n×l-(n-1)]次功能函数运算就可以获得很精确的结果(其中m表示输出模式个数,n表示输入变量个数,l表示求解降维后单变量函数前几阶矩时所需的高斯积分点数)。这相比于SDP方法和Monte Carlo数字模拟法显得更加高效,极大地节约了模型运算成本。

考虑如下多输出功能函数,

点击浏览公式

式中，x₁,x₂,x₃服从^{[0, 1]}上的均匀分布,且相互独立。计算结果如表 1所示。

表 1中S_i(Y)′表示基于文献[6]作了修改的文献[3]中提出的多输出情况下重要性测度指标;M-MDR表示乘法降维方法,MCS表示Monte Carlo数字模拟法。

从计算结果来看,本文提出的多输出情况下重要性测度新指标S_i^m(Y)得出的输入变量重要性排序与文献[3]中的原指标S_i(Y)以及S_i(Y)′并不一致。S^m_i (Y)给出的重要性排序为X₁>X₃>X₂,而S_i(Y)给出的重要性排序为X₃>X₁>X₂,且输入变量X₂的重要性几乎可以忽略。因为S_i(Y)直接对不同输出的方差进行加和,而利用乘法降维方法解出的y₁、y₂和y₃ 3个输出方差为2.682、0.207、0.119,这使得y₂和y₃ 2个输出情况提供的有效信息几乎完全损失了,最终得出的重要性排序以及输入变量间的相对重要性程度与输出y₁一致。S_i^m(Y)先对原始输出作无量纲化处理,得到的利用乘法降维方法解出的新y₁、y₂和y₃ 3个输出无量纲方差为0.628、0.494、0.540,这样3种输出情况所提供的重要性测度信息都被充分有效利用,X₂对于系统的综合影响得到合理估计,最终得到合理的新重要性排序。因此相比较而言,多输出情况下重要性测度新指标S^m_i (Y)充分考虑了各个输出模式中输入变量的重要性测度信息,更为合理可取。S_i(Y)和S_i(Y)′得出各输入变量间的相对重要性程度差异也表明引入文献[6]的改进措施可以更加完整地反映X_i对每个输出情况下响应量方差所造成的影响,进而使新指标S^m_i (Y)衡量的综合影响更加完整合理。

从计算效率和精度方面来说,乘法降维方法具有很大的优势,只需要3×(3×5-2)=45次模型运算即可得到满意的结果。而运用Monte Carlo数字模拟法求解则需抽取20 000个样本点,进行双层循环调用模型才能得到较好的收敛结果,效率很低,代价昂贵。

复合材料层合结构因其高的比强度、比刚度和低质量等优越性能而在航空航天领域得到了广泛应用。但是,由于复合材料力学性能具有较大的分散性和随机性,因此研究这些不确定性对复合材料结构的输出响应性能影响具有十分重要的意义。为此,本算例将以单层板材料力学性能、铺设角、铺层厚度及加载载荷视为随机输入变量,复合材料结构输出位移和强度比视为多个输出响应量,研究不同输入变量的不确定性对复合材料结构多输出响应量的方差贡献大小。

层合板悬臂梁的上面板受到垂直于面板的压力F=0.1 MPa,梁的截面尺寸和网格形式如图 1,其中a=50 mm,b=40 mm,c=t=2 mm,悬臂梁长度l=600 mm,各单层板厚度为d=0.25 mm,梁腹板的铺层顺序为(45°／-45°／-45°／45°)_S,梁面板的铺层顺序为(45°／-45°／0°／90°)_S。单层板性质数据为E₁₁=135 GPa,E₂₂=8.8 GPa,ν₁₂=0.33,X_t=X_c=1.5 GPa,Y_t=40 MPa,Y_c=246 MPa,S=68 MPa。上述材料性质、铺设角、铺层厚度以及压力作为随机输入变量,且假定相互独立,均服从正态分布,变异系数都为0.06,铺设角的分散性以0.5°的标准差来表示。

图 1 工字梁截面示意图和网格划分

图选项

利用PCL语言控制 Patran／Nastran 有限元分析软件调用上述工字梁模型 17×5-16=69次,求得相应的输出响应最大位移和临界强度比的值;最后运用本文所述方法即可求得基本输入变量对于多输出响应最大位移和临界强度比的总的影响程度。结果如图 2所示。

图 2 复合材料梁的新指标求解结果

图选项

某一随机输入变量对于不同的输出相应模式的重要性有很大不同。影响输出响应最大位移方差的主要输入变量排序为d>F>E₁₁,而影响输出响应强度比方差的主要输入变量排序却为0°>Y_t>d>F。基于响应量方差的输入变量对2种输出方差的综合影响的重要性排序为0°>Y_t>d>F,这个排序与对强度比方差影响排序一致且重要性程度相当,这是此工程实例中由于强度比方差远大于最大位移方差的原因,从而使得原指标完全没有利用到最大位移输出模式提供的信息。因此对于此类的多输出响应的工程问题,文献[3]中提出的多输出情况下重要性测度指标难以提供有效的实际指导意义。而本文提出重要性测度新指标可以有效地避免上述问题,充分

利用每个输出模式提供的重要性信息。对于本工程实例,新指标很好地度量了0°和E₁₁2个输入变量在多输出模式下的综合重要性程度,体现了E₁₁的重要性影响。所以多输出情况下重要性测度新指标合理衡量了输入变量对各个输出响应量方差的综合影响。从计算代价方面看,相比较于提出本例的文献[1]中所采用的SDP法(需1 000次模型调用),新指标求解所需的17×5-16=69次模型调用是十分高效的。

本文主要提出了一种针对结构系统多输出情况下基于方差的重要性测度的新指标,用以衡量各输入变量的变异性对多输出结构系统变异性综合影响的重要程度,指导工程实际按新指标对输入变量进行重要性度量和排序。一方面,新指标中先对原始输出模型作无量纲化处理,从而充分考虑了每个输出模式的作用;另一方面,新指标采用修改的基于方差的主重要性测度指标,从而更完整地保留输入变量对响应量方差的影响信息。本文将乘法降维替代模型引入到多输出情况下基于方差的主重要性测度的新指标,准确高效。使用乘法降维替代模型求解新指标只需要m×[n×l-(n-1)]次功能函数运算就可以获得很精确的结果(其中m表示输出模式个数,n表示输入变量个数,l表示求解降维后单变量函数前几阶矩时所需的高斯积分点数)。这相比于SDP方法和Monte Carlo数字模拟法显得更加高效,极大地节约了模型运算成本。