力学大会集2015MS75航空航天结构与多学科优化设计

荣见华1,2，1,2，肖1,2，卿腾飞1,3， 汽车与机械，长沙, 工程车辆轻量化与可靠性技术湖南省高校，长沙 长沙理工大学工程车辆安全性设计与可靠性技术湖南省,长沙,:本文研究了基于梯度和凸分离的非线性近似函数的优化问题求解算法。基于MMA EvolutionaltopologyoptimizationdesignsbasedonMMAanddualsmoothapproximatefunctionsJianHuaRong1,2,LuoJi1,2，TingTingXiao1,2，TengFeiQing1,3，YuLiaohongSchoolofAutoandMechanicEngineering,ChangshaUniversityofScienceandTechnology,Changsha，410114,；r ;rjh@KeyLaboratoryofLightweightandReliabilityTechnologyforEngineeringVehicle,CollegeofHunanProvince,Changsha,410114,;KeyLaboratoryofSafetyDesignandReliabilityTechnologyforEngineeringVehicle,HunanProvince,：Anoptimizationalgorithmofconstraintnonlinearoptimizationproblemswithsimpleboundsisinvestigated,basedonconvexandseparableapproximatedapproximations.Andanovelevolutionaltopologyoptimizationmethodwithmultiplediscementconstraintsandsimpleoriginalboundsisdeveloped,incorporatedwiththeproposedsmoothdualalgorithmandavarieddiscementconstraintlimitscheme.First,tosolvethenon-differentiableproblemofthedesignvariablesolutionexpressedinaclosedformattheLagrangemultiplierpointsforthedualproblem,amapfunctionisintroduced,thenon-smoothexpressionsofthedesignvariablesolutionsarereformulatedassmoothingfunctionsoftheLagrangemultipliersbyusingthisnovelsmoothingfunction.Aseriesofapproximatemodelsanditsalgorithmaredevelopedandimprovedfortheevolutionarytopologyoptimizationproblemwithmultiplediscementconstraints,basedonthesmoothdualalgorithmandvariedconstraintlimits.Comparisonsandysesofnumericalsimulationresultsobtainedbyadoptingtheproposedevolutionarytopologyproposedsmoothdualalgorithm.Thenumericalsimulationresultsshowthattheproposedevolutionaltopologyoptimizationmethodwithavarieddiscementconstraintlimitscheme,canprovideaseriesofnearclearprofilesoftopologies,andsolvethetopologyoptimizationproblemswithmultipleconstraints.Anditiseasytoextracttheoptimumandaninferioroptimalstructurefromthesolutionsobtainedintheoptimizationprocess.Examplesalsoshowthattheproposedmethodisfeasible,effectiveandefficientforsolvinglarge-scaleoptimizationproblems.:Constrainednonlinearprogramming;MethodofMovingAsymptotes;DualEvolutionaltopologyoptimization;Convexquadraticprogramming;Evolutionalstructural1.引 拓扑优化方法得到了迅速发展和广泛应用。均匀化法的思想是将拓扑优化问题转化成一个大的设计域上材料分布的尺度优化问题。许多众所周知的拓扑优化方法都了 Penalization(SIMP) 在SIMP方法中，设计变量(ρ)是指相对密度，它是在ρmin和1ρ 代表空单元。，Bendsøe提出对设计变量的中间值用幂函数插值方法惩罚。另一方面，ESO方法则采用了二进制设计变量ρ0,1，因为它的二进制特点，所以在ESO方法中不需要对设计变量 现，在优化过获得的结果是一系列清晰的拓扑构型（没有灰域），所以很容条件,日乘子几乎等于零。根据结构体积或者重量最小优化要求,结构所有的单元近年来 18]等人结合有理式RAMP材料模型和对偶理论，提出了变约束限法和一些有不同约束统对偶问题的闭式解中边界限约束导致设计变量xi(λ)对日乘子的不可微问题，本文旨在改进和完善基于近似和对偶理论的多位移约束拓扑优化方法,且能得一系列清晰拓扑构型。首先，基于分离的、凸近似方法和光滑近似设计变量函数，提出了一套新的解决非线性优化问题的对偶方法，目的是提高优化的效率和可行性。其次，结合光滑对偶方法和变位移约束限法，发展了多位移约束的渐进拓扑优化方法。给出了一些仿真算例的比较和分析，验证了所提出方法的可行性和有效性。将第e号单元的拓扑变量设为ρe,用过滤gv(ρe来识别单元体积gk(ρe来识g(ρ)(ρ)αv

v

Vegv(ρe)Ve

Kegk(ρ)K

式中αv,υ是两个给定的经验参数,算例中取αv1.25,υ6。VeKe分别代表着当拓扑变为ρ时相应的单元重量和单元刚度矩阵VeKe分别表示单元固有重量和单 刚度矩阵将1/ρe作为设计变量,位移的导数可写(u)/(1/ρ)ruj(u)/(1/ρ) ruj(jue)TKe(ue r

rru

0ru式中

ru

为第r组实载荷作用的结构第j个自由度方向的位移Ke表示e1时第e号单0的刚度矩阵ju)表示单位载荷Fj引起的位移(rue及jue分别表示(ru)及ju中只与第e号单元所有自由度对应的元素组成的单元位移矢量,故在单元级能计算方程(3)的右0本文中总数。可设计单元的单元编号为nv(v1,2,Q，拓扑变量表示为ρnv，上限为1，下限为ρmin0.00001。不可设计单元编号为

21]，对每一个结构性能函数都可以按照移动渐近线方法（MMA）以式(4)进行近似展开u~(ρ)u(ρ(k))Qp(k) )Qq(k) vr r νj,rU(k) U(k (k νj,r L(k (k) (kv

(k

(k)2

ρ

ρ

)p(k)U

ρ

0.001

nq(k)ρ(k

(l)2

v ν

Lν 1.001

ρmaxρmin

(k

ρmax、ρmin分别为ρ的上限值和下限值

ruj

代表零

的最大值，而

ruj(ρ(k))

ruj(ρ(k)代表零与 中的最大值。U(k)、L(k)为ρn的渐近

、下限[19-21]，其他符号参见文献[15-18]v利用式（4）可以近似得到ruj(ρ)关于设计变xv1/n(v1,2,Q的二阶偏导vv表示式见式（7）（式中ϕr105v

2ruj(xv )2ruj(xv ) ) vj,r v(U(k)ρ(k) (U(k)ρ(k) 2(ρ(k))3q(k 2(ρ(k))4q(k

(k vj,r vj,r

u

) 2ruj第r组载荷作用下，2ruj u~(x)u(x(k))

ruj(x(k)ruj(x(k)

1

(xvx(k) ruruj(x(k)

(1/ρ

r2ruj(xv

。2ruj。 n(1/ρnv

2 (xv 位移约束、体积最小问题的结构拓扑优化基本模型可表示

rj/rujnρminnv

)Uj/ruj nv

(j1,,J;r1,2,, v1,2,,ru0(ρ(0)式中Uj表示结构第j个约束方向的位移约束限值

为第r组实载荷作用的结第j个自由度方向的初始位移值；位移约束点的个数用J表示，V表示优化结构的体积minV

)αvV0P(ρ)αvV i ss u 0 ku ) rjr

u0

(0)

(j1,,J;r1,2,, ρminρj1rj(v1,2,, i式中，V0表示不可设计单元的原始体积值,V0表示可设计单元的原始体积值。其中第i 迭代步的k1为引进的位移约束的变化的序列上限值，可表示为 U

(k)

(k)

Uk1rj

min(βruj

,αUjruj

rujραU (k

(k

(k

rj

min(βruj

,ruj

r 式中,α、β均为变约束限的经验参数，本文算例中取松驰参数α=1.03；β按式(12)值β

sn sn位移约束函数在设计变量x(k)1/ρ(k)(v1,2,,Q)处表示成二阶级数形式，得到式(13 )所示的近似优化模min

V

0

0 u~(x)/u

)

/u

) (j1,,J;r1,2,, r

r

U,k r

xv v1,,Q12~(k

~(k (k) ~(k))

V

(k)

V V(x)V~

v1(k

(xv

(xvxv)ruruj(x(k)

1

(xv 22ruj

(k)

)u(x))

u

)

(xx

(xvx

r r

2 (x V(x(k))

(ρ(k)αv

(k)αv0 n Vn(ρi)

/ p

~(k)V

(k)(α1) αv(xv V/max(V

2V(x(k

)(αv2)V0/max(V0 (xvrujruj(x(k)j

Dv,((r1)Jj)

/ru0(ρ(0))

v1,2,,Q;j1,2,,Q;r1,2,,

(ruj(x(k))Uk1)/ru0(ρ(0)

((r1)J 式中，xU,k1/ρmin，xL,k1,v1,2,,Q。并定义 JL λ参考Svanberg[20方法，引进人工变量w1,2,,mT,w0扩大求解范围，既使当模型(13)不存在非空的可行解集合时,也能获得一个违背约束条件较小的解。大量数值仿真试验表明,大部分情况下,这个采用了变约束限的优化模型(13)的非可行解仍属于原优化问题模型(9)的解。另外,参照文献[21,22],λ利用上一迭代步得到的近似优化模型的日乘子，将约束条件函数的二阶量引入目标函数中。去掉式(13)中目标函数的常数项,将模型(13)转化为式(15)效二次规划优min (k s.t.:Dv,j xv)v wj,j1,2,, L,k

U ,v1,2,, x式中QQ1

2~(kV(x)

)LJ

0 /ruj0

(k)x(k))2

x(k))v1

(x

(x

)(x

r1

(c

1dw2

λ其中(kλ(r1)J

为上一迭代步的相应的日乘子值。

jj

2j2~(k)

LJ(k

~(kT(k)V

2ruj

)/u0(ρ(0) S(k)V (x

(r1)J

(x

) v r1v

vvmax:φ( 式中φ(λ)

L(x,λ);xL,kxxU

λv

Q

L(x,λ)

T(k)(xx(k))2S(k)(xx(k))

(cw

dw2v1

j

2jmλλ

Q (xx(k))k)w

j

Dv v 这里设置dj1,j1,2,mλ,取cj为一个较大的值λ以便优化模型(13)存在可行解时,人工变量ww1w2,wmT对利用模型(15)的解获得λv化模型(13)的可行解没有影响。本文算例中设置cj1000j1,2,mλx(k,l1)(λ)表示在第k外循环的第l1内循环迭代步最优点的这些设计变量和人工变量值，利用Karush-v Kuhn-Tucker(KKT)条件可 x(k,l1)(λ)maxxL,k,minxU,k,x*(k,l1)(λ),v1,2,, λx*(k,l1)(λ)x(k)(S(k)λ

jDvj)/

(k v j dd xiRmλRwjRRλ的连续函数,但由于可能分别出x(k,l(λ)xL,k, x(k,l(λxU,li1,2,Q及λjcj，j1,2,mλxiwjλ 通过函数θ(t) p (k,l z(k,l1)(λ)/ (λ) p1ln(1e 1),v1,2,, λ(k,lU(x*(k,l1)(λ)xU,k)/wj(λj)p1λ(k,lU(x*(k,l1)(λ)xU,k)/

x

p2ln(1

2)v从而x(k,l1)(λ)的导数可表示vx(k,l1) z(k,l1)(λ)/pz(k,l1)v

e z(k,l1)

v

,j1,2,, vz(k,l1)v

*(k,l U,k

e e

(λ) /p2 x*(k,l1) (k Dvj/

wq(wq(λqdj(1

1cj) 1

e(λjcj)/(p1dj),qj,j1,2,, )λ

,q通过式(21)和式(22)近似光滑函数的替代,与式(18)相联的对偶函数φ(λ)相对λ是可微的。参考文献[23]的引理6.1.1,对所有的j1,2,mλq1,2,mλλRmλ,φ(λ) (x*x(k))kw(λ

Dvjv vj 2φ(λ)Q Dv

φ(λ)近m展开成如下二次式 m (λλ(k,l))(λλ(k,l)

λj1λ

s

sλλ(k,l)

j mλmλj1s

mλ

j

λλ (2φ(λ) sλλ(k,l)j sλλ(k,l) j1 j1j式中,λ(k,l)是第k外循环第l内循环的第j个日乘子。采用φ(λ)的二阶近似函数且jminmm

1λTGλHT2λ

hr

1λWW(λW(λλ

2Wλ

(k

λλ(k,l 12W

sλ

(k,l

2Wλ2Wgrs

2

λ

(k,l

λλ(k,l)

对于设计变量个数远大于约束个数的非线性优化问题,式(27)具有对称性。为了利用对称二次矩阵的函数的二次规划求解算法,式(27)至式(28)进的日乘子值λ。由式(20)和式(21-22)可得光滑的设计变量和人工变量值。) 入上述中进行后续计算3.2步骤)0.0p201，ξ0.,1)，设置内循环迭代步l0；采用二次规划法求对偶问题(27)的解λ(k,l1)。计算式（19）及式（21）和式(22)的设计变量x(k,l1)和人工变量w(k,l1)x(k,l1)和人工变量w(k,l1),以及x(k,l1)x(k,l1)的相对误差，如小于ε，则停止，并取x(k1)x(k,l1),转下一个外循环迭代步。否则令ll1p1(l1)ξp1lp2(l1)ξp2l，x(k,l1)w(k,l1)分别作为式(25a)的x*wj(λj)的近似值,转步骤2。3.3人工变量初始迭代值w(1)如果设计变量初始迭代值x(1)给出后，需确定人工变量初始迭代w(1)。如果式（15）所表示的近似二次规划优化模型的设计变λ固定后，用wRn表示的优化问题n个单变量优化问题λ cjwj1dj

vkwjjwjj

j1,2,, j 日乘子的导数j将(λc)/(pd j1,2,,m作为未知量用wmax0.0,vk替代w,并用式(22)反 1 。

j1,2,,mλ。再用式（24b）可得到内循环l1时wj对日乘子的导带有置信区间的序列近似拓扑优化模型和其求为了进一步展示光滑对偶求解算法的特征和带有变位移约束限的模型(13)的特minNVminNV αv Q~rujU (v ρnq

p1 ρ

(k

,J;r1,2,, k 这里q是一个经验参数,本文算例中q在[0.015,0.03类似地,(30)转化为如下的X(k)处的QP子优化问题~(k min(k(k

(x)(cjwj djwj2 2w0,c1,2,, ~(k) Q2

(k L

x ,w22ruj(x(k)(x

~(k)

(k)~(k)(x)

1V2

)

/u0

x(k))2V

(x 2(x

(r1)J

r ) v1

r1 ~(k

(x)

ruj(x)/ruj )Uj/ruj

,c(r1)Jj,j1,2,,J;r1,2,,L ~(k~L,k~U ~(k~L,k~U

],q

其他符号参见第3节。采用文献[11]的方式将式(31)和(32)的目标函数和约束条件中近似函数与所有实体材料单元变量相关的一阶和二阶导数进行过滤修正，然后把过滤后的值重新代入上述中进行后续计算。第3优化模型(31)。算图1为左端固定的一悬Lx0.16mLy0.10m、厚度t=0.001m，作用在自由端中点的铅垂载荷是3kN；弹性模量E＝210GPa，泊桑比ν0.3，材料密度33图1主要的位移约束限变化因子β在优化迭代的第一阶段取0.06；后迭代阶段取0.005结构优化历程，经验参数q取0.015。图2中每一幅图给出了在某一迭代步解的拓扑变从图2可以看到,在所有迭代步,获得的解都包含灰单元(‘grey’ (j)提取的最优拓(05(a)第20迭代步拓(b)第25迭代步拓(c)第40迭代步拓(d)第60迭代步拓(e)第78迭代步拓(f)第95迭代步拓(g)第120迭代步拓(h)第136迭代步拓图4采用本文变约束限的序列近似拓扑优化模型和其算法获得的(的位移约束灵敏度相关联。变约束限的序列近似拓扑优化模型和其算法可以获多位移约束下 型结构的拓扑优化设图6是一个一端固支一端滚支的典型0型示意图。其设计域尺寸、支承和载荷情况如图6所示，将垂直载荷F=1000N加载在A、B和C三点上，板的厚度是2mm，假设载荷作用点处的位移约束均为0.12mm。初始结构域离散成为3200个相同的四节点平面应力有限单元，即划分为8040的有限元网格，其中H=0.5m；假定弹性模量E=210Pa；泊松比ν0.3ρ=7800kg/m3。0体积和位移的优化过程较平稳，最终都达到了最优状态。变约束限的序列近似拓图6Mic型结构设计模 (a)第20迭代步拓 (b)第35迭代步拓 (c)第50迭代步拓 (d)第61迭代步拓(e)第75迭代步拓 (f)第90迭代步拓 (g)第104迭代步拓 (h)最优拓图7采用本文变约束限的序列近似拓扑优化模型和其算法获得的 (总体积结(1)本文引进的函数映射θ能近似地将对偶子问题的闭式解xi(λ)转化为日乘子矢法求解结果。算例表明变约束限的序列近似拓扑优化模型和其算法可以获取参考M.P.Bendsøe,N.Kikuchi,Generatingoptimaltopologiesinstructuraldesignusingahomogenizationmethod,Comput.MethodsAppl.Mech.Eng.,1988,71：197–224.M.P.Bendsøe,O.Sigmund,TopologyOptimization:Theory,MethodsandApplications,Springer,Berlin,Heidelberg,2003.M.P.Bendsøe,Optimalshapedesignasamaterialdistributionproblem,Struct.Optim.，1989，1：Y.M.Xie,G.P.Steven,Asimpleevolutionaryprocedureforstructuraloptimization,Comput.Struct.49(1993)885–886.Y.M.Xie,G.P.Steven,Evolutionarystructuraloptimizationfordynamicproblems,Comput.Struct.,1996，58：1067–1073.J.H.Rong,Y.M.Xie,X.Y.YangandQ.Q.Liang,Topologyoptimizationofstructuresunderdynamicresponseconstraints.JournalofSoundandVibration，2000，234(2)：177-189.J.H.Rong,Y.M.XieandX.Y.Yang,Animprovedmethodforevolutionarystructuraloptimizationagainstbuckling.Computers&Structures，2001，79：253-263.X.Y.Yang,Y.M.Xie,G.P.Steven,O.M.Querin,Bidirectionalevolutionarymethodforstiffnessoptimization,AIAAJ.，1999，37(11)：1483–1488.P.Tanskanen,Theevolutionarystructuraloptimizationmethod:theoreticalaspects,Comput.MethodsAppl.Mech.Eng.，2002，191：5485–5498.X.Huang,Y.M.XieConvergentandmesh-independentsolutionsforthebi-directionalevolutionarystructuraloptimizationmethod.FiniteElementsin ysisandDesign43(14)(2007):1039-1049.HuangX,XieYM.Evolutionarytopologyoptimizationofcontinuumstructureswithanadditionaldiscementconstraint.Struct.Multidiscip.Optim.2010,40:409–16.[12]X.Huang,Y.M.Xie,Bi-directionalevolutionarytopologyoptimizationofcontinuumstructureswithoneormultiplematerials,Comput.Mech.43(3)(2009)393–401.X.Huang,S.W.Zhou,Y.M.Xie,Q.Li,Topologyoptimizationofmicrostructuresofcellularmaterialsandcompositesformacrostructures,Comput.Mater.Sci.67(2013)397–407.JianhuaRong,XiaohuaLiu,