CN116205106A - 一种基于simp法的边界光滑处理方法 - Google Patents

Abstract

本发明涉及优化算法领域，且公开了一种基于SIMP法的边界光滑处理方法，步骤一：定义边界条件、材料参数等并确定目标函数，步骤二：SIMP法结构模型初始化，步骤三：有限元分析并计算目标函数、敏度等，步骤四：使用合适的方法进行设计变量更新，步骤五：对设计变量进行敏度过滤，步骤六：引入高一维度的水平集曲面，并寻找合适的与水平集曲面相交的等值面，并投影在xoy面，步骤七：检查是否满足优化条件，若满足则输出结果，若不满足则返回步骤三，该基于SIMP法的边界光滑处理方法，经过该方法后处理的拓扑构型不再有灰度单元，并且具有光滑边界，便于制造。

Description

一种基于SIMP法的边界光滑处理方法

技术领域

本发明涉及优化算法领域，具体为一种基于SIMP法的边界光滑处理方法。

背景技术

SIMP(变密度)法是当前连续体拓扑优化设计的主流方法之一，具有编程简便，收敛稳定等优点，得到了广泛的工程应用。

虽然SIMP有着各种优点，但是其缺点同样明显。SIMP法是一种使用0-1之间的变密度单元对单元材料状态进行描述的一种方法。当材料单元使用0表达时，该单元为空域单元，当材料单元使用1表达时，该单元表示实体单元。但当单元使用0到1之间的数值进行表达时则为灰度单元，通常在SIMP法拓扑优化结果的边界处都会出现灰度单元，灰度单元是一种虚拟的情况，并无实际的结构物理状态与之对应。另一方面，由于SIMP法也是基于有限元网格划分，并将每个单元作为确定优化设计变量的对象，因此SIMP法所得到的的拓扑构型会出现边界锯齿现象，从而降低了拓扑优化构型的可制造性。

理论上将单元无限细的划分能避免这种情况，但是对于单元无限细的划分将会极大增加设计变量个数和有限元分析的计算量，将会大大延长拓扑优化所需要的计算机硬件成本和计算耗时，为此我们提出了一种基于SIMP法的边界光滑处理方法。

发明内容

(一)解决的技术问题

针对现有技术的不足，本发明提供了一种基于SIMP法的边界光滑处理方法，在每一个SIMP法的优化迭代步中，首先将SIMP的拓扑优化结果的设计变量过滤将SIMP法的灰度单元部分变成一个沿着边界梯度下降的变密度区域。再引入更高一维度的水平集曲面，水平集曲面的x、y坐标为原设计域平面的各单元中心点处的x、y坐标，水平集曲面的z坐标为经过设计变量过滤后的单元密度值，所有的单元密度值都处于[0,1]域上。

通过一个z＝待定常数的与xoy平面相平行的等值平面与水平集曲面相交，两曲面的交线在xoy平面上的投影即为拓扑优化的边界。

通过以上方法，不仅很好的解决了边界灰度单元问题，并且得到了光滑的结构边界，大大提升了结构的可制造性。

(二)技术方案

为实现上述所述目的，本发明提供如下技术方案：一种基于SIMP法的边界光滑处理方法，包括以下步骤：

第一步：定义边界条件以及材料参数确定目标函数c；

第二步：SIMP法结构模型初始化；

第三步：结合目标函数c和敏度，进行有限元分析并计算；

第四步：基于目标函数c，选用优化准则法对设计变量x_e进行更新；

第五步：对设计变量x_e进行敏度过滤，首先确定过滤半径r_min，将原始单元灵敏度平均分配给节点，然后将根据距离加权后的敏度相加得到修正后的单元灵敏度；

第六步：引入高一维度的水平集曲面，并寻找合适的与水平集曲面相交的等值面，并投影在xoy面；

第七步：检查是否满足优化条件，若满足则输出结果，若不满足则返回第三步。

优选的，所述第一步包括以下内容：

S1：将一个待处理物体进行有限元分析划分成等边长为1的等尺寸矩形单元，右端施加一个F＝1的向下的力；

S2：设定体积分数约束、实体部分的弹性模量E₀、空相部分的弹性模量Emin、材料的泊松比μ、惩罚因子p以及过滤半径r_min，最小柔顺度作为目标函数c。

优选的，所述第二步中SIMP插值函数表达式：

式中，p为惩罚因子，E_e为第e个单元的弹性模量，x_e为设计变量。

优选的，所述第三步包括以下内容：

目标函数c和敏度由以下公式进行计算：

式中，c是目标函数，U是全局位移矢量，K是全局刚度矩阵,V是体积分数，u_e是单位位移矢量，k₀是单元刚度矩阵；

U使用以下公式求解：

KU＝F；

式中，F是力矢量。

优选的，所述第四步中的优化准则法的具体更新方案如下：

m为移动界限，取0.1到0.3，η为阻尼系数，取0.3到0.5，B_e为启发式更新因子；

启发式更新因子B_e通过以下方案求出：

式中，λ是拉格朗日乘子。

优选的，所述第五步中设计变量x_e的灵敏度过滤方案如下：

式中，H_ei为卷积因子，N_e为元素到过滤中心的距离Δ(e，i)小于过滤半径r_min的合集；

卷积因子H_ei用以下方案确定：

H_ei＝(0，r_min-Δ(e，i))。

优选的，所述第六步中设计变量x_e的值作为z轴参数、变量位置作为x、y轴参数。

优选的，所述第七步的优化条件为体积分数约束以及引入的change≤0.01；

change≤0.01时认为设计变量x_e已经稳定，并且已经达到该算法的柔顺度最小；

change＝max(|x_e(i)-x_e(i-1)|)；

式中，x_e(i)表示第i次迭代的设计变量，x_e(i-1)表示第i次迭代前一次迭代的设计变量。

(三)有益效果

与现有技术相比，本发明提供了一种基于SIMP法的边界光滑处理方法，具备以下有益效果：

1、该基于SIMP法的边界光滑处理方法，经过该方法后处理的拓扑构型不再有灰度单元，并且具有光滑边界，便于制造。

2、该基于SIMP法的边界光滑处理方法，具有编程简便，收敛稳定等优点，能够具备有广泛的工程应用。

附图说明

图1为算法流程示意图；

图2为设计变量x_e初始布置部分细节示意图；

图3为敏度过滤器示意图；

图4为未进行灵敏度过滤示意图；

图5为灵敏度过滤后示意图；

图6为常规SIMP法拓扑优化拓扑构型示意图；

图7为边界光滑时的水平集曲面示意图；

图8为边界光滑后的拓扑构型示意图；

图9为柔顺度收敛示意图；

图10为change函数收敛曲线示意图。

具体实施方式

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

请参阅图1-10，一种基于SIMP法的边界光滑处理方法，包括以下步骤：

步骤一：定义边界条件、材料参数等并确定目标函数。

将一个尺寸为60×30悬臂梁进行有限元分析划分成60×30个等边长为1的等尺寸矩形单元，右端施加一个F＝1的向下的力，体积分数约束(优化后实体部分与设计域面积之比)设为0.5，实体部分的弹性模量E₀＝1，为防止发生矩阵奇异性问题空相部分的弹性模量为Emin＝1*e-9,材料的泊松比μ＝0.3，惩罚因子p＝3，过滤半径r_min＝1.5(1.5倍的单元边长)，最小柔顺度作为目标函数。

步骤二：SIMP法结构模型初始化，参阅图2。

SIMP插值函数表达式：

式中，p为惩罚因子，E_e为第e个单元的弹性模量，x_e为设计变量。

优化前的初始设计域为满域，每个单元的均为E₀的单元，结合体积分数的约束，初始每个单元的x_e均设计为体积分数约束的0.5。

步骤三：有限元分析并计算目标函数、敏度等。

结合步骤一的初始化条件目标函数和敏度由以下公式进行计算：

式中，c是目标函数，U是全局位移矢量，K是全局刚度矩阵,V是体积分数。u_e是单位位移矢量，k₀是单元刚度矩阵。

其中U使用以下公式求解：

KU＝F

(5)

式中，F是力矢量。

步骤四：使用合适的方法进行设计变量更新。

如公式(2)可知，目标函数c可以直接通过更新设计变量x_e进行改变更新。这里选用优化准则法对设计变量x_e进行更新，优化准则法的具体更新方案如下：

式中，m为移动界限，通常取0.1到0.3，η为阻尼系数，用来确保算法的收敛性和稳定性，通常取0.3到0.5。B_e为启发式更新因子。

其中启发式更新因子B_e通过以下方案求出：

式中，λ是拉格朗日乘子，通过二分法求得。

步骤五：对设计变量进行敏度过滤。

灵敏度过滤时，首先确定过滤半径(这里取1.5倍单元边长)，将原始单元灵敏度平均分配给节点，然后将根据距离加权后的敏度相加得到修正后的单元灵敏度。对于设计变量的灵敏度过滤方案如下：

式中，H_ei为卷积因子，N_e为元素到过滤中心的距离Δ(e，i)小于过滤半径r_min的合集。

卷积因子H_ei用以下方案确定：

Hei＝(0，rmin-Δ(e，i))

(9)

参阅图4与图5对比灵敏度过滤前后，没有灵敏度过滤的图形有明显的棋盘格现象，经过灵敏度过滤后很好的解决了棋盘格现象，并且从图5可看出，边界灰度单元的值处理成了一组梯度减少的数值。

步骤六：引入高一维度的水平集曲面，并寻找合适的与水平集曲面相交的等值面，并投影在xoy面。

该算例将设计变量的值作为z轴参数、变量位置作为x、y轴参数，取z＝0.5的等值面为相交曲面。对比图6与图8，用本方法进行边界光滑后不仅很好地处理了边界的灰度单元，并且得到了边界光滑的拓扑构型。

步骤七：检查是否满足优化条件，若满足则输出结果，若不满足则返回步骤三。

该算例的体积分数约束为0.5，且需满足柔顺度最小。该算例的体积分数在SIMP函数初始化时已经满足，因此只需要满足柔顺度最小即可，此时引入一个change函数，该函数代表每次迭代设计变量的变化程度。当change≤0.01时，即可认为设计变量已经稳定，并且已经达到该算法的柔顺度最小。change由以下方法求出：

change＝max(|xe(i)-xe(i-1)|)

(10)

式中，x_e(i)表示第i次迭代的设计变量。

由图9、图10可知，该算例在第67步收敛，收敛时体积分数为13.3069，并且此时change值为0.01，满足收敛条件change≤0.01。

尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。

Claims (7)

1.一种基于SIMP法的边界光滑处理方法，其特征在于，包括以下步骤：

第一步：将一个待处理物体进行有限元分析划分成等边长为1的等尺寸矩形单元，右端施加一个F的向下的力，设定体积分数约束、实体部分的弹性模量E₀、空相部分的弹性模量E_min、材料的泊松比μ、惩罚因子p以及过滤半径r_min，最小柔顺度作为目标函数c；

第二步：SIMP法结构模型初始化；

第三步：进行有限元分析并计算目标函数和单元灵敏度；

第四步：基于目标函数c，选用优化准则法对设计变量x_e进行更新；

第五步：对设计变量x_e进行敏度过滤，首先确定过滤半径r_min，将原始单元灵敏度平均分配给节点，然后将根据距离加权后的敏度相加得到修正后的单元灵敏度；

第六步：引入高一维度的水平集曲面，并寻找合适的与水平集曲面相交的等值面，并投影在xoy面；

第七步：检查是否满足优化条件，若满足则输出结果，若不满足则返回第三步。

2.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第二步中SIMP插值函数表达式：

式中，p为惩罚因子，E_e为第e个单元的弹性模量，x_e为设计变量。

3.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第三步包括以下内容：

目标函数c和敏度由以下公式进行计算：

式中，c是目标函数，U是全局位移矢量，K是全局刚度矩阵,V是体积分数，u_e是单位位移矢量，k₀是单元刚度矩阵；

U使用以下公式求解：

KU＝F；

式中，F是力矢量。

4.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第四步中的优化准则法的具体更新方案如下：

m为移动界限，取0.1到0.3，η为阻尼系数，取0.3到0.5，B_e为启发式更新因子；

启发式更新因子B_e通过以下方案求出：

式中，λ是拉格朗日乘子。

5.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第五步中设计变量x_e的灵敏度过滤方案如下：

式中，H_ei为卷积因子，N_e为元素到过滤中心的距离Δ(e，i)小于过滤半径r_min的合集；

卷积因子H_ei用以下方案确定：

H_ei＝(0，r_min-Δ(e，i))。

6.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第六步中设计变量x_e的值作为z轴参数、变量位置作为x、y轴参数。

7.根据权利要求1所述的一种基于SIMP法的边界光滑处理方法，其特征在于：所述第七步的优化条件为体积分数约束以及引入的change≤0.01；

change≤0.01时认为设计变量x_e已经稳定，并且已经达到该算法的柔顺度最小；

change＝max(|x_e(i)-x_e(i-1)|)；

式中，x_e(i)表示第i次迭代的设计变量，x_e(i-1)表示第i次迭代前一次迭代的设计变量。

Images