Введение

Глобальное изменение климата нельзя остановить или предотвратить. Это естественный циклический процесс развития планеты Земля, который не подвластен человеку, к нему надо приспосабливаться.

Горные породы, длительное время (от нескольких лет до многих тысячелетий) находящиеся при температурах ниже 0°С и сцементированные льдом, получили название многолетней или вечной мерзлоты. Зона многолетней мерзлоты занимает более 65 % территории России или около 10,7 млн.кв.км. и охватывает значительную часть Сибири и Дальнего Востока. В пределах криолитозоны России сосредоточено более 30 % разведанных запасов всей нефти страны, около 60 % природного газа. Изменение состояния пород коры безусловно повлияет на инженерно-технические сооружения, окружающую среду и хозяйственную деятельность человека.

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

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

Система уравнений

Систему (ВГ) можно записать в виде сопряжённой трехмерной математической модели тепломассопереноса в вязких и пористых средах, которая в терминах функции тока, давления, температуры и концентрации выглядит так (рис. 1).

Рис. 1. Система уравнений тепломассопереноса, основанной на теории вихревой гидродинамики ВГ

Это система уравнений шестого и более порядка, в зависимости от количества компонентов (m). Первые четыре уравнения описывают гидродинамику, причём первые три, записанные через функцию тока, векторный потенциал, свободную конвекцию, вызываемую внутренними силами (градиентом плотности), а четвёртое, записанное через скалярный потенциал, вынужденную конвекцию, вызываемую любыми внешними силами. Пятое уравнение – теплопроводности, шестое – массопереноса. Прямая связь между гидродинамическим блоком и уравнениями тепломассопереноса осуществляется через скорости конвекции (фильтрации) уравнение (7) – обратная через плотность среды (флюида) уравнение состояния (8). Для решения данной системы задаются начальные и граничные условия (9). Любые начальные состояния, значения с предыдущего момента времени, внешние воздействия (транзитные потоки и т.д.) или взаимодействия разных веществ между собой и их фазовые переходы задаются с помощью функций источника, записанных в правой части для каждого уравнения.

Для решения системы уравнений, описывающей процессы конвективного тепломассопереноса в трехмерных областях, необходимо использовать численные методы, позволяющие в каждом узле дискретной сетки, покрывающей модельную область, рассчитать значения искомых функций. Для численной реализации краевых задач используется абсолютно устойчивый локально-одномерный метод решения [4]. Суть его заключается в расщеплении уравнения по координатам и последующем решении одномерных уравнений методом прогонки, в сумме являющихся решением многомерного уравнения. Эта неявная схема решения имеет суммарную аппроксимацию порядка O(τ+ ћ²), где τ — шаг по времени, ћ — шаг сетки. Для аппроксимации конвективных членов в уравнениях (10,11) используется схема с направленными разностями [5]. Кроме того, эти уравнения расщепляются по процессам, т. е. в каждый момент времени вначале учитывается вклад только кондуктивного и дисперсионного членов, рассчитанных по неявной схеме, а затем к ним добавлялись значения конвективных составляющих, рассчитанных по явной схеме.

Подробное описание модели и методики расчета сделано в работе [2]. На основе данной модели был разработан пакет программ на языке Fortran для персональных компьютеров, позволяющий решать широкий круг задач тепломассопереноса в пористых средах.

Постановка задачи и схематизация

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

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

Процесс промерзания – оттаивания моделировался с заданным шагом понижения температуры, а затем оттаивания с тем же шагом её повышения. При промерзании в каждый момент времени на верхней границе температура менялась с шагом понижения до момента времени установления границы с нулевой температурой, а при оттаивании менялась с таким же шагом повышения до момента времени полного оттаивания разреза

Решались задачи для разных разрезов. 1) Сплошные слабопроницаемые породы. 2) Имеющие проницаемые пласты с разгрузкой в речную долину, расположенную в правой части расчётного блока в виде разлома с повышенной проницаемостью. 3) По середине блока имелось возвышение высотой до 200 метров.

Рис. 2. Сплошные слабопроницаемые породы. Прямоугольник – изолированная площадь

Рис. 3. Разрез пород, имеющий проницаемые пласты с разгрузкой в речную долину. Левая температура, правая коэффициент сопротивления пластов. Стрелками показаны скорости фильтрации подземных вод

Рис. 4. Разрез пород с возвышение высотой до 200 м. Левая температура, правая концентрация вещества среды (плотность). Голубым отрисован воздух

Для расчета был выбран блок в виде прямоугольного параллелепипеда, размером в горизонтальной плоскости от 3000 до 5000 м, вертикальный размер принимался от 300 до 2500 м. Блок был разбит объемной, неравномерной сеткой с шагом по горизонтали 100-500 м. Вертикальный шаг выбирался в зависимости от шага по времени. Для шага по времени 1 сутки равнялся 1-10 метров, 1 месяц —  5-10 метров. 1 год 20-50 метров. Коэффициент температуропроводности пород принимался одинаковый для всего блока и равнялся 1×10^-6 м²/сек.

В начальных условиях задавался температурный разрез с градиентом температуры 2 — 4°С/100м, на нижней границе задавалась постоянная температура, соответствующая её глубине определённой по заданному вертикальному шагу сетки. На вертикальных границах задавалась изоляция, и они удалялись на расстояние, на котором не влияли на работу системы. Для блока имеющего речную долину в разрезе задавались слабопроницаемые (1×10^-4, 1×10^-5дарси) пласты и высокопроницаемые 1-2 дарси (1дарси = 10^-12м² – К_ф = 0.86 м/сут) через которые шла разгрузка подземных вод в речную долину, при замерзании пород проницаемость их уменьшалась до уровня слабопроницаемых.

Результаты вычислений

Промерзание пород

Результаты расчёта показывают, что при расчёте варианта с постоянной на верхней границе средней за год отрицательной температурой (-1) – (-3)°С, а на нижней +10°С, для блока со сплошной слабопроницаемой породой мощностью до 400 м, глубина промерзание составляет 50 – 120 метров. Промерзание идёт неравномерно. Впервые годы скорость промерзания составляет 5-6 м/год, затем замедляется и к 400 — 1200 годам почти полностью затухает. При температурах до -7°С температура затухает к 1800 году, глубина промерзания увеличивается до 190 метров. При снижении температуры на верхней границе глубина промерзания увеличивается и при температуре -15°С составляет 230 метров, а при температуре -25°С 265 метров. К 2000 — 2400 годам в этих случаях устанавливается равновесная температура пород с нулевой границей (рис 5). Нулевая граница между мёрзлыми и талыми породами не постоянна. При повышении температуры на поверхности она смещается к верху, а при понижении к низу.

Рис. 5. Распределение температуры при промерзании с 1 года до 3000 года, при температуре на верхней границе -8°С на нижней +10°С. Точки расположены на вертикальной координате Z в центральной части блока. Шаг по глубине 6 метров по времени 1 год

При расчётах с шагом по времени один месяц на верхней границе задавалась средняя температура для каждого месяца, которая менялась в течении года от +20°С в июле до -40°С в январе. Такие колебания значений температуры влияли на процесс промерзания пород, который зависел от величины интервала времени с положительными температурами и разницы между положительными и отрицательными. За интервал месяцев с положительными температурами, верхняя часть разреза протаивала, а затем при наступлении интервала месяцев с отрицательными температурами снова начинала промерзать. Такие колебания температуры снижали скорость промерзания пород разреза. При интервале времени с положительной температурой 3-5 месяцев и средне по месяцу максимальной температурой в июле +15°С, глубина оттаивания составляла 1-1.5 метров, при повышении до +20°С, максимум до 2 метров. На глубине более 2 метров тоже были колебания значений температуры, но в диапазоне отрицательных температур, которые затухали с глубиной (рис. 6).

Рис. 6. Распределение температуры при промерзании за 20000 месяцев (~1700 лет) при колебаниях температуры на верхней границе от +17°С до -40°С на нижней постоянной исходной +10°С. Точки расположены на вертикальной координате Z в центральной части блока. Шаг по глубине 2 метра, по времени 1 месяц

При мощности разреза до 1000 – 3000 метров и температуре на нижней границе от +30°С до +100°С, на верхней от -10°С до -30°С процесс промерзания, до установления постоянной нулевой границы может продолжаться десятки тысяч лет до 90000 лет, а глубина промерзания достигать 500 — 1600 метров (рис. 7).

Рис. 7. Распределение температуры при промерзании за 90000 лет при температуре на верхней границе -30°С на нижней +60°С. Шаг по глубине 20-100 метра, по времени 1 год. Глубина промерзания 1400 метров

Проведены расчёты по оценке глубины промерзания пород для разных значений температур на верхней и нижней границах. На верхней границе блока задавались данные, взятые с карты распределение среднемесячных и среднегодовых температур в районах России. На нижней границе задавались значения температуры исходя из оценки их роста на 500 метров глубины, это примерно 10, 12 и 15°С, полученных по результатам буровых работ в северо-восточной части страны [1]. Результаты расчётов представлены в таблице.

Таблица. Глубина промерзания пород при разных граничных условиях

Промерзание – оттаивание пород

В расчётах процесса промерзания – оттаивания в каждый момент времени задавалось снижении с положительной или нулевой температуры на верхней границе на 0.02- 0.04°С в год, что составляло в среднем 3-3.5°С за 100 лет. Породы постепенно начинают промерзать и граница с нулевой температурой погружаться с начальной скоростью около 5 м/год с замедлением и полным затуханием. После установление нулевой границы на определённой глубине, на верхней границе блока задавалось повышение температуры с таким же шагом, что и понижение. Температура разреза начинала повышаться, породы оттаивать, а нулевая граница подниматься. Деградация мерзлоты шла снизу (рис. 8).

Рис. 8. Температуры за 6000 лет. При промерзании с 1 го до 2400 года и оттаивание при повышении, а затем понижении температуры на верхней границе на 0.04°С в год. Точки расположены на вертикальной координате Z в центральной части блока. Шаг по глубине 6 метров по времени 1 год

При расчётах с шагом по времени один месяц на верхней границе задавалась средняя температура для каждого месяца, которая менялась в течении года от +20°С в июле до -40°С в январе. За интервал месяцев с положительными температурами, верхняя часть разреза протаивала, а затем при наступлении интервала месяцев с отрицательными температурами снова начинала промерзать. Такие колебания температуры снижали скорость промерзания пород разреза. При промерзании пород, температура в первой верхней точке (до 2 метров) в интервал месяцев с положительными температурами может подниматься до плюсовых значений (рис. 9).

Рис. 9. Распределение температуры за 60000 месяцев (5000 лет). При промерзании за 20000 месяцев (~1700 лет) при заданных на верхней границе в каждый момент времени средних месячных температурах, а затем оттаивания при повышении значений этих температур на 0.003°С в каждом месяце. Точки расположены на вертикальной координате Z в центральной части блока. Шаг по глубине 2 метра, по времени 1 месяц

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

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

Рис. 10. Распределение температуры за 8000 лет, при заданной на верхней границе 0-ой температуре и изолированного участка, а) за 3600 лет промерзания с уменьшением её на 0.03°С в год, б) затем до 8000 лет оттаивания с увеличением её на 0.03°С в год

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

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

Рис. 11. Распределение температуры (слева) и коэффициента сопротивления пластов (справа) при заданной на верхней границе температуре -25°С. за 9000 лет и глубине промерзания 400 метров

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

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

Рис. 12. Распределение температуры, левая первый момент времени, правая 9000 момент времени, шаг по времени год, глубина промерзания 600 метров

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

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

Рис. 13. Температуры за 10000 лет. При промерзании с 1 го до 4000 года и оттаивание при повышении, а затем понижении температуры на верхней границе на 0.04°С в год. Точки расположены на вертикальной координате Z в центральной части блока. Шаг по глубине 6 метров по времени 1 год

Выводы

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

Анализ результатов численных расчётов по влияние глобального потепления на вечную мерзлоту показал.

Вечная мерзлота формируется за счёт понижения температуры на поверхности коры. Это процесс взаимодействия двух потоков, теплового потока, идущего из недр Земли и холодного потока, идущего с поверхности коры.

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

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

В районах России, где существует вечная мерзлота, средняя годовая температура равна (-8… -20) град. Вечная мерзлота сохраниться пока эта температура не повысится до 0 — +1 град, даже при повышении этой температуры на 4 град. за 100 лет, для этого потребуется 200-500 лет.

Определён критерий, по которому можно сделать оценку времени (промежутка времени) полного окончания оттаивания (деградации) мерзлоты на рассматриваемом участке.