基于少样本不确定性的结构可靠性拓扑优化设计

何佳琦，贾晓璇，钟杰华，许元男，罗阳军

（1. 大连理工大学航空航天学院，大连，116024；2. 中国运载火箭技术研究院，北京，100076）

0 引 言

实际航空航天等重要工程中，广泛存在大量的不确定性信息，如载荷不恒定、材料属性波动以及结构尺寸误差等。为了避免由不确定性因素导致的结构性能偏差，工程上在结构设计阶段常采用安全系数法来得到一种较为保守的结构设计形式。然而，传统的安全系数法可能设置过大的安全系数导致设计过于保守，也可能因未全面考虑多种不确定性因素的耦合而使得优化结构面临失效的风险。

结构可靠性优化设计是以结构可靠性作为结构的安全裕度度量开展结构优化设计。该类优化设计的统一思路是：首先基于先验知识与样本数据将不确定信息量化为数学模型；进而根据不确定性变量与结构安全性能指标之间的函数关系建立不确定性传播模型；最后以结构安全性指标可靠性作为目标函数或约束函数，采用梯度方法或智能优化算法开展结构参数优化或拓扑优化工作。大量学者采用概率模型对不确定信息进行量化，进而以结构失效概率或可靠度作为结构可靠性度量开展结构可靠性优化设计[1,2]。然而，基于概率模型的结构优化设计的合理性极大依赖于概率模型的准确性。通过结构优化设计得出的结构形式往往逼近力学性能极限，而用来量化制造误差、载荷等不确定性的概率模型往往会由于样本数量的局限难以保证其准确性。针对少样本的不确定信息量化问题，目前已有一些学者提出了根据样本信息构建凸集模型的有效方法[3～5]。进而，一些学者以非概率可靠性指标作为结构可靠性度量开展结构优化设计工作[6,7]。

对于结构优化设计方法，在考虑不确定性因素后，很难通过经验判断最优结构拓扑形式，因此，合理的方式是直接采用考虑不确定性因素的拓扑优化方法开展结构设计工作。对于包含较大数目有限元单元数量的可靠性拓扑优化问题，已有的可靠性拓扑优化算法往往局限于传统的拓扑优化方法[8～10]而导致计算效率较低。罗阳军等提出了一种基于材料场级数展开（Material-Field Series-Expansion，MFSE）表征结构拓扑演变的新思路[11]，使结构拓扑优化问题的设计变量完全摆脱有限元网格依赖性，并实现了拓扑优化设计变量的显著降低。该拓扑优化方法已经成功应用于力学问题[12]与光学等。

本文提供了一种基于少样本的结构可靠性拓扑优化设计方法，寻求在满足结构可靠性要求下，质量最轻的结构拓扑优化形式。首先根据样本数据，基于非概率凸集模型对不确定性变量进行量化，进而采用功能度量法将可靠性指标约束转化为功能度量以提高优化过程的稳健性。其次，采用一种序列单循环的方法将原嵌套优化问题解耦，使用MFSE方法大幅度减少设计变量数目提高优化计算效率。分别提供了二维平面应力单元算例、三维实体单元算例来验证本文提供的方法的合理性与适用性。

1 基于样本数据构建非概率凸集模型方法与非概率可靠性指标定义

对于凸集模型，大体可以分为区间模型和椭球模型两类。n维区间模型从几何形貌上来讲是一个n维超立方盒，其内在地假定了所有的不确定性参数互不相关且各自独立。而在实际工程中，对于如载荷一类不确定性变量之间一定存在某种相关性，极少存在所有不确定性变量同时达到界限的情况。椭球模型可以在一定程度上反映不确定性变量之间的相关性，因此对于样本数量不充足且相互间存在相关性的不确定性变量，适合用多椭球模型来进行描述。构建多椭球模型，首先需要根据相关与否将不确定性变量分成若干组，进而对每一组变量构建单椭球模型。多椭球模型本质上是多个单椭球模型的集合，其数学表达式为

1.1 基于样本数据的多椭球凸集模型构建方法

对于不确定性变量X，其单椭球模型Ω的数学表达式为

将矩阵W做特征值分解：

式中Σ为矩阵W的特征值组成的对角阵，对角线元素代表各个维度的半轴长；Q为对应的特征向量矩阵，其正交向量由Gramm-schmidt Orthogonalization公式产生[8]，是主轴方向的函数，主轴方向向量包含Ni-1个方向角分量。

对于不确定性变量X，通过采样可获得样本集｛xi,i= 1,2,…,M｝，用多维椭球模型量化不确定性变量本质上就是找到体积最小的能包络所有不确定性变量可能取值范围的超椭球模型。当椭球模型方向角θ给定时，可将每一个样本点xk在多椭球主轴方向与其它与主轴垂直方向上投影，可得到样本点在第l轴上投影的区间长度dl，即：

式中为xi在椭球模型第l个轴上的投影。在得到样本点在椭球模型各轴方向上的投影范围区间长度dk后，包络上述多维区间模型的最小多维椭球半径可表示为

多维椭球的体积正比于eV：

上式是参数θk（k= 1,…,N-1）的函数，可采用基于Kriging代理模型的EI序列加点算法得到使Ve最小的一组参数θk（k= 1,…,N-1）。

1.2 非概率可靠性指标定义

概率可靠性和非概率可靠性作为结构可靠性分析的两个分支，都是用来衡量结构对各种不确定性因素的抵抗能力。概率模型强调可接受行为的概率，而非概率可靠性强调的是可接受行为的范围。因此，非概率可靠性度量可以理解为结构性能指标可接受的波动范围。

首先将不确定性变量X无量纲化，即进而将多椭球模型投影为单位圆，即超曲面（失效面）g（q）=0将整个设计空间划分为安全域与失效域，可将坐标原点到失效面g（q）=0的最短距离定义为非概率可靠性指标β：

式（7）中如果 1β=，则表示失效面与不确定性变量范围相切，结构处于临界失效状态；如果 1β＜，则表示不确定性变量范围部分或全部处于失效域，结构处于失效状态；如果 1β＞，则表示不确定性变量范围全部处于安全域，结构非概率可靠性指标β越大，表示不确定性范围距离失效域越远，结构越安全。

2 非概率可靠性拓扑优化定义与求解方法

2.1 基于样本数据的多椭球凸集模型构建方法

非概率可靠度拓扑优化问题可归纳为下述优化列式：

式中ρ为单元伪密度向量；q为标准化后的不确定参数；g j（ρ,q）为结构的功能函数。

该优化模型中，可靠性指标约束条件的处理较为困难。为增加优化迭代过程的稳定性和提高收敛效率，通常可采用功能度量方法PMA将可靠性指标约束转换为其更为稳健的等效形式。功能度量方法并不直接比较原优化问题中的可靠性指标与其目标下限值，而是在不确定性允许的分布范围内找出使最可能失效功能函数值最小的点，也可称MPP点，要求该最小功能函数值非负，即

上述可靠度拓扑优化问题为一个嵌套的优化模型。两层循环方法的计算量非常巨大，并不适合于大型复杂问题的求解。为解决直接求解嵌套优化问题计算花费较大的困难，采用序列单循环方法，通过构造一系列近似的确定性子优化问题，将内外层优化连续求解，该序列确定性子优化问题与求解MPP点的优化列式分别为式（10）与式（11）：

对于优化列式（11），可将第j个目标函数在q（k）处做一阶泰勒展开，即：

2.2 基于MFSE的拓扑优化方法

针对式（10）的优化问题，可采用传统SIMP方法求解，即用每一个有限元单元的有无表征结构拓扑形式。考虑到实际航空航天结构有限元网格数量巨大，由此会造成优化问题设计变量ρ数目巨大。本文采用采用MSFE方法求解确定性拓扑优化问题，即用一个材料场（空间有界场）表征结构拓扑形式，见图1。

图1 材料场表征结构拓扑形式示意 Fig.1 Schematic illustration of MFSE Method

材料场级数展开表达式为

式中Npt为观察点数量；η为设计变量；λ与ψ分别为相关矩阵C的特征向量与特征值，相关矩阵C由式（15）、式（16）定义：

式中lc为材料场相关长度。

图2为不同相关长度下的材料场，其中图2b材料场相关长度大于图2a材料场中的相关长度。从几何上来看，大于中相关长度越大，材料场越平滑。在拓扑优化过程中，相关长度可起到控制拓扑结构细节程度的作用，即，相关长度越小，拓扑结构细节程度越高、细节构件尺度越小。

图2 不同相关长度下的材料场 Fig.2 Material Fields with Different Correlation Lengths

式（14）中的材料场级数展开本质上是用个场函数的线性叠加来构成材料场。考虑到较小的特征值与其对应的特征向量对材料场影响较小，因此可给定一小值ε，经过判断若则将材料场函数截断为M项，即：

考虑到φ（x）的值域是（ -∞,+∞），可基于Sigmoid函数将φ（x）投影到区间［0,1］，即：

式中 参数β可用来调整随φ（x）增大或减小（x）趋近于0或1的速率。β取值越大（x）趋近于0的速率快，同时Sigmoid函数非线性程度越高。β取值可在迭代过程中逐步增加使最终的（x）取值趋近于0或1。当空间位置x处，（x）=0认为该处为空相，反之，（x）=1认为该处结构为实相。

室内的体验效果比模型和图纸的真实感和沉浸感更为强烈，不过其弊端是不能像在真实样板房中感受各种材质的质感。可以通过提升渲染的质量来弥补虚拟环境中造成的色差。

3 非概率可靠性拓扑优化求解流程

可靠性拓扑优化算法实施流程如图3所示。a）根据第1.1节阐述的方法将样本数据量化为多椭球凸集模型；b）根据结构设计域形式划分独立于有限元网格的观察点，并根据1.2节阐述的方法构建材料场函数；c）根据式（13）更新不确定性变量；d）基于移动渐进线法更新设计变量；e）重复步骤c与步骤d，直至满 足 收 敛 条 件 ：与，迭代结束，输出满足可靠性指标的最优结构拓扑形式。

4 算 例

4.1 三维算例

算例设计域为长方体结构，如图4所示，长度l= 48 mm，宽度w=24 mm，高度h=24 mm。该结构在左端面固定，承受3个面内不确定性载荷F1、F2与F3，载荷间存在相关性。结构材料弹性模量存在不确定性，其变化区间为E∈[1.2×1011，2.1×1011] Pa。

图4 三维算例结构与载荷形式 Fig.4 Design Domain for 3-dimensional Structure with Multiple Loads

根据不确定性样本信息（如图5中所示样本点），可将载荷不确定性与材料弹性模量不确定性构建为多椭球模型，即一个由用来表征载荷变量的三维椭球与表征材料弹性模量变量的区间所组成的集合，三维椭球模型如图5所示。基于MFSE算法可将设计变量从110 592下降到1420。给定可靠性指标设计值为2，功能函数为结构最大位移不超过1 mm。基于本文提出的非概率可靠性拓扑优化算法得到结构拓扑构型见图6，最优结构体积分数为18.6%，迭代过程如图7所示。

图5 样本点与三维椭球模型 Fig.5 3-dimensional Ellipsoid Model with Samples

图6 三维算例最优拓扑结构设计 Fig.6 The Optimized Topology Design of 3-dimensional Structure

图7 三维算例迭代过程 Fig.7 Iteration History for 3-dimensional Structure with Multiple Loads

4.2 圆柱壳正交加筋算例

圆柱壳直径d=500 mm，高度h=300 mm、厚度t=1.5 mm。针对壳内壁开展正交加筋形式拓扑优化设计，蒙皮与加强筋材料弹性模量皆为E=2.1×1011Pa。蒙皮厚度为0.5 mm，加强筋高度为1.5 mm。圆柱壳结构下端面固定，上端面承受2个不确定性剖面集中载荷，分别为相互垂直的剪力Px和Py。载荷作用点距离上端面l=300 mm。该算例结构与载荷如图8所示，样本点与量化其不确定性的椭圆模型如图9所示。可靠性指标设计值为2，功能函数为结构最大位移不超过 2 mm。

图8 圆柱壳加筋结构与剖面载荷示意 Fig.8 Schematic Illustration of Stiffened Cylindrical Shell Structure with Profile Loads

为实现工程中常见的正交加筋形式，可采用两个单向材料场叠加的方式。对于单向材料场的表征，只需将式（16）中正交方向的相关长度设置为无穷大即可。对于两个材料场与相互覆盖区域的其插值函数可取为

基于本文提出的非概率可靠性拓扑优化算法得到圆柱壳加筋结构最优正交加筋拓扑构型如图10所示，其体积分数为37.4%。

图10 最优正交加筋拓扑构型 Fig.10 The Optimized Topology Design of the Orthogonal Stiffened Cylindrical Shell Structure

采用传统设计方法往往忽略载荷间的关联性，采用极值载荷（或在载荷平均值基础上乘以可包络极值的安全系数）开展设计，即不确定性载荷分别取为：Px=14 kN，Py=9 kN；而采用非概率可靠性设计，在给定非概率可靠性指标范围内，最危险工况对应载荷为：Px=13.93 kN，Py=8.27 kN。相比较于传统设计方法，非概率可靠性设计方法对应的工况可以与可靠性指标一一对应起来，而传统设计方法选取的载荷工况往往缺乏理论依据且容易造成过于保守的情况。此外，基于材料场表征拓扑构型，在大幅度缩减设计变量的同时，易于实现工程上常用的正交加筋构型或其它规则加筋形态的拓扑优化设计。

5 结束语

本文针对实际航空航天结构设计过程中，广泛存在的不确定性样本少且结构有限元模型规模巨大的情况，提出了一种高效的可靠性拓扑优化设计方法。首先将不确定性样本量化为多椭球模型，并根据结构设计域构建材料场函数，进而结合功能度量法与序列单循环策略求解优化问题，最后将最优材料场映射回原设计空间得到最优结构拓扑形式。算例验证了基于MFSE法减少设计变量规模的效率以及本文提出的可靠性拓扑优化设计方法的适用性与合理性。