冲击矿压的优化有限元分析.pdf

第3 0 卷第1 期 2 0 0 1 年1 月 中国矿业大学学报 J o u r n a lo fC h i n aU n i v e r s i l yo fM i n i n g ＆T e c h n o l o g y V 0 1 ．3 0N o ．1 J a n ．2 0 0 1 文章编号，1 0 0 01 9 6 4 2 0 0 1 O l 一0 0 5 9 0 4 冲击矿压的优化有限元分析 刘卫群，缪协兴，杨静 中国矿业大学理学院，江苏徐州2 2 1 0 0 8 摘要利用优化非协调元和杂交元分析了三河尖矿冲击矿压，并与A N S Y S 程序的计算结果进 行了比较．认为优化有限元在冲击矿压乃至岩石材料的数值模拟中能给出更为精确的结果，也可 减小大型计算中由于人力、机时所限粗分网格产生的误差．分析中，采用了广义平均应力率准则 来判断冲击点位置和冲击时刻，所得结论与实际情况吻合． 关键词优化有限元；冲击矿压预测；广义平均应力率准则 中图分类号T D3 2 4文献标识码A 目前有限元商用软件中多是采用协调的位移 单元．这种单元过分强调了单元边界位移的连续， 实际结果显示当单元划分较少时会因为过份强调 单元边界连续而丧失部分柔度；而庞大的单元划分 从人力和机时考虑，特别是三维问题，有时又是不 现实的，这就需要人们从单元上去寻求另外的解决 办法．非协调单元无疑是一种有效的办法，因为它 放松了单元边界的连续，从而可以弥补前面所说的 协调元损失的那部分柔度，不过，由于非协调元受 收敛性的影响，不宜滥用．在国际上，能通过小片试 验一般作为非协调元收敛的标准[ 1 ] ，下文中的 N H l l 单元和基于此的1 8 B 杂交元均能通过小片 试验，在文献[ 2 ，3 ] 中被称为优化单元．本文应用这 两种单元分析了煤矿的冲击矿压 冲击地压 现象． 近年来，随着开采深度和强度的加大，冲击矿 压发生的频率和危害程度都有上升的趋势．各国学 者在研究冲击矿压的过程中总结了一些规律““] ， 但各地发生冲击矿压的情况差异决定了这些规律 仍然存在局限性．下面对徐州三河尖矿9 1 0 8 工作 面发生的冲击矿压用不同单元进行模拟分析． 1 非协调元N i t l l 和1 8 口杂交应力元 8 节点六面体等参元是协调元，收敛性好，所 以在各种商用软件中都有采用，这里只需将这种单 元稍加改造，可得到性能更好的N H l l 非协调元． 记8 节点六面体等参元 图1 的线性插值位 图18 节点六面体元 F i g ．1 8n o d eh e x a h e d r a le l e m e n t 对u 。附加非协调内位移 h 1 “a 一[ P 矿P ] { A 2 净N - A ， 2 l J 则单元位移“一u 。 u 。对单元自然坐标是二次完 备的．但是对非规则形状的三维体元，式 1 不满足 分片检验条件 P a t c hT e s tC o n d i t i o n ，i ．e ．P T C ， 而只有满足了P T C 才能保证解的收敛性，现加以 修正．定义函数 N ’一[ e _ 翱． 3 相应的内参为A 。，此时 “一‰ “- “ 一[ 札舢㈥ Ⅳ⋯。 4 式中”为修正项；口， 为协调和非协调节点位 移． 分片检验条件的强形式能表示成 收藕日期l2 0 0 0 0 2 2 9 基盒项目。国家自然科学基金重点鹰助项目 5 9 7 3 4 0 9 0 ，煤炭高校青年基金资助项目 作者筒介。刘卫群 1 9 7 0 一 ，男，上海市人，中国矿业大学讲师．博士研究生，从事计算力学、采动岩体力学方面的研究． K D矗 秒玑 o ￡ 1 8 量。∑H 鼢 驴 万方数据 中国矿业大学学报 第2 9 卷 矧础一o ， ㈣ 积分号表示沿单元边界积分，n f z ，m ，n 是单 元边界的外法线方向余弦．为了满足这一条件．必 觚一f 长} Ⅳ强非奇异．这是可以办舭 复记P f { ≥} c 札Ⅳ。，出 c 只只，， c s ， P 黔蹦‘一D ， 7 “ 眦蜘一J 】v E ㈣㈥． 8 孵 0 0 0 0 0 为满足虚功方程，单元边界高阶应力在非协调 位移上所做的功必须等于零，即 f 矿妇h 巩出 a a ’肘{ } 一o ， c - ， 式中啦一{ 兰。} [ Ⅳ3M Ⅳ] { } ， 肘一f ，I Ⅳ1 Ⅳ，．v 。lr J ] c 一；，，。。。出 [ M ，肘，] ， 1 2 这里西i9 ] 前6 行与[ 9 t9 。] 相同，7 ～9 行与 脚。p - ] 的4 ～6 行相同． 因为鼢1 ≠0 ，则 肘㈥一㈨引引_ o Ⅲ， 若行列式{ M 。1 ≠0 ，则p ，可由芦。表示，并可从式 9 中消去，最终应力试解可写成 州‰] ㈥一卵， 1 4 其中钆一9 。一9 。J l f i l M 。． 至此，应力的单独插值已完成，位移插值形式 不变．因为用到的单元应力参数共有1 8 个分量，故 称1 8 B 杂交应力元．1 8 口元应力是单独插值计算的， 这样，就得到了满足分片检验条件的非协调位移形 式．通过这种方法构造出来的8 节点六面体非协调 元称为N H l l 元． 杂交元法是采用多变量插值，但在最终的矩阵 方程中只包含节点位移未知量的方法．如果应力单 独插值，生成的单元就是杂交应力元． 在N H l l 元基础上．现对应力单独进行插值， 设 a 一| o 。 a ，d } T } 。z 。t 。P 1 1 一 以十叽一晟巾z 州㈥． 9 这里，巩为常应力；巩为高阶应力{ 屡，既，p 。为单元 应力参数 扉一{ 主} ，卢。一{ 竺。} ，卢x 一{ 竺} ； 000 7 00 000 0f0 00} 00 0 1 0 无须通过位移得到，故解的精确性更高． 2 冲击矿压的非协调元分析及广义平均应 力率准则 三河尖煤矿9 1 0 8 ／作面位于该矿西一采区，采 深6 0 0 m ；上邻9 1 0 6 已采工作面，下邻9 1 1o 未掘进工 作面．该面走向长3 9 8 m ，倾斜长7 6 ～1 2 5 m ，平均 9 5 m ．煤层倾角8 ～2 5 。，平均1 8 。；煤层厚0 ．9 ～ 3 ．8 5 m ，平均2 ．8 m ，并有o ～0 ．9 m 的夹矸．伪顶为 厚0 ．5 m 的粉砂质泥岩，直接顶和老顶合二为一． 为厚9 ．3 m 的中砂岩，直接底为厚0 ．9 m 的泥岩，老 底为厚0 ．5 m 的粉砂岩和细砂岩互层．发生冲击矿 压是当回采工作面已推进1 0 m 时．其位置在巷道 与工作面交接处 见图2 ． v ，m 2 5 0 2 1 0 8 0 6 0寓 E } j i i i i ；；；j i E } E ；i 自 0 9 03 0 0 ，，m 图29 1 0 8 T 作面推进至9 0 I T * 时模型剖面 F 塘．2 M o d e lp r o f i l eo fc o a l [ a c e9 1 0 8 w h e na d v a n c i n g9 0 m 0 O 0 0 叩0 O O 0 拿0 O 0 O O 0 O f 0 O 却O 0 0 0 O 叩0 0 0 O 0车0 O O 0 弦O O 0 O 0 e 0 0 O O一旧旧旧旧㈨∞ 一 ]尹尹[ 万方数据 第1 期 刘卫群等冲击矿压的优化有限元分析 2 ．1 冲击矿压的非协调元分析 分析模型取到老顶，煤层作为底层．四周边界 是两个沿y 方向的边限制z 位移，两个沿, 7 2 方向的 边限制Y 位移，底面限制z 位移．载荷的选取通过 下面的计算确定假设老顶以上各岩层平均密度2 4 k N ／m 3 ，老顶所在深度为6 0 0 m ，则上表面均布压力 g 一1 4 ．1 1 2 M P a ． 岩层材料参数 表1 ．计算采用中国科学技术 大学计算力学研究室自编有限元程序，所选单元为 N H l l ，因程序容量所限，共划分4 5 7 个单元，8 3 7 个 节点 图2 ． 表19 1 0 8 “ r 作面岩性参数 T a b l e1 p r o p e r t i e so rr o c ks t r a t af o rc o a lf a c e9 1 0 8 烘2 ．830 ．3 01 3 ．2 3 泥岩‘0 ．51 00 ．2 42 54 8 中砂岩～ 9 ．32 50 2 22 6 ．4 6 * 伪顶；* * 直接顶 老顶． 选取工作面位于T 一9 0 ．1 0 0 ，1 1 0 ，1 2 0 ，1 3 0 和 1 4 0 m6 个不同位置工况进行研究．不同工况下所 取样点坐标和数值结果列于表2 和表3 中． 表29 1 0 8 w 作面所选样点的坐标 T a 啪e2C o o r d l 呻t e so fc h o 钟n0 0 m t sf o rc o a l f a c e9 1 0 8 序号x ／m y ／m ／m 18 02 1 02 ．8 29 02 1 02 ．8 310 02 1 02 8 4 1 1 0 2 1 02 ．8 5 1 2 02 1 0 2 ．8 裹3 口一的非协调元数值结果 3 a b l e3N m e r i e 州r e s l l I bo t 口- ，埘 N H l le l e m e n tm e t h y lf o rc o R If a c e9 1 [ 1 8 M P h 2 ．2 广义平均应力率准则 为了分析冲击矿压发生过程，现引入广义应力 率的概念，记为凼．下标z 指与工作面推进长度“ 相关，以区别于绝对变化参量时问t ． 广义平均应力率指工作面每推进单位长度应 力的变化量，极限意义下为 匆一慨玺 面d o ． 1 5 为了在符号上与极限意义下的广义应力率击 区别开来．这里将广义平均应力率记为孑。 表4 给出了工作面推进过程中各点的广义平均 应力率． 裹4 广义平均应力率的非协调元分析结果 T a b l e4 A v e r a g e dg 嘶盯A I ‘z e ds t r e 脯r a t 器 c a l c u l a t e db yN H ”e l e m e n tm e t h y lk P a ／m 推进距／m ，芋{ } 一 9 0 ～1 0 01 0 0 ～1 1 0 1 1 0 ～1 2 01 2 0 ～1 3 01 3 0 ～1 4 0 一6 4 1 3 9 3 4 - 17 1 2 9 1 9 2 7 2 2 13 5 广义平均应力率准则认为，即将发生冲击矿压 时，冲击点附近广义平均应力率最大．具体地说，就 是冲击点附近的应力变化在即将发生冲击矿压时 要比冲击前和冲击后来得剧烈；冲击时，冲击点附 近的应力变化要比其他位置来得剧烈．另外，在冲 击发生之前，随着工作面的推进，冲击点附近广义 平均应力率由小变大，之后则由大变小． 利用这一准则和表4 的分析结果可以看出， 9 1 0 8 工作面由9 0 m 推进至1 0 0 m 时，各点应力变化 较其它过程都剧烈，因此这个过程是最危险时段． 工作面推进过程中样点2 ，坐标为 9 0 ，2 1 0 ，2 ．8 点的广义应力变化率始终较其它位置更大一些，因 此样点2 附近是危险区域，这与实际情况相符． 3 冲击矿压的杂交应力元分析 与非协调元分析采用相同模型，计算仍用中国 科学技术大学自编程序完成，单元选取18 口杂交应 力元，工作面不同位置的数据结果以及应力变化情 况见表5 和表6 ． 衰5 9 1 0 8 工作面口。，的杂交元数值结果 T a b l e5N a m e r l e a lr e s u l t so f 口_ ，f o rc o a lf a c e 如蛆垃】业盟 婴塑 坐曼 鲤 丛旦 出。 推进距／m ”8 丽而而函面丽一 12 7 ．1 92 6 ．6 12 6 ．4 22 6 ．3 12 62 7 2 62 6 22 9 ．1 12 7 ．8 4Z 7 ．5 62 7 ．4 42 7 ．3 82 7 3 5 32 6 ．3 32 5 ．9 92 5 ．7 92 5 7 l2 5 ．6 72 5 ．6 6 42 1 ．0 82 0 ．9 02 0 ．8 02 0 ．7 42 0 ．7 42 0 ．7 4 51 9 ．7 619 ．6 61 9 ．6 11 9 ．5 6 1 9 ．5 6 1 9 ．5 6 Z Z 盟 Z 坠 Z Z 表6 广义平均应力辜的杂交元分析结果 T a b l e6A v e r a g e dg e n e r a l i z e ds t r e 镕F a t e s ￡e ￡ 鲤丝 盟盟 罂型巴 坠鲤蔓￡ 也 串导 堡丝里 里 ⋯9 0 ～1 0 01 0 0 ～1 1 01 1 0 ～1 2 01 2 0 ～1 3 01 3 0 ～1 4 0 15 81 91 14～1 212 72 81 2 6 3 33 42 0一841 41 8一l O一6 00 51055 00 万方数据 中国矿业大学学报 第2 9 卷 由广义平均应力率准则可得到与非协调元判 断相同的结论． 4冲击矿压的A N S Y S 分析及与优化单元 结果的比较 为了便于比较，用A N S Y S 程序计算9 1 0 8 I 作 面时，模型仍选取非协调元和杂交元所使用的模 型，只是网格划分得更密一些．计算中将模型划分 为4 9 1 2 个单元。7 0 2 0 个节点．应力结果列于表7 中， 具体的冲击矿压判断过程与前面两种单元时相同， 预测结果也一致． 衰79 1 0 8 “ 1 - 作面口～的A N S Y S 数值计算结果 T a b l e7N u m e r i c a lr e s u l t so f 口‘。f o r 型 望 ] 竖 垡塾坠坠 序号 苎堂堕 婴 塑 塑 塑 塑 l2 6 ．2 32 5 ．5 22 5 ．2 32 5 ．0 62 49 62 4 ．9 3 22 8 ．3 72 6 ．7 82 6 ．3 82 6 ．1 72 6 ．0 42 5 ．9 8 3 2 5 ．1 5 2 4 ．6 4 2 4 ．3 62 4 ．2 12 4 ．1 32 4 ．1 1 4 1 9 ．8 4 1 9 ．6 31 9 ．4 31 9 ．3 1 1 9 ．z 21 9 ，2 1 51 7 ．9 11 7 ．6 81 7 ．5 3i 74 6 i 7 ．4 2 1 7 ．4 2 61 6 ．0 61 5 ．9 11 5 ．8 51 5 ．7 91 5 7 7i 5 7 7 将非协调元、杂交应力元和A N S Y S 用到的8 节点六面体元对应结果作一比较，不难发现结果基 本一致．但在非协调元和杂交元程序中用了比 A N S Y S 分析少得多的单元和节点，计算花费的时 间自然也较A N S Y S 少，如果单元和结点划分与 A N S Y S 相同，优化单元将得到更为精确的结果． 5 结束语 本文利用有限元程序分析了三河尖矿9 1 0 8 工 作面冲击矿压发生时的应力分布情况，其中非协调 元N H l l 及杂交元1 8 p 都是首次用于此类岩土工程 问题的分析，在稀疏网格情况下，得到与A N S Y S 加密网格同样精确的结果．两种单元都放松边界位 移连续的要求，且能通过分片检验条件，不存在任 何收敛性问题，相信它们在岩土工程中定会有广阔 的应用前景．利用对数值结果的分析，本文提出了 广义平均应力率的判别方法．该方法认为，冲击矿 压发生时．冲击点附近的应力变化较该处其他时刻 的变化和其他位置的变化都剧烈．利用这一判别方 法对9 1 0 8 I 作面冲击矿压的发生地点和时闻的分 析是合理的． 参考文献； E 1 ] I r o n sB ，R a z z a q u eA ．E x p e r i e n c ew i t hp a t c ht e s tf o r c o n v e r g e n c eo ff i n i t ee l e m e n tm e t h o d [ A ] ．I n A z i z A ．M a t h e m a t i c a lf o u n d a t i o n so ft h ef i n i t ee l e m e n t m e t h o d [ c ] ．U K A c a d e m i cP r e s s ，1 9 7 2 ．5 5 7 5 8 7 ． [ 2 ] 吴长春，B u f f e rH ．多变量有限元相容性与模式优 化[ J ] ．中国科学 A 辑 ．1 9 9 0 ，9 9 4 6 9 5 6 ． [ 3 ] 吴长春．卞学} 黄非协调原理与杂交元方法[ M ] ．北 京{ 科学出版社，19 9 7 ．1 3 91 4 5 ． E 4 3胡光伟．冲击矿压机理的数值分析[ D ] ．徐州中国 矿业大学数力系，I9 9 9 ． [ 5 ] M O i l e rW ．N u m e r i c a ls i m u l a t i o no fr o c kb u r s t s [ J ] ． M i n i n gS c i e n c ea n dT e c h n o l o g y ，1 9 9 1 ．1 2 2 7 4 2 ． [ 6 3 佩图霍夫．冲击矿压和突出的力学计算方法[ M ] ．段 克信泽．北京煤炭工业出版社，1 9 9 4 ．1 1 7 - 1 4 1 ． A n a l y s i so fR o c kB u r s tb yO p t i m i z e dF E M s L l UW e i q u n ，M I A OX i e - x i n g ，Y A N GJ i n g C o l l e g eo fS c i e n c e s ．C U M T ．X u z h o u ，J i a n g s u2 2 1 0 0 8 ，C h i n a A b s t r a c t I n c o m p a t i b l ee l e m e n tN H l1a n dh y b r i ds t r e s se l e m e n t1 8 t 9w e r ee m p l o y e dt oa n a l y z er o c kb u r s to c c u r r i n gi nc o a lf a c e9 1 0 8o fS a n h e j i a nC o a lM i n e ．C o m p a r e dw i t ht h er e s u l t sf r o mA N S Y S ，t h et w oo p t i m i z e dF E M sc a ne x e c u t em o r ea c c u r a t en u m e r i e a ls i m u l a t i o nj nt h er e s e a r c ho fr o c kb u r s ta sw e l 】a so t h e r r o c kp r o b l e m sa n dl e s s e nt h ee r r o r sf r o mi n a d e q u a t e l ym e s h i n gi nh u g ec a l c u l a t i o n sd u et ot h el i m i t a t i o no fl a b o ra n dc o m p u t e rt i m e ．I na d d i t i o n ，t h ec r i t e r i o no fa v e r a g e dg e n e r a l i z e ds t r e s sr a t ew a sp r o p o s e da n du s e dt o d e t e r m i n et h es i t ea n dm o m e n to fr o c kb u r s t ．T h eg i v e nc o n c l u s i o ni si na c c o r d a n c ew i t ht h ef a c to fc o a lf a c e 9 1 0 8v e r yw e l l ． K e yw o r d s lo p t i m i z e dF E M s } f o r e c a s t i n gr o c kb u r s t ；c r i t e r i o no fa v e r a g e dg e n e r a l i z e ds t r e s sr a t e 万方数据