技术领域
[0001] 本
发明属于结构优化相关技术领域，更具体地，涉及一种含有最大尺寸约束的渐进式结构拓扑优化方法。
背景技术
[0002] 结构拓扑优化设计方法是新型的数字化结构设计方式，通过建立包含目标函数与约束方程的数学模型，对结构进行有限元数值分析，在设计域中按照优化准则或者数学规划方法
迭代出满足目标要求的材料分布，结构拓扑优化设计方法广泛应用于
车身设计、航空航天、通信
电子等结构设计。
[0003] 目前，本领域相关技术人员已经做出了一些研究，如文献1“James K.Guest(2009)Imposing maximum length scale in topology optimization 37:463-473”提出一种基于
密度变量法的设计方法，通过定义以单元为中心的圆形或球形影响区域，将最大尺寸约束问题转化为约束每个影响区域的材料体积占比约束问题，实现最大尺寸约束下的结构拓扑优化设计。然而文献1采用的拓扑优化方法基于密度变量法，设计变量可连续在0至1之间取值，结构中存在灰度单元，难以精准评估结构的当前尺寸；且其他已经提出的方法大多计算效率低下，尺寸约束难以被严格满足。
发明内容
[0004] 针对
现有技术的以上
缺陷或改进需求，本发明提供了一种含有最大尺寸约束的渐进式结构拓扑优化方法，其基于现有拓扑优化设计方法的特点，研究及设计了一种兼顾效率及
稳定性的含有最大尺寸约束的渐进式结构拓扑优化方法。所述渐进式结构拓扑优化方法在设计域经过离散后得到离散单元，离散单元的设计标量服从0-1分布；在每一步迭代得到拓扑结构之后，在每个离散单元周围换分区域，将材料分布体积限制在某个
阈值，阈值由最大半径决定以作为尺寸控制的参数，具体是将每个区域的变量通过P范数凝聚为单个变量，再通过引入拉格朗日乘子在原目标函数添加凝聚后的测度。同时，将最大的局部体积作为拉格朗日乘子更新判据，使得整体结构严格满足尺寸约束，如此在结构拓扑优化中引入明确的尺寸约束机制，优化结果显示此方法可以有效地控制整体结构的尺寸，并具有较强的鲁棒性。
[0005] 为实现上述目的，按照本发明的一个方面，提供了一种含有最大尺寸约束的渐进式结构拓扑优化方法，所述渐进式结构拓扑优化方法主要包括以下步骤：
[0006] (1)对待优化结构的整个设计域进行有限元离散，离散得到的每一个离散单元对应一个离散的拓扑设计变量，并以离散单元的形心为圆心将整个设计域划分成多个圆形影响区域，该圆形影响区域为子区域；对离散后的结构在单元
节点上施加
载荷和约束并进行
有限元分析计算得到节点位移矩阵和结构柔顺度值；同时，计算整体结构的局部尺寸测度，即每一个子区域的面积值，并对得到的所有局部尺寸测度进行P范数凝聚以得到全局测度，进而更新拉格朗日乘子λ以求得伴随解，即更新后的全局测度；
[0007] (2)计算每个离散单元的灵敏度数值，并对得到的灵敏度数值进行过滤处理；
[0008] (3)将得到的离散单元的灵敏度值按照大小进行排序，并根据得到的排序结果对待优化结构的材料进行增删和更新；
[0009] (4)对步骤(3)得到的离散单元的灵敏度值进行平均，并计算得到每迭代五步对应的柔顺度之和的变化量Δc，将得到的变化量Δc与1％进行对比，若Δc小于等于1％且待优化结构的许用材料用量分数达到设定值，则结束，由此完成待优化结构的优化；否则转至步骤(1)。
[0010] 进一步地，柔顺度按照下式进行计算：
[0011] c＝UTKU (1)
[0012] 式中，K为整体
刚度矩阵；U为节点位移矩阵；c为柔顺度值。
[0013] 进一步地，步骤(1)中，以每个离散单元的形心为圆心，rmax为半径将设计域划分为圆形影响区域，该圆形影响区域为子区域，并确定形心包含在圆形影响区域内的离散单元是属于该圆形影响区域的，通过计算每个子区域中的实体单元个数来作为约束的判据，若某个子区域被实体单元填满，则认为这种现象是违反最大尺寸约束的。
[0014] 进一步地，位于子区域中心的单元为中心单元，对于中心单元x，计算其子区域的面积φ(x，dmax)，当且仅当：
[0015] max【φ(x，dmax)】＜φcr
[0016] 时，则判定在整个结构中，无违反最大尺寸约束的现象发生，整体结构满足最大尺寸约束；其中，φcr为子区域包含的离散单元个数；每个中心单元的子区域包含的单元数量相等；dmax＝2rmax为直径。
[0017] 进一步地，全局测度φpn由以下公式计算获得：
[0018]
[0019] 式中，nel为离散单元的个数；p为范数
指定的指数值；φj为每个子区域的面积。
[0020] 进一步地，将拉格朗日乘子按照公式(3)和公式(4)进行更新，直至公式(5)满足或wt-wmin≤10-6：
[0021]
[0022]
[0023] max【φ(x，dmax)】＜φcr (5)
[0024] 式中，λt为第t次更新后的拉格朗日乘子；wt为第t次按照二分法寻找最优的拉格朗日乘子所设定的设计变量，w1设置为1，即λ1的初始值为0；max【φ(x，dmax)】为各子区域的面积值的最大值，φcr为子区域包含的离散单元个数，也即最大尺寸测度值；wmin为设定wt的最小阈值。
[0025] 进一步地，离散单元的灵敏度值按照以下公式进行计算：
[0026]
[0027] 式中，xi为第i个离散单元对应的拓扑设计变量；αi、ui分别为第i个有离散单元的灵敏度数值、位移向量； 为第i个有限元离散单元为实体材料时的刚度矩阵；p为范数指定的指数值；φpn为全局测度。
[0028] 进一步地，采用公式(7)对所述灵敏度数值进行过滤：
[0029]
[0030] 式中，灵敏度数值过滤权重因子 由过滤半径rmin和第i个离散单元与第j个人离散单元之间的中心距离Δ(i，j)计算而得；αj为第j个离散单元的灵敏度数值。
[0031] 进一步地，对待优化结构的材料进行更新时采用的更新准则为：
[0032]
[0033] 式中，αth为材料增长和去除的阈值；l为迭代步数； 为第l步迭代时的第i个离散单元对应的拓扑设计变量。
[0034] 进一步地，变化量Δc的计算公式为：
[0035]
[0036] 式中，cl为第l步迭代时的柔顺度。
[0037] 总体而言，通过本发明所构思的以上技术方案与现有技术相比，本发明提供的含有最大尺寸约束的渐进式结构拓扑优化方法主要具有以下有益效果：
[0038] 1.无中间密度单元，由于采用渐进结构优化法，优化结果自始至终保持变量的离散性，相较于SIMP法，物理意义明确，优化出的结果便于制造和加工。
[0039] 2.提出基于区域测度的严格尺寸约束的同时，将多个最大尺寸测度凝聚为单一全局测度，兼具有设计效率和稳健性，结合更新
算法高效快速，因而在划分单元较多的情况下，仍然能保持较高的结构更新效率，且可较快生成满足最大尺寸约束的优化结构。
[0040] 3.过滤半径除可以避免优化结构出现棋盘格现象外，也可以起到控制结构最小局部尺寸的作用；同时调整过滤半径rmin和最大尺寸控制半径rmax，即可达到协同控制优化尺寸的目的，使得整个结构满足理想中的均布尺寸要求，同时得到更多种的构型以满足制造要求。
附图说明
[0041] 图1是本发明提供的含有最大尺寸约束的渐进式结构拓扑优化方法的流程示意图；
[0042] 图2是图1中的含有最大尺寸约束的渐进式结构拓扑优化方法涉及的测度示意图；
[0043] 图3是本发明
实施例提供的中点加载
悬臂梁工况的示意图；
[0044] 图4是无最大尺寸约束的悬臂梁优化结果示意图；
[0045] 图5是含最大尺寸约束的悬臂梁优化结果示意图；
[0046] 图6是引入最大尺寸约束后、构型在优化迭代中的变化示意图，图6中的(a)、(b)、(c)、(d)、(e)、(f)、(g)、(h)、(i)分别是图3中的L型悬臂梁采用图1中的含最大尺寸约束的渐进式结构拓扑优化方法进行优化设计时的第2步、第7步、第15步、第28步、第56步、第64步、第71步、第82步及第113步的拓扑优化设计示意图。
具体实施方式
[0047] 为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
[0048] 请参阅图1及图2，本发明提供的含有最大尺寸约束的渐进式结构拓扑优化方法，所述渐进式结构拓扑优化方法采用离散的拓扑设计变量，拓扑设计变量的取值1或者0对应材料的有无，具有明确的物理意义。同时，所述方法采用快速高效的渐进式结构拓扑优化方法能够保证在有限优化迭代步数内使得优化所得结果满足最大尺寸约束，且兼具高效性及稳定性。
[0049] 所述的含有最大尺寸约束的渐进式结构拓扑优化方法主要包括以下步骤：
[0050] 步骤一，对待优化结构的整个设计域进行有限元离散，离散得到的每一个离散单元对应一个离散的拓扑设计变量，并以离散单元的形心为圆心将整个设计域划分成多个圆形影响区域，该圆形影响区域为子区域；对离散后的结构在单元节点上施加载荷和约束并进行有限元分析计算得到节点位移矩阵和结构柔顺度值；同时，计算整体结构的局部尺寸测度，即每一个子区域的面积值，并对得到的所有局部尺寸测度进行P范数凝聚以得到全局测度，进而更新拉格朗日乘子λ以求得伴随解，即更新后的全局测度。
[0051] 具体地，对待优化结构的整个设计域进行有限元离散，离散得到的每一个离散单元都对应着一个离散的拓扑设计变量(xi＝0/1)，0表示该单元无材料分布，1代表该单元有材料分布，按照待优化结构的实际工况确定载荷和约束，并在编号后的单元节点
自由度矩阵上施加预定的边界条件(即载荷和约束)，并计算整体刚度矩阵及节点位移矩阵，以此计算待优化结构的柔顺度，柔顺度按照下式进行计算：
[0052] c＝UTKU (1)
[0053] 式中，K为整体刚度矩阵；U为节点位移矩阵；c为柔顺度值。
[0054] 同时，以每个离散单元的形心为圆心，rmax为半径(dmax＝2rmax为直径)将设计域划分为圆形影响区域，该圆形影响区域为子区域，并确定形心包含在圆形影响区域内的离散单元是属于该圆形影响区域的，通过计算每个子区域中的实体单元个数，也即区域面积来作为约束的判据，若某个子区域被实体单元填满，则认为这种现象是违反最大尺寸约束的。具体地，位于子区域中心的单元为中心单元，对于中心单元x，计算其子区域的面积φ(x，dmax)，当且仅当：
[0055] max【φ(x，dmax)】＜φcr
[0056] 时，则判定在整个结构中，无违反尺寸的现象发生，整体结构满足最大尺寸约束。其中，φcr为子区域包含的离散单元个数；每个中心单元的子区域包含的单元数量相等。
[0057] 此外，需要计算整体结构的局部尺寸测度，即每一个子区域的面积值，并将得到的所有局部尺寸测度进行P范围凝聚以得到全局测度，进而更新拉格朗日乘子λ以求得伴随解，即更新后的全局测度。具体地，将多个局部尺寸测度凝聚为一个全局测度φpn，全局测度φpn由以下公式计算获得：
[0058]
[0059] 式中，nel为离散单元的个数，p为范数指定的指数值，φj为每个子区域的面积。
[0060] 将拉格朗日乘子按照公式(3)和公式(4)进行更新，直至公式(5)满足或wt-wmin≤10-6：
[0061]
[0062]
[0063] max【φ(x，dmax)】＜φcr (5)
[0064] 式中，λt为第t次更新后的拉格朗日乘子；wt为第t次按照二分法寻找最优的拉格朗日乘子所设定的设计变量，w1设置为1，即λ1的初始值为0；max【φ(x，dmax)】为各子区域的面积值的最大值，φcr为子区域包含的离散单元个数，也即最大尺寸测度值；wmin为设定wt的最小阈值。
[0065] 其中，需要预先确定基本参数，如材料
弹性模量、泊松比、目标体积分数、材料去除率等；准备进行有限元分析时，具体包含有限单元划分、节点编号、刚度矩阵组装、B矩阵组装等；工况加载时，具体为按照节点自由度施加约束和方向
力，为进行的有限元分析求解做准备，同时考虑过滤矩阵和区域面积计算矩阵，为算法迭代更新做准备。
[0066] 步骤二，计算每个离散单元的灵敏度数值，并对得到的灵敏度数值进行过滤处理。
[0067] 具体地，离散单元的灵敏度值按照以下公式进行计算：
[0068]
[0069] 式中，xi为第i个离散单元对应的拓扑设计变量；αi、ui分别为第i个离散单元的灵敏度数值、位移向量； 为第i个有限元离散单元为实体材料时的刚度矩阵。
[0070] 采用公式(7)对所述灵敏度数值进行过滤：
[0071]
[0072] 式中，灵敏度数值过滤权重因子 由过滤半径rmin和第i个离散单元与第j个离散单元之间的中心距离Δ(i，j)计算而得；αj为第j个离散单元的灵敏度数值。
[0073] 步骤三，将得到的离散单元的灵敏度值按照大小进行排序，并根据得到的排序结果对待优化结构的材料进行增删和更新。
[0074] 具体地，对待优化结构的材料进行更新时采用的更新准则为：
[0075]
[0076] 式中，αth为材料增长和去除的阈值；l为迭代步数； 为第l步迭代时的第i个离散单元对应的拓扑设计变量。
[0077] 步骤四，对步骤三得到的离散单元的灵敏度值进行平均，并计算得到每迭代五步对应的柔顺度之和的变化量Δc，将得到的变化量Δc与1％进行对比，若Δc小于等于1％且待优化结构的许用材料用量分数达到设定值，则结束，由此完成待优化结构的优化；否则转至步骤一。
[0078] 具体地，将得到的灵敏度按照公式(9)进行平均：
[0079]
[0080] 其中，n为引入的平均步数，本实施方式中n取2。
[0081] 迭代终止标准为：
[0082]
[0083] 也即计算得到每五步柔顺度之和的变化量Δc，并将Δc与1％进行对比。若大于1％或/及待优化结构的许用材料用量分数未达到设定值，则转至步骤一，否则结束，由此完成待优化结构的优化。其中，cl为第l步迭代时的柔顺度。
[0084] 以下以具体实施例来对本发明进行进一步地详细说明。
[0085] 请参阅图3、图4、图5及图6，以某悬臂梁构型为例，拓扑结构中的区域I和区域II并未完全被材料填充，则判定这两个区域未违反最大尺寸约束，而区域II充满材料，违反最大尺寸约束。
[0086] 如图3为中点加载悬臂梁工况，若未在优化过程中设置最大尺寸约束，只设置过滤半径rmin为15倍网格高度，则最终优化得到图4中所示构型，材料分布不均，整体结构分支粗细不一。
[0087] 在引入最大尺寸约束后得到满足最大尺寸约束的全新优化构型，具体方法包括以下步骤：
[0088] (1)许用材料用量分数及材料去除率分别设定为50％及1％；
[0089] (2)初始设计域为遍布材料的长宽比为1.6的长方形区域，将该设计域均匀离散成320x200共计64000个四节点正方形单元，每一个单元对应一个独立的设计变量(xi＝0/1)，将材料的
杨氏模量和泊松比分别设为1MPa和0.3，悬臂梁左端固定，将集中载荷施加悬臂梁右端中心节点处，对结构进行有限元分析，计算位移和结构整体柔顺度。
[0090] (3)范数参数设为8，通过P范数法凝聚子区域面积得到全局尺寸测度φPN。
[0091] (4)设置最大尺寸半径rmax为20倍网格高度(直径dmax为40倍网格高度)，计算最大尺寸测度，并据此求解伴随问题，更新获得伴随解λ。
[0092] (5)计算每个离散拓扑变量对应的灵敏度数值。
[0093] (6)灵敏度数值的过滤半径rmin设为15倍网格高度，对灵敏度数值进行过滤处理。
[0094] (7)根据定义的材料去除率和体积分数，结合灵敏度数值，计算当前优化迭代步的许用材料用量并更新设计变量。
[0095] (8)根据设计变量值更新有限元模型的单元材料模型，重复上述过程直至许用材料用量分数得到满足且拓扑构型稳定收敛。
[0096] 最终在119步时迭代稳定，得到图5中所示构型，整体结构的材料分布较为均匀，结构分支粗细较为接近。在引入最大尺寸约束后，第2、7、15、28、56、64、71、82、113步的优化结果如图6所示，整体结构在满足尺寸约束的同时，逐渐将边框去除，使得柔顺度值最小化、整体结构的刚度最大化。
[0097] 表1为有无约束时，优化过程中max【φ(x，dmax)】、柔顺度值的变化。对于rmax＝20，尺寸标准设定值φcr＝1245。
[0098] 表1有无约束情况下的迭代结果
[0099]
[0100] 表1中的数据显示，在引入最大尺寸约束后，测度值max【φ(x，dmax)】在第28步时已降到标定的φcr值之下，也即满足最大尺寸约束；在无最大尺寸约束的优化过程中，测度值则一直未降到标定值之下，这表明本方法是有效的。
[0101] 对于整个优化过程中的柔度值，满足最大尺寸约束的结构柔顺度值则大于无约束的结构柔顺度，整体结构柔顺度值也会稳定在一个相对合理的数值，表明本方法是合理的。
[0102] 优化结果显示，采用本发明所提供的方法，以结最大化结构刚度作为优化目标，在引入本发明所述的最大尺寸约束后获得的设计构型在结构上明显不同于纯刚度优化获得的构型；在结构上避免了材料分布的过度集中，结合过滤半径亦可防止材料分布的过度分散，从而形成对结构的协同尺寸控制。与此同时，对比其余已提出的控制尺寸方法，未引入中间密度单元，因而本优化方法从根本上避免了灰度单元的出现，且由于本方法运用灵敏度排序即可决定材料的增删，因而极大地提高了计算效率，对于较细网格的划分计算优势尤其明显。
[0103] 本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何
修改、等同替换和改进等，均应包含在本发明的保护范围之内。
