LS-DYNA软件开展爆炸冲击波计算时需考虑的问题

李晓勇, 崔村燕, 陈景鹏, 李 渊, 赵蓓蕾

(1.装备学院研究生管理大队,北京101416; 2.装备学院航天装备系,北京101416)

LS-DYNA软件开展爆炸冲击波计算时需考虑的问题

李晓勇1, 崔村燕2, 陈景鹏2, 李 渊1, 赵蓓蕾1

(1.装备学院研究生管理大队,北京101416; 2.装备学院航天装备系,北京101416)

针对LS-DYNA软件模拟空气中爆炸冲击波的传播过程时边界条件设置对计算结果的影响,验证了采用部分建模的方法对具有对称性的模型进行数值模拟的可行性;研究了无反射边界面对冲击波超压值的影响;提出把计算区域厚度与炸药厚度之比d作为判断三维模拟时模型厚度合理性的标准,并指出在进行冲击波三维数值模拟时d应大于2。

爆炸冲击波;无反射边界;三维模拟;LS-DYNA软件

爆炸冲击波是由于爆炸产生的一种强压缩波[1],它是爆炸所释放能量的扩散和传播的一种形式,是产生破坏作用的主要因素之一,研究其传播规律对于有效防护及利用具有很重要的意义。

由于爆炸实验的危险性高、耗费大,而且实验条件严格且重复性差,致使利用实验来研究冲击波的难度较大。随着计算机数值仿真技术的发展,数值模拟已成为研究爆炸冲击波的主要手段[2]。

LS-DYNA软件由J.O.Hallquist博士1976年在美国劳伦斯利弗莫尔国家实验室主持开发完成,可用于分析爆炸与高速冲击等过程中的大变形动力响应问题,能够对爆炸冲击波的生成、传播与耦合过程进行模拟仿真和数值计算,已经广泛应用于国防、民用等领域[3]。

为提高数值模拟的精度,国内一些研究人员针对模拟过程中如何设置算法做了专门研究,例如,北京工业大学石磊[4]对爆炸冲击波数值计算网格划分方法进行了研究;卢红琴[5]105-108针对网格密度对计算结果的影响进行了分析;绪梅[6]等以爆炸相似率为基础,提出了计算机仿真计算爆炸冲击波对目标作用的等效模型,对利用LSDYNA进行爆炸模型的等效条件做了研究。

在阅读文献时,发现不同研究者在利用LSDYNA软件进行爆炸冲击波相关研究时,对于炸药等效模型和计算区域边界等计算初始条件的设置依据没有给出清晰的解释,从而出现不同研究者在研究相同或者类似问题时,给出的量化值有较大差别,同时也给初涉爆炸冲击动力学计算研究的学生造成很大困扰。针对这种现象,本文针对上述初始条件的设置对于计算结果的影响展开讨论,以期为从事基于LS-DYNA软件爆炸冲击动力计算的研究者提供有益的参考。

1 材料模型和等效模型

数值模拟采用cm-g-μs单位制,基本材料为TNT炸药、空气。TNT炸药采用MAT-HIGHEXPLOSLVE-BURN材料模型[7],爆轰压力P和单位体积内能及相对体积V的关系采用JWL状态方程进行描述。JWL状态方程能精确描述凝聚炸药的实验过程,且具有明确的物理意义,可以对爆炸过程中产生的压力做出与实验结果非常相近的预测。

式中：A、B、R1、R2、ω为JWL状态方程的参数,其值由实验确定;E为炸药的内能;V为当前相对体积。采用等效TNT装药进行数值仿真时,参数设置如表1所示。

表1 TNT炸药参数设置表

空气采用MAT_NULL材料模型和线性多项式状态方程EOS_LNIEAR_POLYNOMIAL进行描述[8]

式中：μ=(ρ/ρ0)-1;ρ为当前密度;ρ0为初始密度;E为炸药的内能;C0、C1、C2、C3、C4、C5、C6为状态方程参数。仿真中参数设置如表2所示。

表2 空气材料参数设置表

本文拟对模型具有对称性的问题,采用部分建模的方法,获得等效的模型来进行数值模拟。例如对于立方体炸药在空气域内爆炸的问题,可以利用炸药体的对称性,只建立其1/4个模型,通过在对称面上施加约束,获得整个模型的等效模型来进行计算。

2 数值模拟

2.1 炸药部分建模等效模型可行性的验证

在利用LS-DYNA进行数值模拟时,针对能否利用模型的对称性,采用部分建模的方法进行数值模拟这一问题,对以下2种模型数值模拟结果的差异性进行了探究。

模型I模拟过程中不考虑对象模型的对称性,采用整体建模方法,进行数值模拟。

模型II模拟过程中考虑炸药及空气模型的对称性,采用部分建模的方法,对整体模型取2次对称后,建立其1/4模型,在对称面上施加对称面法线方向的位移约束。模型I中,炸药为边长25 cm的正方体,炸点为炸药体中心,空气域长度为800 cm,整个炸药体浸没在空气中。模型I中空气域所有边界面设置为无反射边界面;模型II中炸药和空气为模型I沿对称面2次取对称后的1/4,在空气及炸药体的2个对称面上施加法线方向上的位移约束,其余各面设置为无反射边界面,炸点为取对称前炸药体的中心位置。2次模拟中选取与炸药相对位置相同,各点间距都为50 cm的A、B 2组测量点。

最终测得2个模型中各测量点冲击波超压波形如图1所示。

可以看出,2种模拟方法得到的空气冲击波超压波形基本一致。各点处超压值、到达时间及超压作用时间基本相同。

由此可得,在对爆炸冲击波进行数值模拟时,对于对称性模型,可以采用部分建模的方法,在对称面上施加法线方向位移约束后进行数值模拟,所得结果与对模型进行整体建模的结果基本相同,在同样的精度要求下,可节省大量计算时间。

图1 2次模拟中测量点的冲击波超压波形

2.2 计算区域边界设置对模拟结果的影响

理论上,在计算区域内,与爆炸点间距离相同的点处超压值应相同。为探究边界设置的影响,采用2.1中部分建模的方法,建立空气中圆柱体炸药爆炸的1/4模型,进行2次模拟,如图2所示。

图2 数值模拟模型

1/4圆柱体高2.92 cm,炸点等效圆柱体的中心。2次模拟中,空气域中M面、N面及与xoy面平行的2个面都设置为无反射边界面,炸药体S1、S2设置法线方向的位移约束。空气域Q、R面,模拟I中设置为无反射边界面,模拟II中施加法线方向上的位移约束。设置A、B 2组等间隔15.8 cm的检测点。

2次模拟结果1 004μs时刻超压峰值云图对比如图3所示。2种情况下,炸药体的边界条件设置一样,产生的爆炸效果相同,但由于Q、R边界条件不同,模拟I中,空气域Q、R附近,超压值明显小于同时刻空气域内部的超压;模拟II中, Q、R边界附近超压值则与空气域内部无明显区别。

图3 2次模拟冲击波超压云图对比

由爆炸相似率知,爆炸冲击波超压P是比例距离Z=R/W1/3的函数,其中R为测点与爆心之间的距离(m);W为等效TNT药量(kg)。2次模拟中A、B 2组测点处P随Z的变化情况如图4所示。模拟I中,B组超压值小于A组中相对应点处超压值,并呈现出越靠近爆炸点,差异越大的趋势,随Z增加,二者差距逐渐减小,当Z=3时, A、B超压值基本一致。与之相比,模拟II中,Z较小时,A、B监测点超压值相差较大,当Z＞1.0时,超压值基本一致。总体上看,2次模拟中的A组监测点超压值基本相同,但边界处B组监测点超压值差异较大,说明无反射边界面对其附近区域模拟结果有一定影响。

图4 2次模拟中对角线处与底边上超压对比

为进一步研究这种影响与比例距离之间的关系,以模拟I的计算设置为例,选取Z分别为1.5、2.0、2.5、3.0的测量点,改变检测点与x轴间的垂直距离(y),得各监测点超压如图5所示。

图5 不同比例距离下超压值随距边界面距离的变化

由图5可见,爆炸源附近区域,边界面上超压值同空气域内超压值相差较大,随Z增加,差距逐渐减小,边界面对P的影响逐渐减弱。

利用数值模拟研究冲击波传播或绕流的过程中经常设置无反射边界面,由于Z不同情况下,边界面对超压峰值的影响不同,因此模拟过程中要针对研究对象所处的比例距离范围,合理的建立模型并选取测量点,保证所关注区域及测量点与无反射边界面的距离足够远,以免造成模拟结果有较大偏差,得出错误结论。如文献[9]中给出,当P分别为0.101 35～0.202 70 MPa和0.202 70～0.101 35 MPa时,空气冲击波对建筑结构的破坏程度分别为防震钢筋混凝土破坏、小房屋倒塌和大型钢架结构破坏。根据本文中的模拟结果,在比例距离为1.5处,靠近边界面上的超压值为0.17 MPa,而空气域内部得出的超压值则大于0.27 MPa,即距离相同,得到的破坏程度不同,应用在实际工程中,导致防护失当。

2.3 三维模拟中计算区域边界与炸药间的距离对模拟结果的影响

冲击波三维数值模拟中无反射边界面常被应用于设置边界条件。为探究在三维数值模拟中可采用的最小厚度值,建立z轴方向厚度l分别为37.5、50、75、125和150 cm的5个模型,空气域长度1 000 cm。采用2.1中所述部分建模的方法进行数值模拟。炸药为12.5 cm×12.5 cm× 25 cm的长方体,炸点在等效立方体的中心,爆炸效果相当于边长q为25 cm的正方体爆炸。5组模型的其余设置都相同,定义模型的相对厚度d=l/q,则相对厚度分别为1.5、2.0、3.0、5.0、6.0 cm。

每个模型均在炸药上方水平高度150 cm处取11个测量点,各点相邻点间距50 cm。得不同厚度模型各测量点超压值随Z的变化如图6所示。其中P经验公式值根据式(3)得来,此公式为1979年亨里奇(Henrych)得出,适用于各种化学炸药[10],许多学者将其作为参考值[5,11]。

图6 不同模型厚度模拟结果

由图6可见,随着d值增加,测量点受无反射边界面的影响减弱,测得超压值逐渐增大,越靠近爆炸点,与经验公式值差异越大。d=1.5时,冲击波超压明显要偏小;d＞2.0的各个模型中,模拟结果逐渐趋近于经验公式值,并在Z＞1.5的区域能够与经验公式吻合的较好,当d＞6.0时,在Z＞1.5的区域测量值略高于经验公式值。

结果表明,计算区域边界与炸药的距离对模拟结果具有一定的影响,主要体现在无反射边界面附近区域超压值偏小,因此在利用LS-DYNA对冲击的传播进行数值模拟时,应充分考虑这一因素。在对冲击波进行三维数值模拟时,其模型的相对厚度d至少应大于2,同时在数值模拟中需要针对所研究的问题,在模型建立和测量点的选取上要注意合理性,最大程度的减小无反射边界面对模拟结果的影响。

2.4 数值模拟结果与实验结果比较进行验证

文献[12]12中对于空气中爆炸冲击波超压的传播做过如图7所示实验,实验中用炸药柱直径和高度都为30 mm,炸药质量为32 g,3个传感器测量点如图7中所示,距离炸药中心距离分别为80、90、100 cm,相应比例距离为2.52、2.38和3.15。

图7 空中爆炸实验

为验证2.3结论的正确性,根据这一实验建立模型相对厚度d分别为5.0和1.5的A、B 2模型进行模拟,空气域所有界面均设置为无反射边界面,a、b、c 3点与实验中测量点相对应。由于数值模拟中炸药模型的密度与实验中炸药密度略有差异,根据炸药质量相等的原则,建立直径和高度都为29.2 mm,质量为31.86 g的炸药圆柱,炸点为炸药中心。炸药中心距离空气域左、右界面距离分别为95 cm、115 cm,距离上、下界面距离为127.2 cm、15 cm。

文献[5]106中指出,炸药体网格划分密度对模拟结果影响较大,炸药体边长与单元格边长之比大于8时,模拟结果精确度较高,本次模拟设置为10。最终得2组数值模拟结果如图8所示,b、c测量点超压值分布与空气域内部超压处于同一颜色水平,测量点的选取达到了预期效果。模型B与模型A相比,受与xoy面平行的无反射边界面的影响,其超压值偏小。

为验证结果的准确性,在炸药上方取一组间隔为15.87 cm(等效比例距离为0.5)的监测点,并将模拟结果同亨里奇提出的经验公式值及根据文献[12]12中的数据作对比,如图9所示。

由图9可以看出,Z＜1.5时,模拟结果比经验公式和文献数据值小,Z＞1.5时,模拟A结果与根据亨里奇经验公式得出的结果较为吻合,在Z＞2的区域内,其结果同文献[12]12及经验公式的结果都很接近,但模型B中各比例距离处,其超压值仅为模型A的50%左右。

图8 模拟结果

图9 模拟结果同经验公式及文献[12]12中数据对比

图10 模型A中与实验对应测量点处模拟冲击波超压波形

测得A模型a、b、c的冲击波波形如图10所示。图7实验中的A、B、C 3测量点超压分别为0.074、0.117、0.092 MPa,与之对应数值模拟结果为0.072、0.111、0.088 MPa,二者最大误差均不超过6%,在模型厚度d=5时,模拟结果与实验值符合的很好,得出的结果具有一定的说服力,很好的验证了2.3节中关于三维模拟时模型厚度选取的结论。

3 结 论

改变LS-DYNA对爆炸冲击波进行数值模拟的边界条件,得到以下结论：

一是采用部分建模方法进行数值模拟,得到的结果与采用整体建模的方法差别不大,因此实际计算时可采用部分建模的方式,以节省计算时间。

二数值模拟的结果会受到无反射边界面的影响,这种影响在爆炸源附近尤为明显,因此在进行数值模拟时,要针对研究对象特点,合理选取测量点,最大程度削弱无反射边界面的影响。

三是进行三维数值模拟时,模型厚度最小应不低于3倍炸药体的厚度,在进行有障碍物情况的绕流研究时,模型厚度还应继续增加。

References)

[1]恽寿榕,赵衡阳.爆炸力学[M].北京：国防工业出版社, 2005：2-3.

[2]刘科种,徐更光,辛春亮,等.近水面水下爆炸对结构冲击的数值模拟[J].兵工学报,2010,31(1)：64-68.

[3]尚晓江,苏建宇,王化峰.基于LS-DYNA动力分析方法与工程实例[M].北京：中国水利水电出版社,2008：1-3.

[4]石磊,杜修力,樊鑫.爆炸冲击波数值计算网格划分方法研究[J].北京工业大学学报,2006,36(11)：1465-1470.

[5]卢红琴,刘伟庆.空中爆炸冲击波的数值模拟研究[J].武汉理工大学学报,2009,31(19)：105-108.

[6]绪梅,王良厚.爆炸冲击波对目标作用的等效计算模型[J].防化研究,2004(3)：17-23.

[7]Livermore Software Technology Corporation.LS-DYNA keyword user's manual[M].California：Livermore Software Technology Corporation,2003：2042-2044.

[8]时党勇,李裕春,张胜民.基于ANSYS/LS-DYNA 8.1进行显式动力分析[M].北京：清华大学出版社,2005：244-247.

[9]陈新华,聂万胜.液体推进剂爆炸危害性评估方法及应用[M].北京：国防工业出版社,2005：126.

[10]亨利奇.爆炸动力学及其应用[M].熊建国,译.北京：科学出版社,1987：127.

[11]马云玲,赵丽君,聂建新.异型防爆墙抗空气冲击波的数值模拟[J].爆破,2010,27(1)：26-30.

[12]宋浦,肖川,梁安定,等.炸药空中与水中爆炸冲击波超压的换算关系[J].火炸药学报,2008,31(4)：10-13.

(编辑：田丽韫)

Several Considerations When Calcutating the Blast Wave Using LS-DYNA

LI Xiaoyong1, CUI Cunyan2, CHEN Jingpeng2, LI Yuan1, ZHAO Beilei1

(1.Department of Graduate Management,Equipment Academy,Beijing 101416,China; 2.Department of Space Equipment,Equipment Academy,Beijing 101416,China)

Boundary conditions have great effects on numerical results when LS-DYNA is used in simulating the blast wave.The practicability of using the part model instead of the whole model is verified in this paper.The influence of the non-reflecting boundary on overpressure is analyzed.The ratio of the calculation area thickness and the explosive thickness,which is named d,is appointed as a standard to judge whether the calculation area thickness is reasonable.And the value of d should be more than 2 when three-dimensional calculation is carried out.

air blast wave;non-reflecting boundary;three-dimensional simulation;LS-DYNA soft

O 38

2095-3828(2014)04-0079-06

ADOI10.3783/j.issn.2095-3828.2014.04.018

2014-01-18

部委级资助项目

李晓勇(1988-),男,硕士研究生.主要研究方向：发射场安全.崔村燕,女,副教授,硕士生导师.