光滑粒子法与有限元的耦合算法及其在冲击动力学中
第27卷第6期
2007年11月爆炸与冲击EXPLOSl0NANDSHOCKV01.27.No.6Nov..2007WAVES
文章编号:lOOl一1455(2007)06一0522一07
光滑粒子法与有限元的耦合算法
及其在冲击动力学中的应用。
王吉,王肖钧,卞梁
(中国科学技术大学力学与机械工程系,安徽合肥230026)
摘要:扼要讨论了光滑粒子法的离散思想,充分利用光滑粒子法和有限元方法各自的优点,提出了一种
初始时刻用有限元建模,计算过程中大变形单元自动转换为光滑粒子的耦合算法。高速碰撞的系列算例说
明,耦合算法不但适宜于计算大变形冲击动力学问题,而且由于集两种方法的优点于一身,可以更高效地模拟
一些高速碰撞问题,提高计算效率。
关键词:固体力学;耦合算法;光滑粒子法;高速碰撞;有限元
中图分类号:0347;0413.1国标学科代码:130・1570文献标志码:A
1引言
高速碰撞是典型的冲击力学问题,高速碰撞数值模拟计算一直是冲击工程数值计算的难点。当采用Lagrange有限元或有限差分方法计算高速碰撞问题时,必然遇到网格大畸变和滑移面处理等一系列数值计算中的关键问题,他们常常成为Lagrange有限元或有限差分法在工程应用中的制约因素。虽然Euler方法不存在计算大畸变问题的困难,但是难以准确描述各类界面。因此寻找适宜于高速碰撞数值模拟的新的计算方法,一直是计算力学工作者的努力目标。光滑粒子法(smoothed
namicsparticlehydrody—method)便是其中的一种典型方法。光滑粒子法是一种纯Lagrange无网格算法,通过带质量的粒子离散计算域,粒子本身便代表材料,不同材料的粒子自然地构成界面,不同材料粒子的相对运动便形成所谓界面的滑移。因此从理论上说,光滑粒子法可以比较“自然地”模拟高速碰撞、侵彻贯穿等物理现象,因此是冲击力学中一种比较理想和有前途的数值方法,但其计算效率要低于有限元法。为了提高计算效率,G.R.Johnson等[1屯]将光滑粒子法与有限元方法耦合起来,发挥他们各自的长处,在变形大的地方采用光滑粒子方法,其他区域使用有限元算法,开展侵彻贯穿方面的数值计算,取得了一些有意义的结果。
本文中在G.R.Johnson等[1吨]的工作基础上,提出一种初始时刻用有限元方法建模、随着变形的增大、大变形区域的有限元单元自动转换成光滑粒子的耦合算法。并详细讨论耦合算法中的几个关键问题,如粒子生成的条件与方式、滑移计算与界面计算等。通过经典Taylor碰撞、高速和超高速碰撞等问题的计算实例,分析耦合算法的长处和特点。
2光滑粒子法
J.K.Chen等口]从Taylor展开式出发,引入核函数的概念,建立了光滑粒子法的基本算式。设函数,(工’),对粒子i所在位置工;进行Taylor展开,并引入核函数w;=Ⅳ(一一藏I,^),其中^为光滑长度。在解域n中积分,便可得到函数,(x’)一阶偏导数核估计的求和形式
・收稿日期:2006一05—22;修回日期:2006-09一04
基金项目:国家自然科学基金项目(10572134)作者简介:王吉(1979一),男,博士研究生。
第6期王吉等:光滑粒子法与有限元的耦合算法及其在冲击动力学中的应用523
怠§竺!厶二五!鲨幽一型§竺』!茎二兰!!里g!£lD』az?台f1、”7
10』
式中:N表示粒子数,由于核函数w是强尖峰函数,因此实际计算时,只需考虑w某一邻域里的粒子,mJ是粒子zJ的质量,WⅢ表示离散核函数Wd对空间坐标z4的偏导数。(1)式是以空间坐标z}为方向的代数方程组,求解该方程组,便可获得函数厂(x’)沿三个方向的一阶偏导数。可以证明由此获得的一阶偏导数的核估计,对于内点精度为二阶,而对于边界点则为一阶。采用(1)式对守恒方程中的空间导数进行离散,便可得到光滑粒子法的基本算式。
3光滑粒子法和有限元法的耦合
在高速和超高速碰撞的数值模拟计算中,Lagrange有限元常常遭遇局部区域一些网格畸变过大,以及侵彻面推进过程中两侧单元滑移计算的困难。虽然光滑粒子法可以自然地避开这些问题,但是其计算效率远低于有限元方法,因此其工程应用受到了很大限制。考虑到网格畸变过大和滑移计算只发生在局部区域的一些个别单元上,倘若能将出现计算困难的单元自动转换为光滑粒子,按光滑粒子法计算,其余区域依然用有限元方法计算,则可充分发挥两者的优势,切实解决高速碰撞中的计算难题。然而耦合算法也提出了新的课题:(1)如何实现有限元单元向光滑粒子的自动转换;(2)光滑粒子与有限元单元间的界面计算问题。
3.1有限元单元向光滑粒子的转换
Lagrange有限元在高速侵彻数值计算中,滑移面两侧单元常常是大畸变单元,因此需要依据网格具体变形情况设计光滑粒子生成的条件。在耦合算法中,我们借用有限元计算中处理大畸变单元的思想,以等效塑性应变作为单元转换的判据,即当位于滑移面两侧Lagrange单元的等效塑性应变达到某一特定值时,将此单元转换成光滑粒子。此时,光滑粒子的质量、体积和应力都等同于被置换Lagrange单元相应的物理量,而粒子的几何量,如坐标和速度则取为该Lagrange单元节点的平均值,再利用光滑粒子的速度分布,求出该时刻粒子的应变率张量和旋转率张量。光滑粒子的光滑长度^可以近似通过质量守恒条件获得
矗=忍o(风/lD)17”
间的维数,轴对称情况下取u一3。
图1为位于滑移面上的Lagrange单元向光滑
粒子转换的示意图,虚线N为初始滑移面,实线M
为单元转换以后的滑移面。以该时刻为例,滑移面
上有4个单元满足转换条件而转换成光滑粒子(黑
色的圆)。
应该指出,有限元算法里,等效塑性应变是作为
单元因畸变过大而破坏所设立的判据,但在耦合算
法里,等效塑性应变是单元转换为粒子的阈值条件,
两者的物理意义并不相同。耦合算法里的等效塑性
应变直接关系到光滑粒子的个数。由于粒子数的多
寡对计算结果的影响并不是本文研究的重点,因此
在下面的算例里将参考G.R.Johnson的工作‘¨,将Fig.1
阈值取为0.5,以说明本文方法的有效性。
3.2光滑粒子与有限元交界面的滑移计算
界面处理也是高速碰撞计算中的关键问题,在耦合算法中,同样也存在着光滑粒子与有限元界面
(滑移面)之间的接触问题。如图2所示,图中每一个圆都代表一个光滑粒子,当粒子侵入滑移面M时图1(2)式中:^。表示初始设定的光滑长度,它的选取带有经验性,并和初始的Lagrange单元尺寸相关,u为空L。g,。ng。单元与光滑粒子的转换示意图DescriptionofthesPHnodegenerationalgorithm
524爆炸与冲击第27卷
rml口l+m2耽+mf口f=mlul+m2钝+m∥f
l优l可1r1+m2忱r2+mf可ir‘=ml西r1+m2谠r2+m∥:ri
lu:一西=(u:一西)(r{一,.。)/(r。一r。)(3)
Areaofsmoothedparticles
SlidinginterfAce
AreaofLagrangefjniteelemems
图2Lagrange单元与光滑粒子的界面模型
Fig.2图3交界面附近粒子的计算示意图Fig.3ComputationsoftheSPHnodes
nearSPHnodeslidingonastandardgridtheslidinginterface
4计算实例和讨论
4.1Taylor碰撞
Taylor碰撞是典型的高速碰撞问题。1948年,G.I.Taylor研究了一平头柱形杆弹法向撞击刚性
靶板问题。从一维弹塑性波理论出发,采用刚塑性本构模型,解析获得了杆弹撞击后柱体变形与材料动态屈服强度间的近似关系。由于Taylor碰撞体现了材料大变形高应变率等典型的冲击力学特性,因此通过Taylor碰撞研究材料力学性能一直受到关注。w.H.Gust[4]和G.R.Johnson等[23都先后通过数值方法研究了Taylor碰撞问题,通过与实验结果的比较,验证本构模型和计算方法的合理性和准确性。为此我们采用光滑粒子法及耦合算法对Taylor碰撞开展数值模拟,以便与已有的理论研究和实验资料相比较,确定本文所述方法和所编程序的可靠性。
设直径Do=7.6mm、长度Lo=25.4mm的钢柱,以速度矾一221m/s法向撞击刚性平板,材料参
数分别为[5]:密度胁一7.890T/m3,冲击压缩实验中测得的D“线性关系中的常数c。一3.6km/s,s=1。80,剪切模量户=80GPa,简单拉伸条件下的屈服应力y0=o.5GPa,Gruneisen系数r=1.81。设材
第6期王吉等:光滑粒子法与有限元的耦合算法及其在冲击动力学中的应用525料为理想弹塑性体。图4给出了碰撞后的变形示意图,定义变形后柱体的长度为L,撞击端直径为D,离固壁o.2L0处的直径为W。张刚明等‘63给出了同样算例的实验结果,L=19.8mm,D一13.7mm,W
、
=8.6mm。
分别采用纯光滑粒子法和耦合算法对Taylor碰撞开展数值模
拟计算。初始时刻,耦合算法以三角形单元建模,共658个单元,单
元数与纯光滑粒子法模型的630个粒子相近,计算过程中,随着有限
元单元变形的发展,部分单元自动按上述方法转换成光滑粒子。
表1计算和实验结果的比较
TabIe1C伽叩n幅ionbetweencOmputatiO璐andexperiment
图4TayIor碰撞示意图
Fig.4Taylorimpact
图5~6分别给出了纯光滑粒子法和耦合算法在lops和40ps时的变形图,耦合算法中杆弹上部仍为有限元,下部已转换为光滑粒子。表1给出了计算值与实验结果的比较。虽然几个特征尺寸耦合算法要比光滑粒子法更接近实验值,但是从撞击端的总体构形上看,两者仍有差异,这说明耦合算法中粒子的生成和界面算法对结果有一定影响。通过与实验结果的比较,进一步改进相关算法,将是本文后续工作的重点。然而Taylor碰撞的初步计算结果表明,采用本文给出的思想建立耦合算法是成功的。
图510肛s时变形图
Fig.5Deformationatz—lOps
图640ps时变形图
Fig.6Deformationat‘一40
ps
爆炸与冲击第27卷
4.2有限元、光滑粒子法和耦合算法的对}匕计算
钨合金平头杆弹,弹长20mm、直径10mm,以初速‰=1.毒km/s正撞击厚为10mm的钢靶,分别采用有限元法、光滑粒子法和耦合算法进行对比计算,靶板和杆弹视为理想弹塑性体。表2给出了计算中用到的两种材料的参数。
表2两种材料的参数
TabIe2Pa甩mete伟0ftw0materials
采用有限元程序计算时,单元取为三角形单元,共2048个单元。采用光滑粒子法计算时,共1024个粒子。而在耦合计算中,初始模型只取了800个三角形单元。
图7~8分别给出了三种方法在不同时刻(3弘s和10肛s)的计算结果。对于畸变过大的单元,有限元方法用消蚀法处理,即当单元等效塑性应变超过某一阈值时,该单元被认为破坏,只能承受静水压力。计算中钨合金和钢靶的极限等效应变分别取为1.o和o.9。光滑粒子法由于不存在单元,等效塑性应变增大,相当于粒子距离减小,因此不会出现有限元方法中因单元变形太大导致计算无法进行下去的情形。而在耦合算法中,由于只有部分单元转换成粒子(如3ps时为98个粒子,10弘s时为310个粒子),因此既能有效处理畸变过大的计算问题,又保留了有限元计算效率高的优势。从图中可以看到,三种方法的变形图像相差不是很大。由于在有限元方法中,破坏的单元对子弹不再存在阻力作用,因此相同时刻杆弹的速度要比光滑粒子法计算的结果大。有限元方法中子弹前部速度比较小的单元由于达到了破坏条件,成为破坏单元,在图上没有画出。光滑粒子法中,由于子弹前部的变形较大,形成的蘑菇头要比有限元方法的大,从而开坑半径要大一些。而耦合算法计算结果更趋近于光滑粒子法。
图73肚s时变形图
Fig.7Deformationat£=3ps
图810ps时变形图
Fig.8Deformationatt=lO肛s
4.3超高速碰撞计算
考虑高速平头弹丸侵彻半无限厚靶板的情况,子弹长20mm、弹径10mm,初始弹速4km/s,子弹和靶板的材料与算例2相同。光滑粒子法共用了10496个粒子。耦合算法中共用了8200个单元(粒
子)。图9~11分别给出了两种算法在几个不同时刻的变形图。图上颜色淡一点的是子弹的粒子,而颜
第6期王吉等:光滑粒子法与有限元的耦合算法及其在冲击动力学中的应用527
图9
Fig.95弘s时变形图Deformationat£一5肛s
图1010肛s时变形图
Fig.10Deformationatf=10_c上s
图1120弘s时变形图
Fig.11Deformationat£一20ps
色深一点的是靶板的粒子,耦合算法中有大量未转换成粒子的有限元单元。从图上可以清晰地看出靶板的成坑过程。刚开始时,靶板表面形成一个比弹径略大的坑。随着子弹的不断侵彻,头部逐渐形成蘑
菇头,阻碍子弹的侵彻,而径向速度不断增大,使得开坑半径变大,同时在耦合算法中随着变形的增大,
528爆炸与冲击第27卷有更多的单元转换成粒子。随着粒子径向速度的消耗,到一定时候坑的径向扩展就停止了,而坑的轴向仍在发展直至弹体粒子的动能消耗殆尽。两种算法的计算结果相差不大,耦合算法所得的坑深和开坑半径略小于光滑粒子法。总体说来,用光滑粒子法或耦合算法模拟高速碰撞贯穿现象是合适的,变形图像与实际比较相似。用有限元方法计算超高速碰撞时会遇到很大困难,即便采用消蚀法,也不可能象光滑粒子法一样模拟出开坑过程。不过光滑粒子法中的粒子搜索很费机时,计算效率要明显低于有限元法,而耦合算法集两者的优点于一身,其计算效率要远高于光滑粒子法,这一点在计算模型的单元较多时显得尤为明显。以此算例为例,采用光滑粒子法计算需数个小时,而耦合算法中,由于靶板的绝大部分计算域始终是有限元单元,光滑粒子只集中在弹靶接触的局部区域,粒子数并不多(20弘s时为1个粒子),因此整个计算过程只需十分钟多,计算效率得到了大幅度的提升。612
参考文献:
[1]JohnsonGR.LinkingofLagrangianparticlemethodst0standardfiniteelementmethodsforhighvelocityimpact
computations[J].Nuclear
[2]Er培ineeringandDesign,1994.50:265—274.testJohnsonGR.HolmquistTJ.Evaluationofcylinder_impact
ofdataforconstitutivemodelconstants[J].JournalAppliedPhysics,1988,64(8)l3901—3910.
JK,BeraunJE.JihcJ.Animprovementfortensileinstabilityinsm00thedparticlehydrodynamics[J].[3]Chen
ComputationalMechanics,1999,23:279—287.
[4]GustwH.Highimpactdeformationofmetalcylindersatelevatedtemperatures[J].JournalofAppliedPhysics,
1982,53(5)13566—3575.
[5]LiberskyLD,PetschekAG.HighstrainLagrangianhydrodynamicslAthree—dimensionalSPHcodefordynamic
materialresponse[J].JournalofComputationalPhysics,1993,109:67—75.
[6]张刚明,王肖钧.高速碰撞数值计算中的光滑粒子法[J].计算物理,2003,20(5):447-454.
ZHANGGarIg-mirlg,WANG
velocityXia州un.Sm∞thedparticlehydrodynamicsmethodtonumericalsimulationofhyper‘impact[J].ChineseJournalofComputatiorlalPhysics,2003,20(5):447—454.
Linking0fsm00thedparticlehydrOdynamicsmethodtOstandard
finiteelementmethodanditsapplicationinimpactdynamics
WANGJi+,WANGXiao-jun,BIANLiang
(Lhf铆Pr5f£yo厂SciP行ce口,ld:心c^nDZogyD,CI^fnn,HP^f230026,A扎^Mf,C冼f行口)
Abstract:Thispapergivesabriefdescriptionofsmoothedparticlehydrodynamics(SPH)method.BytakingtheadvantagesofSPHandstandardfiniteelementmethods,acoupledmethodisproposed,inwhichthecomputationalmodelismadeupoffiniteelementsininitialtimeandthelargedistortionalelementsareconvertedautomaticallyintosmoothedparticles.Numericalexamplesforhypervelocity
aimpactbythecoupIedmethodshowthatthecoupledmethodis
oflargedeformationproblemsinimpactdynamics.moreefficientmethodtosimulation
Keywords:solidmechanics;coupledmethod;SPH;hypervelocityimpact;finiteelement
・CorrespondingauthorlWANGJi
E-mailaddressljily456@163.comTelephone:86—551—3606734
光滑粒子法与有限元的耦合算法及其在冲击动力学中的应用
作者:
作者单位:
刊名:
英文刊名:
年,卷(期):
被引用次数:王吉, 王肖钧, 卞梁, WANG Ji, WANG Xiao-jun, BIAN Liang中国科学技术大学力学与机械工程系,安徽,合肥,230026爆炸与冲击EXPLOSION AND SHOCK WAVES2007,27(6)12次
参考文献(6条)
1.Johnson G R Linking of Lagrangian particle methods to standard finite element methods for highvelocity impact computations 1994
2.Johnson G R;Holmquist T J Evaluation of cylinder-impact test data for constitutive model constants
[外文期刊] 1988(08)
3.Chen J K;Beraun J E;Jih C J An improvement for tensile instability in smoothed particlehydrodynamics[外文期刊] 1999
4.Gust W H High impact deformation of metal cylinders at elevated temperatures 1982(05)
5.Libersky L D;Petschek A G High strain Lagrangian hydrodynamics:A three-dimensional SPH code fordynamic material response 1993
6.张刚明;王肖钧 高速碰撞数值计算中的光滑粒子法[期刊论文]-计算物理 2003(05)
本文读者也读过(2条)
1. 陈敏.隋允康.阳志光 爆炸力学数值模拟中有限元模型讨论[会议论文]-2007
2. 卞梁.王肖钧.章杰.BIAN Liang.WANG Xiao-Jun.ZHANG Jie SPH/FEM耦合算法在陶瓷复合靶抗侵彻数值模拟中的应用[期刊论文]-高压物理学报2010,24(3)
引证文献(12条)
1.刘天生.张晋红.李长顺 三种算法在侵彻模拟中的对比研究[期刊论文]-弹箭与制导学报 2009(3)
2.顾卫平.骆卫东.童文.史佑民 复合材料蜂窝夹层板低速冲击损伤模拟及实验研究[期刊论文]-机械科学与技术2013(7)
3.卞梁.王肖钧.章杰 SPH/FEM耦合算法在陶瓷复合靶抗侵彻数值模拟中的应用[期刊论文]-高压物理学报 2010(3)
4.张晓天.贾光辉.黄海 基于节点分离Lagrange有限元方法的超高速碰撞碎片云数值模拟[期刊论文]-爆炸与冲击2010(5)
5.初文华.张阿漫.明付仁.杨文山 SPH-FEM耦合算法在爆炸螺栓解锁分离过程中的应用[期刊论文]-振动与冲击2012(23)
6.崔溦.宋慧芳.张社荣 土中爆炸作用下箱涵动力响应的SPH-FE耦合分析[期刊论文]-爆炸与冲击 2012(5)
7.肖毅华.胡德安.韩旭.杨刚 一种自适应轴对称FEM-SPH耦合算法及其在高速冲击模拟中的应用[期刊论文]-爆炸与冲击 2012(4)
8.史晓宁.温垚珂.王亚平 基于SPH算法的钢球高速撞击肥皂数值仿真[期刊论文]-计算机技术与发展 2012(1)
9.张晓天.贾光辉.黄海 基于FE重构方法的冲击破碎仿真[期刊论文]-计算力学学报 2011(5)
10.肖毅华.胡德安.韩旭 一种有限元-光滑粒子流体动力学耦合算法[期刊论文]-计算物理 2011(2)
11.姜峰.唐晓成.任冰.陈洪洲.张莉 基于SPH方法沉箱式防波堤的数值研究[期刊论文]-沈阳建筑大学学报(自然科
学版) 2013(5)
12.胡德安.韩旭.肖毅华.杨刚 光滑粒子法及其与有限元耦合算法的研究进展[期刊论文]-力学学报 2013(5)
引用本文格式:王吉.王肖钧.卞梁.WANG Ji.WANG Xiao-jun.BIAN Liang 光滑粒子法与有限元的耦合算法及其在冲击动力学中的应用[期刊论文]-爆炸与冲击 2007(6)