1、非连续变形分析及其在爆破数值模拟中的应用第 57 卷第 3 期有色金属(矿山部分)2005 年 5 月非连续变形分析及其在爆破数值模拟中的应用王根旺刘红岩刘晓琳(河南煤炭建设集团有限责任公司机械化工程公司,北京理工大学爆炸灾害预防控制国家重点实验室,中国兵器工业第二.八研究所)摘要:本文首先介绍了非连续变形方法的基本原理和求解平衡方程 ,然后对该方法的计算程序流程作了详细的介绍,分析了非连续变形方法在岩土工程计算中的特点.最后通过对一个台阶爆破过程的模拟说明了该数值分析方法在岩土爆破工程分析中的适用性.关键词:非连续变形分析程序流程图台阶爆破模拟非连续变形分析(discontinuousdef
2、ormationanalysis,简称 DDA)是石根华博士提出的一种新型的数值分析方法.DDA 方法是在块体运动学基础上,并部分吸收离散单元法关于接触形式的描述方法和刚体弹簧模型中位移函数的构造方法等优点而发展起来的一种与有限元法并行的数值分析方法.它是以被有限的不连续面切割而成的离散块体为基本单元,通过块体间的接触和作用在各块体上的位移约束条件,把若干个单独的块体连接成一个块体系统,并由系统最小势能原理建立总体平衡方程组,通过在块体接触面上施加或去掉刚性弹簧的办法,使块体间不存在侵入和张拉现象,由此也使边界条件得以满足.其刚体位移和块体变形采用统一的有限元格式求解,不仅允许块体本身有位移和
3、变形,而且还允许块体间有滑动,转动,张开等运动形式,能够得到大变形和大位移解.不连续变形独特的性质是:完备的运动学理论及数值实现;全一阶(或高阶) 多项式位移模式 ;严格的平衡假定;正确的能量消耗和高效率的数值计算.因而,非连续变形分析方法自提出以来受到了国内外的广泛关注和深入研究“.本文就在介绍非连续变形分析方法基本原理的基础上,对其计算程序的流程及结构进行了分析,最后通过对一个工程实例的计算说明了该程序在岩土工程爆破数值模拟计算中的有效性和适用性.成的离散块体为基本单元,通过块体间的接触和作用在各块体上的位移约束条件,把若干个单独的块体连接成一个块体系统,并由系统最小势能原理建立总体平衡方
4、程组.如果待分析的块体系统中有 n个块体,则这个块体系统的平衡方程为:KllKl2Kl3K2lK22K23KlKl12KrI3K厂 DK2nlID:;II;KJLDFlF2:FDDA 是对每一载荷或时间按增量方法分时步进行计算的,大位移和大变形是由小位移和小变形累加而成的,在每一步中,所有点的位移都是小的,因而在每一时间步中都满足小位移和小变形条件.对块体之间接触问题的处理是 DDA 方法得以实施的关键也是 DDA 方法的精华部分.在 DDA 方法中,它允许块体间发生滑动和界面间的张开或闭合.从总体上来讲,DDA 中的运动学条件主要包括:接触准则,进入线定义,嵌入判断,嵌入防止及进入位置等.利
5、用 DDA 方法在计算块体运动和变形时,它要求块体之间必须满足无拉伸和无嵌入的条件.这是通过在块体间添加或移去接触弹簧来实现的,如果块体在接触处发生了侵入,则施加刚度很大的弹簧将其沿原路推回;如果两块体间有接触拉力,则应移去两块体间的弹簧,所以这一过程也是通过反复的开合迭代来实现的.1DDA 方法的基本原理非连续变形分析是以被有限的不连续面切割而2DDA.程序计算流程王根旺助理工程师河南郑州 450009DDA 分析方法也是一种数值分析方法,它象有限元一样最终也是通过程序来解决问题的,所以能第 3 期王根旺等:非连续变形分析及其在爆破数值模拟中的应用 41够正确地理解和使用该程序是利用这种数值
6、方法分析问题的关键.该方法和其他数值分析方法的基本步骤是一样的,也是首先输入计算的几何参数和物理参数,绘出所要分析对象的几何图形,然后施加相应的约束和载荷进行求解,最后对处理结果进行相应的分析.所不同的是,由于该方法的计算程序主要是针对块体系统进行分析的,所以搜索块体之间的接触和运动情况成为该程序的一个主要内容,而且也是该算法的精华部分.所以从图 1 所示的计算程序流程中可以看出,只有当每个加载步或时间步内的循环次数大于 6 次且最大位移比小于或等于 3时,这一当前步的计算才算完成,否则均是在该时间步内的循环.对每一个时间步都经过这样的循环,当总的计算时间步数循环结束后,该计算程序才算计算最终
7、完成.所以同相应的有限元计算程序相比,计算量要至少增加 6 倍以上,对机时的耗费较大.图 1DDA 程序计算流程图其中:nn=1,n5 表示循环从第一步一直到最后一步全部完成.从该程序的计算流程图还可以看出,该程序求解问题的方法和有限元相类似,也是首先对单个块体进行分析,然后组集成总体刚度矩阵,只是在求解方程组的时候,它采用的是非零存储法,这比有限元程序采用的方程组解法要有更高的效率.3 算例3.1 计算模型本文中所建立的计算模型是在一二维分析域内的一圆形装药台阶爆破计算模型.在该半无限域内的岩体被三组节理所切割而形成块体.计算模型的尺寸为如图 2 所示,图中的圆形代表平面内的圆形装药.本计算
8、模型采用理想的装药结构,即不考虑装药堵塞部分与整体部分的差异,而把整个台阶内的岩体作为性质完全相同的同一种材料,所有的节理的性质也完全相同.YX图 2 台阶爆破计算模型当炸药在岩石中爆炸后,产生的爆炸载荷是非常复杂的,为了能够利用 DDA 对该爆破过程进行模拟,先要必须对其进行合理的简化.假设爆炸载荷以压力形式的均布载荷方向垂直作用于炮孔壁上,根据爆炸载荷的实际特点,将载荷曲线简化为三角波,其最大的峰值压力为 50MPa,升压时间为 80s,卸载时间为 300s.在本文中,为了方便载荷的施加,专门开发了载荷施加的前处理程序,一方面是为了使载荷的施加简单化,另一方面也是为了保证载荷施加的准确性.
9、3.2 计算参数的选取计算参数包括岩石的材料参数,节理面参数和计算时步参数.在爆炸载荷作用下,材料的力学性能表现出很大的动态性能,这里应取岩石材料的动态参数:动态弹性模量 E=6.010“Pa,动泊松比=0.20,密度 p=2.710kg/m;岩体内节理面42 有色金属(矿山部分) 第 57 卷的摩擦角 0=60.,粘结力 C=1.610Pa;计算时步参数为:总计算时步为 700 步,时间间隔由系统自动选取.3.3DDA 模拟结果分析根据以上所建立的力学模型,采用 DDA 方法对台阶爆破的全过程进行了模拟.图 3(a),(b),(C)分别为在 420 时步,580 时步和 700 时步时爆破过
10、程中的 3 个瞬时状态.一一(a)420 时步(b)580 时步(c)700 时步图 3 爆破过程的模拟示意图从模拟结果可以看出,爆破漏斗的形成过程基间隔的确定等,都无严格的理论选择依据,而需要根本上符合试验和理论分析结果,即在炮孑 L 中产生的据经验进行选取 .同时从本文中的算例中还可以看爆炸冲击波作用于炮孑 L 周围的岩石上,推动岩石向出 ,DDA 还不能用来模拟单元块体在爆炸载荷作用着约束力比较小的方向发展,也即是沿着最小抵抗下的破坏过程,而只能把被节理裂隙切割而成的岩线的方向发展.对于台阶爆破来说,存在着上方和块当作一个最小的不可分隔的单元,仅能发生整体前方两个自由面,岩石也即沿着这两
11、个自由面的方变形而不能被再次破坏,这对于爆破中远区来说还向向外运动.这也说明了自由面对岩石破坏程度及是比较符合的,但是对于爆破近区而言,这一点却是破坏方式的影响规律.不尽符合实际,这也是现有 DDA 方法在模拟爆炸作同时还应该指出,尽管这种数值分析方法很适合用方面的一个不尽完美的地方.不过相信随着这种于用来模拟块体在外力作用下的破坏过程,但是它也方法在理论上的不断完善和在实践中的不断推广应只是提供了一个大致的规律性的结果,如果需要深人用,这些不足将会最终得到改进和完善.还应该综合考虑与理论分析和试验研究相结合.4 结语本文在介绍 DDA 基本理论和其程序结构的基础上对二维平面问题,台阶爆破过程
12、进行了模拟,成功地再现了台阶爆破条件下圆形药包爆破破岩的过程.这充分地说明了 DDA 这种非连续性的数值分析方法在模拟被节理所分割的岩体在爆炸等外载荷作用下的运动规律所具有的独特优越性.然而,尽管 DDA 在非连续介质的计算中有很多的优点,但是由于其出现的时间还不是很长,其在很多方面还不是很完善,如弹簧罚值的选取,合理时间参考文献lShiG.H.andGoodmnR.E.DiscontinuousDeformationAnalysis.Proceedingsofthe25U.S.SymposiumofRockMechanics.Evanston.USA,1984:2692772 石根华.着,裴
13、觉民译.数值流形方法与非连续变形分析M.北京:清华大学出版社,19973 刘军,李仲奎.非连续变形分析方法研究现状及发展趋势J.岩石力学与工程,2003,23(5):8398454 刘红岩,杨军,陈鹏万.爆破漏斗形成过程的 DDA 模拟分析J.工程爆破,2004,10(2):17 2O5 戴晨,朱传云,舒大强,等.DDA 及其在爆破过程仿真中的应用J.爆破,2001,18(增刊):46U(上接第 2l 页)的产状,也是吻合的,早期矿脉不与变质砂岩沟通,极少部分矿脉,赋存在侵人体附近的外接触带中.矿体内广泛发育热液蚀变,围岩蚀变类型与规模直接或间接地指示矿化作用的强度.从白钨矿,辉钼矿,绿柱石,黑钨矿等矿物组合推测,矿床内主要金属矿物形成于高温阶段,气成与中温阶段仅是矿化的前奏与延续.因此认为本区细脉型钨矿床应是岩浆期后高温热液矿床.参考文献1 福建省地质局三 0 六地质队编制.福建清流行洛坑钨钼矿储量报告.1996,3口