（地质工程专业论文）大断层附近煤层开采的防水煤柱留设离散元分析.pdf

大断层附近煤层开采的防水煤柱留设离散元分析 摘要：讨论了离散元法的理论基础及离散元模型分类与特点。针对煤矿开 采中的采空区围岩冒落、破裂、变形现象，用u d e c 4 0 软件对桃园矿煤层 采动后f 2 断层附近位移变化量和应力分布进行数值计算，并利用其模拟结 果对f 2 断层煤柱留设宽度适合性进行了深入研究。主要内容有： 分析了煤系地层与相关含水层的含水隔水性质，对f 2 断层带的物 质组成与物理力学性质进行了评价分析。 采用规程中公式，计算出对应于不同安全系数的1 0 煤与8 l 煤的断层防 水煤柱留设宽度。采用基于m o h r - - c o u l u m b 本构连续介质与j o i n ta r e a c o n t a c t - - c o u l o m bs l i pw i t hr e s i d u a ls t r e n g t h 节理介质的离散元本构模型， 针对桃园煤矿地层及f 2 断层实际，模拟了1 0 煤和8 1 煤各种煤柱留设条件 下的煤层采动围岩变形破坏特性，研究了断层防水煤柱宽度和采厚的非线 性关系，给出了l0 煤和8 i 煤煤柱宽度的合理值。 对断层防水煤柱公式实用中的缺陷进行了深入分析，指出公式的力 学涵义与应用特点，指出公式中煤柱宽度与煤层采厚线性关系的不合理 性，以及与“三下 理论及裂隙发育实际规律间的矛盾性。 对岩层移动角引入防水煤柱留设计算的适用条件与缺陷也进行了 探讨。 关键词：离散元法：断层；防水煤柱；采动裂隙 d i s c r e t ee l e m e n ta n a l y s i so nt h ew a t e r p r o o fp i l l a rt o c o a ls e a mm i n i n gg o a fb e s i d el a r g ef a u l t a b s t r a c t n o to n l yt h et h e o r e t i c a lb a s i so fd i s c r e t ee l e m e n tm e t h o d ( d e m ) b u ta l s o t h ec l a s s i f i c a t i o na n dc h a r a c t e r i s t i c so fd e mw e r ed i s c u s s e d f o rc o a lm i n i n gi n t h em i n e d o u ta r e as u r r o u n d i n gr o c kc a v i n g ，c r a c k i n g ，d e f o r m a t i o np h e n o m e n a ， u d e c 4 0s o f t w a r ew a su s e di nt h en u m e r i c a lc a l c u l a t i o nt ot h ed i s t r i b u t i o no f d i s p l a c e m e n t sa n ds t r e s s e ss u r r o u d i n gm i n i n gg o a v e so fs e a m1o 撑a n ds e a m 8 1a b o u tf a u l tf 2i nt a o y u a nc o a lm i n e ，a n dp r o p e rw i d t ho fw a t e r p r o o fp i l l a rt o t h ef a u l t m a i nc o n t e n t sa r e ： w a t e r b e a r i n ga n dw a t e r r e s i s t i n gf e a t u r e so fc o a l b e a r i n gb e d sa n d r e l a t e ds t r a t aw e r ea n a l y z e d e v a l u a t i o na n da n a l y s i sa b o u tt h em a t e r i a l c o m p o s i t i o na n dp h y s i c sm e c h a n i c so ft h ef a u l tf 2w e r ea n a l y z e d ( 室) o r d e rf o r m u l aw a su s e dt oc a l c u l a t et h ew i d t ho fw a t e r p r o o fp i l l a ro f s e a m s10 拌c o a la n d81 # c o r r e s p o n d i n gt od i f f e r e n ts a f e t yf a c t o r b a s e do n t e c t o n i cm o d e lo fd e ma b o u tm o h r c o u l u m bc o n s t i t u t i v ec o n t i n u o u sm e d i u m a n dj o i n ta r e ac o n t a c t - c o u l o m bs l 妒w i t hr e s i d u a ls t r e n g t hj o i n t sm e d i a ， r e f e r r i n gt ot a o y u a nc o a lm i n es t r a t i g r a p h ya n df a u l tf 2 ，av a r i e t yo fc o a l p i l l a r so fs e a m 10 # a n ds e a m81 # u n d e rt h ec o n d i t i o n so fd e f o r m a t i o na n d f a i l u r ec h a r a c t e r i s t i c so fr o c k sw e r es i m u l a t e d a n d n o n l i n e a rr e l a t i o n s h i p b e t w e e nt h ef a u l tw a t e r p r o o fp i l l a rw i d t ha n dt h i c k n e s so ft h es e a m sw a s r e s e a r c h e d ，a n dt h er e a s o n a b l ev a l u eo fc o a lp i l l a rw i d t ha b o u ts e a m1o 撑a n d s e a m8i # i sa l s ot r i e d b yai n d e p t ha n a l y s i st op r a c t i c a ld e f i c i e n c i e so ff a u l tw a t e r p r o o f p i l l a rf o r m u l ai nt h ep r o f e s s i o n a lr e g u l a t i o n ，t h e m e c h a n i c sm e a n i n ga n d a p p l i c a t i o n c h a r a c t e r i s t i c so ft h ef o r m u l a ，t h ei r r a t i o n a l i t yo fal i n e a r r e l a t i o n s h i p b e t w e e ns e a mt h i c k n e s sa n dp i l l a rw i d t h ，a sw e l la s t h e c o n t r a d i c t i o nb e t w e e nt h a ta n dt h ea c t u a ll a w so ff r a c t u r e d ，w e r ep o i n t e do u t d i s p l a c e m e n ta n g l et h e o r yw a sl e dt oa n a l y z et h ea p p l i c a b l ec o n d i t i o n s a n dd e f e c t so ft h ed e s i g nc a l c u l a t i o no nw a t e r p r o o fp i l l a rt of a u l t s k e y w o r d s ：d i s c r e t e e l e m e n t m e t h o d ( d e m ) ；f a u l t ；w a t e r - p r o o fp i l l a r ； m i n i n g i n d u c e df r a c t u r e 图2 1 图2 2 图3 1 图3 2 图3 3 图3 4 图4 1 图4 2 图4 3 图4 4 图4 5 图4 6 图5 1 图5 2 图5 3 图5 - 4 图5 5 图5 - 6 图5 7 图5 8 图5 - 9 图5 1 0 图5 1 1 图6 1 图6 2 图6 3 图6 4 图6 5 图6 6 图6 7 图6 8 图6 - 9 图6 1 0 图6 1 1 图6 1 2 插图清单 块体模型叠合示意图7 块体模型同步循环计算过程9 矿区交通图1 4 矿区地质构造示意图1 8 小断层走向玫瑰花图2 0 f 2 断层与奥灰强含水层对接示意图2 4 煤层与强含水层接触示意图。2 8 煤层与强含水层接触且局部被覆盖2 9 煤层与强含水层接触且煤层位于含水层上方意图3 0 “三带”分布示意图3 l 应力重分布示意图3 2 岩层移动角示意图3 3 研究区模型示意图3 8 单轴抗压实验照片4 0 剪切试验照片4 0 离散元初次划分计算模型4 2 模型边界应力赋值示意图4 3 初始化测试点位移平衡图4 4 初始化位移变化图4 4 模型初始化z 方向应力云图4 5 模型初始化y 方向应力云图4 5 模型初始化y 方向应力云图4 5 回采面示意图4 6 推进3 0 m 裂隙分布5 1 推进4 5 m 裂隙分布5 2 推进6 0 m 裂隙分布5 2 推进9 0 m 裂隙分布5 2 推进1 2 0 m 裂隙分布5 2 煤柱3 5 m 的位移及裂隙分布5 4 煤柱4 6 m 的位移及裂隙分布5 5 煤柱5 8 m 的位移及裂隙分布5 6 煤柱5 8 m 应力分布5 8 煤柱3 5 m 应力分布。5 9 8 l 煤位移及裂隙分布6 1 8 l 煤应力变化情况6 2 图6 1 3 图6 1 4 图7 1 图7 2 图7 3 8 l 煤3 5 m 煤柱裂隙发育图6 3 8 l 煤4 1 m 煤柱裂隙发育图6 3 公式理模型示意图6 4 无断层裂隙带高度6 7 岩层移动角对断层影响示意图6 8 表格清单 表3 1 煤层物理力学参数1 7 表3 2 断层情况一览表1 9 表3 3 小断层分类表1 9 表3 4 钻孔揭露的f 2 断层信息2 3 表3 5f 2 断层土体颗粒组成表2 5 表3 - 6f 2 断层带土工试验成果表一2 6 表3 7 断层带岩石物理、力学性质试验成果表2 6 表4 1 煤层开采的冒落带高度计算公式一3 2 表4 2 煤层开采的导水裂隙带高度计算公式一3 2 表4 3 各煤层防水煤柱宽度计算值3 4 表4 4 各煤层顶底板力学参数3 4 表5 1 模型几何参数3 7 表5 2 岩体力学参数4 0 表5 3 节理力学参数4 1 表6 11 0 煤防水煤柱宽度与围岩变形指标对比5 3 表6 28 l 煤与1 0 煤防水煤柱宽度与围岩变形指标对比6 0 表6 38 l 煤断层防水煤柱宽度值6 3 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究丁作及取得的研究成果。据我所 知，除了文中特别加以标注和致谢的地方外，论文中不包含其他人已经发表或撰写过的研究成果， 也不包含为获得 盒匿兰些太堂或其他教育机构的学位或证书而使用过的材料。与我一同 工作的同志对本研究所做的任何贡献均已在论文中作了明确的说明并表示谢意。 学位论文作者签名：j 多理签字日期：易。和年彩月刁日 学位论文版权使用授权书 本学位论文作者完全了解金月曼工些厶堂有关保留、使用学位论文的规定，有权保留并向国 家有关部门或机构送交论文的复印件和磁盘，允许论文被查阅和借阅。本人授权金照些态堂可 以将学位论文的全部或部分内容编入有关数据库进行检索，可以采用影印、缩印或扫描等复制手 段保存、汇编学位论文。 ( 保密的学位论文在解密后适用本授权书) 学位论文作者签名： 蝴期：劫缈年乒月刁日 学位论文作者毕业后去向： 工作单位： 通讯地址： 翩摊。孙，涉钢叮纱 签字日期：2 一明口年仁月7 日 电话： 邮编： 致谢 本文是在导师葛晓光教授的精心指导下完成的。导师在完成自己繁重的科 研和管理工作的同时，在本人论文选题、调研、撰写以及定稿过程中给予了大 量的指导，并在生活上关怀备至，使我能够顺利完成学业。在我攻读硕士学位 期间，导师严谨求实的治学态度、渊博的知识以及严肃认真的生活态度给我留 下了深刻的印象，我不仅从导师身上学到了专业知识，而且学到了很多为人处 世的道理，这些对我以后的学习和工作将受益终生。在论文完成之际，我谨向 尊敬的导师表示忠心的感谢和深深的敬意。 衷心感谢葛晓光教授、王国强教授、崔可锐教授、刘东甲教授、谭晓慧副 教授等老师在我攻读硕士学位期间给与我学习和生活上的帮助，感谢在论文撰 写和答辩过程中，您们给予的教导与帮助。 最后，特别感谢家人给予我精神上和物质上的支持和鼓励，没有你们的支 持，就不会有我的今天。 作者：王耀 2 0 10 年3 月2 6 同 第一章绪论 1 1 研究课题的提出 在煤矿生产和建设中经常遇到各种各样的断层问题；其中大断层因落差大、 断层带内力学性质复杂而对一些井巷工程影响很大【l 】。断层作为岩体中重要的 软弱结构面，对岩体的完整性产生严重的破坏，大大改变了岩体的力学强度： 此外开采后引起的扰动使得岩体内的应力重新分布，应力状态发生一定的变化， 可能使得断层内部产生裂隙或使老的裂隙活化，大量的资料显示，8 0 以上的 底板突水是发生在断裂构造附近1 2 】。在煤层工作面采动过程中，采空区上部覆 岩在张裂作用下发生移动、垮落，表现出上部覆岩移动变化的“三带”性，即冒 落带、导水裂隙带和弯曲带( 整体移动带) 。倘若导水裂隙带与断层连通，将有 可能对煤层工作面形成突水性威胁。对于解决此类问题现最常用的方法是留设 防水煤柱1 3 j 。 淮北桃园矿位于宿州市境内，是个年产量超过1 2 0 万吨的国有大型煤矿， 其中f 2 断层为该矿六采区与八采区的分界断层，其落差大于4 0 0 m ，断层带宽 约8 0 m ，是煤矿矿井内少见的大断层。f 2 断层带内角砾大小悬殊，大块角砾主 要由细砂岩、粉砂岩、泥岩组成，细颗粒主要由强度较低的泥岩类构成。1 0 煤 层、8 l 煤层作为八采区的主要开采煤层通过断层带与奥陶系灰岩强含水层对接， 为了在确保生产安全的基础上合理确定断层保护煤岩柱宽度，必须掌握采动冒 落引起的覆岩塌陷移动以及断层带及其周边岩层的物理力学变化情况，分析断 层带的变形程度可否达到引发突水的程度，为防水煤柱设计提供理论依据。为 此本论文以f 2 断层带附近煤层开采中的岩层破坏为研究对象，以非连续介质力 学数值计算方法为研究手段，模拟断层在采动过程变形规律，推测断层活化的 可能性以及断层水突出的可能性，为防水煤柱留设宽度的取值提供一种合理的 力学评价依据。 1 2 国内外研究现状 煤层在开采过程中，采动对覆岩的破坏是一个极其复杂的力学问题，主要 表现为岩层结构的强烈非均质性、应力变化的复杂性、采动岩层变形和破坏的 非线性【4 j 。为了获取合理的防水煤柱留设宽度值，减小采动破坏可能引发的突 水问题，有关学者做了大量研究。一些学者通过力学解析方法提出防水煤柱留 设宽度；一些学者则是通过数值模拟，分析采动后断层带附近的变形破坏程度 和上部覆岩形成的“三带 高度问题，进一步分析防水煤柱留设宽度值【5 】。 1 2 1 煤柱留设宽度力学解析法研究现状 彭文庆【6 i 等在研究不同断层倾角下防水煤柱合理宽度中应用弹塑性极限平 衡理论求解采空区煤柱屈服宽度，认为屈服宽度可以作为保护煤柱宽度的留设 值，并且运用弹塑性极限平衡理论、流体力学及土力学推导出断层倾角大小不 同时，防水煤柱宽度的计算公式。在计算公式中作者认为：煤层埋深、煤层厚 度、煤层倾角、煤的硬度、巷道断面大小、断层倾角和内摩擦角等作为煤柱宽 度的考虑因素。 王其芳1 7 】在断层保护煤柱留设方法研究中，以岩层移动角作为对象，分析 断层对开采沉陷的规律，认为现有井田内断层附近工作面留取煤柱宽度存在不 合理性；提出应以开采引起的地表移动变形边界角作为留设保护煤柱的主要参 数，或在适当条件下以边界角和岩层移动角的平均值来设计保护煤柱。 朱刘娟等j 在深部开采条件下岩层移动角的研究中指出，规程中给出 的移动角，在深部开采时不为定值，在留设保护煤柱时，不能使用规程中给出 的移动角，同时给出了深部开采中移动角的计算公式。 施龙青【9 】等在采场底板断层防水煤柱留设研究中将采场附近煤体划分为 弹性区和非弹性区，弹性区具有隔水性；根据采场支承压力分布特征，建立采 场底板突水力学模型，认为底板突水通道由附近断层和被矿山压力破坏的底板 岩层联合组成，推导出采场断层防水煤柱留设计算公式。 1 2 2 运用数值模拟对断层破坏和三带高度的研究 李晓昭等【lo 】利用数值模拟方法对原始地质条件下的不导水断层在采动影 响下透水机理进行分析，发现断层破碎带会成为采动后围岩变形和应力传播的 屏障，导致断层带与采动空间之间变形和应力集中加剧，断裂带内产生较大的 指向采动空间的位移梯度。 左建平等【ll 】在对深部煤矿采动诱发断层活动规律实验研究中通过二维相 似模型实验运用数字散斑相关方法( d s c m ) 和不连续变形方法( d d a ) 对实 验进行分析。结果表明随着开采进度的推进，采动多形成的岩体裂隙网络的分 布区域逐步向工作面和上覆岩层方向扩展，断层的存在阻碍了裂隙向断层下盘 的扩展；当采空区距断层15 4 0 m 时，断层开始活化，断层附近的水平位移呈马 鞍形，进而表明断层附近发生了剪切作用。 谭志祥【l2 】等运用a d i n a 有限元分析程序模拟中硬覆岩下裂隙发育情况， 发现推进由断层上盘向下盘前行时，比反方向推进产生的裂隙高度增长率大； 同时断层倾角不同时，裂隙增长率有所差异；重复采动条件下断层对裂隙带高 度增长率比初次采动减小一半多。通过不同情况研究建立导水裂隙带高度增长 率计算公式。 崔芳鹏等【1 3 l 在对断层防水煤柱安全宽度的计算分析中认为：对于断层防水 煤柱的安全宽度计算方法要考虑顺层水压的计算方法，导水裂隙带和保护带高 度等间接计算法，取两者较大值为保护煤柱宽度值。根据f l a e 3 d 数值模拟， 发现随着防水煤柱宽度减小，断层面上位移增大，当到达临界防水煤柱值时， 岩层破坏的垂向位移将发展到断层面上。 2 b 万奎【1 4 】等通过建立力学模型，根据弹性理论中半平面体边界上受法向分 布力的应力分析，推导出作用在采场底板的支承压力，分析了断层倾角对断层 面上剪切应力、法向应力及断层活化的影响；通过r f p a 2 d f l o w 模拟表明，采 动过程中小倾角正断层底板易产生贯通煤层底板、断层和含水层的裂隙，产生 渗流且渗流主要集中在工作面下方。 k i m i l5 j 在研究断层长度与位移量之间的关系中发现：位移的最大量与断层 带的长度之间存在如下关系：d 。= c l ”。 i s l a m 等 1 6 1 运用边界单元法，依据摩尔库伦准则对煤层开采导致的断层活 化进行了研究，发现：断层带应力重分布和围岩的形变在很大程度上受到断层 岩性的影响。 d o n n e l l y 【l 7 1 8 j 对于开采引起附近断层活化现象进行了研究，并讨论了断层 带不同岩性对煤田范围内存在断层时地面沉降的延迟速度作用。 g a l y b i n 坶j 研究了二维弹性问题，认为在煤层开采后当莫尔库伦准则满足 时，可运用解析法求出滑动增长量，并推测断层局部的滑动位移将会较大。 张均峰等【2 0 】通过流固耦合理论，利用有限元数值方法模拟承压水体上采 煤过程中局部断层参数对煤层、围岩的应力、变形和水压力变化的影响，得出 防水煤柱的留设距离随断距增大而增大；在断层倾角较大时，断层更易产生张 裂隙。因此，防水煤柱留设宽度随倾角加大而增加。 h a o 2 l 】运用离散元方法模拟了各种不同断层参数对采动后围岩稳定性的 影响，通过实验模拟找出了最大位移量、塑性变形区域与断层参数之间的关系。 1 2 3 国内外离散元应用现状 自从王泳嘉等将离散单元法的引入我国并在岩土工程方面得到了迅速的 推广【2 2 1 。现在离散元已经被大量运用在力学计算方面、岩土力学性质研究、地 下工程和煤矿开采分析等方面。 在力学计算方面，成名【2 3 1 和f e r r e r o 等【2 4 】分别将离散单元法运用于有关高 速冲击问题、地球物理和地质力学中。在高速冲击模拟中成名以能量等效为原 则，构建新的求解弹塑性轴对称问题的离散元模型，通过引入断裂准则，模拟 了钢板受冲击荷载产生层裂的过程。 岩土力学性质分析中主要运用p f c 软件对循环荷载下砂土液化特性颗粒流 进行数值模拟【2 5 1 和砂性岩土的性质进行模拟【26 1 。由于p f c 采用球体模型，因 此，能够很好的模拟实际情况。 在地下工程中应用，伍永田1 2 7 1 运用u d e c 对采矿区的稳定性进行了数值 模拟分析并预测了塌陷而导致的可能地面沉降范围，发现模拟结果与实际情况 较吻合。 魏龙海【2 8 1 对碎石土隧道自稳性进行了三维离散分析研究，探讨了粒径与隧 道自稳性之间的关系，同时模拟分析了碎石土隧道蹋落量与跨度之间的关系， 提出不同粒径碎石土的自稳跨度。 王吉亮【2 9 j 等采用离散元法对隧道塌方运动全过程进行模拟，对塌方过程中 块体形心处的位移、运动速度、加速度、应力的变化规律进行分析，证实节理 岩体隧道塌方与节理面强度及应力释放量有关。 陈陆望等1 3 0 j 运用u d e c 离散元分析软件，模拟煤矿开采中厚松散层及超薄 覆岩条件下防水煤柱覆岩突水可能性，发现u d e c 能够模拟出覆岩采动后实际 形成的“三带”问题，发现随着放采比的增大，煤层顶板冒落带高度与裂隙带高 度逐渐增大，裂隙带中的破坏裂隙从彼此独立的封闭系统逐渐发育成相互沟通 的开放系统。 胡戈等1 3 l j 对叠加开采顶板变形破坏运用离散元进行模拟研究发现：叠加开 采煤层工作面，煤层顶板破坏高度较单煤层开采明显增加；“三下”采煤规程 和矿井水文地质规程中的导裂带高度计算中存在不合适的地方。 1 3 课题研究意义 在传统的分析研究中，对于断层的力学模拟多基于连续介质的性质予以分 析。当前的计算方法主要有：有限元( f e m ) 、有限差分( f d m ) 、边界元( b e m ) 、 无网格法( e f m ) 1 3 2 j ，常用分析方法主要为有限元法。有限元法的特点是【3 3 】： 将结构离散化，用有限个容易分析的单元来表示复杂的对象，单元之间通过有 限个节点相互连接，然后根据变形协调条件综合求解。可供岩土力学分析的有 限元常用软件有f l a c 3 d 、a d i n a 、a n s y s 等。这些基于连续性介质模型的 软件在开展实际模拟时，主要缺点是无法解释和分析模拟材料在破坏后出现的 非连续现象，即当煤层工作面回采之后，上部覆岩发生变形、破坏，形成冒落 带( 存在岩体的移动、脱落) 、裂隙带( 岩体中出现新的裂隙，块体的整体性被 破坏) ，和弯曲带( 整体发生的塑性变形区) 。由于连续性介质方法无法分析采 动破坏后岩体产生裂隙的行为；不能考虑岩体破坏后各个块体之间的相互作用， 也无法考虑节理和结构面在模拟中的作用。因此，采用连续性介质模型研究煤 层开采中的岩体破坏时，理论基础上存在根本性缺陷。模拟结果和实际情况之 间往往存在较大出入。 鉴于这些原因，本次研究采用了非连续介质力学分析法。目前非连续介质 力学分析法中主要的计算方法是离散元( d e m ) 和不连续变形分析( d d a ) d 4 3 5 】。其中运用最多的是离散单元法。根据分析散粒体力学行为的方法它又可以 分为( 1 ) 统计力学模型、( 2 ) 经典牛顿力学模型、( 3 ) 混合模型1 36 l 。其中经 典牛顿力学模型是最常用的，它是将岩体划分为有限个刚性、弹性或弹塑性块 体，以牛顿第二运动定律为基础，结合不同本构关系，考虑块体受力后的运动及由 此导致的受力状态和块体运动随时间的变化。它允许块体间发生平动、转动， 甚至脱离母体下落。相对传统的连续介质分析方法，非连续介质分析的优点主 4 要是模拟非连续性介质受力后的力学性质和运动状态，使得模拟结果与物体不 连续变形情况相吻合，可真实反映出“三带”发育情况和回采后断层带附近的岩 体变形、破坏情况，从而为煤层开采设计提供准确可靠的依据。 1 4 课题研究主要内容 本次研究针对淮北桃园矿f 2 断层防水煤柱合理留设问题，开展了综合分 析。对于防水煤柱留设宽度我们通过数值模拟进行讨论研究。内容包括： 分析研究f 2 断层的产状、工程力学性质、1 0 煤与8 l 煤层通过f 2 断层与 对盘奥灰强含水层的对接情况。 分析煤层采动过程中，围岩的应力应变规律、覆岩破坏的“三带”表现形 式和相关高度。运用u d e c 4 0 数值分析程序对1 0 煤工作面围岩动态破坏进行 模拟，研究不同煤层及不同煤柱留设宽度条件下的断层带内位移变化情况，推 测断层带可否产生活化现象，同时分析煤层回采过程中断层带内位移和应力变 化的动态过程。 通过比较断层带附近的变形、破坏，讨论现行规范中防水煤柱计算公式 的合理性；讨论将岩层移动角作为煤柱留设的一项考虑因素引入计算时，对实 际问题的适应性。 1 5 课题创新点 本研究有如下特点：尽管近几年离散元软件被普遍运用于岩土工程分析 中，但只在岩质边坡稳定性、地下洞室开挖和煤层放项等方面有应用实例。本 次研究首次采用离散单元法模拟煤层与奥灰强含水层通过大型断层带对接条件 下，煤层采动引起断层带的形变和应力分布特点。 依据数值模拟的计算结果验证现行规范中导水断层煤柱解析解或经验公 式解的合理性。 结合实例分析对规范中有关断层煤柱公式在一定断层产状条件下的适应 性提出质疑和意见。 第二章离散元理论及软件 2 1 离散单元法的提出及发展 3 7 - 3 9 】 离散单元法是由p a c u n d a l l 于1 9 7 1 年提出的一种不连续介质的数值计算 方法，该方法最早用于研究准静态和动力条件下的节理系统或块体集合的力学 问题。1 9 7 9 年c u n d a l l 和s t r a c k 又提出了用于土力学的颗粒离散元法，并推出 了二维圆盘程序b a l l 和三维圆球程序t r u b a l ，即现在的应用软件p f c 2 d 3 d 的 原型，被称之为软颗粒模型。1 9 8 5 年，c u n d a l l 完成了可以考虑块体本身变形 的离散单元法，并编制了著名的u d e c ( u n i v e r s a ld i s c r e t ee l e m e n tc o d e ) 程序， 经过不断的修改与扩展，该程序已广泛应用于岩土力学及相关领域中。1 9 8 6 年 用于进行岩石块体三维离散模拟的软件3 d e c 问世。后来又经v o g e l 、l o r i g 、 b r a d y 等人的发展，离散元很快被用于边坡、基础、巷道稳定、放矿力学等方 面。在第一届全国岩石力学数值计算及模型试验讨论会上，王泳嘉教授首次介 绍了离散单元法的基本理论和应用例子，从而将离散元引入到我国。 2 2 离散元的理论基础及离散元模型 2 2 1 离散元的理论基础 离散单元法的理论以牛顿第二运动定律为依据【4 0 4 。作为处理不连续体或 像煤层采动后造成岩体连续性破坏问题的离散单元法，其基本思想是，把不连 续体分割为刚性元素的几何体，使各个刚性元素满足运动方程，当划分块体在 自身重力或者外力作用下产生加速度和相应的位移时，块体的空间状态发生变 化，位移后的块体与所接触的块体产生叠合。然后根据力一位移关系及运动方 程，使用时步迭代的方法求解各刚性元素的运动方程，进一步求不连续体的整 体运动形态。 2 2 2 离散元模型 常用的离散单元模型主要有块体模型和球体模型，块体模型主要用于模拟 分析岩石变形、破碎而造成的岩石连续性破坏现象；球体模型则主要被用于模 拟松散介质的不连续性特征。在此，简单介绍以下块体模型和球体模型。 2 2 2 1 块体模型 在块体模型模拟时，为达到更好的模拟效果常作如下假设【4 2 “4 】： 计算的各个步骤中，块体始终被当作是理想刚性块体，各个块体只是发 生空间上的平动、转动或脱离母体，自身并不发生形变。 所有块体的接触，当为面接触时，假定只有两个点接触，若非面接触则 认定只有一个点接触。 块体模型在进行实际模拟的情况下需要从三个方面进行主要的考虑，即： 块体的表示方法、块体的基本运动方程、块体的接触模型。 6 ( 1 ) 块体表示方法【4 5 l 块体的变形在实际工程中有着不可忽视的作用。因此，通常需要对块体变 形进行最简单的表示：把块体处理成为一个假定有应变分布的简单单元【4 9 5 们。 由于该方法是将块体整体化处理，因此对材料体系的变形会有一定的影响。目 前经过改进后常用的有两种表示变形块体变形性的方法，即“叠加法”和“直接 分割法”。叠加法是通过几个模型形状的叠加得到整个块体的复杂变形形式。缺 点在于叠加后不易体现出材料的非线性；而分割法是把块体划分为内部单元和 边界单元，以增加自由度的数量，划分的单元数目决定了可能出现的变形复杂 程度。 ( 2 ) 块体运动基本方程【3 8 , 4 2 , 4 5 】 a 力与位移的关系 离散单元中采用的是最简单的力和位移关系，并且认为在整个计算过程中 这个关系是绝对可逆的；两块体相互作用的力与接触块体之间的“叠合”量有关。 这里所谓的“叠合”是计算时假定的一个量，将它乘上一个比例系数作为力的一 种度量。 假设块体之间的法向力凡正比于他们之间法向“叠合”u 。( 见图2 1 ( a ) ) ， 即 e = 后。“。 ( 2 1 ) 式中k 一方向刚度系数。 f 。 ( b ) 图2 1 块体模型叠合示意图 倘若离散元块体的边界相互“叠合”( 见图2 1 ( b ) ) ，则有两个角点与界面 接触，可用界面两端的作用力来代替界面上的力。但是，界面接触的实际情况 要比图2 1 ( b ) 这种情况复杂的多。由于无法确定究竟是什么样的点相接触， 7 因此，一般采用最为简单的两个角点相接触的“界面叠合”模式。 对比块体法向力r ，块体所受的剪切，与块体运动和加载的历史和途径有 关，所以剪切力可以用增量表示为： z x f , = k ，坑 ( 2 - 2 ) 式中茂一界面的切向刚度系数。 式( 2 1 ) 和( 2 - 2 ) 所表示的是力和位移的弹性关系情况，当直接顶中岩 土发生破裂脱离母体的时，该弹性关系将不成立。因此，对于塑性剪切破坏情 况，需要在每次迭代时检查剪切力凡是否大于c + e t a n c p 。若大于则表示产生 了滑动，剪切力取极限值et a n 妒。 b 相邻块体间运动方程 依据块体的几何形状及邻近块体的接触关系，运用( 2 1 ) 与( 2 2 ) 依次计算出 某一特定块体上的一组力，根据计算出的这一组力可以算出它们的合力和合力 矩，并根据牛顿第二运动定律确定块体质心的线加速度和角加速度，进而确定 缸内的速度、位移及转动量。 坐：一f ( 2 - 3 ) o t所 上式左边用向前差分格式对于时间f 可表达为： 丝：坐兰型 ( 2 4 ) 将上面两式整理可得出： z i ( ，+ 出= f t 一出+ f ( t ) m a t ( 2 - 5 ) 式( 2 5 ) 中，时步点的速度用位移的形式写成： z 2 + 出= z 2 ( ，) a t + u 7 ( 2 - 6 ) 依据图2 1 ，力的产生依赖于位移，所以力和位移的计算是在同一时步内进行 的，这个循环前进的过程可用图2 - 2 的网络表示。 对于块体同时受到多个力及重力作用时，其速度方程将改写为： 扩加= t i + ( m + g 小血 ( 2 - 7 ) 矿= + ( m 。1 ) x a t ( 2 8 ) 8 图2 - 2 块体模型同步循环计算过程 式中矽一块体对于其形心的角加速度； l 块体的惯性距； g ，块体的重力加速度分量： 矽，一块体形心的速度分量；( i = x ，y ) 最后将得出新的速度可根据下式确定出块体的新位置： 妒+ 扪= 出+ ( 2 9 ) 谚坩= 谚7 a t + o d ( 2 1 0 ) 式中矽一块体绕其形心的转动量： 爵一块体形心的坐标量；( i = x ，y ) 上述的整个运动方程不包含使问题收敛于稳定解所需要的粘性阻尼，为此， 考虑引入粘性阻尼吸收动能，使系统收敛于稳定状态时的基本运动方程。 离散单元法现在最常用的粘性阻尼主要有：质量阻尼和刚度阻尼。当包含 有阻尼器时，离散元法的基本运动方程为： m i i ( t ) + c z 2 ( ，) + k u ( t ) = ( ，) ( 2 1 1 ) 式中所一单元质量； “一位移； 卜时间； c 一阻尼系数。 其中： c = c d m + c a k 式中c 口一质量阻尼比例系数； 9 ( 2 1 2 ) c 。一刚度阻尼比例系数； 枉刚度系数； 厂一单元所收的外荷载。 ( 3 ) 块体接触模型【4 4 】 8 接触的类型 在离散单元法中，常常把岩体中的节理面表示为两个相邻块体之间的接触 面。由于块体可能运动或者旋转，两个块体之间的接触就包括了角一角接触、 角一边接触和边一边接触三种情况。 离散单元中的接触按照接触方式分又可分为两种，即硬接触和软接触。所 谓的硬接触方式不允许两个块体之间相互嵌入，软接触则允许发生可测量到的 变形，软接触变形由一个在接触部位分别在法向和切向的有限法向和切向刚度 表示。对比两种接触方式，不难发现在岩石工程中模拟节理岩体形状时软接触 更加符合现实情况。在u d e c 中仅仅引入软接触，节理的变形源于在法向的挤 压和沿节理面方向的剪切作用。 b 接触的简化表示 当块体相互作用发生在相邻块体的角或两个相对的块体的角时，块体可能 竖直的锁住或悬空，因为假定块体的角是尖的或者具有无限的刚度。在这种情 况下，常常应用角圆化的方法使得两个相对的角相遇时块体能相互光滑的滑过 去，不至于相邻块体发生悬空现象。在块体模型分析软件u d e c 中，解决此类 问题常用方法是为每一个块体设定一个圆形弧来圆化。圆弧由实际尖点到相邻 边的切点的距离而确定。 c 节理的本构模型 节理的本构模型表示两个块体之间的相互作用力的关系。假定在一个步时 内法向和切向应力位移关系呈线性的，并且分别由法向刚度和切向刚度决定。 2 2 2 2 球体模型 对于松散介质而言，内部颗粒位移是相互独立的。它既不是理想的弹性材 料，也不是理想的塑性材料【4 6 1 。这就决定了对于这种介质在加载、卸载过程中 所表现出的复杂特征不能应用块体模型表示。为此，19 7 9 年c u n d a l l 和s t r a c k 又提出了球体模型。在该模型中，与块体模型的相似点在于单元之间也是靠相 邻单元之间接触的“重叠”相互作用产生接触力。两个单元中心距离与两个单元 半径之和大小的比较作为单位是否接触的判断依据。如果前者大于后者则认为 两个单元不接触；反之两个单元接触。在运用球体模型进行分析时，常作如下 假设【47 j ： 颗粒为刚性体，颗粒系统的变形是这些颗粒接触点变形的总和： 颗粒之间的接触发生在很小的区域内，即为点接触； 颗粒接触特性为软接触，但颗粒间的重叠量与颗粒尺寸相比很小，颗粒 l o 本身的变形相对于颗粒的平移和转动没有影响； 在每个时步内，扰动不能从任一颗粒同时传播到相邻颗粒。在所有时间 内，任一颗粒的合力可以由与其接触的颗粒之间的相互作用唯一确定。 ( 1 ) 力和位移关系 对于两个球体单元而言，接触只发生在焦点位置，并且颗粒间的接触力只是 颗粒处于相互挤压状态下才存在，颗粒间的拉应力和粘滞力需要通过定义颗粒 间的接触模型来表现。球体模型中的接触模型主要有“球一球”接触和“球一墙” 接触两种。 ( 2 ) 单元运动规律【4 8 4 9 1 在离散单元中，颗粒单元的位置随着时间的改变而运动，其中单元的位置 和接触力决定着颗粒的运动情况。在球体模型中，单元运动同样以牛顿第二定 律为基础。根据牛顿第二运动定律，可以得出离散单元法的基本运动方程为： f = m ( 2 ；一) 平移运动方程( 2 1 3 ) m 3 = 如旋转运动方程( 2 1 4 ) 式中f 一球体所受到的合力； m 一球体质量； z 一球体加速度； 岛一体积力加速度矢量； m 3 一球体所受到的合力矩； ，一球体的转动惯量； 破一球体的角加速度。 对于半径为r 的球形单元，倘若球体的质量是均匀分布的，那么球体的质心将 和球体的几何中心向重合，球体的主轴转动惯量将相同，即，= i ，。因此，球 体的旋转运动方程将简化为： m ，= 地 ( 2 15 ) 式中i = 一1 册r z 涸盘单元) 2 ，、 z _ ，_ m r 3 ( 圆球单元) 利用显式有限中心差分法对式( 2 13 ) 和( 2 15 ) 表示的离散元运动方程求解， 取时间步长为，。 根据& + ，时刻的速度，计算t 时刻的单元移动加速度和转动加速度，为： 弘扣m 掣删h ( 2 - 1 6 ) 桫，= 壶( 妒驯孙一科卜驯2 ) 将式( 2 1 6 ) 带入式( 2 1 4 ) 和( 2 1 5 ) 中计算出出2 + t 时刻的速度为： 最后运用上式中计算得出速度来更新单元的中心位置坐标： x ! 件出) = x ? + 戈! “f ，2 ) a t ( 2 1 8 ) ( 3 ) 单元计算循环 对于球体模型的单元计算循环与块体模型的计算循环类似，具体见图2 2 ， 实质上都是力一位移定律和牛顿第二运动定律进行遍历计算的过程。在计算过 程中，每个单元单独循环使用牛顿第二运动定律进行计算，每个接触点也对应 运用力一位移定律进行计算。应用运动定律给出一定时间内在受力情况下，颗 粒单元所产生的运动，最后运用力和位移方程计算出运动完成后颗粒单元位移 及其所产生的接触力分布。所有单元不停地重复整个计算过程最终将完成整个 离散单元法的动态模拟。 2 3 离散元软件 从p a c u n d a l l 于1 9 71 年提出离散元分析方法后，以此为基础的应用软 件应运而生，主要有u d e c 、3 d e c 、p f c 、p f c 3 d 、e d e m 和t r u d e c 等。 现在比