CN109190233B - 一种结构拓扑优化方法 - Google Patents

一种结构拓扑优化方法 Download PDF

Info

Publication number
CN109190233B
CN109190233B CN201810979982.1A CN201810979982A CN109190233B CN 109190233 B CN109190233 B CN 109190233B CN 201810979982 A CN201810979982 A CN 201810979982A CN 109190233 B CN109190233 B CN 109190233B
Authority
CN
China
Prior art keywords
optimization
level set
delta
topology
phi
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810979982.1A
Other languages
English (en)
Other versions
CN109190233A (zh
Inventor
魏鹏
汪雯雯
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
South China University of Technology SCUT
Original Assignee
South China University of Technology SCUT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by South China University of Technology SCUT filed Critical South China University of Technology SCUT
Priority to CN201810979982.1A priority Critical patent/CN109190233B/zh
Publication of CN109190233A publication Critical patent/CN109190233A/zh
Application granted granted Critical
Publication of CN109190233B publication Critical patent/CN109190233B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/18—Network design, e.g. design based on topological or interconnect aspects of utility systems, piping, heating ventilation air conditioning [HVAC] or cabling

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Geometry (AREA)
  • Computer Hardware Design (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • Computational Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Image Generation (AREA)
  • Design And Manufacture Of Integrated Circuits (AREA)

Abstract

本发明公开了一种结构拓扑优化方法，包括如下步骤：定义初始结构和相应的水平集函数；将当前的节点水平集函数值φN平均插值得到各单元初始水平集函数值φ；引入可调节宽度的水平集带概念和用于控制中间密度范围的参数Δ，将单元水平集函数值φ通过参数为Δ的映射函数处理得到分布在0‑1之间的单元密度值H(φ)，即得到含有中间密度的拓扑结构；对拓扑结构进行有限元分析，更新节点水平集函数值φN；迭代并进行收敛判断，得到具有清晰边界的结构设计；输出优化结果。本发明可实现连续的拓扑演化，解决边界相关的优化问题，确保优化结果的物理可制造性，可得到合理的结构拓扑，最终收敛时可以得到通过零水平集表达的具有清晰边界的结构设计。

Description

一种结构拓扑优化方法
技术领域
本发明涉及结构优化设计相关技术领域，特别涉及一种基于可调节宽度水平集带的结构拓扑优化方法。
背景技术
基于材料插值模型的变密度法(以SIMP方法为代表)和水平集方法是目前拓扑优化的两种常用方法。
变密度法是引入假想的可在0-1之间连续取值的单元伪密度作为设计变量，然后假定材料物理属性值如弹性模量与伪密度之间存在函数关系的拓扑优化方法。其中，SIMP模型是变密度法中最常用的模型，其单元弹性模量E取为单元伪密度的函数：
E(ρ)＝ρpE0
其中，p为惩罚因子，ρ为材料伪密度，E0为实体材料的弹性模量。
在实际的应用中，当p取值较小时，优化结构含有大部分的灰色单元即中间密度单元，而当p取值较大时又往往会出现收敛过快从而陷入局部最优的问题。因此，目前实际应用中一般取p＝3。该方法通过惩罚因子p的作用实现材料伪密度趋近于0或1从而逐渐去除结构中的中间密度单元，但该方法有两个问题：1.对于三维复杂结构最终结果往往仍然存在大量中间密度单元无法消除；2.对于同时涉及质量与刚度的优化问题如结构自重及动力学问题，需要质量与刚度同时进行惩罚，但由于质量与刚度的物理属性不同，同样的惩罚因子会因为惩罚效果不匹配导致收敛问题，往往需要进行额外处理。
基于水平集的拓扑优化方法，是利用高维标量水平集函数场φ的零等值线(二维问题)或零等值面(三维问题)，即φ＝0，来隐式描述结构的几何轮廓或不同材料的交界面，并利用特定的速度场驱动结构的边界演化，得到非0即1的离散材料分布和明确的结构边界，但水平集方法由于不存在中间密度，拓扑演化过程不连续，往往存在初始设计的依赖性问题，即初始设计对最终设计的影响较大。
发明内容
本发明的目的在于克服现有技术的缺点与不足，提供一种结构拓扑优化方法，此方法能充分发挥中间密度的拓扑表现潜力，并得到收敛的具有明确边界的优化结果，以解决边界相关的优化问题和确保优化结果的物理可制造性。
本发明的目的通过以下的技术方案实现：一种结构拓扑优化方法，包括如下步骤：
S1、选择一个设计域并进行离散和有限元网格划分，导出节点信息、单元信息；
S2、选定需要优化的目标函数，根据实际工作状况施加位移约束条件和荷载；
S3、给定初始水平集带宽度2Δinit及各节点的初始水平集函数值φNinit，使材料在设计域范围内满布；
S4、将当前的节点水平集函数值φN平均插值得到各单元初始水平集函数值φ；
S5、将单元水平集函数值φ通过参数为Δ的映射函数处理得到分布在0-1之间的单元密度值H(φ)，即得到含有中间密度的拓扑结构，Δ表示用于控制中间密度范围的参数；
S6、对上一步得到的拓扑结构进行有限元分析，得到当前用于驱动速度场演化的响应量，进而得到新的节点水平集函数值φN；
S7、迭代收敛判断，当优化目标达到标准时，优化迭代结束，执行步骤S8，否则按式Δ＝Δ-Δt更新参数Δ值，并重复步骤S4至步骤S7，Δt表示预设的每步减小幅度；
S8、输出优化结果，从而得到通过零水平集表达的具有清晰边界的结构设计。
优选的，所述步骤S1中通过编译代码对设计域进行离散和有限元网格划分。
优选的，所述步骤S1中在有限元建模、分析软件中对设计域进行离散和有限元网格划分。
优选的，所述步骤S1中的节点信息包括所有节点的编号及坐标，单元信息包括单元编号以及组成每个单元的节点编号。
优选的，所述步骤S2中目标函数的优化目标包括结构应变能最小化、结构刚度最大化、动态振动频率优化。
优选的，所述步骤S2中，位移约束条件的约束信息包括约束点编号及坐标、被约束的自由度，荷载信息包括：受力点编号及坐标、受力方向对应的自由度和受力大小。
优选的，所述步骤S5中的映射函数为Heaviside函数，其表达式如下：
Figure BDA0001776638430000031
其中，a为接近0的正数，Δ是控制中间密度范围的参数。
优选的，所述步骤S7中设定每步减小幅度为Δt，初始值为Δinit，下限为Δmin；那么在每一步迭代中，可调节的水平集带宽度为2Δ。当Δmin取非常接近0的正数时，可调节的水平集带宽度为接近0的2Δmin，拓扑结构的中间密度范围极小，当Δmin取为0时，则完全去除中间密度，得到具有明确边界的拓扑优化结果。
本发明与现有技术相比，具有如下优点和有益效果：
1、本发明利用使单元密度可以取0-1之间的中间值这一策略，实现连续的拓扑演化，从而解决边界相关的优化问题和确保优化结果的物理可制造性。
2、本发明方法，开始取一个比较大的水平集带宽度(例如10倍单元尺寸)，然后逐渐减小，从而通过水平集带宽度这一个参数的改变实现从密度法到水平集法的平滑过渡，在开始阶段可以得到合理的结构拓扑，最终收敛时可以得到通过零水平集表达的具有清晰边界的结构设计。
3、本发明在映射函数处理阶段引入Heaviside函数，可以使不同区域的材料性质光滑过度，从而得到合理的结构拓扑，最终收敛时可以得到通过零水平集表达的具有明确边界的结构设计。
附图说明
图1是本发明单元初始水平集函数值φinit示意图。
图2是本发明对单元初始水平集函数值φinit进行Heaviside处理，所得单元密度值H(φ)示意图。
图3是本发明迭代过程中，单元水平集函数值φ示意图。
图4是本发明迭代过程中，水平集带宽度为2Δ，对上一步水平集函数进行Heaviside处理，所得单元密度值H(φ)的示意图。
图5是本发明迭代终止时，单元水平集函数值φ示意图。
图6是本发明迭代终止时，水平集带宽度为2Δmin，对上一步水平集函数进行Heaviside处理，所得单元密度值H(φ)示意图。
图7是本发明一种结构拓扑优化方法的流程图。
图8是本发明的实施例一个矩形初始设计域的示意图。
图9是本发明的实施例取一个比较大的水平集带初始宽度Δinit＝5，然后逐渐减小至Δmin＝0.5获得的实施例的拓扑结果，(a)初始拓扑图，(b)迭代过程中的拓扑图、(c)最终拓扑图。
图10是本发明的实施例取一个固定的比较小的水平集带初始宽度Δ＝0.5获得的实施例的拓扑结果，(a)初始拓扑图，(b)迭代过程中的拓扑图、(c)最终拓扑图。
具体实施方式
为了更好的理解本发明的技术方案，下面结合附图详细描述本发明提供的实施例，但本发明的实施方式不限于此。
实施例1
如图7所示，本实施例一种结构拓扑优化方法，包括步骤如下：
S1、本实施例中以图8所示的二维平面应力结构为例来进行结构拓扑优化的说明，该结构中矩形初始设计域左端为固定约束，长0.8m，宽0.4m，厚度为0.01m，实体材料的弹性模量为2.1×1011Pa，空白材料的弹性模量为2.1×102pa，泊松比为0.3，结构在右端中点处受到竖直方向集荷载100N的作用，以结构应变能最小化为优化目标，施加体积等式约束，体积分数取50％。
S2、如图9-10所示采用四节点平面应力单元，通过有限元建模，利用Abaqus将设计域划分为80×40＝3200个单元的网格，并设置单元各节点初始水平集函数值φNinit，使材料在设计域内均匀分布，给定初始水平集带宽度2Ainit。
S3、将当前的节点水平集函数值φN平均插值得到单元水平集函数值φ，如图1所示。
S4、将单元水平集函数值φ通过参数为Δ的Heaviside函数处理得到分布在0-1之间的单元密度值H(φ)，如图2所示，即得到含有中间密度的拓扑结构；
Heaviside函数表达式如下：
Figure BDA0001776638430000061
其中，a取接近0的正数，Δ是控制中间密度范围的参数。
S5、针对步骤S4得到的含有中间密度的拓扑结构，采用基于水平集的拓扑优化方法进行有限元分析，该方法可采用现有技术中已有的方法进行计算，并更新节点水平集函数值φN，以现有技术中基于径向基函数的参数化水平集拓扑优化方法为例，借助一应用实例来对φN的更新进行说明；
总体刚度矩阵K计算公式如下：
Figure BDA0001776638430000062
其中，m为设计域内单元总数，Emin为弱材料弹性模量，E0为实体材料弹性模量，
Figure BDA0001776638430000063
为单元刚度矩阵；
求解平衡方程KU＝F，计算各单元的应变能值Ce＝UTKU，其中，T为转置运算符号，U为位移场，根据各节点所在单元的应变能值Ce加权得到各节点应变能值Cn，以各节点应变能值Cn作为驱动速度场演化的响应量，并求解各节点水平集函数更新值；
此处参考基于径向基函数的参数化水平集拓扑优化方法的水平集函数更新策略，通过更新系数矩阵α(t)来更新φN：
Figure BDA0001776638430000064
Vn＝Cn-λ
φN(x)＝G(x)α(t) (1)
Figure BDA0001776638430000071
Figure BDA0001776638430000072
其中，G(x)为全局径向基MQ函数，x为设计域中节点坐标，xi为基于径向基函数的参数化水平集拓扑优化方法中选取的knot点坐标，λ为处理体积约束的拉格朗日乘子，ci为自由形状参数，
Figure BDA0001776638430000073
为向量微分算子，公式(3)由式(1)代入公式(2)中得到，以更新系数矩阵α(t)。
S6、迭代收敛判断，比较结构应变能在两次连续迭代的相对变化量，当变化量小于一个事先给定的数值时，优化迭代结束，否则按式Δ＝Δ-Δt减小参数Δ值，并重复S3-S6；
具体的，在迭代开始时，预设每步减小幅度为Δt，初始值为Δinit，下限为Δmin。开始时水平集带宽度较宽，在迭代过程中，水平集带宽度逐渐减小，随着迭代的进行，到某一步时，水平集函数值φ如图3所示，所得单元密度值H(φ)如图4所示，可看到拓扑结构的中间密度范围已经得到缩小，持续迭代，在迭代满足条件终止时，水平集函数值φ如图5所示，此时的单元密度值H(φ)如图6所示，中间密度已基本完全去除，得到具有明确边界的拓扑优化结果。
S7、优化结果输出。
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (8)

1.一种基于可调节宽度水平集带的结构拓扑优化方法，其特征在于，包括如下步骤：
S1、对待优化结构的设计域进行离散和有限元网格划分，导出节点信息、单元信息；
S2、选定待优化结构需要优化的目标函数，根据实际工作状况对待优化结构施加位移约束条件和荷载；
S3、给定初始水平集带宽度2Δinit及各节点的初始水平集函数值φNinit，使材料在设计域范围内满布；
S4、将当前的节点水平集函数值φN平均插值得到各单元初始水平集函数值φ；
S5、将单元水平集函数值φ通过参数为Δ的映射函数处理得到分布在0-1之间的单元密度值H(φ)，即得到含有中间密度的拓扑结构，Δ表示用于控制中间密度范围的参数；
S6、对上一步得到的拓扑结构采用基于水平集的拓扑优化方法进行有限元分析，得到当前用于驱动速度场演化的响应量，进而得到新的节点水平集函数值φN；
S7、迭代收敛判断，当优化目标达到标准时，优化迭代结束，执行步骤S8，否则按式Δ＝Δ-Δt更新参数Δ值，并重复步骤S4至步骤S7，Δt表示预设的每步减小幅度；
S8、输出优化的结构拓扑。
2.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S1中通过编译代码对设计域进行离散和有限元网格划分。
3.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S1中在有限元建模、分析软件中对设计域进行离散和有限元网格划分。
4.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S1中的节点信息包括所有节点的编号及坐标，单元信息包括单元编号以及组成每个单元的节点编号。
5.根据权利要求1所述的结构拓扑优化方法，其特征在于：所述步骤S2中目标函数的优化目标包括结构应变能最小化、结构刚度最大化、动态振动频率优化。
6.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S2中，位移约束条件的约束信息包括约束点编号及坐标、被约束的自由度，荷载信息包括：受力点编号及坐标、受力方向对应的自由度和受力大小。
7.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S5中的映射函数为Heaviside函数，其表达式如下：
Figure FDA0002480412490000021
其中，a为接近0的正数，Δ是控制中间密度范围的参数。
8.根据权利要求1所述的结构拓扑优化方法，其特征在于，所述步骤S7中设定每步减小幅度为Δt，初始值为Δinit，下限为Δmin；那么在每一步迭代中，可调节的水平集带宽度为2Δ；当Δmin取非常接近0的正数时，可调节的水平集带宽度为接近0的2Δmin，拓扑结构的中间密度范围极小，当Δmin取为0时，则完全去除中间密度，得到具有明确边界的拓扑优化结果。
CN201810979982.1A 2018-08-24 2018-08-24 一种结构拓扑优化方法 Active CN109190233B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810979982.1A CN109190233B (zh) 2018-08-24 2018-08-24 一种结构拓扑优化方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810979982.1A CN109190233B (zh) 2018-08-24 2018-08-24 一种结构拓扑优化方法

Publications (2)

Publication Number Publication Date
CN109190233A CN109190233A (zh) 2019-01-11
CN109190233B true CN109190233B (zh) 2020-11-24

Family

ID=64916038

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810979982.1A Active CN109190233B (zh) 2018-08-24 2018-08-24 一种结构拓扑优化方法

Country Status (1)

Country Link
CN (1) CN109190233B (zh)

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110321611A (zh) * 2019-06-24 2019-10-11 华中科技大学 一种多材料结构拓扑优化方法
CN110555263B (zh) * 2019-08-30 2023-04-07 华南理工大学 一种用于曲壳结构优化设计的水平集拓扑优化方法
CN110837709B (zh) * 2019-11-05 2022-05-13 广西艾盛创制科技有限公司 一种用于榫卯连接结构设计的归选式拓扑优化方法
CN111027110B (zh) * 2019-11-27 2023-06-30 中国科学院光电技术研究所 一种连续体结构拓扑与形状尺寸综合优化方法
CN112100877B (zh) * 2020-08-10 2022-05-24 华南理工大学 一种结构刚度高效拓扑优化方法及系统
CN112100882B (zh) * 2020-08-27 2024-03-15 华南理工大学 一种具有光滑边界的连续体结构密度演化拓扑优化方法
CN113191044B (zh) * 2021-04-13 2023-03-28 华中科技大学 一种单材料多孔结构的拓扑优化设计方法
CN113673123B (zh) * 2021-07-02 2023-09-01 华南理工大学 一种实体-壳耦合结构拓扑优化方法
CN113705060B (zh) * 2021-10-21 2022-04-15 中南大学 考虑边界优化的拓扑优化方法、系统及存储介质

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR101628823B1 (ko) * 2015-02-24 2016-06-09 한양대학교 산학협력단 조화 탐색법을 이용한 구조물의 위상 및 형상 동시 최적화 장치 및 방법
CN105956292A (zh) * 2016-05-05 2016-09-21 河北工业大学 一种进化水平集结构拓扑优化方法
CN107025340A (zh) * 2017-03-30 2017-08-08 华中科技大学 一种适用于增材制造的自支撑网状结构拓扑优化设计方法
CN107315872A (zh) * 2017-06-23 2017-11-03 华中科技大学 一种高效的结构频率响应拓扑优化方法
CN107766624A (zh) * 2017-09-28 2018-03-06 华中科技大学 一种基于多模板快速推进算法的结构拓扑优化方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR101628823B1 (ko) * 2015-02-24 2016-06-09 한양대학교 산학협력단 조화 탐색법을 이용한 구조물의 위상 및 형상 동시 최적화 장치 및 방법
CN105956292A (zh) * 2016-05-05 2016-09-21 河北工业大学 一种进化水平集结构拓扑优化方法
CN107025340A (zh) * 2017-03-30 2017-08-08 华中科技大学 一种适用于增材制造的自支撑网状结构拓扑优化设计方法
CN107315872A (zh) * 2017-06-23 2017-11-03 华中科技大学 一种高效的结构频率响应拓扑优化方法
CN107766624A (zh) * 2017-09-28 2018-03-06 华中科技大学 一种基于多模板快速推进算法的结构拓扑优化方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
A level set method for structural topology optimization and its applications;Mei Yulin 等;《Advances in Engineering Software》;20040731;第35卷(第7期);415-441 *
On projection methods，convergence and robust formulations in topology optimization;Fengwen wang 等;《Structural and Multidisciplinary optimization》;20110630;第43卷(第6期);767-784 *
面向建筑结构设计的水平集拓扑优化方法研究;李祖宇;《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》;20180615(第6期);C038-345 *

Also Published As

Publication number Publication date
CN109190233A (zh) 2019-01-11

Similar Documents

Publication Publication Date Title
CN109190233B (zh) 一种结构拓扑优化方法
CN110069800B (zh) 具有光滑边界表达的三维结构拓扑优化设计方法及设备
CN107341316B (zh) 设计相关压力载荷作用下的结构形状-拓扑联合优化方法
Zhang et al. Optimal topology design of internal stiffeners for machine pedestal structures using biological branching phenomena
CN110414127B (zh) 一种面向增材制造的支撑体积约束拓扑优化方法
CN108763658A (zh) 基于等几何方法的组合薄壁结构固有频率设计方法
CN110555263B (zh) 一种用于曲壳结构优化设计的水平集拓扑优化方法
CN110263384B (zh) 基于Laplacian微分域变形的三维网格曲面变厚度偏置造型方法
CN111814246B (zh) 一种基于生成对抗网络的翼型反设计方法
CN113569360B (zh) 一种风力机叶片抗颤振翼型簇设计方法
CN110210079B (zh) 一种面向整机动态特性的机床支承件质量匹配方法
CN107025335A (zh) 基于状态变量离散化的仿真计算方法和仿真系统
CN111859763A (zh) 有限元模拟方法、系统及介质
CN113313322B (zh) Moea/d挤压工艺参数多目标优化方法及装置
CN113536623B (zh) 一种材料不确定性结构稳健性拓扑优化设计方法
CN113239584B (zh) 一种优化增材制造方法及系统
CN113779802A (zh) 基于无网格efgm和等几何分析耦合方法的结构拓扑优化技术
CN111597724B (zh) 考虑频带约束的结构动力学拓扑优化方法、系统
CN116756851A (zh) 基于nffd背景网格的参数化网格变形方法及系统
CN113036769B (zh) 一种电力系统静态电压稳定分析方法及系统
CN115272594A (zh) 一种基于geotools的等值面生成方法
CN110826182B (zh) 基于顶点法和序贯优化策略的飞行器结构气动弹性设计法
CN115455754A (zh) 一种基于数字孪生的矿山液压支架设计方法
CN114329320A (zh) 一种基于启发式训练数据采样的偏微分方程数值求解方法
CN110555267B (zh) 一种基于隐式b-样条的参数化水平集结构拓扑优化方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant