高性能计算应用21

21
高性能计算应用
大规模柴油机动力学抗冲击并行仿真计算
●丁峻宏 王 惠 刘 波
上海超级计算中心 上海 ２０１２０３ ｊｈｄｉｎｇ＠ｓｓｃ．ｎｅｔ．ｃｎ
●宋雅丽 任建军
沪东重机有限公司 上海 ２００１２９
摘要：
针对柴油机动力冲击仿真问题，采用显式动力算法，对柴油机复杂结构进行了详细建模和
精细仿真，并对处于不同工作状态下的柴油机开展了机构运动和动力冲击问题的单独和综合模
拟。

为满足计算调试和多组工况带来的大规模计算需求，问题求解借助了高性能计算资源；通
过不同计算平台上多核并行计算和性能数据比较，分析了独特并行加速现象，并探求了合理的
并行加速策略。

计算表明，动态仿真方法有利于柴油机运动和动力问题的一体化研究，而借助
高性能计算资源方能为问题求解和深入研究提供必要保障，文中结论可对柴油机类似复杂冲击
问题研究提供一定参考。

关键词：柴油机，显式算法，机构运动，动力冲击，并行计算
柴油机是目前船舶上使用最广的动力装置，在船只受到不可预测的外部冲击载荷作用下，作为主动力源的柴油机抗冲击性能将直接影响到船舶的续航，因此对其开展研究具有重要意义。

包括柴油机在内的船舶设备结构力学性能相对复杂，研究难度大，此外由于一些应用比较特殊，因而关于其冲击响应特性研究的国外文献资料并不太多。

虽然国内这方面工作起步较晚，但近年来相关研究逐渐涌现，随着计算数学和计算力学的发展，基于仿真计算并结合少量实验的方法对船用设备开展抗冲击研究正不断走向深入。

船舶设备冲击响应计算分析方法包括模态叠加法、直接积分法、冲击因子法、有限元法以及多体动力学方法等［１］。

随着计算机软硬件技术的飞速发展以及对船舶设备研究的深入，结合有限元法和多体动力学法特点的刚柔混合模拟实际运用越来越多并逐渐成为主流。

比如从有限元模型子结构缩减出发，转入多体动力学软件中分析求得力学边界条件，最后恢复有限元精细模型计算动应力［２，３］。

再比如借助有限元软件进行模态分析得到柔性体中性输出文件，转入多体动力学软件计算系统动力学特性［４－６］。

计算过程涉及跨软件操作和中间数据传递，需要由专业研究人员才能完成。

对于柴油机等复杂机械设备来说，无论是有限元法还是多体动力学法，最关键的是对各零部件之间的联接处理，比如运动联接、螺栓联接，焊铆联接等。

多体动力学法多用于模拟分析复杂机械系统运动规律和部件受力环境，但是刚体定义过多必然失真，而定义过多柔体将增加模型和处理过程的复杂性。

有限元法则可以更好地模拟船用设备柔体冲击力学特性，网格细化和细节模拟也将有助于深化研究，但大量零部件机构运动模拟以及计算网格规模扩大对普通计算机的能力提出了挑战。

更重要的是，对于处于运行状态下的柴油机承受冲击问题来说，仿真模型应该要同时考虑系统机构运动和动力冲击两方面特性的影响。

在针对大型复杂柴油机结构开展研究分析时，如若在同一模型中能同时直接模拟这两个问题［７，８］，减少中间数据传递，则有利于提高工业企业仿真计算的一体化、便捷性甚至准确性。

本文针对某型号１６缸柴油机，利用有限元前处理软件Ｈｙｐｅｒｍｅｓｈ以及ＡＢＡＱＵＳ，对其所包含大量零部件进行了合理建模和一体化组装，建立了大规模柴油机动力学冲击仿真模型，计算则基于超级计算
注：　原文发表于《振动与冲击》２０１４年第３３卷第２期第１６３－１６７页，　此处略有改动．
22
《高性能计算发展与应用》 ２０１４年第一期 总第四十六期
机、显式动力算法和区域分解并行求解方法，研究了冲击载荷作用下柴油机中重要零部件受到的冲击响应和安全性能，并就如何有效使用高性能计算资源来解决此类大规模仿真分析进行了探讨。

１．　原理与方法
显式算法非常适合解决柴油机动态响应求解时所遇到的机构运动大位移以及多零部件瞬态冲击问题，它由前一增量步向前推算动态行为，通过节点逐个求解而不必迭代，通过小时间增量来获得高精度求解。

不仅能有效避免隐式算法求解复杂动力学问题时由于反复迭代造成的收敛困难、以及硬盘空间存储问题，而且此算法具有很好的并行加速比，有助于利用硬件环境的改善来提高所求解问题的计算规模和复杂程度。

１．１　显式算法基本方程
三维有限元体系的运动方程可表达为：
（１）
其中，
，
和
分别为体系节点的加速
度、速度和位移矢量，［Ｍ］为总体质量矩阵，［Ｋ］为总体刚度矩阵，［Ｃ］为总体阻尼矩阵，｛Ｐ｝为总体载荷矢量。

显式算法基于ｔｎ时刻及之前的信息来求解ｔｎ＋１时刻的位移
ｎ＋１
，结合位移Ｔａｙｌｏｒ表达式和中心差分法，
可分别得到速度和加速度的表达式为：
（２）
（３）
将公式（２）和（３）代入ｔｎ时刻的公式（１），可以得到下一时刻
ｎ＋１
的求解表达式：
（４）
在整个时域范围内，可由上述积分递推公式求得各个离散时间点处的位移、速度和加速度。

１．２　多ＣＰＵ／核并行求解
ＡＢＡＱＵＳ显式并行算法可分为Ｄｏｍａｉｎ和Ｌｏｏｐ两种不同方式，通常来说基于ＭＰＩ的Ｄｏｍａｉｎ并行求解方式效率更高。

如图１所示，采用并行处理模式进行大规模问题求解时，采用ＭＥＴＩＳ区域分解方法（Ｄｏｍａｉｎ　Ｄｅｃｏｍｐｏｓｉｔｉｏｎ　Ｍｅｔｈｏｄ，简称ＤＤＭ）来处理分析对象，即把整个模拟问题区域划分成很多相对小的求解子区域，　把每个子区域分配给不同的处理器核心分别进行求解，区域分解允许每个处理器核心独立于其他
处理单元来求解自己的那部分问题，而处理器及其核心相互之间通过交互机制进行数据交换，最后将每个子区域的解综合起来，得到整个模型区域的全局解。

图１ 并行求解的模型分区和ＣＰＵ布局
２．　应用实例
２．１　柴油机装配体有限元建模
真实的柴油机及其附属结构非常复杂，包含上千个零部件，需要进行合理和必要的建模和简化，本文主要针对柴油机内部运动件和主要固定件开展分析。

借助前处理软件，首先建立了包括曲轴、平衡块、连杆、活塞和机架等零部件在内的柴油机全三维详细有限元模型，在对柴油机各零部件进行网格划分时，需要对几何模型进行必要的结构特征简化，如局部过小倒角等。

但是重要零部件如整体式曲轴保留了内部油孔油道等细部特征，以及曲柄臂与连杆轴（曲柄销）和主轴颈交接处的过渡圆角。

然而，各零件非重要几何特征细节的保留不便于六面体单元的划分，因此本文通过局部网格加密（如曲轴油孔油道和轴颈过渡区域）、增加四面体单元数量（如仅曲轴就划分了２８万单元）并保证整体网格质量（如严格执行翘曲度、长宽比和扭曲度等检查）等手段以求进一步提高模型计算精度。

柴油机的零部件基本上可以分为固定件和运动件两大类，在本文模型中，机架、主轴承座和气缸套等为固定件，而曲轴组件、连杆组件和飞轮等为运动件，需要通过多种连接单元和仿真手段模拟了各零部件之间的机构连接关系和装配关系。

对于本文所建柴油机计算模型而言，考虑到的零部件之间的主要装配关系包括：１）曲轴和连杆之间、曲轴与主轴承座（主轴瓦）、连杆和活塞之间均定义为圆柱副；２）曲轴与飞轮及齿轮之间定义为固定副；３）活塞和气缸套之间定义为平动副；４）气缸套与机架之
23
高性能计算应用
间、机架与主轴承座均设有绑定约束。

包括四面体和六面体在内的整个柴油机网格模型单元总数超过２１６万，节点总数近６６万。

图２显示
的是隐蔽部分机架网格后的柴油机整体计算模型。

图２ 柴油机整体有限元模型
２．２　相关边界、载荷、参数和工况设定
柴油机工作时，曲柄连杆机构主要外力为最高燃烧压力和惯性力，并对外输出扭矩。

本项目中针对柴油机整体计算模型力学边界条件，以及各种载荷和参数因素的设定，主要考虑了以下几个方面：
（１）选取机架底部对于螺栓连接的部分区域的节点集合，在其上施加固定边界；同时，沿冲击方向，施加如图３所示的历时５０ｍｓ峰值为３９ｇ的短程加速度冲击曲线。

（２）在柴油机运行状态下，针对曲轴等运动部件施加工作转速１１０ｒａｄ／ｓ；将计算得到的输出扭矩作为阻力矩施加在曲轴输出端，以使曲轴在整体计算中保持转速平稳。

（３）柴油机一个完整发火周期为０．１１４ｓ。

根据不同气缸发火顺序和发火时间间隔，在活塞盖的顶端施加各自对应的压力时程曲线（见图４）。

（４）曲轴连杆等运动件采用４２ＣｒＭｏＡ，屈服强度９３０ＭＰａ；机架采用铸铁，屈服强度为２５０ＭＰａ。

（５）　当柴油机静止受冲击时，直接施加冲击加速度时程曲线，计算历时５０ｍｓ；当柴油机正常运行受冲击时，设计工况使得冲击波第一个加速度峰值恰
好与Ａ４缸发火瞬间重合，总计算时间历时７５ｍｓ。

图３ 冲击加速度时程曲线
图４ 气缸内燃烧压力时程曲线
３．　计算和分析
３．１　柴油机抗冲击性能评估
针对此计算模型进行了大量冲击模拟，限于篇幅，本文在此仅以柴油机受垂向冲击为例，介绍重要零部件处于停机和运行两种不同状态下的结构响应分析。

为作相互比较，也另外进行了运行但不受冲击时的柴油机机构运动仿真。

根据相关行业规范，要求设备上由冲击载荷引起的应力不得超过静态屈服极限。

对此，可首先得到各重要零部件各时刻的应力云图，再遍历寻求零部件在受冲击时间历程中最大应力数值，见表１。

通过进一步分析比较发现：
（１）柴油机静止受垂向冲击时，曲轴紧邻飞轮的１＃连杆轴颈上将产生最大７８ＭＰａ的等效应力，机架底部固支表面靠飞轮端边缘也有超过８０ＭＰａ的数值出现；主轴承座和活塞上有近７０ＭＰａ的冲击应力产生，而连杆和活塞销上产生的冲击应力相对较小。

（２）柴油机运行但不受冲击时，曲轴上最低应力数值也在１００ＭＰａ以上，最高近３００ＭＰａ；各连杆上最大应力分布在２５０￣３３０ＭＰａ之间，Ａ４连杆上有超过３００ＭＰａ的数值（见图５）；各活塞和活塞销最大应力分布在１００￣２００ＭＰａ之间。

固定件如主轴承座最大值出现在紧邻飞轮端的１＃主轴承座上。

表１ 垂向冲击时不同工况下零部件最大应力对比
计算工况部件名称
静止受冲击（ＭＰａ）运行无冲击（ＭＰａ）运行受冲击（ＭＰａ）曲轴７８．２３２９４．６３２２．６连杆４４．２８３２３．１３７８．３活塞销２９．９３２０２．４２２２．１活塞６７．３０１６５．０１９８．２　主轴承座６６．５４１６６．６１８０．３机架
８１．７９
８０．１９
９７．６０
（３）柴油机运行受垂向冲击时，各零部件最大应
力数值均有数十ＭＰａ的增加，曲轴应力极值出现位置（见图６）和“运行无冲击”工况相比没有变化，但
24
《高性能计算发展与应用》 ２０１４年第一期 总第四十六期
是连杆和活塞等则出现在不同位置和不同时刻。

主轴承座和机架上最大应力出现位置基本没有变化，
应力增幅也不大。

图５ 连杆上最大等效应力出现时刻云图
图６ 曲轴上最大等效应力出现时刻云图
（４）各零部件受冲击后产生的最大冲击应力数值均小于各自屈服极限，柴油机总体上能够承受峰值峰值为３９ｇ的加速度垂向冲击。

采用４２ＣｒＭｏＡ材料的曲轴连杆等运动件安全裕度比较大，但是针对采用铸铁材料的机架开展进一步研究发现，如果考虑气缸内燃烧压力沿气缸盖内侧和连接螺栓传递到机架的影响，那么需要注意机架螺栓孔中连接区域的局部应力。

（５）综合比较发现，“柴油机运行受冲击”工况并非“静止受冲击”和“运行无冲击”两工况的简单叠加，系统冲击响应具有非线性特性。

以曲轴为例（见图７），当柴油机静止时，曲轴最大应力出现在冲击加速度到达第一个峰值时刻附近，而当柴油机运行时，不仅最大应力出现位置改变，而且应力变化趋势也主要受机构运动特性影响，即便受冲击
应力增幅也不太大。

（１）　静止状态
（２）　运行状态
图７　不同状态下曲轴最大等效应力时变曲线
３．２　高性能计算平台并行求解
本文建立的柴油机抗冲击计算模型规模大，历时较长，如果使用普通工作站（１￣４个ＣＰＵ核心；内
存８￣１６ＧＢ），计算一个工况至少需要５￣１０天，且需预留硬盘空间５０ＧＢ以上。

如果计算工况多，则对计算时间和计算效率提出了更高要求，必须要借助高性能计算资源。

上海超级计算中心目前有两套计算平台，一是运行四年多的“魔方”超级计算机，其中用于工程计算的Ａ区由８２个计算节点构成，每个节点含八个ＡＭＤ四核处理器（主频１．９　ＧＨｚ，共享１２８ＧＢ内存）；二是今年新建“蜂鸟”集群，该集群包括６５台ＨＳ２３刀片计算节点，每个节点含两个Ｉｎｔｅｌ八核Ｉｎｔｅｌ处理器（主频２．６ＧＨｚ，　共享６４Ｇ内存）。

针对柴油机抗冲击计算模型的单个计算工况，在这两个高性能计算平台上分别进行了多核并行计算测试和性能比较，数据汇列于表２和表３中，并据此得到了以下一些认识：
（１）在使用“魔方”计算机进行求解时，随着ＣＰＵ计算核心数的增加，单个工况所需计算时间迅速降低，在使用６４核时，可能将单工况求解时间压缩至近一个工作日。

８核以内多核并行求解时，出现超线性加速比现象；１６核以上多核并行计算时，并行效率才回落到１００％以下，并随着计算核心数的比例增加呈明显下降趋势。

对此进行分析认为：１）由于本文建立的柴油机冲击计算模型规模较大，在使用单核计算时，高速缓存（Ｃａｃｈｅ）存入数据量小，因此需要更多调用内存，运行时间因而变长；２）随着核数增加，每个进程分发问题规模随之变小，高速缓存可容纳必要处理数据相对增加，Ｃａｃｈｅ命中率提高，因此总体运行性能大幅提高。

整体上高速缓存效应优势大于通信
25
高性能计算应用
开销导致的性能损失，从而造成超线性加速现象；
３）针对本文问题，使用１６核以上甚至更多核数时，
由于进程间开销增加以及负载均衡性降低的影响逐
渐加大，因此才会出现加速比性能和并行效率数值
回落至正常数据水平，并使并行效率逐渐降低。

（２）在使用“蜂鸟”计算机进行计算时，随着
ＣＰＵ计算核心数的增加，单个工况所需计算时间也
随之不断降低，在使用１６核时，就可以将单工况求
解时间迅速降低至一个工作日，与“魔方”相比计
算性能优势比较明显。

随着计算核心数的进一步增
加，加速比和计算效率呈有规律的正常下降趋势，
无超线性加速现象出现。

对此进行分析认为：１）由于蜂鸟ＣＰＵ性能比魔
方ＣＰＵ不仅体现在主频的提高，缓存尤其是三级缓
存更增加了十倍，单核数据处理能力相对有了很大
提升，因此没有出现超线性加速现象；２）区域分解
属于粗粒度算法，意味着计算开销要远大于通信开
销，但当计算核数增加到一定数目使得通信开销影
响加剧时也必须予以重视；３）虽然这两台不同计算
集群的计算节点之间均采用ＩｎｆｉｎｉＢａｎｄ网高速互联，
但数据表明“蜂鸟”跨节点计算时能够相对有更佳
表现，数十核求解问题时仍可使并行效率保持在高
位。

表２ 魔方计算机上的单工况并行计算性能比较
性能比较计算核数计算时间
（ｈ）
加速比
并行效率
（％）
１１３７９１．００１００２５７８．４２．３８１１９４３２７．８４．２１１０５８１５８．３８．２０１０２１６８６．７３１５．９９９．９３２４６．３８２９．７９２．９６４２７．３１５０．５７９．０
表３ 蜂鸟计算机上的单工况并行计算性能比较
性能比较
计算核数
计算时间
（ｈ）
加速比
并行效率
（％）
１４０１．８１．００１００
２２００．９１．９９９９．５
４１０１．９３．９４９８．５
８５１．７３７．７７９７．１
１６２６．５１１５．２９４．６
３２１３．７９２９．１９１．０
６４７．５６７５３．１８３．０４．　结论
本文详细阐述了在高性能计算平台基础上，如何将显式有限元算法、并行计算以及加速策略有效结合起来，对大型柴油机结构动力冲击问题进行数值仿真，给出了具体实现方法和结果分析，并有以下认识：
（１）以高性能计算资源为支撑，可全力对柴油机动力冲击问题进行细致建模和有限元模拟、同时能真实再现柴油机处于不同工作状态下时机构运动和动力冲击问题的一体化分析，并将此大规模复杂计算仿真模型的问题求解和深入研究付诸实现。

（２）基于区域分解法的并行策略和集群架构高性能计算机，显式算法计算量开销大于通讯量，非常适合柴油机等设备大规模冲击问题的求解。

数据表明，对本文百万级计算模型而言，６４核时加速比和并行效率依然不低，仍表现出较好的可扩展性。

（３）不同体系结构超级计算机的并行计算性能差别较大，工业用户进行同类问题仿真时，可据此进行合理的资源搭建或配置。

此外，在进行超线性加速和并行性能数据差异比较时，要和不同ＣＰＵ具体性能参数以及网络连接技术参数结合起来分析。

参考文献
［１］　贺少华，吴新跃．　舰艇机械设备冲击响应仿真建模计算方法综述［Ｊ］．　船舶科学技术，　２０１１，　３３（１）：　１４－３３．
［２］　Ｐｏｌｏｎｉ　Ｃ，　Ｇｅｒｅｍｉａ　Ｐ，　Ｃｌａｒｉｃｈ　Ａ．　Ｍｕｌｔｉ－ｏｂｊｅｃｔｉｖｅ　ｒｏｂｕｓｔ　ｄｅｓｉｇｎ　ｏｐｔｉｍｉｚａｔｉｏｎ　ｏｆ　ａｎ　ｅｎｇｉｎｅ　ｃｒａｎｋｓｈａｆｔ［Ｃ］．　ＩＩＩ　Ｅｕｒｏｐｅａｎ　Ｃｏｎｆｅｒｅｎｃｅ　ｏｎ　Ｃｏｍｐｕｔａｔｉｏｎａｌ　Ｍｅｃｈａｎｉｃｓ：　Ｓｏｌｉｄｓ，　Ｓｔｒｕｃｔｕｒｅｓ　ａｎｄ　Ｃｏｕｐｌｅｄ　Ｐｒｏｂｌｅｍｓ　ｉｎ　Ｅｎｇｉｎｅｅｒｉｎｇ，　Ｌｉｓｂｏｎ，　Ｊｕｎｅ　５－８，　２００６．
［３］　陈亮，　宋希庚，　明章杰等．　基于有限元和多体动力学的柴油机曲轴动态强度与冲击响应分析［Ｊ］．　振动与冲击，　２００８，　２７（１１）：　１８６－１８９．
［４］　Ｙｉｌｍａｚ　Ｙ，　Ａｎｌａｓ　Ｇ．　Ａｎ　ｉｎｖｅｓｔｉｇａｔｉｏｎ　ｏｆ　ｔｈｅ　ｅｆｆｅｃｔ　ｏｆ　ｃｏｕｎｔｅｒｗｅｉｇｈｔ　ｃｏｎｆｉｇｕｒａｔｉｏｎ　ｏｎ　ｍａｉｎ　ｂｅａｒｉｎｇ　ｌｏａｄ　ａｎｄ　ｃｒａｎｋｓｈａｆｔ　ｂｅｎｄｉｎｇ　ｓｔｒｅｓｓ　［Ｊ］．　Ａｄｖａｎｃｅｓ　ｉｎ　Ｅｎｇｉｎｅｅｒｉｎｇ　Ｓｏｆｔｗａｒｅ，　２００９，　４０：　９５－１０４．
［５］　郭磊，　郝志勇，　林琼．　柴油机曲轴与气缸体系统动力学仿真研究［Ｊ］．浙江大学学报（工学版），２００７，４１（５）：　７８０－７８４．
［６］　吴问鲍．　柴油机抗冲击性能分析［Ｊ］．　船海工程，　２０１０，　３９（１）：　６３－６５．
［７］　汪玉，　赵建华，　杜俭业等．　基于多体动力学有限元计算的一种舰用柴油机抗冲击性能仿真［Ｊ］．　振动与冲击，　２００９，　２８（１１）：　８７－９０．
［８］　王承禹，　常汉宝，　郑仁伯等．　基于ＤＹＴＲＡＮ的发动机曲轴系冲击动力学问题［Ｊ］．　内燃机工程，２００９，３０（４）：　８７－９２．。