刘志林, 孙巍巍, 王晓鸣
(1.南京理工大学 智能弹药技术国防重点学科实验室,南京 210094; 2.南京理工大学 理学院土木工程系,南京 210094)
基于颗粒流离散元模型的弹丸侵彻细观混凝土数值模拟方法研究
刘志林1, 孙巍巍2, 王晓鸣1
(1.南京理工大学 智能弹药技术国防重点学科实验室,南京210094; 2.南京理工大学 理学院土木工程系,南京210094)
摘要:利用颗粒离散单元法,研究弹丸侵彻细观混凝土模型中弹丸受到介质的阻应力与侵彻速度的关系。采用蒙特卡罗法随机生成并投放混凝土骨料且骨料的粒径分布满足级配曲线。通过对混凝土颗粒离散元细观力学模型进行单轴压缩实验、巴西劈裂实验和双轴压缩实验的参数反演,确定细观模型参数,能使细观混凝土模型具有和一般混凝土等效的力学性能。分析了骨料、过渡层和砂浆三相材料各细观参数对混凝土单轴压缩应力应变关系影响以及锥形弹和平头弹弹丸直径对侵彻阻应力的影响。将颗粒离散元细观力学模型方法计算的弹丸阻应力与空腔膨胀理论计算模型相比较,表明计算离散元方法具有良好的精度和实用性。
关键词:离散元;颗粒流;细观模型;混凝土;侵彻
混凝土是国民经济建设和国防建设中广泛使用的建筑材料,其可被看作一种由碎石、沙子、水泥、添加剂等组成的复合材料,其力学性能和内部结构非常复杂。以混凝土为目标靶的侵彻问题是研究钻地弹的主要问题。弹丸侵彻混凝土是一个伴随着混凝土的大变形、高应变率和高压的过程,响应十分复杂。因此,找到准确描述混凝土介质力学性能的模型和算法显得尤为重要。
颗粒离散元(Particle Flow Code)是Cundall[1]基于分子动力学原理提出的通过圆形颗粒单元运动与相互作用来模拟颗粒材料力学性能的一种非连续介质力学数值方法,颗粒单元的运动与相互关系满足牛顿第二定律和力-位移定律,材料的宏观力学性能主要由细观颗粒体的强度、颗粒黏结强度、颗粒尺寸、形状等影响。由于混凝土是一种典型的非均质各向异性材料,较之连续介质力学计算方法,颗粒离散元在模拟混凝土材料细观破坏过程方面有着独特的优势。Ghazvinian 等[2]用颗粒离散元研究了平面非持久性张开节理的失效机制;徐文杰等[3-4]基于数字图像建模技术,研究了土石混合体等非均质材料的细观力学特征;Qin等[5]用颗粒离散元在细观尺度研究了三相混凝土(石子、砂浆和过渡层)的动态失效行为,对单轴压缩和动态劈拉实验的动态增强因子进行了验证模拟,取得较好结果。混凝土的基本的力学参数(如:单轴压缩强度,拉伸强度,剪切强度)是弹丸侵彻混凝土研究中的重要参数,也是弹丸侵彻混凝土的阻力研究中的空腔膨胀理论方法中的重要参数[6],在弹丸侵彻混凝土数值计算研究中,需对数值计算中混凝土的模型的基本力学进行必要的标定。
本文主要利用颗粒离散元数值计算方法,建立考虑骨料、过渡层和砂浆的三相细观混凝土模型,通过对混凝土单轴压缩实验、双轴压缩实验、巴西实验的模拟,确定颗粒离散元细观力学参数,并在此基础上研究了弹丸直径对阻应力的影响。
1混凝土颗粒离散元细观力学模型
混凝土细观尺度的颗粒离散元模型的建立主要分两个步骤:① 随机骨料的生成和投放;由于混凝土中骨料分布的随机性,本研究采用蒙特卡罗法,确定生成骨料的直径与投放位置,生成的骨料不重叠,循环迭代至骨料含量达到设计要求,停止生成骨料,保存骨料的直径和位置信息。② 细观成分(骨料、砂浆和界面)赋予细观力学属性。将生成的骨料信息投射到颗粒流离散元背景网格上,本文通过文献[2]中颗粒排布规律生成背景网格。在被骨料投射到的单元赋予骨料细观力学属性,骨料外一层单元赋予过渡层细观力学属性,其它单元则赋予砂浆离散元属性,如图1所示。
骨料的粒径分布可以用级配来表示,Walraven[6]基于Fuller公式,将三维级配曲线转化为试件内二维截面上任意一点具有D (1) 式中:D0为筛孔直径,Pk为骨料占混凝土总体积分数,Dmax为最大骨料直径。 图1 单个骨料投射到背景网格示意图Fig.1 Schematic of a single aggregate projected onto the background grid 本研究中粗骨料粒径为5~30 mm,素混凝土的各成分的配比:水泥∶砂∶石子∶水=1∶2.65∶6.2∶0.61。骨料体积分数可近似为40%,根据Walraven提出的式(1)计算级配,各组粒径骨料的含量如表1,根据上述方法建立的细观离散元模型如图2所示。 表1 骨料不同粒径的体积分数 图2 骨料随机投放的细观离散元模型Fig.2 Meso-scale discrete element model of aggregate random delivery 2混凝土力学性能离散元细观参数研究 目前在颗粒离散元计算中并没有完善的理论可以直接从微观特性来预见宏观特性,为使模拟结果与实测结果相吻合,需要反复调整微观参数。本研究通过对细观模型进行单轴压缩、双轴压缩、巴西劈裂实验标定符合要求的细观参数,确定细观模型能很好的描述C45等级的混凝土的力学性能,为研究侵彻阻力的模拟研究提供可信的细观模型参数。经过大量的反复试算,最终模型参数见表2,表中kn、ks、σn、σs、μ、d分别为颗粒法向刚度、切向刚度、方向黏结强度、切向黏结强度以及摩擦因数。 表2 模型参数 2.1颗粒流离散元计算原理 颗粒离散元是利用显示差分算法和离散元理论开发的微/细观力学程序,其通过离散单元法来模拟圆形颗粒介质的运动及其相互作用。计算模型中的所有颗粒都假设成刚体,颗粒之间的接触力的大小与颗粒间接触重叠的大小相关,颗粒之间的运动与接触力的计算主要通过反复迭代牛顿第二运动定律和力-位移定律进行,直至满足预定条件,停止计算。 力-位移定律采用线弹性接触本构模型,接触力的计算公式如下: (2) (3) (4) 用接触黏结模型将整个模型的颗粒单元在接触处黏结起来,黏结强度由法向黏结强度和切向黏结强度组成。当接触黏结存在时颗粒不会发生滑动,即有恒定法向刚度与切向刚的弹簧作用于颗粒接触位置。当颗粒间重合量小于零时,颗粒间存在法向拉力,如果法向接触力大于等于法向黏结强度时,法向黏结消失,颗粒法向方向将不再承受拉力;当切向接触力达到切向黏结强度时,切向黏结消失,切向运动满足滑移模型。滑移发生的条件为: (5) (6) 具体计算模型如图3所示。 图3 接触黏结模型和滑移模型[7]Fig.3 Contact bond model and slip model 2.2混凝土单轴压缩模拟 混凝土单轴压缩加载条件下,将会出现受力后变形、内部微裂纹的发展、损伤积累、达到强度极限和最终破坏等一系列变化过程,其应力应变关系是研究和分析混凝土结构承载力和变形的主要依据,在弹丸侵彻混凝土研究中,混凝土的单轴压缩强度是一个非常重要的靶体指标。 本研究采用线性刚度模型、接触黏结模型与滑移摩擦模型进行混凝土的细观离散元数值仿真,需要确定的细观参数主要有:kn、ks、σn、σs、μ。混凝土单轴抗压强度实验选用150×150×150标准试件进行模拟,细观模型如图4所示。 对试件的恒应变率加载方式是通过对试件两端的刚性墙体Wall施加速度载荷来实现的,不同的墙体速度对应不同的加载应变率,为了降低骨料随机分布带来的离散性对结果的影响,需对图2所示的三种模型进行平行压缩试验,取三者结果的平均值作为实验结果,三种模型应力应变见图4。三种模型应力应变曲线在峰值前接近重合,峰值变化跳动不大(见表3),软化段整体的趋势完全一致。 图4 三种随机离散模型的计算结果Fig.4 The results of three randomized discrete model 图5 本文计算结果与经典Hognestad模型比较 Fig.5 The results compare this paper with the classical model of Hognestad 三个模型对应的单轴压缩强度的峰值应力以及峰值应变列入表3中。三种随机投放模型模拟峰值应力的平均值为49.1 MPa,峰值应变的平均值为2.02×10-3。如图5所示,细观离散元模型得到的应力应变曲线与经典美国Hognestad[8]模型相吻合。图6为单轴压缩仿真模型的加载方式示意图,图7与图8显示了试件加载后的破坏情况,单轴压缩条件下,立方体试件加载后出现贯穿裂缝,裂纹均沿骨料与骨料间薄弱的过渡层发展至试件端面,骨料间破坏的砂浆单元多是法向黏结破坏,少数的切向黏结破坏多是发生在过渡层单元上。单轴压缩条件下,单元的法向黏结破坏数明显多于切向黏结破坏数。 表3 三种随机投放模型单轴压缩模拟结果 图6 加载方式Fig.6Loadingmode图7 压缩破坏图Fig.7Specimenfailureundercompression图8 压缩黏结破坏图Fig.8Bondfailureundercompression 本研究为了能够了解混凝土中骨料、过渡层和砂浆各细观参数对混凝土整体性能的影响,对骨料、过渡层和砂浆的刚度、强度和摩擦细观参数做了分析计算。对刚度和强度参数分别作了放大十倍和缩小十倍的对比,结果见图9~图11。比较可知: (1) 三相材料(骨料、过渡层和砂浆)的刚度是由kn和ks控制,整体刚度随着各三相材料各材料的刚度的增加而增加,刚度值的改变会影响整体强度值。 (2) 混凝土中骨料的强度值最大,过渡层的强度值最小。骨料的强度增加时,将不会对整体的强度峰值产生影响,对软化段的性能会产生一定影响(增加了断裂能)。减小至一定强度时峰值则会减小;过渡层的强度减小时,将不会对整体的强度峰值产生影响,增加时峰值则会增加;砂浆的强度值的增加或减小,整体强度值则会增加或减小。三相材料的强度值的改变都不会改变整体的刚度。本研究对骨料、过渡层和砂浆强度的影响,只针对本研究中的细观参数放大或缩小十倍情况下的定性研究,混凝土整体强度与骨料、过渡层和砂浆强度的关系的定量关系还需要进一步深入研究。 (3) 骨料、过渡层和砂浆的摩擦因数的改变可以改变整体强度峰值,而不影响整体刚度。对于强度最大的骨料,摩擦因数的增加对强度的贡献有限;相比于强度最低的过渡层,过渡层摩擦因数的增加对整体强度的增加幅度明显高于骨料和砂浆。 图9 骨料的刚度、强度以及摩擦因数的影响Fig.9 The influence of the stiffness, strength and friction parameters of aggregates 图10 过渡层的刚度、强度以及摩擦因数的影响Fig.10 The influence of the stiffness, strength and friction parameters of interfacial transition layer 图11 砂浆的刚度、强度以及摩擦因数的影响Fig.11 The influence of the stiffness, strength and friction parameters of mortar 2.3混凝土双轴抗压模拟 在弹丸侵彻混凝土介质过程中,混凝土在其与弹丸作用近区[9]处于高静水压力的复杂受力环境,混凝土在高静水压力下的特性对侵彻混凝土介质问题的研究十分必要,本研究模拟了混凝土模型在双轴压缩条件下的围压对强度的影响。 在混凝土单轴抗压模拟模型的基础上,在细观离散元模型试件两侧加上由伺服器控制的墙体Wall,提供恒定的侧向压力,计算结果如图12所示,围压增大,峰值应力和峰值应变随着增大,对应值见表4。根据表4中测压与峰值应力绘制莫尔圆包络线,可以得到静水压力与剪切应力关系,线性拟合可以得到混凝土基本材料参数λ=1.466,τ0=13.04 MPa(见图13和图14)。 表4 不同围压下应力应变曲线的应力峰值及其对应的应变 图12 不同围压下的轴向应力应变曲线Fig.12 Axial stress-strain curves under different confining pressures 图13 莫尔圆包络线Fig.13 Mohr envelope 图14 静水压力与剪应力关系Fig.14 The relationship of Hydrostatic pressure and shear stress 2.4混凝土巴西劈拉实验模拟 混凝土介质的抗拉强度远小于其抗压性能,很多混凝土材料结构发生破坏都是由于抗拉强度不足引起的,混凝土的抗拉强度在力学性能指标中十分重要。秦川[7]在细观层次上分析和研究了混凝土在霍普金森杆上冲击劈拉实验,采用基于背景网格的混凝土细观力学预处理方法[8]和细观参数反演技术[10],建立颗粒离散元细观力学模型,对混凝土冲击劈拉实验进行了细观数值模拟,仿真结果与实验吻合较好。 与文献[11]动态冲击研究不同,本文主要研究准静态条件下试件的劈拉性能。巴西劈拉模拟的加载方式如图15,在试件上下端面施加刚性墙体Wall,墙体Wall速度恒定直至试件破坏停止加载,墙体上的受力载荷的历史信息将存储在history文件中,巴西劈拉试件直径为150 mm。当墙体恒定加载速度为0.015 m/s时,三个模型计算结果如图16,平均最大破坏载荷为1.411×106N。根据巴西试样拉伸强度计算式(7)可知试件的拉伸强度为5.99 MPa。巴西劈裂试件破坏图见图17和图18,破坏模式与单轴压缩类似,裂纹组要沿薄弱的过渡层发展,直至形成贯穿裂纹发展到试件端面。 (7) 图15 劈裂实验细观离散元模型及加载方式Fig.15 The meso-scale discrete element model and Loading mode of split tests 图16 位移载荷曲线Fig.16 Load-deformation curves 图17 黏结破坏分布图(浅色破坏为法向黏结破坏,深色为切向破坏)Fig.17 Bond failure distribution (light colour:normal bond failure,dark colour:tangential bond failure) 图18 试件破坏图Fig.18 Specimen failure under split 模型最大破坏载荷/N最大破坏位移/11.416×1061.99×10-421.441×1062.22×10-431.376×1062.18×10-4平均值1.411×1062.13×10-4 3侵彻阻力仿真计算 弹丸侵彻混凝土介质的研究中,混凝土常被当作均质各项同性体来处理,弹丸在侵彻过程中受到的混凝土介质的阻力也是在此基础上计算而得,比较著名的侵彻阻力理论计算模型如空腔膨胀理论、微分面力法、Amini-Anderson模型等,无一例外都是没有考虑到混凝土材料各向异性的非均质材料特性。国内有部分研究是通过有限元软件建立骨料砂浆以及过渡层的三相有限元模型来模拟计算混凝土介质力学性质以及侵彻性能,但其材料模型都没有脱离宏观本构,很难全面准确的描述混凝土断裂损伤性能。与连续介质力学方法不同的是,颗粒流计算方法试图从微观角度研究介质的力学性能和行为,适合求解大位移和非线性问题。 以上基本力学实验模拟验证了模型各相细观参数能很好的描述C45混凝土的力学性能,通过建立弹丸侵彻混凝土靶体的细观模型对侵彻阻力特性进行研究。为了研究弹丸侵彻速度与侵彻阻力之间的关系,由墙体Wall以恒定的速度侵彻混凝土,侵彻速度分别为400 m/s、600 m/s、800 m/s、1 000 m/s、1 200 m/s、1 400 m/s、1 600 m/s。混凝土靶体细观模型如图19(a)所示,靶体宽0.9 m,厚度为1 m。图中黑色颗粒为骨料,骨料间的空白为砂浆。墙体侵彻靶体过程中,靶体破坏情况见图19(c),靶体表面撞击点位置有靶体崩落,形成开坑区,墙体侵彻过后形成隧道区,隧道近区骨料破碎,在墙体前端有贯穿裂缝,延伸至靶体侧面。力链分布不连续,在裂纹处间断。黏结破坏分布见图19(b),弹丸隧道近区的黏结破坏密度最大,切向黏结破坏主要发生在近区附近,在远离近区的破坏黏结以法向黏结破坏为主。 (黑色颗粒为骨料) (黑色为切向黏结破坏)图19 细观离散元混凝土靶及其侵彻响应图Fig.19 Meso-scale discrete element concrete target and penetration response 为了消除自由表面对阻力特性的影响,取墙体Wall侵彻深度大于2D后的弹丸阻力的平均值为此速度对应的侵彻阻力,D为弹丸直径。为使墙体Wall侵彻行程相同,即在靶体中撞击到的骨料砂浆和过渡层相同,不同侵彻速度侵彻条件下侵彻时间将随着速度的增大而减小。模拟结果见图20。随着侵彻速度的增加,侵彻阻力随着增加,且震荡幅值也在增大。取其平均值作为对应速度下侵彻阻力,将阻应力与空腔膨胀理论计算值作对比。图中空腔膨胀模型依次来自文献[13-16],比较结果显示,本文计算阻应力在1 000 m/s以下,阻力值在现有空腔膨胀模型计算范围内;当速度大于1 000 m/s时,本文计算阻应力比现有空腔膨胀阻应力小。在实际情况下,多数研究均发现在速度1 000 m/s以上时,现有的空腔膨胀模型计算阻力都偏大,文献[15]中几种空腔膨胀理论计算的侵彻深度值比实际实验侵彻深度值小,说明理论计算的阻力偏大。在速度大于1 000 m/s后偏离的越明显。可见本文离散元计算阻力趋势符合现实情况。 图20 不同侵彻速度下轴向阻力10-4s内的时程曲线Fig.20 Curves of axial penetration resistance under different velocity within 10-4s 图21 本文计算阻应力与经典空腔膨胀理论模型比较Fig.21 The resistance stresses calculated of this article compared with the classical model of cavity expansion theory 文献[14]中提到,弹丸侵彻阻应力可能与弹丸直径有关。本研究做了锥形弹和平头弹不同弹径不同速度下侵彻阻力的计算模拟,计算结果如图22。图22(a)中显示了锥形弹丸D/d在4~32范围内对应速度为400、800和1 600 m/s时的侵彻阻应力,结果显示三个速度下,随着D/d的增大,阻应力逐渐减小,减小的趋势逐渐平缓。平头弹规律与锥形弹类似,减小的趋势在D/d=10左右就趋近于平缓。对于速度为1 600 m/s,D/d=32对应的阻应力相对较低的原因可能为自由面边界的影响已经显现,故不能在此尺寸的靶上再做D/d>32侵彻阻力研究。 (a) 锥形弹 (b) 平头弹图22 弹丸直径与颗粒直径比对阻应力的影响Fig.22 The influence of projectile diameter and particle diameter ratio to resistance stress 4结论 本研究主要利用颗粒流离散元计算方法,建立含骨料、砂浆和过渡层的混凝土细观力学模型, 对模型进行单轴压缩、巴西劈裂和双轴实验的参数反演,确定混凝土细观模型参数,为侵彻阻应力研究提供可信材料细观参数。对侵彻阻应力研究发现: (1) 本研究计算的弹丸侵彻阻应力与空腔模型理论计算比较发现,速度在1 000 m/s以下时,本文计算阻应力在现有空腔膨胀模型计算范围内;速度大于1 000 m/s时,本文计算阻应力则比现有空腔膨胀理论计算小。 (2) 弹丸直径对侵彻阻应力有影响,在D/d>4范围内,锥形弹丸阻应力随着D/d的增大逐渐减小,趋势逐渐趋于平缓;平头弹在D/d=10以后变化很小,即弹丸阻应力在弹径某一阈值后将不会随着弹丸直径变化,小于直径阈值时,弹丸阻应力随着弹径的增大而减小。 参 考 文 献 [ 1 ] Cundall P A, Strack O D L. A discrete numerical model for graunlar assemblies [J]. Geotechnique,1979, 29(1): 47-65. [ 2 ] Ghazvinian A, Sarfarazi V, Schubert W, et al. A study of the failure mechanism of planar non-persistent open joints using PFC2D[J]. Rock Mechanics Rock Engineering, 2012,45:677-693. [ 3 ] 徐文杰,胡瑞林,王艳萍. 基于数字图像的非均质岩土材料细观结构PFC2D模型[J]. 煤炭学报,2007,32(4),385-389. XU Wen-jie, HU Rui-lin, WANG Yan-ping. PFC2D model form esostructure of inhomogeneous geomaterial based on digital image processing[J].Journal of China Coal Society, 2007,32(4):385-389. [ 4 ] 丁秀丽,李耀旭,王新.基于数字图像的土石混合体力学性质的颗粒流模拟[J].岩石力学与工程学报,2010(3):477-484. DING Xiu-li, LI Yao-xu, WANG Xin. Particle flow modeling mechanical properties of soil and rock mixtures based on digital image[J]. Chinese Journal of Rock Mechanics and Engineering, 2010(3):477-484. [ 5 ] Qin Chuan, Zhang Chu-han. Numerical study of dynamic behavior of concrete by meso-scale particle element modeling[J]. International Journal of Impact Engineering, 2011(38):1011-1021. [ 6 ] Walraven J C, Reinhardt H W. Theory and experiments on the mechanical behavior of cracks in plain and reinforced concrete subjected to shear loading[J].Heron,1991,26(1A):26-35. [ 7 ] Itasca. PFC2D-Particle flow code in 2-dimensions,version 4.0,user’s manual [M].Minneapolis, MN:Itasca Consulting Group,Inc.,2008. [ 8 ] 东南大学,同济大学,天津大学.混凝土结构设计原理[M].北京:中国建筑工业出版社,2011. [ 9 ] 葛涛,王明洋. 坚硬岩石在强冲击载荷作用近区的性状研究[J].爆炸与冲击,2007,27(4):306-310. GE Tao, WANG Ming-yang. Character near strong impact loading zone in hard rock [J]. Explosion and Shock Waves,2007,27(4):0306-0311. [10] 秦川,武明鑫,张楚汉.混凝土冲击劈拉实验与细观离散元数值仿真[J].水力发电学报,2013,32(1):196-205. QIN Chuan, WU Ming-xin, ZHANG Chu-han. Impact splitting tensile experiments of concrete and numerical modeling by meso-scale discrete elements[J]. Journal of Hydroelectric Engineering, 2013,32(1):196-205. [11] 秦川,郭长青,张楚汉. 基于背景网格的混凝土细观力学预处理方法[J]. 水利学报,2011,42(8): 941-948. QIN Chuan,GUO Chang-qing,ZHANG Chu-han. A pre-processing scheme based on background grid approach for meso-concrete mechanics [J]. Journal of Hydraulic Engineering,2011,42(8): 941-948. [12] 秦川. 基于细观力学的混凝土动态特性研究[D]. 北京: 清华大学,2011. [13] Luk V K, Forrestal M J. An empirical equation for penetration depth of ogive-nose projectile into concrete targets[J]. International Journal of Impact Engineering 1994,15(4): 395-405. [14] Forrestal M J,Tzou D Y. A spherical cavity-expansion penetration model for concrete targets [J]. International Journal of Solids and Structures, 1997, 34:4127-4146. [15] 李志康,黄风雷.考虑混凝土孔隙压实效应的球形空腔膨胀理论[J].岩土力学,2010(5):1481-1485. LI Zhi-kang, HUANG Feng-lei. A spherical cavity expansion theory of concrete considering voids compacted effects[J]. Rock and Soil Mechanics, 2010(5):1481-1485. [16] 王明洋,郑大亮,钱七虎. 弹体对混凝土介质侵彻、贯穿的比例换算关系[J].爆炸与冲击, 2004, 24(2):108-115. WANG Ming-yang, ZHENG Da-liang, QIAN Qi-hu. The scaling problems of penetration and perforationfor projectile into concrete media[J].Explosion and Shock Wave, 2004, 24(2):108-115. Numerical simulation for projectile penetrating meso-scale concrete based on particle flow discrete element model LIUZhi-lin1,SUNWei-wei2,WANGXiao-ming1 (1. National Defense Key Laboratory of ZNDY, Nanjing University of Science and Technology, Nanjing 210094, China;2. Civil Engineering Department, Nanjing University of Science and Technology, Nanjing 210094, China) Abstract:The relationship between penetration resistance and velocity of penetration in a projectile penetrating meso-scale concrete model was studied here using the particle flow discrete element method. Concrete aggregates whose size distribution satisfied a grading curve were randomly generated and put in using Monte Carlo method. In order to get the same effective mechanical properties as those of the macro-scale model, the parameters of the micro-scale model were determined by applying the parameter inversion technique in uniaxial compression, splitting tensile, and biaxial compression tests using the particle flow discrete element micromechanical model. The influences of the aggregate, transition layer and mortar 3-phase material microscopic parameters on concrete uniaxial compression stress-strain relation and the effects of diameter of projectiles with flat and conical nose on penetration resistance stress were analyzed. The comparison of the numerical results of resistance stress calculated with the particle flow discrete element micromechanical model and the analytical results based on the cavity expansion theory showed that the discrete element model has good applicability and accuracy. Key words:discrete element; particle flow; meso-scale model; concrete; penetration 中图分类号:O385 文献标志码:A DOI:10.13465/j.cnki.jvs.2016.04.026 通信作者孙巍巍 男,副教授,1978年生 收稿日期:2014-11-10修改稿收到日期:2015-03-16 基金项目:国家重点基础研究发展计划(973);国家自然科学基金(51308297) 第一作者 刘志林 男,博士生,1988年生 E-mail:sww717@163.com