CN111814383B - 一种基于b样条密度法的自支撑结构拓扑优化设计方法 - Google Patents

Abstract

本发明提供了一种基于B样条密度法的自支撑结构拓扑优化设计方法。首先采用B样条函数表示结构的密度场，通过求导得到密度场的梯度，并根据梯度模以及曲率半径的大小确定结构的边界区域。在边界区域内通过梯度方向计算出边界的倾斜角，并且构造倾角约束，从而控制边界处倾斜角大于设备规定的临界角度。此外，根据不同的成型方向，构造包含结构域的矩形区域，在区域内布置若干均匀分布的检测点，并通过控制这些检测点处密度值之间的关系实现结构内V型区的消除。通过以上两个约束实现结构的自支撑设计。

Description

一种基于B样条密度法的自支撑结构拓扑优化设计方法

技术领域

本发明涉及一种结构拓扑优化设计方法，特别涉及一种基于B样条密度法的自支撑结构拓扑优化设计方法。

背景技术

增材制造技术采用逐层制造的方式成型零件，因此在制造时具有较高的自由度，可以成型形状复杂的零件。然而，增材制造技术依然存在许多制造约束。例如，若结构中存在悬空部位，常见的增材制造技术需要在这些部位下添加额外的支撑避免结构坍塌，且在打印完成后需要去除这些支撑结构。这不仅会造成原料和成本的浪费，延长制造周期，且容易在去除支撑时损伤结构表面。因此在结构拓扑优化设计阶段将自支撑性考虑在内，在满足结构自支撑的同时得到其最优构型，具有非常重要的工程应用价值。

文献“Zhang K,Cheng G,Xu L.Topology optimization considering overhangconstraint in additive manufacturing[J].Computers&Structures,2019,212:86-100.”公开了一种自支撑结构构型的拓扑优化方法。文献提出的方法以密度法为基础，引入两个悬空角约束实现自支撑结构的优化设计。首先，根据有限元模型中单元的局部密度分布估算结构边界的法向向量，从而计算出结构边界的倾斜角，并约束该倾斜角使其大于打印设备规定的临界角度；其次，通过使单元同一水平线上某侧单元的平均密度值大于该单元的密度值，控制结构的最小尺寸，实现结构中V型区的消除。

文献提出的方法虽然可以得到自支撑结构，但是其结构边界处的法向向量通过周围若干单元的密度线性插值得出，计算方式复杂且不够精确。此外，该方法提出的悬空角约束基于规则网格与规则形状设计域施加，并未给出单元形状及排布方式不规则时的解决方案，不利于工程实际中的应用。

发明内容

为了克服现有自支撑结构的拓扑优化设计方法实用性差的不足，本发明提供了一种基于B样条密度法的自支撑结构拓扑优化设计方法。首先采用B样条函数表示结构的密度场，通过求导得到密度场的梯度，并根据梯度模以及曲率半径的大小确定结构的边界区域。在边界区域内通过梯度方向计算出边界的倾斜角，并且构造倾角约束，从而控制边界处倾斜角大于设备规定的临界角度。此外，根据不同的成型方向，构造包含结构域的矩形区域，在区域内布置若干均匀分布的检测点，并通过控制这些检测点处密度值之间的关系实现结构内V型区的消除。通过以上两个约束实现结构的自支撑设计。

相比背景技术的设计方法，本发明中结构边界处梯度直接通过求导计算得出，因此悬空角度的计算更为方便与准确。此外，引入的悬空角约束与分析单元无关，不依赖局部单元的位置与形状，因此可以采用不规则形状的单元进行有限元分析，对复杂形状设计域进行优化设计，具有较强的工程实用性。

本发明的技术方案为：

一种基于B样条密度法的自支撑结构拓扑优化设计方法，包括以下步骤：

步骤1：建立自支撑结构的结构设计域Ω_d，并构造出表示结构设计域的密度场ρ_d；

根据结构设计域Ω_d的形状，构造能够将结构设计域Ω_d包含在内的矩形区域Ω；在矩形区域Ω内部构造B样条曲面作为该区域内的B样条密度场ρ；对B样条密度场ρ进行Heaviside投影得到投影后的B样条密度场

将投影后的B样条密度场

与表示结构设计域的密度场ρ_d进行布尔交集运算，得到结构的最终密度场

步骤2：通过对结构设计域Ω_d内布置的若干积分点求导得到结构密度场

的梯度

并继而得到梯度模，依据结构设计域Ω_d内各积分点梯度模以及曲率半径，识别得到结构的边界区域Ω_b；

步骤3：考虑到要在结构边界上施加倾角约束，这里在边界区域内通过梯度方向计算出边界的倾斜角，并且构造倾角约束，从而控制边界处倾斜角大于设备规定的临界角度。

对于步骤2得到的边界区域Ω_b，计算边界区域Ω_b内各个积分点处，B样条密度场ρ的梯度

进而根据公式

计算边界区域Ω_b内各个积分点处的倾角α余弦值cosα，其中b为方向指向增材制造成型方向的单位向量，

表示B样条密度场的梯度模；

根据倾角α与设备规定的临界角度α₀的余弦值之差，在结构的边界区域Ω_b内进行积分，得到倾角约束值q：

其中L(δ)是关于变量δ的分段函数，当δ<0时，L(δ)＝0；当δ>0时，L(δ)＝δ²。

利用倾角约束值q可以建立构造倾角约束函数：

ε_q是为了松弛该倾角约束而引入的极小的数，通过施加该倾角约束，使得结构边界处任意点的倾角大于设备规定的临界角度。

步骤4：根据增材制造的成型方向，构造将区域Ω包含在内的新矩形区域Ω_t，Ω_t的一条边平行于增材制造成型方向，另一条边垂直于增材制造成型方向；在矩形区域Ω_t内沿自身两条边的方向均匀布置若干检测点，并取位于自支撑结构底边正方向一侧的检测点为有效检测点，所述正方向为增材制造成型方向；

根据结构最终密度场

计算有效检测点处的密度值，若有效检测点位于结构设计域Ω_d外，则该有效检测点的密度取0，得到有效检测点的密度值向量

在有效检测点中，将最底层检测点排除在外后，剩余的n_t个检测点为参与计算的检测点。这里的最底层就是指沿增材制造成型方向的最下层。

对于第i个需参与计算的检测点，约束该点处的密度值

应小于位于其下方三角区域内n_c个点处密度值中的最大值，用公式表示为：

这里所谓下方三角区域，也是在增材制造成型方向上所定义的下方，n_c个点其实就是第i个需参与计算的检测点下一层中的n_c个点，优选包括第i个需参与计算的检测点正下方的点，以及所述正下方的点两侧的若干点，共同组成这n_c个点。

分别表示第i个检测点下方三角区域内n_c个点的密度值，下标中第一项i表示三角区域顶部检测点对应的序号，第二项1～n_c分别表示三角区域内部的n_c个检测点。为了便于使用梯度算法进行优化，采用KS函数近似上式中的max函数，并计算整体V型区约束值：

w是KS函数中的参数，ε_g为了松弛该V型区约束而引入的极小的数，h(δ)是关于变量δ的函数，其形式为：

其中参数μ的大小决定该函数在δ＝0处变化的剧烈程度。

利用V型区约束值g可以建立整体V型区约束函数

通过施加该V型区约束，实现优化结果中V型区的消除。

步骤5：将矩形区域Ω离散为有限元单元，同时定义载荷和边界条件；

步骤6：定义拓扑优化问题为：

min C+λg

其中P表示矩形区域Ω内所构造B样条曲面的控制点的值，n_x和n_y分别表示矩形区域Ω两个方向上的控制点的数目，P_i中下标i表示控制点的序号，所有控制点取值介于0和1之间。C表示结构的柔顺度。g表示V型区约束，在优化问题中作为惩罚项加入到目标函数中，λ是惩罚参数。K、F和U分别表示结构的总体刚度矩阵、总体载荷向量和位移向量。V、V₀和

分别表示结构的当前体积分数、总体积和体积分数上限。q表示结构边界区域的倾角约束。

步骤7：采用优化算法对步骤6的拓扑优化问题进行求解，得到最优结果。

进一步的，步骤2中识别结构的边界区域的具体过程为：

对于结构设计域Ω_d内的积分点(x,y)，在积分点(x,y)处，结构密度场

的梯度

以及梯度模

为：

其中x和y分别为笛卡尔坐标系中的坐标。由于梯度模的大小反映了密度场变化的剧烈程度，而密度场在结构边界附近变化更为剧烈，因此可以通过梯度模的大小确定结构边界。预先给定一个梯度模下限

则根据梯度模确定的结构的边界区域Ω_b可以表示为：

为了排除结构内孔洞上方的小圆角区域对要施加的倾角约束的影响，计算所有区域的曲率半径，并且将曲率半径较小的区域排除在外。根据隐函数曲率半径的计算公式，点(x,y)处的曲率半径r表示为：

通过预先给定一个曲率半径下限r_min，则去掉小圆角之后的边界区域Ω_b可以重新表示为：

有益效果

本发明的有益效果是：该方法通过使用B样条函数构造结构密度场，从而根据结构边界处的梯度求得边界倾斜角度，控制悬空角的大小。与此同时，通过在结构域内布置均匀的检测点，控制检测点处密度值之间的关系，实现结构内部V型区的消除。相比背景技术的设计方法，本发明中结构边界处的梯度直接通过求导计算，倾角控制更为方便与准确。此外，本发明提出的悬空角约束与分析单元无关，不依赖局部单元的位置与形状，可对复杂形状设计域进行优化设计，具有较强的工程实用性。

本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。

附图说明

本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得明显和容易理解，其中：

图1是本发明实施例中模型几何尺寸与边界条件示意图。

图2是本发明实施例中结构域内检测点分布的示意图。

图3是本发明实施例的设计结果图。

具体实施方式

下面详细描述本发明的实施例，所述实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。

参照图1-3，以固定载荷下简支梁结构为例说明。由于该简支梁的对称性，仅取其右半部分进行优化设计。其杨氏模量为10，泊松比为0.3。规定该结构的成型方向用向量b＝(cos 45°,sin 45°)表示，打印设备的临界成型角度为α₀＝45°。选取设计该简支梁中的材料分布使得其柔顺度最小，总材料用量体积分数最大为50％。具体方法步骤为：

1、构造结构密度场。首先构造结构设计域的密度场ρ_d。构造正方形区域Ω将设计域Ω_d包含在内，采用左右重节点3次B样条基函数进行插值得到正方形区域Ω内的B样条密度场ρ。对该B样条密度场进行Heaviside投影，并与设计域密度场ρ_d进行布尔交集运算，得到结构最终密度场

2、识别结构边界。该实施例中预先给定的梯度模下限为

曲率半径下限为r_min＝1.5。因此边界区域Ω_b表示为

3、施加倾角约束。在边界区域Ω_b内选取积分点，通过求导计算积分点处的梯度并进一步计算其倾角的余弦值，计算公式为

根据该倾斜角余弦值与设备规定临界角度的余弦值之差，在区域Ω_b内进行积分，构造倾角约束函数：

该实施例中允许区域Ω_b内的倾斜角度存在0.5°的容差，因此ε_q计算为

4、V型区的消除。沿增材制造的成型方向及其垂直方向，构造将正方形区域Ω包含在内的新正方形区域。沿着该新正方形区域的两条边规则且均匀地布置若干检测点，并根据结构最终密度场

计算得出这些检测点处的密度值

该实施例中总检测点数目为113×114，有效检测点数目为113×80，去除有效检测点中最底层检测点之后剩余应参与计算检测点数目为113×79。整体V型区约束函数为：

式中n_t＝8927，w＝10，n_c＝5，ε_g＝0.15。h(δ)函数中参数μ的取值为μ＝40。

5、采用数目为120×120的矩形网格将Ω离散。定义边界条件为：简支梁左边界约束x方向位移，右下角约束y方向位移，在简支梁左上角施加沿y轴负方向的载荷，大小为1。

6、定义拓扑优化问题：

min C+λg

式中P表示所有控制点的值，其初始值均为0.5，上限和下限分别是1和0。n_x和n_y分别表示x方向与y方向上控制点的数目，n_x＝43，n_y＝43。惩罚参数λ在优化前100步取值为6×10⁵，然后每100步乘以系数0.5直到达到其最小值1×10⁵。

7、有限元分析和优化求解。对以上建立的模型进行有限元分析，并分别对目标函数和约束函数进行灵敏度分析。选取全局收敛移动渐近线优化算法(GCMMA)进行优化设计，得到最优化结果。

本实施例方法经过实施例500步迭代后得到优化设计结果。初始结构的体积分数为50％，柔顺度为50.6803，优化结果中结构体积分数为50％，柔顺度为9.9733，柔顺度减少了80.32％。在不考虑自支撑的情况下，用相同初始结构优化得到的结构柔顺度为9.0152，因此实现结构自支撑牺牲了10.63％的柔顺度。相比背景技术的设计方法，本发明中结构边界处的梯度直接通过求导计算，倾角控制更为方便与准确。此外，本发明提出的悬空角约束与分析单元无关，不依赖局部单元的位置与形状，可对复杂形状设计域进行优化设计，具有较强的工程实用性。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在不脱离本发明的原理和宗旨的情况下在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (2)

1.一种基于B样条密度法的自支撑结构拓扑优化设计方法，其特征在于：包括以下步骤：

步骤1：建立自支撑结构的结构设计域Ω_d，并构造结构设计域的密度场ρ_d；

构造能够将结构设计域Ω_d包含在内的矩形区域Ω；在矩形区域Ω内部构造B样条曲面作为该区域内的B样条密度场ρ；对B样条密度场ρ进行投影得到投影后的B样条密度场

将投影后的B样条密度场

与表示结构设计域的密度场ρ_d进行布尔交集运算，得到结构的最终密度场

步骤2：通过对结构设计域Ω_d内布置的若干积分点求导得到结构密度场

的梯度

并继而得到梯度模，依据结构设计域Ω_d内各积分点梯度模以及曲率半径，识别得到结构的边界区域Ω_b：

对于结构设计域Ω_d内的积分点(x,y)，在积分点(x,y)处，结构密度场

的梯度

以及梯度模

为：

并根据公式

得到点(x,y)处的曲率半径r_(x,y)；利用设定的梯度模下限

和曲率半径下限r_min，得到边界区域Ω_b为

步骤3：在边界区域Ω_b中，对各个积分点的倾角余弦值与设备规定的临界角度余弦值之差进行积分，得到倾角约束值q：

对于步骤2得到的边界区域Ω_b，计算边界区域Ω_b内各个积分点处，B样条密度场ρ的梯度

进而根据公式

计算边界区域Ω_b内各个积分点处的倾角α余弦值cosα，其中b为方向指向增材制造成型方向的单位向量，

表示B样条密度场的梯度模；

根据倾角α与设备规定的临界角度α₀的余弦值之差，在结构的边界区域Ω_b内进行积分，得到倾角约束值q：

其中L(δ)是关于变量δ的分段函数，当δ<0时，L(δ)＝0；当δ>0时，L(δ)＝δ²；

步骤4：根据增材制造的成型方向，构造将区域Ω包含在内的新矩形区域Ω_t，Ω_t的一条边平行于增材制造成型方向，另一条边垂直于增材制造成型方向；在矩形区域Ω_t内沿自身两条边的方向均匀布置若干检测点，并取位于自支撑结构底边正方向一侧的检测点为有效检测点，所述正方向为增材制造成型方向；

根据结构最终密度场

计算有效检测点处的密度值，若有效检测点位于结构设计域Ω_d外，则该有效检测点的密度取0，得到有效检测点的密度值向量

在有效检测点中，将最底层检测点排除在外后，剩余的n_t个检测点为参与计算的检测点；对于第i个参与计算的检测点，找到其下方的n_c个点，

分别表示第i个检测点下方n_c个点的密度值；

根据公式

计算整体V型区约束值g，其中w是KS函数参数，ε_g为引入的极小数，h(δ)是关于变量δ的函数，其形式为：

μ为设定参数；

步骤5：将矩形区域Ω离散为有限元单元，同时定义载荷和边界条件；

步骤6：定义拓扑优化问题为：

minC+λg

其中P表示矩形区域Ω内所构造B样条曲面的控制点的值，n_x和n_y分别表示矩形区域Ω两个方向上的控制点的数目，P_i中下标i表示控制点的序号，所有控制点取值介于0和1之间；C表示结构的柔顺度；g表示V型区约束，在优化问题中作为惩罚项加入到目标函数中，λ是惩罚参数；K、F和U分别表示结构的总体刚度矩阵、总体载荷向量和位移向量；V、V₀和

分别表示结构的当前体积分数、总体积和体积分数上限；q表示结构边界区域的倾角约束；

步骤7：采用优化算法对步骤6的拓扑优化问题进行求解，得到最优结果。

2.根据权利要求1所述一种基于B样条密度法的自支撑结构拓扑优化设计方法，其特征在于：步骤4中第i个参与计算的检测点下方的n_c个点为第i个参与计算的检测点下方一层检测点中的n_c个点，包括第i个参与计算的检测点正下方的点，以及所述正下方的点两侧的若干点。

Images