Доклады Российской академии наук. Науки о Земле, 2023, T. 508, № 2, стр. 245-252

Разработка методики мониторинга состояния газогидратных залежей восточно-сибирского шельфа

В. А. Чеверда 1*, Д. С. Братчиков 1, К. Г. Гадыльшин 2, Е. Н. Голубева 3, В. В. Малахова 3, Г. В. Решетова 3

1 Институт математики им. С.Л. Соболева Сибирского отделения Российской академии наук
Новосибирск, Россия

2 Институт нефтегазовой геологии и геофизики им. А.А. Трофимука Сибирского отделения Российской академии наук
Новосибирск, Россия

3 Институт вычислительной математики и математической геофизики Сибирского отделения Российской академии наук
Новосибирск, Россия

* E-mail: vova_chev@mail.ru

Поступила в редакцию 21.09.2022
После доработки 06.10.2022
Принята к публикации 13.10.2022

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

Аннотация

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

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

ВВЕДЕНИЕ

Арктический шельф является одним из крупнейших источников выбросов метана в атмосферу [1]. Одним из возможных физических механизмов, ответственных за этот процесс, является разложение газогидратов, законсервированных в криолитозоне континентальных отложений на шельфе арктических морей [2]. Несмотря на отсутствие надежных геологических данных о наличии скоплений газогидратов на северном шельфе России, существуют косвенные признаки, подтверждающие возможность их образования, а именно, наличие термобарических условий, благоприятных для формирования и сохранения залежей газогидратов. Так, существование в настоящее время субаквальной мерзлоты, сформировавшейся в ледниковую эпоху вследствие понижения уровня моря, в донных отложениях Восточно-Сибирского шельфа подтверждается результатами буровых работ [3]. Присутствие в толще донных отложений многолетнемерзлых пород (ММП) обеспечивает условия для формирования в них зон стабильности гидратов метана (ЗСГМ) при глубине воды менее 100 м [4].

Реликтовые гидраты, связанные с существованием многолетнемерзлых пород, находясь на небольшой глубине, особенно чувствительны к изменениям климата. Повышение уровня моря и последующее затопление шельфа привели к повышению температуры осадочного слоя и заглублению верхней границы этих пород. Современные процессы, связанные с глобальным потеплением, способствуют проникновению дополнительного тепла в осадочный слой, что усиливает процесс деградации многолетнемерзлых пород и может приводить к дестабилизации газовых гидратов. Данные наблюдений для 1920–2009 гг. в мелководной части шельфа и прибрежной зоне моря Лаптевых и Восточно-Сибирского моря отражают значительное повышение придонной температуры (аномалии до 2.1°С), начавшееся с середины 1980-х годов [5]. Второе десятилетие 21 столетия оказалось для сибирских арктических морей значительно теплее первого, а летом 2019 и 2020 г. поверхностная температура в Карском море, в море Лаптевых и прилегающих глубоководных акваториях достигала максимальных значений, превышающих порог в 90% всех среднемесячных значений за период с 1981 по 2010 г. [6]. Накопленный морскими водами в летний период запас тепла реализуется как потоками в атмосферу, так и способствует прогреву придонного слоя шельфовых морей за счет осенне-зимнего перемешивания. Устойчивое повышение температуры придонных вод и передача тепла в осадочный слой могут приводить к дестабилизации газогидратных залежей с последующим поднятием огромного количества газа.

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

1. ЗОНЫ МНОГОЛЕТНЕМЕРЗЛЫХ ПОРОД И СТАБИЛЬНОСТЬ ГАЗОГИДРАТОВ

Численное моделирование процессов формирования многолетнемерзлых пород проводилось применительно к шельфовой области моря Лаптевых в районе Новосибирских островов с современными глубинами порядка 20 м. Выбор области исследования определялся следующими факторами:

– наличием данных по геологическому строению, а также геотермических наблюдений до глубины 200 м;

– возможностью верификации палеогеографического сценария и численной модели для данной области на основе имеющихся данных по скважинам в этом районе [7];

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

Для исследования процессов формирования многолетнемерзлых пород была использована одномерная модель теплофизических процессов в донных отложениях с учетом фазовых переходов между мерзлым и талым грунтом [7, 9]. Граничное условие на поверхности донных отложений определялось периодами трансгрессий-регрессий с учетом изменения уровня моря за последние 200 тыс. лет [7]. В восточном секторе Арктики отсутствуют палеотемпературные данные по целому ряду временных интервалов. Это обстоятельство не позволяет построить сценарий динамики температуры воздуха или пород, используя только региональные данные. В этих целях для его построения применяется методика использования изотопных палеотемпературных кривых, полученных по ледниковым кернам Восточной Антарктиды и Гренландии, которые отражают ход глобальных колебаний климата. Применение изотопных палеотемпературных кривых обеспечивает построение непрерывных во времени региональных кривых динамики температуры воздуха и пород при использовании дискретных региональных данных [7] (рис. 1).

Рис. 1.

Результаты численного моделирования на основе модели теплопереноса. На верхней панели – палеотемпературный сценарий для периода 20 тыс. лет назад – 1950 г. и 1950–2200 гг.(выделен розовым цветом). На нижней – эволюция ММП и ЗСГМ за последние 20 тыс. лет: нижняя граница ММП (PERMB) – синяя линия, верхняя граница ММП (PERMT) – голубая линия, верхняя граница ЗСГМ (HSZT) – зеленая линия. Область возможного существования ГГ выделена серым цветом.

Для задания геологического строения и физических свойств пород использовались буровые данные по строению Новосибирских островов [7]. Интенсивность геотермального потока задавалась 60 мВт/м2 в соответствии с данными для рассматриваемой области, представленными в работе [10]. Для обеспечения отсутствия влияния нижней границы расчетной области на динамику подошвы мерзлоты задаваемая по вертикали область геологической модели составила 1.5 км. Для решения задачи Стефана использовался метод ловли фронта в узел пространственной сетки. Численная реализация модели основана на методе прогонки на дискретной вычислительной сетке с вертикальным шагом 0.5 м и неявной схеме по времени с шагом 1 мес. Одновременно с вычислением термического состояния донных отложений вычислялись термодинамические границы зоны стабильности гидратов метана [9].

В работе [7] была проведена верификация данной модели с использованием реальных геологических разрезов по двум скважинам на Новосибирских островах. Геологический разрез, в верхних 80 и 200 м был задан в соответствии с результатами бурения, а ниже – в соответствии с данными геологических исследований. Коррекция свойств пород позволила получить соответствие натурных и расчетных данных.

С использованием модели теплопереноса в донных отложениях при изменениях граничных условий на верхней границе по построенному палеогеографическому сценарию [7] мы оценили изменения мощности многолетнемерзлых пород и зоны стабильности гидратов метана до 1950 г. (рис. 1). В численных расчетах для периода 1950–2200 гг. учитывался рост среднегодовой температуры придонной воды в 0.02°C за год. Изменение мощности многолетнемерзлых пород за последние 20 тыс. лет, по результатам численного эксперимента, представлено на рис. 1. В период последнего ледникового максимума (около 20 тыс. лет назад) мощность мерзлых пород достигала 440 м. Согласно полученным результатам, максимальные различия в состоянии многолетнемерзлых пород происходят после затопления шельфа семь тысяч лет назад и перехода континентальных многолетнемерзлых пород в субаквальное состояние. Последующее потепление и затопление шельфа морской водой привели к росту температуры донных пород и оттаиванию мерзлого слоя. По модельным оценкам мощность реликтовых многолетнемерзлых пород в донных отложениях рассматриваемой области шельфа в современных условиях составила 300 м (рис. 1). Дальнейшее потепление до 2200 г. приводит к ускорению их деградации со стороны верхней границы, которая заглубляется к 2200 г. до глубины 75 м (рис. 1).

Эволюция зоны стабильности гидратов метана и конфигурация ее верхней границы связаны с динамикой термобарических условий, которые определяются динамикой многолетнемерзлых пород и изменения уровня моря. За последние 20 тыс. лет моделирования происходит рост среднегодовой температуры атмосферы на 15°С. Это приводит к заглублению верхней границы ЗСГМ на 20 м (см. рис. 1). Дальнейший рост температуры поверхности донных отложений в период затопления шельфа до –1.8°С усиливает этот процесс. Нарушение условий для формирования зоны стабильности гидратов метана за последние несколько тыс. лет, а следовательно, возможная деградация газовых гидратов на глубине до 160 м ниже дна приводит к скоплению метана в мерзлых горизонтах донных отложений. Существование ММП может также способствовать переходу гидратных залежей в метастабильное состояние за счет эффекта самоконсервации в условиях отрицательной температуры пород [11]. Последующее разрушение мерзлоты, которое происходит в современный период и в будущем, может стать причиной разрушения таких залежей. Верхняя граница ЗСГМ в современный период находится на глубине 160 м ниже морского дна и постепенно заглубляется, что также способствует разрушению газовых гидратов при их существовании на этих глубинах (см. рис. 1).

2. МОДЕЛИРОВАНИЕ СЕЙСМИЧЕСКИХ ВОЛН В СРЕДАХ С ПОГЛОЩЕНИЕМ

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

• скорости распространения P- и S-волн;

• добротности P- и S-волн.

Единственность решения динамической обратной задачи сейсмики для вязкоупругих сред доказана в работе [13], однако в ней не были рассмотрены вопросы чувствительности результатов к помехам в данных и точности проведенных вычислений.

Основной параметр вязкоупругой среды, характеризующий поглощение сейсмической энергии и называемый добротностъю и обозначаемый Q, определяется формулой:

(1)
${{Q}^{{ - 1}}} = \frac{1}{{2\pi }}\frac{{\Delta E}}{E}$
Здесь E задает сейсмическую энергию в единице объема в единицу времени, а ∆E характеризует ее потерю в единице объема за единицу времени. В результате многочисленных натурных экспериментов установлен замечательный факт: добротность в сейсмическом диапазоне частот для очень широкого спектра геологических сред постоянна [14].

В данной работе использовалась наиболее распространенная в настоящее время модель вязкоупругой среды: Обобщенная Стандартная Линейная Модель Твердого Тела (Generalized Standard Linear Solid model или GSLS) [15]. Достоинство этой модели заключается в возможности введения поглощения независимо для продольных и поперечных волн путем использования следующих частотно-зависимых коэффициентов Ламэ:

• Для поперечных волн

(2)
$\mu (\omega ) = {{\mu }_{r}}\left( {1 + \frac{{i\omega {{\tau }^{S}}}}{L}\sum\limits_{l = 1}^L {\frac{{{{\tau }_{l}}}}{{1 + i\omega {{\tau }_{l}}}}} } \right) = {{\mu }_{r}}\left( {1 + {{\tau }^{S}}S(\omega )} \right).$

• Для продольных волн

(3)
$\begin{gathered} \lambda \left( \omega \right) + 2\mu \left( \omega \right) = \left( {{{\lambda }_{r}} + 2{{\mu }_{r}}} \right)\left( {1 + \frac{{i\omega {{\tau }^{P}}}}{L}\sum\limits_{l = 1}^L {\frac{{{{\tau }_{l}}}}{{1 + i\omega {{\tau }_{l}}}}} } \right) = \\ = \left( {{{\lambda }_{r}} + 2{{\mu }_{r}}} \right)\left( {1 + {{\tau }^{P}}S(\omega )} \right). \\ \end{gathered} $

Здесь ${{\mu }_{r}}$ и ${{\lambda }_{r}}$ обозначают коэффициенты Ламэ, соответствующие нулевой частоте, т.е. статической теории упругости, ${{\tau }_{l}}$ задают времена релаксаций для напряжений, обеспечивающие постоянство добротности в заданном частотном диапазоне, а ${{\tau }^{{S,P}}}$ – времена релаксаций для смещений, которые определяют изменчивость добротности по пространству, L – количество релаксационных механизмов. Символ r здесь означает, что данная величина берется на нулевой частоте. Путем прямолинейных вычислений получаем, что добротности P- и S-волн в рамках GSLS выражаются следующим образом:

(4)
${{Q}_{{P,S}}} = \frac{{\sum\limits_{l = 1}^L {\left( {1 + \omega {{\tau }^{{P,S}}}\frac{{{{\omega }^{2}}\tau _{l}^{2}}}{{1 + {{\omega }^{2}}\tau _{l}^{2}}}} \right)} }}{{\omega {{\tau }^{{P,S}}}\sum\limits_{l - 1}^L {\frac{{\omega \tau }}{{1 + {{\omega }^{2}}\tau _{l}^{2}}}} }}$

Параметры ${{\tau }_{l}}$ определяют времена релаксации напряжений и имеют размерности времени. Величины τP,S определяют значения добротности Q для P- и S-волн. Таким образом, построение математической модели вязкоупругой среды для заданного распределения добротностей ${{Q}_{{P,S}}}$ производится в два шага:

1. Вычисление времен релаксаций τl, обеспечивающих постоянство добротности путем минимизации функционала11:

${{\tau }_{l}} = \arg \min {{\left\| {\operatorname{const} - {{Q}_{{P,S}}}} \right\|}_{2}}.$

2. Определение пространственного распределения величин ${{\tau }^{P}}$ и ${{\tau }^{S}}$, обеспечивающих заданное значение добротностей в целевой области.

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

(5)
$\left\{ \begin{gathered} {{\omega }^{2}}\rho {{u}_{x}} + \frac{\partial }{{\partial x}}\left( {\left( {\lambda + 2\mu } \right)\left( {1 + {{S}^{P}}{{\tau }^{P}}} \right)\operatorname{div} \vec {u} - 2\mu \left( {1 + {{S}^{s}}{{\tau }^{S}}} \right)\frac{{\partial {{u}_{z}}}}{{\partial z}}} \right) + \frac{\partial }{{\partial z}}\left( {\mu \left( {1 + {{S}^{S}}{{\tau }^{S}}} \right)\left( {\frac{{\partial {{u}_{x}}}}{{\partial z}} + \frac{{\partial {{u}_{z}}}}{{\partial x}}} \right)} \right) = \hfill \\ = F\left( \omega \right){{G}_{x}}\left( {\vec {x} - {{{\vec {x}}}_{0}}} \right); \hfill \\ {{\omega }^{2}}\rho {{u}_{z}} + \frac{\partial }{{\partial z}}\left( {\left( {\lambda + 2\mu } \right)\left( {1 + {{S}^{P}}{{\tau }^{P}}} \right)\operatorname{div} \vec {u} - 2\mu \left( {1 + {{S}^{S}}{{\tau }^{S}}} \right)\frac{{\partial {{u}_{z}}}}{{\partial z}}} \right) + \frac{\partial }{{\partial x}}\left( {\mu \left( {1 + {{S}^{S}}{{\tau }^{S}}} \right)\left( {\frac{{\partial {{u}_{x}}}}{{\partial z}} + \frac{{\partial {{u}_{z}}}}{{\partial x}}} \right)} \right) = \hfill \\ = F(\omega ){{G}_{2}}(\vec {x} - {{{\vec {x}}}_{0}}), \hfill \\ \end{gathered} \right.$
где ${{u}_{x}}$, ${{u}_{z}}$ – компоненты скорости смещений. Параметры Ламэ здесь зависят от частоты в соответствии с соотношениями (2) и (3). Форма волны, излучаемой источником, определяется функцией $F(\omega )$, а его тип (источник типа сосредоточенной силы, объемного расширения, тензор сейсмического момента и др.) задается функциями ${{G}_{1}}$, ${{G}_{2}}$.

3. ВЛИЯНИЕ ИЗМЕНЧИВОСТИ ДОБРОТНОСТИ НА СЕЙСМИЧЕСКИЕ ВОЛНОВЫЕ ПОЛЯ: ЧИСЛЕННЫЕ ЭКСПЕРИМЕНТЫ

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

Используемое при моделировании распределение скоростей продольных волн можно видеть на рис. 2. Ввиду отсутствия данных о распределении скоростей поперечных волн в этой области, был применен стандартный прием, предполагающий скорость поперечных волн равной скорости продольных, деленной на $\sqrt 3 $. Точно также нам недоступны данные прямых измерений добротности, однако по данными бурения, слой, залегающий на глубине 160–210 м, содержит залежи газогидратов. В исходном замороженном состоянии газогидраты по своим механическим параметрам очень близки к идеально-упругой среде, однако при повышении температуры они становятся все менее и менее консолидированными, пока не произойдет их полное разложение на газ и воду. Естественно, что при этом существенно изменяются и механические свойства породы. В частности, уменьшается добротность упругой среды. Именно это свойство газогидратов, проявляющееся в процессе их растепления, и представляется нам перспективным в качестве прогностического признака приближения высвобождения значительных объемов метана. Для того, чтобы оценить, как сильно уменьшение добротности влияет на сейсмические волновые поля, была выполнена серия численных экспериментов, описывающих поведение сейсмического волнового поля в геологической среде с изменяющейся добротностью в слое на глубине 160–210 м. Расчеты проводились для трех значений добротности: соответствующей идеально-упругой среде (добротность 1000) и двум значениям температуры, найденным в результате численного моделирования (см. раздел 2) для задачи формирования областей с термобарическими условиями, обеспечивающими стабильность газогидратов. Изменчивость добротности мы брали, исходя из полученного в результате численного моделирования прогноза повышения температуры группы слоев с соответствующими термобарическими условиями.

Рис. 2.

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

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

${{Q}_{P}}$, ${{Q}_{S}} \gg 1$, соответствует идеально-упругой среде;

${{Q}_{P}} = 30$, ${{Q}_{S}} = 25$;

${{Q}_{P}} = 20$, ${{Q}_{S}} = 10$.

На рис. 3 приведены синтетические сейсмические трассы, рассчитанные для точечного источника типа вертикальной силы в точке, расположенной на расстоянии 250 м от источника. Как видно, изменчивость добротности проявляется и в динамике, и в кинематике волнового поля. Проведенные нами прямые расчеты показали, что относительное среднеквадратичное уклонение волновых полей при переходе от идеально-упругой среды к вязкоупругой с указанными выше добротностями достигает 10%.

Рис. 3.

Сейсмические трассы для трех значений добротности. Черная – идеально-упругая среда, красная вязкоупругая среда с добротностями ${{Q}_{P}} = 20$, ${{Q}_{S}} = 15$, синяя – вязкоупругая среда с добротностями ${{Q}_{P}} = 15$, ${{Q}_{S}} = 10$.

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

4. СВЯЗАННОСТЬ ИЗМЕНЧИВОСТИ ДОБРОТНОСТИ С ИЗМЕНЧИВОСТЬЮ ДРУГИХ УПРУГИХ ПАРАМЕТРОВ

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

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

Напомним кратко его основные моменты. Заметим, прежде всего, что система уравнений (5) определяет нелинейный оператор, устанавливающий отображение из пространства коэффициентов Ламэ $\lambda $, $\mu $ и времен релаксации ${{\tau }^{P}}$, ${{\tau }^{S}}$ в пространство данных, полученных для заданных источников $({{x}_{s}},{{z}_{s}})$ и приемников $({{x}_{r}},{{z}_{r}})$:

(6)
$\vec {u}({{x}_{r}},{{z}_{r}};{{x}_{s}},{{z}_{s}};\omega ) = B\left[ {\lambda ,\mu ,{{\tau }^{P}},{{\tau }^{S}}} \right]\quad {{\omega }_{1}} \leqslant \omega \leqslant {{\omega }_{2}}.$

Решение уравнения (6) ищется нелинейным методом наименьших квадратов, заключающемся в отыскании минимума целевого функционала:

(7)
$\begin{gathered} \Phi \left[ {\lambda ,\mu ,{{\tau }^{P}},{{\tau }^{S}}} \right] = \\ = {{\left\| {\vec {u}\left( {{{x}_{r}},{{z}_{r}};{{x}_{s}},{{z}_{s}};\omega } \right) - B\left[ {\lambda ,\mu ,{{\tau }^{P}},{{\tau }^{S}}} \right]} \right\|}^{2}}, \\ \end{gathered} $
где в правой части берется ${{L}_{2}}$-норма по всем переменным $({{x}_{r}},{{z}_{r}};{{x}_{s}},{{z}_{s}};\omega )$.

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

(8)
$\begin{gathered} \left[ {{{\lambda }_{{k + 1}}},{{\mu }_{{k + 1}}},\tau _{{k + 1}}^{P},\tau _{{k + 1}}^{S}} \right] = \left[ {{{\lambda }_{k}},{{\mu }_{k}},\tau _{k}^{P},\tau _{k}^{S}} \right] - \\ - \;{{\alpha }_{k}}\nabla \Phi \left[ {{{\lambda }_{k}},{{\mu }_{k}},\tau _{k}^{P},\tau _{k}^{S}} \right], \\ \end{gathered} $
(9)
${{\alpha }_{k}} = \arg \min \Phi \left[ {\left[ {{{\lambda }_{k}},{{\mu }_{k}},\tau _{k}^{P},\tau _{k}^{S}} \right] - \alpha \nabla \Phi } \right].$

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

Для проведения численных экспериментов использовалась модель среды, представленной на рис. 2, в которой будем определять изменчивость добротности в слое на глубине 160–210 м.

В качестве начального приближения мы взяли идеально-упругую среду, т.е. среду без поглощения. Соответствующая трасса на рис. 3 обозначена красным. Далее была инициирована минимизация целевого функционала (7) методом сопряженных градиентов. В качестве критерия остановки использовалась стабилизация процесса, т.е., когда относительное изменение параметров на очередном шаге оказывалось меньше наперед заданной величины, в данном случае ${{10}^{{ - 6}}}$. Как видно из приведенных на рис. 4 результатов, добротность и продольных, и поперечных волн восстановлена с приемлемой точностью. Необходимо подчеркнуть, что весьма важным фактором здесь является отсутствие артефактов связанности параметров. Другими словами, мы видим проявление изменчивости упругих параметров только там, где они и должны быть. В частности, интервалы вариации скоростей распространения волн никоим образом не накладываются на интервалы изменчивости добротности и наоборот. Этот факт имеет принципиальное значение для корректного решения многопараметрической обратной задачи при осуществлении мониторинга, позволяя корректно локализовать изменения каждого из упругих параметров. Причем, как мы видим, корректно определяется не только локализация такой изменчивости, но и ее численные значения. Применительно к мониторингу состояния газогидратов – это именно добротности. Таким образом, применение метода обращения полного волнового поля дает надежную информацию об изменчивости упругих параметров залежей газогидратов.

Рис. 4.

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

5. ЗАКЛЮЧЕНИЕ

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

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

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

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

  1. Shakhova N., Semiletov I., Chuvilin E. Understanding the permafrost–hydrate system and associated methane releases in the East Siberian Arctic Shelf // Geosciences. 2019. 9. 251.

  2. Collett T.S., Lee M.W., Agena W.F., Miller J.J., Lewis K.A., Zyrianova M.V., et al. Permafrost-associated natural gas hydrate occurrences on the Alaska North Slope // Marine and Petroleum Geology. 2011. 28 (2). P. 279–294.

  3. Rekant P., Bauch H.A., Schwenk T., Portnov A., Gusev E., Spiess V., et al. Evolution of subsea permafrost landscapes in Arctic Siberia since the Late Pleistocene: A synoptic insight from acoustic data of the Laptev Sea // Arktos. 2015. 1 (1). 11.

  4. Malakhova V. The response of the Arctic Ocean gas hydrate associated with subsea permafrost to natural and anthropogenic climate changes, IOP Conf. Ser. // Earth and Environmental Sci. 2020. 606. 012035.

  5. Dmitrenko I., Kirillov S., Tremblay L., Kassens H., Anisimov O., Lavrov S., Razumov S., Grigoriev M. Recent changes in shelf hydrography in the Siberian Arctic: Potential for subsea permafrost instability // J. Geophys. Res. 2011. V. 116: C10027.

  6. Golubeva E., Kraineva M., Platov G., Iakshina D., Tarkhanova M. Marine Heatwaves in Siberian Arctic Seas and Adjacent Region // Remote Sens. 2021. 13. 4436.

  7. Gavrilov A., Malakhova V., Pizhankova E., Popova A. Permafrost and gas hydrate stability zone of the glacial part of the East Siberian shelf // Geosciences. 2020. 10. P. 484.

  8. Юрганов Л.Н., Лейфер А. Оценки эмиссии метана от некоторых арктических и приарктических районов по данным орбитального интерферометра IASI // Современные проблемы дистанционного зондирования Земли из космоса. 2016. Т. 13. № 3. С. 173–183.

  9. Malakhova V., Eliseev A. Salt diffusion effect on the submarine permafrost state and distribution as well as on the stability zone of methane hydrates on the Laptev Sea shelf // Ice and Snow. 2020. 60. P. 533–546.

  10. Davies J.H. Global map of Solid Earth surface heat flow // Geochem. Geophyst. Geosyst. 2013. V. 14 (10). P. 4608–4622.

  11. Chuvilin E., Bukhanov B., Davletshina D., Grebenkin S., Istomin V. Dissociation and self‑preservation of gas hydrates in permafrost // Geosciences. 2018. V. 8. P. 431–442.

  12. Sun Y.F., Goldberg D. Hydrocarbon signatures from high-resolution attenuation profiles, SEG Technical Program Expanded Abstracts. 1998. P. 996–999.

  13. Romanov V.G. The two-dimensional inverse problem for the equation of viscoelasticity // Siberian mathematical journal. 2012. V. 53 (6). P. 1401–1412.

  14. Liu H.-P., Anderson D.L., Kanamori H. Velocity dispersion due to anelasticity; implications for seismology and mantle composition // Geophysics. 1976. V. 47. P. 41–58.

  15. Hao Q., Greenhalgh S. The generalized standard-linear-solid model and the corresponding viscoacoustic wave equations revisited // Geophysical Journal International. 2019. 219. P. 1939–1947.

  16. Knyazev A.V., Lashuk I. Steepest Descent and Conjugate Gradient Methods with Variable Preconditioning // SIAM Journal on Matrix Analysis and Applications. 2008. 29 (4). 1267.

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