基于等效静态载荷法的B柱多材料拓扑优化方法

技术领域

本发明涉及一种车身B柱设计方法，更确切地说，涉及一种基于等效静态载荷法的B柱多材料拓扑优化方法。

背景技术

近年来，随着车身制造工艺和材料性能的提升，对汽车轻量化的要求也有着进一步的提高，设计良好的车身结构能够在保证汽车安全性的同时最大程度地提高轻量化效果。B柱是汽车车身较为重要的抗撞结构，合理的B柱结构能够提高整车的侧面碰撞安全性能和轻量化程度。

拓扑优化技术为车身结构提供了高效的设计手段，目前，针对单材料拓扑优化技术的研究和应用较为成熟，但随着车身性能和轻量化需求的提高，单材料结构已不能充分满足工程应用的需要。为了最大化地挖掘结构的轻量化潜力，多材料车身结构被提出，并以其优异的性能和更高的设计自由度成为车身结构的发展趋势。通过对B柱结构进行多材料拓扑优化，可以找到多种材料在B柱上的最佳分布位置，为后续的尺寸优化等详细设计提供理论依据，优化后的B柱结构可以提高材料利用率，减轻重量，提高车辆侧面碰撞安全性能。

传统的拓扑优化方法多基于线性静态工况，对于汽车碰撞等高度非线性的动态工况下的优化问题，不能准确的模拟动态过程中结构的受力情况，可能导致拓扑优化结果不能满足实际工程的需要；而对于直接处理动态工况的拓扑优化来说，较高的灵敏度信息计算成本难以驱动拓扑优化的进行。等效静态载荷法是目前较为流行的处理动态载荷的动态拓扑优化方法。通过计算动态工况下结构的等效静态载荷并将其作为静态拓扑优化模型中的边界载荷，实现了动态拓扑优化问题向静态拓扑优化问题的转化。

发明内容

基于背景所述，本发明将多材料拓扑优化方法和等效静态载荷法进行有效结合并应用在车身B柱结构上，提供一种基于等效静态载荷法的B柱多材料拓扑优化方法。本发明针对的是现有利用单材料拓扑优化方法设计的B柱结构设计自由度较小，不能最大化B柱结构轻量化潜力的问题，以及现有静态拓扑优化方法不能很好的模拟侧面碰撞工况下B柱动态过程的问题。

为解决上述技术问题，本发明是采用如下技术方案实现的：

一种基于等效静态载荷法的B柱多材料拓扑优化方法，包括如下步骤：

1)B柱多材料插值模型的建立

基于变密度法中的SIMP法建立多材料插值模型，模型中自变量包含单元的拓扑变量和材料选择变量，因变量为单元的弹性模量；

2)B柱等效静态载荷计算模型的建立，包含

2.1)B柱初始碰撞分析模型

在Hypermesh软件中的LS-Dyna模块下，首先建立侧面碰撞工况下的B柱有限元模型，然后在模型中添加*CONTROL_TERMINATION、*DATABASE_BINARY_D3PLOT、*DATABASE_OPTION、*DATABASE_HISTORY_NODE四个控制卡片；

2.2)B柱初始刚度求解模型

利用Hypermesh软件中的Convert功能对建立的B柱碰撞分析模型进行转化得到Optistruct模块下的B柱模型，然后在模型中添加刚度求解工况；

3)B柱等效静态载荷计算模型的批处理

利用MATLAB平台在初始碰撞分析模型和刚度求解模型对应的k文件和fem文件中写入组件信息、材料信息、属性信息以及单元与组件的对应信息；

4)基于等效静态载荷法的B柱结构多材料拓扑优化模型的建立

q个候选材料的基于等效静态载荷法的B柱结构多材料拓扑优化模型如下：

模型中包含设计变量、约束条件、优化目标；

式中：

ρ为单元密度即拓扑变量，表示每个单元的有无；m

C

为第j个时间步的等效静态载荷，K(b)为刚度矩阵，z

M

5)等效静态载荷的计算及输入

5.1)等效静态载荷的计算

在结构动态非线性平衡方程中响应相同位移的线性静态载荷集，计算式为：

5.2)等效静态载荷的输入

利用MATLAB平台中的fprintf函数计算得到力和力矩两种形式的等效静态载荷，输入到基于等效静态载荷法的B柱结构多材料拓扑优化模型中；

6)模型的收敛判定与设计变量的更新

利用MATLAB平台中的system函数将等效静态载荷输入后的多材料拓扑优化模型提交到Optistruct求解器进行计算，计算后对碰撞分析模型和刚度求解模型中的设计变量进行更新，依次循环，直至满足收敛条件，优化结束。

其中，所述的多材料插值模型为：

式中：E(ρ,m

其中，在建立B柱初始碰撞分析模型时：

约束B柱两端节点的全部自由度，

用*CONTROL_TERMINATION控制碰撞仿真的总时间，

用*DATABASE_BINARY_D3PLOT控制输出动画间隔，二者的比值加1即为时间步的数量，

用*DATABASE_OPTION控制节点位移值、截面峰值力、总能量输出的时间间隔，

用*DATABASE_HISTORY_NODE输出节点位移信息；

用Tool面板下的renumber功能将模型中B柱结构的全部节点ID和单元ID从1开始重新排列。

其中，在利用B柱初始刚度求解模型求解中，对B柱结构的所有节点都施加ASET约束，通过添加PARAM控制卡片，勾选EXTOUT选项并将类型设置为DMIGPCH，输出B柱结构所有节点的静态刚度值。

其中，在B柱等效静态载荷计算模型的批处理中，创建的组件数目与单元数目相同，每个单元的ID与其所属的组件ID、材料ID和属性ID统一。

其中，基于等效静态载荷法的B柱结构多材料拓扑优化模型中，工况的数量根据碰撞分析模型的总时间设置，每个工况下的载荷力和载荷力矩为空，只设置相应节点的位移约束。

其中，有限元方法中，结构的动态非线性平衡方程为：

式中：M(b)为质量矩阵；C(b)为阻尼矩阵；K

其中，在等效静态载荷的计算中，利用MATLAB平台对各个节点的动态位移值和静态刚度值进行提取，得到位移矩阵和刚度矩阵。

其中，模型的收敛判定是：当满足以下两个条件中任意一个时，认为优化结束；

第一个条件是，当循环次数达到最大循环次数时；

第二个条件是，通过计算相邻两次循环最后一个迭代步的单元相对密度变化量Δ

式中：

其中，模型设计变量的更新是：基于MATLAB平台对各个单元材料选择变量进行提取，然后按照公式(7)对初始碰撞分析模型和刚度求解模型各个单元的材料弹性模量和材料密度进行更新，更新后开始下一次循环，直至满足收敛条件：

式中：

与现有技术相比，本发明的有益效果是：

本发明所述的基于等效静态载荷法的B柱多材料拓扑优化方法与现有利用静态单材料拓扑优化方法进行B柱结构设计相比，可以准确的模拟侧面碰撞工况下B柱的动态过程，得到更满足工程实际的多材料B柱拓扑优化结果，同时提高了拓扑优化所带来的经济效益以及优化方案的可行性。

附图说明

图1是基于等效静态载荷法的多材料拓扑优化方法的流程图。

图2是单元弹性模量与两个设计变量之间的关系图。

图3是B柱初始刚度求解模型中待输入的信息截图。

图4(a)是B柱初始刚度求解模型中输入后的组件信息截图。

图4(b)是B柱初始刚度求解模型中输入后的材料信息截图。

图4(c)是B柱初始刚度求解模型中输入后的属性信息截图。

图4(d)是B柱初始刚度求解模型中输入后的单元与组件对应信息截图。

图5是B柱结构动态位移信息截图。

图6是B柱结构静态刚度信息截图。

图7(a)是提取后的位移矩阵图。

图7(b)是提取后的刚度矩阵图。

图8是等效静态载荷矩阵图。

图9是基于等效静态载荷法的B柱结构多材料拓扑优化模型中待输入的载荷信息图。

图10是等效静态载荷输入后的B柱结构多材料拓扑优化模型中的载荷信息截图。

图11是多材料拓扑优化过程中的单元相对密度信息截图。

图12是实施例提供的B柱碰撞分析模型图。

图13是实施例提供的B柱刚度求解模型图。

图14是不同视角下基于等效静态载荷法的B柱多材料拓扑优化结果展示。

图15(a)是实施例提供的该乘用车材料替换前的B柱结构图。

图15(b)是实施例提供的该乘用车材料替换后的B柱结构图。

具体实施方式

下面结合附图和实施例对本发明进行详细的描述，本领域的技术人员应该知道，以下实施例并不是对本发明技术方案作的唯一限定，凡是在本发明技术方案精神实质下所做的任何等同变换或改动，均应视为属于本发明的保护范围。

目前B柱结构拓扑优化方法往往是利用相关商业化软件得到人为简化的静态工况下的单材料结构，这种方法并没有充分地挖掘B柱结构的轻量化潜力以及抗撞性能。因此为了发挥多材料B柱结构的优势并充分考虑侧面碰撞中的非线性问题，本发明提供一种基于等效静态载荷法的B柱多材料拓扑优化方法。

基于等效静态载荷法的B柱多材料拓扑优化方法的流程如图1所示，包括如下步骤：

1.B柱多材料插值模型的建立

本发明选用的拓扑优化方法是基于变密度法中的SIMP法，SIMP法又称为各项正交惩罚材料密度法。该方法以每个单元的密度作为设计变量，在优化过程中，通过改变每个单元的密度，使结构中每个单元的杨氏模量发生变化，最终达到材料在设计空间中的最优分布。在实现多材料B柱结构的过程中，通过引入新的材料选择变量来实现对材料种类的选择。

以两种候选材料(下文称为材料1和材料2)的拓扑优化为例，其材料插值模型为：

式中：E(ρ,m

式中：m

2.B柱等效静态载荷计算模型的建立

本发明将等效静态载荷法与多材料拓扑优化方法相结合，进行基于等效静态载荷法的B柱多材料拓扑优化。首先建立B柱等效静态载荷计算模型，即B柱初始碰撞分析模型和B柱初始刚度求解模型用于等效静态载荷的计算。

2.1)B柱初始碰撞分析模型

碰撞分析模型即B柱结构在碰撞等动态工况下的有限元仿真模型，用于B柱结构节点位移的求解。

该模型在前处理软件Hypermesh中的LS-Dyna模块下建立，首先建立侧面碰撞工况下的B柱有限元模型，即带有质量的刚性体以一定速度撞击B柱。设置好接触条件与边界条件，将刚性体作为主动部件，将B柱作为从动部件，约束B柱两端节点的全部自由度。然后在模型中添加*CONTROL_TERMINATION、*DATABASE_BINARY_D3PLOT、*DATABASE_OPTION、*DATABASE_HISTORY_NODE四个控制卡片。其中，*CONTROL_TERMINATION卡片用来控制碰撞仿真的总时间；*DATABASE_BINARY_D3PLOT卡片用来控制输出动画间隔，二者的比值加1即为时间步的数量，例如，碰撞总时间设置为100ms，输出动画间隔设置为1ms，加上初始0时刻的时间步，则总共有101个时间步；*DATABASE_OPTION卡片用来控制节点位移值、截面峰值力、总能量等输出的时间间隔；*DATABASE_HISTORY_NODE控制卡片用来输出节点位移信息，本发明选择B柱结构的全部节点作为位移输出点。

为了方便后续对节点位移值、节点刚度值以及各个单元的相对密度进行读取，本发明利用Tool面板下的renumber功能将模型初始设计空间中B柱结构的全部节点ID和单元ID从1开始进行重新排列。

2.2)B柱初始刚度求解模型

刚度求解模型即B柱结构的静态有限元分析模型，用于B柱结构节点刚度的求解。

本发明利用Hypermesh软件中的Convert功能对上一步建立的初始B柱碰撞分析模型进行转化得到Optistruct模块下的B柱模型，然后在模型中添加刚度求解工况。为了在B柱刚度求解模型的计算结果中生成B柱的静态刚度信息，对B柱结构的所有节点都施加了ASET约束，通过添加PARAM控制卡片，勾选EXTOUT选项并将类型设置为DMIGPCH以输出B柱结构所有节点的静态刚度值。

3.B柱等效静态载荷计算模型的批处理

本发明在上一步建立的B柱初始碰撞分析模型和B柱初始刚度求解模型中，并不包括组件信息、材料信息、属性信息以及单元与组件的对应信息，需要利用MATLAB平台直接在二者对应的k文件和fem文件中自动写入。k文件和fem文件本质上都属于ASCII编码文件，具有非常高的可读性和可操作性，二者的区别主要在于文件格式的宽度不同，k文件信息保存10个字符，而fem文件信息保存8个字符。以B柱初始刚度求解模型为例，模型中待输入的信息如图3所示。

本发明针对碰撞分析模型和刚度求解模型初始设计空间中的每个单元都创建一个组件与之对应，因此创建的组件数目与单元数目相同，即如果待优化的B柱结构含有H个单元，则创建H个组件，同时创建H个材料及属性与单元相对应。一般来说，在Hypermesh软件中，单元必须存在于组件中，即单元信息中必须包含每个单元所对应的组件ID。

综上所述，为了方便在多材料拓扑优化后根据单元ID更新设计变量，本发明将初始碰撞分析模型和初始刚度求解模型中每个单元的ID与其所属的组件ID、材料ID和属性ID进行统一，即第h个单元对应第h个组件、第h个材料以及第h个属性。本发明按照k文件和fem文件的格式标准令ID从1开始创建组件、材料以及属性，创建的数量与单元数量相同，并将单元ID与组件ID进行对应。作为对照，以批处理后的刚度求解模型为例，模型中输入后的信息如图4(a)～图4(d)所示。对比图3和图4(a)～图4(d)可以看出，组件信息、材料信息、属性信息以及单元与组件的对应信息输入到了初始刚度求解模型中的对应位置。

4.基于等效静态载荷法的B柱结构多材料拓扑优化模型的建立

基于等效静态载荷法的B柱结构多材料拓扑优化模型在前处理软件Hypermesh中的Optistruct模块下建立。首先设置好优化三要素，即设计变量、设计约束和优化目标，该模型的设计变量为B柱结构中每个单元的拓扑变量和材料选择变量，设计约束为质量分数小于给定值M

式中：C

n个候选材料的B柱结构基于等效静态载荷法的多材料拓扑优化数学表达式如下：

式中，M

5.等效静态载荷的计算及输入

等效静态载荷法的基本原理是通过动态非线性分析将动态工况划分成若干个时间步，时间步的数量与工况的数量一致均为J,且每个时间步间隔与步骤4种工况设置一致，计算出每个时间步下结构的动态位移值，然后将动态位移值与结构的静态刚度值相乘即可得到若干个时间步下的等效静态载荷，将等效静态载荷作为静态优化模型中的若干个工况进行加载，这样即可使静态优化模型在每个时间步产生的位移值与动态非线性分析中各个时间步的位移值相等，以此实现动态载荷的处理。

5.1)等效静态载荷的计算

在有限元方法中，结构的动态非线性平衡方程为：

式中：M(b)为质量矩阵；C(b)为阻尼矩阵；以第j个时间步为例，K

通过式(2)的动态非线性分析计算得到z

其中，K(b)为结构的刚度矩阵通过线性静态分析也就是对刚度求解模型提交计算后得到。

本发明利用MATLAB平台中的system函数自动将批处理后的B柱碰撞分析模型和刚度求解模型分别提交到LS-Dyna和Optistruct求解器进行计算。在碰撞分析模型和刚度求解模型计算完成后，会对应的生成B柱结构动态位移文件和B柱结构静态刚度文件，文件中分别保存着B柱结构中全部节点的位移信息和全部节点的刚度信息，如图5、6所示。本发明利用MATLAB平台对各个节点的动态位移值和静态刚度值进行提取，得到位移矩阵和刚度矩阵，如图7(a)、7(b)所示。在刚度矩阵和位移矩阵提取完成后，根据公式(3)将二者相乘即可得到等效静态载荷矩阵，如图8所示，等效静态载荷包括三个方向的力和三个方向的力矩两种形式。

5.2)等效静态载荷的输入

本发明建立的基于等效静态载荷法的B柱结构多材料拓扑优化模型中每个工况只包含节点位移约束信息，不包括载荷力信息和载荷力矩信息，模型中的待输入信息如图9所示。

本发明利用MATLAB平台中的fprintf函数自动将上一步计算得到的力和力矩两种形式的等效静态载荷输入到初始多材料拓扑优化模型中(就是公式(5))。等效静态载荷输入后的B柱结构多材料拓扑优化模型中的载荷信息如图10所示。对比图9和图10可以看出，等效静态载荷作为边界载荷以力和力矩两种形式输入到了B柱结构多材料拓扑优化模型中。

6.模型的收敛判定与设计变量的更新

本发明利用MATLAB平台中的system函数自动将等效静态载荷输入后的多材料拓扑优化模型提交到Optistruct求解器进行计算。计算后对碰撞分析模型和刚度求解模型中的设计变量进行更新，依次循环，直至满足收敛条件，优化结束。

6.1)收敛判定

本发明的收敛判定满足两个条件中任意一个时，认为优化结果收敛，优化结束。第一个条件是，当循环次数达到最大循环次数时优化结束；第二个条件是，通过计算相邻两次循环最后一个迭代步的单元相对密度变化量Δ

式中：

6.2)设计变量更新

在循环过程中需要对碰撞分析模型和刚度求解模型中的设计变量进行更新。多材料拓扑优化模型计算完成后会生成SH文件，文件中保存着优化后材料选择变量(m

式中：

实施例：以下给出一个算例验证本方法的有效性

基于上述的设计方法，本发明在MATLAB中编写了一个基于等效静态载荷法的B柱多材料拓扑优化方法的程序代码。用本发明对侧面碰撞工况下的某乘用车B柱结构进行求解，得到该乘用车的多材料B柱结构形式。

首先在Hypermesh前处理软件中的LS-Dyna模块下建立B柱结构的碰撞分析模型，参照该乘用车原始B柱结构的尺寸和形状建立B柱结构的初始设计空间，该B柱有限元模型的尺寸：长×宽×高为1000mm×200mm×100mm，采用四边形壳单元划分有限元网格，壳单元厚度为2mm，网格尺寸为12mm。该B柱模型(不包括刚性体)共含有4474个节点，4472个单元。建立的B柱碰撞分析模型如图12所示，约束B柱两端节点的全部(六个)自由度，用一个质量为1000kg的刚性体以50km/h的速度沿Y轴正方向撞击B柱。仿真终止时间设置为20ms，动画输出时间间隔和节点位移输出时间间隔都设置为2ms，除去初始0时刻对应的时间步，该模型总共有10个时间步。

然后通过Convert功能将B柱碰撞分析模型转换为Optistruct模块下的刚度求解模型。只保留B柱结构的有限元模型，删除刚性体有限元网格及自由度约束等信息，建立generic工况并对全部节点施加ASET约束，添加控制卡片DMIGPCH对B柱结构静态刚度进行输出。建立的B柱刚度求解模型如图13所示。

然后建立B柱多材料拓扑优化模型，设置好多材料拓扑优化三要素并创建10个linear static工况。该B柱结构的设计变量为全部单元的相对密度，设计约束为mf<0.6，优化目标为十个工况下的加权柔度最小(即刚度最大)，候选材料包括钢、铝两种材料，并将最小成员尺寸约束设置为120mm。

将最大循环次数设置为10，运行MATLAB程序，经过8次循环后优化结束，程序运行时间2610s。不同视角下基于等效静态载荷法的B柱多材料拓扑优化结果如图14所示。

拓扑结果显示，该B柱碰撞侧大部分都是材料钢，仅在上面一小部分应用材料铝，而非碰撞侧恰好相反，仅在下面一小部分应用材料钢，大部分都是材料铝，体现了在保证整车侧面碰撞安全性能的同时进行轻量化的思想。

根据上述拓扑优化结果对该乘用车的B柱结构进行材料替换，如图15(a)所示，该乘用车的B柱结构由内板、外板和加强板构成，其中内板分为上下两部分，通过焊点连接在一起，整个B柱结构选用钢这一种材料。参照图14的拓扑优化结果，如图15(b)所示，将外板也分为上下两部分，通过共节点的方式连接在一起，并将内外板的上部分结构材料替换成铝，材料替换后的B柱质量由10.48kg减轻到9.08kg。

本发明针对动态工况下的B柱结构提出了基于等效静态载荷法的多材料拓扑优化方法。基于MATLAB平台将多材料拓扑优化方法和等效静态载荷法相结合，进一步挖掘B柱结构的轻量化潜力。本发明所给算例表明，该方法能够得到动态工况下的多材料B柱结构形式，其轻量化效果相对于单材料结构有着显著提高。

以上所述仅为本发明的优选实例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。反对本发明所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。