国产有限元软件iSolver案例分享 | 模态及响应谱分析

1. 引言

iSolver为一个完全自主的面向工程应用的通用结构有限元软件，对标Nastran、Ansys、Abaqus设计和实现，具备结构有限元常用分析类型和单元、材料、载荷等基础算法组件，精度和商软误差<0.1%。以指示牌模态及响应谱分析为例，演示iSolver的分析流程，并将iSolver和Abaqus及Ansys workbench计算结果进行对比。

2. 有限元模型

指示牌尺寸为1900mm×200mm×4000mm。有限元模型如下图所示：

 

采用壳单元进行模拟，圆管和牌面分别采用钢和铝两种材料，厚度为：10mm钢、6mm钢、4mm铝。材料及界面设置如下图所示：

 

边界条件为指示牌底部设置固定约束，如下图所示：

 

加速度谱设置如下：

 

3. 分析难点

在地震等激励已知的情况下，为了寻找结构中一个或多个位置的最大位移、应力或加速度，不管响应历程中最大值出现的时刻，精确的方法是采用瞬态动力学。

不过地震载荷等随时间一般无规律变化，瞬态分析为了捕捉载荷，时间步长必须取得很小，因而费时且昂贵，所以很多时候采用响应谱来近似得到，这种方法虽然不够精确，但分析速度快，而且由于分析过程简单快捷。

对同一类对象，每个人用这种分析方法得到的结果差异性不会像瞬态分析那么大，也就是更容易得到统一的仿真结果。

以前也有行业外的完全不了解仿真的人问我们：工程复杂问题仿真分析都是理想化模型，仿真本身和试验有很大误差，那么仿真的意义是什么？

按我们的理解只要大家都这么校核结果，慢慢改进校核方法，很容易形成行业的仿真流程和规范，只要按这种规范来做仿真模型就可以认为合格，当通过这种仿真校核规范的产品在实际过程中依然出事故时，那么就改规范校核标准来逐步规避这种事故。

仿真和规范迭代发展，这也是工程上仿真的意义，而不是单纯的追求和实际物理现象的完全一致。

响应谱分析在iSolver内部分三步执行：

1. 首先采用模态分析Frequency步计算每个结构的模态i的振型{      }和模态各个方向的参与因子      , {      }代表模态      的最大响应。
2. 将实测的地震等随时间无规律变化载荷曲线通用FFT改为频域曲线      。
3. 在上述曲线中找到单个模态频率fi对应的载荷幅值      ，      计算每一个模态的最大位移矢量       。
4. 对于给定激励的最大响应，最大值通过估算结构本身最大响应的方式组合而得，譬如采用ABS方法，那么采用各阶模态效应绝对值相加的方法得到最大响应，但这种结果过于保守，因为一般情况，结构在外界激励作用下各阶模态响应峰值不可能同步出现。因此还可选择NRL等其它更加合理的模态效应组合方式，以准确估计结构在外界激励作用下的总体响应峰值。至于具体是哪种模态效应，每个行业不同，这就是规范，有些行业的规范不止是简单的模态组合，还有对于譬如多个模态频率相近的情况的处理等等。

从上面响应谱分析的流程来说，算法极其简单，真正的难点主要是软件实现时的响应谱分析步流程如何融入到整体的分析步流程中。iSolver整体框架如下图所示：包括静力分析、模态分析、稳态、瞬态、声学、非线性等六大分析类型，同时每个大的分析类型下又各自有细节的小的分析类型，每个小的分析类型都是通过增量迭代步的统一流程调用基础的数学、单元、全局的工具库函数，只是做的操作不同。

一般有两种方式来实现这种根据不同分析类型实现不同具体操作的流程：

1. 只有一个整体的实际分析流程，每次要做一个具体操作时，用if else先判断是哪个分析类型，来判断具体要做哪个操作，一个简单的例子就是if 模态分析/稳态分析/响应谱分析/基于模态的瞬态分析时，那么组装质量阵，else不组装质量阵。这种类型相对搞理论的人来说非常直观，需要每个编程人员对所有分析的流程都要熟悉。但扩展性相对较差，加一个譬如本文中的响应谱分析需要从分析流程的头到尾加一堆的if判断。很多网上的开源代码都是这种整体分析流程结构，如果是Fortran写的上古代码，你还会看到很多的Goto语句，这些很难维护。
2. 有很多个分析流程，针对每个分析类型都有一个不同的分析流程，同时各个分析流程中的如果能公用的函数都直接复用，这就已经不仅仅是专业理论和普通的编程能力了，而是要深入了解面向对象的编程技术了。缺点是专业和软件两个方面要求都比较高的编程人员难找，优势是容易扩展，如果架构足够好，当有一个新的分析类型时可以不影响以前的流程。

iSolver的Step类为一个纯虚类，仅仅涉及了任何一个有限元分析流程的共性部分，各种不同的分析的个性化部分则交给从Step类继承的派生类之相应的方法去处理，在主函数中遍历Step，对每个Step采用工厂模式根据用户选择的分析问题类型名称譬如“StaticLinearPerturbationStep”而生成相应的从Step类继承的派生类的实例，并返回相应的指针。

在每个Step派生类主要完成三方面的工作：

1. Init（）函数，分析步前的数据初始化，譬如对模态分析一般节点坐标直接是从求解输入文件中直接读取就行，但如果前面分析步是静力分析，那么将静力分析的最终位移场量叠加到初始位置作为各个节点坐标。
2. Solve（）函数，具体每个分析步对应的个性化的过程，如组装刚度矩阵，迭代控制，调用矩阵求解，求解分析的核心算法一般在此。
3. Update（）函数，每个分析步结束时对输出结果或者中间状态物理量的更新，譬如静力分析输出应力、应变、位移（包括平动和转动）等结果，塑性应变则作为状态量输出等。而模态分析则不需要输出塑性应变，但需要输出各个模态的参与因子，这样下一步的响应谱能直接调用。

4. 结果对比

将相同模型分别通过iSolver、Abaqus和Ansys Workbench进行模态和响应谱分析，并将三者的计算结果进行对比。

模态频率值如下表所示：

一阶模态位移云图如下：

 

二阶模态位移云图如下：

 

三阶模态位移云图如下：

 

四阶模态位移云图如下：

 

五阶模态位移云图如下：

 

六阶模态位移云图如下：

 

iSolver谱分析Mises应力云图如下：

Abaqus谱分析Mises应力云图如下：

Ansys谱分析Mises应力云图如下：

结论：基于iSolver软件对指示牌进行了模态及响应谱分析，并将计算结果与Abaqus和Ansys进行了对比。三者模态振型均变现出一致性，iSolver与Abaqus模态频率值一致，略大于Ansys计算结果；iSolver与Abaqus的响应谱分析的最大应力值一致，均为1.62MPa，略小于Ansys计算结果。