CN112784489B - 一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法 - Google Patents

Abstract

本发明公开了一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法。该方法首先根据连续体结构的服役特点，基于非概率集合模型考虑有限样本条件参数的不确定性效应，以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以多元二次曲面样条函数的系数作为设计变量，建立考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型，基于优化准则法求解优化模型，通过迭代逐渐获得在给定外载和边界条件下的水平集函数，进而确定结构的构型。本发明针对拓扑优化中考虑不确定性参数影响的动力学响应分析的计算量巨大的问题，在三方面采用有效地手段提高计算效率，最终完成一定目标的减重并且使得结构整体响应优化。

Description

一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法

技术领域

本发明涉及结构拓扑优化设计技术领域，特别涉及一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法，该方法考虑材料弹性模量、外部载荷存在一定的范围的不确定性，建立了以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型。并且使用优化准则法(OC：Optimality Criteria method)求解了该优化模型。

背景技术

拓扑优化对应结构优化的概念设计阶段，是在结构设计最开始的阶段，在给定荷载条件和边界条件的情况下，找到设计域内满足约束条件并且能达到设计目标最优化的材料分布情况从而确定结构在满足设计要求下的最佳传力路径的设计手段。借助于拓扑优化，设计者可以更容易提出初步的设计方案，缩短设计周期。相较于尺寸优化而言，拓扑优化需要谋求设计域内的满足设计要求并且达到性能最优的最佳的材料分布方案，从而获得结构设计所达到的目标函数的最值，设计变量更多，设计的可能性更多，这些特点使得拓扑优化问题求解也更加困难。拓扑优化成为结构优化领域中难度最大、水平更高的研究课题。在诸多工业领域，拓扑优化的应用成效已经日渐显著，逐渐体现出其在结构设计阶段的优越性。例如卫星以及火箭的部分结构、民机以及军机的部分结构、汽车和船舶的某些结构的设计、微机电系统都已经引入拓扑优化的方法，以谋求经济效益和工程材料损耗的“同步极限”。经过几十年的发展，连续体结构拓扑优化获得极大地发展，研究人员提出了很多拓扑优化的方法，例如均匀化方法、渐进结构法(ESO)、SIMP法和水平集法等多种方法。其中参数化水平集方法由于其突出的优点，近年来经过不断的发展已经成为拓扑优化的主流方法之一。

由于结构服役环境的时变特性，设计一个结构在不同时间范围内保持良好的功能是结构动力学的关键因素。动力学性能设计是目前结构设计技术中极具挑战性的研究方向，是各类航空航天飞行器(飞机、运载火箭、卫星、飞船、空间站等)以及其它高端机械装备性能提升的关键技术之一。大量事实表明，结构的动力学性能设计直接影响航空航天飞行器的重量、刚度分布及其它性能，不合理设计甚至导致各类严重事故，因此研究结构动力学拓扑优化具有重要的理论意义和工程应用价值。

目前工程结构的优化设计方法大多基于确定性数学模型的拓扑优化方法，为了研究的方便，往往会忽略实际工程中存在诸多不确定性因素的现实情况，这就使得优化的结果与实际应用情况存在一定的差距，这些不确定性因素的累积可能会对结构产生显著的影响，结构甚至会因此失效。某型直升机油箱原采用复合材料，由于诸多不确定因素的影响最终更改为金属材料。结构设计方法经过几十年的发展，研究人员发现鲁棒性优化方法是分析工程中不确定性因素的影响的实用且有效手段。不确定性表征方法主要可以分为概率模型、模糊模型和非概率模型。基于非概率理论的不确定性分析方法针对贫信息，少样本的不确定参数具有很好的处理效果。

将非概率理论与动力学拓扑优化相结合是进行拓扑优化时更加合理也更加方便的考虑，这样所得到的拓扑优化设计结果考虑了一定参数的不确定性，使得拓扑优化设计更加合理更加贴近现实。

发明内容

本发明要解决的技术问题是：克服现有技术的不足，提供一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法。本发明采用参数化水平集法模拟拓扑优化设计中结构的构型的改变，采用等效静态载荷法处理多次求解动力学平衡方程的问题，采用基于区间模型的顶点法表征和分析不确定性的影响，采用优化准则法求解优化模型，建立考虑材料弹性模量、外部载荷存在一定范围的不确定性的动力学鲁棒性拓扑优化方法，在动力学鲁棒性拓扑优化的不确定性表征分析、动力学平衡方程求解、拓扑优化中结构模拟和优化模型求解阶段采用高效手段实现针对连续体结构的动力学鲁棒性拓扑优化设计。本发明充分考虑实际工程问题中普遍存在的不确定性因素，以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ(MultiQuadric spline)样条函数的系数作为设计变量，建立考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型，在动力学鲁棒性拓扑优化的不确定性表征分析、动力学平衡方程求解、拓扑优化中结构模拟和优化模型求解阶段采用高效手段实现针对连续体结构的动力学鲁棒性拓扑优化设计。所得到的设计结果更加符合真实情况，工程适用性、安全性和经济性得到平衡和兼顾。

本发明采用的技术方案为：一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法，用于连续体结构的拓扑优化设计中，实现不确定性影响下的平面连续体结构受时变载荷作用以及固支、简支边界约束的鲁棒性拓扑优化设计，在拓扑优化问题中基于参数化水平集法和顶点法求解区间形式的等效静态平衡方程，并且考虑了材料的弹性模量和载荷数值存在的有界不确定性，实现步骤如下：

步骤一：给定平面连续体结构、时变载荷以及约束条件和设计域与非设计域，基于拓扑优化数学模型，以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，建立考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型；

步骤二：考虑材料弹性模量，外部载荷存在的不确定性，采用区间变量来表示整体刚度区间矩阵、整体阻尼区间矩阵、整体质量区间矩阵和载荷区间向量，采用区间向量来表示位移、速度和加速度区间向量，根据有限元控制方程的区间形式，使用顶点法，求出初始结构的时变载荷影响下的位移、速度和加速度在有界不确定参数影响下的上下界；

步骤三：根据等效静态载荷法，计算各个离散时间点处的结构的等效静态载荷，建立与动力学平衡方程等效的多个离散时间点处的静态平衡方程，将拓扑优化中动力学平衡的状态约束条件转化为多个离散时间点处的等效静力学平衡约束；根据顶点法，求解考虑不确定性参数影响的区间形式等效静态平衡方程，获得结构在离散时间处的各个结点位移区间；

步骤四：依据形状导数计算优化模型中以平均动柔度表示的目标函数关于水平集函数演化伪时间的灵敏度和优化模型中约束条件关于演化伪时间的灵敏度；

步骤五：采用优化准则法，求解以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，所建立的考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型，构造水平集函数的演化速度场；

步骤六：基于演化速度场，求解水平集函数的演化方程，通过更新设计变量以更新水平集函数模拟结构构型的不断改变；

步骤七：如果当前设计满足材料占比约束条件并且满足迭代稳定性条件，则迭代结束，将当前拓扑优化的结果作为最终的优化结果；如果不满足材料占比约束条件或者迭代不稳定，则重复步骤三至步骤六；

所述步骤一：基于一般的拓扑优化模型，建立基于参数化水平集法的动力学鲁棒性拓扑优化的数学模型：

这表示优化模型以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量。针对线弹性结构以材料占比为约束条件，通过MQ样条函数对水平集函数插值进而以MQ样条函数的系数作为设计变量，考虑一般力学问题中载荷的时变特性以及工程中广泛存在的多源不确定性对连续体结构进行拓扑优化设计。

其中，

表示第i个控制点处的MQ样条基函数的系数，下标N表示控制点的个数，

是为了保证计算稳定性而附加的一阶多项式系数，x是控制设计域内某点的位置变量

是优化的目标函数，

是多源不确定性影响下的结构在时变载荷作用下的平均动柔度的区间中心值，

是多源不确定性影响下的结构在时变载荷作用下的平均动柔度的区间半径，ξ₁,ξ₂分别代表平均动柔度区间半径和区间中心值的权重，可以通过调整两个参数的数值大小以反映对不同部分的关注程度，

是随着伪时间

不断演化的水平集函数，上标I表示该变量为区间形式，上标c表示区间的中心值，上标r表示区间的半径，b表示有界但不确定性的不确定性变量，u表示结构的位移响应，K表示结构基于有限元离散的整体刚度矩阵，F_eq表示基于等效静态载荷法的离散时间t_i处的结构等效静态载荷列向量，

表示离散时间t_i处的弹性结构平衡状态方程；D表示全部的结构设计域，

是描述结构材料占比的约束的函数，V_targ表示拓扑优化设计的材料占比约束条件。

表示水平集函数

的阶跃函数，形式为

水平集法是将结构空间Rⁿ中的边界隐式表示到高一维空间Rⁿ⁺¹的水平集函数中去，将比结构空间高一维空间Rⁿ⁺¹中水平集函数的零水平集作为结构空间Rⁿ的边界描述，对二维问题而言，界面曲线

表示为一个三维标量函数

的零水平集：

其中x表示二维空间中某点的坐标(x,y)，Ω表示结构设计域内存在材料的区域，

表示结构设计域内不存在材料的区域，

表示结构设计域内的边界。水平集函数演化方程为：

这是一个哈密顿-雅克比(H-J：Hamilton-Jacobi)偏微分方程，V_n表示水平集函数零水平集的法向演化速度，

表示水平集函数梯度，通过构造合适的法向演化速度场，可以驱动水平集函数的演化，使得结构边界向特定的方向演化，以满足拓扑优化设计的需要。该H-J偏微分方程的求解需要满足一系列条件不利于求解，通过引入MQ样条函数对水平集函数进行插值，可以将H-J偏微分方程转化为常微分方程求解，提高计算效率。利用MQ样条基函数对水平集函数插值表示：

其中，g_i(x)为第i个控制点处的MQ样条函数，MQ样条函数以及其在X方向的梯度如图2所示，

为保证计算稳定性的一阶多项式，将MQ样条基函数插值表示的水平集函数代入的水平集函数的演化方程，以上面的H-J偏微分方程将转化为常微分方程：

其中：

所述步骤二：针对某个确定的水平集函数，考虑材料弹性模量，外部载荷存在一定范围的不确定性，采用区间变量M^I，K^I，C^I和F^I来表示整体刚度区间矩阵,整体阻尼区间矩阵和载荷区间向量，采用区间向量

来表示离散时间t_i处的结点位移区间向量，根据有限元的位移控制方程有：

其中

为结构结点的加速度区间向量，

为结构结点的速度区间向量，基于不确定性分析方法中的顶点法以及动力学响应分析方法中的Newmark-β算法，首先计算加速度

在有界但不确定参数影响下的上下界：

其中

分别代表遍历Nd个不确定性参数上下界的2^Nd次动力学响应分析中的第j个结点处的加速度最小值和最大值，以此确定加速度区间向量。基于Newmark-β算法和顶点法可以确定性速度和位移的区间向量。

所述步骤三：根据等效静态载荷法，有以下等式成立：

其中F_eq由结构的外部载荷、阻尼载荷、惯性载荷三部分组成，M,C,K,F分别为确定性条件下的整体质量、阻尼、考虑诸多不确定性的影响，基于等效静态载荷法，将拓扑优化中动力学平衡方程的求解转化为了有限个静力学平衡方程的求解，这样的转化提高了计算的效率。上述方程将变为区间形式的平衡方程：

整体刚度区间矩阵K^I直接由不确定性参数b^I影响，由等效静态载荷法，各个离散时间t_i处的等效载荷区间向量具有以下表达式：

其中阻尼采用瑞利阻尼C^I(b^I)＝α_mM^I(b^I)+α_kK^I(b^I)，α_m,α_k代表刚度和质量对阻尼的影响系数，

和

分别代表初始结构的加速度和速度区间向量，两者可以将初始结构的水平集函数代入步骤二求解获得初始结构在时变载荷影响下的响应区间向量。整体阻尼区间矩阵的上下界受整体质量区间矩阵和整体刚度区间矩阵的影响：

C^U＝α_mM^U(b^I)+α_kK^U(b^I)

C^L＝α_mM^L(b^I)+α_kK^L(b^I)

进而可以确定各个离散时间处的等效静态载荷：

整体位移区间向量表示为：

其中：

通过以上分析，得出如下结论:1).每次拓扑优化的计算代价与参数不确定性参数Nd数量呈指数关系。2).每次拓扑优化的计算成本与单元数N和离散次数Nt之间存在乘积关系。因此，采用简单的不确定性分析方法来降低计算成本是非常重要的，这也是采用顶点法的原因。阐述了在某拓扑优化过程中以区间形式求解等效静力平衡方程的思想。通过顶点法得到了结构位移响应区间的上界和下界。

所述步骤四：根据形状导数的概念推导步骤一所建立的拓扑优化模型关于水平集演化伪时间

的灵敏度信息，在获得优化模型灵敏度信息的基础上可以建立合适的演化速度场，推动水平集函数的演化使得满足约束条件的同时获得目标函数的极小值。优化模型中目标函数的表达式：

其中：

目标函数由两部分构成，第一项是结构在时变载荷和不确定性影响下的平均动柔度区间的中心值，第二项是结构在时变载荷和不确定性影响下的平均动柔度区间的半径，结构的平均动柔度由全时段柔度积分转化为基于有限元动力学时间离散的各个离散时间点处的结构总柔度求和构成，其中离散时间点处的结构总柔度由有限元动力学空间离散的各个单元的柔度求和代替积分，这样将积分转化为求和的计算方法提高了计算效率。其中单元柔度的表达式为：

目标函数关于伪时间的灵敏度为：

目标函数的灵敏度分为两部分，一部分是结构在时变载荷和不确定性影响下的平均动柔度区间的中心值的灵敏度，第二部分是结构在时变载荷和不确定性影响下的平均动柔度区间的半径的灵敏度。第一部分可以详写为：

依据形状导数的原理：

将上式代入第一部分的灵敏度表达式有：

这样就推导了结构在时变载荷和不确定性影响下的平均动柔度区间的中心值的灵敏度，同样的，结构在时变载荷和不确定性影响下的平均动柔度区间的半径的灵敏度表达为：

同样的基于形状导数原理可以推导材料占比约束条件的灵敏度：

到目前为止，已经推导出了动力学鲁棒性拓扑优化模型对伪时间的灵敏度。在每次等效静力学拓扑优化过程中，根据已经推导的优化模型的灵敏度信息，构造合适的速度场，然后经历水平集函数随着伪时间的演化，使目标函数向下并满足约束条件。

所述步骤五：基于梯度算法求解优化模型相比于遗传算法等搜索算法具有更高的效率，本发明所采用的是优化准则法，优化准则法适用于设计变量多，约束条件少的优化问题。本发明所建立的动力学鲁棒性拓扑优化模型设计变量数量与有限元空间离散的结点数量成线性正相关，约束条件只有材料占比约束，这适合采用优化准则法求解。基于优化准则法，将带约束的优化模型转化为无约束优化模型：

其中，

无约束优化目标函数，λ为约束条件的权重系数。为了使得无约束目标函数随着水平集函数的演化不断下降，需要求解无约束优化目标函数关于演化伪时间的灵敏度信息，基于灵敏度信息求解方法，可以求得以下结果：

从无约束优化目标函数关于演化伪时间的灵敏度表达式可以建立法向演化速度场V_n：

如此建立的速度场可以使得无约束目标函数随着水平集函数的演化不断下降，其中参数λ在多次拓扑优化中的迭代公式为：

其中k为迭代次数，μ为迭代参数λ的权重系数，n_c为加速参数，可以使得拓扑优化在前n_c次迭代中对材料占比约束条件有一个加速过程，更快的满足约束条件，参数η的迭代公式为：

η^k+1＝min(η^k+Δη,η_max)k＞n_c

其中Δη为参数变化步长，η_max为参数上界。

所述步骤六：水平集函数的演化是通过H-J偏微分方程的求解来驱动的，利用MQ样条基函数对水平集函数插值表示：

其中，g_i(x)为第i个控制点处的MQ样条函数，

为保证计算稳定性的一阶多项式，将MQ样条基函数插值表示的水平集函数代入演化方程，上面的H-J偏微分方程将转化为常微分方程：

其中：

该耦合非线性常微分方程可以用一阶正演欧拉法求解，近似解为：

这是拓扑优化设计变量的更新方程，通过求解该代数方程可以实现设计变量的更新，进而实现水平集函数的更新，代表结构构型的改变。其中i为拓扑优化的循环次数，

表示水平集函数演化的时间步长，G是控制点处的MQ基函数数值矩阵，B是控制点处演化速度数值矩阵，

通过求解设计变量更新方程，同时也实现了水平集函数随着伪时间的不断演化，代表着结构的不断改变，直到迭代收敛。

所述步骤七：拓扑优化的收敛条件为：

这表示当优化模型中的目标函数连续6次都满足波动程度小于0.2％并且材料占比约束条件连续6次都满足波动程度小于0.2％。从拓扑优化的收敛条件看出，收敛的思路是获得稳定的满足约束条件的同时地目标函数的极小值。当拓扑优化迭代过程中满足收敛条件时，迭代终止，输出结构最终的拓扑优化设计；当拓扑优化迭代不满足收敛条件，则重复步骤三到步骤六，继续驱动结构改变指导迭代终止；

本发明的一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法，其高效性体现在：所述步骤一中使用MQ样条函数对水平集函数进行插值，将H-J偏微分方程的求解转化为常微分方程的求解，提高了计算效率；所述步骤三中使用等效静态载荷法将动力学平衡状态方程转化为多个等效静态平衡状态方程，提高了计算效率；所述步骤二、三中不确定性量化和不确定性传播方法的顶点法的适用；所述步骤五中，采用优化准则法求解所建立的动力学鲁棒性拓扑优化模型。

本发明实现了顶点法、参数化水平集法和等效静态载荷法的结合，其中创新之处在于基于顶点法和参数化水平集法的区间形式的等效静态平衡方程的求解，主要贡献在于在动力学鲁棒性拓扑优化的多个阶段提高了计算效率，最终实现了一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法。

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

本发明提供了一种针对连续体结构的高效的动力学鲁棒性拓扑优化设计的新思路，在连续体结构的拓扑优化设计中考虑了一定参数的不确定性和服役环境的时变特性，使得拓扑优化设计更加合理更加贴近现实。本文采用参数化水平集法模拟拓扑优化设计中结构的构型的改变，采用等效静态载荷法处理多次求解动力学平衡方程的问题，采用基于区间模型的顶点法表征和分析不确定性的影响，采用优化准则法求解优化模型，建立考虑材料弹性模量、外部载荷存在一定范围的不确定性的动力学鲁棒性拓扑优化方法。并且在动力学鲁棒性拓扑优化的不确定性表征分析、动力学平衡方程求解、拓扑优化中结构模拟和优化模型求解阶段采用高效手段实现针对连续体结构的动力学鲁棒性拓扑优化设计，在结构的安全性、经济性两方面得以兼顾和平衡。

附图说明

图1是本发明提出的动力学鲁棒性拓扑优化整体流程；

图2是MQ样条函数及其X方向的梯度；(a)MQ样条函数(b)MQ样条函数在X方向的梯度；

图3是实施例的设计域、边界条件和载荷；

图4是实施例在不同拓扑优化策略下的最终设计结果；(a)静力学确定性(V_targ＝0.4)；(b)静力学鲁棒性(V_targ＝0.4)；(c)动力学确定性(V_targ＝0.4)；(d)动力学鲁棒性(V_targ＝0.4)；(e)动力学鲁棒性(V_targ＝0.35)；(f)动力学鲁棒性(V_targ＝0.45)；

图5是实施例在不同拓扑优化策略下的迭代历程；

(a):静力学确定性拓扑优化，材料占比约束0.4

(b):静力学鲁棒性拓扑优化，材料占比约束0.4

(c):动力学确定性拓扑优化，材料占比约束0.4

(d):动力学鲁棒性拓扑优化，材料占比约束0.4

(e):动力学鲁棒性拓扑优化，材料占比约束0.35

(f):动力学鲁棒性拓扑优化，材料占比约束0.45。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整的描述，显然，所描述的实施例仅为本发明的一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域的普通技术人员在不付出创造性劳动的前提下所获得的所有其他实施例，都属于本发明的保护范围。

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

步骤一：给定平面连续体结构、时变载荷以及约束条件和设计域与非设计域，基于一般的拓扑优化数学模型，以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，建立考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型；

步骤二：考虑材料弹性模量，外部载荷存在一定范围的不确定性，采用区间变量来表示整体刚度区间矩阵、整体阻尼区间矩阵、整体质量区间矩阵和载荷区间向量，采用区间向量来表示位移、速度和加速度区间向量，根据有限元控制方程的区间形式，使用顶点法，求出初始结构的时变载荷影响下的位移、速度和加速度在有界不确定参数影响下的上下界；

步骤三：根据等效静态载荷法，计算各个离散时间点处的结构的等效静态载荷，建立与动力学平衡方程等效的多个离散时间点处的静态平衡方程，将拓扑优化中动力学平衡的状态约束条件转化为多个离散时间点处的等效静力学平衡约束。根据顶点法，求解考虑不确定性参数影响的区间形式等效静态平衡方程，获得结构在离散时间处的各个结点位移区间。

步骤四：依据形状导数的思想计算优化模型中目标函数关于水平集函数演化伪时间的灵敏度和优化模型中约束条件关于演化伪时间的灵敏度；

步骤五：采用优化准则法，求解以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，所建立的考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型，构造水平集函数的演化速度场；

步骤六：基于演化速度场，求解水平集函数的演化方程，通过更新设计变量以更新水平集函数模拟结构构型的不断改变；

步骤七：如果当前设计满足材料占比约束条件并且满足迭代稳定性条件，则迭代结束，将当前拓扑优化的结果作为最终的优化结果。如果不满足材料占比约束条件或者迭代不稳定，则重复步骤三至步骤六。

所述步骤一：基于一般的拓扑优化模型，建立基于参数化水平集法的动力学鲁棒性拓扑优化的数学模型：

这表示优化模型以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量。针对线弹性结构以材料占比为约束条件，通过MQ样条函数对水平集函数插值进而以MQ样条函数的系数作为设计变量，考虑一般力学问题中载荷的时变特性以及工程中广泛存在的多源不确定性对连续体结构进行拓扑优化设计。

其中，

表示第i个控制点处的MQ样条基函数的系数，下标N表示控制点的个数，

是为了保证计算稳定性而附加的一阶多项式系数，x是控制设计域内某点的位置变量

是优化的目标函数，

是多源不确定性影响下的结构在时变载荷作用下的平均动柔度的区间中心值，

是多源不确定性影响下的结构在时变载荷作用下的平均动柔度的区间半径，ξ₁,ξ₂分别代表平均动柔度区间半径和区间中心值的权重，可以通过调整两个参数的数值大小以反映对不同部分的关注程度，

是随着伪时间

不断演化的水平集函数，上标I表示该变量为区间形式，上标c表示区间的中心值，上标r表示区间的半径，b^I表示有界但不确定性的不确定性变量，u表示结构的位移响应，K表示结构基于有限元离散的整体刚度矩阵，F_eq表示基于等效静态载荷法的离散时间t_i处的结构等效静态载荷列向量，

表示离散时间t_i处的弹性结构平衡状态方程；D表示全部的结构设计域，

是描述结构材料占比的约束的函数，V_targ表示拓扑优化设计的材料占比约束条件。

表示水平集函数

的阶跃函数，形式为：

水平集法是将结构空间Rⁿ中的边界隐式表示到高一维空间Rⁿ⁺¹的水平集函数中去，将比结构空间高一维空间Rⁿ⁺¹中水平集函数的零水平集作为结构空间Rⁿ的边界描述，对二维问题而言，界面曲线

表示为一个三维标量函数

的零水平集：

其中x表示二维空间中某点的坐标(x,y)，Ω表示结构设计域内存在材料的区域，

表示结构设计域内不存在材料的区域，

表示结构设计域内的边界。水平集函数演化方程为：

这是一个哈密顿-雅克比(H-J：Hamilton-Jacobi)偏微分方程，V_n表示水平集函数零水平集的法向演化速度，

表示水平集函数梯度，通过构造合适的法向演化速度场，可以驱动水平集函数的演化，使得结构边界向特定的方向演化，以满足拓扑优化设计的需要。该H-J偏微分方程的求解需要满足一系列条件不利于求解，通过引入MQ样条函数对水平集函数进行插值，可以将H-J偏微分方程转化为常微分方程求解，提高计算效率。利用MQ样条基函数对水平集函数插值表示：

其中，g_i(x)为第i个控制点处的MQ样条函数，MQ样条函数以及其在X方向的梯度如图2(a)、(b)所示，

为保证计算稳定性的一阶多项式，将MQ样条基函数插值表示的水平集函数代入的水平集函数的演化方程，以上面的H-J偏微分方程将转化为常微分方程：

其中：

所述步骤二：针对某个确定的水平集函数，考虑材料弹性模量，外部载荷存在一定范围的不确定性，采用区间变量M^I，K^I，C^I和F^I来表示整体刚度区间矩阵,整体阻尼区间矩阵和载荷区间向量，采用区间向量

来表示离散时间t_i处的结点位移区间向量，根据有限元的位移控制方程有：

其中

为结构结点的加速度区间向量，

为结构结点的速度区间向量，基于不确定性分析方法中的顶点法以及动力学响应分析方法中的Newmark-β算法，首先计算加速度

在有界但不确定参数影响下的上下界：

其中

分别代表遍历Nd个不确定性参数上下界的2^Nd次动力学响应分析中的第j个结点处的加速度最小值和最大值，以此确定加速度区间向量。基于Newmark-β算法和顶点法可以确定性速度和位移的区间向量。

所述步骤三：根据等效静态载荷法，有以下等式成立：

其中F_eq由结构的外部载荷、阻尼载荷、惯性载荷三部分组成，M,C,K,F分别为确定性条件下的整体质量、阻尼、考虑诸多不确定性的影响，基于等效静态载荷法，将拓扑优化中动力学平衡方程的求解转化为了有限个静力学平衡方程的求解，这样的转化提高了计算的效率。上述方程将变为区间形式的平衡方程：

整体刚度区间矩阵K^I直接由不确定性参数b^I影响，由等效静态载荷法，各个离散时间t_i处的等效载荷区间向量具有以下表达式：

其中阻尼采用瑞利阻尼

α_m,α_k代表刚度和质量对阻尼的影响系数，

和

分别代表初始结构的加速度和速度区间向量，两者可以将初始结构的水平集函数代入步骤二求解获得初始结构在时变载荷影响下的响应区间向量。整体阻尼区间矩阵的上下界受整体质量区间矩阵和整体刚度区间矩阵的影响：

C^U＝α_mM^U(b^I)+α_kK^U(b^I)

C^L＝α_mM^L(b^I)+α_kK^L(b^I)

进而可以确定各个离散时间处的等效静态载荷：

整体位移区间向量表示为：

其中：

通过以上分析，得出如下结论:1).每次拓扑优化的计算代价与参数不确定性参数Nd数量呈指数关系。2).每次拓扑优化的计算成本与单元数N和离散次数Nt之间存在乘积关系。因此，采用简单的不确定性分析方法来降低计算成本是非常重要的，这也是采用顶点法的原因。阐述了在某拓扑优化过程中以区间形式求解等效静力平衡方程的思想。通过顶点法得到了结构位移响应区间的上界和下界。

所述步骤四：根据形状导数的概念推导步骤一所建立的拓扑优化模型关于水平集演化伪时间

的灵敏度信息，在获得优化模型灵敏度信息的基础上可以建立合适的演化速度场，推动水平集函数的演化使得满足约束条件的同时获得目标函数的极小值。优化模型中目标函数的表达式：

其中：

目标函数由两部分构成，第一项是结构在时变载荷和不确定性影响下的平均动柔度区间的中心值，第二项是结构在时变载荷和不确定性影响下的平均动柔度区间的半径，结构的平均动柔度由全时段柔度积分转化为基于有限元动力学时间离散的各个离散时间点处的结构总柔度求和构成，其中离散时间点处的结构总柔度由有限元动力学空间离散的各个单元的柔度求和代替积分，这样将积分转化为求和的计算方法提高了计算效率。其中单元柔度的表达式为：

目标函数关于伪时间的灵敏度为：

目标函数的灵敏度分为两部分，一部分是结构在时变载荷和不确定性影响下的平均动柔度区间的中心值的灵敏度，第二部分是结构在时变载荷和不确定性影响下的平均动柔度区间的半径的灵敏度。第一部分可以详写为：

依据形状导数的原理：

将上式代入第一部分的灵敏度表达式有：

这样就推导了结构在时变载荷和不确定性影响下的平均动柔度区间的中心值的灵敏度，同样的，结构在时变载荷和不确定性影响下的平均动柔度区间的半径的灵敏度表达为：

同样的基于形状导数原理可以推导材料占比约束条件的灵敏度：

到目前为止，已经推导出了动力学鲁棒性拓扑优化模型对伪时间的灵敏度。在每次等效静力学拓扑优化过程中，根据已经推导的优化模型的灵敏度信息，构造合适的速度场，然后经历水平集函数随着伪时间的演化，使目标函数向下并满足约束条件。

所述步骤五：基于梯度算法求解优化模型相比于遗传算法等搜索算法具有更高的效率，本专利所采用的是优化准则法，优化准则法适用于设计变量多，约束条件少的优化问题。本专利所建立的动力学鲁棒性拓扑优化模型设计变量数量与有限元空间离散的结点数量成线性正相关，约束条件只有材料占比约束，这适合采用优化准则法求解。基于优化准则法，将带约束的优化模型转化为无约束优化模型：

其中，

无约束优化目标函数，λ为约束条件的权重系数。为了使得无约束目标函数随着水平集函数的演化不断下降，需要求解无约束优化目标函数关于演化伪时间的灵敏度信息，基于权利5要求的灵敏度信息求解方法，可以求得以下结果：

从无约束优化目标函数关于演化伪时间的灵敏度表达式可以建立法向演化速度场V_n：

如此建立的速度场可以使得无约束目标函数随着水平集函数的演化不断下降，其中参数λ在多次拓扑优化中的迭代公式为：

其中k为迭代次数，μ为迭代参数λ的权重系数，n_c为加速参数，可以使得拓扑优化在前n_c次迭代中对材料占比约束条件有一个加速过程，更快的满足约束条件，参数η的迭代公式为：

η^k+1＝min(η^k+Δη,η_max)k＞n_c

其中Δη为参数变化步长，η_max为参数上界。

所述步骤六：水平集函数的演化是通过H-J偏微分方程的求解来驱动的，利用MQ样条基函数对水平集函数插值表示：

其中，g_i(x)为第i个控制点处的MQ样条函数，

为保证计算稳定性的一阶多项式，将MQ样条基函数插值表示的水平集函数代入演化方程，上面的H-J偏微分方程将转化为常微分方程：

其中：

该耦合非线性常微分方程可以用一阶正演欧拉法求解，近似解为：

这是拓扑优化设计变量的更新方程，通过求解该代数方程可以实现设计变量的更新，进而实现水平集函数的更新，代表结构构型的改变。其中i为拓扑优化的循环次数，

表示水平集函数演化的时间步长，G是控制点处的MQ基函数数值矩阵，B是控制点处演化速度数值矩阵，

通过求解设计变量更新方程，同时也实现了水平集函数随着伪时间的不断演化，代表着结构的不断改变，直到迭代收敛。

所述步骤七：拓扑优化的收敛条件为：

这表示当优化模型中的目标函数连续6次都满足波动程度小于0.2％并且材料占比约束条件连续6次都满足波动程度小于0.2％。从拓扑优化的收敛条件看出，收敛的思路是获得稳定的满足约束条件的同时地目标函数的极小值。当拓扑优化迭代过程中满足收敛条件时，迭代终止，输出结构最终的拓扑优化设计；当拓扑优化迭代不满足收敛条件，则重复步骤三到步骤六，继续驱动结构改变指导迭代终止；

本发明的一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法，其高效性体现在：所述步骤一中使用MQ样条函数对水平集函数进行插值，将H-J偏微分方程的求解转化为常微分方程的求解，提高了计算效率；所述步骤三中使用等效静态载荷法将动力学平衡状态方程转化为多个等效静态平衡状态方程，提高了计算效率；所述步骤二、三中不确定性量化和不确定性传播方法的顶点法的适用；所述步骤五中，采用优化准则法求解所建立的动力学鲁棒性拓扑优化模型。

本发明实现了顶点法、参数化水平集法和等效静态载荷法的结合，其中创新之处在于基于顶点法和参数化水平集法的区间形式的等效静态平衡方程的求解，主要贡献在于在动力学鲁棒性拓扑优化的多个阶段提高了计算效率，最终实现了一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法。

实施例：

本发明所公开的动力学鲁棒性拓扑优化方法针对任意平面结构的动力学拓扑优化设计都适用，如图3所示的一般矩形平板的动力学拓扑优化问题，采用本发明所提出的一种高效的动力学鲁棒性拓扑优化方法可以说明动力学和鲁棒性拓扑优化设计的必要性。设计区域为0.6m×0.3m的矩形区域，划分为60×30个单元。材料弹性模量70GPa，泊松比μ＝0.3。矩形区域上边界的左半部分施加固支约束，下边界右顶点施加时变集中载荷，其中心值的满足F^c＝500sin(4πt)e^-3t+500，材料占比约束条件设置为0.35，0.4和0.45。设弹性模量E、载荷F相对名义值均有5％的波动。不同设计策略下拓扑优化的构型如图4(a)-(f)所示，不同设计策略下拓扑优化的迭代历程曲线如图5(a)-(f)所示。

表1实施例在不同拓扑优化策略下的结果

结果说明：

1)从实施例工况(a)/(c)和(b)/(d)的结果可以看出，动力学拓扑优化得到的最终设计结果与静力学拓扑优化得到的最终设计结果有很大的不同。此外，动力学拓扑优化策略下的平均动柔度区间的中心值和半径值都小于相同条件下静力学拓扑优化的结果。这说明有必要考虑结构服役环境的时变效应也就是说明动力学拓扑优化的必要性。

2)从实施例工况(a)/(b)和(c)/(d)也可以看出，鲁棒性拓扑优化得到的最终结果与确定性拓扑优化得到的最终结果有显著差异。鲁棒性拓扑优化策略下的平均动柔度区间的中心值和半径值都小于相同条件下确定性拓扑优化的结果。结果表明，与确定性拓扑优化相比，鲁棒性拓扑优化能有效降低不确定性对结构的影响，这说明有必要考虑结构服役环境的不确定性参数的影响也就是说明鲁棒性拓扑优化的必要性。

3)实施例工况(d)/(e)/(f)的结果说明了动力学鲁棒性拓扑优化方法针对不同材料占比约束下的有效性。当材料占比较大时，结构的平均动柔度较小，即结构在时变荷载作用下的响应较小。

本发明提供了一种针对连续体结构的高效的动力学鲁棒性拓扑优化设计的新思路，在连续体结构的动力学拓扑优化设计中考虑了一定参数的不确定性，使得拓扑优化设计更加合理更加贴近现实。本文采用参数化水平集法模拟拓扑优化设计中结构的构型的改变，采用等效静态载荷法处理多次求解动力学平衡方程的问题，采用基于区间模型的顶点法表征和分析不确定性的影响，采用优化准则法求解优化模型，建立考虑材料弹性模量、外部载荷存在一定范围的不确定性的动力学鲁棒性拓扑优化方法，在动力学鲁棒性拓扑优化的不确定性表征分析、动力学平衡方程求解、拓扑优化中结构模拟和优化模型求解阶段采用高效手段实现针对连续体结构的动力学鲁棒性拓扑优化设计，在结构的安全性、经济性两方面得以兼顾和平衡。

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

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

Claims (1)

1.一种针对连续体结构的高效的动力学鲁棒性拓扑优化方法，用于连续体结构的拓扑优化设计中，实现不确定性影响下的平面连续体结构受时变载荷作用以及固支、简支边界约束的鲁棒性拓扑优化设计，其特征在于：在拓扑优化问题中基于参数化水平集法和顶点法求解区间形式的等效静态平衡方程，并且考虑了材料的弹性模量和载荷数值存在的有界不确定性，实现步骤如下：

步骤一：给定平面连续体结构、时变载荷以及约束条件和设计域与非设计域，基于拓扑优化数学模型，以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，建立考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型；

步骤二：考虑材料弹性模量，外部载荷存在的不确定性，采用区间变量来表示整体刚度区间矩阵、整体阻尼区间矩阵、整体质量区间矩阵和载荷区间向量，采用区间向量来表示位移、速度和加速度区间向量，根据有限元控制方程的区间形式，使用顶点法，求出初始结构的位移、速度和加速度在有界不确定参数影响下的上下界；

步骤三：根据等效静态载荷法，计算各个离散时间点处的结构的等效静态载荷，建立与动力学平衡方程等效的多个离散时间点处的静态平衡方程，将拓扑优化中动力学平衡的状态约束条件转化为多个离散时间点处的等效静力学平衡约束；根据顶点法，求解考虑不确定性参数影响的区间形式等效静态平衡方程，获得结构在离散时间处的各个结点位移区间；

步骤四：依据形状导数计算优化模型中以平均动柔度表示的目标函数关于水平集函数演化伪时间的灵敏度和优化模型中约束条件关于演化伪时间的灵敏度；

步骤五：采用优化准则法，求解以材料占比为约束条件，以不确定性参数影响下受载荷全时段作用的结构平均动柔度的区间中心值与区间半径为优化目标，以MQ样条函数的系数作为设计变量，所建立的考虑有界不确定性参数影响的动力学鲁棒性拓扑优化模型，构造水平集函数的演化速度场；

步骤六：基于演化速度场，求解水平集函数的演化方程，通过更新设计变量以更新水平集函数模拟结构构型的不断改变；

步骤七：如果当前设计满足材料占比约束条件并且满足迭代稳定性条件，则迭代结束，将当前拓扑优化的结果作为最终的优化结果；如果不满足材料占比约束条件或者迭代不稳定，则重复步骤三至步骤六。

Images