CN112446163A - 基于参数化水平集的能量有限元拓扑优化方法 - Google Patents

Abstract

基于参数化水平集的能量有限元拓扑优化方法，通过使用参数化水平集进行结构描述，从而减小变量数量，使得结构的描述不再依赖背景网格；采用能量有限元法对振动响应进行分析，计算结果在稀疏网格中仍具有较高的精度，进一步减小计算压力；本发明可解决大型高频振动结构优化设计中所遇到的计算精度不高、计算量大等问题，适应于大尺寸高频振动结构的优化设计问题。

Description

基于参数化水平集的能量有限元拓扑优化方法

技术领域

本发明属于结构的动态性能优化设计技术领域，具体涉及基于参数化水平集的能量有限元拓扑优化方法。

技术背景

结构的振动与噪声问题一直是机械设备设计制造中亟待解决的关键问题之一，设计者往往在机械设备的设计研究阶段，利用计算模拟方法对结构系统在外界激励下的振动响应情况进行分析研究，从而对机械结构的设计进行优化改进，提高机械结构的性能；经典有限元方法是研究结构振动问题时应用最广泛的一种传统结构动力学分析方法，该方法需要将结构划分为若干网格，从而计算网格节点上的动态位移；为了保证网格能够捕捉到结构系统在振动时的变形情况，网格尺寸应小于波长的六分之一，因此在研究高频振动时需要加密网格或提高单元插值阶次；同时传统的拓扑优化方法例如变密度法(SIMP)、进化算法(ESO)、RAMP法等基于像素描述，这使得在研究高频段大型结构的振动问题时需要大量的网格，计算量极大，会因计算机计算能力等的限制而受到制约；对于高频范围内大尺寸三维结构系统的动力学分析，由于网格尺寸的要求与分析质量密切相关，采用经典有限元法和传统拓扑优化方法求解结构的运动响应并不合适。

发明内容

为了克服上述现有技术的缺点，本发明的目的在于提供基于参数化水平集的能量有限元拓扑优化方法，使用参数化水平集方法(MMC)描述拓扑结构并投影到能量有限元模型，采用能量有限元方法分析机械结构系统的振动响应情况，在减少计算网格数量的同时仍可以获得较为准确的振动分析结果，降低了运算量。

为了达到上述目标，本发明采取的技术方案为：

基于参数化水平集的能量有限元拓扑优化方法，包括以下步骤：

1)定义设计工况：

设计域的结构为矩形加筋板，加筋板由基板和加强筋组成，基板尺寸参数为l×w×h_plate，基板由厚度为h_stiffener的加强筋加强，定义加筋板的四边固定，加筋板上表面受均匀分布的单位动态载荷；

2)确定设计变量：

用组件来描述加强筋，取矩形作为组件的形状，每个组件包含矩形中心坐标(x₀,y₀)、长度L、宽度T以及倾斜角θ共计5个变量，在基板上均匀布置n个组件，将其作为初始布局，此时共有5n个变量，将这些变量有序的储存到向量v中；

3)基于参数化水平集的结构几何描述：

3.1)n个组件均匀分布基板上表面，构成初始结构；

3.2)使用几何参数表示第q个组件的水平集函数φ_q，q＝1,2…n:

3.3)使用Heaviside公式对每个组件的水平集函数φ_q标准化：

3.4)组装所有组件水平集：

3.5)使用Heaviside公式对整个结构进行第二次标准化，得到整个结构的最终水平集φ_sum：

4)将结构的几何描述投影到能量有限元模型:

每个单元的厚度取决于水平集的覆盖范围：

上式中，h^j为第j个单元的厚度，h_stiffener和h_plate分别为加强筋和基板的厚度，

为第j个单元的第i个节点的水平集值，n_nod为每个单元的总节点数；

5)能量有限元计算：

5.1)构建单元能量矩阵：

每个单元中弹性波能量平衡的控制微分方程为：

上式中，e为能量密度，c_g为弹性波群速度，η为阻尼系数，ω为角频率，∏为输入功率；

控制方程的弱形式为：

上式中，

为加强筋边界的法向量，N为形函数，{e^j}为第j个单元的节点能量密度向量；

弱形式控制方程的矩阵形式为：

[K^j]{e^j}＝{F^j}+{Q^j} (8)

其中：

上式中，K^j为单元能量矩阵，F^j为单元输入功率，∏(x,y)为节点(x,y)上的输入功率，Q^j为单元边界上的能流；

5.2)能量的反射和折射：

能流在厚度不同的相邻单元边界上发生反射和折射，能量传递系数τ₁₂和反射系数r₁₁分别为：

上式中，θ₁为入射角，θ₂为折射角，A_f1、C_f1、D_f1、A_f2、B_f2分别为入射波、反射波、近场消散波、折射波、消散波的振幅；

5.3)耦合单元分析：

5.3.1)不同厚度的单元间耦合：

不同厚度的相邻单元边界增加新节点，此时能量有限元表达式为：

上式中，K为未耦合的全局能量矩阵，K_q为相邻单元的耦合矩阵；

5.3.2)相同厚度的单元间耦合：

所有单元边界增加新节点得到能量有限元新网格——非耦合网格，非耦合网格中全局能量矩阵装配步骤如下：

5.3.2.1)非耦合网格中无公共节点耦合的全局能量矩阵:

5.3.2.2)耦合水平边界单元节点：定义i1、i2、i3、i4为一个单元的四个节点，定义K(i,...)表示K(i1,i1)、K(i1,i2)、K(i1,i3)、K(i1,i4)，当i1节点和i2节点需要耦合时，首先将K(i1,...)添加到第i2行，再将K(i2,...)添加到第i1行；

5.3.2.3)按照上述规则耦合垂直边界单元节点；

5.3.2.4)按照上述规则耦合对角线边界单元节点；

5.3.2.5)将组装所得全局能量矩阵中相互耦合的节点i_coup1、节点i_coup2进行操作：第i_coup2行的K_{icoup2,icoup2}直接移动至第i_coup1列，且第i_coup1列保持不变；得到全局能量矩阵的最终形式：

5.4)确定优化模型：

以结构的动态性能最佳为优化目标，将结构所储存的能量设定为衡量结构动态性能的指标，命名为能量柔度，能量柔度即为关于设计变量的目标函数；设定优化结构材料用量不得超过设计许用材料用量，将材料用量作为约束函数；

优化数学模型如下：

上式中，变量v为所有组件的几何参数，J(v)为目标函数，M(v)和M_upp分别为加强筋材料用量和最大设计许用材料用量，n_e为总单元数，h^j和S^j分别为第j个单元的厚度和面积；

5.5)敏度分析：

5.5.1)目标函数敏度：

目标函数的敏度计算公式如下：

上式中，

为全局能量矩阵的敏度，

为耦合矩阵的敏度；

5.5.2)约束函数敏度：

约束函数的敏度计算公式如下：

6)迭代优化：

将能量有限元计算结果以及敏度带入移动渐近线优化算法(MMA)中，迭代更新变量，直至目标函数在满足约束条件的情况下收敛为止，此时获得满足材料用量约束条件下的加筋板的最优结构布局；

7)适应性处理：

按照生产工艺要求圆整加筋板的结构布局，从而获得加筋板的结构最终布局。

为适应不同设计需求，使用时并不局限于所述的约束及优化目标，设计者能够加入质量评价、强度评价、刚度评价、疲劳寿命评价，评价方法通过能量有限元计算获得。

本发明的有益效果为：

由于本发明方法基于参数化水平集来对拓扑结构进行描述，与传统基于像素的结构描述相比，其优化中所需的变量数量大大减小，且结构的描述不再依赖背景网格，在稀疏网格中应用拓扑优化时更符合能量有限元的要求；这使得设计的灵活性不再取决于网格的分辨率，而是更依赖于几何参数和组件的数量，能够在稀疏网格内实现大尺寸高振动结构的拓扑优化；因为本方法采用能量有限元法对振动响应进行分析，所以在稀疏网格中仍具有较高的精度，其计算结果可以为一般设计者所接受，这大大减小了计算压力，可以解决大型高频振动结构优化设计中所遇到的计算精度不高、计算量大等问题；从传统设计角度来看，本发明方法可以进一步扩展，增加质量、强度、刚度、疲劳寿命、等设计目标或约束，这些扩展帮助本方法适应传统设计需求，两者互为表里、相得益彰，对大尺寸高频振动结构的优化设计而言意义非常。

附图说明

图1为本发明实施例的边界条件示意图。

图2为本发明实施例中组件初始布局图。

图3为本发明实施例中非耦合网格示意图。

图4为本发明实施例中满足约束条件的最终优化结果。

具体实施方式

下面结合附图和实施例对本发明作进一步说明，实施例采用大尺寸矩形加筋板。

基于参数化水平集的能量有限元拓扑优化方法，包括以下步骤：

1)定义设计工况：

本实施例设计域的结构为大尺寸矩形加筋板，加筋板由基板和加强筋组成，基板尺寸参数为6m×6m×0.05m，基板由厚度为0.05m的加强筋加强，定义加筋板的四边固定，加筋板上表面受均匀分布的动态载荷Fe^jωt,其中F＝6400N，角频率ω＝1000，如图1所示；

2)确定设计变量：

用组件来描述加强筋，本实施例中取矩形作为组件的形状，每个组件包含矩形中心坐标(x₀,y₀)、长度L、宽度T以及倾斜角θ共计5个变量，本实施例在基板上均匀布置32个加强筋，如图2所示，将其作为初始布局，此时共有160个变量，将这些变量有序的储存到向量v中；

3)基于参数化水平集的结构几何描述：

3.1)本实施例将32个组件均匀分布基板上表面，构成初始结构；

3.2)使用几何参数表示第q个组件的水平集函数φ_q，q＝1,2…n：

3.3)使用Heaviside公式对每个组件的水平集函数φ_q标准化：

3.4)组装所有组件水平集：

3.5)本实施例使用Heaviside公式对整个结构进行第二次标准化，其中Heaviside公式的大数b取10000，得到整个结构的最终水平集φ_sum：

4)将结构的几何描述投影到能量有限元模型:

每个单元的厚度取决于水平集的覆盖范围：

上式中，h^j为第j个单元的厚度，h_stiffener和h_plate分别为加强筋和基板的厚度，本实施例中两者取相同值0.05m，

为第j个单元的第i个节点的水平集值，n_nod为每个单元的总节点数，本实施例中为4；

5)能量有限元计算：

5.1)构建单元能量矩阵：

本实施例中单元的弹性波能量平衡的控制微分方程为：

上式中，e为能量密度，c_g为弹性波群速度，η为阻尼系数，ω为角频率，∏为输入功率；

控制方程的弱形式为：

上式中，

为加强筋边界的法向量，N为形函数，{e^j}为第j个单元的节点能量密度向量；

弱形式控制方程的矩阵形式为：

[K^j]{e^j}＝{F^j}+{Q^j} (8)

其中：

上式中，K^j为单元能量矩阵，F^j为单元输入功率，∏(x,y)为节点(x,y)上的输入功率，Q^j为单元边界上的能流；

5.2)能量的反射和折射：

能流在厚度不同的相邻单元边界上发生反射和折射，能量传递系数τ₁₂和反射系数r₁₁分别为：

上式中，θ₁为入射角，θ₂为折射角，A_f1、C_f1、D_f1、A_f2、B_f2分别为入射波、反射波、近场消散波、折射波、消散波的振幅；

5.3)耦合单元分析：

5.3.1)不同厚度的单元间耦合：

不同厚度的相邻单元边界增加新节点，此时能量有限元表达式为：

上式中，K为未耦合的全局能量矩阵，K_q为相邻单元的耦合矩阵；

5.3.2)相同厚度单元间耦合：

所有单元边界增加新节点得到能量有限元新网格——非耦合网格，如图3所示，非耦合网格中全局能量矩阵装配步骤如下：

5.3.2.1)非耦合网格中无公共节点耦合的全局能量矩阵：

5.3.2.2)耦合水平边界单元节点：定义i1、i2、i3、i4为一个单元的四个节点，定义K(i,...)表示K(i1,i1)、K(i1,i2)、K(i1,i3)、K(i1,i4)，当i1节点和i2节点需要耦合时，首先将K(i1,...)添加到第i2行，再将K(i2,...)添加到第i1行；

5.3.2.3)按照上述规则耦合垂直边界单元节点；

5.3.2.4)按照上述规则耦合对角线边界单元节点；

5.3.2.5)将组装所得全局能量矩阵中相互耦合的节点i_coup1、节点i_coup2进行操作：第i_coup2行的

直接移动至第i_coup1列，且第i_coup1列保持不变；得到全局能量矩阵的最终形式：

5.4)确定优化模型：

本实施例以大尺寸矩形加筋板的动态性能最佳为优化目标，将其所储存的能量设定为衡量结构动态性能的指标，命名为能量柔度，能量柔度即为关于设计变量的目标函数；本实施例的约束条件是加强筋材料用量不得超过基板材料用量的43.5％；

优化数学模型如下：

本实施例中，变量v为所有组件的几何参数，J(x)为目标函数，M(x)和M_upp分别为加强筋材料用量和最大设计许用材料用量，总单元数n_e为1600，h^j和S^j分别为第j个单元的厚度和面积；

5.5)敏度分析：

5.5.1)目标函数敏度：

本实施例的目标函数的敏度计算公式如下：

上式中，

为全局能量矩阵的敏度，

为耦合矩阵的敏度；

5.5.2)约束函数敏度：

本实施例的约束函数的敏度计算公式如下：

6)迭代优化：

将步骤5)能量有限元计算得到的目标函数值、约束函数值以及它们的敏度带入优化算法——移动渐近线法中，迭代更新变量，直至目标函数在满足约束条件的情况下收敛为止，此时获得满足材料用量约束条件下的加筋板的最优结构布局；

7)适应性处理：

按照生产工艺要求圆整加筋板的结构布局，从而获得加筋板的结构最终布局，如图4所示。

为适应不同设计需求，使用时并不局限于所述的约束及优化目标，设计者能够加入质量评价、强度评价、刚度评价、疲劳寿命评价等等；本方法旨在提供基于参数化水平集的能量有限元拓扑优化的设计思路，其它评价方法可通过能量有限元计算获得。

Claims (2)

1.基于参数化水平集的能量有限元拓扑优化方法，其特征在于，包括以下步骤：

1)定义设计工况：

设计域的结构为矩形加筋板，加筋板由基板和加强筋组成，基板尺寸参数为l×w×h_plate，基板由厚度为h_stiffener的加强筋加强，定义加筋板的四边固定，加筋板上表面受均匀分布的单位动态载荷；

2)确定设计变量：

用组件来描述加筋板，取矩形作为组件的形状，每个组件包含矩形中心坐标(x₀,y₀)、长度L、宽度T以及倾斜角θ共计5个变量，在基板上均匀布置n个组件，将其作为初始布局，此时共有5n个变量，将这些变量有序的储存到向量v中；

3)基于参数化水平集的结构几何描述：

3.1)n个组件均匀分布基板上表面，构成初始结构；

3.2)使用几何参数表示第q个组件的水平集函数φ_q，q＝1,2…n：

3.3)使用Heaviside公式对每个组件的水平集函数φ_q标准化：

3.4)组装所有组件水平集：

3.5)使用Heaviside公式对整个结构进行第二次标准化，得到整个结构的最终水平集φ_sum：

4)将结构的几何描述投影到能量有限元模型:

每个单元的厚度取决于水平集的覆盖范围：

上式中，h^j为第j个单元的厚度，h_stiffener和h_plate分别为加强筋和基板的厚度，

为第j个单元的第i个节点的水平集值，n_nod为每个单元的总节点数；

5)能量有限元计算：

5.1)构建单元能量矩阵：

单元中弹性波能量平衡的控制微分方程为：

上式中，e为能量密度，c_g为弹性波群速度，η为阻尼系数，ω为角频率，Π为输入功率；

控制方程的弱形式为：

上式中，

为加强筋边界的法向量，N为形函数，{e^j}为第j个单元的节点能量密度向量；

弱形式控制方程的矩阵形式为：

[K^j]{e^j}＝{F^j}+{Q^j} (8)

其中：

上式中，K^j为单元能量矩阵，F^j为单元输入功率，∏(x,y)为节点(x,y)上的输入功率，Q^j为单元边界上的能流；

5.2)能量的反射和折射：

能流在厚度不同的相邻单元边界上发生反射和折射，能量传递系数τ₁₂和反射系数r₁₁分别为：

上式中，θ₁为入射角，θ₂为折射角，A_f1、C_f1、D_f1、A_f2、B_f2分别为入射波、反射波、近场消散波、折射波、消散波的振幅；

5.3)耦合单元分析：

5.3.1)不同厚度单元间耦合：

不同厚度的相邻单元边界增加新节点，此时能量有限元表达式为：

上式中，K为未耦合的全局能量矩阵，K_q为相邻单元的耦合矩阵；

5.3.2)相同厚度单元间耦合：

所有单元边界增加新节点得到能量有限元新网格——非耦合网格，非耦合网格中全局能量矩阵装配步骤如下：

5.3.2.1)非耦合网格中无公共节点耦合的全局能量矩阵:

5.3.2.2)耦合水平边界单元节点：定义i1、i2、i3、i4为一个单元的四个节点，定义K(i,...)表示K(i1,i1)、K(i1,i2)、K(i1,i3)、K(i1,i4)，当i1节点和i2节点需要耦合时，首先将K(i1,...)添加到第i2行，再将K(i2,...)添加到第i1行；

5.3.2.3)按照上述规则耦合垂直边界单元节点；

5.3.2.4)按照上述规则耦合对角线边界单元节点；

5.3.2.5)将组装所得全局能量矩阵中相互耦合的节点i_coup1、节点i_coup2进行操作：第i_coup2行的

直接移动至第i_coup1列，且第i_coup1列保持不变；得到全局能量矩阵的最终形式：

5.4)确定优化模型：

以结构的动态性能最佳为优化目标，将结构所储存的能量设定为衡量结构动态性能的指标，命名为能量柔度，能量柔度即为关于设计变量的目标函数；设定优化结构材料用量不得超过设计许用材料用量，将材料用量作为约束函数；

优化数学模型如下：

上式中，变量v为所有组件的几何参数，J(v)为目标函数，M(v)和M_upp分别为加强筋材料用量和最大设计许用材料用量，n_e为总单元数，h^j和S^j分别为第j个单元的厚度和面积；

5.5)敏度分析：

5.5.1)目标函数敏度：

目标函数的敏度计算公式如下：

上式中，

为全局能量矩阵的敏度，

为耦合矩阵的敏度；

5.5.2)约束函数敏度：

约束函数的敏度计算公式如下：

6)迭代优化：

将能量有限元计算结果以及敏度带入移动渐近线优化算法(MMA)中，迭代更新变量，直至目标函数在满足约束条件的情况下收敛为止，此时获得满足材料用量约束条件下加筋板的最优结构布局；

7)适应性处理：

按照生产工艺要求圆整加筋板的结构布局，从而获得加筋板的结构最终布局。

2.根据权利要求1所述的基于参数化水平集的能量有限元拓扑优化方法，其特征在于，为适应不同设计需求，使用时并不局限于所述的约束及优化目标，设计者能够加入质量评价、强度评价、刚度评价、疲劳寿命评价，评价方法通过能量有限元计算获得。

Images