“旋转的硬币”——不同求解器对比

前言    

随着各大CAE公司逐个被收购，才认识到原本的CAE巨头在更大的巨头面前也只是小鱼，感叹资本的力量。

但目前看来，貌似这些收购暂时未影响CAE行业的格局。

在显式求解器领域，主要仍然是 Ls-dyna，Pam-Crash，abaqus，Radioss(以及OpenRadioss)，这4大显式求解器。

本文将通过使用这4个不同的求解器，来模拟同一现象——“旋转的硬币”及其倒下的过程。

在保证每个模型设置一致的前提条件下。判断能否得到相同的计算结果。包括：硬币质量，边界条件，旋转中心，摩擦系数等影响因素。

模型前处理  

其中，硬币采用1元硬币，模型如下图所示，网格模型采用与之前《抛硬币》的文章中所用模型一致，直径30mm，厚度1.9mm，质量为6.1g。为了区分正反面，在硬币正面设置数字标识。

硬币模型

硬币质量信息

下面，本文将使用上述提到的4个不同的求解器进行模拟，由于各求解器设置方法略有不同，小编将会在不同的前处理软件中实现前处理的设置。最终，会确保各模型的旋转速度、旋转中心、摩擦系数等参数一致。

1  Abaqus求解器    

Abaqus的inp求解文件可以在Abaqus CAE软件，Hypermesh或者Ansa中创建，第三方软件与Abaqus CAE相比各有优劣，由于是导入其他求解器的网格进行转换，本文使用了Hypermesh进行设置。

1-1 边界条件设置   

Abaqus中没有刚性墙，所以需要建立一个壳单元的part作为地面，将其抓取为刚体，并约束其主节点。如下图所示：

abaqus中地面及其约束设置

1-2 载荷设置   

对硬币施加初始旋转速度，设置如下图所示：

abaqus中硬币旋转初速度设置

重力载荷施加如下图所示：

abaqus中重力载荷设置

1-3 接触设置   

Abaqus中，硬币与地面接触采用了通用接触的设置方法，参与接触的对象选择模型中所有的部件；通过接触属性定义摩擦系数，摩擦系数定义为0.2。接触设置以及相关卡片如下图所示。

abaqus中硬币与地面接触设置

2 Ls-dyna求解器    

由于Ls-dyna求解文件卡片的格式固定，逻辑一致性较好，因此主流前处理软件对Ls-dyna的支持都相对较好，且相对优于LS-PrePost。本文的求解K文件在Hypermesh前处理软件新界面中进行设置完成。新界面Hypermesh软件与经典界面相比，GUI设计更加现代化，增加了更多的快捷键，释放了左手，老CAEr们趁早上手吧，毕竟老的要淘汰 。

2-1 边界条件设置   

Ls-dyna中可以通过刚性墙关键字*RIGIDWALL定义地面。在关键字卡片中设置刚性墙的位置，和从属关系以及摩擦系数。当然，如果想使用Abaqus中的设置方法，用壳单元作为刚体的方式也可以，Ls-dyna中刚体的定义方式可以通过材料实现。

刚性墙的设置如下图所示：

 

Ls-dyna中刚性墙设置

2-2 载荷设置   

在Ls-dyna中，使用*INITIAL_VELOCITY_GENERATION关键字定义旋转初速度，关于卡片的解释如下：

-OMEGA - 关于转动轴的角速度；

-VX、VY、VZ - 指定X、Y、Z方向的初始平移速度；

-XC、YC、ZC - 全局坐标系下旋转中心点的X、Y、Z坐标；

-NX、NY、NZ - 转动轴的X、Y、Z方向余弦；当NX=-999时，NY、NZ定义两个点，决定旋转轴的方向。

定义如下图所示：

Ls-dyna中硬币旋转初速度设置

在Ls-dyna中，定义重力载荷，使用关键字*LOAD_BODY，需要定义曲线和缩放值。需要注意的是，重力方向为定义方向的反方向。如下图所示：

Ls-dyna中重力载荷设置卡片

3 PamCrash求解器    

由于PamCrash求解文件卡片的格式并不固定，而且随着求解器版本的更新，新增关键字的同时，有些关键字会被淘汰，甚至某些关键字的默认选项的意义会发生变化。因此，相对于第三方前处理软件，PamCrash原生的前处理软件对其关键字的支持是最好的，即ESI Group公司的VCP软件。本文将使用VCP软件对PamCrash的求解文件进行设置。

3-1边界条件设置   

PamCrash中刚性墙关键字已经被淘汰，因此，设置刚性地面，需要设置单独的Shell Part，将其抓取为刚体，并约束刚体中心。刚体抓取可以使用传统Rbody关键字创建刚体，也可以使用MTOCO关键字，本文使用Rbody创建刚体。需要注意的是，由于PamCrash模型的接触有接触厚度的概念，在定义地面位置时，需要将该距离预留出来。如下图所示：

PamCrash中刚性墙设置方法

PamCrash中约束设置

3-2 载荷设置   

在PamCrash中定义初速度使用关键字INVEL▫/▫，在卡片中需要定义速度值和局部坐标系，两者共同决定最终硬币旋转的绕轴和速度。确保与前文旋转中心一致的关键在于局部坐标系的定义，局部坐标系的原点决定其轴所在的位置。具体设置如下图所示。

PamCrash中旋转初速度的设置

在PamCrash中设置重力载荷，使用关键字ACFLD▫/▫，同样由曲线和缩放值，共同组合决定最后重力载荷的大小和持续时间。

PamCrash中重力载荷的设置

3-3 接触设置   

在PamCrash中，硬币与地面接触采用了33号接触，主、从关系分别为地面和硬币；在接触卡片中定义摩擦系数，类型为常数，值为0.2。接触卡片如下图所示。

PamCrash中硬币和地面的接触卡片设置

4 Radioss求解器    

Radioss求解器是Altair旗下的显式动力学求解器，功能强大，2022年9月，Altair发布了其开源版本 OpenRadioss，这一操作在当时引起了广泛的讨论。目前，OpenRadioss社区仍在不断完善其功能。本文在Hypermesh前处理软件中完成Radioss求解文件的设置。

4-1 边界条件设置   

Radioss与Ls-dyna类似，可以通过刚性墙关键字 /RWALL定义地面。在关键字卡片中设置刚性墙的位置，和从属关系以及摩擦系数。如下图所示。

Radioss中刚性墙设置

4-2 载荷设置   

在Radioss中，使用 /INIVEL/AXIS 关键字定义旋转初速度，类似于PamCrash的设置方法，需要定义速度值和局部坐标系，两者共同决定最终硬币旋转的绕轴和速度。如下图所示。

Radioss中旋转初速度的设置

在Radioss中设置重力载荷，使用关键字 /GRAV，与前两者类似，同样由曲线和缩放值，共同组合决定最后重力载荷的大小和持续时间。如下图所示。

Radioss中重力载荷的设置

结果后处理  

为了方便对比，各求解文件有相同的计算时间和一致的输出频率，各求解器的计算结果如下图所示，动画结果表明，4者略有差异，倒下的时间不同，但最后都得到了正面朝上的结果 。

不同求解器的计算结果

视 频    

小编将PamCrash前处理的操作过程制作成了视频，

结束语

本文通过使用不同的求解器，对“旋转的硬币”直至倒下这一过程进行了模拟。并未证明或对比任何精度问题，仅作有趣的探索。发现 4者在前处理设置一致的前提条件下，动画结果基本一致，且均得到了正面朝上的结果。

观察结果动画，可以发现，4个结果还是有些差异，硬币倒下的时刻不相同。大致可以说明同为显式求解器，大方向的基础算法，大同小异，但某些细节之处仍有差异。针对该差异，笔者并未做更深入的探讨研究，仅做抛砖引玉，读者若有兴趣，可以做下一步更深入的研究。

希望本文能为你带来一些收获，感谢阅读！

文字 · TechManLxs

校核 · Nicle

视频 · TechManLxs

奇妙仿真CAE

准确的仿真比比皆是

有趣的仿真万里挑一

长按扫码，关注我们