CN107273580B - 一种确定多相双模量材料布局问题体积约束的方法 - Google Patents

Abstract

本发明公开的一种确定多相双模量材料布局问题体积约束的方法，具体为：步骤1、建立设计域的有限元模型，初始化算法中的参数，令外迭代编号N_out＝1；步骤2、经步骤1后，定义优化模型，按照当前设计变量及约束条件初始化优化问题；步骤3、经步骤2后，对第i次优化问题寻找当前体积约束下多相双模量材料最优布局；步骤4、经步骤3后，判断拉、压区域材料分布，若所有材料的体积约束满足拉、压区域分离，则继续进行下一步，即转入步骤5；反之，则令N_out＝N_out+1，调整双模量材料体积率，更新结构，返回到步骤2；步骤5、停止并保存，完成优化。本发明的方法主要能针对材料相数不超过3的多相双模量布局优化。

Description

一种确定多相双模量材料布局问题体积约束的方法

技术领域

本发明属于多相材料体积率获取方法技术领域，具体涉及一种确定多相双模量材料布局问题体积约束的方法。

背景技术

多相材料通常由多种材料组合而成，且能超越其中任何一种材料性能。

对多相材料布局优化的数学模型可描述为如下形式：

在上式中：ρ_e为设计域中单元e的设计变量；v_e为单元e的体积；ρ_i,e表示在单元e中前i相材料的总体积率，对于i＝1的情况，ρ₀为0；F和U分别表示总体节点力向量和总体节点位移向量；

为最终结构中i相材料的体积率；K为刚度矩阵，K和F可通过有限元方程求解。

通常，对于以材料体积率为约束的多相材料拓扑优化问题，每一相材料的体积率是根据结构需要人为给定。当材料为双模量材料时，采用材料替换法解决双模量材料结构的拓扑优化问题，以避免由于材料非线性导致计算量巨大的结构再分析过程，采用两个各向同性材料近似原有的分段线性双模量材料，基于局部应力状态进行设计变量更新。

对于设计域内的多种材料，需要进行插值，采用的是传统的变密度法。变密度法为结构拓扑优化的有效方法，其插值模型将离散的设计问题转化为连续设计变量的求解；为了消除中间密度的数量，引入人工惩罚方案进行计算。

当设计域中的材料总数为“k”时，等效弹性模量表示为如下形式：

式中，E_k为第k相材料的实际弹性模量，其中的设计变量为ρ₁,ρ₂,…,ρ_k-1，ρ_k-1，其表示单元中，前k-1相材料占前k相材料在单元中的比例，设计变量的个数为(k-1)×N_e。

对于材料i，在单元e中所占的实际比例可如图1所示，可用ζ_i,j表示为如下形式：

为了描述材料在抗拉性能和抗压性能之间的差异，定义拉压模量比，具体如下：

其中，E_T和E_C分别为材料的拉模量和压模量。

为了避免多相双模量材料布局问题中拉压材料体积率定义不合理的现象，需采用一种合理的方法，在优化循环外部加入一层循环，在内层循环中进行多相双模量材料布局优化，在外循环中调整拉压特性不同的双模量材料比例。

由于基于传统变密度法求解多相材料布局问题难于收敛，当求解体积率未知的多相双模量布局优化将会出现运算量过大无法求解的情况，难以解决大规模工程实际问题，并且为了避免材料相数过多引起的制造困难，需针对材料相数不超过3的多相双模量布局优化。

发明内容

本发明的目的在于提供一种确定多相双模量材料布局问题体积约束的方法，能针对材料相数不超过3的多相双模量布局优化。

本发明所采用的技术方案是，一种确定多相双模量材料布局问题体积约束的方法，具体按照以下步骤实施：

步骤1、建立设计域的有限元模型，初始化算法中的参数，令外迭代编号N_out＝1；

步骤2、经步骤1后，定义优化模型，按照当前设计变量及约束条件初始化优化问题；

步骤3、经步骤2后，对第i次优化问题寻找当前体积约束下多相双模量材料最优布局；

步骤4、经步骤3后，判断拉、压区域材料分布，若所有材料的体积约束满足拉、压区域分离，则继续进行下一步，即转入步骤5；反之，则令N_out＝N_out+1，调整双模量材料体积率，更新结构，返回到步骤2；

步骤5、停止并保存，完成优化。

本发明的特点还在于：

步骤1具体按照以下步骤实施：

步骤1.1、建立设计域的有限元模型，初始化算法中的参数，对材料按照拉模量的大小进行降序排列，使其弹性模量满足如下关系：

E₁＞E₂；

步骤1.2、经步骤1.1后，按照工程给定的材料体积率范围定义双模量材料体积率初值，令拉模量较强材料体积率为V_T，压模量较强材料体积率为V_C。

步骤3具体按照以下步骤实施：

步骤3.1、按照优化模型，对结构进行有限元求解；

步骤3.2、经步骤3.1后，根据有限元求解得到的结果，使用材料置换方案；

材料置换方案的目的是：将结构重分析拆分并融入优化中的设计变量更新过程；尽管双模量材料被替换为两相各向同性材料，即：拉模量和压模量，实际在每一次结构变形分析中只使用其中的一个模量，但考虑到双模量材料的应力依赖特性，拉压模量的选择要根据局部应力状态确定；对于给定的设计域，σ_s和ε_s表示单元三个方向的主应力和主应变，其中s＝1,2,3；

若σ₁＞σ₂＞σ₃＞0，单元受拉；

若0≥σ₁＞σ₂＞σ₃，单元受压；

当单元处于复杂应力状态下时，单元的拉压状态能通过计算其拉应变能密度SED_T和压应变能密度SED_C来获得，具体如下：

式中，拉应变能密度和压应变能密度定义如下形式：

式(1)和式(2)中：N_G表示等参元中一维高斯积分点的个数；σ_s和ε_s、表示单元三个方向的主应力和主应变；w_G为高斯积分点的加权系数，其中s＝1,2,3；

当σ₁＞σ₂＞σ₃＞0时，第r相材料的拉模量被选择，并根据式(1)计算拉应变能密度；

当0≥σ₁＞σ₂＞σ₃时，第r相材料的压模量被选择，并根据式(2)计算压应变能密度；

当单元处于复杂应力状态时，若σ₁·σ₃＜0，材料弹性模量的选择取决于单元的拉压应变能密度大小，拉压应变能密度的大小取决于第二主应力的值；

步骤3.3、经步骤3.2后，修正局部刚度：

步骤3.4、使用MMA算法进行设计变量更新，计算并更新每个单元对应的弹性模量；

步骤3.5、经步骤3.4后，判断收敛：

当不满足收敛准则或迭代次数未达到上限时，转至步骤3.1，反之保存并终止；

收敛准则定义为结构柔度的收敛容差：

整数kn＝4表示4个连续的迭代步都满足收敛条件；

当收敛准则满足时，迭代终止；

η为较小的正数，取为0.01；同时，迭代次数上限设为100步。

步骤3.3具体按照以下步骤实施：

设计变量在第N_out次迭代中，单元在复杂应力状态下的描述如下：

另一个为“旧”SED，也叫做有效SED，即单元在原始双模量材料填充时的相同应力状态的描述，定义为如下形式：

其中，ψ_r为r相材料在单元中所占的比例，具体算法如下：

若单元在第N_out-1次迭代中呈现压模量且在第N_out次迭代中为拉模量，则有如下关系式：

式中，

表示第r相双模量材料的拉压比；

若单元在第N_out-1次迭代中呈现拉模量且在第N_out次迭代中为压模量，则有如下关系式：

当单元处于纯拉和纯压状态下两个应变能密度值相等，当处于复杂应力状态下，两个应变能密度值不同，为描述其差异大小定义如下形式：

单元e中的双模量材料被“新的”各向同性材料替代后的刚度矩阵，具体表示为如下形式：

式中，B_e为应变位移算子，D_e为替换为多相各向同性材料后单元的弹性矩阵；

修正后的刚度矩阵表示为如下形式：

修正后的局部刚度为对单元e的设计变量的二次调整，并用于计算目标方程，目标函数则修正为如下形式：

步骤4具体按照以下步骤实施：

步骤4.1、对优化完成的结构进行后处理，按照当前体积率及每个单元中所含各材料的百分比，过滤掉每个单元中的次要材料，使每个单元仅填充一种实体材料；

步骤4.2、经步骤4.1后，在设计域内对所有单元进行搜索，判断单元处于受拉状态、受压状态、还是复杂应力状态，并进行标记；

步骤4.3、经步骤4.2后，对受拉区域所有单元判断，若存在抗压模量较强材料，则标记，并计算出受拉区抗压材料总体积V_T-C；

对受压区域所有单元判断，若存在抗拉模量较强材料，则标记，并计算出受压区抗拉材料总体积V_C-T；

步骤4.4、经步骤4.3后，判断若V_T-C＝V_C-T＝0，则进入下一步，即进入步骤5；否则，判断若V_T-C>0，则减少抗压材料体积率V_T-C/2，使抗压材料体积率为V_T＝V_T-V_T-C/2，若V_C-T>0，则减少抗拉材料体积率V_C-T/2，使抗拉材料体积率为V_T＝V_T-V_T-C/2。

本发明的有益效果在于：

本发明一种确定多相双模量材料布局问题体积约束的方法，该方法在实际应用中主要针对材料相数不超过3的多相双模量布局优化，克服了现有技术中的难题。

本发明一种确定多相双模量材料布局问题体积约束的方法，能够更大程度的发挥材料拉压特性，在计算中得到材料体积率，相比于预先设定材料体积约束的方法更加准确布置材料。

附图说明

图1是材料i在单元e中所占的实际比例分布图；

图2是材料第一种应力状态示意图；

图3是材料第二种应力状态示意图；

图4是材料第三种应力状态示意图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

本发明一种确定多相双模量材料布局问题体积约束的方法，具体按照以下步骤实施：

步骤1、建立设计域的有限元模型，初始化算法中的参数，令外迭代编号N_out＝1，具体按照以下步骤实施：

步骤1.1、建立设计域的有限元模型，初始化算法中的参数，对材料按照拉模量的大小进行降序排列，使其弹性模量满足如下关系：

E₁＞E₂；

步骤1.2、经步骤1.1后，按照工程给定的材料体积率范围(实际需要)定义双模量材料体积率初值，令拉模量较强材料体积率为V_T，压模量较强材料体积率为V_C。

步骤2、经步骤1后，定义优化模型，按照当前设计变量及约束条件初始化优化问题。

步骤3、经步骤2后，对第i次优化问题寻找当前体积约束下多相双模量材料最优布局，如图1所示，具体按照以下步骤实施：

步骤3.1、按照优化模型，对结构进行有限元求解；

步骤3.2、经步骤3.1后，根据有限元求解得到的结果，使用材料置换方案；

材料置换方案的目的是：将结构重分析拆分并融入优化中的设计变量更新过程；尽管双模量材料被替换为两相各向同性材料(拉模量和压模量)，实际在每一次结构变形分析中只使用其中的一个模量，但考虑到双模量材料的应力依赖特性，拉压模量的选择要根据局部应力状态确定；对于给定的设计域，σ_s和ε_s(s＝1,2,3)表示单元三个方向的主应力和主应变；

若σ₁＞σ₂＞σ₃＞0，如图2所示，单元受拉；

若0≥σ₁＞σ₂＞σ₃，如图3所示，单元受压；

当单元处于复杂应力状态下时，如图4中的(a)和(b)所示，单元的拉压状态能通过计算其拉应变能密度SED_T和压应变能密度SED_C来获得，具体如下：

式中，拉应变能密度(StrainEnergyDensities：SED)和压应变能密度(CompressionSED：CSED)定义如下形式：

式(1)和式(2)中：N_G表示等参元中一维高斯积分点的个数；σ_s和ε_s(s＝1,2,3)表示单元三个方向的主应力和主应变；w_G为高斯积分点的加权系数；

当σ₁＞σ₂＞σ₃＞0时，第r相材料的拉模量被选择，并根据式(1)计算拉应变能密度；

当0≥σ₁＞σ₂＞σ₃时，第r相材料的压模量被选择，并根据式(2)计算压应变能密度；

当单元处于复杂应力状态时，若σ₁·σ₃＜0，材料弹性模量的选择取决于单元的拉压应变能密度大小，拉压应变能密度的大小取决于第二主应力的值；

步骤3.3、经步骤3.2后，修正局部刚度：

设计变量在第N_out次迭代中，单元在复杂应力状态下的描述如下：

另一个为“旧”SED，也叫做有效SED，即单元在原始双模量材料填充时的相同应力状态的描述，定义为如下形式：

其中，ψ_r为r相材料在单元中所占的比例，具体算法如下：

若单元在第N_out-1次迭代中呈现压模量且在第N_out次迭代中为拉模量，则有如下关系式：

式中，

表示第r相双模量材料的拉压比；

若单元在第N_out-1次迭代中呈现拉模量且在第N_out次迭代中为压模量，则有如下关系式：

当单元处于纯拉和纯压状态下两个应变能密度值相等，当处于复杂应力状态下，两个应变能密度值不同，为描述其差异大小定义如下形式：

单元e中的双模量材料被“新的”各向同性材料替代后的刚度矩阵，具体表示为如下形式：

式中，B_e为应变位移算子，D_e为替换为多相各向同性材料后单元的弹性矩阵；

修正后的刚度矩阵表示为如下形式：

修正后的局部刚度为对单元e的设计变量的二次调整，并用于计算目标方程，目标函数则修正为如下形式：

步骤3.4、使用MMA算法进行设计变量更新，计算并更新每个单元对应的弹性模量；

步骤3.5、经步骤3.4后，判断收敛：

当不满足收敛准则或迭代次数未达到上限时，转至步骤3.1，反之保存并终止；

收敛准则定义为结构柔度的收敛容差：

整数kn＝4表示4个连续的迭代步都满足收敛条件；

当收敛准则满足时，迭代终止；

η为较小的正数，取为0.01；同时，迭代次数上限设为100步。

步骤4、经步骤3后，判断拉、压区域材料分布，若所有材料的体积约束满足拉、压区域分离，则继续进行下一步，即转入步骤5；反之，则令N_out＝N_out+1，调整双模量材料体积率，更新结构，返回到步骤2，具体按照以下步骤实施：

步骤4.1、对优化完成的结构进行后处理，按照当前体积率及每个单元中所含各材料的百分比，过滤掉每个单元中的次要材料，使每个单元仅填充一种实体材料；

步骤4.2、经步骤4.1后，在设计域内对所有单元进行搜索，判断单元处于受拉状态、受压状态、还是复杂应力状态，并进行标记；

步骤4.3、经步骤4.2后，对受拉区域所有单元判断，若存在抗压模量较强材料，则标记，标出受拉区抗压材料总体积V_T-C；

对受压区域所有单元判断，若存在抗拉模量较强材料，则标记，标出受压区抗拉材料总体积V_C-T；

步骤4.4、经步骤4.3后，判断若V_T-C＝V_C-T＝0，则进入下一步，即进入步骤5；否则，判断若V_T-C>0，则减少抗压材料体积率V_T-C/2，使抗压材料体积率为V_T＝V_T-V_T-C/2，若V_C-T>0，则减少抗拉材料体积率V_C-T/2，使抗拉材料体积率为V_T＝V_T-V_T-C/2。

步骤5、停止并保存，完成优化。

本发明一种确定多相双模量材料布局问题体积约束的方法，该方法在实际应用中主要针对材料相数不超过3的多相双模量布局优化，克服了现有技术中的难题。采用本发明一种确定多相双模量材料布局问题体积约束的方法，可以避免预先设定各相双模量材料体积率而造成的拉压区域分布材料不符合实际需求的情况，在优化过程中调整所需的拉压材料比例，物尽其用，更大限度的发挥各相材料性能。

Claims (3)

1.一种确定多相双模量材料布局问题体积约束的方法，其特征在于，具体按照以下步骤实施：

步骤1、建立设计域的有限元模型，初始化算法中的参数，令外迭代编号N_out＝1；

步骤2、经步骤1后，定义优化模型，按照当前设计变量及约束条件初始化优化问题；

步骤3、经步骤2后，对第i次优化问题寻找当前体积约束下多相双模量材料最优布局；

步骤4、经步骤3后，判断拉、压区域材料分布，若所有材料的体积约束满足拉、压区域分离，则继续进行下一步，即转入步骤5；反之，则令N_out＝N_out+1，调整双模量材料体积率，更新结构，返回到步骤2；

步骤5、停止并保存，完成优化；

所述步骤1具体按照以下步骤实施：

步骤1.1、建立设计域的有限元模型，初始化算法中的参数，对材料按照拉模量的大小进行降序排列，使其弹性模量满足如下关系：

E₁＞E₂；

步骤1.2、经步骤1.1后，按照工程给定的材料体积率范围定义双模量材料体积率初值，令拉模量较强材料体积率为V_T，压模量较强材料体积率为V_C；

所述步骤4具体按照以下步骤实施：

步骤4.1、对优化完成的结构进行后处理，按照当前体积率及每个单元中所含各材料的百分比，过滤掉每个单元中的次要材料，使每个单元仅填充一种实体材料；

步骤4.2、经步骤4.1后，在设计域内对所有单元进行搜索，判断单元处于受拉状态、受压状态、还是复杂应力状态，并进行标记；

步骤4.3、经步骤4.2后，对受拉区域所有单元判断，若存在抗压模量较强材料，则标记，标出受拉区抗压材料总体积V_T-C；

对受压区域所有单元判断，若存在抗拉模量较强材料，则标记，标出受压区抗拉材料总体积V_C-T；

步骤4.4、经步骤4.3后，判断若V_T-C＝V_C-T＝0，则进入下一步，即进入步骤5；否则，判断若V_T-C>0，则减少抗压材料体积率V_T-C/2，使抗压材料体积率为V_T＝V_T-V_T-C/2，若V_C-T>0，则减少抗拉材料体积率V_C-T/2，使抗拉材料体积率为V_T＝V_T-V_T-C/2。

2.根据权利要求1所述的一种确定多相双模量材料布局问题体积约束的方法，其特征在于，所述步骤3具体按照以下步骤实施：

步骤3.1、按照优化模型，对结构进行有限元求解；

步骤3.2、经步骤3.1后，根据有限元求解得到的结果，使用材料置换方案；

材料置换方案的目的是：将结构重分析拆分并融入优化中的设计变量更新过程；尽管双模量材料被替换为两相各向同性材料，即：拉模量和压模量，实际在每一次结构变形分析中只使用其中的一个模量，但考虑到双模量材料的应力依赖特性，拉压模量的选择要根据局部应力状态确定；对于给定的设计域，σ_s和ε_s表示单元三个方向的主应力和主应变，其中s＝1,2,3

若σ₁＞σ₂＞σ₃＞0，单元受拉；

若0≥σ₁＞σ₂＞σ₃，单元受压；

当单元处于复杂应力状态下时，单元的拉压状态能通过计算其拉应变能密度SED_T和压应变能密度SED_C来获得，具体如下：

式中，拉应变能密度和压应变能密度定义如下形式：

式(1)和式(2)中：N_G表示等参元中一维高斯积分点的个数；σ_s和ε_s表示单元三个方向的主应力和主应变，其中s＝1,2,3；w_G为高斯积分点的加权系数；

当σ₁＞σ₂＞σ₃＞0时，第r相材料的拉模量被选择，并根据式(1)计算拉应变能密度；

当0≥σ₁＞σ₂＞σ₃时，第r相材料的压模量被选择，并根据式(2)计算压应变能密度；

当单元处于复杂应力状态时，若σ₁·σ₃＜0，材料弹性模量的选择取决于单元的拉压应变能密度大小，拉压应变能密度的大小取决于第二主应力的值；

步骤3.3、经步骤3.2后，修正局部刚度：

步骤3.4、使用MMA算法进行设计变量更新，计算并更新每个单元对应的弹性模量；

步骤3.5、经步骤3.4后，判断收敛：

当不满足收敛准则或迭代次数未达到上限时，转至步骤3.1，反之保存并终止；

收敛准则定义为结构柔度的收敛容差：

整数kn＝4表示4个连续的迭代步都满足收敛条件；

当收敛准则满足时，迭代终止；

η为较小的正数，取为0.01；同时，迭代次数上限设为100步。

3.根据权利要求2所述的一种确定多相双模量材料布局问题体积约束的方法，其特征在于，所述步骤3.3具体按照以下步骤实施：

设计变量在第N_out次迭代中，单元在复杂应力状态下的描述如下：

另一个为“旧”SED，也叫做有效SED，即单元在原始双模量材料填充时的相同应力状态的描述，定义为如下形式：

其中，ψ_r为r相材料在单元中所占的比例，具体算法如下：

若单元在第N_out-1次迭代中呈现压模量且在第N_out次迭代中为拉模量，则有如下关系式：

式中，

表示第r相双模量材料的拉压比；

若单元在第N_out-1次迭代中呈现拉模量且在第N_out次迭代中为压模量，则有如下关系式：

当单元处于纯拉和纯压状态下两个应变能密度值相等，当处于复杂应力状态下，两个应变能密度值不同，为描述其差异大小定义如下形式：

单元e中的双模量材料被“新的”各向同性材料替代后的刚度矩阵，具体表示为如下形式：

式中，B_e为应变位移算子，D_e为替换为多相各向同性材料后单元的弹性矩阵；

修正后的刚度矩阵表示为如下形式：

修正后的局部刚度为对单元e的设计变量的二次调整，并用于计算目标方程，目标函数则修正为如下形式：

Images