CN114999591A - 一种多构型点阵结构的拓扑优化方法 - Google Patents

一种多构型点阵结构的拓扑优化方法 Download PDF

Info

Publication number
CN114999591A
CN114999591A CN202210507811.5A CN202210507811A CN114999591A CN 114999591 A CN114999591 A CN 114999591A CN 202210507811 A CN202210507811 A CN 202210507811A CN 114999591 A CN114999591 A CN 114999591A
Authority
CN
China
Prior art keywords
lattice
configuration
lattice structure
matrix
unit cells
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.)
Granted
Application number
CN202210507811.5A
Other languages
English (en)
Other versions
CN114999591B (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.)
Huazhong University of Science and Technology
Original Assignee
Huazhong University of Science and Technology
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 Huazhong University of Science and Technology filed Critical Huazhong University of Science and Technology
Priority to CN202210507811.5A priority Critical patent/CN114999591B/zh
Priority claimed from CN202210507811.5A external-priority patent/CN114999591B/zh
Publication of CN114999591A publication Critical patent/CN114999591A/zh
Application granted granted Critical
Publication of CN114999591B publication Critical patent/CN114999591B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16CCOMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C60/00Computational materials science, i.e. ICT specially adapted for investigating the physical or chemical properties of materials or phenomena associated with their design, synthesis, processing, characterisation or utilisation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/25Design optimisation, verification or simulation using particle-based methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Computing Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Complex Calculations (AREA)

Abstract

本发明公开了一种多构型点阵结构的拓扑优化方法,属于拓扑结构优化技术领域,所述方法包括:基于水平集函数描述多种构型点阵的单胞,并在相同等效密度区间中获取各个所述构型基础点阵单胞对应的一系列点阵单胞,并计算其等效弹性张量;建立多构型点阵结构每个有限单元对应的弹性张量关于单元相对密度的插值函数,利用所述插值函数构建多构型点阵结构的拓扑优化模型;采用交替活动相算法将所述拓扑优化问题分解为多个子问题并对其求解,从而对点阵结构的宏观布局进行优化,得到目标多构型点阵结构。本发明实现了多种拓扑构型不同的点阵单胞在宏观设计域中的非均匀分布,充分发挥了材料潜力,提升了点阵结构的力学性能。

Description

一种多构型点阵结构的拓扑优化方法
技术领域
本发明属于拓扑结构优化技术领域,更具体地,涉及一种多构型点阵 结构的拓扑优化方法。
背景技术
点阵结构具有高的比刚度/强度、减振吸能、隔热防热、抗冲击能力强 等特性。此外,点阵结构的微结构以杆件为主,可设计性强,承载性能优 越,因此广泛应用于航空航天、生物医学、自动化建筑工程等领域。在点 阵结构设计中考虑多种不同构型的点阵单胞,可以更大程度地发挥点阵结 构的设计潜能,使得点阵结构在给定的载荷条件下实现更好的力学性能。
针对多构型点阵结构的拓扑优化,本领域相关技术人员已做了一些研 究,如文献:“C Wang,J H Zhu,W H Zhang,et al.Concurrent topology optimization design ofstructures and non-uniform parameterized lattice microstructures.Structuraland Multidisciplinary Optimization,2018,58(1): 35-50.”通过引入额外的参数来定义微结构单胞构型,从而避免拓扑优化过 程中高昂的数值均匀化迭代计算成本,也使得结构设计的建模变得更容易, 最终可以得到多种微结构非均匀分布的宏观构型。
然而,该方法的微结构单胞的拓扑构型及其数量受到限制,并不能耦 合更多种拓扑构型差异较大的点阵单胞。
发明内容
针对现有技术的以上缺陷或改进需求,本发明提供了一种多构型点阵 结构的拓扑优化方法,其目的在于,通过基于水平集函数构建不同构型的 点阵单胞,针对每一种构型的基础点阵单胞获得相同等效密度区间中的一 系列点阵单胞及其对应的等效弹性张量,建立多构型点阵结构每个有限单 元对应的弹性张量关于单元相对密度的插值函数以构建多构型点阵结构的 拓扑优化模型,采用交替活动相算法将拓扑优化问题划分为多个子问题并 进行求解,实现拓扑优化过程,由此解决现有微结构单胞的拓扑构型及其 数量受到限制导致点阵结构设计空间受限的技术问题。
为实现上述目的,按照本发明的一个方面,提供了一种多构型点阵结 构的拓扑优化方法,包括:
S1:基于水平集函数描述多种构型点阵的单胞,并在相同等效密度区 间中获取各个所述构型基础点阵单胞对应的一系列点阵单胞,并计算其等 效弹性张量;
S2:建立多构型点阵结构中每个有限单元对应的弹性张量关于单元相 对密度的插值函数,利用所述插值函数构建多构型点阵结构的拓扑优化模 型;
S3:采用交替活动相算法将所述拓扑优化问题分解为多个子问题并对 其求解,从而对点阵结构的宏观布局进行优化,得到目标多构型点阵结构。
在其中一个实施例中,所述S1包括:
S11:利用一根桁架的水平集函数
Figure BDA0003636715630000021
描述各种 构型点阵的单胞;φ3D(x)=max(φ3D,c(x,y,z),φ3D,s1(x,y,z),φ3D,s2(x,y,z));D表 示一个固定的欧拉参考空间,x表示在空间D内点的坐标,
Figure BDA0003636715630000022
表示所述一 根桁架的结构边界,Ω是所述一根桁架所占据的空间,φ3D,c(x,y,z)、 φ3D,s1(x,y,z)和φ3D,s2(x,y,z)分布表示所述一根桁架所包含的圆柱体、第一球体 和第二球体的水平集函数;所述第一球体的球心和所述第二球体的球心分 别与所述圆柱体的两个底面的圆心相互重合;
S12:针对每一种构型的点阵单胞,利用形状插值技术获得相同等效密 度区间中的其他点阵单胞;
S13:计算所有所述点阵单胞的等效弹性张量。
在其中一个实施例中,所述S12包括:针对每一种构型的基础点阵单 胞,利用公式
Figure BDA0003636715630000031
计算相同等效密度区间中的其他点阵单胞的水 平集函数;φpro是每种构型点阵单胞的初始水平集函数,对应等效密度为0.01 的点阵单胞;
Figure BDA0003636715630000032
是插值系数矩阵,
Figure BDA0003636715630000033
在其中一个实施例中,所述S13包括:利用公式
Figure BDA0003636715630000034
计算各个点阵单胞的等效弹性张量; 其中,
Figure BDA0003636715630000035
表示等效弹性张量矩阵中的分量,Y表示点阵单胞的体积,i,j,k,l 是索引向量;
Figure BDA0003636715630000036
代表的是对应于单元测试应变场
Figure BDA0003636715630000037
的单元位移解,ke为单元刚度矩阵;DH表示点阵单胞的等效弹性张量矩阵,三维情况下的DH的表现形式为:
Figure BDA0003636715630000038
在其中一个实施例中,所述S2包括:
S21:建立多构型点阵结构每个有限单元对应的弹性张量关于单元相对 密度的插值函数;
S22:利用所述插值函数计算宏观结构整体的刚度矩阵K,利用所述刚 度矩阵K构建多构型点阵结构的拓扑优化模型,其表达式为:
Figure BDA0003636715630000041
其中,
Figure BDA0003636715630000042
表示设计变量,是指宏观设计域内第i种构型在第j个有限单 元中的相对密度;C为宏观结构整体的柔度值,F是外部施加的载荷矩阵, K为宏观结构整体的刚度矩阵,U是施加载荷后引起的位移矩阵,Ω表示整 个宏观设计域,V*表示宏观结构整体的体积分数约束,li和ui分别表示设计 变量的上下限,为xmin到1之间的某个值,xmin取0.001。
在其中一个实施例中,所述S21包括:定义单元相对密度的插值函数 为:
Figure BDA0003636715630000043
(Dj)k表示宏观设计域中每个有限单元惩罚后的弹 性张量矩阵,Di表示第i种构型的等效弹性张量矩阵,p为惩罚因子。
在其中一个实施例中,所述S22中:体积分数约束为:V*=1,即所述 目标多构型点阵结构的点阵单胞完全填充整个宏观设计域。
在其中一个实施例中,所述S3包括:采用交替活动相算法开展拓扑优 化设计,将多种构型点阵结构的拓扑优化问题转化为两种点阵单胞构型的 拓扑优化子问题,迭代过程分为内部迭代和外部迭代两个部分;其中,内 部迭代的拓扑优化模型可以表示为:
Figure BDA0003636715630000051
其中,
Figure BDA0003636715630000052
表示的是第a种构型在第j个有限单元内的相对密度,即内部 迭代的主设计变量,va表示的是第a种构型在整个宏观设计域内的体积分数, 为预先给定的参数。
在其中一个实施例中,所述多种构型的点阵属于桁架点阵和/或TPMS 点阵。
总体而言,通过本发明所构思的以上技术方案与现有技术相比,能够 取得下列有益效果:
(1)本发明提供了一种多构型点阵结构的拓扑优化方法,通过基于水 平集函数构建不同构型的点阵单胞,针对每一种构型的基础点阵单胞获得 相同等效密度区间中的其他一系列点阵单胞及其对应的等效弹性张量,建 立多构型点阵结构每个有限单元对应的弹性张量关于单元相对密度的插值 函数以构建多构型点阵结构的拓扑优化模型,采用交替活动相算法将拓扑 优化问题划分为多个子问题并进行求解,实现拓扑优化过程,由此解决现 有微结构单胞的拓扑构型及其数量受到限制导致点阵结构设计空间受限的 技术问题。本发明实现了多种拓扑构型不同的点阵单胞在宏观设计域中的 非均匀分布,充分发挥了材料潜力,提升了点阵结构的力学性能。
(2)本发明由于基于多相材料的插值模型进行拓扑优化设计,可以实 现多种拓扑构型差异较大的点阵结构在宏观设计域中的合理分布。
(3)本发明由于采用交替活动相算法进行求解,在子问题中涉及的设 计变量和体积约束更少,可以有效地降低拓扑优化计算成本。
(4)本发明在预定义点阵单胞时已充分考虑了单胞之间的连接性,通 过预设连接点进行微结构之间的连接,保证了点阵微结构单胞之间良好的 连接性。
附图说明
图1是本发明一实施例提供的一种多构型点阵结构的拓扑优化方法的 流程图;
图2是本发明一实施例提供的棱边立方点阵、面心立方点阵、棱边体 心立方点阵和体心立方点阵的构型示意图;
图3是本发明一实施例提供的任意一种点阵结构中一根桁架所包含的 一个圆柱体、两个球体之间的几何关系示意图;
图4a是本发明一实施例提供的棱边立方点阵的等效密度及其对应的等 效弹性张量示意图;
图4b是本发明一实施例提供的面心立方点阵的等效密度及其对应的等 效弹性张量示意图;
图4c是本发明一实施例提供的棱边体心立方点阵的等效密度及其对应 的等效弹性张量示意图;
图4d是本发明一实施例提供的体心立方点阵的等效密度及其对应的等 效弹性张量示意图;
图5是本发明一实施例提供的四种点阵结构的等效杨氏模量对比图;
图6是本发明一实施例中构建的宏观点阵结构设计域、载荷及边界条 件示意图;
图7是本发明一实施例中构建的图6中点阵结构优化后的多构型点阵 结构示意图;
图8是本发明一实施例中构建的图6中点阵结构优化后四种点阵单胞 独立分布的示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图 及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体 实施例仅仅用以解释本发明,并不用于限定本发明。此外,下面所描述的 本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可 以相互组合。
本发明提供了一种多构型点阵结构的拓扑优化方法,包括:
S1:基于水平集函数描述多种构型点阵的单胞,并在相同等效密度区 间中获取各个所述构型基础点阵单胞对应的一系列点阵单胞,并计算其等 效弹性张量;
S2:建立多构型点阵结构中每个有限单元对应的弹性张量关于单元相 对密度的插值函数,利用所述插值函数构建多构型点阵结构的拓扑优化模 型;
S3:采用交替活动相算法将所述拓扑优化问题分解为多个子问题并对 其求解,从而对点阵结构的宏观布局进行优化,得到目标多构型点阵结构。
在其中一个实施例中,所述S1包括:
S11:利用一根桁架的水平集函数
Figure BDA0003636715630000071
描述各种 构型点阵的单胞;φ3D(x)=max(φ3D,c(x,y,z),φ3D,s1(x,y,z),φ3D,s2(x,y,z));D表 示一个固定的欧拉参考空间,x表示在空间D内点的坐标,
Figure BDA0003636715630000072
表示所述一 根桁架的结构边界,Ω是所述一根桁架所占据的空间,φ3D,c(x,y,z)、 φ3D,s1(x,y,z)和φ3D,s2(x,y,z)分布表示所述一根桁架所包含的圆柱体、第一球体 和第二球体的水平集函数;所述第一球体的球心和所述第二球体的球心分 别与所述圆柱体的两个底面的圆心相互重合;
S12:针对每一种构型的点阵单胞,利用形状插值技术获得相同等效密 度区间中的其他点阵单胞;
S13:计算所有所述点阵单胞的等效弹性张量。
步骤1、使用水平集函数描述点阵单胞的拓扑构型,构建多种构型不同 的点阵单胞,针对每一种构型的基础点阵单胞使用形状插值技术获得相同 等效密度区间中的一系列点阵单胞,并计算所有获得的点阵单胞的等效弹 性张量。
具体的,使用水平集函数描述点阵单胞的拓扑构型,并将等效密度值 为0.01的点阵作为基础点阵,并对该基础点阵使用形状插值技术获得一系 列的样本点阵,使用数值均匀化法求解所有点阵单胞的等效弹性张量。
其中,所述点阵单胞的构型为棱边立方点阵、面心立方点阵、棱边体 心立方点阵和体心立方点阵,以体心立方点阵为例,其包含4根桁架,所 述体心立方点阵的水平集函数为:
Figure BDA0003636715630000081
其中,φs(x)=max(φi),φi=φ3D,i(x),i=1,2,3,4.φ3D,i是所述四根桁架的水平集函数,x表示在空间D内的点的坐标,D是一个固定的欧拉参考空间,
Figure BDA0003636715630000082
表示所述体心立方点阵的结构边界,Ωs是所述体心立方点阵所占据的空 间,所述Ωs满足Ωs=Ω1∪Ω2∪Ω3∪Ω4,Ω1234分别是四根桁架所占据的 空间。
所述任意一种点阵中,一根桁架的水平集函数为:
Figure BDA0003636715630000083
φ3D(x)=max(φ3D,c(x,y,z),φ3D,s1(x,y,z),φ3D,s2(x,y,z)),
φ3D,c(x,y,z)=min(φ3D,c1(x,y,z),φ3D,c2(x,y,z)),
φ3D,c1(x,y,z)=(L3D/2)2-(cosθ3D·Ld)23D,c2(x,y,z)=(t3D/2)2-(sinθ3D·Ld)2,
Figure BDA0003636715630000091
Figure BDA0003636715630000092
dx=x-x0,dy=y-y0,dz=z-z0,
Figure BDA0003636715630000093
φ3D,s1(x,y,z)=(t3D/2)2-(x-x1)2+(y-y1)2+(z-z1)2,
φ3D,s2(x,y,z)=(t3D/2)2-(x-x2)2+(y-y2)2+(z-z2)2
其中,D是一个固定的欧拉参考空间,x表示在空间D内的点的坐标,
Figure BDA0003636715630000099
表示所述一根桁架的结构边界,Ω是所述一根桁架所占据的空间, φ3D,c(x,y,z)、φ3D,s1(x,y,z)和φ3D,s2(x,y,z)分别表示一根桁架所包含的一个圆 柱体、两个球体的水平集函数、所述两个球体的球心与所述圆柱体的两个 底面的圆心相互重合,(x1,y1,z1)和(x2,y2,z2)分别表示所述两个球体的球心 坐标,(x0,y0,z0)表示所述两个球心连线的中点坐标,t3D和L3D分别表示所述 圆柱体底面的直径和圆柱体的长度。
本实施方式中,所述样本点阵的数量为50个,且样本点阵的等效密度 呈等差数列,所述样本点阵的密度值范围为[0.01,1],形状插值技术对应的模 型为:
Figure BDA0003636715630000094
其中,φe是希望获得的特定等效密度下的每种构型 点阵单胞的水平集函数,φpro是所述每种构型点阵单胞的初始水平集函数, 对应等效密度为0.01的点阵单胞,
Figure BDA0003636715630000095
是插值系数矩阵,
Figure BDA0003636715630000096
的取值范围为
Figure BDA0003636715630000097
的值可以通过二分法计算得到。
在其中一个实施例中,所述S12包括:针对每一种构型的点阵单胞, 利用公式
Figure BDA0003636715630000098
计算相同等效密度区间中的点阵单胞的水平集函 数;φpro是每种构型点阵单胞的初始水平集函数,对应等效密度为0.01的点 阵单胞;
Figure BDA0003636715630000101
是插值系数矩阵,
Figure BDA0003636715630000102
在其中一个实施例中,所述S13包括:
利用公式
Figure BDA0003636715630000103
计算各个点阵单胞的 等效弹性张量;其中,
Figure BDA0003636715630000104
表示等效弹性张量矩阵中的分量,Y表示点阵单 胞的体积,i,j,k,l是索引向量;
Figure BDA0003636715630000105
代表的是对应于单元测试应变场
Figure BDA0003636715630000106
的单元位移解,ke为单元刚度矩阵;DH表示点阵单胞的等效弹性张量矩阵, 三维情况下的DH的表现形式为:
Figure BDA0003636715630000107
在其中一个实施例中,所述S2包括:S21:建立多构型点阵结构每个 有限单元对应的弹性张量关于单元相对密度的插值函数;S22:利用所述插 值函数计算宏观结构整体的刚度矩阵K,利用所述刚度矩阵K构建多构型 点阵结构的拓扑优化模型,其表达式为:
Figure BDA0003636715630000108
其中,
Figure BDA0003636715630000111
表示设计变量,是指宏观设计域内第i种构型在第j个有限单 元中的相对密度;C为宏观结构整体的柔度值,F是外部施加的载荷矩阵, K为宏观结构整体的刚度矩阵,U是施加载荷后引起的位移矩阵,Ω表示整 个宏观设计域,V*表示宏观结构整体的体积分数约束,li和ui分别表示设计 变量的上下限,为xmin到1之间的某个值,xmin通常取0.001。
在其中一个实施例中,所述S21包括:定义单元相对密度的插值函数 为:
Figure BDA0003636715630000112
(Dj)k表示宏观设计域中每个有限单元惩罚后的弹 性张量矩阵,Di表示第i种构型的等效弹性张量矩阵,p为惩罚因子。
在其中一个实施例中,所述S22中:体积分数约束为:V*=1,即所述 目标多构型点阵结构的点阵单胞完全填充整个宏观设计域。
在其中一个实施例中,所述S3包括:采用交替活动相算法开展拓扑优 化设计,将多种构型点阵结构的拓扑优化问题转化为两种点阵单胞构型的 拓扑优化子问题,迭代过程分为内部迭代和外部迭代两个部分;其中,内 部迭代的拓扑优化模型可以表示为:
Figure BDA0003636715630000113
其中,
Figure BDA0003636715630000114
表示的是第a种构型在第j个有限单元内的相对密度,即内部 迭代的主设计变量,va表示的是第a种构型在整个宏观设计域内的体积分数, 为预先给定的参数。
在其中一个实施例中,所述多种构型的点阵至少包括两种以上的点阵 单胞,所述多种构型的点阵可以属于桁架点阵和/或TPMS点阵,举例来说, 棱边立方点阵、面心立方点阵、体心立方点阵等桁架点阵。需要说明的是, 还可以为其他点阵,此处不做限制。
举例来说,请参阅图1,本发明提供的多构型点阵结构拓扑优化方法, 所述优化设计方法主要包括以下步骤:
本实施例中的待优化的点阵结构,其设计域、载荷及边界条件如图6 所示,本发明中的密度是指点阵单胞中实体部分体积所占该点阵单胞单元 体积的体积率,本实施例中的优化目标设置为点阵结构的柔度值最小,优 化约束为等效密度为0.2的体心立方点阵、等效密度为0.4的棱边体心立方 点阵、等效密度为0.6的面心立方点阵以及等效密度为0.8的棱边立方点阵 的体积分数分别为0.4、0.1、0.1和0.4,总体积约束为0.5。所有构型的点阵单胞均采用同一种材料,弹性模量为E=2750MPa,泊松比μ=0.38。
如图1所示,本发明一种多构型点阵结构拓扑优化方法包括如下步骤:
步骤一,基于水平集函数构建不同构型的点阵单胞,并采用数值均匀 化法计算点阵单胞的等效弹性张量,具体包括以下子步骤:
(1.1)使用水平集函数描述点阵单胞的拓扑构型,构建多种构型不同 的点阵单胞,此处点阵单胞的类型为棱边立方点阵、面心立方点阵、棱边 体心立方点阵和体心立方点阵。以体心立方点阵为例,其包含4根桁架, 所述体心立方点阵的水平集函数为:
Figure BDA0003636715630000121
其中,φs(x)=max(φi),φi=φ3D,i(x),i=1,2,3,4.φ3D,i是所述四根桁架的水平集函数,x表示在空间D内的点的坐标,D是一个固定的欧拉参考空间,
Figure BDA0003636715630000122
表示所述体心立方点阵的结构边界,Ωs是所述体心立方点阵所占据的空 间,所述Ωs满足Ωs=Ω1∪Ω2∪Ω3∪Ω4,Ω1234分别是四根桁架所占据的 空间。
所述任意一种点阵中,一根桁架的水平集函数如下:
Figure BDA0003636715630000131
φ3D(x)=max(φ3D,c(x,y,z),φ3D,s1(x,y,z),φ3D,s2(x,y,z)),
φ3D,c(x,y,z)=min(φ3D,c1(x,y,z),φ3D,c2(x,y,z)),
φ3D,c1(x,y,z)=(L3D/2)2-(cosθ3D·Ld)23D,c2(x,y,z)=(t3D/2)2-(sinθ3D·Ld)2,
Figure BDA0003636715630000132
Figure BDA0003636715630000133
dx=x-x0,dy=y-y0,dz=z-z0,
Figure BDA0003636715630000134
φ3D,s1(x,y,z)=(t3D/2)2-(x-x1)2+(y-y1)2+(z-z1)2,
φ3D,s2(x,y,z)=(t3D/2)2-(x-x2)2+(y-y2)2+(z-z2)2
其中,D是一个固定的欧拉参考空间,x表示在空间D内的点的坐标,
Figure BDA0003636715630000135
表示所述一根桁架的结构边界,Ω是所述一根桁架所占据的空间, φ3D,c(x,y,z)、φ3D,s1(x,y,z)和φ3D,s2(x,y,z)分别表示一根桁架所包含的一个圆 柱体、两个球体的水平集函数、所述两个球体的球心与所述圆柱体的两个 底面的圆心相互重合,(x1,y1,z1)和(x2,y2,z2)分别表示所述两个球体的球心 坐标,(x0,y0,z0)表示所述两个球心连线的中点坐标,t3D和L3D分别表示所述 圆柱体的底面直径和圆柱体的长度。
(1.2)将等效密度值为0.01的点阵作为基础点阵,其水平集函数为 φpro(x),并对该基础点阵使用形状插值技术,以获得四个系列的样本点阵, 样本点阵的等效密度呈等差数列排列,样本点阵的密度值范围为[0.01,1],使 用数值均匀化法计算所有获得的点阵的等效弹性张量。
(1.3)若等效弹性张量满足以下形式:
Figure BDA0003636715630000141
则点阵单胞的等效杨氏模量可依据下式进行计算:
Figure BDA0003636715630000142
其中,DH表示点阵单胞的等效弹性张量矩阵,D11、D12和D44是矩阵中 三个互相独立的等效弹性常数,均可通过数值均匀化法计算得到。可根据 此计算式对点阵单胞的力学性能进行对比,确定强弱材料。
步骤二,建立多构型点阵结构等效弹性张量关于单元相对密度的插值 函数,构建多构型点阵结构拓扑优化模型,具体包括以下子步骤:
(2.1)建立多构型点阵结构每个有限单元对应的弹性张量关于单元相 对密度的插值函数:
Figure BDA0003636715630000143
其中,(Dj)k表示宏观设计 域中每个有限单元惩罚后的弹性张量矩阵,
Figure BDA0003636715630000144
表示的是设计变量,是指宏 观设计域内第i种构型在第j个有限单元中的相对密度,Di表示第i种构型 的等效弹性张量矩阵,p为惩罚因子,一般取p=3。
(2.2)构建多构型点阵结构拓扑优化模型:
Figure BDA0003636715630000151
其中,C为宏观结构整体的柔度值,F是外部施加的载荷矩阵,K为宏 观结构整体的刚度矩阵,U是施加载荷后引起的位移矩阵,Ω表示整个宏观 设计域,V*表示宏观结构整体的体积分数约束,li和ui分别表示设计变量的 上下限,为xmin到1之间的某个值,xmin是为了避免刚度矩阵出现奇异现象, 通常取0.001。
步骤三,采用交替活动相算法将问题分解为一系列子问题,并进行求 解,得到宏观结构中每个有限单元的最终点阵构型,从而实现拓扑优化过 程,获得多构型点阵结构,具体包括以下子步骤:
(3.1)使用有限元的思想将点阵结构离散化为60×12×10=7200个空间八 节点单元,施加载荷并设置边界条件,以点阵结构柔度最小为目标,采用 交替活动相算法开展拓扑优化设计,将多种构型点阵结构的拓扑优化问题 转化为一系列两种点阵单胞构型的拓扑优化子问题,迭代过程分为内部迭 代和外部迭代两个部分。其中,内部迭代的拓扑优化模型可以表示为:
Figure BDA0003636715630000152
其中,
Figure BDA0003636715630000161
表示的是第a种构型在第j个有限单元内的相对密度,即内部 迭代的主设计变量,va表示的是第a种构型在整个宏观设计域内的体积分数, 为预先给定的参数。
(3.2)假设内部迭代中更新的设计变量为第a种构型和第b种构型在每 个有限单元中的相对密度,由于固定了其他m-2种构型的分布保持不变, 则内部迭代的总设计变量可表示为:
Figure BDA0003636715630000162
在内部迭代中,对设计变量
Figure BDA0003636715630000163
进行迭代更新之后,设计变量
Figure BDA0003636715630000164
也可通过 下式进行更新:
Figure BDA0003636715630000165
计算目标函数及约束条件对设计变量的敏度信息,目标函数的敏度计 算公式如下:
Figure BDA0003636715630000166
其中,
Figure BDA0003636715630000167
是宏观点阵结构的柔度对设计变量
Figure BDA0003636715630000168
的灵敏度,B代表的 是应变-位移矩阵,也可称为几何矩阵,刚度矩阵K可通过下式进行计算:
Figure BDA0003636715630000169
约束条件的敏度计算公式如下:
Figure BDA00036367156300001610
(3.3)在计算敏度信息时利用邻近有限元单元敏度信息对当前有限元 单元的敏度信息进行过滤,避免棋盘格、网格依赖性等数值不稳定现象的 发生。
(3.4)采用优化准则法(Optimality Criteria,OC)更新设计变量,根据 优化结果判断目标函数是否满足设定收敛条件,若满足收敛条件,则输出 当前点阵结构每个有限单元的构型信息,否则继续执行步骤(3.2)。
以下以三维简支梁约束的点阵结构的设计来进一步说明本发明。如图2 所示是本发明所构建的棱边立方点阵、面心立方点阵、棱边体心立方点阵 和体心立方点阵的构型示意图,图3所示是本发明所构建的任意一种点阵 结构中一根桁架所包含的一个圆柱体、两个球体之间的几何关系示意图。
如图4a、图4b、图4c和图4d分别是本发明所构建的棱边立方点阵、 面心立方点阵、棱边体心立方点阵和体心立方点阵的等效密度及其对应的 等效弹性张量示意图。
图5所示是本发明所构建的四种点阵结构的等效杨氏模量对比图,可 以看到,在设定的点阵单胞等效密度的范围内,低密度区域棱边立方单胞 的等效杨氏模量最大,面心立方与棱边体心立方单胞的相差不大,而体心 立方单胞的最小;在中密度区域,四种单胞的等效杨氏模量值差别更为明 显;而在高密度区域,四种单胞的等效杨氏模量值非常接近。在整体趋势 上,等效杨氏模量由大到小的顺序应为:棱边立方单胞、面心立方单胞、 棱边体心立方单胞、体心立方单胞。所以在构型预定义时,本实施例将以 等效密度为0.2的体心立方单胞、等效密度为0.4的棱边体心立方单胞、等 效密度为0.6的面心立方单胞以及等效密度为0.8的棱边立方单胞作为后续 拓扑优化设计的微结构材料。
如图6所示,简支约束的点阵结构设计域长L=0.3m,宽W=0.05m,高 H=0.06m。分布载荷施加在结构的上表面中线处,大小为q=4×105N/m,方 向为竖直向下。设计域的下表面左边界线处由固定铰支座进行约束,所有 方向上的自由度被完全约束;下表面右边界线处由滚动铰支座进行约束, 该点竖直方向的自由度被约束。宏观设计域被离散化为60×10×12个空间八 节点单元。优化约束为等效密度为0.2的体心立方单胞、等效密度为0.4的棱边体心立方单胞、等效密度为0.6的面心立方单胞以及等效密度为0.8的 棱边立方单胞的体积分数分别为0.4、0.1、0.1和0.4,总体积约束为0.5, 优化目标是点阵结构的柔度值最小。
如图7所示是优化后的多构型点阵结构示意图,图8所示是优化后四 种点阵单胞独立分布的示意图。可以明显看到,红色主材料是分布在主要 的传力路径上,对整个结构起到了很好的支撑作用,而剩余的三种材料则 主要分布在宏观结构的内部或未承载/未约束的边界处。所有点阵材料的分 布情况与其对应点阵单胞的力学性能表现是一致的。
本领域的技术人员容易理解,以上所述仅为本发明的较佳实施例而已, 并不用以限制本发明,凡在本发明的精神和原则之内所作的任何修改、等 同替换和改进等,均应包含在本发明的保护范围之内。

Claims (9)

1.一种多构型点阵结构的拓扑优化方法,其特征在于,包括:
S1:基于水平集函数描述多种构型点阵的单胞,并在相同等效密度区间中获取各个所述构型基础点阵单胞对应的一系列点阵单胞,并计算其等效弹性张量;
S2:建立多构型点阵结构中每个有限单元对应的弹性张量关于单元相对密度的插值函数,利用所述插值函数构建多构型点阵结构的拓扑优化模型;
S3:采用交替活动相算法将所述拓扑优化问题分解为多个子问题并对其求解,从而对点阵结构的宏观布局进行优化,得到目标多构型点阵结构。
2.如权利要求1所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S1包括:
S11:利用一根桁架的水平集函数
Figure FDA0003636715620000011
描述各种构型点阵的单胞;φ3D(x)=max(φ3D,c(x,y,z),φ3D,s1(x,y,z),φ3D,s2(x,y,z));D表示一个固定的欧拉参考空间,x表示在空间D内点的坐标,
Figure FDA0003636715620000012
表示所述一根桁架的结构边界,Ω是所述一根桁架所占据的空间,φ3D,c(x,y,z)、φ3D,s1(x,y,z)和φ3D,s2(x,y,z)分布表示所述一根桁架所包含的圆柱体、第一球体和第二球体的水平集函数;所述第一球体的球心和所述第二球体的球心分别与所述圆柱体的两个底面的圆心相互重合;
S12:针对每一种构型的点阵单胞,利用形状插值技术获得相同等效密度区间中的其他点阵单胞;
S13:计算所有所述点阵单胞的等效弹性张量。
3.如权利要求2所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S12包括:针对每一种构型的基础点阵单胞,利用公式
Figure FDA0003636715620000013
计算相同等效密度区间中其他的点阵单胞的水平集函数;φpro是每种构型点阵单胞的初始水平集函数,对应等效密度为0.01的点阵单胞;
Figure FDA0003636715620000021
是插值系数矩阵,
Figure FDA0003636715620000022
4.如权利要求2所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S13包括:利用公式
Figure FDA0003636715620000023
计算各个点阵单胞的等效弹性张量;其中,
Figure FDA0003636715620000024
表示等效弹性张量矩阵中的分量,Y表示点阵单胞的体积,i,j,k,l是索引向量;
Figure FDA0003636715620000025
代表的是对应于单元测试应变场
Figure FDA0003636715620000026
的单元位移解,ke为单元刚度矩阵;DH表示点阵单胞的等效弹性张量矩阵,三维情况下的DH的表现形式为:
Figure FDA0003636715620000027
5.如权利要求1所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S2包括:
S21:建立多构型点阵结构每个有限单元对应的弹性张量关于单元相对密度的插值函数;
S22:利用所述插值函数计算宏观结构整体的刚度矩阵K,利用所述刚度矩阵K构建多构型点阵结构的拓扑优化模型,其表达式为:
Figure FDA0003636715620000031
Figure FDA0003636715620000032
Subject to:F=KU
Figure FDA0003636715620000033
其中,
Figure FDA0003636715620000034
表示设计变量,是指宏观设计域内第i种构型在第j个有限单元中的相对密度;C为宏观结构整体的柔度值,F是外部施加的载荷矩阵,K为宏观结构整体的刚度矩阵,U是施加载荷后引起的位移矩阵,Ω表示整个宏观设计域,V*表示宏观结构整体的体积分数约束,li和ui分别表示设计变量的上下限,为xmin到1之间的某个值,xmin取0.001。
6.如权利要求5所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S21包括:定义单元相对密度的插值函数为:
Figure FDA0003636715620000035
(Dj)k表示宏观设计域中每个有限单元惩罚后的弹性张量矩阵,Di表示第i种构型的等效弹性张量矩阵,p为惩罚因子。
7.如权利要求5所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S22中:体积分数约束为:V*=1,即所述目标多构型点阵结构的点阵单胞完全填充整个宏观设计域。
8.如权利要求5所述的多构型点阵结构的拓扑优化方法,其特征在于,所述S3包括:采用交替活动相算法开展拓扑优化设计,将多种构型点阵结构的拓扑优化问题转化为两种点阵单胞构型的拓扑优化子问题,迭代过程分为内部迭代和外部迭代两个部分;其中,内部迭代的拓扑优化模型可以表示为:
Figure FDA0003636715620000041
Figure FDA0003636715620000042
Subject to:F=KU
Figure FDA0003636715620000043
Figure FDA0003636715620000044
0<xmin≤la≤ua≤1
其中,
Figure FDA0003636715620000045
表示的是第a种构型在第j个有限单元内的相对密度,即内部迭代的主设计变量,va表示的是第a种构型在整个宏观设计域内的体积分数,为预先给定的参数。
9.如权利要求1-8任一项所述的多构型点阵结构的拓扑优化方法,其特征在于,所述多种构型的点阵属于桁架点阵和/或TPMS点阵。
CN202210507811.5A 2022-05-10 一种多构型点阵结构的拓扑优化方法 Active CN114999591B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210507811.5A CN114999591B (zh) 2022-05-10 一种多构型点阵结构的拓扑优化方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210507811.5A CN114999591B (zh) 2022-05-10 一种多构型点阵结构的拓扑优化方法

Publications (2)

Publication Number Publication Date
CN114999591A true CN114999591A (zh) 2022-09-02
CN114999591B CN114999591B (zh) 2024-07-16

Family

ID=

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116579151A (zh) * 2023-05-05 2023-08-11 大连理工大学 一种基于mmc框架的非均匀点阵结构优化设计方法
GB2619567A (en) * 2022-06-09 2023-12-13 Ocado Innovation Ltd 3-D lattice optimisation
CN117473836A (zh) * 2023-11-16 2024-01-30 北京理工大学 一种薄壁-多类点阵填充结构一体化设计方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107025340A (zh) * 2017-03-30 2017-08-08 华中科技大学 一种适用于增材制造的自支撑网状结构拓扑优化设计方法
CN112182929A (zh) * 2020-09-18 2021-01-05 北京航空航天大学 一种考虑尺寸控制的多孔材料跨尺度可靠性拓扑优化方法
CN113345536A (zh) * 2021-05-31 2021-09-03 山东大学 一种基于极限各向异性点阵材料的结构拓扑优化方法
CN114254408A (zh) * 2021-12-17 2022-03-29 华中科技大学 一种基于代理模型的梯度点阵等几何拓扑优化方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107025340A (zh) * 2017-03-30 2017-08-08 华中科技大学 一种适用于增材制造的自支撑网状结构拓扑优化设计方法
CN112182929A (zh) * 2020-09-18 2021-01-05 北京航空航天大学 一种考虑尺寸控制的多孔材料跨尺度可靠性拓扑优化方法
CN113345536A (zh) * 2021-05-31 2021-09-03 山东大学 一种基于极限各向异性点阵材料的结构拓扑优化方法
CN114254408A (zh) * 2021-12-17 2022-03-29 华中科技大学 一种基于代理模型的梯度点阵等几何拓扑优化方法

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2619567A (en) * 2022-06-09 2023-12-13 Ocado Innovation Ltd 3-D lattice optimisation
GB2619567B (en) * 2022-06-09 2024-07-17 Ocado Innovation Ltd 3-D lattice optimisation
CN116579151A (zh) * 2023-05-05 2023-08-11 大连理工大学 一种基于mmc框架的非均匀点阵结构优化设计方法
CN116579151B (zh) * 2023-05-05 2024-01-30 大连理工大学 一种基于mmc框架的非均匀点阵结构优化设计方法
CN117473836A (zh) * 2023-11-16 2024-01-30 北京理工大学 一种薄壁-多类点阵填充结构一体化设计方法

Similar Documents

Publication Publication Date Title
Zhang et al. Topology optimization for concurrent design of layer-wise graded lattice materials and structures
CN110795873B (zh) 一种考虑尺寸控制的跨尺度拓扑优化方法
Kaveh et al. Size optimization of space trusses using Big Bang–Big Crunch algorithm
Zhang et al. Concurrent topology optimization for cellular structures with nonuniform microstructures based on the kriging metamodel
Do et al. A modified symbiotic organisms search (mSOS) algorithm for optimization of pin-jointed structures
Zhou et al. Design of graded two-phase microstructures for tailored elasticity gradients
CN112836411B (zh) 加筋板壳结构的优化方法、装置、计算机设备和存储介质
CN112395700B (zh) 一种代理模型驱动的梯度点阵夹层结构优化方法
CN112182929A (zh) 一种考虑尺寸控制的多孔材料跨尺度可靠性拓扑优化方法
CN110941924B (zh) 一种多组件系统集成一体化的多尺度拓扑优化设计方法
CN111241738A (zh) 一种考虑破损-安全条件的连续体位移与频率约束拓扑优化设计方法
CN111950149A (zh) 基于参数化水平集法的连续体结构非概率拓扑优化方法
CN111027110A (zh) 一种连续体结构拓扑与形状尺寸综合优化方法
CN108647405A (zh) 多层级点阵结构拓扑优化设计的子结构插值模型建模方法
CN113345536B (zh) 一种基于极限各向异性点阵材料的结构拓扑优化方法
CN110955938B (zh) 一种具有梯度多孔夹芯的夹层结构拓扑优化方法
CN114254408A (zh) 一种基于代理模型的梯度点阵等几何拓扑优化方法
CN112446163A (zh) 基于参数化水平集的能量有限元拓扑优化方法
CN115203997A (zh) 一种基于多变量设计的点阵-实体复合结构拓扑优化方法
CN110751729A (zh) 基于腐蚀-扩散算子的准周期层级结构拓扑优化方法
CN116341179B (zh) 五模超材料骨支架的多目标等几何多尺度拓扑优化方法
CN110717208B (zh) 一种基于连续梯度微结构的多尺度频率响应拓扑优化方法
CN114999591A (zh) 一种多构型点阵结构的拓扑优化方法
CN114999591B (zh) 一种多构型点阵结构的拓扑优化方法
CN109299499B (zh) 考虑修正因子的多步骤结构优化设计方法及飞行器

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