基于影响区域节点密度插值的结构拓扑优化方法

吉先磊，张美艳，王　皓

(复旦大学 航空航天系，上海 200433)



基于影响区域节点密度插值的结构拓扑优化方法

吉先磊，张美艳，王皓

(复旦大学 航空航天系，上海 200433)

摘要：变密度法是结构拓扑优化中一种常用的处理方法，但在优化过程中会出现棋盘格、网格依赖性、灰度单元等数值不稳定现象.为解决单元变密度法中出现的棋盘格问题，在优化过程中引入节点密度作为设计变量.在提出基于影响区域节点密度插值优化方法的同时，给出了一种基于等密度线的后处理方法.计算结果表明该方法在消除棋盘格现象的同时，处理了网格依赖性问题，并消除了结果中的灰度单元现象.

关键词：拓扑优化； 节点密度； 区域插值

连续体结构拓扑优化是在一个特定的设计区域内，在边界、载荷等约束条件已知的条件下，为达到结构刚度最大、基频最大等目标而寻求材料最优分布的一种结构设计方法.拓扑优化设计方法是目前连续体结构优化设计研究的热点，取得了迅速发展，成为结构设计的一种有效方法.

在拓扑优化众多的求解方法中，变密度法是一种常见的处理方法.在变密度法中通常设定单元密度为设计变量，通过单元密度的迭代获得最优的拓扑结构.然而由于有限元单元密度相互独立，普遍存在棋盘格和网格依赖性等问题[1]，在单元法中需要借助额外的方法如灵敏度过滤[2]、周长约束[3]、局部梯度[4]等进行处理.为了从根本上避免棋盘格问题的出现，可以将设计变量定义在节点位置，使得单元密度场本身具有连续性.基于节点变量处理棋盘格问题的思想可以追溯到Youn等[5]提出的基于节点的密度再分配后处理方法，该方法在单元密度迭代之后，定义节点密度为环绕单元密度的平均值，再由节点密度插值出新的单元密度，通过保证相邻单元密度的连续性避免棋盘格的出现.之后Rahmatalla等[6]，Guest等[7]在此基础上提出了将有限元网格节点密度作为设计变量的优化方法.Kang等[8]，占金青等[9]，Yi等[10]等推广了该方法，将设计变量定义在非网格节点位置进行处理.龙凯等[11]进一步提出了基于物质点描述的优化方法.节点密度法可以避免棋盘格的出现，然而目前对于节点密度下处理网格依赖性和灰度单元问题的研究还较为有限.

本文将设计变量定义在网格节点位置，通过影响区域节点密度插值的方法在消除棋盘格的同时，避免了网格依赖性的出现.

1影响区域节点密度插值方法

变密度法拓扑优化以有限元单元作为优化的基本单位，通过单元伪密度的引入将0，1离散优化问题转换为[0，1]上的连续优化问题.变密度法假定单元的弹性模量与单元伪密度之间的关系满足

(1)

其中：ρe为单元的伪密度，ρe在[0，1]上连续变化，ρe为0时表示该单元物质不存在，ρe为1时表示该单元为实体密度；E0为材料固体部分弹性模量；Ee为单元的等效的弹性模量；p为惩罚参数，文中取为3.

如图1(a)所示(见第312页)，在目前常用的节点密度法中，通常将单元的密度定义为该单元包含节点密度的插值函数.然而随着设计区域网格划分的不同，往往会产生严重的网格依赖性问题，主要表现在随着网格的细分，优化结果中结构最小尺寸会随之变小.在优化过程中为获得更精确的结果往往需要细化网格，而实际应用中又需要对结构中的最小尺寸做出限制.为了达到这样的目标，这里引入影响区域内节点密度插值的方法.如图1(b)所示(见第312页)，在节点密度作为设计变量的基础上，定义单元e的伪密度值为该单元插值区域上节点伪密度的加权平均，即：

ρe=ωe，iρi，

(2)

其中：ρe为加权平均的单元伪密度；ρi为节点伪密度值，是优化过程中的设计变量；ωe，i为对应节点变量的权函数，且

(3)

其中：rmin为限制的结构最小尺寸；dist(e，i)表示单元e中心到节点i的距离；Ωr为距e单元中心距离小于rmin的子区域.

同时为便于获得最终0，1分布的结构，在优化的结果中引入“密度等值线”作为参考.如图2(见第312页)所示的区域，其4个端点位置密度分别为1.0，1.0，1.0，0.001，定义其0.5密度等值线为图中虚线所示线段，定义大于等密度线的部分为实体结构，小于等密度线的部分为孔洞结构，帮助获得最终的优化结构.

2优化模型及灵敏度

基于节点变密度法，以结构刚度最大化为目标，体积为约束，连续体结构拓扑优化的模型可以表示为

(4)

其中： C为目标结构的总体应变能；ue为单元位移向量；U为整体位移向量；F为外在载荷；K为整体刚度矩阵；n为有限元单元数目；f为目标体积分数；V0为设计区域总体积；ρe为单元伪密度；ρi为节点伪密度变量.

基于变密度法，由(1)式，目标函数对有限单元e伪密度的灵敏度表示为

(5)

通过链式求导法则，进而可以表示出目标函数对节点密度变量的灵敏度为

(6)

(7)

3数值算例

下面采用2维平面拓扑优化算例对本文方法进行讨论.

3.1算例1

考虑如图3所示的梁结构刚度最大化问题，结构设计区域的尺寸为60cm×10cm×0.1cm，载荷在上边界中点位置，大小为F=3500N，材料弹性模量为E=210GPa，泊松比为v=0.3，目标体积为整个设计区域体积的50%.

采用平面四边形单元对设计区域进行离散，网格数量为120×20，网格尺寸为0.5cm×0.5cm.分别采用单元节点密度插值、本文的影响区域节点密度插值方法进行优化，两种节点密度方法下的拓扑优化结果对比如图4所示.

如图4(a)所示，在采用单元节点密度插值时，得到结构的最小尺寸与网格的粗细相关联，不可做进一步处理.图4(b-1)、图4(b-2)、图4(b-3)为采用本文影响区域节点密度插值方法获得的结果，通过不同插值区域的选取，在网格固定的条件下，可以明显改变结果中的最小结构尺寸，从而避免网格依赖性对优化结果的影响.

由于节点密度定义的变量值为节点密度值，这里将图4中采用节点密度法得到的结构中的节点密度进行分类，讨论最终结构中节点密度的0，1分布问题.

从表1可以看出，采用影响区域节点密度插值法获得的结果中，随着插值半径的扩大，“1”密度节点个数逐渐增多，即在高密度区域材料分布更为集中.而在低密度区域，“0”密度节点变少，中间密度节点变多，即在这些区域材料分布更为分散化.为此，可以在较大影响半径下采用较高等密度线的方法获得最终的优化结构.如对图4中影响区域半径为1.5cm的优化结果，采用等密度线获得的最终结构如图5所示.

从图5中可以看出，在采用影响区域插值方法时，较大插值半径下，材料更加集中在一些核心区域，在此意义上，等密度线可以在结构设计后期对最终结构提供参考.

表1　结果中密度值分布

3.2算例2

为考察本文方法在三角形单元中的应用，考虑如图6圆盘结构刚度最大化优化问题，圆盘中间固定，半径为50cm，载荷在圆周角30°、90°、150°、210°、270°、330°位置，方向为切向方向，大小为F=100N，材料弹性模量为E=210GPa，泊松比为v=0.3.

以目标体积不超过设计区域体积的60%作为优化目标，采用三角形单元进行离散，采用上文的方法获得的优化结果如图7所示.

从图7可以看出，本文提出的方法同样适用于三角形网格单元，可以有效地处理非规则设计区域，并通过等密度值线获得清晰的拓扑结构.

4结语

针对变密度法拓扑优化中出现的网格依赖性的现象，本文讨论了基于影响区域节点密度插值的优化方法.将设计变量定义在有限元网格节点位置，采用影响区域插值的方法获得有限元单元的伪密度值，通过插值半径的选取控制结构中的最小尺寸，从而处理网格依赖性问题.并提出了基于密度等值线的清晰化方法，得到了0，1分布的优化结果.经过具体的数值算例，证明了本文方法可以在四边形单元和三角形单元中有效地处理网格依赖性问题，获得清晰的拓扑结构.

参考文献：

[1]GEORGE R， TOMASZ L. Topology optimization in structural and continuum mechanics [M]. Germany： Springer， 2014.

[2]SIGMUND O. Design of material structures using topology optimization [D]. Denmark： Department of Solid Mechanics， Technical University of Denmark，1994.

[3]HABER R B， BENDSE M P， JOG C. A new approach to variable-topology shape design using a constraint on the perimeter [J].StructuralOptimization， 1996，11(1)： 1-12.

[4]PETERSSON J， SIGMUND O. Slope constrained topology optimization [J].InternationalJournalforNumericalMethodsinEngineering， 1998，41(8)： 1417-1434.

[5]SUNG K Y， SANG H P. A study on the shape extraction process in the structural topology optimization using homogenized material [J].ComputersandStructures， 1997，62(3)： 527-538.

[6]RAHMATALLA S F， SWAN C C. A Q4/Q4 continuum structural topology optimization implementation [J].StructMultidiscOptim， 2004，27(1)： 130-135.

[7]GUEST J K， PREVOST J H， BELYTSCHKO T. Achieving minimum length scale in topology optimization using nodal design variables and projection functions [J].InternationalJournalforNumericalMethodsinEngineering， 2004，61(2)： 238-254.

[8]KANG Z， WANG Y Q. Structural topology optimization based on non-local shepard interpolation of density field [J].ComputerMethodsinAppliedMechanicsandEngineering， 2011，200(49)： 3515-3525.

[9]占金青，杨康，黄志超.基于节点变量法的连续体结构拓扑优化设计 [J].农业机械学报，2014，45(9)： 329-333.

[10]YI J J， ZENG T， RONG J H，etal. A topology optimization method based on element independent nodal density [J].JCentSouthUniv， 2014，21(2)： 558-566.

[11]龙凯，陈广华.基于物质点描述的双向渐进式拓扑优化方法 [J].工程力学，2012，29(8)： 308-312.

[12]SIGMUND O. A 99 line topology optimization code written in matlab [J].StructMultidiscOptim， 2001，21(2)： 120-127.

[13]BENDSOE M P. Optimization of structural topology， shape and material [M]. Germany： Springer， 1995.

文章编号：0427-7104(2016)03-0310-05

收稿日期：2016-01-05

基金项目：国家自然科学基金(11572089)

作者简介：吉先磊(1990—)，男，硕士研究生；王皓(1963—)，男，副教授，通讯联系人，E-mail： wanghao@fudan.edu.cn.

中图分类号：O 342

文献标志码：A

Topology Optimization Based on Regional Node Density Interpolation Method

JI Xianlei， ZHANG Meiyan， WANG Hao

(DepartmentofAeronauticsandAstronautics，FudanUniversity，Shanghai200433，China)

Abstract：Variable density method is a common approach in topology optimization， while there will be problems like checkerboard， grid dependence， gray cells and other instability occurs in the optimization procedure. In order to solve the checkerboard problem in topology optimization using variable density method， the node density is introduced in the optimization process as a design variable. A new method based on regional node density interpolation is proposed， as well as a node density based on clarity approach is introduced. The results show that the presented node density method can effectively eliminate the grid dependence problem as well as the checkerboard phenomenon in topology optimization.

Keywords：topology optimization； nodal density； regional interpolation