УДК 622.253
ПРИМЕНЕНИЕ ПРОГРАММНОГО КОМПЛЕКСА «РЯОгЕ^АЬЬ» ДЛЯ РАСЧЕТА ИСКУССТВЕННОГО ЗАМОРАЖИВАНИЯ ПОРОД
Л.Ю. Левин,А.В. Богомягков, О.С. Паршаков, М.А. Семин
Представлены общие сведения о разработанной программе «Гго2епЖа11», предназначенной для моделирования ледопородного ограждения вертикальных шахтных стволов. Особенностью разработанной программы является возможность калибровки математической модели на предмет наилучшего соответствия с данными экспериментальных измерений температуры путем численного решения обратной задачи Стефана. Кроме того, предусмотрен учет ряда дополнительных физических процессов, влияющих на состояние ледопородного ограждения. На основании выполняемых в программе расчетов можно разрабатывать технические мероприятия, направленные на обеспечение эффективности строительства шахтных стволов в сложных гидрогеологических условиях, с учетом влияния технологических факторов и возможных аварийных ситуаций.
Ключевые слова: численное моделирование, мониторинг, ледопородное ограждение, породный массив, искусственное замораживание, обратная задача Стефана.
В настоящее время, искусственное замораживание пород является необходимой составляющей процесса строительства практически любого вертикального шахтного ствола в сложных гидрогеологических условиях [1]. Суть методаискусственного замораживания заключается в создании защитного ограждения из мерзлой водонасыщенной породы вокруг строящейся горной выработки для восприятия горного и гидростатического давления, предотвращения попадания грунтовых вод в выработку. Для этого по контуру, окружающему горную выработку, бурится серия замораживающих скважин, в которые монтируются колонки, по которым в дальнейшем циркулирует хладоноситель с температурой ниже нуля.
Защитное ограждение из мерзлой водонасыщенной породы, формирующееся в окрестности замораживающих скважин, в литературе называется ледопородным ограждением (ЛПО), его схематичное изображение приведено на рис. 1. Эффективность функционирования ЛПО в конечном счете зависит от его толщины и сплошности по периметру и глубине контура замораживающих скважин. При достижении достаточной толщины и сплошности ЛПО, начинается проходка горной выработки.
Для оценки текущего состояния ЛПО и во избежание возникновения аварийных ситуаций во время проходки выработки, необходимо организовать систематический термометрический контроль состояния замораживаемого породного массива [2, 3, 4]. Систематический контроль подразумевает проведение измерений температуры породного массива в контрольных скважинах и параметров работы замораживающих станций (расхода, температуры рассола) с заданной периодичностью, последующая оперативная обработ-
ка и интерпретация измеренных данных, формирование отчетов о состоянии ЛПО. Под интерпретацией в данном случае понимается восстановление поля температур и границ ледопородного ограждения по всему объему породного массива путем проведения численного моделирования тепло - и массоперено-са в породном массиве с учетом влияния технологических факторов.
Рис. 1. Принципиальная схема ЛПО: 1 - зона охлаждения, 2 - зона льда, 3 - шахтный ствол, 4 - замораживающие колонки, 5 - фронт фазового перехода
Математическое моделирование процесса замораживания породного массива вокруг строящихся шахтных стволов является одной из основных частей контроля формирования и состояния ледопородного ограждения [512]. Моделирование позволяет рассчитывать распределение термодинамических и гидравлических параметров породного массива по всему объему, подверженному термическому воздействию в результате работы системы замораживания. На основании результатов моделирования устанавливаются фактические и прогнозные параметры ЛПО.
На данный моментв мире существует ряд программ для численного теплового расчета грунтов, такие как Frost 3D, ANSYS, Abaqus, Comsol. Они способны моделировать распределения температур в расчетной области с учетом фазовых превращений, конвективного и кондуктивного переноса теплоты и содержат инструменты для визуализации тепловых полей в сечении цветовым распределением и изолиниями. Эти программы имеют широкий диапазон областей применения и рассчитаны на построение сложных 3 D моделей. Вследствие этого они имеют высокий порог вхождения для специалиста.
Помимо этого, при моделировании во всех существующих программных комплексах имеется следующая проблема. Как показывает практика [13], чаще всего, модельные значения температуры, вычисленные посредством решения прямой задачи тепло- и массопереноса, со временем начинают значительно расходиться со значениями фактических температур пород в скважинах. Несоответствие теоретических и экспериментальных температур объясняется недостаточностью выборки образцов керна и объемов работ по проведению разведочного бурения перед началом строительства шахтного ствола, а также неоднородностью и анизотропностью реального обводненного породного массива.Для уменьшения расхождения модельных и экспериментальных значений температур, требуется калибровать теплофизические параметры горных пород для наилучшего соответствия.
Одним из способов калибровки параметров модели является численное решение обратной задачи Стефана [13]. Численное решение этой задачи может быть найдено с помощью метода естественной регуляризации, заключающегося в минимизации функционала рассогласования теоретически рассчитанных и экспериментально измеренных температур породного массива в контрольных скважинах.
Цель настоящего исследования заключалась в разработке специализированного и простого в освоении программного обеспечения, с помощью которого можно осуществлять более точный и быстрый прогноз параметров ЛПО строящихся шахтных стволов с возможностью калибровки теплофизи-ческих параметров горных пород путем численного решения обратной задачи Стефана.
В результате исследования была разработана специализированная программа «Frozenwall», описанию которой посвящена настоящая статья.
Структура программы «FrozenWall». Программа«FrozenWall» осуществляет расчет теплопереноса и фазовых превращений влаги в замораживаемом породном массиве и предназначена для проведения контроля за массивом, построения отчетов, информирующих о состоянии ЛПО.Она имеет ряд инструментов для решения различных прикладныхзадач, связанных с прогнозированием термодинамических процессов, происходящих в породном массиве в условиях его замораживания. Укрупненный алгоритм работы программы представлен на рис. 2.
модельных и фактических значений температуры замораживаемых горных пород и в случае, когда погрешность превышает заданную точность
^тах, проводить численную минимизацию функционала невязки I, полученного в результате регуляризации обратной задачи Стефана.
На рис. 3 представлено главное окно программы «FrozenWall» с тематическими вкладками, содержащими различные команды.
Вкладка «Исходные данные» включает в себя команды, предназначенные для формирования геологической слоистой модели участка породно-
го массива и расчетных исходных данных в готовые таблицы, «Расчет» — для выполнения расчетов теплораспределения в замораживаемом породном массивеи проведения калибровки расчетных параметров на основании термометрии контрольных скважин, «Анализ результатов» — для визуализации термодинамических процессов и составления отчета о состоянии ЛПО.
/Ввод начальных исходных данных (геология и гидрогеология, теплофизические и физико-механические свойства горных пород)
да
> t
Оформление результатов
(отчет о состоянии ледопородного ограждения)
Рис. 2. Блок-схема, иллюстрирующая алгоритм работы программы
«Frozen Wall»
Разработанный итерационный алгоритм позволяет определять невязку Исходные данные. В качестве исходных данных для расчета температурного поля замораживаемого породного массива принимаются физические свойства горных пород, график замораживания, включающий данные о температуре хладоносителя, и его объемном расходе, отклонения устьев и забоев скважин от проектного направления, измеренные посредством инкли-нометрии, параметры воздуха, поступающего на проветривания горных работ, фактическая отметка забоя и другие. Данные заносятся в программу путем заполнения таблиц во вкладке «Исходные данные». Помимо это-
го,дополнительные технологические параметры задаются в специальном окне, открывающемся при помощи команды «Файл — Расчетные параметры». Они включают в себя параметры бетона, хладоносителя, замораживающих колонок, ствола, и прочие параметры. Они сгруппированы по соответствующим вкладкам окна.
Рис. 3. Главное окно программы «FrozenWall»
Расчет. Проведение расчетов теплораспределения в замораживаемом породном массиве выполняется прежде всего посредством численного решения прямой задачи теплопроводности с движущейся границей фазового перехода (прямой задачи Стефана). Прямая задача Стефана в декартовых координатах задается следующей системой уравнений, записанной в энтальпий-ной постановке:
дИ (Т) _ д
д Л дТ Л д + —
дх
X—
у дх у
ду
дТ X— ду
Х = Х1д (1 -ф/се ) + Xsd Ф
И (Т ) =
Р!дс!д (Т - Тр1 ) + Р1д™Ь>
Р1д^ (1 -ф/се ),
psdcsd (Т - Тр 2) >
Хп -УТ -а(Т-ь - Т)
Тр1 < Т Тр 2 < Т < Тр1
Т < Т
р2
о
-X
т\
То (к);
ские координаты, м; г — физическое время, сек; Хд, X ^ — теплопроводность
среды в зонах охлаждения и льда соответственно, Вт/(м-К); сд, с^ — удельная теплоемкость среды в зонах охлаждения и льда соответственно, Дж/(кг-К); рд, р ^ — плотность среды в зонах охлаждения и льда соответ-
(полной кристаллизации) соответственно, °С; фгсе — концентрация твердой фазы подземных вод; Ь — удельная теплота плавления льда, Дж/кг; w — влаго-содержание массива, кг/кг; Тр — температура стенок замораживающей скважины, °С; То — температура непотревоженного породного массива на удалении, °С; а — коэффициент теплообмена на границе с замораживающими скважинами, Вт/(м -°С); к — средняя глубина слоя горной породы, м; = - границы со всеми замораживающими скважинами /=/,...,Л';
О.оиг — внешняя граница области моделирования, п — координата вдоль нормали к границе Ор, м.
Численное решение прямой задачи осуществляется с помощью метода конечных разностей. Вблизи границы с замораживающими колонками задается сгущение конечно-разностной сетки для корректного моделирования высоких градиентов температуры в начальный период замораживания, а также при изменении режимов работы замораживающих колонок (рис. 4).
Для уточнения теплофизических параметров математической модели массо- и теплообменных процессов, происходящих в обводненном массиве горных пород, проводится их калибровка на основании фактических температур, измеренных в контрольных скважинах.
ная температура в 1-й скважине.
Посредством условия (1) прямая задача Стефана переопределяется и становится обратной, после чего она регуляризуется и численно решается методом градиентного спуска [11, 14]. То есть происходить минимизация функционала:
ственно, кг/м3; Тр1,Тр2 — температура начала и конца фазового перехода
— модель-
(1)
Кроме того, в программе «FrozenWall» предусмотрен учет ряда дополнительных физических процессов, влияющих на состояние параметров ледопородного ограждения: теплообмен горных пород с воздухом в шахтном стволе, влияние экзотермической реакции твердения бетона при возведении крепи, термическое сопротивление замораживающих скважин, выход из строя замораживающих колонок.
Учет теплообмена горных пород с воздухом в шахтном стволе. В программном обеспечении реализована возможность задания переменных во времени температуры и скорости воздуха, поступающего для проветривания шахтного ствола в процессе его проходки. Зависимости температуры и скорости от времени задаются как дискретный набор точек. Значения для промежуточных значений времени получаются посредством линейной интерполяции.
При расчете теплообмена между породным массивом, крепью ствола и поступающим на проветривание атмосферным воздухом, тепловой поток рассчитывается согласно закону Ньютона:
д = аУа (г)[ Та (г)- Тт ],
где q — тепловой поток на единице поверхности, Вт/м ; а — коэффициент теплоотдачи, зависящий от скорости Va движения воздуха в стволе,
Вт/(м2-°С); Та — температура воздуха в стволе, °С; Тт — температура породного массива (крепи) на границе с замораживающей колонкой, °С.
Влияние экзотермической реакция твердения бетона при возведении крепи. В момент появления шахтного ствола в исследуемом слое горных пород и, соответственно, передовой бетонной крепи на последней задается равномерно распределенный мгновенный источник выделения тепла.
Суммарная выделившаяся теплота рассчитывается исходя из толщины ^ бетонной крепи для рассматриваемого слоя горных пород и удельной теплоты фазового перехода Ьс твердения на 1 кг бетона. Считается, что после затвердевания бетона, в нем не остается влаги.
Уравнение теплопереноса для бетонной крепи имеет следующий вид:
ди
— = ХАТт + - )РсЬс >
где р — плотность бетона, кг/м ; су — теплоемкость бетона, Дж/(кг-°С); И — удельная энтальпия, Дж/м ; г —время, с; X — теплопроводность бетона, Вт/(м-°С); Тт — температура породного массива (бетонной крепи), °С; —
момент проходки ствола, с; рс — плотность бетона, кг/м3.
Учет термического сопротивления замораживающих скважин. Если в межтрубном пространстве между замораживающей колонкой и кондуктором, имеются пустоты, которые обусловлены некачественным производством тампонажных работ, то термическое сопротивление некоторых замораживающих скважин значительно возрастает, вследствие чего замораживание слоев горных пород происходит крайне неравномерно, в этом случае в программе предусмотрено задание термического сопротивления для каждой замораживающей колонки.
Термическое сопротивление Я при идеальных прямых, центрированных кондукторе и колонке рассчитывается по формуле
Я = —1—1п °к - 2 ак , 2 пХ1 Озк
где X — коэффициент теплопроводности воздуха при отрицательной температуре, Вт/(м-°С); Вк и Взк — диаметры труб кондуктора и замораживающей
колонки по внешней стенке, м; ак — толщина стенки трубы кондуктора, м; I — длина участка трубы, м.
Выход из строя замораживающих колонок. В программе также реализована возможность моделировать ситуацию, когда нарушена работа выборочного набора замораживающих колонок или их аварийный выход из строя. Для этого для каждой из замораживающих колонок введена функция состояния S (t), значение которой варьируется от 0 до 1, где 1 соответствует случаю, когда колонка работает в штатном режиме, а значение 0 - случаю ее аварийного отключения.
Граничное условие теплообмена на стенке замораживающей колонки имеет вид
Чг = h Si (t) Va(t)\Ta(t)-Tm], где Si (t) — функция состояния замораживающей колонки; q — тепловой поток на единице поверхности, Вт/м ; h — коэффициент теплоотдачи, зависящий от скорости Va движения воздуха в стволе, Вт/(м -°С); Ta — температура
воздуха в стволе, °С; Tm — температура породного массива (крепи) на границе с замораживающей колонкой, °С.
Вывод результатов. По результатам моделирования программа формирует отчет о состоянии ЛПО в формате PortableDocumentFormat (PDF). Отчет содержит следующую информацию:
графики температуры прямого и обратного потока хладоносителя, и его объемного расхода за весь период, включая текущую дату;
графики среднесуточных температур в каждой контрольно-термической скважине с отметками глубин, наложенные на укрупненную литологию модели;
карты распределения температур в горизонтальной плоскости для каждого слоя породного массива, учитываемого в модели;
карта распределения температур в вертикальном сечении через область ледопородного ограждения с нанесением контура ствола в проходке;
карта распределения температур в главной вертикальной плоскости, проходящей через контур замораживающих колонок;
таблица фактических толщин ледопородного ограждения в слоях на текущие и предыдущие сутки.
Примеры некоторых из графиков представлены ниже.
Рис. 5. Графики температуры прямого и обратного потока
хладоносителя
15 -10 -5 0 5 10 15
Рис. 6. Карта распределения температур в горизонтальной плоскости
На рис. 6 представлено температурное поле на 343-и сутки замораживания горного массива в исследуемом слое песчаника, через который осуществлялась проходка Скипового ствола № 1 рудника Талицкого ГОК.
По итогам проведенной работы сотрудниками Горного института УрО РАН представлена простая в освоении программа, которая позволяет оперативно проводить расчеты состояния ЛПО и получать отчеты, включающие в себя данные о фактическом состоянии ЛПО, графики среднесуточных температур в контрольно-термических скважинах, карты распределения температур в разных плоскостях.
Также для поддержания необходимой точности проводимых расчетов, предусмотрена возможность калибровки параметров модели путем решения обратной задачи Стефана. Программа успешно прошла испытания при ведении термометрического мониторинга состояния ЛПО на рудниках Талицкого ГОК и Нежинского ГОК.
Исследование выполнено при финансовой поддержке УрО РАН в рамках научного проекта № 18-5-5-5 и при финансовой поддержке Совета по грантам Президента РФ в рамках научного проекта № МК-6244.2018.5.
Список литературы
1. Роменский А.А. Обоснование параметров проходнического цикла и ледопородного ограждения при строительстве вертикальных стволов: дис. ... канд. техн. наук. М., 1983. 227 с.
2. ВСН 189-78 Ведомственные строительные нормы. Инструкция по проектированию и производству работ по искусственному замораживанию грунтов при строительстве метрополитенов и тоннелей, утв. Приказом Главного технического управления Министерства транспортного строительства СССР 10.05.1978. № 21. 79 с.
3. ПБ 03-428-02 Правила безопасности при строительстве подземных сооружений, утв. постановлением Госгортехнадзора России 01.11.2001 г. № 49. 167 с.
4. СП 103.13330.2012 Защита горных выработок от подземных и поверхностных вод, утв. приказом Министерства регионального развития Российской Федерации (Минрегион России) 30.06.2012 г. 68 с.
5. Амосов П.В., Лукичев С.В., Наговицын О.В. Влияние пористости породного массива и температуры хладоносителя на скорость создания сплошного ледопородного ограждения // Вестник Кольского научного центра РАН. 2016. № 4 (27). С. 43-50.
6. Математическое моделирование искусственного замораживания грунтов / П.Н. Вабищевич, М.В. Васильева, В.Ф. Горнов, Н.В. Павлова // Вычислительные технологии. 2014. Т. 19, № 4. С. 19-31.
7. Левин Л.Ю., Семин М.А., Паршаков О.С. Математический метод прогнозирования толщины ледопородного ограждения при проходке стволов // ФТПРПИ. 2017. № 5. С. 154-161.
8. Николаев П.В., Шуплик М.Н. Математическое моделирование процесса искусственного замораживания грунтов с применением твердого диок-
сида углерода // Горный информационно-аналитический бюллетень. 2015. № 11. С. 243-251.
9. Хайден Т.В., Вегнер Б. Современная технология замораживания пород на примере двух объектов, находящихся на стадии строительства, состоящих из пяти вертикальных шахтных стволов // Горный журнал. 2014. № 9. С. 65-69.
10. Левин Л.Ю., Семин М.А., Зайцев А.В. Контроль и прогноз формирования ледопородного ограждения с использованием оптоволоконных технологий // Сб. науч. трудов. Стратегия и процессы освоения георесурсов. Пермь, 2016. С. 236-238.
11. Schmall P., Brinckerhoff P., Ott K., Madsen P. ESA's SEM challenge. // North American Tunneling Journal. 2013. P. 24-26.
12. Sopko J. Coupled Heat Transfer and Groundwater Flow Models for Ground Freezing // ASCE GeoFlorida. 2017. P. 719-729.
13. Левин Л.Ю., Семин М.А., Зайцев А.В. Калибровка теплофизиче-ских свойств породного массива при моделировании формирования ледопо-родного ограждения строящихся шахтных стволов // Физико-технические проблемы разработки месторождений полезных ископаемых. 2019. № 1. С. 172-184.
14. Jarny Y., Ozisik M.N., Bardon J.P. A general optimization method using adjoint equation for solving multidimensional inverse heat conduction // International Journal of Heat and Mass Transfer. 1991. Vol. 34, № 11. P. 2911-2919.
Левин Лев Юрьевич, д-р техн. наук, зав. отд. аэрологии и теплофизики, aerolog_lev@mail. ru, Россия, Пермь, Горный институт УрО РАН,
Сёмин Михаил Александрович, канд. техн. наук, науч. сотр., seminma@outlook.com, Россия, Пермь, Горный институт УрО РАН,
Богомягков Александр Васильевич, млад. науч. сотр., bavaerolog@,gmail.com, Россия, Пермь, Горный институт УрО РАН,
Паршаков Олег Сергеевич, млад. науч. сотр., olegparshakov@,gmail.com, Россия, Пермь, Горный институт УрО РАН
THE APPLICA TION OF "FROZEN WALL " SOFTWARE IN SIMULA TION OF ARTIFICIAL
GROUND FREEZING
L.Y.Levin, A.V. Bogomyagkov, O.S. Parshakov, M.A. Semin
The paper presents general information about the software application "Frozen Wall", which was designed to simulate frozen wall formation around constructed vertical shafts. The main feature of the developed application is the possibility of calibrating the mathematical model for the best fit with the experimental temperature measurements by numerical solution of the inverse Stefan problem. In addition, it takes into account a number of technological processes that affect the state of the frozen wall. Based on calculations performed in the application, it is
possible to develop technical measures aimed at ensuring the efficiency of mine shafts construction in difficult hydro-geological conditions.
Key words: numerical simulation, monitoring, frozen wall, rock mass, artificial ground freezing, inverse Stefan problem.
Levin Lev Yurievich, doctor of technical sciences, head of the mine ventilation and thermal physics department, aerolog_lev@mail. ru, Russia, Perm, Mining Institute of the Ural Branch of the Russian Academy of Sciences,
Semin Mikhail Aleksandrovich, candidate of technical sciences, researcher, seminma@,outlook. com, Russia, Perm, Mining Institute of the Ural Branch of the Russian Academy of Sciences,
Bogomyagkov Aleksandr Vasilievich, junior researcher, bavaerolog@gmail. com, Russia, Perm, Mining Institute of the Ural Branch of the Russian Academy of Sciences,
Parshakov Oleg Sergeevich, junior researcher, olegparshakov@gmail. com, Russia, Perm, Mining Institute of the Ural Branch of the Russian Academy of Sciences
Reference
1. Romensky A. A. Substantiation of the parameters of the tunneling cycle and ice-bearing fencing in the construction of vertical shafts: Diss. ... Cand. tech. sciences'. Moscow, 1983. 227 PP.
2. VSN 189-78 Departmental building codes. Instructions for the design and production of works on artificial freezing of soils in the construction of subways and tunnels, UTV. Order of the Main technical Department of the Ministry of transport construction of the USSR 10.05.1978. No. 21. 79 PP.
3. PB 03-428-02 safety Rules for the construction of underground structures, UTV. by the resolution of Gosgortekhnadzor of Russia 01.11.2001 No. 49. 167 PP.
4. SP 103.13330.2012 Protection of mine workings from underground and surface waters, UTV. order of the Ministry of regional development of the Russian Federation (Minregion of Russia) 30.06.2012 68 p.
5. Amosov P. V., Lukichev S. V., Nagovitsyn O. V. Influence of porosity of rock mass and coolant temperature on the speed of creation of continuous ice-bearing fence.Vestnik Kola scientific center RAS. 2016. No. 4 (27). Pp. 43-50.
6. Mathematical modeling of artificial soil freezing / P. N. vabischevich, M. V. Vasi-lyeva, V. F. Gornov, N. V. Pavlova // Computational technologies. 2014. Vol. 19, No. 4. Pp. 1931.
7. Levin L. Yu., Semin M. A., Parshakov O. S. Mathematical method of predicting the thickness of the ice-bearing fence when sinking trunks / / fftprpi. 2017. No. 5. Pp. 154-161.
8. Nikolaev P. V., Shuplik M. N. Mathematical modeling of the process of artificial freezing of soils with the use of solid carbon dioxide / / Mining information and analytical Bulletin. 2015. No. 11. Pp. 243-251.
9. Hayden T. V., Wegner B. Modern technology of rock freezing on the example of two objects under construction, consisting of five vertical mine shafts. Gorny Zhurnal. 2014. No. 9. Pp. 65-69.
10. Levin L. Yu., Semin M. A., Zaitsev A.V. Control and forecast of the formation of ice-bearing fencing with the use of fiber-optic technologies. labours'. Strategy and processes of development of georesources. Perm, 2016. Pp. 236-238.
11. Schmall P., Brinckerhoff P., Ott K., Madsen P. ESA's SEM challenge. // North American Tunneling Journal. 2013. P. 24-26.
12. Sopko J. Coupled Heat Transfer and Groundwater Flow Models for Ground Freezing // ASCE GeoFlorida. 2017. P. 719-729.
13. Levin L. Yu., Semin M. A., Zaitsev A.V. Calibration of thermal and physical properties of rock mass in modeling the formation of ice-bearing fencing of mine shafts under construction. Physico-technical problems of mineral deposits development. 2019. No. 1. Pp. 172184.
14. Jarny Y., Ozisik M. N., Bardon J. P. A general optimization method using adjoint equation for solving multidimensional inverse heat conduction // International Journal of Heat and Mass Transfer. 1991. Vol. 34, No. 11. P. 2911-2919.