CN106529044B - 一种柔顺机构0-1变构型拓扑图提取方法 - Google Patents
一种柔顺机构0-1变构型拓扑图提取方法 Download PDFInfo
- Publication number
- CN106529044B CN106529044B CN201611002663.2A CN201611002663A CN106529044B CN 106529044 B CN106529044 B CN 106529044B CN 201611002663 A CN201611002663 A CN 201611002663A CN 106529044 B CN106529044 B CN 106529044B
- Authority
- CN
- China
- Prior art keywords
- compliant mechanism
- density
- topological
- optimization model
- cell
- 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
Links
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
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/02—CAD in a network environment, e.g. collaborative CAD or distributed simulation
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)
- Complex Calculations (AREA)
- Design And Manufacture Of Integrated Circuits (AREA)
Abstract
本发明公开了一种柔顺机构0‑1变构型拓扑图提取方法,包括以下步骤:首先建立柔顺机构第一层拓扑优化模型;接着,对拓扑优化模型进行优化,提取出含有中间单元的拓扑图;然后,对含有中间单元的柔顺机构拓扑优化模型做进一步优化,提取出柔顺机构0‑1拓扑图;之后,建立下一层柔顺机构拓扑优化模型,并对其进行优化,通过不断增大拓扑优化模型体积约束值的方式逐渐逼近最佳解,最终提取出满足目标体积的理想的柔顺机构0‑1变构型拓扑图。本发明采用体积约束延拓方法,通过不断增大拓扑优化模型体积约束值的方式逐渐逼近最佳解,能最大限度地保留所有可能单元,提取出柔顺机构0‑1变构型拓扑图。
Description
技术领域
本发明属于柔顺机构拓扑优化设计中的拓扑图提取技术领域,特别涉及一种柔顺机构0-1变构型拓扑图提取方法。
背景技术
随着微机电系统、微加工和微操作以及新材料等领域的迅速发展,柔顺机构的设计已经成为目前国内外机构学领域的研究热点。
采用拓扑优化方法研究柔顺机构的设计只需给定设计域和指定输入输出位置,无须从一个已知的刚性机构出发,且所得的机构具有分布式柔性的优越性能,因而引起了极大的重视。
这种方法通常是以有限元分析为基础,在拓扑优化的初始阶段,首先将设计域离散成一定数量的有限元网格,再利用优化方法确定单元材料的保留与删除,以满足预定的目标和约束条件。
在优化结果中,单元密度的理想取值为0或1,当单元密度取值为0时,表示该单元被删除,单元密度取值为1时,该单元被保留,于是,由高密度单元连接构成机构的拓扑图。但是,这类整数规划问题往往很难求解。因此,人们常采用松弛法,如均匀化方法和相对密度法,用连续变量的优化模型逼近原离散变量的优化模型,于是,拓扑设计变量可以在区间[0,1]内连续取值,使优化结果中出现中间密度单元。这一现象与工程实际不符。因此,在柔顺机构拓扑优化设计中,必须发展一种有效的方法,以克服中间单元问题,提取出清晰的柔顺机构0-1拓扑图。
目前,在包括柔顺机构在内的拓扑优化领域中,主要包括以下几种方法:
一是采用形如ρP,(P>1)的惩罚形式,以减小中间密度单元的刚度影响,然而,这种方法难以彻底去除中间单元。另一种方法是阈值法,但由于其阈值选择是启发式的,因而必然影响到方法的执行效果。此外,又发展了在优化模型的目标函数或约束条件中追加显式的惩罚函数方法,但这种额外的惩罚项也可能给优化软件寻找可行解造成困难。类似的惩罚函数法还有混合的SINH方法和凹约束方法,不过,这两种方法也依然不能提取出满意的0-1拓扑图。另外,基于图像的过滤技术也被尝试来克服中间单元问题,然而,由于其灵敏度计算复杂,使得计算效率大大降低,且所得的拓扑图具有网格依赖性。
近来,又提出了一种修改最佳准则表达形式的启发式方法,但这种方法的最佳拓扑图具有参数敏感性。除了上述方法之外,渐进方法也是一种常用的去除中间单元的方法。它是通过逐渐地去除无效材料来获得最佳解,其优点在于不需引入太多的数学理论,不过由于该方法属于整数规划范畴,灵敏度连续性差,且被去除的单元不能再被利用,因而计算效率低,并易于导致结果不收敛。在此基础上,又发展了双向渐进方法,这种方法旨在去除无效材料的同时,也能添加有效材料,以改进收敛效果,但是这一方法本质上仍是整数规划,且无法精确评估所添加的材料单元的灵敏度。另外,还有模拟退火和SIMP相结合的方法以及单元连结参数化方法等,从它们的优化结果来看,也依然无法提取出理想的0-1拓扑图。并且,上述方法有着易于导致棋盘格的共同缺点。
此外,基于水平集的拓扑优化方法也是解决中间单元问题的一种有效方法,但是这种方法也具有初始敏感性、不能生成新孔、计算效率低和难以收敛到不光滑的角点等缺陷,虽然目前已提出一些改进的方法,但这些问题尚未完全得到很好的解决。
发明内容
本发明的目的在于克服现有技术的缺点与不足,提供一种柔顺机构0-1变构型拓扑图提取方法,能够提取出满足目标体积的理想的柔顺机构0-1变构型拓扑图。
本发明的目的通过以下的技术方案实现:一种柔顺机构0-1变构型拓扑图提取方法,包括以下步骤:首先建立柔顺机构第一层拓扑优化模型;接着,对拓扑优化模型进行优化,提取出含有中间单元的拓扑图;然后,对含有中间单元的柔顺机构拓扑优化模型做进一步优化,提取出柔顺机构0-1拓扑图;之后,建立下一层柔顺机构拓扑优化模型,并对其进行优化,通过不断增大拓扑优化模型体积约束值的方式逐渐逼近最佳解,最终提取出满足目标体积的理想的柔顺机构0-1变构型拓扑图。
具体步骤如下:
S1:建立柔顺机构第一层拓扑优化模型;
S2:对拓扑优化模型进行优化,步骤如下:
S2-1:对拓扑优化模型进行有限元分析,并计算灵敏度;
S2-2:对灵敏度进行过滤,以消除拓扑图中的棋盘格;
S2-3:基于最佳准则更新设计变量;
S2-4:重复步骤S2-1至S2-3,直到优化迭代收敛,提取出含有中间单元的柔顺机构拓扑图;
S3:对含有中间单元的柔顺机构拓扑优化模型做进一步优化,步骤如下:
S3-1:对拓扑优化模型进行有限元分析,计算灵敏度;
S3-2:对灵敏度进行过滤;
S3-3:基于最佳准则更新设计变量,保留已有的实心单元,抑制密度低于预设值的单元;
S3-4:重复步骤S3-1至S3-3,直到优化迭代收敛,提取出柔顺机构0-1拓扑图;
S4:判断柔顺机构拓扑优化迭代次数是否大于预设阈值,以及前后两次迭代的单元密度变化最大值是否小于一个阈值,若上述任一条件满足则终止循环并输出结果,否则继续执行以下步骤;
S5:以柔顺机构0-1拓扑图为基准,建立下一层柔顺机构拓扑优化模型;
S6:重复步骤S2至S5,如此不断增大拓扑优化模型的体积约束值,直至提取出满足目标体积V*的柔顺机构0-1变构型拓扑图。
优选的,步骤S1中建立柔顺机构第一层拓扑优化模型的具体步骤如下:
以Ω1为设计域,其为柔顺机构初始拓扑优化模型可利用材料域;Pi和Po分别为柔顺机构载荷输入点和位移输出点;Fin和Fd分别为柔顺机构输入载荷和沿输出位移方向的虚拟单位载荷;kin和kout分别为柔顺机构输入和输出弹簧刚度;该柔顺机构的目标体积比为θ*,初始阶段的拓扑优化层数为L=1,其体积比为将设计域离散为Nx行Ny列;
柔顺机构的应变能和互应变能如下:
Es=∫Ωε(u)TDε(u)dΩ=UTKU
式中,Es是系统的应变能,Ems是系统的互应变能,D为弹性矩阵,K是系统整体刚度矩阵,U是Fin作用下的节点位移向量,Ud是Fd作用下的节点位移向量,ε(u)和u是设计域内任一点在载荷F作用下的应变和弹性变形,ε(ud)和ud是设计域内任一点在载荷Fd作用下的应变和弹性变形;
通过多目标优化得到柔顺机构的应变能和互应变能的关系如下:
Min:式中,Min代表最小值;
首先,基于相对密度法松弛设计变量,使柔顺机构单元密度在0-1范围内取值,即:
0<ρmin≤ρi≤ρmax=1,i=1,2,…,N1;
式中,ρi是单元i的密度,ρmin是单元密度下限,ρmax是单元密度上限,N1是自然数,表明初始阶段可利用材料域的单元数,且N1=NxNy;
柔顺机构拓扑优化模型的整体刚度矩阵是:
式中,Ki是单元i的刚度矩阵,Ve是任一实心单元的材料体积,V是体积,P为密度ρi的指数,且P∈Z,P>1,Z表示整数,K0为任一实心单元的单元刚度矩阵,且B是任一实心单元的应变矩阵;该柔顺机构第一层拓扑优化模型的体积约束如下:
式中,V0是柔顺机构初始拓扑优化模型可利用材料的体积,ρ是由ρi所构成的列向量,i=1,2,…,N;
根据上述目标和约束条件,得到柔顺机构第一层拓扑优化模型如下:
优选的,步骤S2-1中,优化目标的灵敏度为:
体积约束的灵敏度为:i=1,...,N1。
优选的,步骤S2-2中,过滤后的灵敏度如下:
式中,Ne为单元i的邻域,该邻域内的各单元中心到单元i的中心的距离小于或等于过滤半径r,ρj是单元i的邻域内单元j的密度,是卷积因子,如下:
式中,xi是单元i的中心坐标,xj是单元j的中心坐标。
优选的,步骤S2-3中,最佳准则如下:
式中,L≥1,ΩL,2是第L-1层提取的柔顺机构0-1拓扑图中密度为ρmax且恒定不变的单元集合,当L=1时,ΩL,2为空集,Ω1\ΩL,2是步骤S2的可利用材料域中密度可变的单元集合,t为迭代数,η(0<η<1)为松弛因子,ζ为一个移动极限,为一个非负数,即:
式中,由Kuhn–Tucker必要条件导出,如下:
式中,λ为拉格朗日乘子,在设计变量更新过程中采用二分法得到,以使更新后的密度满足体积约束,即:
式中,是步骤S2的体积比,和分别是集合Ω1\ΩL,2中所更新的三类设计变量的和。
进一步的,基于最佳准则更新设计变量的具体过程是:
首先,给定λ的最小值和最大值,分别是λ1和λ2;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式,得到更新后的单元密度;
之后,计算体积约束余量再根据值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。
优选的,步骤S3-3中,最佳准则如下:
式中,β为给定的阈值,用于抑制密度低于预设值的单元,且ρmin<β≤1;在迭代过程中,对β采用延拓方法进行修正;首先判断以下体积约束条件是否成立:
如果不成立则β不变;如果体积约束条件恒成立,再判断以下迭代数条件是否成立:
t-t1-[(t-t1)/t2]t2≡1
式中,符号[]表示取整运算,t1为步骤S2-4中提取出含有中间单元的柔顺机构拓扑图的总迭代数,t2为允许β恒定不变的迭代数;
如果不成立则β不变;如果迭代数条件恒成立,β取为
β=β+min(Δβ,1-β)
式中,Δβ为给定的β的增量,min代表最小值;
更新后的密度应满足如下体积约束:
式中,是步骤S2中的体积比,和分别是集合Ω1\ΩL,2中所更新的各类设计变量的和。
进一步的,步骤S3-3中,
首先,给定λ的最小值和最大值,分别是λ1和λ2;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式,得到更新后的单元密度;
之后,计算体积约束余量再根据值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。
优选的,步骤S5中,建立下一层柔顺机构拓扑优化模型的方法是:首先修改柔顺机构拓扑优化的层数L,使L=L+1;接着,修改体积比:采用体积约束延拓方法,对下一层优化模型,体积比取为:
式中,max代表最大值,α是体积约束延拓因子初值,且0<α<1,Δα是体积约束延拓因子增量,且0<Δα<1;
最后,得到下一层(即第L,(L≥2)层)柔顺机构拓扑优化模型为:
本发明与现有技术相比,具有如下优点和有益效果:
1、本发明采用体积约束延拓方法,通过不断增大拓扑优化模型体积约束值的方式逐渐逼近最佳解,能最大限度地保留所有可能单元,提取出柔顺机构0-1变构型拓扑图。
2、本发明通过给定阈值抑制较低密度单元,同时,采用延拓方法对阈值进行修正,可有效克服基于SIMP方法松弛设计变量时所产生的中间单元问题,提取出理想的柔顺机构0-1拓扑图。
3、本发明具有连续设计变量,灵敏度连续性好,而且概念简单,方法易于执行,不易导致棋盘格现象。
4、本发明克服了传统方法的参数敏感性和网格依赖性的缺陷,具有较好的收敛稳定性。
附图说明
图1为本发明柔顺机构0-1变构型拓扑图提取方法流程图;
图2为实施例柔顺力-位移反向机构设计域和边界条件示意图;
图3为第一层拓扑优化结果示意图;
图4为第二层拓扑优化结果示意图。
具体实施方式
下面结合实施例及附图对本发明作进一步详细的描述,但本发明的实施方式不限于此。
如图1所示,为本发明所提出的一种柔顺机构0-1拓扑图提取方法流程图。一种柔顺机构0-1变构型拓扑图提取方法,包括以下步骤:首先建立柔顺机构第一层拓扑优化模型;接着,对拓扑优化模型进行优化,提取出含有中间单元的拓扑图;然后,对含有中间单元的柔顺机构拓扑优化模型做进一步优化,提取出柔顺机构0-1拓扑图;之后,建立下一层柔顺机构拓扑优化模型,并对其进行优化,通过不断增大拓扑优化模型体积约束值的方式逐渐逼近最佳解,最终提取出满足目标体积的理想的柔顺机构0-1变构型拓扑图。
具体步骤如下:
S1:建立柔顺机构第一层拓扑优化模型;
本发明实施例是一种典型的柔顺力-位移反向机构,其设计域和边界条件如图2所示。其中,设计域Ω1大小为120×120,材料的弹性模量和泊松比分别为E=1和v=0.3,Pi和Po分别为载荷输入点和位移输出点,输入载荷Fin=1,同时,沿输出位移方向作用有一虚拟单位载荷Fd,输入和输出弹簧刚度分别为kin=1和kout=0.001。
柔顺机构的目标体积比为θ*=0.25,灵敏度过滤半径为r=2.5,初始阶段的拓扑优化层数为L=1,其体积比为将设计域离散为Nx行Ny列;柔顺机构的应变能和互应变能如下:
Es=∫Ωε(u)TDε(u)dΩ=UTKU
式中,Es是系统的应变能,应变能越小则表明系统的刚度越大;Ems是系统的互应变能,互应变能越大则表明系统的柔性越大;D为弹性矩阵,K是系统整体刚度矩阵,U是Fin作用下的节点位移向量,Ud是Fd作用下的节点位移向量,ε(u)和u是设计域内任一点在载荷F作用下的应变和弹性变形,ε(ud)和ud是设计域内任一点在载荷Fd作用下的应变和弹性变形;
由于设计域和边界条件具有对称性,因此仅取设计域的下半区域进行优化。采用SIMP方法松弛密度变量,以柔顺机构的最大柔度和最大刚度为优化目标,通过多目标优化而得到柔顺机构的应变能和互应变能的关系如下:
Min:式中,Min代表最小值;
首先,基于相对密度法松弛设计变量,使柔顺机构单元密度在0-1范围内取值,即:
0<ρmin≤ρi≤ρmax=1,i=1,2,…,N1;
式中,ρi是单元i的密度,ρmin是单元密度下限,ρmax是单元密度上限,N1是自然数,表明初始阶段可利用材料域的单元数,且N1=NxNy;
柔顺机构拓扑优化模型的整体刚度矩阵是:
式中,Ki是单元i的刚度矩阵,Ve是任一实心单元的材料体积,V是体积,P为密度ρi的指数,且P∈Z,P>1,Z表示整数,本实施例中,为了加快收敛,取P=3;K0为任一实心单元的单元刚度矩阵,且B是任一实心单元的应变矩阵;该柔顺机构第一层拓扑优化模型的体积约束如下:
式中,V0是柔顺机构初始拓扑优化模型可利用材料的体积,ρ是由ρi,i=1,2,…,N所构成的列向量;
根据上述目标和约束条件,得到柔顺机构第一层拓扑优化模型如下:
S2:对拓扑优化模型进行优化,步骤如下:
S2-1:对拓扑优化模型进行有限元分析,并计算灵敏度;
为了获得优化模型的拓扑优化最佳结果,按以下公式对其进行有限元分析,求出节点位移向量U和Ud:
KU=Fin;
KUd=Fd。
进一步得出优化目标和体积约束的灵敏度,优化目标的灵敏度为:
体积约束的灵敏度为:i=1,...,N1。
S2-2:对灵敏度进行过滤,以消除拓扑图中的棋盘格;
过滤后的灵敏度如下:
式中,Ne为单元i的邻域,该邻域内的各单元中心到单元i的中心的距离小于或等于灵敏度过滤半径r,ρj是单元i的邻域内单元j的密度,是卷积因子,即:式中,xi是单元i的中心坐标,xj是单元j的中心坐标。
S2-3:基于最佳准则更新设计变量,最佳准则公式如下:
式中,L≥1,ΩL,2是第L-1层提取的柔顺机构0-1拓扑图中密度为ρmax且恒定不变的单元集合,当L=1时,ΩL,2为空集,Ω1\ΩL,2是第L层步骤S2的可利用材料域中密度可变的单元集合,t为迭代数,η(0<η<1)为松弛因子,ζ为一个较小的移动极限,为一个非负数,即:
式中,由Kuhn–Tucker必要条件导出,如下:
式中,λ为拉格朗日乘子,在设计变量更新过程中采用二分法得到,以使更新后的密度满足体积约束,即:
式中,是步骤S2的体积比,和分别是集合Ω1\ΩL,2中所更新的三类设计变量的和。
基于最佳准则更新设计变量的具体过程是:
首先,给定λ的最小值和最大值,分别是λ1=0和λ2=100000;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式(最佳准则公式),得到更新后的单元密度;
之后,计算体积约束余量再根据值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。此条件满足后则不再变化的值,即收敛。
下一层优化模型,重新确定λ。
S2-4:重复步骤S2-1至S2-3,直到优化迭代收敛,提取出含有中间单元的柔顺机构拓扑图;
S3:对含有中间单元的柔顺机构拓扑优化模型做进一步优化,步骤如下:
S3-1:对拓扑优化模型进行有限元分析,计算灵敏度;
S3-2:对灵敏度进行过滤;
S3-3:基于最佳准则更新设计变量,保留已有的实心单元,抑制较低密度单元(密度低于预设值);
最佳准则如下:
式中,β为给定的阈值,其作用是抑制较低密度单元,且ρmin<β≤1;
在迭代过程中,对β采用延拓方法进行修正;首先判断以下体积约束条件是否成立:
如果体积约束条件恒成立,再判断以下迭代数条件是否成立:
t-t1-[(t-t1)/t2]t2≡1
式中,符号[]表示取整运算t1为提取出含有中间单元的柔顺机构拓扑图的总迭代数,t2为允许β恒定不变的迭代数;如果迭代数条件恒成立,β取为
β=β+min(Δβ,1-β)
式中,Δβ为给定的β的增量,min代表最小值。更新后的密度应满足如下体积约束:
式中,是步骤S2的体积比,和分别是集合Ω1\ΩL,2中所更新的各类设计变量的和。
S3-4:重复步骤S3-1至S3-3,直到优化迭代收敛,提取出柔顺机构0-1拓扑图;
首先,给定λ的最小值和最大值,分别是λ1和λ2;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式,得到更新后的单元密度;
之后,计算体积约束余量再根据值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。
S4:判断柔顺机构拓扑优化迭代次数是否大于预设阈值,以及前后两次迭代的单元密度变化最大值是否小于一个阈值,若上述任一条件满足则终止循环并输出结果,否则继续执行以下步骤;
S5:以柔顺机构0-1拓扑图为基准,建立下一层柔顺机构拓扑优化模型;
首先修改柔顺机构拓扑优化的层数L,使L=L+1;接着,修改体积比:采用体积约束延拓方法,对下一层优化模型,体积比取为:
式中,max代表最大值,α是体积约束延拓因子初值,且0<α<1,Δα是体积约束延拓因子增量,且0<Δα<1;
最后,得到下一层(即第L,(L≥2)层)柔顺机构拓扑优化模型为:
式中,ΩL,2是第L-1层提取的柔顺机构0-1拓扑图中密度为ρmax且恒定不变的单元集合,Ω1\ΩL,2是第L,(L≥2)层可利用材料域中密度可变的单元集合。
S6:重复步骤S2至S5,如此不断增大拓扑优化模型的体积约束值,直至提取出满足目标体积V*的柔顺机构0-1变构型拓扑图。图3为第一层拓扑优化结果示意图,图4为第二层拓扑优化结果示意图。
上述实施例为本发明较佳的实施方式,但本发明的实施方式并不受上述实施例的限制,其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化,均应为等效的置换方式,都包含在本发明的保护范围之内。
Claims (8)
1.一种柔顺机构0-1变构型拓扑图提取方法,其特征在于,具体步骤如下:
S1:建立柔顺机构第一层拓扑优化模型;
S2:对拓扑优化模型进行优化,步骤如下:
S2-1:对拓扑优化模型进行有限元分析,并计算灵敏度;
S2-2:对灵敏度进行过滤,以消除拓扑图中的棋盘格;
S2-3:基于最佳准则更新设计变量;
S2-4:重复步骤S2-1至S2-3,直到优化迭代收敛,提取出含有中间单元的柔顺机构拓扑图;
S3:对含有中间单元的柔顺机构拓扑优化模型做进一步优化,步骤如下:
S3-1:对拓扑优化模型进行有限元分析,计算灵敏度;
S3-2:对灵敏度进行过滤;
S3-3:基于最佳准则更新设计变量,保留已有的实心单元,抑制密度低于预设值的单元;
步骤S3-3中,最佳准则如下:
式中,β为给定的阈值,用于抑制密度低于预设值的单元,且ρmin<β≤1;
ρi是单元i的密度,ρmin是单元密度下限,ρmax是单元密度上限,L≥1,ΩL,2是第L-1层提取的柔顺机构0-1拓扑图中密度为ρmax且恒定不变的单元集合,当L=1时,ΩL,2为空集,Ω1\ΩL,2是步骤S2的可利用材料域中密度可变的单元集合,t为迭代数,η为松弛因子,0<η<1,ζ为一个移动极限,为一个非负数,即:
在迭代过程中,对β采用延拓方法进行修正;首先判断以下体积约束条件是否成立:
Ve是任一实心单元的材料体积,V0是柔顺机构初始拓扑优化模型可利用材料的体积,是步骤S2中的体积比,如果不成立则β不变;如果体积约束条件恒成立,再判断以下迭代数条件是否成立:
t-t1-[(t-t1)/t2]t2≡1
式中,符号[]表示取整运算,t1为步骤S2-4中提取出含有中间单元的柔顺机构拓扑图的总迭代数,t2为允许β恒定不变的迭代数;
如果不成立则β不变;如果迭代数条件恒成立,β取为
β=β+min(Δβ,1-β)
式中,Δβ为给定的β的增量,min代表最小值;
更新后的密度应满足如下体积约束:
式中,是步骤S2中的体积比,和分别是集合Ω1\ΩL,2中所更新的各类设计变量的和;
S3-4:重复步骤S3-1至S3-3,直到优化迭代收敛,提取出柔顺机构0-1拓扑图;
S4:判断柔顺机构拓扑优化迭代次数是否大于预设阈值,以及前后两次迭代的单元密度变化最大值是否小于一个阈值,若上述任一条件满足则终止循环并输出结果,否则继续执行以下步骤;
S5:以柔顺机构0-1拓扑图为基准,建立下一层柔顺机构拓扑优化模型;
S6:重复步骤S2至S5,如此不断增大拓扑优化模型的体积约束值,直至提取出满足目标体积V*的柔顺机构0-1变构型拓扑图。
2.根据权利要求1所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S1中建立柔顺机构第一层拓扑优化模型的具体步骤如下:
以Ω1为设计域,其为柔顺机构初始拓扑优化模型可利用材料域;Pi和Po分别为柔顺机构载荷输入点和位移输出点;Fin和Fd分别为柔顺机构输入载荷和沿输出位移方向的虚拟单位载荷;kin和kout分别为柔顺机构输入和输出弹簧刚度;该柔顺机构的目标体积比为θ*,初始阶段的拓扑优化层数为L=1,其体积比为 将设计域离散为Nx行Ny列;
柔顺机构的应变能和互应变能如下:
Es=∫Ωε(u)TDε(u)dΩ=UTKU
式中,Es是系统的应变能,Ems是系统的互应变能,D为弹性矩阵,K是系统整体刚度矩阵,U是Fin作用下的节点位移向量,Ud是Fd作用下的节点位移向量,ε(u)和u是设计域内任一点在载荷F作用下的应变和弹性变形,ε(ud)和ud是设计域内任一点在载荷Fd作用下的应变和弹性变形;
通过多目标优化得到柔顺机构的应变能和互应变能的关系如下:
式中,Min代表最小值;
首先,基于相对密度法松弛设计变量,使柔顺机构单元密度在0-1范围内取值,即:
0<ρmin≤ρi≤ρmax=1,i=1,2,…,N1;
式中,ρi是单元i的密度,ρmin是单元密度下限,ρmax是单元密度上限,N1是自然数,表明初始阶段可利用材料域的单元数,且N1=NxNy;
柔顺机构拓扑优化模型的整体刚度矩阵是:
式中,Ki是单元i的刚度矩阵,Ve是任一实心单元的材料体积,V是体积,P为密度ρi的指数,且P∈Z,P>1,Z表示整数,K0为任一实心单元的单元刚度矩阵,且B是任一实心单元的应变矩阵;该柔顺机构第一层拓扑优化模型的体积约束如下:
式中,V0是柔顺机构初始拓扑优化模型可利用材料的体积,ρ是由ρi所构成的列向量,i=1,2,…,N;
根据下述目标和约束条件,得到柔顺机构第一层拓扑优化模型如下:
3.根据权利要求2所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S2-1中,优化目标的灵敏度为:
体积约束的灵敏度为:
4.根据权利要求2所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S2-2中,过滤后的灵敏度如下:
式中,Ne为单元i的邻域,该邻域内的各单元中心到单元i的中心的距离小于或等于灵敏度过滤半径r,ρj是单元i的邻域内单元j的密度,是卷积因子,如下:
式中,xi是单元i的中心坐标,xj是单元j的中心坐标。
5.根据权利要求2所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S2-3中,最佳准则如下:
式中,L≥1,ΩL,2是第L-1层提取的柔顺机构0-1拓扑图中密度为ρmax且恒定不变的单元集合,当L=1时,ΩL,2为空集,Ω1\ΩL,2是步骤S2的可利用材料域中密度可变的单元集合,t为迭代数,η(0<η<1)为松弛因子,ζ为一个移动极限,为一个非负数,即:
式中,由Kuhn–Tucker必要条件导出,如下:
式中,λ为拉格朗日乘子,在设计变量更新过程中采用二分法得到,以使更新后的密度满足体积约束,即:
式中,是步骤S2中的体积比,和分别是集合Ω1\ΩL,2中所更新的三类设计变量的和。
6.根据权利要求4所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,基于最佳准则更新设计变量的具体过程是:
首先,给定λ的最小值和最大值,分别是λ1和λ2;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式,得到更新后的单元密度;
之后,计算体积约束余量再根据▽V值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。
7.根据权利要求2所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S3-3中,
首先,给定λ的最小值和最大值,分别是λ1和λ2;
然后,计算它们的中点值:λ=0.5(λ1+λ2),并代入密度变量更新公式,得到更新后的单元密度;
之后,计算体积约束余量再根据▽V值修改λ1或λ2的值:如果则使λ1=λ,否则如果则使λ2=λ;
重复λ1和λ2的中点值的计算、单元密度更新以及λ1或λ2的修改过程,直至λ2与λ1的差小于一个阈值为止。
8.根据权利要求2所述的柔顺机构0-1变构型拓扑图提取方法,其特征在于,步骤S5中,建立下一层柔顺机构拓扑优化模型的方法是:首先修改柔顺机构拓扑优化的层数L,使L=L+1;接着,修改体积比:采用体积约束延拓方法,对下一层优化模型,体积比取为:
式中,max代表最大值,α是体积约束延拓因子初值,且0<α<1,Δα是体积约束延拓因子增量,且0<Δα<1;
最后,得到下一层即第L(L≥2)层柔顺机构拓扑优化模型为:
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201611002663.2A CN106529044B (zh) | 2016-11-14 | 2016-11-14 | 一种柔顺机构0-1变构型拓扑图提取方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201611002663.2A CN106529044B (zh) | 2016-11-14 | 2016-11-14 | 一种柔顺机构0-1变构型拓扑图提取方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN106529044A CN106529044A (zh) | 2017-03-22 |
CN106529044B true CN106529044B (zh) | 2019-12-10 |
Family
ID=58351890
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201611002663.2A Active CN106529044B (zh) | 2016-11-14 | 2016-11-14 | 一种柔顺机构0-1变构型拓扑图提取方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN106529044B (zh) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109472052A (zh) * | 2018-10-12 | 2019-03-15 | 中国航天空气动力技术研究院 | 一种高速机械手支座结构拓扑优化设计方法 |
CN109583091B (zh) * | 2018-11-30 | 2022-11-29 | 长沙理工大学 | 基于自适应约束的柔性机构拓扑优化设计方法 |
CN110414165B (zh) * | 2019-08-01 | 2022-06-10 | 华东交通大学 | 一种基于全局应力约束的多相材料柔顺机构拓扑优化方法 |
CN111027110B (zh) * | 2019-11-27 | 2023-06-30 | 中国科学院光电技术研究所 | 一种连续体结构拓扑与形状尺寸综合优化方法 |
CN111027253B (zh) * | 2019-12-17 | 2021-12-03 | 华中科技大学 | 一种含有最大尺寸约束的渐进式结构拓扑优化方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103020361A (zh) * | 2012-12-17 | 2013-04-03 | 华南理工大学 | 一种柔顺机构无棋盘格拓扑图提取方法 |
CN103049606A (zh) * | 2012-12-17 | 2013-04-17 | 华南理工大学 | 一种柔顺机构0-1拓扑图提取方法 |
-
2016
- 2016-11-14 CN CN201611002663.2A patent/CN106529044B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103020361A (zh) * | 2012-12-17 | 2013-04-03 | 华南理工大学 | 一种柔顺机构无棋盘格拓扑图提取方法 |
CN103049606A (zh) * | 2012-12-17 | 2013-04-17 | 华南理工大学 | 一种柔顺机构0-1拓扑图提取方法 |
Non-Patent Citations (2)
Title |
---|
多层渐进黑白拓扑优化设计方法;付永清 等;《中国机械工程》;20170606;第1292-1299页 * |
结构及柔顺机构拓扑优化设计中的拓扑图提取;付永清 等;《力学进展》;20060225;第75-84页 * |
Also Published As
Publication number | Publication date |
---|---|
CN106529044A (zh) | 2017-03-22 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106529044B (zh) | 一种柔顺机构0-1变构型拓扑图提取方法 | |
AU2019442319B2 (en) | Structural topology optimization method based on material-field reduction series expansion | |
CN110069800B (zh) | 具有光滑边界表达的三维结构拓扑优化设计方法及设备 | |
CN112836411B (zh) | 加筋板壳结构的优化方法、装置、计算机设备和存储介质 | |
Paulino et al. | A modified Q4/Q4 element for topology optimization | |
CN106372347A (zh) | 改进双向渐进法的等效静载荷法动态响应拓扑优化方法 | |
CN108763827B (zh) | 一种输电塔有限元模型建立方法及装置 | |
Williamschen et al. | Parallel anisotropic block-based adaptive mesh refinement algorithm for three-dimensional flows | |
CN113887095B (zh) | 一种基于等几何分析的渐进式结构拓扑优化方法 | |
CN103020361B (zh) | 一种柔顺机构无棋盘格拓扑图提取方法 | |
CN111523270A (zh) | 一种改进的连续体结构拓扑优化后处理方法 | |
US20170011143A1 (en) | Search and poll method for solving multi-fidelity optimization problems | |
CN110069864A (zh) | 一种结合变密度法的改进双向渐进结构拓扑优化方法 | |
Zhang et al. | Parallel high-order anisotropic block-based adaptive mesh refinement finite-volume scheme | |
CN111079327B (zh) | 一种基于各向异性过滤技术的异形曲面加筋拓扑优化方法 | |
CN106599049B (zh) | 一种决策表数据约简方法 | |
Freret et al. | Enhanced anisotropic block-based adaptive mesh refinement for three-dimensional inviscid and viscous compressible flows | |
CN103049606A (zh) | 一种柔顺机构0-1拓扑图提取方法 | |
Andrews et al. | Deterministic control of randomly-terminated processes | |
Xiao et al. | A hybrid method for density-related topology optimization | |
Mavriplis | Viscous flow analysis using a parallel unstructured multigrid solver | |
CN106874611B (zh) | 一种基于超体积迭代策略的含区间参数结构响应区间的分析方法 | |
CN114444012A (zh) | 一种基于神经网络架构搜索的偏微分方程数值求解方法 | |
Zhu et al. | Explicit Isogeometric Topology Optimization Method with Suitably Graded Truncated Hierarchical B-Spline. | |
JP2008197921A (ja) | シミュレーション装置、シミュレーション方法及びシミュレーションプログラム |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | 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 |