Концепция моделирования гидродинамических процессов на шахтных полях

РЕФЕРАТ

Концепция моделирования гидродинамических процессов на шахтных полях

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

Рис. 1. Принципы классификации моделей гидрогазодинамических (ГГД) процессов в ГТС горнопромышленного региона Созданный комплекс моделей предназначен для составления прогнозов нестационарных физико-химических процессов, протекающих в геотехнических системах, а также долговременных последствий, затрагивающих их наиболее подвижные элементы: гидросферу, атмосферу, агро-экосистемы.

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

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

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

Ряд параметров и расчетных характеристик задаются в виде многомерных случайных величин или случайных полей. Для таких параметров в любой точке массива P в момент времени t выделяются закономерная и случайная составляющие

. (1)

Пусть (t) — случайный многомерный процесс, описывающий пространственно-временные изменения ГГД параметров или расчетных характеристик, _i=(t_i) — фактические значения этих параметров в моменты времени t_i, _i,m — измеренные в системе опорных точек — пунктов наблюдения с некоторой погрешностью в те же моменты, _i,c — значения, рассчитанные на основании моделирования. Задача оптимального прогнозирования на основе неполных исходных данных состоит в минимизации математического ожидания погрешностей прогнозируемых величин на интервале времени :

. (2)

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

В табл. 2 приняты следующие обозначения: H — напор (уровень) подземных вод; H_m — уровень шахтных вод; Q — водоприток в шахту; Q_tr — расход перебрасываемых вод из одной шахты в другую; C_gw, C_m и C_gr — концентрации веществ в подземных, шахтных водах и в поровом растворе верхних слоев грунта; q_gw и q_sw - массопотоки в грунтовые и в поверхностные воды; P_w — давление рудничного газа в выработках; Q_g — его поток через массив и систему вентиляции; C_atm — концентрация газообразных примесей и пыли в атмосфере. Обозначения для параметров, задаваемых для расчета, следующие: K — коэффициент фильтрации; - раскрытие трещин, и — азимут и угол их падения; n и V_f — пористость (трещиноватость) и объем трещин в части массива; m — мощность водоносных и водоупорных слоев; H_b — напор на границе области фильтрации; - инфильтрация; D_f и D_g — коэффициенты диффузии водных растворов и газов в горных породах; n_e и — параметры сорбции и нейтрализации (распада) веществ в горных породах; _g — плотность газа в выработках; - водонасыщенность пород, V_w — объем выработанного пространства, d_w — размер выработок; p_g — газоносность пород; q_ex — темпы добычи углесодержащих пород; - концентрация метана в рудничном газе, C_s, q_s — концентрация веществ в шахтных породах и скорость их выщелачивания; u — скорость ветра, _g — параметр осаждения и поглощения примесей из атмосферы на поверхности земли; d_p — гранулометрический состав пылевых выбросов; V — скорость фильтрации.

Таблица 2. Взаимосвязь основных параметров ГГД процессов

В табл. 2. также обозначено: «!» — достоверно определяемый параметр, «?» — распределенный параметр или случайная величина (функция, процесс), «+» — статистически оцениваемый параметр; «-» — определение не обязательно; СД — возможность сопоставления с исходными данными: «Н» — низкая, «С» — средняя, «В» — высокая На этапе построения гидрогазодинамических моделей решаются следующие основные задачи: 1) разбиение массива на структурные элементы (слои, блоки), 2) построение полей проницаемости, пористости, 3) задание пустотности (пространственно-временного распределения выработок по объему), 4) задание инфильтрационных водопритоков, 5) определение граничных условий и контактов с соседними шахтными полями (в том числе переток через сбойки), 6) определение методики расчета водопритока, как на уровне шахтного поля, так и в масштабе отдельных горизонтов и блоков, 7) учет технологий отработки (особенности водоотлива и дренажа, консервации шахты), 8) учет переменной газоносности пород, различных путей газовыделения и миграции газа, 9) стыковка моделей разных масштабов, 10) стыковка гидродинамических, миграционных и газодинамических моделей на границах твердой, жидкой и газообразной сред, 11) учет взаимного влияния процессов в различных фазах.

Моделирование ГГД процессов в подработанном массиве осложняется слоистостью, наличием зон водопроводящих трещин, зон сдвижения, последовательностью отработки. Адекватный учет этих особенностей предполагает пространственную дискретизацию массива. Размер блока дискретизации должен допускать осреднение параметров и их оценку (измерение) в натурных условиях. С этой точки зрения блоки размером менее 50 м по горизонтали и 15−20 м по вертикали нецелесообразны: при этом размеры пустот (выработок) становятся соизмеримыми с размерами блока. Напротив, увеличение блоков, соответствующее региональным моделям, в большей мере соответствует гипотезе сплошности среды, на которой строится модель фильтрации. При этом напор, рассчитанный на крупномасштабных моделях, особенно в массиве с крутопадающими слоями при наличии осложняющих факторов, следует трактовать как средне вероятный в пределах блока.

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

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

а) б) Рис. 2. Изменение коэффициента фильтрации K (м/сут) с глубиной z на полях шахт: а) Александр-Запад, б) им. Калинина; крайние значения; - средние значения; N — число измерений (столбцы)

Анализ графиков свидетельствует о большой вариации K в пределах одного слоя, что характерно для всех глубин. Это может быть связано: 1) с наличием нескольких подзон различной проницаемости, 2) с макронеоднородностью и масштабными эффектами, 3) с хаотической неоднородностью. Значительный разброс данных измерений подтверждает необходимость использования моделей, в которых проницаемость пород является случайным полем или распределенным параметром.

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

Закономерная неоднородность на уровне отдельных блоков может рассматриваться как стохастическая в пределах слоя (горизонта) или его больших подзон. В разработанных моделях задается пространственно-временная изменчивость параметров, связанная с неодновременным проведением горных работ на различных горизонтах отработки.

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

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

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

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

Стыковка моделей миграции в грунтах и подземных водах осуществляется через параметры массопотока на границах «подземные — поверхностные воды» и на границе зоны аэрации. Модели миграции в верхних слоях грунтов и отвалах шахтных пород и воздушной среде сочленяются посредством параметров осаждения твердых частиц и поглощения газообразных примесей на подстилающей поверхности. Принципы моделирования движения пылевой и газовой компонент в геотехнических системах изложены далее в главе 3.

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

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

С целью адекватного отображения фильтрационной неоднородности и пространственного распределения зон отработки шахтного поля производится его разбивка на блоки в виде прямоугольников. На каждом горизонте толщиной z, где ведется отработка, задается послойно общий объем выработанного пространства V_w,l (l=1,… N_z, N_z — число слоев). При схематизации структуры массива на каждом горизонте отработки выделяется несколько зон с отличающимися значениями фильтрационных параметров и объемами выработанного пространства. Тогда, в соответствии с горизонтальным распределением объема V_w по слоям и с учетом последовательности отработки для каждого блока можно задать объем выработок, время их появления и погашения. Исходя из объема V_w,l определяется площадь горизонтального сечения выработок S_h на заданной глубине z как отношение объема пустот в слое к его толщине:.

Пример такой фильтрационной схематизации для шахтного поля в условиях Центрального района Донбасса показан на рис. 3 и 4. Наличие зон усиленной инфильтрации (отстойников шахтных вод, понижений рельефа) моделируется путем задания неравномерной по площади интенсивности инфильтрационного питания .

шахта гидродинамический массив Рис. 3. Распределение выработанного пространства V_w и среднего коэффициента фильтрации K_f по глубине В условиях нарушения естественных геологических структур при разработке пластов, наличия зон обрушения с повышенной вертикальной проницаемостью расчет изменения уровня подземных вод целесообразно вести на основе уравнения плановой неустановившейся фильтрации [5,6]

. (3)

Здесь Н — уровень подземных вод, T_x и T_y — проводимость вдоль осей Ox и Oy соответственно, n_f — трещинная пористость пород массива, Q_w — интенсивность стоков в выработки, — интенсивность инфильтрации.

В сложных гидрогеологических условиях фильтрационный расчет выполняется с помощью численных моделей. Для решения уравнения (3) использована попеременно-треугольная схема метода конечных разностей, сочетающая преимущества явной схемы с хорошей вычислительной устойчивостью [7,8].

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

Рис. 4. Расположение горизонтов отработки и очагов инфильтрации в пределах шахтного поля:1,…, 5 — контуры горизонтов отработки на отметках 1 — от 75 до 210 м; 2 — от -100 до 75 м; 3 — от -300 до -100 м; 4 — от -500 до -300 м; 5 — от -600 до -500 м; S — шахтный ствол Параметры водопроводимости определяются по формулам

(4)

где K_x и K_y — значения коэффициента фильтрации вдоль осей Ox и Oy соответственно. При послойной аппроксимации параметра проницаемости интегрирование в выражениях (4) заменяется соответствующей суммой.

Особенностью разработанной методики является способ расчета притока в выработки в соответствии с принятой схематизацией подработанного массива. Водоприток в выработки в пределах блока (i, j) определяется суммированием водопритоков по всем дренируемым горизонтам отработки, число которых равно N_ij

. (5)

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

Приток воды из окружающих пород в выработки в блоке на k-м горизонте отработки можно вычислить по линейному закону

(6)

где K_w,ijk — коэффициент фильтрации на данном горизонте в блоке (i, j), S_w,ijk — площадь поверхности «выработки — окружающий массив», H_w,ijk — перепад напора в водонасыщенных породах вокруг выработок и поверхностью стока, L_a,ijk — средняя длина пути фильтрации из произвольной точки массива к поверхности выработок. Параметры, входящие в выражение (6), можно оценить на основе статистических характеристик нарушенного горного массива.

Так, площадь фильтрации S_w (далее индексы i, j, k для простоты опущены) определяется следующим образом. Исходя из приближения к цилиндрической форме выработок, можно записать для объема пустот V_w в пределах некоторого блока-параллелепипеда

(7)

где l_w(r) — длина выработки радиуса r, — ее среднее значение для выработок всех размеров, p_w(r) — плотность распределения выработок по объему блока, [r_min, r_max] - диапазон изменения радиусов выработок.

При аналогичных допущениях

(8)

где p_s(r) — плотность распределения выработок по площади их боковой поверхности в зависимости от радиуса. Сопоставляя (2.7) и (2.8), получим

. (9)

Параметры распределений p_w(r) и p_s(r) можно оценить при наличии информации о размерах выработок в данном блоке или части шахтного поля. При равномерном распределении (p_w(r) = p_s(r) =1/(r_max — r_min)₎

. (10)

Тогда при обычном диапазоне радиусов выработок (r_min=0,5 м, r_max=2,2 м) S_w1,32 V_w. С увеличением доли малых выработок числовой коэффициент перед V_w возрастает, хотя обычно не превышает 2−2,5. Увеличение площади контакта «выработка — окружающий массив» аналогично увеличению удельной поверхности скелета пористой среды при уменьшении размеров частиц.

При осреднении напора по вертикали в принятой расчетной схеме перепад напора H_w,ijk рассчитывается как разность между уровнем подземных вод и средним положением дренируемых (или затапливаемых) выработок на горизонте отработки:

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

Параметр L_a в формуле (6) можно оценить следующим способом. Предположим, что пустотность (доля объема выработок в массиве) в горизонтальном или вертикальном сечении блока равна его объемной пустотности _a= V_w/V_b, где V_b — объем блока. Тогда суммарная площадь сечения выработок через некоторую грань площадью S составляет S_w, = S. Средняя площадь сечения отдельной выработки гранью блока

(11)

где — угол между выработкой и гранью блока. При равномерном распределении p_s(r) и тех же значениях r_max и r_min величина составит около 9 м². Предположим далее, что отношение площади непроницаемой фазы к площади выработок в сечении остается одинаковым в масштабе блока и для более мелкого блока вокруг отдельной выработки:

(12)

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

В пределах сечения площадью S₁ вокруг выработки идет фильтрация к границе «выработка — массив». Среднее расстояние от произвольной точки такого сечения до границы можно определить, заменив его одной из простых фигур (круг, квадрат) той же площади (рис. 5). В случае круга выделим внешнюю границу — окружность радиусом и внутреннюю — окружность радиусом. Среднее расстояние от произвольной внутренней точки такой фигуры до малой окружности составит

. (13)

Рис. 5. Схема к расчету водопритока в выработки из массива

При типичных значениях пустотности и размеров выработок значение r_a составляет несколько метров, уменьшаясь в сильно нарушенных участках и возрастая в частях шахтного поля с малым объемом выработанного пространства. С учетом извилистости пустот следует принять L_a= r_a, где — коэффициент извилистости, принимаемый по аналогии с пористой средой.

Изменение уровня воды в стволе в период осушения и добычи определяется по графику отработки угольных пластов. После отключения водоотлива на каждом шаге по времени величина H_s рассчитывается по формуле [4]:

(14)

где Q_w,(t) — суммарный приток в выработки в момент времени t, S_h(H_s) — площадь горизонтальной поверхности выработок на уровне H_s.

Разработанная методика применена для расчета изменения уровня подземных вод на шахтном поле в период его разработки, и после отключения водоотлива. Исходные данные принимались по горно-гидрогеологическим условиям, типичным для Центрального района Донбасса. На границе шахтного поля размером 56 км установился стационарный уровень подземных вод 200 м. Трещинная пористость пород массива составляет 0,03. При затоплении, для учета влагосодержания в трещинах, вместо этой величины в уравнении (3) задавался меньший параметр водоотдачи. Параметры, характеризующие распределение коэффициента фильтрации и объема выработанного пространства по глубине, приведены на рис. 3 и 4. Расчеты проводились на прямоугольной сетке с размером блоков 5050 м, с максимальным шагом по времени 10 сут. Результаты моделирования представлены на рис. 6 — 8.

Анализ изменения суммарного водопритока в шахту (рис. 6) показывает, что наибольший приток наблюдается при отработке верхних слоев, характеризующихся большей проницаемостью. Этим объясняются, в частности, локальные максимумы величины Q_w в первые 10−20 лет после начала разработки. В этот период уровень подземных вод находится преимущественно в пределах двух более проницаемых верхних горизонтов отработки. По мере заглубления выработок и дренирования нижних горизонтов водоприток стабилизируется.

Рис. 6. Изменение водопритока (1), уровня воды в шахтном стволе (2) и уровня подземных вод вблизи шахтного ствола (3)

Рис. 7. Уровень подземных вод на шахтном поле перед отключением водоотлива а)

б) Рис. 8. Изменение уровня подземных вод в шахтном поле после отключения водоотлива: а) вдоль профиля А-А, б) вдоль профиля В-В (см. рис. 7); 1, 2, 3 — моменты времени.

Период незначительного уменьшения водопритока (28−55 лет после начала отработки) характеризуется практически постоянным уровнем воды в шахтном стволе и сформировавшейся воронкой депрессии подземных вод (рис. 7), которая повторяет контуры зон отработки. Наибольшее понижение приурочено к тем зонам, где, в соответствии с принятыми исходными данными, объем выработанного пространства больше. По мере подъема уровня подземных вод при затоплении отмеченные неравномерности сглаживаются (рис. 8).

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

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

Использование и охрана подземных вод / Н. А. Маринов, А. Е. Орадовская, Е. В. Пиннекер и др. Новосибирск, Наука, 1983. 231 с.

Гонтаревский В.П., Кулешов В. М. Обеспечение эколого-гидробезопасности при ликвидации шахт // Уголь Украины. — 1999. — № 11. — С. 22−25.

ДБН Ф. 2.2−1-95. Состав и содержание материалов оценки воздействий на окружающую среду при проектировании и строительстве предприятий, зданий и сооружений. Основные положения проектирования / Госкомградостроительства Украины, Минэкобезопасности Украины. — К.: Укрархстройинформ, — 1996. — 14 с.

КД 12.12.004−98 «Ликвидация угольных шахт. Защита земной поверхности от затопления горных выработок. Рекомендации» / УкрНИМИ, 1998. — 46 с.

Инструкции по безопасному ведению горных работ у затопленных выработок / ВНИМИ, 1984. — 66 с.

Бошенятов Е. В. Влияние метаморфизма горных пород на величину безопасной глубины разработки под водными объектами // Горная геомеханика и маркшейдерское дело: Сб. науч тр. ВНИМИ. — СПб, 1999. — С. 150−153.

Панасенко Г. П. О методике определения мощности зоны водопроводящих трещиноватых пород над очистными выработками угольных шахт Восточного Донбасса // В сб. «Результаты геологических исследований на территории Нижнего Дона и Нижней Волги». — Ростов-на-Дону: Изд-во РГУ, 1967. — С. 54−61.

Мохов А. В. Определение высоты зоны водопроводящих трещин на каменноугольных месторождениях в натурных условиях. — М., 1988. — 25 с.

Беседа Н. И. Особенности гидрогеологиченских условий вскрытия и разработки углей на больших глубинах в юго-западной части Донбасса / УкрНИИНТИ. — К., 1970. -49 с.

Методические указания по изучению и прогнозированию гидрогеологических условий при разработке угольных месторождений Донецкого бассейна. Днепропетровск: ДО ИМР, 1979. — 100 с.

Беседа Н.И., Сляднев В. А., Яковлев Е. А. Региональные техногенные изменения геологической среды Донбасса под влиянием горных работ. — К. — 1977. — 66 с.

Цабут И.И., Улицкий О. А., Шворак Н. В. Анализ многолетней динамики обводнения угольных шахт Центрального углепромышленного района Донбасса // Труды ДонГТУ, Серия горно-геологическая, вып. 11, Донецк, ДонГТУ, 2000. — С. 124−128.

Ермаков В.Н., Улицкий О. А., СпожакинА.И. Изменение гидродинамического режима шахт при затоплении // Уголь Украины. — 1998. — № 6. — С. 11−13.

Ермаков В.Н., Улицкий О. А., Питаленко Е. И. и др. Изменение гидрогеологических и геомеханических условий при закрытии шахт // Наук. Праці ДонНТУ, Серія гірничо-геологічна, вип. 321, Донецьк, ДонНТУ, 2001. — С. 69−73.

Ермаков В.Н., Семенов А. Н., Улицкий О. А. и др. Развитие процессов подтопления земной поверхности под влиянием закрывающихся шахт // Уголь Украины. — 2001. — № 6. — С. 12−14

ДНАОП 1.1.30 — 1.01.96. Правила безпеки у вугільних шахтах. — Донецьк: Донеччина, 2001. — 495 с.

Ермаков В.Н., Петренко С. Я., Касимов О. И., Кочерга В. Н. О предотвращении выделения газов из ликвидируемых шахт Стахановского региона // Уголь Украины, — № 5. — 1999. — С. 15−17.

Василянский Н.П., Кочерга В. Н. Предотвращение газовыделения на земную поверхность путем дегазации выработанных пространств / Сб. науч. трудов МакНИИ. — 1987. — С. 27−32.

Айкер Х., Хесбрюгге С. Газовыделение из закрытых шахт и погашенных выработок // Глюкауф. — № 23. — 1984. — С. 29−34.

Robinson R. Mine gas hazards in the surface environment // Mining Technology, Section A. — Vol. 109. — 2000. — P. A228−236.