基于离散单元法的岩体破碎关键块体定位方法与流程

本发明属于岩土工程领域，特别涉及一种基于离散单元法的岩体破碎关键块体定位方法。

背景技术：

岩体破碎工程广泛存在和应用于滑坡、山崩等地质灾害、矿山开采、油气井钻进、地质勘探、隧道掘进以及国防建设等领域中。深入研究岩体破碎过程中产生的岩块的运动对揭示岩石破碎机理具有十分重要的意义。按照石根华博士与GOODMAN R.E提出的块体理论，岩体破碎过程中产生的块体可分为关键块体和非关键块体，其中，关键块体是指边坡开挖或岩石破碎过程中，在工程作用力和自重作用下首先产生失稳滑动的块体，这些块体的变形很大程度上决定了岩体工程的稳定性。关键块体的定位有助于正确认识岩体的破碎机理，科学地评价岩体的长期稳定性，对公路、铁路、水利、矿山等边坡防治工程中的防护网安装等具有指导意义。此外，通过对关键块体进行定位，预先对关键块体进行加固、防护等处理，有助于预防工程地质灾害。

目前，关键块体的定位方法主要有矢量分析法、赤平投影法及赤平投影解析法，其中，矢量分析法难度较大，计算过程复杂，且所得结果不够直观；赤平投影法虽然较直观，但由于作图误差，难免会丢失或人为增加一些关键块体，其结果的精度难以保证；而赤平投影解析法则只考虑了结构面的形状，未考虑其具体空间位置，结果只能定性地分析出某个块体是否为关键块体，但无法确定该关键块体的具体位置。基于上述技术现状，开发出更简单和准确地确定岩体破碎关键块体的方法，对于深入了解岩体破碎机理和直到实际岩体破碎工程都将产生重要的意义。

技术实现要素：

本发明的目的在于克服现有技术的不足，提供基于离散单元法的岩体岩石破碎关键块体定位方法，以提高岩体破碎关键块体定位的准确性并简化岩体破碎关键块体定位方法。

本发明提供的基于离散单元法的岩体破碎关键块体定位方法，步骤如下：

(1)采用离散单元法将目标岩体表示为由N个离散颗粒组合成的集合体，赋予每个颗粒不同的ID编号并记录各颗粒在目标岩体中的三维坐标，所述ID编号为0～N-1的正整数，记作ID0，ID2，…，IDN-1，然后采用离散单元法对由离散颗粒组合成的集合体的破碎过程进行数值模拟，数值模拟结束后，得到M对胶结颗粒对，将胶结颗粒对用[IDm,IDn]表示，将M对胶结颗粒对的ID编号列于表A中，表A表示为{[IDm₁,IDn₁]，[[IDm₂,IDn₂]，…，[IDm_M,IDn_M]}；

(2)对表A中的M对胶结颗粒对进行循环检索，筛选出步骤(1)对破碎过程数值模拟完成后产生的关键块体；关键块体的筛选方法如下：

①将表A中的第1对胶结颗粒对中各颗粒的ID编号依次与表A中的其他胶结颗粒对中各颗粒的ID编号进行比较，比较的结果为：第1对胶结颗粒对与表A中的其他胶结颗粒对中没有相同ID编号的颗粒，或者是第1对胶结颗粒对与表A中的其他胶结颗粒对中存在ID编号相同的颗粒；

若第1对胶结颗粒对与表A中的其他胶结颗粒对中没有相同ID编号的颗粒，则将第1对胶结颗粒对从列表A中删除；

若第1对胶结颗粒对与表A中的其他胶结颗粒对中存在ID编号相同的颗粒，则将具有相同ID编号的颗粒的胶结颗粒对从表A中删除并将其ID编号列于表B中；将表B中各颗粒的ID编号与表A中的胶结颗粒对中各颗粒的ID编号进行比较，将具有与表B中各颗粒ID编号相同的胶结颗粒对从表A中删除并将其ID编号列于表B中，重复前述操作，直到表B中各颗粒与表A中的胶结颗粒对中没有相同ID编号的颗粒时停止比较，然后删除表B中重复的ID编号得到表C，若表C中的ID编号的数量g≥h个，h为不小于10的整数，则这g个ID编号的颗粒组成的岩体块体即为岩体破碎关键块体；

②按照步骤①的操作，依次对经过上一步骤处理的表A进行处理，直到表A中无剩余胶结颗粒对时结束筛选，即筛选出所有的关键块体；

③根据步骤②筛选出的关键块体中各颗粒的ID编号和在目标岩体中的三维坐标，确定出各关键块体在目标岩体中的位置，即完成岩体破碎关键块体的定位。

上述方法中，h的值的大小根据目标岩体的的大小进行确定，通常，h为10～1000的整数。

上述方法中，离散单元法是模拟固体颗粒动态模型的数值方法，该数值方法被广泛应用于岩土工程中(参考文献[1]Cundall P.A.,Strack O.D.L.,Geotechnique,1979,29(1):47-65.参考文献[2]Cundall P.A..A comprter model for simulating progressive:large-scale movements in blocky rock systems,in Proceeding of the International Symposium on Rock Fracture.1971:Nancy,France.)，该方法采用胶结颗粒模型代表岩体试样，即将岩体视为刚性球体颗粒的集合体，并在颗粒之间的接触处形成具有一定强度的胶结(参考文献[3]刘连峰,王泳嘉.有色金属.1995,3:16-19.)，胶结的断裂即代表岩体的破坏。上述方法的步骤(1)中，可采用离散元软件PFC、UDEC，以及开源代码ESyS-Particle等工具将目标岩体表示为由N个离散颗粒组合成的集合体并对由离散颗粒组合成的集合体的破碎过程进行数值模拟，具体模拟方法可参考文献[4]～[6](参考文献[4]Jiang MJ,Chen H,Crosta GB.,Int J Rock Mech Min Sci.2015,78:175-189.参考文献[5]Hazzard JF,Young RP,Maxwell SC.,J Geophys Res.2000,105(B7):16683-16697.参考文献[6]Abe S.,Place D.,Mora P.,Pure and Applied Geophysics,2004,161(11-12):2265-2277.)。

与现有技术相比，本发明具有以下有益效果：

1.本发明提供的基于离散单元法的岩体破碎关键块体定位方法是一种岩体破碎关键块体定位的新方法，该方法在离散单元法对目标岩体的破碎过程进行数值模拟的基础上进行岩体破碎关键块体的筛选定位，具有定位过程简单、定位精度高的优势，克服了现有矢量分析法难度较大，计算过程复杂、赤平投影法精度难以保证、以及赤平投影解析法无法确定关键块体位置的不足。

附图说明

图1为将目标岩体表示为由离散颗粒组合成的集合体的示意图；

图2为采用离散单元法模拟岩体破碎过程的示意图，其中，图(a)为模拟前的示意图，图(b)为模拟完成后的示意图；

图3为实施例中得到的11个关键块体的形态示意图；

图4为实施例中关键块体在目标岩体中的分布示意图，图中a为关键块体1、b为关键块体2、c为关键块体3、d为关键块体4、f为关键块体6、g为关键块体7、h为关键块体8、i为关键块体9、k为关键块体11。

具体实施方式

下面通过实施例对本发明所述基于离散单元法的岩体破碎关键块体定位方法作进一步说明。有必要指出的是，以下实施例只用于对本发明作进一步说明，不能理解为对本发明保护范围的限制，所属领域技术人员根据上述发明内容，对本发明做出一些非本质的改进和调整进行具体实施，仍属于发明的保护范围。

实施例

本实施例中，采用本发明所述方法对岩体岩石在自重作用下破碎过程中产生的关键块体进行定位，具体步骤如下：

(1)选定的目标岩体是尺寸为100mm*50mm*25mm的长方体岩体，采用离散单元法将该目标岩体表示为由10880个离散颗粒组合成的集合体，如图1所示，赋予每个颗粒不同的ID编号，所述ID编号为0～10879的正整数，记作ID0，ID2，…，ID10879，每个ID编号代表的颗粒都对应了一个其在目标岩体中的三维坐标(x,y,z)，由于数量巨大，此处不列出每个ID编号代表的颗粒的三维坐标；然后采用离散单元法对由离散颗粒组合成的集合体的破碎过程进行数值模拟，数值模拟结束后，由离散颗粒组合成的集合体中会有部分颗粒胶结在一起形成块体，图2为采用离散单元法模拟岩体破碎过程的示意图，其中，图(a)为模拟前的示意图，图(b)为模拟完成后的示意图。

该步骤具体采用的是离散元开源代码ESyS-Particle，基于Pyhon计算机语言编写的程序将目标岩体表示为由离散颗粒组合成的集合体并进行数值模拟的，数值模拟结束后得到11325对胶结颗粒对，将胶结颗粒对用[IDm,IDn]表示，将11325对胶结颗粒对的ID编号列于表A中，表A表示为{[373,436]，[468,469]，[468,1214]，[469,1214]，[868,1590]，…，[10872,10873]，[10874,10875]，[10875,10876]，[10876,10877]，[10877,10878]}，由于数据量巨大，此处不一一列出表A中11325对胶结颗粒对的ID编号。

(2)对表A中的11325对对胶结颗粒对进行循环检索，筛选出步骤(1)对破碎过程数值模拟完成后产生的关键块体；

关键块体的筛选方法如下：

①将表A中的第1对胶结颗粒对中各颗粒的ID编号依次与表A中的其他胶结颗粒对中各颗粒的ID编号进行比较，比较的结果为：第1对胶结颗粒对与表A中的其他胶结颗粒对中没有相同ID编号的颗粒，或者是第1对胶结颗粒对与表A中的其他胶结颗粒对中存在ID编号相同的颗粒；

若第1对胶结颗粒对与表A中的其他胶结颗粒对中没有相同ID编号的颗粒，则将第1对胶结颗粒对从列表A中删除；

若第1对胶结颗粒对与表A中的其他胶结颗粒对中存在ID编号相同的颗粒，则将具有相同ID编号的颗粒的胶结颗粒对从表A中删除并将其ID编号列于表B中；将表B中各颗粒的ID编号与表A中的胶结颗粒对中各颗粒的ID编号进行比较，将具有与表B中各颗粒ID编号相同的胶结颗粒对从表A中删除并将其ID编号列于表B中，重复前述操作，直到表B中各颗粒与表A中的胶结颗粒对中没有相同ID编号的颗粒时停止比较，然后删除表B中重复的ID编号得到表C，若表C中的ID编号的数量g≥h个，h＝10，则这g个ID编号的颗粒组成的岩体块体即为岩体破碎关键块体；

以下述三类情况为例，具体说明关键块体的筛选方法：

(Ⅰ)对表A中的第1对胶结颗粒对[373,436]中两个颗粒的ID编号依次与表A中的第2对胶结颗粒对[468,469]、第3对胶结颗粒对[468,1214]，…，第11325对胶结颗粒对[10877,10878]中各颗粒的ID编号进行比较，比较结束发现第2对～第11325对胶结颗粒对中没有任何一对胶结颗粒对中的颗粒与第1对胶结颗粒对中的颗粒存在相同的ID编号，因此，将第1对胶结颗粒对[373,436]从表A中删除；

(Ⅱ)经过步骤(a)的筛选后，表A中还剩余11314对胶结颗粒对，此时表A变为{[468,469]，[468,1214]，[469,1214]，[868,1590]，…，[10872,10873]，[10874,10875]，[10875,10876]，[10876,10877]，[10877,10878]}。对表A中的第1对胶结颗粒对[468,469]中两个颗粒的ID编号依次与表A中的第2对胶结颗粒对[468,1214]、第3对胶结颗粒对[469,1214]，…，第11324对胶结颗粒对[10877,10878]中各颗粒的ID编号进行比较，比较结束发现第1对与第2对胶结颗粒对[468,1214]存在ID编号相同的颗粒ID468、第1对与第3对胶结颗粒对[469,1214]存在ID编号相同的颗粒ID469，将第1对胶结颗粒对[468,469]与上述第2对胶结颗粒对[468,1214]、第3对胶结颗粒对[469,1214]从表A中删除并将它们的ID编号按由小到大的顺序列于表B中，表B为[468,468,469,469,1214,1214]，将表B中各颗粒的ID编号与表A中的胶结颗粒对中各颗粒的ID编号进行比较，结果发现表B中各颗粒与表A中的胶结颗粒对中没有相同ID编号的颗粒，此时停止比较，删除表B中重复的ID编号得到表C，表C为[468,469,1214]，由于表C中ID编号的数量只有3个，即表C中只有3个颗粒，未达到10个，不属于岩体破碎关键块体。

(Ⅲ)经过多次筛选后，表A中的第1对胶结颗粒变为[3709,3710]，将该胶结颗粒对[3709,3710]中两个颗粒的ID编号依次与表A中的第2对～最后一对胶结颗粒对中各颗粒的ID编号进行比较，比较结束发现第1对胶结颗粒对与以下5对胶结颗粒对存在ID编号相同的颗粒，即[3709,3804]，[3710,3711]，[3710,3808]，[3710,4399]，[3710,4401]，将第1对胶结颗粒对[3709,3710]与上述5对胶结颗粒对从表A中删除并将它们的ID编号按由小到大的顺序列于表B中，表B为[3709,3709,3710,3710,3710,3710,3710,3711,3804,3808,4399,4401]；

将表B中各颗粒的ID编号依次与表A中的胶结颗粒对中各颗粒的ID编号进行比较，结果发现表A中有以下9对胶结颗粒对[3711,3807]，[3711,3808]，[3807,3808]，[3808,4402]，[3808,4497]，[3808,4499]，[4399,4402]，[4399,4404]，[4401,4404]具有与表B中各颗粒ID编号相同的胶结颗粒，将前述9对胶结颗粒对从表A中删除并将其ID编号列于表B中并将表B中所有ID编号按由小到大的顺序排序，此时表B变为[3709,3709,3710,3710,3710,3710,3710,3711,3711,3711,3804,3807,3807,3808,3808,3808,3808,3808,3808,4399,4399,4399,4401,4401,4402,4402,4404,4404,4497,4499]；

将表B中各颗粒的ID编号依次与表A中的胶结颗粒对中各颗粒的ID编号进行比较，结果发现表A中有以下5对胶结颗粒对[3807,3811]，[4402,4498]，[4407,4497]，[4497,4498]，[4497,4499]具有与表B中各颗粒ID编号相同的胶结颗粒，将前述5对胶结颗粒对从表A中删除并将其ID编号列于表B中并将表B中所有ID编号按由小到大的顺序排序，此时表B变为[3709,3709,3710,3710,3710,3710,3710,3711,3711,3711,3804,3807,3807,3807,3808,3808,3808,3808,3808,3808,3811，4399,4399,4399,4401,4401,4402,4402,4402,4404,4404,4407,4497,4497,4497,4497,4498,4498,4499,4499]；

将表B中各颗粒的ID编号依次与表A中的胶结颗粒对中各颗粒的ID编号进行比较，结果发现表A中有以下1对胶结颗粒对[3811,4499]具有与表B中各颗粒ID编号相同的胶结颗粒，将前述1对胶结颗粒对从表A中删除并将其ID编号列于表B中并将表B中所有ID编号按由小到大的顺序排序，此时表B变为[3709,3709,3710,3710,3710,3710,3710,3711,3711,3711,3804,3807,3807,3807,3808,3808,3808,3808,3808,3808,3811,3811,4399,4399,4399,4401,4401,4402,4402,4402,4404,4404,4407,4497,4497,4497,4497,4498,4498,4499,4499,4499]；

将表B中各颗粒的ID编号依次与表A中的胶结颗粒对中各颗粒的ID编号进行比较，结果发现表B中各颗粒与表A中的胶结颗粒对中没有相同ID编号的颗粒，此时停止比较，删除表B中重复的ID编号得到表C，表C为[3709,3710,3711,3804,3807,3808,3811,4399,4401,4402,4404,4407,4497,4498,4499]，表C中ID编号的数量为15个，即表C中有15个颗粒，属于岩体破碎关键块体。

②按照步骤①的操作，依次对经过上一步骤处理的表A进行处理，直到表A中无剩余胶结颗粒对时结束筛选，共筛选出11个关键块体，分别为：

关键块体1[1204,1718,1722,1810,1811,1813,1814,1815,1816,1820,1822,1824,1895,1896,1897,1898,1901,2419,2420,2421,2423,2426,2427,2428,2429,2519,2521,2522,2523,2524,2525,2526,2527,2528,2529,2530,2531,2532,2533,2534,2535,2605,2606,2607,2608,2609,2610,2611,3035,3037,3128,3129,3130,3132,3133,3137,3138,3142,3234,3235,3236,3237,3238,3239,3240,3241,3242,3244,3312,3314,3315,3316,3317,3318,3319,3321,3946,4025]；

关键块体2[1733,1838,1839,1840,1841,1842,1843,1850,1851,1906,1909,1915,2360,2367,2368,2369,2370,2436,2439,2440,2441,2443,2445,2447,2450,2452,2453,2454,2455,2456,2457,2458,2459,2460,2542,2543,2544,2545,2546,2547,2548,2549,2550]；

关键块体3[1885,2300,2319,2394,2398,2399,2401,2402,2403,2404,2406,2407,2408,2410,2411,2412,2414,2480,2481,2488,2490,2491,2494,2496,2497,2501,2502,2503,2504,2505,2506,2507,2508,2509,2512,2514,2515,2516,2586,2588,2589]；

关键块体4[3080,3082,3182,3281,3707,3791,3792,3793,3794,3795,3796,3798,3799,3800,3801,3891,3892,3894,3990,3992,3994,4051,4052,4487,4488,4491,4492,4493,4494,4582,4585,4588,4589,4590,4591,4592,4593,4594,4595,4687,4688,4689,4690,4691,4693,4694,4695,4753,4755,4756]；

关键块体5[3709,3710,3711,3804,3807,3808,3811,4399,4401,4402,4404,4407,4497,4498,4499]；

关键块体6[3850,3936,4532,4533,4534,4633,4635,4636,4637,4638,4640,4721,5342,5343,5345,5346,5347,5435]；

关键块体7[3856,4445,4459,4460,4461,4463,4464,4537,4538,4539,4541,4544,4545,4554,4641,4642,4643,4644,4645,4646,4648,4653,4654,4657,4724,4726,5081,5082,5084,5086,5087,5091,5097,5099,5101,5166,5170,5171,5172,5173,5175,5176,5178,5179,5181,5183,5184,5185,5186,5187,5188,5189,5190,5197]；

关键块体8[4921,4922,4923,4924,4925,4927,4928,4933,5017,5018,5020,5115,5116,5117,5120,5127,5131,5133,5228,5540,5542,5543,5544,5545,5548,5556,5628,5629,5630,5631,5632,5633,5634,5635,5636,5637,5638,5639,5640,5641,5642,5643,5644,5721,5722,5723,5724,5725,5726,5727,5728,5729,5730,5731,5732,5733]；

关键块体9[5757,5763,5854,5855,5856,5857,5858,5860,5861,5862,5863,5864,5866,5867,5868,5870,5942,5943,5944,5945,5946,5947,5948,5949,5950,5951,5952,5953,5954,5955,5956,5957,5958,5959,6052,6058,6059,6060,6062,6063,6143,6460,6461]；

关键块体10[6278,6291,6293,6381,6382,6383,6384,6386,6388,6389,6394,6888,6891,6899,6959,6965,6966,6968,6971,6973,6974,6975,6978,6979,6980,6981,6982,6983,6985,6987,7054,7060,7061,7062,7064,7066,7067,7068,7069,7070,7071,7072,7073,7074,7075,7076,7077,7078,7079,7080,7081,7082]；

关键块体11[6706,6707,6708,6709,6799,6800,6801,6802,6865,6866,6867,7402,7404,7406,7502,7504,7505,7506,7508,7567,7568,7569]。

③根据筛选出的上述11个关键块体中各颗粒的ID编号和在目标岩体中的三维坐标，确定出各关键块体在目标岩体中的位置，即完成岩体破碎关键块体的定位。

根据筛选出的上述11个关键块体中各颗粒的ID编号和空间坐标，可得到上述关键块体1～关键块体11的形态分别依次如图3(a)～(k)所示，关键块体1～关键块体11在目标岩体中分分布示意图如图4所示，其中关键块体5和关键块体10位于目标岩体内部，从图4中无法看到，图4中，a为关键块体1、b为关键块体2、c为关键块体3、d为关键块体4、f为关键块体6、g为关键块体7、h为关键块体8、i为关键块体9、k为关键块体11。