Численное моделирование процесса сдвижения породных массивов применительно к горно-геологическим условиям шахты имени рубана в Кузбассе

УДК 622.831

Ю.Ю. Кутепов, Е.Б. Боргер

ЧИСЛЕННОЕ

МОДЕЛИРОВАНИЕ ПРОЦЕССА СДВИЖЕНИЯ ПОРОДНЫХ МАССИВОВ ПРИМЕНИТЕЛЬНО К ГОРНО-ГЕОЛОГИЧЕСКИМ УСЛОВИЯМ ШАХТЫ ИМЕНИ РУБАНА В КУЗБАССЕ

Выполнен анализ различных методических подходов прогнозирования процессов сдвижения при подземной разработке угольных месторождений. Рассмотрены горно-геологические условия участка шахты им. Рубана в Кузбассе, характеризующиеся последовательной отработкой свиты мощных наклонных пластов при наличии в разрезе значительной по мощности толщи неоген-четвертичных пород. Предложено для прогноза деформаций подрабатываемого массива горных пород и земной поверхности использовать численное моделирование МКЭ и программный комплекс Abaqus. Разработана методика численного моделирования процесса сдвижения, а также описана последовательность построения модели. Произведены геомеханические расчеты напряженно-деформированного состояния подработанного массива и оседания поверхности при последовательной отработке двух лав по одному пласту, а также пяти и трех лав по двум пластам. С целью проверки полученных при моделировании результатов был выполнен расчет параметров сдвижения по «Правилам охраны...». Выполнено сравнение результатов сдвижения при использовании эмпирического метода и моделирования. Для прогнозирования мест провалов на подработанной поверхности земли применено численное моделирование с заданием дополнительного параметра «Tension Cutoff» при описании пластических свойств в модели Кулона-Мора. Ключевые слова: геомеханика, прогноз, сдвижение, подземная разработка, угольные месторождения, численное моделирование, метод конечных элементов, провалы, Кузбасс.

Введение

В настоящее время прогноз сдвижения породного массива при подземной разработке угольных месторождений осуществляется с использованием методических подходов, в которых угловые параметры процесса определяются методом

ISSN 0236-1493. Горный информационно-аналитический бюллетень. 2017. № 5. С. 66-75. © 2017. Ю.Ю. Кутепов, Е.Б. Боргер.

аналогий или на основании эмпирических зависимостей, полученных путем обобщения опыта инструментальных наблюдений на шахтах с различными горно-геологическими условиями. Данные зависимости приведены в «Правилах охраны сооружений и природных объектов от вредного влияния подземных разработках на угольных месторождений» — отраслевом нормативно-методическом документе РФ [1]. В зарубежной практике широкое применение получили полуэмпирические методы прогноза, основанные на различных «функциях влияния», каждая из которых получена, исходя из известных теоретических зависимостей, а коэффициенты определяются по данным натурных наблюдений процессов сдвижения. Зная значения коэффициентов функции, можно получить оседания в любой точке мульды сдвижения путем интегрирования функции по площади. Наибольшее распространение получили тригонометрическая функция Бальса [2], степенная функция Бейера [3], экспоненциальная функция Занна и Эрхардта-Зауэра [4].

Численное моделирование методом конечных элементов (МКЭ) для решения краевых задач механики сплошной среды начало применяться в горном деле с середины 50-х годов прошлого столетия. Использование МКЭ для прогноза сдвижений массива горных пород с целью охраны сооружений от подработки отмечаются в работах Л. Мюллера [5], Г. Кратча [6], А.Б. Фадеева [7] и др. Однако широкое применение данного метода для прогноза сдвижения пород сдерживалось возможностями вычислительной техники. Технологический прогресс и появление современных программных комплексов, таких, например, как Abaqus, существенно расширило возможности использования МКЭ. В настоящее время МКЭ широко применяется для решения различных задач геомеханики [8—10], и в том числе для прогноза процесса сдвижения [11—14]. В данной работе выполнены прогнозы сдвижения подрабатываемых массивов и поверхности земли на шахте им. Рубана в Кузбассе.

Характеристика объекта исследований

Объектом исследования являются горно-геологические условия участка Красноярский шахты им. А.Д. Рубана в г. Ленинск-Кузнецке, где к настоящему моменту отработанны два пласта — Байкаимский, мощностью 2,6 м на глубинах 200—330 м, и По-лысаевский-2, мощностью 4,9 м на глубинах 100—180 м. При этом участок можно условно разделить на зоны однократной и зона двукратной подработки, в каждой из которых выполнено

моделирование процесса сдвижения от отработки: двух смежных лав по пл. Байкаимский; восьми лав по пластам Байкаим-ский и Полысаевский-2. Мощность четвертичных отложений (наносов) в этом районе достигает 50 м.

Описание модели

Моделирование проводилось в программном комплексе Aba-qus. Задача решалась в плоско-деформационной постановке в рамках механики сплошной среды. Использовался явный метод решения дифференциальных уравнений, реализованный в программе в решателе «Abaqus Dynamic Explicit». Для описания физико-механического поведения вмещающих пород, а также наносов использовалась упругопластическая модель Кулона-Мора.

Моделирование процесса сдвижения лав проводилось последовательно в несколько расчетных шагов. На первом шаге создавалось поле начальных напряжений и вводилась гравитация, при этом все возможные деформации запрещались граничными условиями. На втором и последующих шагах поочередно формировались полости, моделирующие отработку лавы. Под действием гравитационных сил массив начинал деформироваться, заполняя образовавшееся пространство, деформации достигали через какое-то время земной поверхности, образовывая на ней характерную мульду сдвижения. Время каждого шага подбиралось таким образом, чтобы к его концу все деформации в массиве реализовались, и система полностью переходила в равновесное состояние, что соответствовало периоду полных деформаций сдвижения. В конечном итоге получалась картина совместных деформаций породного массива от последовательной отработки всех очистных лав.

Численная модель принципиально разделялась на две области. В середине находилась основная расчетная область прямоугольной формы, в которой происходили все деформационные процессы. Расчетная область разбивалась конечно-элементной сеткой из четырехугольных элементов первого порядка. Вокруг нее создавалась область из бесконечных элементов. Эти элементы исполняли роль граничных условий, основным назначением которых являлось предотвращение отражения деформационных волн от границ модели. Использование такой техники необходимо при решении задач в динамической постановке, для того чтобы избежать колебательных процессов в системе, которые значительно искажают результат моделирования. По-

мимо этого, в свойства материала закладывается коэффициент демпфирования для погашения энергии удара в момент опускания кровли выработки, что так же способствовало повышению точности расчетов.

Моделирования сдвижения от двух смежных лав

по одному пласту

Выбранное расчетное сечение расположено в зоне однократной подработки и включает в себя две очистные лавы (№ 1306 и № 1308) по пласту Байкаимский (т = 2,6 м). Лавы были отработаны с интервалом в один год в восходящем порядке. Средняя глубина первой лавы (№ 1308) составляет 180 м (с учетом наносов), второй (№ 1306) — 140 м.

Маркшейдерских станций на этом участке на период отработки лав не закладывалось, поэтому с целью проверки полученных при моделировании результатов был выполнен расчет параметров сдвижения по рекомендованной в РФ методике [1]. Согласно данному прогнозу величина максимального оседания составили для лав № 1308 и № 1306 соответственно 1,94 и 2,13 м. По результатам расчетов построена деформационная картина от одной и двух лав (соответственно серая и черная пунктирные линии) (рис. 1). При моделировании также получены деформационные картины от одной и двух лав, на рис. 1 обозначенных пунктирными линиями соответственно серого и черного цветов.

Величина максимального оседания для лавы № 1308 по результатам моделирования составила 1,67 м, что на 14% отлича-

Координата по оси X, м

Мульда сдвижения расчитанная по ПБ 07-269-98: --для первой лавы > - -для двух лав Мульда сдвижения полученная в АЬадш: —для первой лавы ; —для двух лав

Рис. 1. Сравнение результатов, полученных при моделировании и по правилам охраны

ется от рассчитанной по правилам безопасности (1,94 м). На мульде сдвижения земной поверхности от двух лав получены две точки, соответствующие локальным максимальным оседаниям. Для лавы № 1308 значение максимального оседания составило 1,7 м, а для лавы № 1306 — 2,2 м. При этом следует отметить, что точка максимального оседания от первой лавы при отработке второй сместилась в ее сторону. Величина максимального оседания от лавы № 1306 согласно расчетам по [1] составляет 2,13, что всего на 3% отличается от значения, полученного при моделировании. Отсутствие плавности кривых по результатам моделирования свиде-тельствует о том, что полностью влияния динамических процессов на результаты моделировании при использовании бесконечных элементов и демпфирования избежать не удалось. На количественный результат данное обстоятельство не сказывается, о чем свидетельствуют расположение соответствующих кривых на графиках, которые находятся достаточно близко друг к другу, поэтому полученные результаты моделирования соответствуют реальным представлениям о механизме процесса сдвижения.

Положительным моментом численного моделирования МКЭ является то, что оно позволяет прогнозировать изменение напряженно-деформированного состояния в ранее подработанном массиве, т.е. отследить трансформацию мульды сдвижения при повторной и последующих подработках массива. Такой подход обеспечивает решение более сложных геомеханических задач, например, прогнозировать напряженно-деформированного состояния пород и земной поверхности при отработке свиты угольных пластов различными лавами в любой последовательности.

Моделирование сдвижения по двум пластам

На втором расчетном сечении рассматривалось отработка пяти лав по пл. Байкаимский (нижележащий, т = 2,6 м) и трех лав по пл. Полысаевский-2 (т = 4,9 м). Отработка лав при этом осуществлялась сначала по нижнему пласту в восходящем порядке, а затем по верхнему пласту по такому же принципу. Особенностью условий ведения горных работ по нижнему и верхнему пластам является расположение лав строго друг под другом.

Идеология моделирования сдвижения при последовательной отработке нескольких пластов принципиально не отличается от случая работы по одному пласту. Однако методически можно отметить ряд отличительных особенностей. Во-первых,

Рис. 2. Схема расположения очистных выработок, а также полученные при моделировании мульды сдвижения при отработке пл. Байкаимский и при отработке двух пластов

в данном случае используется большее количество расчетных шагов, а именно девять. Во-вторых, в связи с тем, что сначала осуществляется отработка нижележащего пласта, то при моделировании отработки вышележащего пласта деформации должны складываться, поэтому формирование полостей в модели должно происходить уже после того, как произошел процесс сдвижения от нижних выработок. Это было достигнуто путем удаления из модели элементов, находящихся на месте выработок, на определенном расчетном шаге. Таким образом, полость формируется в уже подработанном массиве.

На рис. 2 приведены две мульды сдвижения, полученные соответственно при отработке только нижнего пласта и последовательной отработки обоих. Величина максимального оседания по результатам моделирования получилась 6,52 м, над серединами лав № 804 и № 1310.

Прогнозирование провалов

Специфической особенностью рассматриваемого шахтного поля являлось образование на поверхности земли негативных горногеологических процессов — провалов, которые развивались здесь практически повсеместно, но в большинстве своем на подрабатываемых территориях при интенсификации отработки и переходу от выемки угля по одному к двум пластам. Опыт подземной отработки угольных месторождений свидетельствует, что провалы происходят либо при крутом падении пласта, либо при небольшой глубине залегания угольных пластов (Н), когда выполняется соотношение Н < 20m (т — вынимаемая мощность пласта). Причина данных явлений, по на-

шему мнению, кроется в особенностях генезиса, строения, состояния и свойств пород мощной неоген-четвертичной толщи.

Инженерно-геологические исследования покровных неоген-четвертичных пород в районе г. Ленинск-Кузнецка показывают их повышенную мощность относительно других районов Кузбасса, превышающую в отдельных случаях 70 м. Кроме того, они характеризуются неоднородным строением, выраженным в чередовании слоев различного генезиса — делювиально-про-лювиального и аллювиального. Верхняя часть разреза сложена лессовидными отложениями. Выделенные генетические разновидности отличаются по вещественному составу, физическому состоянию, физико-механическим свойствам и, как следствие, деформационному поведению под влиянием процессов сдвижения. Особенно контрастируют с другими разновидностями лессовидные суглинки, характеризующиеся наличием рыхлого скелета с жесткими структурными связями, которые при замачивании водой и нарушении каркаса дают резкую просадку, а под действием растягивающих и сдвиговых деформаций образуют разрывы и блоки. Поэтому в отличии классического представления о механизме сдвижения горных пород, когда верхние слои — наносы плавно повторяет деформации, возникшие при развитии процесса сдвижения в коренных углевмещающих породах. По всей видимости, данный постулат не совсем верен для слоистой анизотропной толщи неоген-четвертичных пород повышенной мощности.

Для конкретных условий шахты осуществлена попытка выполнить численное моделирование напряженно-деформированного состояния подрабатываемого массива, при этом в верхней части неоген-четвертичной толщи лессовидные породы харак-

РЕЕС2Т

+2.907е-02

Зоны возможного провалообразования

--+7.091е-03

--+5.б73е-03

--+4.254е-03

--+2.836е-03

--+1.418е-03

—I- +0.000е+00

Лава №1308

Лава

у

♦ СЮВ: ЯиЬап_1_тр.оаЬ АЬачиэ/ЕхрНсК 30ЕХРЕМЕ1ЧСЕ 1*2016х НйНх 2 Рг1 1<оу 18 15:10:36 СМТ+03:00 X б1ер: 1лу1306

1псгетег^ 14023: Step 71 те = 5.000 _Штагу Уаг: РЕЕОТ_

1.

Рис. 3. Зоны возможного провалообразования

теризовались дополнительным параметром «Tension Cutoff», позволяющим описывать их пластические свойства в модели Кулона-Мора. Это необходимо для того, чтобы ограничить прочность материала в зонах растягивающих напряжений. Значение данного параметра принималось соответствующие структурной прочности, которая для лессовидных пород принималась 5-15 кПа.

На рис. 3 приведены результаты расчетов относительных пластических деформаций растяжения в верхней части неоген-четвертичных отложениях. Более насыщенным цветом выделены зоны повышенных деформаций — провалов. Можно отметить, что в основном эти зоны формируются в местах с наибольшей кривизной мульды, что позволяет нам сделать вывод о том, что данное явление связанно с процессом сдвижения и особенностью инженерно-геологического строения покровных отложений.

Заключение

Рассмотрены процессы сдвижения горных пород при подземной разработки угля и выполнен анализ различных методических подходов для их прогнозирования. Применительно к условиям шахты им. Рубана в Кузбассе выполнен прогноз оседания земной поверхности при отработки двух лав по одному пласту, а также пяти и трех лав по двум пластам. Установлена достаточно хорошая сходимость результатов сдвижения при использовании эмпирического метода и моделирования. Предложено для прогнозирования мест провалов на участках применять численное моделирование с использованием дополнительного параметра «Tension Cutoff» для получения зон повышенных растягивающих напряжений.

СПИСОК ЛИТЕРАТУРЫ

1. Правила охраны сооружений и природных объектов от вредного влияния подземных горных разработок на угольных месторождениях. — М.: Недра, 1981. — 288 с.

2. Bals R. Beitrag zur Frage der Vorausberechnung bergbaulicher Senkungen // Markscheidew. Vol. 42. 1931. pp. 98—111.

3. Beyer F. On predicting ground deformation due to mining flat seams. Habilitation Thesis. 1945. Tech. Univ., Berlin.

4. Ehrhardt W, Sauer A. Die Vorausberechnung von Senkung, Schieflage und Kriimmung iiber dem Abbau in flacher Lagerung // Bergb.-Wiss. 1961.

5. MullerL. Der Felsbau, Vol. 1, Stuttgart: Enke 1963.

6. Кратч Г. Сдвижение горных пород и защита подрабатываемых сооружений. — М.: Недра, 1978. — 494 с.

7. Фадеев А. Б. Метод конечных элементов в геомеханике. — М.: Недра, 1987. - 221 с.

8. Кутепов Ю. И., Кутепова Н. А., Карасев М. А., Кутепов Ю. Ю. Прогноз формоизменения намывных массивов гидроотвалов при складировании на них отвальных насыпей // Горный журнал. — 2016. — № 12. — С. 23—27.

9. Кутепов Ю. Ю., Карасев М. А. Оценка устойчивости кровли подземного сооружения // Маркшейдерия и недропользование. — 2015. — № 5 (79). — С. 43—46.

10. Кутепов Ю. И., Кутепова Н. А., Карасев М. А., Фоменко Н. Г. Геомеханическое обоснование отсыпки отвалов «сухих» пород на гидроотвалах // Геоэкология, инженерная геология, гидрогеология, геокриология. — 2015. — № 3. — С. 220—225.

11. Сахно И.Г., Грищенков Н.Н., Голубев Ф.М. Математическое моделирование сдвижений горного массива над одиночной лавой // Науковi пращ УкрНДМ1 НАН Украши. — 2013. — № 13. — С. 209—219.

12. Vyazmensky A., Elmo D, Stead D, Rance J. Combined finite-discrete element modelling of surface subsidence associated with block caving mining // Proceedings of 1st Canada—US rock mechanics symposium. 2007. — pp. 467—475.

13. Nuric A., Nuric S, KricakL., Lapandic I., Husagic R. Numerical Modeling and Computer Simulation of Ground Movement Above Underground Mine // International Journal of Environmental, Chemical, Ecological, Geological and Geophysical Engineering. — 2012. — Vol. 6. — no 2. — pp. 579—587.

14. Dasa R., Singh P. K. Numerical analysis of surface subsidence in asymmetric parallel highway tunnels // Journal of Rock Mechanics and Geotech-nical Engineering. 2017.

КОРОТКО ОБ АВТОРАХ

Кутепов Юрий Юрьевич — аспирант, e-mail: Jura992@gmail.com, Санкт-Петербургский горный университет,

Боргер Елена Борисовна — главный маркшейдер шахтоуправления, АО «СУЭК-Кузбасс», e-mail: BorgerEB@suek.ru.

UDC 622.831

Gornyy informatsionno-analiticheskiy byulleten'. 2017. No. 5, pp. 66-75. Yu.Yu. Kutepov, E.B. Borger

NUMERICAL MODELING OF THE ROCK MASS SUBSIDENCE APPLIED TO GEOLOGICAL CONDITIONS OF THE MINE NAMED AFTER RUBAN IN KUZBASS

The analysis of different prediction methods of subsidence due to longwall coal mining has been done. The analysis of the geotechnical conditions has been done for the site of the coal mine named after Ruban in Kuzbass where there is a suite of thick coal seams and the productive aeries of sandstone and siltstone of suite Kolchuginskaya.Top layers of Neogene-Quaternary overburden soils represent itself subsiding structurally unstable loess loam. It is

proposed to use numerical simulation with FEM for the prediction of deformations of undermined rock mass and the ground surface. Abaqus software has been used for this purpose. The geomechanical calculations of subsidence have been done for two longwalls in one coal seam as well as for five and three longwalls in two coal seams. In order to verify the results obtained in the simulation the calculation of subsidence parameters on the basis of «Rules of protection ...» has been done also. The comparison of subsidence results gained by empirical methods and modeling is provided. For the prediction of failure zones on the surface it is proposed to use finite element modeling by using an additional parameter «Tension Cutoff» describing properties in Mohr-Coulomb plasticity model.

Key words: geomechanics, prediction, subsidence, underground mining, coal deposits, numerical simulation, finite element method, surface failures, Kuzbass.

AUTHORS

Kutepov Yu.Yu., Graduate Student, e-mail: Jura992@gmail.com, Saint Petersburg Mining University, 199106, Saint-Petersburg, Russia, Borger E.B., Chief Surveyor of Mine Administration, e-mail: BorgerEB@suek.ru, JSC «SUEK-Kuzbass», Leninsk-Kuznetsk, Russia.

REFERENCES

1. Pravila okhrany sooruzheniy iprirodnykh ob"ektov ot vrednogo vliyaniya podzemnykh gornykh razrabotok na ugol'nykh mestorozhdeniyakh (Protection Regulations of constructions and natural objects from harmful influence of underground mining in coal deposits), Moscow, Nedra, 1981, 288 p.

2. Bals R. Beitrag zur Frage der Vforausberechnung bergbaulicher Senkungen. Marksc-heidew. Vol. 42. 1931. pp. 98-111.

3. Beyer F. On predicting ground deformation due to mining flat seams. Habilitation Thesis. 1945. Tech. Univ., Berlin.

4. Ehrhardt W., Sauer A. Die Vorausberechnung von Senkung, Schieflage und Kriimmung iiber dem Abbau in flacher Lagerung. Bergb.-Wiss. 1961.

5. Muller L. DerFelsbau, Vol. 1, Stuttgart: Enke 1963.

6. Kratch G. Sdvizhenie gornykh porod i zashchita podrabatyvaemykh sooruzheniy (Mining subsidence engineering), Moscow, Nedra, 1978, 494 p.

7. Fadeev A. B. Metod konechnykh elementov v geomekhanike (Finite element method in geomechanics), Moscow, Nedra, 1987, 221 p.

8. Kutepov Yu. I., Kutepova N. A., Karasev M. A., Kutepov Yu. Yu. Gornyy zhurnal. 2016, no 12, pp. 23-27.

9. Kutepov Yu. Yu., Karasev M. A. Marksheyderiya inedropol'zovanie. 2015, no 5 (79), pp. 43-46.

10. Kutepov Yu. I., Kutepova N. A., Karasev M. A., Fomenko N. G. Geoekologiya, in-zhenernaya geologiya, gidrogeologiya, geokriologiya. 2015, no 3, pp. 220-225.

11. Sakhno I. G., Grishchenkov N. N., Golubev F. M. Naukovipratsi UkrNDMINAN Ukraini. 2013, no 13, pp. 209-219.

12. Vyazmensky A., Elmo D., Stead D., Rance J. Combined finite-discrete element modelling of surface subsidence associated with block caving mining. Proceedings of 1st Canada—US rock mechanics symposium. 2007. pp. 467-475.

13. Nuric A., Nuric S., Kricak L., Lapandic I., Husagic R. Numerical Modeling and Computer Simulation of Ground Movement Above Underground Mine. International Journal of Environmental, Chemical, Ecological, Geological and Geophysical Engineering. 2012. Vol. 6. no 2. pp. 579-587.

14. Dasa R., Singh P. K. Numerical analysis of surface subsidence in asymmetric parallel highway tunnels. Journal of Rock Mechanics and Geotechnical Engineering. 2017.