Геоэкология. Инженерная геология, гидрогеология, геокриология, 2022, № 1, стр. 69-76

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

Т. Г. Макеева^1, *, В. А. Трофимов^2, **

¹ Национальный исследовательский Московский государственный строительный университет (НИУ МГСУ)
129337 Москва, Ярославское шоссе, д. 26, Россия

² Институт проблем комплексного освоения недр им. Н.В. Мельникова РАН
111020 Москва, Крюковский тупик, д. 4, Россия

^* E-mail: makeeva13new@yandex.ru
^** E-mail: asas_2001@mail.ru

Поступила в редакцию 19.10.2021
После доработки 19.11.2021
Принята к публикации 22.11.2021

EDN: LGIJSZ
DOI: 10.31857/S086978092201009X

Полный текст (PDF)

Аннотация

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

Ключевые слова: мульда оседания, метод граничных элементов, основная кровля, выработанное пространство, двухслойная геомеханическая модель

ВВЕДЕНИЕ

Создание подземных выработок всегда сопровождается оседанием земной поверхности, величина которого может достигать 1000 и более миллиметров, в зависимости от размеров выработанного пространства, строения налегающей толщи горных пород и их деформационно-прочностных свойств. Для обеспечения безопасности и безаварийной работы многих наземных сооружений и коммуникаций важно знать и уметь прогнозировать возникновение и протекание различных деформационных процессов на земной поверхности, в зависимости от развития подземных горных работ, связанных с расширением выработанного пространства. Это весьма актуальная проблема как с научной, так и практической точек зрения, хотя усилия по ее решению предпринимаются на протяжении многих десятилетий. До сих пор ежегодно появляются публикации на указанную тему практически во всех высокорейтинговых специализированных журналах [3–6, 9–12]. Востребованность решения этой проблемы связана с практическими нуждами даже в крупных мегаполисах при создании значительных пустот под землей (скажем при постройке метрополитена), которые являются причиной потери устойчивости и разрушения зданий.

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

Предпринимались многочисленные попытки применения разработанных подходов для анализа состояния и поведения массива при наличии выработок неглубокого заложения. Однако при сопоставлении результатов расчетов с реально замеренными величинами оседания на конкретных горнотехнических объектах были выявлены определенные особенности.

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

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

В связи с этим вопросы описания деформирования земной поверхности рассматриваются в рамках другого научного направления – маркшейдерии. В этом направлении развиты многочисленные подходы, достаточно приемлемо описывающие особенности оседания земной поверхности, ее наклонов, растяжения–сжатия и др. Отметим, однако, что все полученные результаты в этом направлении основаны на обобщении многочисленных результатов наблюдений на шахтах и рудниках [7, 8].

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

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

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

МЕТОДЫ И РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЯ

Рассмотрим типичный пример (характерный для угольных месторождений) геологического строения массива горных пород, когда в нем существует монолитный слой достаточно прочных песчаников мощностью 10 м и более с модулем упругости в E_П ~0.5 ⋅ 10⁵ МПа. При этом усредненный “эффективный” модуль упругости остальной слоистой толщи налегающих пород значительно меньше (на порядок и более). Если угольный пласт оказывается непосредственно под этим слоем, то последний будет выполнять роль кровли при отработке пласта и, в определенной степени, защищать его выработки от разрушения. Опыт отработки угольных пластов показывает, что такая ситуация весьма широко распространена в практике.

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

Двухслойная модель массива

Для моделирования состояния и поведения геомеханической системы, включающей мощный слой песчаника и расположенные под ним выработки, будем использовать двухслойную модель массива горных пород. Упомянутый слой песчаника мощностью h_сл образует нижний упругий слой модели с модулем упругости E_П, который фактически несет всю нагрузку от веса пород. В то же время второй слой модели – верхний (от кровли нижнего слоя и вплоть до поверхности), как уже говорилось, имеет малый “усредненный” модуль упругости E $ \ll $ E_П. Физически это означает, что он подвергся выветриванию и по большей части разрушен, в связи с чем играет незначительную роль в процессе деформирования. Фактически им можно пренебречь, а его вес считать пассивной нагрузкой, действующей на нижний слой. Иными словами, это слой несжимаемого песка, деформация поверхности которого (фактически это земная поверхность) совпадает с деформацией кровли нижнего слоя. Постановка задачи для численного расчета методом граничных элементов в этом случае иллюстрируется рис. 1 [1].

Рис. 1.

Схематическое изображение постановки задачи для решения методом граничных элементов.

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

На земной поверхности необходимо смоделировать нулевые нормальные напряжения, в связи с чем задается симметричная относительно некоторой произвольной горизонтальной прямой Y₀ постановка.

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

В полных напряжениях на поверхности выработки будем иметь:

(1)

${{\sigma }_{n}} = 0,\quad {{\tau }_{n}} = 0,\quad {\text{|}}x{\text{|}} < a,\quad y = H.$

На верхней же поверхности введенного слоя должны выполняться соотношения:

(2)

$\begin{gathered} {{\sigma }_{n}} = - \gamma (H - {{h}_{{{\text{сл}}}}}),\quad {{\tau }_{n}} = 0, \\ {\text{|}}x{\text{|}} < \infty ,\quad y = - (H - {{h}_{{{\text{сл}}}}}), \\ \end{gathered} $

где σ_n, τ_n – нормальные и касательные напряжения на границах расчетной области, γ – средний удельный вес пород кровли, h_сл – мощность основной кровли.

Задача решается методом граничных элементов с использованием алгоритма разрывных смещений.

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

Обобщенная зависимость оседания земной поверхности

При этом представляет интерес не проведение отдельных расчетов с заданными значениями параметров H, a, h_сл, а построение некоторых обобщенных зависимостей – конечных соотношений, которые можно было бы использовать для элементарных расчетов и проведения анализа получаемых закономерностей. Это достигается варьированием приведенных параметров и проведением многочисленных расчетов для различных комбинаций этих величин.

В качестве характерного параметра, описывающего оседание земной поверхности, будем считать максимальное оседание в мульде W₀, приуроченное к ее центру. Рисунок 2а иллюстрирует характер зависимости этой величины от безразмерной мощности слоя ${{\bar {h}}_{{{\text{сл}}}}}$ = h_сл/a для разных величин безразмерной глубины $\bar {H}$ = H/a.

Рис. 2.

Максимальное оседание в мульде W₀ в зависимости от безразмерной мощности слоя ${{\bar {h}}_{{{\text{сл}}}}}$ для разных величин безразмерной глубины $\bar {H}$. (а) – линейные, и (б) – логарифмические координаты.

При этом в расчетах величина полупротяженности выработанного пространства принималась равной a = 100 м, а глубины 250, 500, 750 м (кривые 3, 2, 1).

Эти же зависимости в логарифмических координатах отражены на рис. 2б, из которого видно, что они имеют практически прямолинейный характер.

Нормировка всех кривых на рис. 2б на величину $2\alpha {{a}^{2}}\bar {H}$ приводит к единой кривой в координатах $({{\bar {W}}_{0}} = {{W}_{0}}{\text{/}}2\alpha {{a}^{2}}\bar {H},{{\bar {h}}_{{{\text{сл}}}}})$, которая и отражает зависимость безразмерного оседания в центре мульды от безразмерной толщины слоя (рис. 3, сплошная линия А). При этом деформационные свойства пород кровли отображаются комплексным параметром $\alpha = 2(1 - {{\nu }^{2}})\gamma {\text{/}}E$.

Рис. 3.

Обобщенная зависимость максимального безразмерного оседания в мульде ${{\bar {W}}_{0}}$ от безразмерной мощности слоя ${{\bar {h}}_{{{\text{сл}}}}}$.

Здесь же показана возможная аппроксимация полученной кривой посредством прямой линии в двойных логарифмических координатах (пунктирная линия Б) для всего диапазона 0.25 < ${{\bar {h}}_{{{\text{сл}}}}}$ < 5.0.

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

Полученный график А на рис. 3 (либо прямая Б) являются конечным результатом проведенного исследования и позволяет определить максимальное оседание в мульде для тех или иных параметров $\bar {H}$, ${{\bar {h}}_{{{\text{сл}}}}}$, α. В аналитическом виде прямая Б записывается соотношением:

(3)

${{\bar {W}}_{0}} = \frac{{1.498}}{{\bar {h}_{{{\text{сл}}}}^{{1.417}}}}.$

Более точно эту кривую можно аппроксимировать двумя прямыми для диапазонов 0.25 < ${{\bar {h}}_{{{\text{сл}}}}}$ < < 1.5 и 1.5 < ${{\bar {h}}_{{{\text{сл}}}}}$ < 5.0, а именно соотношениями:

(4)

${{\bar {W}}_{0}} = \frac{{0.796}}{{\bar {h}_{{{\text{сл}}}}^{{2.165}}}}\,\,\,{\text{и}}\,\,\,{{\bar {W}}_{0}} = \frac{{1.157}}{{\bar {h}_{{{\text{сл}}}}^{{1.104}}}}.$

Примеры расчета

Приведенные теоретические построения и полученные результаты должны быть в той или иной степени подтверждены экспериментальными данными, полученными in situ. В качестве таковых использовались данные замеров оседаний земной поверхности при разработке угольных месторождений Донбасса¹ ¹, практически вся территория которых покрыта реперными линиями.

Отметим, что использованная выше модель деформирования массива имеет достаточно узкий диапазон применения, а именно:

− должно быть выполнено условие плоской деформации массива;

− показания реперов обусловлены проходкой единственной выработки;

− пологое падение пласта.

Тем не менее, такие условия выполняются на многих шахтах, и для них приведенная модель может быть использована в полной мере.

Рассмотрим несколько конкретных вариантов применения модели.

Угольный пласт мощностью m ≈ 0.82 м с углом падения 10° отрабатывался шахтой Дельта-2 Треста “Ворошиловуголь” одиночной панелью шириной a = 85 м и протяженностью более 300 м на глубине в среднем H = 115 м. При этом в кровле выработки лежит мощный слой песчаников (h_сл = = 25 м). В качестве E_П для слоя была взята средняя табличная величина 0.25 ⋅ 10⁵ МПа. Удельный вес пород кровли –2.5 ⋅ 10⁴ Н/м³.

Используя приведенные значения горнотехнических параметров, вычислим величины характерных параметров модели: $\bar {H}$ = 1.35; ${{\bar {h}}_{{{\text{сл}}}}}$ = = 0.29. Для этих значений в соответствии с графиком рис. 3 будет ${{\bar {W}}_{0}}$ = 11.61, и величина максимального оседания равна ${{W}_{0}} = {{\bar {W}}_{0}}2\alpha {{a}^{2}}\bar {H} \cong 0.434$ м. При этом максимальное оседание поверхности, полученное посредством замеров соответствующей реперной станцией, составило 489 мм. Корректировка в расчетах величины E_П примерно на 10% (E_П = 0.22 ⋅ 10⁵ МПа) приводит к практическому совпадению замеренных и вычисленных результатов. На общем графике рис. 3 этому случаю соответствует точка 1.

Угольный пласт l₁ мощностью m ≈ 1.1 м с углом падения 9° отрабатывался шахтой 5–6 им. Димитрова Треста “Красноармейскуголь” одиночной панелью шириной a = 85 м на глубине в среднем H = 145 м. При этом в кровле выработки лежит мощный слой песчаников (h_сл = 21 м). В качестве E_П для слоя была взята средняя табличная величина 25 ⋅ 10⁵ МПа. Удельный вес пород кровли – 2.5 ⋅ 10⁴ Н/м³.

Используя приведенные значения горнотехнических параметров, вычислим величины характерных параметров модели: $\bar {H}$ = 1.71; ${{\bar {h}}_{{{\text{сл}}}}}$ = 0.25. Для этих значений в соответствии с графиком рис. 3 будет ${{\bar {W}}_{0}}$ = 15.7, и величина максимального оседания равна ${{W}_{0}} = {{\bar {W}}_{0}}2\alpha {{a}^{2}}\bar {H} \cong 0.743$ м.

При этом максимальное оседание поверхности, полученное посредством замеров соответствующей реперной станцией, составило 747 мм (рис. 4). В данном случае наблюдается практически полное совпадение расчетных и экспериментальных данных. На общем графике рис. 3 этому случаю соответствует точка 2.

Рис. 4.

Пример экспериментальных данных (форма мульды оседания) для шахты 5–6 им. Димитрова Треста “Красноармейскуголь”.

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

Аппроксимация формы мульды оседания

Для примера рассмотрим пример из практики, результаты для которого приведены на рис. 4. Измеренные величины на графике показаны точками. При этом точка А в силу тех или иных обстоятельств выпадает из общей закономерности и ее не следует учитывать при дальнейшем анализе. Здесь же прямоугольником показано расположение выработанного пространства (270–350 м). При этом центр мульды оседания расположен со сдвигом в 40 м относительно центра выработанного пространства. Он локализуется при x = 270 м и на рисунке отмечен вертикальным отрезком.

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

Применим разработанный алгоритм к решению описанной задачи, для которой выше была дана оценка максимального оседания в центре мульды, и отобразим форму полученной левой полумульды на рис. 5 (кривая 1).

Рис. 5.

Расчетная аппроксимация оседаний земной поверхности для шахты 5–6 им. Димитрова Треста “Красноармейскуголь”.

При этом для удобства изменим начало отсчета координаты x и будем отсчитывать ее от центра выработанного пространства (на рисунке отображена лишь его левая половина). По-прежнему точками на рисунке обозначены экспериментальные данные, которые лежат на некотором удалении от полученной расчетной кривой 1. Очевидно, что смещение этой расчетной кривой параллельно самой себе на 40 м в положение 2, приводит к вполне приемлемой аппроксимации экспериментальных данных. Подобный сдвиг обусловлен тем, что упругая задача решалась для строго горизонтального пласта, а экспериментальные данные, как уже упоминалось выше, получены для наклонного и, следовательно, имеют сдвиг. Таким образом, модель слоя дает неплохую аппроксимацию не только по величине ${{\bar {W}}_{0}}$, но и по всей мульде оседания.

Аналогичные результаты могут быть получены и по другим случаям оседания земной поверхности, некоторые из которых отражены на рис. 3 точками.

ЗАКЛЮЧЕНИЕ

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

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

При этом величина максимального оседания, как и форма мульды, определяются безразмерными величинами: глубиной заложения выработки, мощностью слоя песчаника в кровле и материальным параметром, соответственно – $\bar {H}$, ${{\bar {h}}_{{{\text{сл}}}}}$, α.

Список литературы

Крауч С., Старфилд А. Методы граничных элементов в механике твердого тела. М.: Мир, 1987. 328 с.
Трубецкой К.Н., Иофис М.А., Кузнецов С.В., Трофимов В.А. Основные закономерности оседания подрабатываемой толщи горных пород и прогиба зависающей кровли на малых и больших глубинах // Физико-технические проблемы разработки полезных ископаемых. 1999. № 3. С. 5–11.
Alejano L.R., Ramirez-Oyanguren P., Taboada J. FDM predictive methodology for subsidence due to flat and inclined coal seam mining // Int. J. of Rock Mechanics and Mining Sciences. 1999. 36 (4). P. 475–491. https://doi.org/10.1016/S0148-9062(99)00022-4
Chuang Liu, Huamin Li & Hani Mitri. Effect of Strata Conditions on Shield Pressure and Surface Subsidence at a Longwall Top Coal Caving Working Face // Rock Mechanics and Rock Engineering. 2019. Vol. 52 (3). P. 1523–1537. https://doi.org/10.1007/s00603-018-1601-3
Hamdi Pooya, Stead Doug, Elmo Davide, Toyra Jimmy. Use of an integrated finite/discrete element method-discrete fracture network approach to characterize surface subsidence associated with sub-level caving // Int. J. of Rock Mechanics and Mining Sciences. 2018. V. 103. P. 55–67. https://doi.org/10.1016/j.ijrmms.2018.01.02
Lizarraga Jose J., Buscamera G. A geospatial model for the analysis of time-dependent land subsidence induced by reservoir depletion // Int. J. of Rock Mechanics and Mining Sciences. 2020. V. 129 (2): 104272. https://doi.org/10.1016/j.ijrmms.2020.104272
Makeeva T., Trofimov V. Forecast of deformations of the land surface from the separate clearing development, displacement and deformation in the main sections of the trough // E3S Web Conf. XXII Int. Scientific Conf. “Construction the Formation of Living Environment” (FORM-2019). 2019. V. 97: 04016. https://doi.org/10.1051/e3sconf/20199704016.
Makeeva T., Trofimov V. Regularities of the day surface deformation during layer mining by consecutive lavas // MATEC Web of Conf. VI Int. Scientific Conf. “Integration, Partnership and Innovation in Construction Science and Education” (IPICSE-2018), 2018. V. 2: 5102013.https://doi.org/10.1051/matecconf/201825102013.
Riesgo PP., Rodríguez F.G., Krzemienb G.A., et al. Subsidence versus natural landslides when dealing with property damage liabilities in underground coal mines // Int. J. of Rock Mechanics and Mining Sciences. 2020. V. 126 (5): 104175. https://doi.org/10.1016/j.ijrmms.2019.104175
Wempen J.M. Application of DInSAR for short period monitoring of initial subsidence due to longwall mining in the mountain west United States // Int. J. of Mining Science and Technology. 2020. V. 30. Is. 1. P. 3–37. https://doi.org/10.1016/j.ijmst.2019.12.011.
Wang Binglong, Xu Jialin, Xuan Davaang. Time function model of dynamic surface subsidence assessment of grout-injected overburden of a coal mine // Int. J. of Rock Mechanics and Mining Sciences. 2018. V. 104. P. 1–8. https://doi.org/10.1016/j.ijrmms.2018.01.044
Yang Xuelin, Wen Guangcai, Dai Linchao, Sun Haitao and Li Xuelong. Ground Subsidence and Surface Cracks Evolution from Shallow-Buried Close-Distance Multi-seam Mining: A Case Study in Bulianta Coal Mine // Rock Mechanics and Rock Engineering. 2019. Vol. 52 (2). P. 2835–2852. https://doi.org/10.1007/s00603-018-1726-4

Дополнительные материалы отсутствуют.

Инструменты

следующая статья выпуска предыдущая статья выпуска содержание выпуска

Геоэкология. Инженерная геология, гидрогеология, геокриология

Архивы выпусков Информация о журнале Отправить рукопись в журнал