CN109446585B - 一种周期性多材料结构的拓扑优化设计方法 - Google Patents

Abstract

本发明公开了一种周期性多材料结构的拓扑优化设计方法，属于结构优化技术领域。该方法基于双向渐近结构优化法构建周期性多材料结构的拓扑优化模型，通过定义初始化参数、设定材料变换规则、求解敏度数并进行结构优化等步骤，以各材料目标体积及最终收敛条件为判断依据进行循环迭代，从而完成周期性多材料结构的拓扑优化设计。该方法能够使所设计的宏观结构不仅具有周期性结构所拥有的规律性、有序性、美观性及良好工艺性等优点，还可具备多材料结构所拥有的高强度、防腐蚀、抗疲劳、高耐磨、强隔热以及轻量化等多种复合特性，可广泛应用于航空航天飞行器、武器装备、高速轨道交通以及工程机械等领域。

Description

一种周期性多材料结构的拓扑优化设计方法

技术领域

本发明属于结构优化技术领域，更具体地，涉及一种周期性多材料结构的拓扑优化设计方法。

背景技术

多材料结构由两种或两种以上的实体材料构成，通常具有高强度、防腐蚀、抗疲劳、高耐磨、强隔热以及轻量化等卓越的复合特性，因而多材料结构有着宽泛的应用背景，尤其适用于航空航天、飞行器与先进武器装备等综合性能要求高的场合。周期性结构则在宏观上呈现一定的规律性、有序性，具有外形美观、制造成本低、工艺性好及便于模块化组装等显著优点，尤其适用于各类起重机主梁、各式桥梁以及三明治结构等具有大长径比的结构。

微观结构通常具有一定的阵列分布属性，所以目前周期性结构拓扑优化的研究工作多数限于微观结构设计，而在针对周期性宏观结构的设计方面还缺乏系统性的设计方法；多材料结构的研究对象在宏观结构与微观结构上都有所报道。而集中了周期性与多材料两个特征的周期性多材料结构在公开出版物中却鲜有说明。不难知道，该结构特点结合了两类结构的工程优点，而这一复合特性是传统的“实-虚”拓扑优化设计无法获得的。本发明关注于周期性多材料结构设计，可为航空航天、飞行器、武器装备与高速轨道交通等国家重点领域提供一种先进的结构优化设计方法。

采用各项同性材料插值法(简称：SIMP法)设计多材料结构时，研究者一般通过不断增加优化过程中的惩罚指数消除中间密度单元，而这一操作不仅无法完全解决该问题，还容易造成拓扑解的收敛困难。其它的水平集法、相场法等拓扑优化方法在设计多材料结构时，存在异质材料间的边界难于描述与表达，连续相间不同材料有交叉重叠情况等问题。虽然引入非线性扩散技术、哈维赛德映射、Cahn–Hilliard模型等额外约束可以显著克服这些弊端，但无疑增加了整个优化程序的计算成本。

发明内容

针对现有技术存在的以上缺陷或改进需求，本发明提供了一种周期性多材料结构的拓扑优化设计方法，其目的在于对宏观结构进行拓扑优化时，同时考虑了结构的多材料与周期性两个特征，所设计的结构无疑也兼备了多材料结构与周期性结构二者的工程优势，不仅可显著提升结构综合性能，还可充分挖掘材料潜力。

为实现上述目的，本发明提供了一种周期性多材料结构的拓扑优化设计方法，包括以下步骤：

(1)给定设计域，基于双向渐近结构优化法构建周期性多材料结构的拓扑优化模型；

(2)定义初始化设计参数，包括：进化率ER、过滤半径r_min、各相材料的目标体积V_j ^*、周期性特征n₁×n₂，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，并建立周期性多材料结构的有限元模型；

(3)初始时，按材料所需属性强弱由大到小进行排序，满设计由属性最强的材料填充，按照属性强弱大小的顺序依次对各材料进行迭代优化；

对于材料j，令当前体积

表示材料j在上次迭代优化时的体积取值，将V_j与目标体积V_j ^*进行比较：

若

则

否则，V_j＝V_j ^*；

(4)基于步骤(3)的材料体积优化结果，根据边界条件与载荷施加的工况，对步骤(2)的有限元模型进行有限元分析，得到单元位移矩阵；

(5)基于步骤(4)得到的单元位移矩阵进行灵敏度分析，求解单元敏度数；

(6)对单元敏度进行如下约束：要求每个子域对应位置单元的敏度数相同，任一子域中某单元的敏度数设定为所有子域相应单元敏度数的平均值；

(7)对单元敏度数按照由大到小的顺序进行排序，并进行如下判断：如果实体单元i满足α_i＜α_th条件，其单元密度从1转变为x_min，表示将材料j转换为材料j+1；如果虚单元i满足α_i＞α_th条件，其单元密度则从x_min转换为1，表示将材料j+1转换为材料j，从而实现多材料结构中任意相邻两种材料的逐步交换；其中，α_th为敏度去除阀值，引入x_min表示单元的最小密度，以避免单元矩阵的奇异；

(8)判断各相材料的目标体积与设定的收敛条件是否同时满足：

若有材料的目标体积不满足，则返回步骤(3)，继续对各材料的体积逐个进行优化；直到所有材料的目标体积均满足，进入收敛条件判断：

若满足收敛条件，则输出当前的设计结果作为最优的周期性多材料结构；否则，返回步骤(4)，继续对结构应变能进行优化。

进一步地，步骤(1)中基于双向渐近结构优化法构建周期性多材料结构的拓扑优化数学模型为：

最小化:

要满足的约束条件:

x_ij＝x_minor 1

x_1b＝x_2b＝…＝x_ab(a＝1,2,...,n₁,b＝1,2,...,n₂)

式中，K和u分别是结构的刚度矩阵和位移矩阵；C为结构的总应变能，

与

分别表示第j种，第m种材料的目标体积；V_i是第i个单元的当前体积；M为设计域中单元的总数量，x_ij表示第j相材料中第i个单元的密度值，

a表示周期性结构中的第a个元胞，b表示单个元胞中的第b个单元，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，x_min表示单元的最小密度；x_ab表示第a个周期性元胞中第b个单元的材料状态。

进一步地，步骤(4)中，选取的材料属性为弹性模量，对于n相材料，将各材料的弹性模量按由大到小的顺序排序为E₁＞E₂＞…＞E_n，分别对应材料1、材料2、…、材料n的弹性模量。

进一步地，步骤(5)中，对于两相材料，单元i的敏度数α_i的计算公式如下所示：

其中，E₁、E₂分别为两相材料设计中材料1和材料2的弹性模量，k_i和u_i分别代表单元i的刚度矩阵和位移矩阵。

进一步地，步骤(5)中，对于n相材料，相邻材料j和j+1的敏度数α_i的计算公式如下所示：

其中，j∈[1,n-1]，n≥2，

和

分别表示由弹性模量为E_j和E_j+1的第j种材料、第j+1种材料计算所得的单元i的刚度矩阵，u_i为单元i的位移矩阵，E_j>E_j+1。

进一步地，步骤(6)中敏度约束如下：

其中，N为结构的总周期数，α_ab表示第a个周期性元胞中的第b个单元敏度数，

表示第a个元胞中第b个单元的原始敏度。

进一步地，步骤(8)中的收敛条件定义为：

其中，C_k-d+1与C_k-X-d+1分别为第k-d+1次、第k-X-d+1次迭代计算所得的结构应变能，k为当前迭代次数，τ是允许的收敛误差，X是整数。

进一步地，步骤(5)中，从第二次循环开始，将基于单元位移矩阵求取的敏度数与前一次迭代时的敏度数的平均值，作为当前敏度数。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，主要有以下有益效果：

(1)该方法能够使所设计的宏观结构不仅具有周期性结构所拥有的规律性、有序性、美观性及良好工艺性等优点，还可具备多材料结构所拥有的高强度、防腐蚀、抗疲劳、高耐磨、强隔热以及轻量化等多种复合特性，在航空航天飞行器、武器装备、高速轨道交通以及工程机械等领域都有着广泛的应用性。

(2)采用双向渐近结构优化法(Bi-directional evolutionary structuraloptimization，以下简称BESO方法)进行多材料结构设计，其独特的阶段设计策略使得计算结果不仅独立于惩罚因子的选取，且不同材料间具有清晰的设计边界，此外，该方法还具备良好的收敛性与高效的计算精度。

BESO方法设计多材料结构时，其独特的阶段设计策略使得计算结果不仅独立于惩罚因子的选取，且不同材料间具有清晰的设计边界。此外，该方法还具备良好的收敛性与高效的计算精度。因而BESO方法在设计多材料结构时具有独特优势。

附图说明

图1是本发明较佳实施提供的一种周期性多材料结构的拓扑优化设计方法流程图；

图2是周期性结构示意图；

图3采用图1中的一种设计周期性多材料结构的拓扑优化方法进行优化的Michell-type结构的设计域边界条件；

图4是图3中的Michell-type结构优化后得到的周期性为4×1的多材料结构示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

图1为本发明实施例提供的一种周期性多材料结构的拓扑优化设计方法流程图。

(1)给定设计域，基于双向渐近结构优化法(简称：BESO法)构建周期性多材料结构的拓扑优化模型；

基于双向渐近结构优化法构建周期性多材料结构的拓扑优化数学模型为：

minimize:

Subject to:

x_ij＝x_minor 1

x_1b＝x_2b＝…＝x_ab(a＝1,2,...,n₁,b＝1,2,...,n₂)

式中，K和u分别是结构的刚度矩阵和位移矩阵。C为宏观结构的总应变能，

与

分别表示第j种，第m种材料的目标体积；V_i是第i种材料的当前体积；M为设计域中单元的总数量，x_ij表示第j相材料中第i个单元的密度值。a表示周期性结构中的第a个元胞，b表示单个元胞中的第b个单元，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，x_min＝0.001表示单元的最小密度，为了避免刚度矩阵的奇异。x_ab表示第a个周期性元胞中第b个单元的材料状态(有或无)，x_ij可进一步表达为：

(2)定义初始化设计参数，主要包括：进化率ER、过滤半径r_min、惩罚因子p、各相材料的目标体积V_j ^*、周期性特征n₁×n₂等，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，并建立周期性多材料结构的有限元模型；

(3)初始时，按材料所需属性强弱由大到小进行排序，初次迭代由满设计开始，满设计由属性最强的材料填充，按照属性强弱大小的顺序对各材料依次进行迭代优化；

后续对未满足目标体积的材料依次进行迭代优化时，对于本次优化的材料j的当前体积V_j，先令

表示材料j在上次迭代优化时的体积取值，再将V_j与目标体积V_j ^*进行比较：

若

则取

否则，取V_j＝V_j ^*；

(4)基于步骤(3)的材料体积优化结果，根据边界条件与载荷施加的工况，对步骤(2)的有限元模型进行有限元分析，得到单元位移矩阵；

(5)基于步骤(4)得到的单元位移矩阵进行灵敏度分析，求解单元敏度数；

对于两相材料设计，敏度数的计算公式如下所示：

其中，E₁、E₂分别为两相材料设计的弹性模量，k_i和u_i分别代表单元的刚度矩阵和位移矩阵。

对于n相多材料设计，其敏度计算公式可表示为：

其中，j∈[1,n-1]，n≥2，

和

分别表示由弹性模量为E_j和E_j+1的第j种材料、第j+1种材料计算所得的单元i的刚度矩阵，u_i为单元i的位移矩阵，E_j>E_j+1；

(6)为了使单元内材料的删除与添加状态保持一致，根据设定的周期性条件对单元敏度进行约束：要求每个子域对应位置单元的敏度数一样大。任一子域中某单元的敏度数设定为所有子域相应单元敏度数的平均值；

(7)基于单元敏度数计算结果对单元敏度数按照由大到小的顺序进行排序，采用优化准则法更新设计变量，逐渐去除对结构总应变能贡献小的单元：如果实体单元i满足α_i＜α_th条件，其单元密度从1转变为x_min，表示将材料j转换为材料j+1；如果虚单元i满足α_i＞α_th条件，其单元密度则从x_min转换为1，表示将材料j+1转换为材料j，从而实现多材料结构中任意相邻两种材料的逐步交换；其中，α_th为敏度去除阀值，x_min表示单元的最小密度，如此实现结构的不断更新与优化；

(8)判断材料的目标体积与设定的收敛条件是否同时满足：

若有材料的目标体积不满足，则返回步骤(3)，继续对各材料的体积逐个进行优化；直到所有材料的目标体积均满足，进入收敛条件判断：

若满足收敛条件，则输出当前的设计结果作为最优的周期性多材料结构；否则，返回步骤(4)，继续对结构应变能进行优化。

收敛条件定义为：

进一步地，步骤(8)中的收敛条件定义为：

其中，C_k-d+1与C_k-X-d+1分别为第k-d+1次、第k-X-d+1次迭代计算所得的结构应变能，k为当前迭代次数，τ是允许的收敛误差，X是整数，为经验值。

下面结合图2～4所示的一个具体实施例来对本发明的上述步骤进行详细说明：

(1)基于双向渐近结构优化法构建周期性多材料结构的拓扑优化模型，表示为：

minimize:

Subject to:

x_ij＝x_minor1

x_1b＝x_2b＝…＝x_ab(a＝1,2,...,n₁,b＝1,2,...,n₂)

式中，K和u分别是结构的刚度矩阵和位移矩阵。C为宏观结构的总应变能，

与

分别表示第j种，第m种材料的目标体积；V_i是第i种材料的当前体积；M为设计域中单元的总数量，x_ij表示第j相材料中第i个单元的密度值。a表示周期性结构中的第a个元胞，b表示单个元胞中的第b个单元，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，x_min表示单元的最小密度，为了避免刚度矩阵的奇异，本实施例取x_min＝0.001。x_ab表示第a个周期性元胞中第b个单元的材料状态(有或无)。

(2)按照图2～4定义设计参数。以设计目标为三相材料为例，体积比设定为：材料1的体积比为0.3，材料2的体积比为0.2，材料3的体积比为0.5；第一阶段进化率为0.01，第二阶段进化率为0.02。惩罚因子p为3，周期性特征n₁×n₂＝4×1。三种材料的弹性模量分别为E₁＝1，E₂＝0.15，E₃＝1E-06，本实施例中E₃取很小的正数是为了避免刚度矩阵的奇异。

(3)首次迭代时，满设计从材料1开始填充；

后续对不满足目标体积的材料逐个进行迭代优化时，对于本次优化的材料j的当前体积V_j，先令

表示材料j在上次迭代优化时的体积取值，再将V_j与目标体积V_j ^*进行比较：

若

则取

否则，取V_j＝V_j ^*；

(4)图3为实施例的设计域边界条件：长30cm，宽10cm的Michell-type结构，下边界中点处有一集中载荷F＝5N，下边界左右两端固定约束。对设计域进行有限元分析，得到单元位移矩阵u_i；

具体地，根据单元刚度＝单元弹性模量*矩阵常量可以求得单元刚度，其中，矩阵常量是8*8的矩阵，当材料的弹性模量为1，泊松比为0.3时矩阵常量内的元素如下：

求解出单元刚度以后，再根据载荷＝刚度×位移求解单元位移矩阵u_i。本实施例中求解单元位移矩阵u_i时，先对相邻材料j和j+1的弹性模量进行插值，两种相邻材料的弹性模量E_j和E_j+1的插值形式E(x_ij)可表示为：

其中，p代表惩罚因子，由上式可以计算出求解单元刚度时的单元弹性模量。上述公式适用于两相以上材料的计算，由于本实施例是以三相材料为例，本实施例取j＝1、2。

(5)基于单元位移矩阵u_i进行单元敏度数计算，计算公式可表示为：

其中，

和

表示由弹性模量为E_j和E_j+1计算所得的单元刚度矩阵；u_i为单元位移矩阵。该公式适用于任意多相材料结构相邻两相材料之间的敏度计算。由于本实施例是以三相材料为例，取j＝1、2。

本实施例中，从第二次迭代开始，将单元的当前敏度数与前一次迭代敏度数进行平均，获得本次迭代的平均敏度数，敏度数直接决定了结构的变化，为了让结构变化平缓一些，不至于出现剧烈的变化，因而对相邻两次迭代的敏度数进行平均；

(6)根据指定周期数进行敏度约束，具体地，由以下公式得到：

其中，N为结构总周期数，α_ab表示第a个元胞中第b个单元周期性约束后的敏度，

表示第a个元胞中第b个单元的原始敏度。如图2所示为水平方向周期数为4，垂直方向周期数为1，结构总周期数N＝4×1的周期性结构。

(7)更新设计变量，得到新的优化结构。

具体地，基于单元敏度数计算结果对单元敏度数进行排序，采用优化准则法逐渐去除对结构总应变能贡献小的单元。如果实体单元i满足α_i＜α_th条件，其单元密度从1转变为x_min，表示将材料j转换为材料j+1；如果虚单元i满足α_i＞α_th条件，其单元密度则从x_min转换为1，表示将材料j+1转换为材料j，从而实现多材料结构中任意相邻两种材料的逐步交换；其中，α_th为敏度去除阀值，x_min表示单元的最小密度；如此实现材料1、材料2、材料3之间的转换及调整，实现结构的不断更新与优化。

(8)判断各相材料的目标体积与设定的收敛条件是否同时满足：

若有材料的目标体积不满足，则返回步骤(3)，继续对各材料的体积逐个进行优化；直到所有材料的目标体积均满足，进入收敛条件判断：

若满足收敛条件，则输出当前的设计结果作为最优的周期性多材料结构；否则，返回步骤(4)，继续对结构应变能进行优化。

具体地，收敛条件公式如下所示：

其中，C_k-d+1与C_k-X-d+1分别为第k-d+1次，第k-X-d+1次迭代计算所得的结构应变能，k为当前迭代数，τ是允许的收敛误差，本实施例取τ＝0.01％，X是整数，本实施例取X＝5。

优化后的周期性为4×1的三相材料拓扑优化结构如图4所示。其优化过程可以简述为：当n＝3时，若弹性模量E₁＞E₂＞E₃，从材料1的满设计开始(满设计指初始设计时，结构全部充满了材料1，因而材料1的体积比为100％，称之为满设计，此后迭代，材料1逐渐减少)，依据设定的进化率，材料1逐渐减少，相应地，材料2逐渐增加直至满足材料2的目标体积；在此过程中，材料1与材料2随着逐步优化实现交换，该优化阶段相当于传统的“实-虚”设计；由于本方法是基于相邻材料j和j+1的敏度进行优化，所以只能在相邻材料间进行交换，当材料2满足目标体积之后，材料1的体积继续减少，材料2的体积继续增加但增加的量立刻转移给材料3，从而保持材料2的指定用量(即目标体积)不变，材料3的体积用量逐步增加直至达到材料3的目标体积，即材料1与材料3的交换需要通过材料1与材料2交换后，材料2立刻与材料3交换来实现。在整个优化过程中，材料1与材料2之间的交换取决于材料1的指定体积用量与敏度数

同样地，材料2与材料3之间的交换取决于材料2的指定体积用量与敏度数

等所有材料都满足体积分数之后，再判断设定的收敛条件是否满足，当所有材料的目标体积与收敛条件同时满足后，整个优化过程结束，输出最优结果。

本方法采用BESO方法，其独特的阶段设计策略使得计算结果不仅独立于惩罚因子的选取，且不同材料间具有清晰的设计边界。此外，还具备良好的收敛性与高效的计算精度。

本发明提供的一种周期性多材料结构的拓扑优化设计方法，其目的在于对宏观结构进行拓扑优化时，同时考虑了结构的多材料与周期性两个特征，所设计的结构无疑也兼备了多材料结构与周期性结构二者的工程优势，不仅可显著提升结构综合性能，还可充分挖掘材料潜力。

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种周期性多材料结构的拓扑优化设计方法，其特征在于，包括以下步骤：

(1)给定设计域，基于双向渐近结构优化法构建周期性多材料结构的拓扑优化模型；

(2)定义初始化设计参数，包括：进化率ER、过滤半径r_min、各相材料的目标体积V_j ^*、周期性特征n₁×n₂，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，并建立周期性多材料结构的有限元模型；

(3)初始时，按材料所需属性强弱由大到小进行排序，满设计由属性最强的材料填充，按照属性强弱大小的顺序依次对各材料进行迭代优化；

对于材料j，令当前体积

表示材料j在上次迭代优化时的体积取值，将V_j与目标体积V_j ^*进行比较：

若

则

否则，V_j＝V_j ^*；

(4)基于步骤(3)的材料体积优化结果，根据边界条件与载荷施加的工况，对步骤(2)的有限元模型进行有限元分析，得到单元位移矩阵；

(5)基于步骤(4)得到的单元位移矩阵进行灵敏度分析，求解单元敏度数；

(6)对单元敏度进行如下约束：要求每个子域对应位置单元的敏度数相同，任一子域中某单元的敏度数设定为所有子域相应单元敏度数的平均值；

(7)对单元敏度数按照由大到小的顺序进行排序，并进行如下判断：如果实体单元i的敏度数α_i满足α_i＜α_th条件，其单元密度从1转变为x_min，表示将材料j转换为材料j+1；如果虚单元i满足α_i＞α_th条件，其单元密度则从x_min转换为1，表示将材料j+1转换为材料j，从而实现多材料结构中任意相邻两种材料的逐步交换；其中，α_th为敏度去除阀值，引入x_min表示单元的最小密度，以避免单元矩阵的奇异；

(8)判断各相材料的目标体积与设定的收敛条件是否同时满足：

若有材料的目标体积不满足，则返回步骤(3)，继续对各材料的体积逐个进行优化；直到所有材料的目标体积均满足，进入收敛条件判断：

若满足收敛条件，则输出当前的设计结果作为最优的周期性多材料结构；否则，返回步骤(4)，继续对结构应变能进行优化；

步骤(1)中基于双向渐近结构优化法构建周期性多材料结构的拓扑优化数学模型为：

最小化:

要满足的约束条件:

x_ij＝x_minor 1

x_1b＝x_2b＝…＝x_ab，a＝1,2,...,n₁,b＝1,2,...,n₂

式中，K和u分别是结构的刚度矩阵和位移矩阵；C为结构的总应变能，

与

分别表示第j种，第m种材料的目标体积；V_i是第i个单元的当前体积；M为设计域中单元的总数量，x_ij表示第j相材料中第i个单元的密度值，

a表示周期性结构中的第a个元胞，b表示单个元胞中的第b个单元，n₁为结构在水平方向的周期数，n₂为结构在垂直方向的周期数，x_min表示单元的最小密度；x_ab表示第a个周期性元胞中第b个单元的材料状态；E_j和E_j+1分别表示第j种材料、第j+1种材料的弹性模量。

2.如权利要求1所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(4)中，选取的材料属性为弹性模量，对于n相材料，将各材料的弹性模量按由大到小的顺序排序为E₁＞E₂＞…＞E_n，分别对应材料1、材料2、…、材料n的弹性模量。

3.如权利要求2所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(5)中，对于n相材料，相邻材料j和j+1的敏度数α_i的计算公式如下所示：

其中，j∈[1,n-1]，n≥2，

和

分别表示由弹性模量为E_j和E_j+1的第j种材料、第j+1种材料计算所得的单元i的刚度矩阵，u_i为单元i的位移矩阵，E_j>E_j+1，p为惩罚因子。

4.如权利要求3所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(5)中，对于两相材料，单元i的敏度数α_i的计算公式如下所示：

其中，E₁、E₂分别为两相材料设计中材料1和材料2的弹性模量，k_i和u_i分别代表单元i的刚度矩阵和位移矩阵。

5.如权利要求3所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(6)中敏度约束如下：

其中，N为结构的总周期数，α_ab表示第a个周期性元胞中的第b个单元敏度数，

表示第a个元胞中第b个单元的原始敏度。

6.如权利要求5所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(8)中的收敛条件定义为：

其中，C_k-d+1与C_k-X-d+1分别为第k-d+1次、第k-X-d+1次迭代计算所得的结构应变能，k为当前迭代次数，τ是允许的收敛误差，X是整数。

7.如权利要求1～6任意一项所述的一种周期性多材料结构的拓扑优化设计方法，其特征在于，步骤(5)中，从第二次循环开始，将基于单元位移矩阵求取的敏度数与前一次迭代时的敏度数的平均值，作为当前敏度数。

Images