学文网摘 要：对不同冲击速度下花岗岩靶板的裂纹扩展分布情况进行了数值模拟分析。采用非线性动力学分析软件Autodyn，Johnson-Homquist II（JH-2模型）损伤本构模型，采用SPH算法，对花岗岩靶板进行冲击模拟，得到不同冲击速度下裂纹扩展分布的情况，研究不同冲击速度对花岗岩靶板裂纹扩展及起裂位置的影响。得到起裂位置、主裂纹垂直长度与不同速度的关系曲线。结果表明：当速度低于15m/s时，起裂位置距中心轴线较远且主裂纹扩展长度较短；当速度大于20m/s时，起裂位置距中心轴线较近且主裂纹扩展长度增加，裂纹分布密度增大，当达到一定速度后裂纹扩展基本不变。
关键词：花岗岩；冲击；数值模拟；SPH算法；裂纹扩展
1 概述
花岗岩由于其硬度高、耐磨性好、耐酸性好，被广泛应用于建筑装修、防护工程、化工等方面，素有“岩石之王”之称。花岗岩等脆性材料在冲击作用下的破坏过程一般是损伤、损伤积累、破碎、失效、主裂纹及裂纹群扩展发生的过程，也就是岩石等脆性材料内部裂纹的萌生、扩展、断裂的过程[1]。
许多学者对花岗岩等脆性材料裂纹扩展方面做了很多研究，秦飞[2]等用边界元法模拟了多裂纹的扩展。胡柳青[3]等对在不同冲击作用下裂纹的动态响应进行了数值模拟分析。王海兵[4]等对不同速度下弹丸冲击花岗岩标靶产生的裂纹分布情况进行了实验和数值分析，并归纳总结了一定速度范围内弹坑深度及最大裂纹长度与弹丸速度之间的关系表达式。Motamedi[5]在网格不发生改变的情况下研究了复合材料中动态裂纹完整性的扩展问题。Daux[6]等建立了用于分析随机分叉和交叉裂纹问题的扩展有限元法。周小平[7]等对压应力状态下的多裂纹扩展问题采用扩展有限元法进行模拟分析。
传统有限元法（Finite Element Method，FEM）以网格划分为基础，在处理高应变率、大变形等问题时会导致网格扭曲，导致计算精度下降甚至计算难以继续等问题。光滑质点流体动力学法（Smooth Partcle Hydrodynamic Method，SPH），属于纯拉格朗日方法，无需进行网格划分，从而能较好的进行高应变率、大变形等不连续问题的处理。随着计算机仿真技术的发展，应用SPH法对塑性金属材料切削过程的仿真得到了成功的应用，但是应用于脆性材料和准脆性材料在冲击侵彻过程中内部裂纹扩展机理研究的文献较少。文章借助于非线性动力学分析软件Autodyn采用SPH算法基于Johnson-Homqu-istII（JH-2模型）损伤本构模型，对锤面为平面的锤头在不同速度下冲击花岗岩靶板进行裂纹扩展数值模拟，并分析在冲击作用下裂纹扩展机理。
2 数值模拟
2.1 SPH方法
光滑离子流体动力（SPH）法可以广泛的对脆性材料和准脆性材料等解体、破碎、固体的层裂、脆性断裂等大变形问题进行数值模拟仿真，不需要重构网格，并且能保证计算的精度。其核心思想是将整个流场的物质离散为一系列具有质量、速度和能量的粒子，每个粒子具有单独的速度、能量、质量特征。是一种纯拉格朗日的粒子方法，即无网格（meshfree）算法，且逻辑简单[8]。标准算法执行步骤流程***（如***1所示）[9]。
SPH法是以插值理论为基础的算法，借助核函数对一点上场变量的值给出核估计，把偏微分形式控制方程转化为积分形式，其粒子运动信息的近似函数定义为：
在SPH算法[10，11]的计算上领域搜索是一个重要的步骤，在计算时，单个指点的影响范围为以2h为半径的球形区域，如***3所示。其目的是在每个时间步列出该区域的所有质点。当质点分离的时候，光滑长度h的距离就会增加，相反，当质点汇聚的时候，光滑长度h的距离就会减小。光滑长度h的可变性使得球形领域内的质点数量保持不变，但需设置最大值与最小值，表示为：
式中：h0为初始光滑长度，若HMIN，HMAX均为1时，则h为固定光滑长度，不随空间和时间变化。
2.2 Johnson-Homquist II（JH-2）损伤本构模型
损伤本构模型常常是用来描述岩石的基本力学行为，过去常用的本构模型有流体模型、包络线模型、Mohr-Coulomb模型、流体-弹塑性模型、RHT模型及Johson-Cook模型等，但是这些模型都不能完整的描述岩石破碎的变形和损伤行为[12]。
Johnson-Homquist II（JH-2）损伤本构模型，假定脆性材料（无论是未损伤还是完全损伤的脆性材料），只考虑其损伤是塑性压缩和剪切造成的。其状态方程主要包括应变率、静水压力以及拉伸应力和多项式等。冲击过程中材料性能最先体现为弹性应变，直到应力大于等于材料的屈服极限时，材料开始产生损伤积累并发生断裂。随着损伤的积累，材料最终完全破碎[13]。
JH-2材料模型（如***4所示）是引入等效应力来表示岩石等材料的一种动态损伤本构模型，其函数与幂函数形式的静水压力与应变率和损伤因子相关[4，14]。
损伤因子D的函数表达形式与等效塑性应变相关，其表达形式为：
式中：K1为体积模量，K2和K3为材料常数；？滋为体积应变；？籽为当前材料密度，？籽0为初始材料密度。
当变形无法恢复时，岩石出现损伤，起裂开始。表示为：
式中：PC为压碎压力；PL为压实压力；？滋C为压碎点体应变；？滋L为压实点体应变。
当起裂开始后，微裂纹逐渐扩展并交叉贯通，岩石损坏破碎开始。表示为：
3 仿真模型建立与结果分析
3.1 仿真模型及材料参数
文章以平面锤头，对花岗岩进行10m/s、15m/s、20m/s、25m/s、30m/s不同速度下的垂直冲击数值模拟仿真，锤头采用Autodyn材料库中TUNG.ALLOY材料，采用该材料的Shock状态方程和Johnson-Cook材料模型。花岗岩采用文献[4]JH-2材料模型参数，JH模型可更好的模拟冲击过程中裂纹的形状，数量和长度。材料模型和材料参数如表1和表2所示。
表1 花岗岩材料模型
3.2 二维模型建立
以锤面为平面的锤头，采用二维算法，锤头和花岗岩靶板均采用SPH算法，锤头宽15mm，高10mm和花岗岩靶板高200mm，宽150mm，如***6所示，计算模型如***7所示。定义粒子尺寸为0.25，采用轴对称的形式。
3.2.1 裂纹形成过程
为研究观察裂纹演化过程，选取平头锤分别以10m/s处的速度垂直冲击花岗岩靶板进行分析，观察裂纹演化过程。如***8所示为10m/s时不同时刻的损伤状态及裂纹分布模拟情况。
冲击速度以10m/s时不同时刻的裂纹扩展情况如***8所示，其中***8（a）和（b）演示了冲击初期花岗岩靶板受冲击时的损伤状态，并有细小裂纹产生，且起裂位置位于靶板顶部轴线中心左右两侧52.5mm处，主裂纹扩展垂直高度为42mm，在30μs时，***8（b）由于受压缩和剪切力的作用使锤面下方出现“Λ”型裂纹。***8（c），（d），（e）和（f）展示了在拉应力下损伤积累后裂纹的扩展，在40μs之后，上部裂纹扩展停止，下部和中部损伤开始积累并逐渐形成裂纹并开始扩展。在50μs后，上部裂纹扩展开始继续扩展且起裂时出现变向，裂纹开始扩展时逐渐趋向于原方向，且靶板的中下部分出裂纹开始扩展。80μs后主裂纹几乎不再扩展，微裂纹裂扩展速度减缓。***8（g）为靶板最终状态下的裂纹扩展情况。
3.2.2 冲击速度对裂纹的影响
为更好的观察研究不同冲击速度下裂纹扩展情况，选取锤头着靶速度分别为10m/s、15m/s、20m/s、25m/s、30m/s五种不同速度进行数值模拟。得结果如***9所示，裂纹扩展与速度的关系如***10，11所示。
由***9可以明显的观察到，冲击速度不同靶板裂纹也不同，冲击速度越大，裂纹的密度随之增大，且主裂纹的延伸长度也随之增大。速度越大，起裂的位置离中心轴线也越近，即所受压缩和剪切作用越大。并且可以明显观察到速度越大，锤头下方裂纹扩展越密集，表现为破碎越彻底。由***10，***11可以看出，速度越大，起裂位置距中心轴线就越近，相应的主裂纹扩展垂直高度就越长。速度大于25m/s时主裂纹垂直高度扩展缓慢。速度大于15m/s时起裂位置向轴线中心靠拢缓慢。分析可知，存在一个最佳速度，使起裂位置更靠近中心轴线，主裂纹扩展垂直长达对大，且裂纹分布密度更均匀。
4 结束语
（1）针对花岗岩的特性，采用SPH法及JH-2损伤模型，在AUTODYN中使用二维轴对称法建立花岗岩靶板冲击破碎模型，形象的反映出裂纹扩展分布情况及不同速度对裂纹扩展的影响。
（2）起裂位置距离中心轴线的距离与速度有关，起裂位置决定了主裂纹的扩展长度，当起裂位置距离中心轴线越近，主裂纹扩展距离越长且裂纹的扩展分布密度越大。
（3）主裂纹扩展垂直高度，并不是随着冲击速度的增加而增加，当冲击速度达到一定时主裂纹扩展速度先增加后基本保持不变，此时微裂纹开始缓慢扩展并逐渐与主裂纹及其他微裂纹相互贯通。
参考文献
[1]席婧仪，辉，朱帝杰，等.岩石不等长裂纹应力强度因子及起裂规律研究[J].岩土工程学报，2015，04：727-733.
[2]秦飞，岑章志，杜庆华.多裂纹扩展分析的边界元方法[J].固体力学学报，2002，04：431-438.
[3]胡柳青，李夕兵，龚声武.冲击载荷作用下裂纹动态响应的数值模拟[J].爆炸与冲击，2006，03：214-221.
[4]王海兵，寿列枫，张建鑫，等.弹丸撞击下花岗岩靶破坏效应实验与数值分析[J].岩石力学与工程学报，2014，02：366-375.
[5]Motamedi D，Mohammadi S.Dynamic Crack Propagation Analysis of Orthotropic Media by the Extended Finite Element Method[J].International Journal of Fracture，2010，161（1）：21-39.
[6]Christophe Daux， Nicolas Mos， John Dolbow， Natara jan Sukumar， Ted Belytschko. Arbitrarybranched and intersecting cracks with the extended finite element method[J]. Int.J.Numer.Meth.Engng.，2000，4812.
[7]周小平，杨海清，董捷.压应力状态下多裂纹扩展过程数值模拟[J].岩土工程学报，2010，02：192-197.
[8]肖毅华.有限元法与光滑粒子法的耦合算法研究[D].湖南大学，2012.
[9]石少卿，汪敏，孙波，等.AUTODYN工程动力分析及应用实例[M].北京：中国建筑工业出版社，2012：16-19.
[10]王泽鹏，胡仁喜，康士廷，等.ANSYS13.0/LS-DYNA非线性有限元分析实例指导教程[M].北京：机械工业出版社，2011：352-362.
[11]李运禄，尹建平，王志***，等.基于SPH算法对混凝土裂纹扩展研究[J].兵器材料科学与工程，2015，05：19-22.
[12]王海兵，张海波，田宙，等.花岗岩介质中爆炸回填堵塞的数值模拟[J].岩石力学与工程学报，2005，24（增1）：5013-5019.
[13]杨震琦，庞宝君，王立闻，等.JH-2模型及其在Al2O3陶瓷低速撞击数值模拟中的应用[J].爆炸与冲击，2010，05：463-471.
[14]石祥超，陶祖文，孟英峰，等.致密砂岩Johnson-Holmquist损伤本构模型参数求取及验证[J].岩石力学与工程学报，2015，S2：3750-
3758.
作者简介：王石安（1990-），男，汉族，河南省三门峡市人，硕士研究生，主要研究领域为机械工程。
*通讯作者：曾海峰（1974-），男，汉族，***昌吉人，副教授，硕士研究生导师，主要研究领域为机械制造及自动化。
作者简介：张雨（1990-），男，汉族，湖北荆州人，硕士研究生，主要研究领域为机电一体化控制技术。
田佳（1992-），男，汉族，陕西蒲城人，硕士研究生，主要研究领域为自动化控制。
转载请注明出处学文网 » 基于SPH方法的花岗岩在冲击作用下裂纹扩展数值模拟