一种基于有界不确定性的连续体结构非概率拓扑优化方法与流程

本发明涉及连续体结构拓扑优化设计技术领域，特别涉及一种基于有界不确定性的连续体结构非概率拓扑优化方法，该方法考虑材料、外载荷和位移许用值的不确定性，在基于优化特征距离d这一非概率可靠性指标的约束下，对连续体结构进行拓扑优化。

背景技术：

近年来，得益于计算机技术的巨大进步，结构优化技术也得到了长足发展，其应用范围涵盖了航空航天、机械、土木、水利、汽车等诸多领域，优化手段从简单的尺寸参数优化扩展到了形状优化以及更具挑战性的拓扑优化。根据研究对象的不同，结构拓扑优化可以分为离散结构拓扑优化和连续体拓扑优化两大类，其中连续体拓扑优化在最近十多年来成为结构优化领域最活跃的研究方向之一。自1988年基于微结构和均匀化思想的拓扑优化方法被提出来之后，已经发展了多种拓扑优化方法，其中代表性的方法有变厚度法，ESO法(evolutionary structural optimization),SIMP法(solid isotropic material with penalization)等。值得注意的是，现有的拓扑优化研究大部分是基于确定性假设的。

随着科技水平的不断进步，工程结构系统的复杂程度在不断增加，不确定性的表现随之也越来越突出。考虑不确定性因素的结构优化设计日益引起重视，在此环境下，基于概率理论的可靠性优化设计(Reliability-Based Design Optimization,RBO)被提了出来，并用于拓扑优化设计中来，即可靠性拓扑优化(Reliability-Based Top Optimization,RBTO)。

然而，在工程结构系统中广泛存在随机、模糊、未知然而有界等多种不确定性信息，而且结构样本数据常常是缺乏的。因而概率可靠性模型和模糊可靠性模型条件往往不能得以满足。对于实际问题，不确定信息的精确统计数据往往不易获得，但是对于不确定信息的不确定界限却比较容易确定，基于这一思想，基于凸集合模型的非概率可靠性的概念被提了出来。然而在非概率可靠性拓扑优化研究方面的成果较少，现有的研究都是间接地利用非概率可靠性指标，而不是直接使用非概率可靠度的灵敏度作为优化准则，此外，在利用非概率集合可靠性模型作为可靠性指标的拓扑优化方法还是空白。

技术实现要素：

本发明要解决的技术问题是：克服现有技术的不足，提供一种基于有界不确定性的连续体结构非概率拓扑优化方法。本发明充分考虑实际工程问题中普遍存在的不确定性因素，以提出的优化特征距离d这一非概率可靠性度量指标作为优化模型的约束条件，所得到的设计结果更加符合真实情况，工程适用性更强。

本发明采用的技术方案为：一种基于有界不确定性的连续体结构非概率拓扑优化方法，其实现步骤如下：

步骤一：基于一般的拓扑优化数学模型，使用非概率可靠性指标作为约束，建立非概率可靠性拓扑优化的数学模型：

其中，V是优化区域的体积，ρ_i和V_i分别为第i个单元的相对密度和体积，N为优化区域划分的单元总数，ρ为单元相对密度的下限，这是为了防止刚度矩阵奇异而设置的一个小值。d_j是第j个约束的可靠性，m为约束的个数。对于SIMP模型，单元的弹性模量是材料相对密度的函数：

其中P＞1是惩罚因子，用于实现对中间密度单元的惩罚。按照经验，一般取P＝3，E₀是完全实心材料的弹性模量；

步骤二：考虑材料弹性模量，载荷大小与位移许可值的不确定性，采用区间变量K^I和F^I来表示整体刚度区间矩阵和载荷区间向量，采用区间向量来表示位移区间向量。根据有限元的位移控制方程有：

K^Iu^I＝F^I

然后使用区间参数顶点法，由位移关于弹性模量和载荷的单调性，求出位移在有界但不确定参数影响下的上下界：

其中其中下标corj表示位移区间向量u^I中对应于第j个位移约束的分量；上标k_i＝1,2，当k_i＝1时表示对应值取下界，当k_i＝2时表示对应值取上界_，即(K^-1)²＝K^-1，，F_i¹＝F_i，

步骤三：采用非概率集合可靠性模型，将实际位移区间和安全位移区间作标准化变换：

其中u_j,a为第j个位移约束的实际位移，u_j,s为第j个位移约束的安全位移，两者都是区间变量，和分别为实际位移u_j,a和安全位移u_j,s的中值，为区间半径。根据结构功能函数

M(u_j,s,u_j,a)＝u_j,s-u_j,a

来判断结构是否安全；

步骤四：定义优化特征距离d这一非概率可靠性指标。优化特征距离d的定义为：原失效平面到目标失效平面的移动距离。其中目标失效平面是与原失效平面平行的平面，并且其可靠度为一给定值。用优化特征距离d这个指标来量化当前设计的非概率可靠度；

步骤五：使用伴随向量法求解位移上下界的灵敏度，然后根据复合函数的求导法则得到优化特征距离d的灵敏度。

步骤六：采用MMA优化算法，以最小化相对体积为目标，以可靠度为约束，利用优化特征距离d和相对体积的灵敏度进行迭代求解，在迭代过程中，如果当前设计不满足可靠度约束d＜0，或者尽管满足可靠度约束，但相较于上一个可行解，目标函数的相对变化百分比大于预设值ε时，将已经完成迭代次数的值增加一，并返回步骤二，否则，进行步骤七；

步骤七：如果当前设计满足可靠度约束d＜0，并相较于上一个可行解，目标函数的相对变化百分比小于预设值ε时，则迭代结束，将当前拓扑优化的结果作为最终的优化结果。

本发明与现有技术相比的优点在于：

本发明提供了一种基于有界不确定性的连续体结构在位移可靠度约束下的非概率拓扑优化设计的新思路，弥补和完善了传统基于概率理论的可靠性设计方法的局限性。所构建的非概率拓扑优化模型，一方面可大幅减小对样本信息的依赖性，另一方面可有效计及并量化不确定性作用下对拓扑结构构型的影响。在对位移约束下的连续体结构进行拓扑优化设计时，可以充分考虑不确定性作用下的结构拓扑变化规律，在确保结构位移满足一定约束条件下可大大降低结构重量，提高性能的同时，降低设计周期和经济成本。

附图说明

图1是本发明针对基于有界不确定性的连续体结构在位移可靠度约束下的非概率拓扑设计流程图；

图2是本发明所用到的非概率集合可靠性模型的一维干涉模型；

图3是本发明所用到的非概率集合可靠性模型示意图；

图4是本发明所用到的非概率集合可靠性模型的六种不同的干涉情况示意图，其中，图4(a)中图4(b)中图4(c)中图4(d)中图4(e)中图4(f)中

图5是本发明用于所提出的优化特征距离d的计算的两种临界斜率示意图；

图6是本发明针对连续体结构拓扑优化的几何模型示意图；

图7是本发明针对连续体结构拓扑优化的优化结果示意图，其中，图7(a)为确定性优化，图7(b)为非概率可靠性优化(R＝0.90)，图7(c)为非概率可靠性优化(R＝0.95)，图7(d)为非概率可靠性优化(R＝0.999)；

图8是本发明针对连续体结构拓扑优化迭代历程曲线，其中，图8(a)为确定性优化，图8(b)为非概率可靠性优化(R＝0.90)，图8(c)为非概率可靠性优化(R＝0.95)，图8(d)为非概率可靠性优化(R＝0.999)。

具体实施方式

下面结合附图以及具体实施例进一步说明本发明。

如图1所示，本发明提出了一种基于有界不确定性的连续体结构非概率拓扑优化方法，包括以下步骤：

(1)考虑材料，外载荷和位移许用值的不确定性，建立包含非概率可靠性约束的优化问题的数学模型：

其中，V是优化区域的体积，ρ_i和V_i分别为第i个单元的相对密度和体积，N为优化区域划分的单元总数，ρ为单元相对密度的下限，这是为了防止刚度矩阵奇异而设置的一个小值。d_j是第j个约束的可靠性，m为约束的个数。对于SIMP模型，单元的弹性模量是材料相对密度的函数：

其中P＞1是惩罚因子，用于实现对中间密度单元的惩罚。按照经验，一般取P＝3，E₀是完全实心材料的弹性模量；

(2)考虑材料，外载荷与位移许用值的不确定性，采用区间变量K^I和F^I来表示整体刚度区间矩阵和载荷区间向量，采用区间向量来表示位移区间向量。根据有限元的位移控制方程有：

K^Iu^I＝F^I

其中u^I为位移区间向量。由于控制方程是线性的，可以采用下面的区间参数顶点法来求解中任意分量的上下界。

区间参数顶点法：如果f(x₁,x₂,…,x_n)对自变量x_i(i＝1,2,…,n)是单调的，将自变量考虑为区间变量时，即：

由函数的单调性可知，f的取值范围为：

其中r为顶点(区间两端点)组合序数，k_i＝1,2，i＝1,2,…,n；r＝1,2,…,2ⁿ。

所以，根据区间参数顶点法，得到第j个约束对应的位移的取值区间为：

其中其中下标corj表示位移区间向量u^I中对应于第j个位移约束的分量；上标k_i＝1,2，当k_i＝1时表示对应值取下界，当k_i＝2时表示对应值取上界，即(K^-1)²＝K^-1，，F_i¹＝F_i，

(3)通过步骤(2)求出位移的上下界之后，就可以基于结构非概率集合可靠性模型，建立以下位移约束下的非概率可靠性模型。

设u_j,a为第j个位移约束的实际位移，u_j,s为第j个位移约束的安全位移，将其取为基本区间变量，即：

将以上两个区间在同一个数轴上表示出来，由于不确定性的存在，两者可能会存在相交的区域，如图2所示。图2中和分别为实际位移u_j,a和安全位移u_j,s的中值。设结构功能函数为：

M(u_j,s,u_j,a)＝u_j,s-u_j,a

其失效平面或极限状态平面为：

M(u_j,s,u_j,a)＝u_j,s-u_j,a＝0

其中，M(u_j,s,u_j,a)＞0表示结构满足约束；M(u_j,s,u_j,a)＜0表示结构不满足约束。对实际位移和安全位移区间变量u_j,a∈u_j,a^I、u_j,s∈u_j,s^I做标准化变换：

其中，则有δu_j,a∈[-1,1]，δ^u_j,s∈[-1,1]。将上式代入失效平面方程，有：

由此可以得到δu_j,s和δu_j,a之间的关系式为：

将上式在直角坐标系中画出，并标示出δu_j,s和δu_j,a的取值范围，如图3所示。

为了求出失效平面分割的两个面积，首先针对一种情况，求解出失效平面与直线δu_j,s＝1的交点。在式(1)中，令δu_j,s＝1，可求解出δu_j,a为：

令可解得

接着求出失效平面与直线δu_j,a＝-1的交点。在式(1)中，令δu_j,a＝-1，可求解出：

令可解得

将满足约束条件的区域面积S_AEF与变量区域的总面积S_ABCD之比定义为结构的非概率可靠度R，则R的表达式为：

代入式(2)和式(3)，得：

考虑到失效平面与变量区域的相交形式的不同，下面分六种情况给出其示意图，如图4所示。六种情况下的R的表达式为：

(4)对于以面积比定义的可靠度，虽然比较直观，但是当失效平面与变量区域处于图4的(e)和(f)两种情况时，非概率可靠度R为常数(0或1)，因此对设计变量的偏导数为零。对于梯度优化算法来说，R在较大范围内梯度都为零，这将会导致优化算法无法找到优化方向，使得优化过程难以收敛。为了解决这一问题，下面将引入优化特征距离d这一非概率可靠性指标。

优化特征距离d的定义为：原失效平面到目标失效平面的移动距离。其中目标失效平面是与原失效平面平行的平面，并且其可靠度为一给定值。如图5，由于可靠度一般接近1，故目标失效平面一般位于变量区域的右下角，图示为目标失效平面与变量区域相交形式的两种临界情况。

首先计算临界情况下失效平面的斜率。对于k₁，有(2×2/k₁×1/2)/4＝1-R，解得k₁＝1/2(1-R)，同理可得k₂＝2(1-R)，针对原失效平面斜率k取值的不同情况，使用直线间的距离公式，并定义原失效平面在目标失效平面上方的距离为正，反之为负，给出优化特征距离d的表达式为：

当d＞0时，失效平面在与目标非概率可靠度R对应的目标失效平面上方，此时由于安全区域的面积小于目标值，对应的非概率可靠度R_d＜0＜R，不满足要求。当d≤0时，失效平面在与目标非概率可靠度R对应的目标失效平面下方，此时由于安全区域的面积大于等于目标值，对应的非概率可靠度R_d＜0≥R，满足设计要求。

(5)本发明采用基于函数梯度的数学规划法——移动渐进线方法(MMA)求解优化问题，因此需要进行目标函数和约束函数(非概率可靠性的性能值)对设计变量(材料相对密度)的灵敏度分析。连续体结构拓扑优化问题一般考虑刚度、频率等全局性约束条件，因而设计变量数目远多于约束条件数目。针对这一特点，本发明采用伴随向量法，实现非概率可靠性拓扑优化中约束函数值对设计变量(材料相对密度)的灵敏度分析。

对于第j(j＝1,2,…,m)个约束，其优化特征距离d_j对单个设计变量x_i(i＝1,2,…,N)的全导数为：

其中：

对于可以直接计算，由于设计变量x_i的数量很多，如果直接采用差分的方式计算和需要进行2N次系统重分析，计算量十分巨大，为了避免直接计算和使计算效率更高，构造如下的约束函数的增广拉格朗日函数：

其中，λ_j(j＝1,2,…,m)为与平衡方程关联的任意的乘子向量，也被称为伴随向量。由于F-Ku＝0，故上式对设计变量v_i求全导数得：

其中：

上式对任意λ均成立，因此可以选取适当的λ使得du/dx_i所在项的系数为零，即令：

利用刚度矩阵的对称性可以将上式改为：

对比有限元的位移控制方程可知，可以仿造有限元位移求解的过程，将视为虚拟载荷，则可将虚拟位移λ求出，这可以利用有限元程序来求解。求解出λ后，约束点位移上下界对设计变量的灵敏度则由下式给出：

其中分别为对应于的伴随向量、总体刚度矩阵和位移向量，λ_j、K_j、u_j分别为对应于u_j,a的伴随向量、总体刚度矩阵和位移向量。在本发明的优化模型中，载荷F不随设计变量变化，即dF/dx_i＝0，则上式可以改写为：

所以约束函数d_j(j＝1,2,…,m)对设计变量的灵敏度为：

此外，目标函数V对设计变量的偏导为：

(6)采用MMA优化算法进行迭代计算，根据当前的单元相对密度、约束函数d对设计变量的灵敏度、目标函数V对设计变量的灵敏度，求解出新的设计变量。在迭代过程中，如果当前设计不满足可靠度约束d＜0，或者尽管满足可靠度约束，但相较于上一个可行解，目标函数的相对变化百分比大于预设值ε时，将已经完成迭代次数的值增加一，并返回步骤二，否则，进行步骤七；

(7)同时考虑可靠度约束和相对变化量，如果当前设计满足可靠度约束d＜0，并相较于上一个可行解，目标函数的相对变化百分比小于预设值ε时，则迭代结束，将当前拓扑优化的结果作为最终的优化结果。

实施例：

为了更充分地了解该发明的特点及其对工程实际的适用性，本发明针对如图6所示的矩形平板进行拓扑优化设计。设计区域为1.4m×0.3m的矩形区域，厚度为0.001m，划分为140×30个单元。材料弹性模量E＝100Mpa，泊松比μ＝0.3。矩形区域的左右两侧固定，施加P₁＝P₃＝2.5N，P₂＝3N三个垂直向下的力，不考虑重力的影响，约束三个加载点的位移，使得u₁<0.8mm，u₁<1.5mm，u₁<0.8mm，选取惩罚因子p＝3。设弹性模量E、载荷F相对名义值均有10％的波动，即E＝[90,110]Mpa，P₁＝[2.25,2.75]N，P₂＝[2.7,3.3]N，P₃＝[2.25,2.75]N；设位移约束u₁、u₂、u₃相对名义值有1％的波动，即u₁＝[0.72,0.88]mm，u₂＝[1.35,1.45]mm，u₁＝[0.72,0.88]mm。

图7分别为确定性拓扑优化结果和非概率可靠性分别为R＝0.90、R＝0.95和R＝0.999时拓扑优化结果对比。可以看到，确定性拓扑优化和不同非概率可靠性拓扑优化出来的结构的构型存在较大的区别，相比于确定性拓扑优化结果，非概率可靠性拓扑优化结果结构更加合理，细长杆较少，结构更为稳定。在使用与非概率可靠性一样的不确定量时，确定性优化结果的三个被约束位移的非概率可靠性仅分别为R₁＝0.4504、R₂＝0.4501、R₃＝0.4501。即确定性优化的结果不足以应对不确定性变量的影响。拓扑优化过程中的迭代历史如图8所示，相较于初始设计，减重效果明显；随着可靠度许用值增加，结构趋于安全，重量有所增加。

综上所述，本发明提出了一种基于有界不确定性的连续体结构非概率拓扑优化方法。首先，建立连续体结构在非概率可靠性约束下的拓扑优化数学模型，然后考虑材料、外载荷和位移许用值不确定性，使用区间参数顶点法计算出位移的上下界；接下来基于非概率集合可靠性模型，建立优化特征距离d这一新的非概率可靠性指标；其次，使用伴随向量法，并结合复合函数的求导法则，求出优化特征距离d的灵敏度；最后使用MMA优化算法，以可靠度为约束，以相对体积为目标，进行优化迭代计算，从而完成了连续体结构的非概率拓扑优化设计。

以上仅是本发明的具体步骤，对本发明的保护范围不构成任何限制；其可扩展应用于含缺陷结构的优化设计领域，凡采用等同变换或者等效替换而形成的技术方案，均落在本发明权利保护范围之内。

本发明未详细阐述部分属于本领域技术人员的公知技术。