王 超 徐 斌,2) 段尊义 荣见华
∗(西北工业大学力学与土木建筑学院,西安 710072)
†(长沙理工大学汽车与机械工程学院,长沙 410114)
拓扑优化是一种重要的结构优化方法,已得到学术界和工业界的广泛关注[1-2].多年来,大量研究致力于如何利用拓扑优化改进设计性能[3-5].随着该方法的进一步推广应用,优化设计的可制造性问题逐渐引起人们的重视[6-7].然而,现有研究关于性能和可制造性的考虑通常是分开处理的,常侧重于某一方面指标.因此,将这两个重要而又不可或缺的设计目标结合起来进行优化设计,将是一项极具实际应用价值且富有挑战性的工作.
早期基于拓扑优化的设计性能研究主要集中在刚度问题上,但在实际工程应用中,多数设计还应考虑强度要求[8-10].因为局部高应力集中易导致构件疲劳损伤和断裂破坏,将严重影响结构的安全使用和耐久性.在这种情况下,考虑应力相关的拓扑优化将发挥重要作用[11-12].此外,研究还表明[13-14]强度最大化设计并不一定等同于刚度最大化设计,这就意味着在刚度最大化准则下获得的优化设计将不一定适用于强度导向下的工程问题.因此,面向强度需求,研究基于应力的拓扑优化问题既符合工程需要也极具研究意义[15-17].
自Duysinx 等[18]提出连续体结构应力拓扑优化研究以来,应力相关的拓扑优化问题受到了越来越多的关注.然而,考虑应力相关的拓扑优化问题有其自身难点,常面临应力奇异性、局部性、非线性及度量准确性等问题[19-20].为此,不少学者提出了不同的解决方法.目前常用的处理方法包括松弛方法[21]、基于P 范数或K-S 函数的凝聚技术[8,14]、分区度量[20]、以及凝聚近似误差修正策略[22]等.此外,隋允康等[23]提出的将局部应力约束向整体应变能约束转化的方法也可有效解决应力局部性问题.尽管上述方法能有效处理常规的拓扑优化问题,但尚未涉及复杂的可制造性拓扑优化问题.
事实上,拓扑优化设计常面临可制造性问题[24].在传统制造技术条件下,优化产生的复杂构型将给产品制备带来巨大挑战.这在一定程度上限制了拓扑优化的推广应用,不利于高性能产品的研发.作为一种先进的制造技术,增材制造能够实现高度复杂结构的有效制备[25].将增材制造与拓扑优化有机结合,能够大幅提高创新产品的研制能力,现已成为重要研究方向[6-7,25].然而,增材制造并非“自由”制造,有其独特的制造约束.要求结构不含封闭孔洞的连通性约束便是其中最重要的制造约束之一[26-27].如果结构含有封闭孔洞,制造后残留在孔洞内的未熔融粉末或支撑结构将难以去除.显然,面向增材制造,考虑可制造连通性的拓扑优化问题将具有重要工程应用价值[27-29].
研究人员已提出一些有效方法.Liu 等[27]利用热传导原理提出了虚拟温度法,在结构孔洞区引入虚拟热源,而在固体区填充隔热材料,然后通过限制结构上的最高温度来抑制封闭孔洞.数学上,该方法可描述为服从泊松方程的标量场约束[29].Zhou 等[28]则提出一种边界约束方法,利用孔洞作为基本设计元素,通过将与所有孔洞中心相关的设计变量限定在设计域之外来确保结构连通性.Xiong 等[30]提出了一种基于BESO 的结构连通性控制方法,通过有选择地生成将孔洞与结构边界连接起来的通道来消除封闭孔洞.在水平集框架下,Xia 等[31]将孔洞抑制转换为对设计速度矢量的约束;Wang 等[32]则提出一种基于水平集函数梯度的约束来描述孔洞抑制条件.然而,这些研究只考虑了传统刚度问题,缺乏对应力问题的必要关注.
综上所述,在拓扑优化领域,应力问题和连通性问题分别已得到了有效处理.然而,鲜有研究同时考虑这两个重要的设计要求.一方面,如前所述,这两个问题均有其特殊困难.另一方面,两种要求的结合涉及相互影响,将加剧优化设计的复杂性.一般而言,对于相对复杂的优化问题,设计在优化过程中往往会出现不可预测的数值响应,这将给优化求解带来各种挑战.尽管研究人员业已提出有效的方法来克服应力或连通性拓扑优化中的诸多困难,但到目前为止,现有应力优化技术与可制造性连通性方法之间的兼容性尚未得到必要验证.因此,同时考虑应力和连通性问题,发展面向增材制造的高强度结构优化方法,不仅具有实际工程应用价值还具有重要的科学研究意义.
基于上述原因,本文将研究协同考虑应力和可制造连通性的拓扑优化问题.建立以结构全局最大应力最小化为目标,以材料体积和结构连通性为约束的拓扑优化模型.通过合理引入不同优化技术来构建有效的优化求解策略,以克服优化过程中的不同数值困难.最后基于典型数值算例,探讨文中模型及方法的可行性及其数值特点,以期为高性能可制造性结构的研发提供有益参考.
针对协同考虑结构强度和可制造连通性的拓扑优化问题,基于经典密度法[33],建立以结构全局最大应力强度准则为目标,以材料用量和结构连通性为约束的拓扑优化模型,其数学列式如下
为有效实现拓扑优化中抑制封闭孔洞的单连通性约束,采用Liu 等[27]提出的基于泊松方程的标量场约束方法[29].为便于该方法的描述,不失一般性,这里采用了静电场物理模型[35]进行表征.
在静电场模型中,假定孔洞单元材料为具有强介电性能的带静电载体,实体单元材料为介电性能较差的不带电绝缘介质,相应材料模型如下
其中,ε 和ρ 分别表示材料介电常数(F/m)和电荷密度(C/m3);εϕ和ρϕ为实常数;µ 为一个较小的正数(避免矩阵奇异),与静电模型的势强度表征有关.此外,该模型中需假定用户指定的域边界接地(即约束边界),以作为零电势参考基准.
当结构含有封闭孔洞时,根据经典电磁学原理可知,带电体(封闭孔洞) 将在其周围介质中激发出强静电场(可由静电标量势ϕ 表征).一旦带电体接地(封闭孔洞与约束边界连通),电场强度将大幅衰减.于是,可通过分析结构静电场的最大势强度ϕmax来评估其连通性.相应地,拓扑优化中抑制封闭孔洞的单连通性约束可等效为对结构最大势强度的约束,即静电势约束ϕEPC.
其中,η·ϕ0为评估结构是否含封闭孔洞的约束阈值;ϕ0为基准阈值,可基于初始化约束模型(仅由孔洞单元构成)求解获得[29];η>0 为约束阈值修正因子,需针对不同问题进行调整,以控制约束效果[35].控制ϕ的二阶微分方程为经典泊松方程[36]
其中,∇为微分符号.静电问题的边界条件如下
其中,ΓD和ΓN分别为狄利克雷边界(即连通性约束边界) 和诺伊曼边界,n为边界外法向单位向量.需要注意的是,基于上述方法获得的优化结果对约束边界条件的选择有依赖性,不同约束边界条件对应不同约束问题,可产生不同优化结果[35].基于式(2)和式(5)可将原优化问题转化为静电场约束问题(式(3))求解分析.
静电场与原方法中的温度场[27]在数学上一致,均为服从泊松方程的标量场.这里采用静电场模型主要有两方面考虑.一方面,静电原理也能较好表征封闭空洞,而不同的物理场具有不同物理规律,可以启发读者进一步发展泊松方法.另一方面,静电理论有助于解析泊松方法中关键的约束关系和边界依赖性,有益于该方法的理解应用[35].
在拓扑优化模型中,可采用如下材料插值模型有效表征式(2)所示材料介电属性
由于问题(1) 为相对复杂的多物理场多约束优化问题,为实现高效稳健求解,需合理引入不同优化技术来克服相应数值困难.
为避免拓扑优化中的数值不稳定性问题[37],并获取清晰的0-1 设计,采取线性密度过滤[38-39]和非线性阈值投影[40]相结合的技术方案.对单元密度变量xi进行线性密度过滤可得过滤后的变量
式中,NEi表示单元中心落在以单元i为中心的半径为rmin的圆形区域内的单元集;加权因子ωik=max(0,rmin−∆(i,k)),∆(i,k)表示单元i和k中心之间的欧式距离.再对过滤后的变量进行非线性阈值投影可得投影后的物理密度
式中,β >0 为控制投影函数非线性的投影参数;x0∈[0,1]为调节投影函数拐点的投影阈值,本文取常用值x0=0.5.为提高应力评估准确性以及优化迭代稳定性,参考文献[19,41],这里采用连续化策略定义β 值:在前10 步迭代中取β=1,此后每步迭代加1 直至β=10.由于物理密度定义了单元属性,在后续研究中将使用物理密度进行结构分析.值得注意的是,过滤和投影会影响连通性约束,考虑这种影响对基于泊松方法的连通性约束的有效实现具有重要意义,感兴趣的读者可以参阅文献[35].
由于最大值函数不可微,不失一般性,本文采用P-norm 全局应力评估策略来近似结构最大应力
式中,ps为应力凝聚参数,其取值对凝聚函数的近似质量和光滑度有较大影响.理论上,当ps→∞时,近似误差将减小,但较大ps取值将加剧凝聚函数的非线性,引发迭代振荡及收敛困难问题,不利于优化问题的高效稳健求解.一般需要在足够光滑度和精确近似之间进行权衡,以确定较合理的ps取值.单元von Mises 应力则为
式中,σi为单元中心处的应力向量.为有效处理应力奇异性问题,这里采用文献[20]的方法,对应力施加不同于刚度的幂指数惩罚
式中,q为应力惩罚因子,这里取常用值q=0.5;D0为实体材料弹性矩阵,Bi为单元i中心处的应变位移矩阵.平面应力状态下,常数矩阵V为
基于最大值函数,表征结构连通性的静电势约束可表示为
式中,φj为节点势;ND为设计域节点总数.为使约束函数可微,同样可采用P-norm 函数来近似全局势约束.需要指出的是,静电约束中的势量级一般较大,以提高结构连通性的可识别性.而较大凝聚项易导致凝聚函数值失准及计算溢出问题.因此,有必要对凝聚函数正则化,得到改进的约束表达
式中,pc为势凝聚参数.如前所述,较大的凝聚参数值不利于高效稳健的优化求解,而适中取值将不可避免的引入近似误差.不同于前述目标函数中的最大近似目标值的趋势控制,约束函数中的P-norm 函数需要有效控制最大近似目标值,而近似误差的存在将极大影响该控制效果.此外,这里的势约束效果对约束阈值有较强依赖性[35],较大的近似误差也不利于设计连通性的有效控制.为此,引入基于稳定控制的近似误差修正技术[22]
为高效求解上述多约束优化问题,采用相适宜的移动渐近线算法[42].下面将利用链式法则和伴随法,推导算法所需的目标和约束函数的一阶灵敏度信息.基于链式法则,任意函数g关于设计变量xi的敏度为
下面将推导目标和约束函数关于物理密度的敏度,其中,体积约束的敏度可直接求导,故而省略.
目标函数式(10)关于单元物理密度的敏度可通过链式法则计算
式中,δei为克罗内克函数;矩阵Li从全局位移矢量中聚集单元i的节点位移矢量.假定载荷与设计无关,将平衡方程KU=F对物理密度求导得到
式中,伴随向量λs为下列伴随方程的解
本节通过经典数值算例来验证所提方法.不失一般性,仅考虑平面应力问题,优化参数均无量纲化.除特别声明外,弹性模量E0=1,泊松比v=0.3,体积分数f=0.4;为避免零密度导致的数值困难,取设计变量下限为0.001;所有域边界均设置为连通性约束边界.为提高收敛效率,算法中的默认移动限设为0.1.当连续迭代之间设计变量的最大变化小于0.01或迭代循环的最大次数达到200 时,优化迭代停止.设计拓扑和应力分布分别以灰度图(白色为空相,黑色为固相)和彩色图(蓝色为低应力,红色为高应力)表示.
考虑经典L 型梁优化问题[11,16],设计条件如图1(a) 所示.结构顶部固支,右端中点受集中载荷作用.为避免应力集中,载荷均匀分布到中点附近的5 个节点上.采用6400 个四节点双线性正方形有限元单元离散设计域.过滤半径rmin=4,连通性约束阈值η=1.00.
为便于比较,首先给出传统体积约束下的刚度最大化设计(图1(b)和图1(c)),以及考虑连通性约束的刚度最大化设计(图1(d) 和图1(e)).从图1(b) 和图1(c) 可以看出,传统刚度设计含有一些明显的封闭孔洞构造(x–y平面内),不利于具有连通性要求的增材制造.此外,设计拐角处存在局部高应力集中,易使结构发生断裂破坏,给结构安全带来较大隐患.
图1 L 型梁问题(a)及其刚度设计(b)(c)和连通性刚度设计(d)(e)Fig.1 L-shape beam problem(a),the stiffness design(b)(c)and the connectivity stiffness design(d)(e)
从图1(d) 和图1(e) 可以看出,引入连通性约束能够有效抑制优化设计中的封闭孔洞,但所得设计出现了严重恶化的局部高应力集中,其最大应力是传统刚度设计的1.87 倍.这表明可制造连通性优化设计中十分有必要考虑应力问题.此外,上述连通性设计的刚度比单纯的刚度设计要弱,这是连通性设计所付出的必要代价[35].
图2 展示了基于本文方法获得的不同应力凝聚参数下的应力最小化连通性设计.可以看出,所得设计均满足连通性要求.这些设计拓扑构型相似,但拐角处结构形状存在明显差异.当凝聚参数较小时(如ps=4),所得设计与连通性刚度设计(图1(d)和图1(e)) 相似.随着凝聚参数增大,拓扑形状发生显著变化,拐角构造逐渐变圆,应力集中效应得到有效缓解,其代价是牺牲设计刚度.当ps=8 和10 时,所得设计几乎接近于满应力设计,但ps=8 时的综合性能(结构刚度和强度) 较优.显然,利用本文方法能够获取应力集中有效缓解的连通性设计,且所得设计不同于传统的刚度最大化连通性设计.此外,上述结果还表明所得设计的性能与凝聚参数取值有关,参数值并非越大越好,合理取值才能有助于获取综合性能较优的设计.
图2 不同应力凝聚参数下应力最小化连通性设计Fig.2 Stress minimization connectivity designs for different stress aggregation parameters
图3 给出了图2 中性能最优设计ps=8 的目标和约束函数的迭代历史.可以看出,本文方法优化迭代历程稳健,最终设计能较好满足各约束条件.
图3 设计ps=8 的优化迭代历史Fig.3 Iteration histories of the design of ps=8
为研究体积分数对应力最小化连通性设计的影响,现考虑体积分数为0.3 和0.5 的优化问题,优化结果如图4 所示.为公平比较,其余优化参数与图2 中ps=8 的设计保持一致.与前述体积分数0.4 的设计相比,设计的刚度和强度在体积分数较小时均显著降低,反之,则明显提高.可见,增加材料用量能够有助于提高所研究设计的强度和刚度性能.
图4 不同体积分数下的应力最小化连通性设计Fig.4 Stress minimization connectivity designs for different volume fractions
为考察连通性约束边界条件对应力最小化连通性设计的影响,考虑如图5(a)∼图5(c) 和图5(d)∼图5(f)所示两种条件下的优化问题.为便于比较,其余优化参数与图2 中ps=8 的设计保持一致.对比图2 和图5 的设计不难发现,不同连通性约束边界条件可产生具有不同拓扑构型且性能各异的设计.该现象遵循了泊松方法的约束边界依赖性规律.特别地,图5(a)∼图5(c) 设计的综合性能不仅比图2 中ps=8 的设计要好(结构强度和刚度均分别提高了11.61%和24.21%),而且比图1(d)和图1(e)中的连通性刚度设计也要强(结构强度和刚度均分别提高了48.17%和5.25%).显然,连通性约束边界条件对所研究设计的拓扑构型及性能均有较大影响.然而,由于边界条件的变化将产生新的优化约束问题,这使得自动寻找最优约束边界条件仍旧是未来极具挑战性的工作[35].
图5 不同连通性约束边界条件下的应力最小化连通性设计,绿色粗线表示约束边界Fig.5 Stress minimization connectivity designs for different connectivity constraint boundary(highlighted by green thick lines)conditions
考虑典型T 型梁优化问题[16,22],设计条件如图6(a)所示,结构左右两侧固支,顶部中点受集中载荷作用.为避免应力集中,载荷均匀分配到中点附近5 个节点上.设计域离散为8800 个四节点双线性正方形有限元单元.过滤半径rmin=3,连通性约束阈值取η=2.01.
图6(b)∼图6(e)分别展示了有/无连通性约束的刚度最大化设计.可以看出,由于封闭孔洞和局部高应力集中问题的存在,这些设计在可制造性和/或结构强度方面存在明显局限性.
图6 T 型梁问题(a)及其刚度设计(b)(c)和连通性刚度设计(d)(e)Fig.6 T-shape beam problem(a),the stiffness design(b)(c)and the connectivity stiffness design(d)(e)
图7 给出了不同应力凝聚参数下的应力最小化连通性设计.类似于L 型梁设计情况,所得设计均能满足连通性约束要求,但设计综合性能随凝聚参数取值不同而有明显差异.应力凝聚参数较小时,设计结果趋近于传统刚度设计,结构刚度最优但强度最弱;随着应力凝聚参数增大,应力集中现象将得到有效缓解,但结构刚度也会削弱.给出的设计中ps=4的综合性能相对较优,当应力凝聚参数值进一步增大时,设计刚度和强度同步开始恶化.上述研究再次论证了本文方法的有效性,同时也揭示了所得设计对应力凝聚参数的依赖性.数值经验表明,对于本文所研究问题,应力凝聚参数取值为4 ∼8 时一般可获得较理想的优化效果.
图8 显示了不同连通性约束边界条件下的应力最小化连通性设计.为公平比较,其余优化参数与图7 中ps=6 的设计保持一致.可以看到,图8(a)∼图8(c) 的设计与图7 中ps=6 的设计相似,但结构强度和刚度均有所提升(分别为6.91%和10.86%),而图8(d)∼图8(f)设计的拓扑构型则明显不同,性能也有很大差异.连通性约束边界条件对所研究设计的重要影响再次得到论证.需要指出的是,数值试验还表明,连通性约束边界条件的变化并不能改变所得设计对应力凝聚参数的依赖.
图7 不同应力凝聚参数下应力最小化连通性设计Fig.7 Stress minimization connectivity designs for different stress aggregation parameters
图8 不同连通性约束边界条件下的应力最小化连通性设计Fig.8 Stress minimization connectivity designs for different connectivity constraint boundary conditions
兼顾性能和可制造性的拓扑优化具有重要工程意义和研究价值.本文面向增材制造,提出了协同考虑结构强度和可制造连通性的拓扑优化方法.构建了体积和连通性标量场约束下的应力最小化拓扑优化模型,开发了有效的优化求解策略,并详细推导了相关灵敏度列式.典型数值算例结果表明,传统的刚度最大化连通性设计不一定能避免局部高应力集中.基于文中模型和方法所获得的应力最小化连通性设计在确保可制造连通性的同时,能够有效减轻上述应力集中效应,且优化迭代历程稳健.优化结果还表明,所得设计与体积分数、连通性约束边界条件和应力凝聚参数有关,其合理配置是获取高性能设计的关键.本文研究论证了可制造连通性拓扑优化中考虑应力问题的必要性和可行性,可为高性能可制造性工程结构的研发提供有益参考.