爆炸与冲击问题Euler型数值方法的并行算法

　　0　引言

　　爆炸与冲击问题是一类涉及高温、高压、高速等极端条件的瞬态动力学问题.在这些极端条件下,爆炸与冲击问题的数值模拟需要对材料的大变形、多种物质的界面以及各种强间断进行处理,比通常的流体力学问题、空气动力学问题及结构动力学问题要复杂得多.

　　现实中的工事错综复杂,工事壁面的约束作用使得冲击波峰值超压将降低得很慢,持续时间很长,与空气介质中的爆炸存在很大的差别.因此对无限空气介质中的爆炸场进行三维数值模拟研究远不能满足现实的要求.对这类复杂模型进行数值模拟显得十分重要,研究爆炸冲击波在工事内的传播过程及其相互作用以及冲击波峰值超压的大小能够有效的判断战斗部对地下工事的毁伤程度.

　　近年来,随着大规模科学与工程计算的需求,许多计算问题已经超出单机所能承受的能力范围,并行计算对于大规模科学与工程计算越来越重要.而三维爆炸与冲击问题的数值模拟的突出问题就是计算规模大及计算时间长.

　　针对这一问题,本文提出了爆炸与冲击问题Euler型数值方法的并行算法,采用国际上流行的MPI标准,设计了三维多物质弹塑性流体动力学并行程序PMMIC-3D程序.MPI(消息传递接口)是消息传递函数库的标准规范,是目前广泛使用的并行编程工具,MPI基于分布式存储,但同样适应于共享存储,具有移植性好、功能强大、效率高等多种优点.

　　通过对两类模型的数值模拟,考察其数值模拟结果是否与实际物理规律相符,验证所采用的并行算法的可行性,并分析数值模拟结果,为以后的进一步的研究工作提供一定的参考作用.

　　1、数学模型及数值方法

　　1.1控制方程

　　爆炸与冲击过程是一个复杂的流体弹塑性过程,计算采用不考虑外力、外源和热传导,非守恒形式的Euler流体弹塑性动力学偏微分方程组:

　　

　　1.2数值计算方法

　　将三个守恒方程按照各项的物理特性,统一写成如下形式:

　　

　　式中,f代表ρ(密度)、e(比内能)和u(速度矢量)等物理量,H为源项.

　　数值计算方法采用算子分裂算法,按物理效应进行分裂,把控制方程中的三组守恒方程按对流项和扩散项分成两步,即Lagrange步和Euler步进行计算.Lagrange步不考虑对流项u·f的影响,仅考虑源项H的作用,即压力和偏应力的梯度效应作用,得到网格各物理量的中间值;Euler步考虑对流项u·f的影响,通过计算网格间的输运量,从而重新确定质量、动量、能量在网格间的分布.

　　1.3界面处理方法

　　PMMIC-3D并行程序界面处理采用模糊界面方法,对于含有多种物质的混合网格,不区分物质间的界面;根据模糊方法计算各物质占有网格的体积比,把体积比作为模糊权重系数;对网格进行分类,不同类网格之间视为物质界面;对介质进行模糊排序,决定输运优先权和制定模糊输运表;根据模糊权重计算输运量并按模糊输运表进行输运;在建模和计算中应用模糊方法,故称为“模糊界面方法”.

　　1.4并行方法

　　PMMIC-3D并行程序采用国际上通用的MPI(消息传递接口)标准,基于域分解的方法,将三维计算区域分成若干子区域,实现了计算域沿着3个方向的任意分割,并创建子区域间的拓扑结构,一个处理器处理一个或多个子区域.每个子区域之间的通信采用捆绑式发送和接收函数,这样可以有效的避免因通信而造成程序的死锁,并设计堵塞和非堵塞两种通信模式去实现计算与通信的重叠.每个子区域边界增加两层虚网格,用于储存临近子区域相关变量.

　　2、三维复杂爆炸流场的数值模拟

　　2.1算例模型一

　　算例模型如图2、图3所示,为一类较为普遍的工事,工事长40m,宽22m,高3m.周围的墙体作为刚性壁处理,其余采用连续边界条件.炸药分别位于A、B、C处,装药量为1吨.计算域尺寸为25.0m×4.50m×45.0,采用等步长网格离散计算域,网格步长为0.05m,共4100万网格.

　　

　　爆轰产物主要参数见表1和表2,其中k0为爆轰产物的多方指数,k1为爆轰产物充分膨胀下的多方指数,b为调节系数.

　　

　　图2显示了冲击波的传播过程及其相互作用,这与实际的物理规律是相吻合的.

　　模型左端的四个分支对称分布,在每个分支设置一个关键点(关键点的分布见图3),它们距离关键点11的距离相同.由图4可以发现关键点14和关键点15的冲击波峰值超压明显高于关键点12和关键点13的冲击波峰值超压,由此可见工事成锐角型的走向能够有效地减小冲击波的毁伤.

　　由数值模拟结果可以得到所有关键点的峰值超压如图4所示.为了便于更清楚看到地下工事内的最小超压值,忽略炸药位置及其附近的关键点峰值超压,在图中不予显示.由图4可知,最小的超压值在0.2MPa左右,而这个压力足以对地下工事的设备以及人员造成严重的毁伤.

　　对比不同位置处爆炸的峰值超压如图5所示.从图中可以明显的看出,炸药在C处爆炸,得到的整体峰值超压明显高于A处和B处,C处为该模型的理想炸点位置.PMMIC-3D并行程序可对工事中爆炸的各种工况进行预估,为寻找最佳炸点位置提供参考.

　　2.2算例模型二

　　算例模型如图6所示,模型长为26m、宽为25m、高为3m.炸药位于模型的中心处,采用球形装药计算.周围的墙体作为刚性壁处理,其余采用连续边界条件.采用等步长离散计算域,计算模型如表3所示.该组算例的爆轰产物主要参数与上述算例所取的参数一致.

　　

　　图7给出了不同时刻冲击波在工事内的传播过程及其相互作用.

　　

　　不同炸药量下各关键点(关键点的分布见图8)处超压的数值模拟结果值如图9所示.

　　由图9中各个位置处的关键点的冲击波峰值超压可以判断对于这样一个模型,当炸药量大于200kg时,大部分位置的峰值超压均达到了5个大气压,而当炸药量为200kg时,各个位置处的冲击波峰值超压比较均衡,但都超过了1个大气压,而这个压力就能使工事内的人员无法生存.

　　3、结论

　　(1)采用MPI标准设计了PMMIC-3D并行程序,算例的数值模拟结果与实际的物理规律相符,验证了所采用的并行算法的合理性以及程序的可靠性;(2)数值模拟结果清晰的显示了冲击波的传播过程及其相互作用,由不同位置处的冲击波峰值超压可以 判断不 同装药量下地下工 事 的 毁 伤程度;(3)通过不同位置处的冲击波峰值超压比较,发现工事成锐角型走向能够有效减小冲击波的毁伤,可为工程设计提供一定的参考.

　　参考文献：
　　[1] 恽寿榕,涂侯杰,梁德寿.爆炸力学计算方法[M].北京:北京 理 工 大 学 出 版 社,1995.
　　[2]Dongarra J,Foster I,Fox G,等.并行计算综论[M].莫则尧,陈 军,曹 小 林,译.北 京:电 子 工 业 出 版 社,2005.
　　[3]Malard J.MPI:A Message-Passing Interface Standard[R].UK:The University of Edinburgh,1994.
　　[4] 宁建国,王成,马天宝.爆炸与冲击动力学[M].北京:国防工业出版社,2010.
　　[5]Fei G L,Ma T B and Hao L.Large-scale high perform-ance computation on 3Dexplosion and shock problems[J].Appl.Math.Mech.-Engl.2011,32:375-382.
　　[6]Ning J G,Chen L W.Fuzzy interface treatment in Eu-lerian method[J].Sci China Ser E,2004,47(5):550-568.