Строение земной коры, тепловой режим и нефтегазоносность Волго-Уральского сегмента Восточно-Европейского кратона тема диссертации и автореферата по ВАК РФ 00.00.00, кандидат наук Огнев Игорь Николаевич

• Огнев Игорь Николаевич
• кандидат наук
• 2023

Оглавление диссертации кандидат наук Огнев Игорь Николаевич

ВВЕДЕНИЕ

ГЛАВА 1. СОВРЕМЕННЫЕ ПРЕДСТАВЛЕНИЯ О СТРОЕНИИ ЗЕМНОЙ КОРЫ И ВЕРХНЕЙ МАНТИИ И ЕГО СВЯЗЬ С ТЕПЛОВЫМ РЕЖИМОМ И НЕФТЕГАЗОНОСНОСТЬЮ

1.1. Текущие представления о строении земной коры и верхней мантии и методы его

изучения

1.1.1. Строение земной коры и верхней мантии

1.1.2. Основные методы изучения строения земной коры и верхней мантии

1.1.3. Метод гравиметрии для изучения строения земной коры и верхней мантии

1.2. Связь строения земной коры с тепловым режимом и нефтегазоносностью

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

осадочного чехла и строение земной коры и верхней мантии

1.2.2. Прочие факторы, связывающие строение земной коры с нефтегазоносностью

ГЛАВА 2. ОБЪЕКТ ИССЛЕДОВАНИЯ - ВОЛГО-УРАЛИЯ

2.1. Тектоническое строение и эволюция Волго-Уралии как части Восточно-Европейского

кратона

2.2. Современные данные о структуре земной коры Волго-Уралии

2.3. Геотермические исследования на территории Волго-Уралии

ГЛАВА 3. МЕТОДИКА ИССЛЕДОВАНИЯ СТРОЕНИЯ ЗЕМНОЙ КОРЫ И ТЕПЛОВОГО РЕЖИМА ВОЛГО-УРАЛИИ

3.1. Методика гравитационного моделирования

3.1.1. Входные данные для гравитационного моделирования

3.1.2. Инверсия поля градиента силы тяжести

3.1.3. Решение прямой задачи гравиметрии

3.2. Методика геотермического моделирования

3.2.1. Входные данные для геотермического моделирования

3.2.2. Расчет распределения температуры в земной коре и верхней мантии в случае

однослойной земной коры

3.2.3. Расчет распределения температуры в земной коре и верхней мантии в случае

многослойной земной коры

3.2.4. Байесовская инверсия для определения геотермической структуры Волго-Уралии ... 77 ГЛАВА 4. РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЯ СТРОЕНИЯ ЗЕМНОЙ КОРЫ И ТЕПЛОВОГО РЕЖИМА ВОЛГО-УРАЛИИ И ИХ СВЯЗИ С НЕФТЕГАЗОНОСНОСТЬЮ

4.1. Результаты гравитационного моделирования

4.1.1. Результаты инверсии поля градиента силы тяжести

4.1.2. Результаты решения прямой задачи гравиметрии

4.1.3. Сравнение полученной модели с другими региональными моделями Мохо

4.2. Результаты геотермического моделирования

4.2.1. Геотермическая модель Волго-Уралии с однослойной земной корой

4.2.2. Геотермическая модель Волго-Уралии с многослойной земной корой

4.3. Связь строения земной коры с нефтегазоносностью Волго-Уралии и перспективы дальнейших исследований

ЗАКЛЮЧЕНИЕ

СПИСОК СОКРАЩЕНИЙ И УСЛОВНЫХ ОБОЗНАЧЕНИЙ

СПИСОК ЛИТЕРАТУРЫ

СПИСОК ИЛЛЮСТРАТИВНОГО МАТЕРИАЛА И ТАБЛИЦ

ВВЕДЕНИЕ

Введение диссертации (часть автореферата) на тему «Строение земной коры, тепловой режим и нефтегазоносность Волго-Уральского сегмента Восточно-Европейского кратона»

Актуальность темы исследования.

Изучение строения земной коры, ее теплового режима и их влияния на нефтегазоносность является актуальной задачей современной геологии и геофизики. Прежде всего, это связано с развитием современных геофизических методов изучения строения земной коры и верхней мантии. Одним из ключевых методов, развитие которого претерпело качественный скачок в последние два десятилетия, является спутниковая гравиметрия. Преимуществом данного метода является полное и равномерное покрытие земного шара спутниковыми измерениями поля силы тяжести с разрешением, достаточным для проведения региональной оценки структуры земной коры. При этом строение земной коры и верхней мантии взаимосвязано с характером теплового поля [Artemieva, 2019; Lösing, Ebbing, Szwillus, 2020], от которого зависят процессы созревания нефтематеринских пород и, как следствие, нефтегазоносность осадочных бассейнов [Beardsmore, Cull, 2001; Hantschel, Kauerauf, 2009], что обуславливает применение данных спутниковой гравиметрии при региональном прогнозе нефтегазоносности.

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

Одним из перспективных регионов Российской Федерации для пилотного исследования строения земной коры по спутниковым гравиметрическим данным, а также изучения теплового режима и их влияния на нефтегазоносность является Волго-Уральская нефтегазоносная провинция (НГП), охватывающая восточную часть Русской плиты и Предуральский краевой прогиб [Лозин, 2002]. Данный регион обладает достаточно высокой степенью геолого-геофизической изученности, выражающейся в том числе в наличии сейсмических данных о глубинном строении земной коры, а также в большой истории скважинных геотермических измерений. Вместе с тем, до сих пор активна дискуссия о процессах созревания и миграции нефти во многих крупных нефтегазовых месторождений данного региона. Таким образом, изучение строения земной коры Волго-Уралии и последующее геотермическое моделирование является важной задачей, способной пролить свет на взаимосвязь пространственного распределения

нефтегазоносности на данной территории со строением земной коры и неоднородностями ее теплофизических параметров.

Степень разработанности темы исследования.

На настоящий момент времени была показана возможность применения спутниковой гравиметрии для изучения строения земной коры как в глобальном, так и в региональном масштабах. Примером глобальной модели структуры границы Мохо является модель GEMMA, построенная на основании данных спутника GOCE [Reguzzoni, Sampietro, 2015]. Примерами региональных моделей земной коры, полученных в результате инверсии спутниковых гравитационных данных, могут служить модели земной коры Египта [Sobh et al., 2019], Гренландии [Steffen, Strykowski, Lund, 2017], Австралии [Aitken, Salmon, Kennett, 2013], Ирана [Eshagh, Ebadi, Tenzer, 2017], Южной Америки [Meijde van der, Juliá, Assump9äo, 2013], Азиатского орогенного пояса [Guy, Holzrichter, Ebbing, 2017] и многих других областей [Braitenberg, Ebbing, 2009; Haas, Ebbing, Szwillus, 2020; Jiang, Jin, McNutt, 2004; Singh, Rao, 2021; Welford et al., 2010].

Строение земной коры Волго-Уральского региона до настоящего времени изучалось при помощи полевых потенциальных и сейсмических методов. Так, Г.Е. Кузнецов внес существенный вклад в изучение структуры границы Мохо Русской плиты в целом и Татарстана в частности по данным гравиметрических исследований [Кузнецов, 2002; Кузнецов, Боровский, 2000]. М.В. Минц и др. разработали модель Восточно-Европейского кратона (ВЕК), включающую в себя подробное описание строения его Волго-Уральского сегмента, базируясь на данных глубинного сейсмического профиля ТАТСЕЙС-2003 [Минц и др., 2010]. И.М. Артемьева и Г. Тюбо создали модель Европы, Гренландии и Североатлантического региона EUNAseis. Данная модель также включает в себя Волго-Уральский регион и построена исключительно по сейсмическим данным [Artemieva, Thybo, 2013]. Отдельно необходимо выделить фундаментальные работы С.В. Богдановой и А.В. Постникова по изучению вещественного состава земной коры Волго-Уралии [Богданова, 1986; Постников, 2002].

Работы по изучению теплового режима Волго-Уралии активно проводились, начиная со второй половины XX века [Булашевич, Щапов, 1978; Непримеров, 1971; Синявский, Непримеров, Николаев, 1978], когда осуществлялась масштабная программа по измерению теплового потока. В результате было опубликовано более 800 измерений поверхностного теплового потока на территории европейской части СССР [Любимова и др., 1973; Глубинный тепловой поток европейской части СССР, 1974]. Одни из первых моделей распределения температуры внутри осадочного чехла Урало-Поволжья и кристаллического фундамента Татарстана были опубликованы в работах Н.Н. Христофоровой на основании многолетней работы [Христофорова, 2002; Христофорова и др., 2004; Христофорова, Христофоров,

Бергеманн, 2008]. При этом, не до конца раскрыт вопрос латеральных вариаций теплофизических параметров Волго-Уралии, за исключением работы А.А. Липаева и др. по теплопроводности [Липаев, Гуревич, Липаев, 2001] и работ Н.С. Боганика, В.Н. Глазнева и др. и Е.А. Любимовой и др. по радиогенному теплообразованию (РТО) [Боганик, 1975; Глазнев и др., 2021; Любимова, Любошиц, Парфенюк, 1983].

Предыдущие исследования связи нефтегазоносности Волго-Уралии со строением земной коры уделяют большое значение разломам земной коры как фактору, влияющему на размещение нефтяных и газовых месторождений [Постников, 2002; Трофимов, Горюнов, Сабиров, 2017]. Часто ведется дискуссия о возможной роли разломов земной коры в качестве зон восполнения углеводородов (УВ) путем дегазации кристаллического фундамента [Муслимов, Плотникова, 2019; Plotnikova, 2008]. Тем не менее, влияние строения земной коры на нефтегазоносность через тепловое поле литосферы практически не рассмотрено, несмотря на значительную работу по изучению распределения температур в скважинах Поволжья [Христофорова, 2002; Христофорова и др., 2004; Христофорова, Христофоров, Бергеманн, 2008].

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

Цель и задачи исследования.

Основной целью диссертационного исследования является определение структуры земной коры, теплового режима и их влияния на нефтегазоносность Волго-Уральского сегмента ВЕК посредством анализа данных спутниковой гравиметрии и геотермических исследований в скважинах.

Для достижения поставленной цели предусмотрено решение следующих задач:

1. Изучить тектоническое строение Волго-Уральского сегмента ВЕК.

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

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

4. Создать трехмерную плотностную модель земной коры и верхней мантии Волго-Уралии в результате решения прямой задачи гравиметрии.

5. Провести статистический анализ геотермической структуры Волго-Уралии, используя полученную ранее модель земной коры и имеющиеся данные о тепловом поле исследуемого региона.

6. Выделить основные закономерности связи нефтегазоносности Волго-Уральской НГП с вариациями структурных и теплофизических параметров земной коры и верхней мантии региона.

7. Определить направления дальнейшей работы для решения задач уточнения строения земной коры Волго-Уралии, вариаций ее теплофизических свойств и оценки нефтегазоносности Волго-Уральского сегмента ВЕК.

Научная новизна исследования.

1. Впервые использованы спутниковые гравитационные данные для моделирования строения земной коры Волго-Уральского сегмента ВЕК.

2. Построена новая модель земной коры Волго-Уралии с использованием спутниковых гравитационных и глубинных сейсмических данных на изучаемой территории.

3. Произведено моделирование геотермической структуры Волго-Уральского сегмента ВЕК с использованием статистического байесовского подхода методом Монте-Карло по схеме марковских цепей (МСМС) на основе ранее созданной структурной модели.

4. Впервые рассчитаны латеральные вариации таких теплофизических параметров Волго-Уралии, как теплопроводность и радиогенное теплообразование земной коры, теплопроводность верхней мантии и мантийный тепловой поток.

5. Выявлены закономерности связи строения земной коры и верхней мантии с тепловым полем, а также строения земной коры, верхней мантии и теплового поля с нефтегазоносностью Волго-Уральского региона.

Теоретическая и практическая значимость работы.

В ходе работы была построена карта глубины поверхности Мохоровичича, согласующаяся с гравитационным полем и сейсмическими данными о строении земной коры на территории Волго-Уралии. В процессе гравитационного моделирования подтвердилась гипотеза о наличии уплотненного материала в нижней коре центральной Волго-Уралии, соответствующей так называемой «Ветлужской синформе» [Минц и др., 2010; Artemieva, ^уЬо, 2013].

Моделирование геотермической структуры Волго-Уралии позволило оценить пространственные вариации теплофизических параметров земной коры и верхней мантии. Выявлена пространственная связь между радиогенным теплообразованием и распространением архейских комплексов метаосадочных горных пород Волго-Уралии, которые, как было показано А.В. Постниковым, также коррелируют с распределением нефтяных месторождений на территории Волго-Уральской НГП [Постников, 2002]. По результатам моделирования выдвинута

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

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

Методология и методы исследования.

Для достижения цели диссертационного исследования и решения поставленных задач использовались следующие методы:

1. Сбор и оцифровка необходимых данных о строении земной коры и структуре теплового поля Волго-Уралии.

2. Расчет гравитационного эффекта осадочного чехла с учетом кривизны Земли при помощи модели тессероидов и исключение данного эффекта из поля градиента силы тяжести с введенной поправкой за рельеф [Uieda, Barbosa, Braitenberg, 2016].

3. Решение обратной задачи гравиметрии с использованием методики П. Хааса и др. для учета латерально-изменяющегося контраста плотностей земная кора - мантия [Haas, Ebbing, Szwillus, 2020].

4. Решение прямой задачи гравиметрии в программном обеспечении IGMAS+ с одновременным учетом гравиметрических измерений и сейсмических данных о структуре земной коры [Götze, Lahmeyer, 1988; Schmidt et al., 2020].

5. Численное моделирование распределения тепловых параметров земной коры и верхней мантии на основе байесовского подхода MCMC по аналогии с исследованием [Lösing, Ebbing, Szwillus, 2020].

6. Проверка статистических гипотез: о нормальности распределения по критерию Шапиро-Уилка [Shapiro, Wilk, 1965], об однородности дисперсий по критерию Левене [Brown, Forsythe, 1974; Levene, 1960], о различии средних по Т-критерию Стьюдента [Гмурман, 2004; Кремер, 2010] для статистической оценки полученных в процессе моделирования структурных и тепловых параметров земной коры и верхней мантии внутри и вне существующих МУВ.

Вся работа по моделированию и визуализации строения земной коры производилась с использованием следующего программного обеспечения: ArcGIS Pro, Surfer, IGMAS+, язык программирования Python с дополнительным модулем для расчета гравитационного эффекта тессероидов «Tesseroids» [Uieda, Barbosa, Braitenberg, 2016], модулем для визуализации пространственных данных PyGMT [Uieda, Leonardo et al., 2021] и модулем для визуализации

математических графиков Matplotlib [Hunter, 2007]. Инверсия спутникового поля градиента силы тяжести для определения глубины Мохо проводилась с помощью кода П. Хааса и др. [Haas, Ebbing, Szwillus, 2020]. Для байесовского анализа теплового поля и визуализации результатов использовался язык программирования Python с дополнительным модулем для визуализации статистической графики Seaborn [Waskom, 2021], а также модули PyGMT и Matplotlib. В качестве основы для байесовского анализа использовался код М. Лезинг и др. [Lösing, Ebbing, Szwillus, 2020]. Для статистического анализа различия значений полученных структурных и теплофизических параметров земной коры и верхней мантии внутри и за пределами МУВ Волго-Уральской НГП был использован модуль SciPy [Virtanen et al., 2020].

Положения, выносимые на защиту.

1. В центральной части Волго-Уралии выделен слой уплотненной нижней коры, подтверждающий гипотезу о магматическом подслаивании в области «Ветлужской синформы», что было показано по результатам построения 3-D модели земной коры Волго-Уральского региона на основании данных спутниковой гравиметрии.

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

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

Степень достоверности и апробация результатов исследования.

Достоверность предлагаемых автором моделей и выводов основывается на контроле полученных моделей по имеющимся априорным данным и другим аналогичным моделям. Для контроля достоверности полученной в результате решения прямой и обратной задач гравиметрии литосферной модели Волго-Уралии использовались как имеющиеся глобальные и региональные модели земной коры, так и данные глубинных сейсмических исследований. Для контроля достоверности полученной геотермической модели использовались данные температурных измерений в девяти глубоких скважинах [Христофорова, 2002; Христофорова и др., 2004], а также карта температур на глубине 1000 м на территории Поволжья [Христофорова и др., 2004].

Основные положения и результаты исследования докладывались на Международных и Всероссийских конференциях: 73-я Международная молодежная научная конференция «Нефть

и газ 2019» (Москва, 2019 г.), 21-я конференция по вопросам геологоразведки и разработки месторождений нефти и газа «Геомодель 2019» (Геленджик, 2019 г.), Международная междисциплинарная научная геоконференция SGEM (International Multidisciplinary Scientific GeoConference: «SGEM») (Вена, Австрия, 2020 г.), 81-я конференция Немецкого геофизического сообщества (Annual meeting of the German Geophysical Society) (Киль, Германия, 2021 г.), Генеральная Ассамблея Европейского геофизического союза (EGU General Assembly) (Вена, Австрия, 2021 г.), Международная научно-практическая конференция «Решение Европейского союза о декарбонизации и новая парадигма развития топливно-энергетического комплекса России» (Казань, 2021), IV Всероссийская с международным участием школа-конференция студентов, аспирантов и молодых ученых «Материалы и технологии XXI века» (Казань, 2021 г.), 82-я конференция Немецкого геофизического сообщества (Мюнхен, Германия, 2022 г.).

Материалы диссертации опубликованы в 7 научных работах: из них 5 статей в рецензируемых научных изданиях и рекомендованных в диссертационном совете КФУ по специальности 1.6.9 «Геофизика». Среди опубликованных работ индексируются в Scopus - 6, в WoS - 4, опубликованы в журналах, рекомендованных ВАК при Министерстве науки и высшего образования - 2.

Структура и объем работы.

Диссертация состоит из введения, четырех глав, заключения, списка литературы, включающего 270 наименований, иллюстрирована 45 рисунками и содержит 11 таблиц. Общий объем работы составляет 146 страниц.

Личный вклад автора в решение поставленных задач.

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

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

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

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

4. Проведена инверсия спутникового поля градиента силы тяжести, по результатам которой произведена предварительная оценка глубины Мохо на изучаемой территории по методике П. Хааса и др. [Haas, Ebbing, Szwillus, 2020].

5. Построена трехмерная плотностная модель земной коры и верхней мантии Волго-Уралии в процессе решения прямой задачи гравиметрии в программном обеспечении IGMAS+.

6. Собраны необходимые данные для геотермического моделирования, куда входят температура и тепловой поток на земной поверхности, температура на границе литосфера -астеносфера (LAB).

7. Дополнен программный код для байесовской инверсии геотермических параметров [Lösing, Ebbing, Szwillus, 2020] с добавлением возможности расчёта геотермических параметров для многослойной модели земной коры.

8. Проведен статистический анализ геотермической структуры Волго-Уралии с использованием полученной ранее модели литосферы и имеющихся данных о тепловом поле исследуемого региона.

9. Рассмотрены статистические различия строения земной коры и верхней мантии и значений теплофизических параметров в пределах Волго-Уральской НГП для нефтегазоносных регионов и регионов без открытых МУВ.

10. Определены направления дальнейшей работы по изучению теплофизических свойств земной коры и анализу нефтегазоносности Волго-Уральского сегмента ВЕК.

Благодарности.

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

Особую благодарность автор выражает профессору, доктору Йоргу Эббингу и его коллегам - Питеру Хаасу и Марин Лезинг из Кильского университета имени Кристиана Альбрехта за совместную работу, ценные комментарии и рекомендации как в процессе выполнения самой работы, так и в процессе публикации результатов исследования.

Отдельная благодарность Г.С. Хамидуллиной, Н.Н. Равиловой, З.М. Слепаку, Д.И. Хасанову, Э.В. Утемову, Е.А. Ячменевой, Э.Р. Зиганшину, П.С. Крылову и всему коллективу кафедры геофизики и геоинформационных технологий института геологии и

нефтегазовых технологий Казанского (Приволжского) федерального университета за поддержку в написании, оформлении и корректировке диссертации.

ГЛАВА 1. СОВРЕМЕННЫЕ ПРЕДСТАВЛЕНИЯ О СТРОЕНИИ ЗЕМНОЙ КОРЫ И ВЕРХНЕЙ МАНТИИ И ЕГО СВЯЗЬ С ТЕПЛОВЫМ РЕЖИМОМ И НЕФТЕГАЗОНОСНОСТЬЮ

1.1. Текущие представления о строении земной коры и верхней мантии и методы его

изучения

1.1.1. Строение земной коры и верхней мантии

На протяжении многих десятилетий считается, что процессы, происходящие в верхней мантии, влияют на формирование земной коры, поэтому изучение земной коры и верхней мантии зачастую происходит в тесной взаимосвязи друг с другом [Hawkesworth et al., 2017; Mooney, 2015; The Earth's Crust and Upper Mantle, 1969]. Данный факт обоснован и тем, что земная кора и верхняя часть верхней мантии имеют схожие механические свойства и обладают кондуктивным механизмом теплопереноса, формируя литосферу - твердую оболочку Земли, не подверженную конвекции [Hawkesworth et al., 2017]. В этой связи, далее кратко рассмотрим особенности строения литосферы в целом и земной коры и литосферной верхней мантии в частности.

Строение литосферы.

Литосфера по своему строению подразделяется на океаническую и континентальную. Океаническая литосфера погружается в мантию в зонах субдукции, в связи с чем цикл ее жизни составляет ~200 млн лет, в то время как континентальная литосфера более стабильна и может достигать возраста до 4.4 млрд лет [Artemieva, 2011; Hawkesworth et al.,2017]. Изучение строения континентальной литосферы исключительно важно как с научной точки зрения для понимания процессов эволюции Земли, тектоники плит и взаимодействия литосферы и мантии, так и с практической точки зрения для понимания процессов минералообразования и геологических угроз [Artemieva, 2011]. В зависимости от применяемых геофизических методов и изучаемых физических свойств, понятие литосферы может определяться по-разному, от чего будет зависеть и мощность литосферы. Различают эластическую, термальную, сейсмическую, электрическую и петрологическую литосферу [Artemieva, 2011]. Наиболее часто используемым подходом к определению литосферы на практике является понятие термальной литосферы.

Термальная литосфера - это верхняя оболочка Земли, в которой кондуктивный теплоперенос преобладает над конвективным, доминирующим в залегающей под ней астеносфере [Allen Р., Allen J., 2013; Artemieva, 2011]. При этом за подошву термальной литосферы может приниматься температура 0.8-0.9 от мантийного солидуса - температуры плавления пород мантии [Pollack, Chapman, 1977]. Также за глубину термальной литосферы некоторые исследователи принимают точку пересечения продолжения кондуктивного

температурного профиля с мантийной адиабатой [Artemieva, 2011; Jaupart, Mareschal, 1999]. Мантийная адиабата является температурным профилем конвективной мантии, показывающим увеличение температуры пород мантии с глубиной под воздействием увеличивающегося горного давления перекрывающих толщ [Turcotte, Schubert, 2014]. Таким образом, за глубину термальной литосферы на практике зачастую принимается изотермическая граница с температурой в диапазоне 1250-1350 °С [Allen Р., Allen J., 2013; Artemieva, 2011] (Рисунок 1).

Рисунок 1. Схематическое изображение температурного профиля (геотермы) в земной коре и

верхней мантии по [Artemieva, 2011; Jaupart, Mareschal, 1999; Turcotte, Schubert, 2014]. Мантийная адиабата и мантийный солидус построены по данным Т. Кацуры [Katsura, 2022] и

П.Ф. Рея [Rey, 2015] соответственно. Температурный профиль в кондуктивной литосфере построен по уравнениям (3) и (13) для земной коры и верхней мантии соответственно (М=40 км, Ç0=60 мВт/м2, Г0=10 °С, ÇM=30 мВт/м2, кс=2.5 Вт/м/К, Ас=1 мкВт/м3)

Глубина термальной литосферы связана с ее возрастом: чем древнее литосфера, тем она холоднее и тем глубже залегает ее подошва [Artemieva, 2006; Artemieva, Mooney, 2001]. Увеличенную мощность древней архей-палеопротерозойской литосферы связывают с активным архейским магматизмом, ведущим к истощению пород литосферной мантии мафическими компонентами, такими как FeO, Al2O3 и CaO [Jordan, 1978; Pollack, 1986]. Истощение верхней мантии железом ведет к формированию более насыщенной магнием вариацией перидотита -гарцбургита, который не только менее плотный, чем неизмененный перидотит, но также и менее

подверженный плавлению, что делает его устойчивым при больших температурах [Allen Р., Allen J., 2013]. Другой механизм, ведущий к увеличению мощности архейской континентальной литосферы - это истощение верхней мантии летучими компонентами или деволатилизация, которая приводит к возрастанию линии солидуса для пород верхней мантии на 300-600 °С, что также затрудняет их плавление [Pollack, 1986]. Оба описанных механизма ведут к формированию регионов архейской литосферы с повышенной мощностью кондуктивного слоя, которая может достигать до 300-400 км по сравнению с 100-200 км в подвижных более молодых зонах за их пределами [Artemieva, 2011]. Такие участки с мощной кондуктивной литосферой играют роль барьеров для конвективного теплопереноса и отводят конвективные мантийные потоки к своей периферии. С учетом того, что конвекция является более эффективным механизмом теплопереноса, чем простой кондуктивный теплоперенос, мобильные зоны на периферии архейской литосферы будут получать повышенный глубинный тепловой поток за счет мантийной конвекции, действующей здесь на меньших глубинах [Ballard, Pollack, 1987]. Отвод конвективного тепла и пониженное радиогенное теплообразование архейской земной коры по сравнению с более молодыми зонами объясняют наблюдаемый контраст поверхностного теплового потока, который равен ~ 41 ± 11 мВт/м2 на архейских областях по сравнению с 46 ± 15 мВт/м2 и 49 ± 16 мВт/м2 в мезо- и неопротерозойских областях соответственно и с >60 мВт/м2 в фанерозойских регионах [Artemieva, 2011; Ballard, Pollack, 1987].

Список литературы диссертационного исследования кандидат наук Огнев Игорь Николаевич, 2023 год

Источники

Сейсмические данные

Глобальный сейсмический каталог USGS Глубинный профиль TATSEIS-2003 Глубинный профиль URSEIS-95

Глубинный профиль UWARS

Chulick, Detweiler, Mooney [2013] Trofimov [2006]

Tryggvason, Brown, Pérez-Estaún [2001], Пучков [2010] Thouvenot et al. [1995]

Гравитационные данные

Гравитационные градиенты GOCE Глобальная гравитационная модель XGM2019e Bouman et al. [2016] Zingerle et al. [2019]

Структурные данные

Рельеф земной поверхности ETOPO1 Мощность осадочного чехла EUNAseis Поверхность LAB Amante, Eakins [2009] Artemieva, Thybo [2013] Artemieva [2019]

Петрофизические данные

Данные о плотности осадков, земной коры и Artemieva [2007]

верхней мантии

Примечание: составлено автором.

Сейсмические данные.

Данные сейсмических измерений являются главным источником информации о мощности земной коры, поэтому они играют ключевую роль в гравитационном моделировании, так как ограничивают возможные оценки глубины Мохо в местах проведения сейсмической съемки. В настоящей работе были использованы сейсмические данные о глубине Мохо в пределах исследуемого региона по глобальному сейсмическому каталогу USGS [Chulick, Detweiler, Mooney, 2013]. В данном каталоге содержится информация о мощности земной коры по результатам глубинных сейсмических исследований, проведенных на территории Русской платформы преимущественно в советское время. В дополнение к используемому каталогу были оцифрованы сейсмические профили, проложенные в конце XX и начале XXI веков на территории Волго-Уральского субкратона. К ним относится геотраверс TATSEIS-2003 [Trofimov, 2006], проходящий через центральную часть субкратона, глубинные профили URSEIS-95, ESRU и UWARS, которые, пересекая Уральские горы, маркируют структуру земной коры восточной границы Волго-Уралии [Пучков, 2010; Рыбалка и др., 2006; Brown et al., 2002a; Thouvenot et al., 1995; Tryggvason, Brown, Perez-Estaun, 2001]. Оценки глубины границы Мохо, согласно используемым сейсмическим данным, представлены на Рисунке 16.

Рисунок 16. Обрамление исследуемого региона с сейсмическими данными о глубине границы Мохо. Рельеф местности приведен согласно модели ETOPO1 [Amante, Eakins, 2009]. Сейсмические оценки глубины Мохо приведены по данным глобального сейсмического каталога USGS и геотраверсов TATSEIS-2003, URSEIS-95, ESRU и UWARS

Гравитационные данные.

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

тяжести был использован грид поля градиента силы тяжести по данным спутника GOCE на высоте 225 км [Bouman et al., 2016]. Решение использовать спутниковые измерения компонент силы тяжести, а именно градиент силы тяжести по данным спутника GOCE обосновано тремя главными доводами. Во-первых, в отличие от наземных измерений, спутниковые измерения поля силы тяжести равномерно покрывают практически всю поверхность Земли, включая изучаемый регион. Во-вторых, на наземные измерения поля силы тяжести большое влияние оказывают внутрикоровые плотностные неоднородности, что усложняет интерпретацию глубинного строения земной коры [Bouman et al., 2015]. В-третьих, как показали Й. Боуман и др. (2015), градиенты поля силы тяжести на высоте полета спутника чувствительны к поверхностям с большим контрастом плотностей, таким как граница Мохо [Bouman et al., 2015]. Они решили прямую задачу гравиметрии для литосферной модели Северо-Восточной окраины Атлантического океана и оценили чувствительность различных компонент тензора градиента силы тяжести к вариациям плотностей как внутри земной коры для слоя модели, расположенного на глубинах от 7.5 до 10 км, так и в переходной зоне между земной корой и верхней мантией для слоя, расположенного на глубине от 30 до 35 км. В результате выяснилось, что компоненты тензора градиента силы тяжести на высоте орбиты спутника практически нечувствительны к внутрикоровым вариациям плотностей и, в то же время, демонстрируют сильные вариации в зависимости от мощности земной коры, что в особенности выражено в поле вертикального градиента силы тяжести.

Помимо гравитационной инверсии, грид Bouman et al. [2016] с введённой поправкой за рельеф был использован и при решении прямой задачи гравиметрии наряду с аномалией поля силы тяжести на земной поверхности в редукции Буге согласно глобальной гравитационной модели XGM2019e [Zingerle et al., 2019].

Структурные данные.

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

Рельеф земной поверхности и мощность осадочного чехла необходимы для внесения поправок за рельеф и за осадочный чехол для подготовки данных к инверсии. В настоящем исследовании была использована глобальная модель рельефа ETOPO1 с разрешением в 1 угловую минуту. Мощность осадочного чехла была получена по сейсмической модели EUNAseis, обладающей разрешением 2 градуса [Amante, Eakins, 2009; Artemieva, Thybo, 2013]. Обе модели были перемасштабированы на географическую сетку с разрешением 1 градус перед проведением инверсии поля градиента силы тяжести.

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

Петрофизические данные.

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

Плотность осадочного чехла описывается экспоненциально возрастающей с глубиной функцией, полученной для отложений ВЕК [Artemieva, 2007]. Плотности верхней и нижней коры были взяты, согласно сейсмическим оценкам, по модели CRUST 1.0 [Laske et al., 2013].

Таблица 6 - Плотностная модель, используемая в диссертации

Слой_Плотность, кг/м3_

Осадочный чехол 2430xz0■045*

Верхняя кора 2750

Нижняя кора 2900

Верхняя мантия 3230

Астеносфера_3225_

*z - У мощности осадочного чехла в км

Примечание: составлено автором.

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

= Л 1е + ТО\

Рт in situ = Рту1 а 2 )' '

где рт - плотность литосферной мантии в стандартных условиях, кг/м3; Тм - температура на границе Мохо, °С; Т0 - температура на границе LAB, °C; а - коэффициент теплового расширения, °C-1.

В настоящем исследовании учитывается тот факт, что архейская литосферная мантия обеднена мафическими компонентами, что приводит к понижению ее плотности [Kaban et al., 2003]. Плотность литосферной мантии ВЕК в стандартных условиях принята за 3340 кг/м3, что соответствует архей-палеопротерозойскому возрасту согласно [Artemieva, 2007]. Температура на границе Мохо принята за 500 °С, что входит в температурный диапазон архей-палеопротерозойской земной коры ВЕК согласно [Artemieva, 2007] и температура на границе LAB - за 1400 °C согласно модели [Artemieva, 2019], которая была использована в настоящей работе. Необходимо отметить, что в природе температура на границе Мохо должна изменяться в зависимости от мощности земной коры [Mareschal, Jaupart, 2013]. Поэтому постоянная температура является упрощением, принятым в данной работе. Коэффициент теплового расширения взят за 3.5 х 10-5 °C-1 [Artemieva, 2007; Artemieva, 2019]. Используя данные параметры, была рассчитана плотность литосферной мантии в естественных условиях, которая составила 3230 кг/м3.

С учетом того, что температура на границе Мохо не влияет на тепловое расширение пород астеносферы, можно рассчитать плотность астеносферной мантии в естественных условиях, приняв температуру в астеносфере за температуру на границе LAB:

ра in situ = ра(1 - аТо), (8)

где ра - плотность астеносферной мантии в стандартных условиях, кг/м3, которая была принята за 3390 кг/м3 [Artemieva, 2007].

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

3.1.2. Инверсия поля градиента силы тяжести

Обработка спутниковых гравитационных данных.

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

тяжести, внести поправку за свободный воздух и поправку за рельеф (например [Sobh et al., 2019; Steffen, Strykowski, Lund, 2017]). В настоящей работе используется вертикальный градиент поля силы тяжести с плотностью промежуточного слоя 2670 кг/м3 и плотностью морской воды 1030 кг/м3. Тем не менее, помимо вышеназванных поправок, следует учитывать эффект осадочного чехла на поле силы тяжести, который тем больше, чем больше мощность осадков, и может вызывать аномалии поля силы тяжести на высоте полета спутника, сравнимые с таковыми от границы Мохо [Steffen, Strykowski, Lund, 2017].

Несмотря на то, что Волго-Уралия обладает неизменной мощностью осадочного чехла в своей кратонной части, осадочный чехол граничащей с ней Прикаспийской низменности и расположенных на востоке структур Камско-Бельского и Серноводско-Абдулинского авлакогенов достигает мощностей до 10-20 км (Рисунок 12) [Шаргородский и др., 2004; Artemieva, Thybo, 2013; Neprochnov, Kosminskaya, Malovitsky, 1970]. По этой причине важно дополнительно внести поправки за осадочный чехол в скорректированное за рельеф поле вертикального градиента силы тяжести. Это позволит получить аномалию градиента силы тяжести, вызванного преимущественно положением границы Мохо:

GM = GTC — Gsed, (9)

где GM - градиент силы тяжести, подготовленный для инверсии, отражающий преимущественно вариации положения границы Мохо, этвеш;

GTC - градиент силы тяжести, скорректированный за рельеф местности, этвеш;

Gsed - градиент силы тяжести, вызванный осадочным чехлом, этвеш.

Учитывая, что моделируемая область занимает географически значительную территорию, большое влияние на результаты решения обратной задачи гравиметрии будет оказывать сферичность Земли [Долгаль и др., 2018]. В текущей работе было принято решение использовать модель тессероидов - сферических призм, для учета эффекта сферичности Земли [Uieda, Barbosa, Braitenberg, 2016]. Во-первых, была рассчитана абсолютная глубина подошвы осадочного чехла на сетке разрешением 1x1° путем вычитания мощности осадочного чехла по модели EUNAseis из рельефа земной поверхности по модели ETOPO1. Во-вторых, вся толща осадочного чехла была разбита на тессероиды с латеральными размерами 1x1° и вертикальной мощностью 1 км. В третьих, каждому тессероиду была назначена своя плотность в зависимости от его глубины согласно зависимости плотности осадочных пород от глубины, полученной для ВЕК [Artemieva, 2007]:

р = 2430 • z0 045, (10)

где z - У глубины середины тессерода, км.

Последним шагом был расчет гравитационного эффекта осадочного чехла на вертикальный градиент силы тяжести с использованием библиотеки Tesseroids в Python. Данный эффект в дальнейшем был вычтен из поля градиента силы тяжести с введенной поправкой за рельеф (Рисунки 17 и 18).

36° Е 42° Е 48° Е 54° Е 60° Е 36° Е 42° Е 48° Е 54° Е 60° Е

(а) 1.........=Е (б) ^^......Щ

0 2 4 6 8 10 12 14 16 18 20 -1.80 -1.55 -1.30 -1.05 -0.80 -0.55 -0.30 -0.05 0.20

Мощность осадочного чехла [ЕиЫА5е1з] (км) Градиент силы тяжести осадочного чехла (этвеш)

Рисунок 17. Мощность осадочного чехла и его влияние на градиент гравитационного поля. (а) Мощность осадочного чехла на территории Волго-Уралии по данным региональной сейсмической модели EUNAseis [Artemieva, ^уЬо, 2013]. (б) Аномалия вертикального градиента силы тяжести, создаваемого осадочным чехлом модели EUNAseis, дискретизированным тессероидами

36° Е 42° Е 48° Е 54° Е 60° Е 36° Е 42° Е 48° Е 54° Е 60° Е

-0.75 -0.50 -0.25 0.00 0.25 -0.50 -0.25 0.00 0.25 0.50 0.75 1.00 1.25 1.50

Градиент силы тяжести с поправкой за рельеф (этвеш) Скорректированный градиент силы тяжести (этвеш)

Рисунок 18. (а) Вертикальный градиент гравитационного поля на территории Волго-Уральского региона, скорректированный за рельеф и (б) градиент гравитационного поля, скорректированный за рельеф и осадочный чехол

Инверсия поля градиента силы тяжести с латерально изменяющимся контрастом плотности.

Для целей инверсии поля градиента силы тяжести был использован новый подход из [Haas, Ebbing, Szwillus, 2020], позволяющий контрасту плотности между земной корой и верхней мантией изменяться латерально согласно с тектоническими регионами, представленными на исследуемой территории. Этот подход заключается в решении обратной задачи с помощью алгоритма Гаусса - Ньютона с тихоновской регуляризацией второго порядка для обеспечения устойчивости решения и требует наличия двух гиперпараметров для инверсии: опорной глубины Мохо zrgf и контраста плотности земная кора - мантия Ар c диапазоном возможных значений. Алгоритм итеративно устанавливает каждый контраст плотности из заданного диапазона с указанным шагом предопределенным тектоническим регионам, сохраняя Ар пространственно неизменным в каждом регионе. Таким образом проверяются все возможные комбинации латерального распределения контрастов плотности. Затем автоматически выбирается комбинация, дающая наименьшую среднеквадратичную ошибку между глубиной Мохо, оцененной с помощью инверсии, и глубиной Мохо, определенной в точках доступных сейсмических измерений.

Хотя в вышеупомянутом алгоритме можно использовать любую производную гравитационного потенциала для инверсии, в настоящем исследовании было принято решение использовать вертикальный градиент силы тяжести на высоте полета спутника, поскольку показано, что он более чувствителен к колебаниям глубины Мохо, чем другие компоненты гравитационного потенциала [Bouman et al., 2016]. Здесь для целей тектонического районирования были использованы главные тектонические структуры Волго-Уралии по Богдановой и др. [Bogdanova, Gorbatschev, Garetsky, 2016], которые включают архейскую кратонную континентальную кору и палеопротерозойские подвижные пояса. Также в отдельную тектоническую область был выделен Уральский ороген из-за его относительно молодого возраста и особого состава земной коры. Еще один, заслуживающий внимания, тектонический регион - Прикаспийский осадочный бассейн. Осадочные толщи в его центральной части достигают 20 км общей мощности и включают слои пермских солей мощностью около 4-5 км [Brunet et al., 1999; Volozh, Talbot, Ismail-Zadeh, 2003], что является уникальным для ВЕК. Соответственно, данный регион был выделен как четвертая тектоническая единица, используемая для гравитационной инверсии.

Для контрастов плотности был выбран диапазон от 350 до 550кг/м3 с шагом 50кг/м3. Такой диапазон был выбран в соответствии с предыдущими выводами для региона на основе спутниковых гравиметрических [Eshagh, Ebadi, Tenzer, 2017] и сейсмических исследований [Rabbel, Kaban, Tesauro, 2013]. По их результатам контраст плотности для данного региона должен составлять от 300 до 600 кг/м3. В итоге алгоритм инверсии был запущен 10 раз с использованием разных опорных глубин Мохо в диапазоне от 41 до 50 км с шагом в 1 км. Такой диапазон был выбран в соответствии со средней глубиной Мохо по имеющимся сейсмическим исследованиям, которая составляет 45 км. Перед выполнением инверсии пространственное разрешение используемых сейсмических глубин Мохо на оцифрованных профилях было уменьшено, чтобы оно соответствовало разрешению сейсмического каталога USGS в районе исследования. Таким образом, разрешение было уменьшено в 4 раза: с 10 до 40 км. Наконец, на основании оценки среднеквадратичных ошибок, была выбрана модель Мохо, которая менее всего расходилась с данными сейсмических исследований. Все расчеты при инверсии гравитационного поля проводились на сетке 1*1°.

3.1.3. Решение прямой задачи гравиметрии

За гравитационной инверсией последовало решение прямой задачи гравиметрии, которое было выполнено с помощью программного обеспечения IGMAS+ [Götze, Lahmeyer, 1988; Schmidt et al., 2020]. IGMAS+ - это геофизический пакет, предназначенный для трехмерного численного моделирования, визуализации и интерпретации потенциальных полей. Данный пакет

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

Перед началом моделирования площадь счёта была увеличена на 2500 км для минимизации краевых эффектов. Это было сделано путем продления мощности смоделированных слоев от краев исследуемой области. Размеры моделируемой области за исключением растянутой части составляют 2672 км в направлении запад - восток х 3236 км в направлении север - юг. Это означает, что проведенное латеральное растяжение модели составляет примерно 95 % ее размера по широте и 75 % по долготе. Глубина модели по вертикали была задана равной 3 00 км, чтобы охватить главные границы, вдоль которых возникают основные контрасты плотности, начиная от подошвы осадочного чехла и заканчивая LAB. Трехмерная модель состоит из многогранников, триангулированных между 67 вертикальными сечениями, ориентированными в направлении запад - восток. Примерное расстояние, разделяющее сечения, составляет 50 км. IGMAS+ позволяет рассчитать поле силы тяжести по имеющейся плотностной модели и провести перекрестное сравнение данного поля с измеренным полем силы тяжести, что позволяет улучшить модель в местах значительного несоответствия рассчитанного и измеренного полей. Полный процесс решения прямой задачи гравиметрии в текущем исследовании можно разделить на пять шагов:

1. Импорт сейсмических, структурных и гравиметрических данных в IGMAS+: (а) граница Мохо, полученная ранее в результате гравитационной инверсии, (б) глубина осадочного чехла по модели EUNAseis [Artemieva, Thybo, 2013], (в) глубина LAB, полученная методом термальной изостазии [Artemieva, 2019], (г) доступные сейсмические измерения глубины Мохо из сейсмического каталога USGS, сейсмических профилей TATSEIS, URSEIS, UWARS и ESRU [Brown et al., 2002b; Chulick, Detweiler, Mooney, 2013; Thouvenot et al., 1995; Trofimov, 2006; Tryggvason, Brown, Pérez-Estaun, 2001], (д) аномалия поля силы тяжести в редукции Буге согласно глобальной модели XGM2019e [Zingerle et al., 2019], (е) аномалия поля градиента силы тяжести с поправкой за рельеф согласно гридам из [Bouman et al., 2016]. Дополнительно земная кора была разделена на верхнюю и нижнюю части с исходно горизонтальной границей раздела. Плотности всех слоев задавались согласно значениям из Таблицы 6. Осадочный слой был дискретизирован на ряд изометрических вокселей мощностью 500 м. Это позволило представить экспоненциальное увеличение плотности осадочного чехла с глубиной.

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

3. Были рассчитаны поля силы тяжести и градиента силы тяжести по текущей модели и была обнаружена область значительного расхождения рассчитанного и измеренного полей, составляющего порядка 95 мГал в центральной части Волго-Уральского субкратона. Это расхождение было связано с телом андерплейтинга относительно большей плотности, расположенного в нижней части коры (см. раздел 4.1.).

4. Был оценен дисбаланс массы (избыток и дефицит) в этом районе с помощью изостатических расчетов, следуя подходу Й. Эббинга [Ebbing, 2007] для Скандинавской горной цепи:

5

psedDsed + pucduc + PlcDlc + PrnDm + PaDa - ^ PrefiDrefi = ALoad/g, (11)

¿=1

где p и D - плотности в кг/м3 и мощности в м осадочного чехла, верхней и нижней земной коры, литосферной мантии и астеносферы геологической модели IGMAS;

prefi и Dreß - плотности в кг/м3 и мощности в м опорной модели, которые равны средним значениям этих параметров, используемых в соответствующих слоях геологической модели; ALoad - изостатическая нагрузка, отражающая избыток и дефицит массы в исследуемом регионе, кг/мс2;

g - нормальное поле силы тяжести, м/с2.

Плотность андерплейтинга была принята равной 3100 кг/м3, что дает разницу между менее плотной нижней корой и предполагаемым телом в -200 кг/м3. Полученный дисбаланс массы из уравнения (11) был разделен на данное значение. Таким образом, была получена мощность высокоплотного нижнего слоя земной коры, связанного с андерплейтингом.

5. Последним шагом было изменение геометрии слоев для достижения хорошей сходимости измеренных и рассчитанных гравитационных полей. Здесь граница Мохо и граница между верхней и нижней корой подверглись дальнейшим модификациям. Граница между верхней и нижней корой была изменена таким образом, чтобы обеспечить лучшую гравитационную сходимость и согласовываться с подошвой «кисло-средней» коры из модели

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

3.2. Методика геотермического моделирования

Восстановление термической структуры литосферы является одной из основных задач, решаемых геофизическими методами. Возможны различные подходы к решению этой задачи в зависимости от имеющихся данных и преследуемых целей. В стабильных кратонных регионах часто принимается допущение о чисто кондуктивном теплопереносе в литосфере (например, [Afonso et al., 2016; Artemieva, Mooney, 2001; Chapman, 1986]). Это позволяет использовать измерения поверхностного теплового потока для расчета геотерм земной коры [Förster et al., 2010], мантийного теплового потока [Cermâk, 1982], мощности термальной литосферы [Artemieva, Mooney, 2001] и других геотермических характеристик земной коры [Perry et al., 2006]. В качестве альтернативы, когда данных о поверхностном тепловом потоке недостаточно или когда теплоперенос в литосфере нельзя считать полностью кондуктивным, можно полагаться на скорости сейсмических волн для определения термической структуры литосферы благодаря связи между упругими и тепловыми свойствами горных пород [Artemieva et al., 2004]. Таким образом Д.Л. Шутт и др. определили температуру на границе Мохо по скоростям продольных волн для западной части США [Schutt, Lowry, Buehler, 2018], а Г. Диаферия и др. исследовали термическую структуру земной коры вдоль итальянского полуострова [Diaferia et al., 2019; Diaferia, Cammarano, Faccenna, 2019]. Кроме того, глубина залегания изотермы Кюри [Li, 2011; Wang, Li, 2015], геотермобарометрия ксенолитов [Kukkonen, Peltonen, 1999], топографические высоты [Artemieva, 2019; Hasterok, Chapman, 2007] и другие геофизические наблюдения (например, [Adâm, 1980; Kukkonen, Jöeleht, 1996; Zeyen, Dérerovâ, Bielik, 2002]) могут быть полезны в качестве дополнительных источников информации для исследования термической структуры литосферы. Вдобавок, химический состав изверженных основных пород может использоваться для определения температуры мантии и глубины термальной литосферы, как это было показано в новом вероятностном подходе, разработанном Б. Оливейра и др. [Oliveira, Afonso, Klöcking, 2021].

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

тепловым потоком и температурой на поверхности LAB согласно уравнениям (3) и (13). Некоторые из перечисленных термических параметров являются неизвестными, что делает невозможным аналитическое нахождение их значений. По этой причине для определения неизвестных термических параметров Волго-Уралии и их латеральных вариаций был использован численный подход на основании байесовского метода MCMC, предложенный М. Лезинг и др. (2020) при анализе геотермической структуры Антарктиды [Lösing, Ebbing, Szwillus, 2020]. Данный подход позволяет найти наиболее вероятные значения термических параметров, задействуя имеющуюся информацию о строении литосферы и тепловом режиме исследуемого региона. В этой связи, проведенная в настоящем исследовании работа по анализу геотермической структуры Волго-Уральского региона может быть охарактеризована тремя основными этапами (Рисунок 19).

Первый этап заключается в сборе имеющихся структурных и термических данных для инверсии методом MCMC. Структурные данные включают в себя глубину границы LAB ZLAB и мощности осадочного чехла, верхней и нижней земной коры (HSED, Нис и HLC). К термическим данным относятся температура на поверхности LAB, температура и тепловой поток на земной поверхности (TLAB, Т0 и Ç0). На втором этапе выполняется алгоритм инверсии методом MCMC для определения возможных латеральных вариаций термических параметров на территории исследования. Список искомых термических параметров будет зависеть от типа принятой геотермической модели. В настоящей работе были рассмотрены две модели: (1) геотермическая модель с однослойной земной корой, (2) геотермическая модель с многослойной земной корой. Таким образом, для модели с однослойной корой искомыми термическими параметрами являются теплопроводность земной коры и верхней мантии и радиогенное теплообразование земной коры (кс, км, Ас). Для модели с многослойной корой в качестве неизвестных термических параметров выступают РТО и теплопроводности осадочного чехла, верхней и нижней земной коры, а также теплопроводность верхней мантии (ASED, Аис, ALC, kSED, кис, kLC, км). На третьем этапе производилась проверка полученных геотермических моделей по имеющимся температурным данным, полученным по прямым измерениям в скважинах.

ВХОДНЫЕ СТРУКТУРНЫЕ ДАННЫЕ

ZlAB> t¡SED> Нцс> Hlc

I ~

ВХОДНЫЕ ТЕРМИЧЕСКИЕ ДАННЫЕ

Tlab> ТО, Qо

Байесовская инверсия Монте-Карло

с марковскими цепями

•...............................»................................

Распределение неизвестных термических параметров для одно- и многослойной моделей земной коры

Проверка геотермических моделей по измеренным температурам земной коры

Рисунок 19. Схематический рабочий процесс геотермического моделирования

3.2.1. Входные данные для геотермического моделирования

Структурные и термические данные, используемые в моделировании, приведены в Таблице 7.

Таблица 7 - Структурные и термические данные, используемые в процессе геотермического

моделирования

Наборы данных Источники

Структурные данные

Глубина LAB Artemieva [2019]

Мощность верхней и нижней коры Ognev, Ebbing, Haas [2022]

Мощность осадочного чехла Artemieva, Thybo [2013]

Термические данные

Тепловой поток на земной Jennings, Hasterok, Lucazeau [2021]

поверхности

Температура на земной поверхности Fick, Hijmans [2017]

Температура на глубине 1000 м Христофорова и др. [2004]

Температурные профили в скважинах Христофорова [2002], Христофорова и др. [2004]

Температура на границе LAB Artemieva [2019]

Примечание: составлено автором.

Структурные данные.

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

коры. Глубина границы LAB заимствована из модели LAB Европы, построенной И.М. Артемьевой (2019) на основании метода термальной изостазии [Artemieva, 2019]. Данная модель хорошо согласуется с измерениями наземного теплового потока и, как показано выше, уже была успешно использована для построения 3-D литосферной модели Волго-Уралии. Мощность осадочного чехла взята из сейсмической модели EUNAseis [Artemieva, Thybo, 2013]. В качестве мощностей верхней и нижней земной коры использованы слои литосферной модели Волго-Уралии, представленной в разделе 4.1.2. [Ognev, Ebbing, Haas, 2022]. Все используемые структурные данные представлены на Рисунке 20.

Рисунок 20. Структурные данные, использованные в геотермическом моделировании. (а) Мощность осадочного чехла по модели EUNAseis [Artemieva, Thybo, 2013], (б) мощность верхней земной коры и (в) мощность нижней земной коры по модели из настоящего исследования [Ognev, Ebbing, Haas, 2022], (г) глубина LAB [Artemieva, 2019]

Термические данные.

В настоящей работе было использовано несколько наборов термических данных. Во-первых, это температура на поверхности LAB, которая была принята равной 1400 °C согласно со значением в модели И.М. Артемьевой [Artemieva, 2019]. Данная температура использовалась в алгоритме MCMC, как значение, которое должно прогнозироваться строящимися геотермическими моделями. Во-вторых, была использована информация о среднегодовой температуре на земной поверхности исследуемого региона по базе данных WorldClim 2 [Fick, Hijmans, 2017], так как при расчете геотерм земной коры и литосферы всегда должна вводиться поправка за значение температуры на поверхности Земли. В-третьих, были использованы измеренные значения теплового потока на поверхности Земли по базе данных Thermoglobe, включающей в себя более 600 измерений для исследуемого региона [Jennmgs, Hasterok, Lucazeau, 2021]. В данном случае тепловой поток на земной поверхности является интегральной характеристикой земных теплопотерь в каждой точке измерений, являясь суммой мантийного теплового потока и радиогенного теплового потока земной коры. Температура и тепловой поток

на земной поверхности были перемасштабированы на сетку с разрешением ~50 км для согласованности со структурной моделью, построенной ранее. Карты температуры и теплового потока на поверхности Земли показаны на Рисунках 21а-б, гистограмма измеренных значений теплового потока приведена на Рисунке 21в. По Рисунку 21б видно, что места измерений теплового потока тяготеют к областям выявленных скоплений УВ в пределах Волго-Уральской НГП на востоке Волго-Уральского субкратона. Данная зависимость объясняется тем фактом, что статически большее число скважин пробурено в пределах нефтегазоносных территорий, позволяя проводить больше геотермических исследований.

Помимо упомянутых наборов термических данных, в настоящей работе также использовались результаты измерений температурных профилей в скважинах Поволжья, проведенные Н.Н. Христофоровой и коллегами. Здесь была использована карта температур на глубине 1000 м [Христофорова и др., 2004] и температурные профили, опубликованные для нескольких глубоких скважин региона [Христофорова, 2002; Христофорова и др., 2004]. Эти данные использовались на последнем этапе геотермического моделирования в процессе проверки полученных моделей.

Температура на поверхности(°С)

15 30 45 60 75

Тепловой поток на земной поверхности (мВт м2)

Рисунок 21. Термические данные, использованные в геотермическом моделировании. (а) Среднегодовая температура на земной поверхности исследуемого региона по базе данных WorldClim 2 [Fick, Hijmans, 2017], (б) тепловой поток на земной поверхности по базе данных Thermoglobe [Jennings, Hasterok, Lucazeau, 2021], (в) гистограмма, основанная на измерениях

теплового потока

3.2.2. Расчет распределения температуры в земной коре и верхней мантии в случае

однослойной земной коры

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

Ас , Qo

2КС кс

Литосферная мантия обладает низким радиогенным теплообразованием, по разным оценкам, составляющим 0.02-0.03 мкВт/м3 [Hasterok, Chapman, 2011; Rudnick, McDonough, O'Connell, 1998]. Во многих геотермических исследованиях РТО литосферной мантии пренебрегают и принимают равным нулю (например, [Afonso et al., 2019; Lösing, Ebbing, Szwillus, 2020]). В случае предположения об отсутствии ТГИ в верхней мантии, ее геотерма будет линейной и будет зависеть от глубины Мохо, теплового потока и температуры на поверхности Мохо:

Tz = TM+^^z-M), (13)

км

где Тг - это температура внутри литосферной мантии на глубине г, °С; Тм - температура на поверхности Мохо, °С;

Qм - тепловой поток на поверхности Мохо или мантийный тепловой поток, Вт/м2; км - теплопроводность литосферной мантии, Вт/м/К; М - глубина поверхности Мохо, м. Теперь, скомбинировав уравнения (12) и (13), получим:

Ю0 - 0.5АСМ\ QM TZ = T0+M С ) +jr(z- М). (14)

\ Кг / Км

Учитывая, что земная кора обладает постоянным радиогенным теплообразованием, тепловой поток на границе Мохоровичича может быть выражен через поверхностный тепловой

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

Ю0 - 0.5АСМ\ Q0 - АСМ TZ = T0 + M Г° с ) + С (Z - М). (15)

\ Кг ) Км

3.2.3. Расчет распределения температуры в земной коре и верхней мантии в случае

многослойной земной коры

Распределение температуры в многослойном полупространстве, состоящем из n слоев, было сформулировано Г.Н. Поллаком [Pollack, 1965]:

п

¿-1

Qo-^ArHt ¿=1

- —•-Н-2 2ki Hl

(16)

где Тп - это температура в основании п-го слоя, к[ и А[ - теплопроводность и РТО каждого из п слоев, НI - мощность каждого из п слоев.

Для земной коры, состоящей из трех слоев: (1) осадочного чехла, (2) верхней ЗК, и (3) нижней ЗК, уравнение (16) может быть переформулировано следующим образом:

(H-i Нп Но\ A-IH-IH? A-IH-IH?

К1 К2 К3' к2 к3 к3

2 2 2ч (17)

<А1Н12 _ А2Н22 _ А3Н32\

k^ к2 к

где к1, к2, к3, и А1, А2, А3 - это теплопроводности и РТО осадочного чехла, верхней и нижней ЗК, Нь Н2, Н3 - мощности этих слоев.

В зависимости от слоя земной коры, для которого производится расчет геотермы, уравнение (16) может быть изменено соответственно. Так, в осадочном чехле температура на глубине 2 с РТО А1 и теплопроводностью к1 может быть получена через:

9 Qo

—К-^ гС^

Температура в верхней ЗК с РТО А2 и теплопроводностью к2 рассчитывается через:

(Нл 2\ АлНл! (АлНл2 Агг2\

+ + + (19)

где 2 = 2 — - это глубина ниже подошвы осадочного чехла.

Также температура нижней ЗК с РТО А3 и теплопроводностью к3 может быть получена как:

/н1 н2 г \

^1+^2+1_\—А1н1н1—А1н1г — А2Н2г — о5

к\ к2 к3/ к2 к3 к3

2 2 (20) А1Н12 + А2Н22 + А^2\ у '

к1 к2 к3 )

где 2 = х — (Н1 + Н2) - это глубина ниже подошвы верхней ЗК.

Для трехслойной ЗК с постоянными РТО в пределах каждого слоя, мантийный тепловой поток может быть выражен через тепловой поток на земной поверхности следующим образом:

Qм = Qo— А1Н1 — А2Н2 — А3Н3. (21)

Комбинируя уравнения (13), (17) и (21), мы можем получить распределение температуры в литосферной мантии для случая многослойной земной коры:

+ ^А _ _ А1Н1Н3 _ А2Н2Н3 _ о ^

ккл ко ко) ко ко ко

(П1 По По\

ц+ц+ц)

(22)

/А1Н12 + А2Н22 + А3Н32\ + д0—А1Н1—А2Н2—А3Н3 ^ _ м)

к^ к2 к3 к

м

3.2.4. Байесовская инверсия для определения геотермической структуры Волго-Уралии

Уравнения (13) и (15) показывают, что для расчета геотерм земной коры и литосферной мантии в случае однослойной земной коры требуется информация о таких тепловых параметрах, как температура на земной поверхности Т0, поверхностный тепловой поток теплопроводность и РТО земной коры, кс и Ас, и теплопроводность верхней мантии км. В случае многослойной земной коры появляются дополнительные неизвестные параметры, такие как мощности, теплопроводности и РТО осадочного чехла, верхней и нижней ЗК (Нь Н2, Н3, к1, к2,

к3, и А1, А2, Аг) согласно уравнениям (17) - (22). В данной работе в пределах исследуемого региона известны вариации теплового потока и температуры на земной поверхности, глубины LAB, а также мощностей отдельных слоев земной коры (Таблица 7). Поэтому необходимо ответить на вопрос, как определить оставшиеся неизвестные термические параметры Волго-Уральской литосферы? Для этой цели был использован байесовский статистический подход [Box, Tiao, 1973; Kitanidis, 2011]. Байесовская инверсия использовалась для решения геофизических обратных задач с конца 1960-х годов [Cary, Chapman, 1988; Keilis-Borok, Yanovskaja, 1967; Mosegaard, Tarantola, 1995; Press, 1968]. Однако ее применимость была усилена современным прогрессом вычислительной техники, что привело к большому количеству геофизических исследований, направленных на решение широкого круга обратных задач с помощью байесовской статистики (например, [Afonso et al., 2016; Bodin et al., 2012; Gallagher, 2012; Minson et al., 2014]). В настоящей работе используется статистический подход байесовской инверсии MCMC, заимствованный и модифицированный по [Lösing, Ebbing, Szwillus, 2020]. Поскольку для территории Волго-Уралии известна глубина LAB и ее температура [Artemieva, 2019], метод MCMC может быть использован для прогнозирования температуры на границе LAB через уравнения (15) и (22) и ее сравнения с известным значением. В результате такой операции будут найдены значения всех неизвестных тепловых параметров для однослойной и многослойной ЗК.

Таким образом, все неизвестные моделируемые тепловые параметры для каждой колонны земной коры и верхней мантии исследуемой области помещаются в вектор т = (км, кс, Ас) для однослойной земной коры и т = (км, klt к2, к3, Alt А2, А-$) для многослойной земной коры. Далее применяется оператор F(m), который предсказывает температуру на глубине LAB, используя байесовский подход MCMC. После этого известные значения температуры на границе LAB помещаются в вектор Т = (TLAB) и производится сравнение прогнозных и измеренных температур LAB в каждой ячейке модели. В конечном итоге целью инверсии является получение такого набора тепловых параметров т, для которого F(m) « Т.

Теорема Байеса позволяет рассчитать апостериорное распределение вероятности тепловых параметров модели с учетом данных по измерениям температур:

L(Tlm) Р(т)

Р(т\Т) = ' , (23)

где Р(т) - это априорное распределение, содержащее априорное знание о значениях моделируемых параметров, L(T\m) - функция правдоподобия, описывающая, насколько вероятны наблюдаемые температуры Т для набора параметров модели т, Р(Т) - это

коэффициент нормализации, не зависящий от модели [Malinverno, Briggs, 2004; Moorkamp et al., 2020].

P(T) также известен как предельная вероятность и остается постоянной для всех возможных моделей. По этой причине Р(Т) сокращается при оценке относительных вероятностей различных моделей [Statistical Inversion of Electromagnetic Logging Data, 2021]. Таким образом, вероятность конкретной геотермической модели P(mlT) будет пропорциональна априорной вероятности Р(т) и функции правдоподобия L(T\m):

Р(т\Т) a L(T\m) P(m). (24)

Здесь для априорной вероятности используется равномерное распределение.

1

lb <m<ub

P(m) = \ub-lb , (25)

0 иначе

где ub и lb - верхняя и нижняя границы диапазонов значений, которые потенциально могут принимать моделируемые тепловые параметры.

Функция правдоподобия L(T\m) может быть получена из нормального распределения (распределения Гаусса):

1 1 (F(m)-T)2

L(T\m)=-==e2 , (26)

где ат - стандартное отклонение, описывающее, насколько хорошо должны совпадать предсказанные и наблюдаемые температуры LAB. Здесь наблюдаемая температура Т принята в качестве среднего, а в качестве стандартного отклонения используется значение 1 °C для того, чтобы обеспечить удовлетворительное соответствие прогнозируемых и измеренных температур.

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

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

Очевидно, что описанная выше методика легко может привести к нереалистичным значениям тепловых параметров после нескольких изменений исходной геотермической модели. В связи с этим в данной работе используется алгоритм Метрополиса-Гастингса для того, чтобы решить, какие из предлагаемых комбинаций тепловых параметров должны быть приняты или отклонены. Упрощенная блок-схема алгоритма показана на Рисунке 22. Технически, после каждого шага изменения модели вычисляются апостериорные вероятности геотермической модели до P(m0idlT) и после ее изменения P(mnewlT). Если P(mnewlT) > P(m0idlT), тогда отношение P(^newlT)/ P(m0i¿lT) будет больше 1 и новая комбинация тепловых параметров будет принята, сохранена и установлена в качестве отправной точки перед следующим изменением модели. Даже в случае, если вероятность модели после ее изменения окажется меньше вероятности модели до изменения, шанс принять mnew все равно остается. Новая предложенная комбинация тепловых параметров mnew будет принята, если отношение вероятностей окажется выше, чем случайно полученное число из равномерного распределения между 0 и 1. Если последнее условие не выполняется, тогда шпет будет отклонена, и старая комбинация mold будет сохранена второй раз и снова использована в качестве отправной точки при следующей итерации. Как видно из уравнений (24) и (25), апостериорная вероятность комбинации тепловых параметров будет равна 0, если значение одного из ее параметров находится вне диапазона априорной модели. При этом, чем дальше прогнозируемые температуры F(m) находятся от фактических Т, тем меньше будет значение функции правдоподобия и, следовательно, меньше окажется значение апостериорной вероятности оцениваемой комбинации тепловых параметров. Такой подход гарантирует, что будут приняты только комбинации с реалистичными значениями тепловых параметров, и в целом цепь Маркова будет иметь тенденцию к большей сходимости предсказанных и измеренных температур LAB.

В результате при запуске описанного алгоритма N раз, будет получена цепь из N комбинаций тепловых параметров. Часто предполагается, что первые несколько членов цепи Маркова не являются репрезентативными и поэтому должны быть отброшены (например, [Lösing, Ebbing, Szwillus, 2020; Roy, Romanowicz, 2017]). Количество отбрасываемых шагов называется периодом приработки. В настоящем исследовании реализован период приработки, равный 10 % от всех итераций, что для данной модели достаточно для перехода цепи Маркова в стационарное состояние.

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

Таблица 8 - Исходная (априорная) геотермическая модель, шаг изменения и возможные

диапазоны тепловых параметров

Тепловой параметр

Исходная Шаг модель изменения mstart_модели

Диапазон априорной модели

Однослойная земная кора

Теплопроводность верхней мантии км (Вт/м/К) Теплопроводность земной коры кс (Вт/м/К) РТО земной коры Ас (мкВт/м3)_

3.8 2.3 0.9

0.05 0.05 0.1

3.5-4.1 1.5-3.0 0.2-1.6

Многослойная земная кора

Теплопроводность верхней мантии км (Вт/м/К)

Теплопроводность осадочного чехла к1 (Вт/м/К)

Теплопроводность верхней ЗК к2 (Вт/м/К)

Теплопроводность нижней ЗК к3 (Вт/м/К)

РТО осадочного чехла А1 (мкВт/м3)

РТО верхней ЗК А2 (мкВт/м3)

РТО нижней ЗК А3 (мкВт/м3)_

3.8 2.5 2.5 2.0 1.5 1.5 0.5

0.05 0.05 0.05 0.05 0.1 0.1 0.1

3.5-4.1 0.5-3.5 1.5-3.5 1.5-2.5 0.2-2.5 0.2-2.5 0.1-1.0

Примечание: составлено автором.

В настоящем исследовании автор предполагает, что теплопроводность Волго-Уральской коры и литосферы должна следовать глобальным трендам распределения теплопроводности. Поэтому для теплопроводности земной коры был взят диапазон кс от 1.5 до 3.0 Вт/м/К при начальном значении 2.3 Вт/м/К [Artemieva, Mooney, 2001; Furlong, Chapman, 1987]. Литосферная мантия состоит преимущественно из перидотита и, как ожидается, будет иметь постоянную теплопроводность около 4.0 Вт/м/К [Artemieva, 2011; Artemieva, Mooney, 2001]. Таким образом, для теплопроводности литосферной мантии км был взят небольшой диапазон 3.5 - 4.1 Вт/м/К при начальном значении 3.8 Вт/м/К. При задании диапазона теплопроводностей для трехслойной

модели земной коры учитывалось, что теплопроводность в осадочных породах может существенно различаться в зависимости от их типа и пористости [Artemieva, 2011; Cermak, Rybach, 1982], что в частности подтверждается для Волго-Уральского региона исследованиями пород нефтяных месторождений Татарстана [Липаев, Гуревич, Липаев, 2001]. В связи с этим для осадочного чехла был установлен широкий априорный диапазон теплопроводности осадочного чехла в 0.5 - 3.5 Вт/м/К при 2.5 Вт/м/К в качестве начального значения. Горные породы высоких стадий метаморфизма слагают кристаллическую верхнюю кору Волго-Уралии [Постников, 2002; Bogdanova et al., 2021]. Для пород такого типа характерны изменения теплопроводности в зависимости от температуры и содержания в них полевого шпата и кварца. В то же время теплопроводность нижней коры почти постоянна (~2.0 Вт/м/К) [Artemieva, 2011; Artemieva, Mooney, 2001]. Таким образом, для теплопроводности верхней земной коры был выбран более широкий диапазон в 1.5 - 3.5 Вт/м/К при исходном значении 2.5 Вт/м/К, в то время как диапазон теплопроводности нижней земной коры составил 1.5 - 2.5 Вт/м/К с 2.0 Вт/м/К в качестве начального значения.

Как указывалось в разделе 2.3., существует недостаток информации о радиогенном теплообразовании ЗК Волго-Уралии. В связи с этим автор полагался по большей части на данные о глобальных особенностях распределения РТО. К. Жаупар и Ж.-К. Марешаль показали, что объемное среднее континентальное РТО земной коры имеет тенденцию к снижению с возрастом земной коры и в среднем должно находиться в диапазоне 0.8 - 1.0 мкВт/м3 [Jaupart, Mareschal, 2014] (Таблица 2). Тем не менее, объемное среднее значение РТО земной коры в зависимости от ее петрологического состава может иметь большие латеральные вариации даже в пределах одной геологической провинции [Jaupart, Mareschal, 2014; Jaupart, Mareschal, Iarotsky, 2016]. Принимая во внимание, что метаморфические породы нижней коры, относящиеся к гранулитовой фации, как сообщается, имеют РТО всего порядка 0.2 - 0.5 мкВт/м3 [Hasterok, Chapman, 2011; Mareschal, Jaupart, 2013] и что среднее объемное РТО верхней коры, оцененное с помощью различных глобальных геохимических моделей, должно быть порядка 1.6 мкВт/м3 [Jaupart, Mareschal, Iarotsky, 2016], автор настоящей работы предполагает, что среднее РТО земной коры Волго-Уралии не может выйти за эти значения. Поэтому для значений РТО земной коры был использован диапазон 0.2 - 1.6 мкВт/м3 в случае однослойной модели земной коры. Здесь начальное значение РТО было взято в середине этого диапазона и равнялось 0.9 мкВт/м3, что соответствует глобальному среднему значению континентального РТО по [Jaupart, Mareschal, 2014]. Несмотря на то, что изменение РТО земной коры с глубиной не является монотонным, в целом, считается, что верхняя часть земной коры обогащена ТГИ по сравненною с ее нижней частью (например, [Artemieva, Mooney, 2001; Mareschal, Jaupart, 2013; Roy, Blackwell, Birch, 1968]). Поэтому при определении диапазонов и априорных значений РТО для трехслойной

модели земной коры были установлены разные диапазоны и априорные значения теплообразования для осадочного чехла, верхней и нижней земной коры. Полученные ранее данные свидетельствуют о том, что РТО осадочных горных пород в разных районах Волго-Уралии может варьировать от 1.0 до 2.3 мкВт/м3 [Любимова, Любошиц, Парфенюк, 1983]. В настоящей модели автор допускает несколько более высокий диапазон РТО осадочных пород: 0.2 - 2.5 мкВт/м3. Поскольку экспериментальные данные о РТО Волго-Уральской кристаллической коры отсутствуют, автором предполагалось, что более кислая и радиоактивная верхняя кора может иметь РТО, изменяющиеся в тех же пределах, что и РТО осадочных пород, а нижнекоровое РТО должно изменяться в пределах меньшего диапазона с целью учета основного состава пород нижней коры. Для нижнекорового РТО был принят диапазон 0.1 - 1.0 мкВт/м3. Здесь в качестве исходного значения было установлено 1.5 мкВт/м3 для верхней коры и осадочного чехла и 0.5 мкВт/м3 для нижней коры. Эти значения взяты согласно петрологическим оценкам РТО в центральных районах Волго-Уралии, произведенными Н.С. Богаником [Боганик, 1975].

Рисунок 22. Упрощенная схема алгоритма Метрополиса-Гастингса, принятая в настоящем исследовании для моделирования вариации теплофизических параметров. На схеме приведены теплофизические параметры для случая однослойной модели земной коры

ГЛАВА 4. РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЯ СТРОЕНИЯ ЗЕМНОЙ КОРЫ И ТЕПЛОВОГО РЕЖИМА ВОЛГО-УРАЛИИ И ИХ СВЯЗИ С НЕФТЕГАЗОНОСНОСТЬЮ

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

4.1. Результаты гравитационного моделирования

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

4.1.1. Результаты инверсии поля градиента силы тяжести

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

Опорная глубина, при которой поверхность Мохо, полученная инверсией поля градиента силы тяжести, лучше всего согласовывалась с сейсмическими данными, равнялась 45 км. Такая относительно глубокая оценка получена благодаря тому, что на большой части сейсмических профилей TATSEIS-2003 и URSEIS-95 глубина границы Мохо составляла более 50 км. Контраст плотности между земной корой и верхней мантией архейской кратонной коры и уральского орогена был принят за 550 кг/м3, тогда как для палеопротерозойских рифтовых систем и Прикаспийского бассейна он был принят за 500 и 350 кг/м3 соответственно (Рисунок 23а). Эти значения близки к предыдущим выводам М. Эшаха и др. [Eshagh et al., 2016], использовавшим гравитационные градиенты GOCE и определившим, что контраст плотности земной коры и мантии на территории Евразии должен находиться в диапазоне 400-600 кг/м3. В то же время, другие сейсмические исследования предполагают несколько меньший контраст плотности -около 300-400 кг/м3 для тектонических условий, аналогичных условиям моделируемого региона [Bassin, Laske, Masters, 2000; Rabbel, Kaban, Tesauro, 2013]. Это несоответствие могло возникнуть ввиду того, что гравитационные методы усредняют плотности земной коры и подкоровой мантии и выражают их разницу в интегрированном сигнале. В то же время методы, основанные на

сейсморазведке, восстанавливают плотности конкретных слоев коры и литосферы и дают более точную оценку контраста плотностей между нижней корой и литосферной мантией. Плотностная модель, используемая при ООГМ с решением прямой задачи гравиметрии в настоящем исследовании, обладает контрастом плотности между земной корой и верхней мантией около 330 кг/м3 (Таблица 6), что ближе к значениям, полученным по сейсмическим оценкам.

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

Рисунок 23. (а) Контраст плотности, определенный с использованием алгоритма Хааса и др. [Haas, Ebbing, Szwillus, 2020] и (б) глубина Мохо, полученная с помощью инверсии поля градиента силы тяжести. Опорная глубина Мохо равна 45 км, контраст плотности земной коры и верхней мантии в 550 кг/м3 отнесен к архейским кратонам и уральскому орогену, в 500 кг/м3 -к палеопротерозойским рифтам, в 350 кг/м3 - к Прикаспийскому бассейну

4.1.2. Результаты решения прямой задачи гравиметрии

Конечным результатом гравитационного моделирования с решением прямой задачи гравиметрии является трехмерная модель строения земной коры Волго-Уральского региона, построенная в IGMAS+. Она включает обновленную модель Мохо вместе с основными слоями земной коры и верхней мантии. Построенная в IGMAS+ модель имеет стандартное отклонение

измеренной и рассчитанной силы тяжести, равное 8.0 мГал, при этом коэффициент корреляции между измеренной и рассчитанной силой тяжести составил 0.91 (Рисунок 24а-в). Для вертикального градиента силы тяжести стандартное отклонение равно 0.13 этвеш, а коэффициент корреляции равен 0.81 (Рисунок 25а-в). Это можно считать приемлемой сходимостью для региональных исследований земной коры (например, [Sobh et al., 2019]). Как видно по Рисунку 25в, в пределах изучаемого региона наблюдается остаточное длинноволновое поле градиента силы тяжести. Потенциально это может свидетельствовать о латеральной плотностной неоднородности, связанной с изменением состава коры и литосферной мантии Волго-Уралии, что не учитывалось в настоящей модели. Степень подобных вариаций плотности, а также возможные причины таких вариаций рассматриваются на территории ВЕК И.М. Артемьевой [Artemieva, 2003] и более детально проанализированы в [Ognev, Utemov, Nurgaliev, 2021]. Общий вид 3D-модели IGMAS+ с расположением вертикальных разрезов представлен на Рисунке 26.

Измеренная силы тяжести (мГал) Рассчитанная сила тяжести (мГал) Разница сил тяжести [измеренная - рассчитанная] (мГал)

Рисунок 24. Сравнение измеренного и рассчитанного поля силы тяжести. (а) Измеренная аномалия силы тяжести в редукции Буге по модели XGM2019e ^^ег1е et а1., 2019], (б) Рассчитанная аномалия силы тяжести в редукции Буге по 3D модели, построенной в IGMAS+, (в) Разница между измеренной и рассчитанной аномалией силы тяжести

-0.90-0.75-0.60-0.45-0.30-0.15 0.00 0.15 0.30 0.45 -0.90-0.75-0.60-0.45-0.30-0.15 0.00 0.15 0.30 0.45 -0.60 -0.45 -0.30 -0.15 0.00 0.15 0.30 0.45

Измеренный градиент силы тяжести (этвеш) Рассчитанный градиент силы тяжести (этвеш) Разница градиентов [измеренный - рассчитанный] (этвеш)

Рисунок 25. Сравнение измеренного и рассчитанного вертикальных градиентов силы тяжести. (а) Измеренный Вертикальный градиент силы тяжести спутника GOCE с поправкой за рельеф на высоте 225 км [Bouman et al., 2016], (б) Рассчитанный по 3D модели, построенной в IGMAS+ вертикальный градиент силы тяжести с поправкой за рельеф, (в) Разница между измеренным и

рассчитанным градиентами силы тяжести

Рисунок 26. Трехмерная модель литосферы Волго-Уралии, построенная в программе IGMAS+. Она состоит из 67 вертикальных разрезов, что дает пространственное разрешение порядка 50 км. Модель включает 5 основных слоев: осадочный чехол, верхнюю земную кору, нижнюю земную кору, литосферную мантию, астеносферу и дополнительный 6-й слой андерплейтинга. Аномалия силы тяжести в редукции Буге, рассчитанная по данной модели, показана сверху

модели

Прежде чем был достигнут вышеобозначенный уровень сходимости полей силы тяжести и градиента силы тяжести, на начальном этапе моделирования в IGMAS+ была выявлена область со значительным несоответствием измеренной и рассчитанной силы тяжести. Эта область была проинтерпретирована и смоделирована как регион с уплотненным слоем нижней коры, образовавшимся в результате андерплейтинга (см. раздел 3.1.3.). Такое расхождение измеренного и рассчитанного гравитационных полей возникло после подстройки инвертированной глубины Мохо к данным сейсмических измерений в северо-западной части сейсмического профиля TATSEIS-2003. Разница между сейсмическими оценками глубины Мохо и моделью Мохо, полученной в результате гравитационной инверсии, показана на Рисунке 27а. На Рисунке 27б показана разница между глубиной Мохо, полученной при решении прямой задачи гравиметрии, и глубиной Мохо, полученной в результате инверсии. Из Рисунка 27 видно, что граница Мохо находится значительно глубже в центре Волго-Уралии в области интерпретируемого андерплейтинга. Эта особенность изначально не видна на карте глубины Мохо, полученной инверсией поля градиента силы тяжести.

Рисунок 27. (а) Сравнение глубины Мохо, полученной в результате гравитационной инверсии с данными сейсмических измерений, где основные расхождения сосредоточены в сейсмическом профиле TATSEIS 2003 в центре модели, (б) Разница между глубиной поверхности Мохо, полученной в IGMAS+ (Рисунок 29) с ее оценкой по результатам гравитационной инверсии

Гипотеза андерплейтинга в центре Волго-Уралии уже обсуждалась различными исследователями. В частности, наличие андерплейтинга в этой области было предложено Г. Тюбо и И.М. Артемьевой [^уЬо, Artemieva, 2013] и рассматривалось в другой литературе [Минц и др., 2010; Bogdanova е! а1., 2010; Bogdanova, Gorbatschev, Garetsky, 2016]. Обнаруженное тело андерплейтинга находится на севере Токмовского мегаблока под Окским блоком (Рисунок 11). Это тело определяется на сейсмическом профиле TATSEIS-2003 как акустически прозрачная область между 650 и 900 км профиля [Trofimov, 2006] (Рисунок 9). Минц и др. интерпретировали эту особенность сейсмической картины как область однородных основных пород, частично метаморфизованных в высокоплотные гранулиты или эклогиты в основании так называемой «ветлужской» синформы [Минц и др., 2010]. Изостатические расчеты по уравнению (11) также показывают тело высокой плотности со средней толщиной около 10 км, которое четко очерчено зоной изостатического дисбаланса в центре Волго-Уралии (Рисунок 28). Другие районы с большим дефицитом массы расположены на юго-востоке исследуемого региона и связаны с Прикаспийской впадиной и Южно-Уральским орогеном. Однако они не соответствуют сколько-нибудь значительной разности измеренного и рассчитанного гравитационного сигналов и обусловлены большим отклонением мощности осадочных пород и земной коры от средних значений на этой территории, что приводит к более высоким значениям дисбаланса массы.

Рисунок 28. Мощность высокоплотного нижнего слоя земной коры согласно изостатическим

расчетам

Глубина границы Мохо построенной в IGMAS+ модели хорошо согласуется с данными сейсмических измерений (Рисунок 29). Средняя разница сейсмической и смоделированной глубины Мохо составляет 0.75 км, стандартное отклонение - 3.31км. Это можно считать удовлетворительным результатом, так как считается, что сейсмические оценки Мохо имеют погрешность не менее 2 км [Ebbing et al., 2012]. Таким образом, построенная модель литосферы Волго-Уралии соответствует имеющимся сейсмическим и гравиметрическим данным c достаточным уровнем сходимости измеренного и рассчитанного гравитационных полей и высокой согласованностью с сейсмическими оценками глубины Мохо.

Рисунок 29. Модель Мохо Волго-Уральского субкратона, полученная в результате гравитационной инверсии с латерально переменными контрастами плотности и последующего ООГМ с решением прямой задачи гравиметрии в IGMAS+ и использованием данных

сейсмических измерений

Большая часть расхождений между сейсмическими данными и моделью Мохо, разработанной в IGMAS+, связана с сейсмическим профилем TATSEIS-2003. Сейсмические глубины Мохо по профилям TATSEIS-2003 и URSEIS-95 в сравнении с моделью показаны на рисунке 30. Как видно, в пределах профиля TATSEIS сейсмическая граница Мохо имеет

несколько крутых впадин [Artemieva, ^уЬо, 2013; Tшfimov, 2006], которые не отражаются в спутниковом поле силы тяжести. Эта особенность привела автора к компромиссному решению: итоговая глубина Мохо построена с учетом трендов, заданных сейсмическими профилями, однако резкие перепады при этом были сглажены для лучшей сходимости измеренных и рассчитанных полей силы тяжести и вертикального градиента силы тяжести. В результате видно, что Мохо опускается на наиболее значительные глубины в центральной кратонной области и поднимается к периферийной зоне Волго-Уралии, где достигает минимальной глубины около 40 км под мощным осадочным разрезом Серноводско-Абдулинского авлакогена. Затем граница раздела земной коры и мантии резко опускается до глубин ниже 50 км после достижения уральской складчатой системы.

Сила тяжести и градиент силы тяжести по сейсмическим профилям ТАТЭЕ13 и иРЭЕ^

-0.4

Измеренное поле силы тяжести в редукции Буге [модель Х6М_2019е] (мГал) Рассчитанное поле силы тяжести в редукции Буге [модель ЮМАЭ+] (мГал)

Измеренный градиент силы тяжести с поправкой за рельеф на высоте 225 км (этвеш)

Рассчитанный градиент силы тяжести с поправкой за рельеф на высоте 225 км (этвеш) [модель ЮМАЭ+]

600 800 Расстояние (км)

Волго-Уральский субкратон

Серноводско-Абдулинские авлакогены

Уралиды

Ыт Рельеф

Верхняя земная кора, р=2750 кг м-Нижняя земная кора, р=2900 кг

Осадочный чехол

ШБЕ^-Эб

Андерплейтинг, р=3100 кг м 3 /

Верхняя мантия, р=3230 кг м-3

Глубина Мохо по профилям ТАТЗЕ18 и иР5Е15 Плотность разреза (кг м 3)

2058 32

Рисунок 30. (а) Измеренные и рассчитанные аномалии силы тяжести в редукции Буге и аномалии вертикального градиента силы тяжести с поправкой за рельеф. (б) Разрез модели IGMAS+ по глубинным профилям TATSEIS-2003 и URSEIS-95 - см. Рисунок 16 для привязки к карте. Геологический разрез увеличен по вертикали в 10 раз, а топография - в 100 раз по

отношению к горизонтальному масштабу

Аналогичную тенденцию к изменению геометрии поверхности границы Мохоровичича можно наблюдать и при рассмотрении разреза, пересекающего Волго-Уралию с севера на юг (Рисунок 31). Наиболее мощная земная кора наблюдается в северной кратонной области, удаленной от палеорифтов. Тогда как под рифтовыми структурами и в примыкающей к ним

кратонной области кора утоняется. Особенно маломощная ЗК наблюдается в южной части разреза под Прикаспийским осадочным бассейном.

- Измеренное поле силы тяжести в редукции Буге [модель Х6М_2019е] (мГал)

---- Рассчитанное поле силы тяжести в редукции Буге [модель ЮМАЗ+] (мГал)

- Измеренный градиент силы тяжести с поправкой за рельеф на высоте 225 км (этвеш)

---- Рассчитанный градиент силы тяжести с поправкой за рельеф на высоте 225 км (этвеш) [модель ЮМА5+]

200 400 600 800 1000 1200

Расстояние (км)

Волго-Ураль-ский субкратон

Камско-Бельские рифты

Прикаспийский бассейн

Верхняя земная кора, р=2750 кг м-

Осадочный чехол

Нижняя земная кора, р=2900 кг М"

• •• о о

Верхняя мантия, р=3230 кг м~:

^^убинаМохопосейс^рофилям Плотность разреза (кг м 3)

2058

Рисунок 31. (а) Измеренные и рассчитанные аномалии силы тяжести в редукции Буге и аномалии вертикального градиента силы тяжести с поправкой за рельеф. (б) Разрез модели IGMAS+ в направлении с севера на юг - см. Рисунок 16 для привязки к карте. Геологический разрез увеличен по вертикали в 10 раз, а топография - в 100 раз по отношению к горизонтальному масштабу. Здесь ВУСК - это Волго-Уральский субкратон, САА -

Серноводско-Абдулинский авлакоген

Модель Волго-Уральской литосферы, построенная в IGMAS+, показала изменение мощности земной коры в диапазоне от 32 до более чем 55 км в некоторых районах. Наиболее маломощная кора мощностью менее 40 км сформировалась на Прикаспийской впадине и Предуральском прогибе, что соответствует наиболее мощным осадочным толщам. Относительно маломощная кора обнаружена также вдоль Центрально-Русской рифтовой системы, а также в северной и южной частях Пачелмского рифта. В осевых частях обоих рифтов мощность земной коры сокращается до 40-42 км, а на окружающей территории кора снова набирает мощность до 44-46 км. Мощная земная кора расположена также под Уральскими горами и в центре Волго-Уральского субкратона. В этих регионах мощность земной коры превышает 50 км. В целом, созданная модель показывает, что в пределах ВЕК архейские кратонные блоки связаны с утолщением земной коры, а палеопротерозойские рифты - с уменьшением ее мощности.

4.1.3. Сравнение полученной модели с другими региональными моделями Мохо

Модель Мохо, построенная в IGMAS+, была сравнена с существующими глобальными и региональными моделями земной коры, охватывающими изучаемый регион. Для сравнения были выбраны глобальная модель CRUST 1.0, построенная преимущественно по сейсмическим данным с использованием гравитационных данных в областях с отсутствующими сейсмическими измерениями [Laske et al., 2013], глобальная модель GEMMA, построенная в результате инверсии гравитационных данных спутника GOCE [Reguzzoni, Sampietro, 2015], и региональная сейсмическая модель EUNAseis [Artemieva, Thybo, 2013]. Разница между глубиной Мохо в модели, представленной в настоящем исследовании, с моделями, упомянутыми выше, показана на Рисунке 32. На нем хорошо видно, что модель Мохо, построенная в настоящей работе, гораздо глубже, чем модель GEMMA, и имеет больше сходства по глубине с моделями CRUST 1.0 и EUNAseis (Рисунок 32а). Это объясняется тем фактом, что как модель CRUST 1.0 и особенно сейсмическая модель EUNAseis, так и настоящая модель ЗК построены с учетом доступных сейсмических данных о строении земной коры по сравнению моделью GEMMA, построенной исключительно по гравитационным данным.

Сравнивая модель Мохо настоящего исследования с моделью CRUST 1.0, видно, что построенная в IGMAS+ модель отличается более мощной корой в центре Волго-Уралии в области андерплейтинга (Рисунок 32б). Скорее всего, такое различие выявлено ввиду того, что при создании модели CRUST 1.0 не были использованы самые последние сейсмические исследования на территории Русской платформы, в том числе был упущен профиль TATSEIS. Видно, что модель EUNAseis, созданная на основе более обширной базы сейсмических данных для Русской платформы, имеет менее выраженное расхождение с моделью из настоящей работы в центре Волго-Уралии в зоне андерплейтинга (Рисунок 32в).

Последней очевидной особенностью является менее глубокое положение границы Мохо по построенной модели на юго-востоке Волго-Уралии в отличие от модели EUNAseis. Такая аномалия возникает из-за того, что сейсмический каталог Геологической службы США и сейсмическая база данных EUNAseis были построены независимо друг от друга и имеют определенные различия в сейсмических оценках Мохо в этом регионе. Модель, построенная в настоящем исследовании, больше соответствует сейсмическим оценкам глубины Мохо по данным каталога Геологической службы США на юго-востоке Волго-Уралии (Рисунок 29), но расходится с оценками Мохо по модели EUNAseis, которая показывает глубину Мохо, меньшую на 3-9 км в юго-восточной части Волго-Уралии и на Южном Урале.

Рисунок 32. Разница в глубинах Мохо между (а) моделью, построенной в настоящем исследовании с помощью IGMAS+, и моделями CRUST 1.0, (б) EUNAseis, (в) GEMMA. На верхней панели показаны карты разницы глубин Мохо в км между моделью из настоящего исследования и моделями, выбранными для сравнения. На нижней панели показаны

гистограммы разницы глубин Мохо в км

4.2. Результаты геотермического моделирования

В результате геотермического моделирования было получено распределение тепловых параметров Волго-Уралии для однослойной и многослойной земной коры. При моделировании использовалось 100 000 итераций для обеспечения стабильности решения и были отброшены

первые 10 % комбинаций тепловых параметров в качестве фазы приработки. Теоретически, чем больше количество итераций, тем стабильнее полученное решение. Однако проведенный синтетический тест показал, что нет значительного улучшения конечного результата после увеличения количества итераций свыше 100 000. Поэтому было принято решение сохранить данное количество итераций для построения геотермических моделей. Представленные окончательные геотермические модели были рассчитаны как средние арифметические последних 90 % комбинаций термических параметров, полученных по методу МСМС.

4.2.1. Геотермическая модель Волго-Уралии с однослойной земной корой

Распределение теплового потока на земной поверхности, мантийной и коровой теплопроводностей, мантийного теплового потока и РТО земной коры для однослойной модели земной коры показаны на Рисунке 33. Соотношения между парами различных моделируемых параметров приведены на Рисунке 34.

Как видно из Рисунка 33б-в, байесовская инверсия предлагает несколько более высокую изменчивость теплопроводности земной коры со средним значением 2.06 Вт/м/К и стандартным отклонением 0.27 Вт/м/К по сравнению с теплопроводностью верхней мантии, среднее значение которой составило 3.74 Вт/м/К и стандартное отклонение 0.14 Вт/м/К. Это ожидаемый результат, учитывая более узкий допущенный диапазон изменения теплопроводности верхней мантии. Видимые на Рисунке 33б вариации теплопроводности земной коры показывают ее уменьшение до ~1.5 - 2.0 Вт/м/К в области Уральских гор и восточной Фенноскандии. Если Уральский регион в целом характеризуется низкими значениями геотермических параметров [Голованова, 2005] и ожидаемо имеет относительно низкую теплопроводность, то аномалия в восточной части Фенноскандии не подтверждается более ранними исследованиями. Оба района также характеризуются наименьшими аномалиями поверхностного теплового потока (Рисунок 33в). Эта связь подтверждается заметной положительной корреляцией между поверхностным тепловым потоком и теплопроводностью земной коры (Рисунок 34).

Как видно по Рисунку 33г, радиогенное теплообразование земной коры Волго-Уралии сильно варьирует в пределах исследуемого региона. Среднее значение РТО для моделируемого региона довольно низкое (0.54 мкВт/м3) и находится в пределах типичных значений РТО для архейской коры, полученных К. Жаупаром и Ж.-К. Марешалем ^аираг!, Mareschal, 2014] (Рисунок 33г). РТО ЗК в пределах кратонной коры Волго-Уралии также имеет низкие значения со средним, равным 0.58 мкВт/м3. Однако в восточных и южных районах субкратона наблюдаются зоны с локальным повышением РТО, где оно может превышать среднее значение более чем в два раза. Эти повышенные уровни РТО могут возникать из-за спорадического присутствия кислых гранитоидов и особенно по причине обилия высокоглиноземистых сланцев

и гнейсов, которые обнаруживаются в архейских коровых комплексах Волго-Уралии [Богданова, 1986; Постников, 2002; Bibikova et al., 2009]. Последние могут показывать повышенное значение РТО, поскольку оно имеет тенденцию коррелировать с содержанием Al2O3 в метаосадочных породах [Hasterok, Gard, Webb, 2018].

Еще одним районом с повышенным РТО является Тимано-Печорский осадочный бассейн, протянувшийся вдоль неопротерозойского Тиманского орогена к северу от Волго-Уралии. Одной из причин такой аномалии может быть более молодая радиогенная протерозойская верхняя кора. Другой причиной может служить утонение литосферы, которое могло привести к увеличению мантийного теплового потока и, как следствие, теплового потока на земной поверхности [Artyushkov, Baer, 1986]. В этом случае аномалия с повышенным РТО будет менее выражена. Однако такое утонение не отражено в модели литосферы И.М. Артемьевой [Artemieva, 2019], используемой в настоящем исследовании, хотя оно и было предложено ранее в модели TC1 ([Artemieva, 2006]).

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

Уралиды, северо-западная часть Волго-Уралии и восточная часть Фенноскандии являются районами низкого теплообразования. Уралиды, согласно результатам моделирования, представляют собой район с наименьшим РТО земной коры, имеющим на всей своей территории значения в диапазоне 0.2-0.4 мкВт/м3. Этот результат согласуется с опубликованными низкими значениями измеренного РТО и тем фактом, что петрологически земная кора под Уральскими горами является более мафической, чем кора ВЕК [Brown et al., 2003; Kukkonen et al., 1997]. Для северо-западной части Волго-Уралии характерно отсутствие метаосадочных пород с предположительно повышенной радиоактивностью [Постников, 2002]. Низкие значения РТО в восточной Фенноскандии объясняются низким поверхностным тепловым потоком, составляющим ~20 мВт/м2, зарегистрированным в архейском карельском блоке Фенноскандии [Cermak, 1993]. Несмотря на то, что архейская кора восточной части Фенноскандии имеет относительно низкие значения РТО, как было показано в исследованиях на территории Финляндии [Jöeleht, Kukkonen, 1998; Kukkonen, Lahtinen, 2001], крайне низкие значения РТО в 0.2-0.3 мкВт/м3, полученные для северо-востока моделируемой области следует рассматривать с осторожностью. Для этой конкретной части изучаемого региона нет опубликованных данных по

РТО и на данной территории присутствует значительное количество кислых пород, которое теоретически должно повышать общую радиоактивность ЗК [Svetov, Slabunov, Kulikov, 2018].

В целом, как видно по Рисунку 34, вариации поверхностного теплового потока сильно коррелируют с распределением РТО, что делает РТО земной коры наиболее значимым фактором, контролирующим поверхностный тепловой поток Волго-Уралии. Такая взаимосвязь согласуется с общими представлениями о вариациях поверхностного теплового потока на континентальной коре (например, [Jaupart, Mareschal, Iarotsky, 2016; Roy, Blackwell, Birch, 1968]).

Мантийный тепловой поток Волго-Уралии характеризуется значительной пространственной изменчивостью (Рисунок 33д). Большая его часть определяется глубиной LAB, что отражается в довольно сильной отрицательной зависимости: чем меньше глубина LAB, тем больше мантийный тепловой поток (Рисунок 34). Такая взаимосвязь ожидается из уравнения (15) и наблюдается за счет значительного изменения мощности термальной литосферы, которое преобладает над влиянием вариаций РТО на мантийный тепловой поток (Рисунок 33д). В результате QM составляет в среднем около 20 мВт/м2 на территории Волго-Уральского субкратона и возрастает до ~25-30 мВт/м2 в Мезенских рифтах к востоку от Балтийского щита и даже до ~35-40 мВт/м2 в Прикаспийской впадине за счет маломощной термальной литосферы (Рисунок 33д). Эти результаты в целом соответствуют работе, проведенной ранее В. Чермаком [Cermák, 1982], который получил карту мантийного теплового потока Европы путем вычитания предполагаемой радиогенной коровой составляющей теплового потока из измерений поверхностного теплового потока. При этом карта мантийного теплового потока В. Чермака не показывает положительной аномалии на востоке Фенноскандии, полученной в настоящей модели. Такая аномалия является результатом неглубокой термальной литосферы в этом регионе, согласно модели И.М. Артемьевой [Artemieva, 2019], используемой в настоящей работе.

После получения однослойной геотермической модели Волго-Уралии, она была проверена с помощью двух доступных наборов данных по температуре земной коры: (1) карта измеренной температуры на глубине 1000 м и (2) профили температуры в нескольких глубоких скважинах [Христофорова, 2002; Христофорова и др., 2004]. В процессе проверки модели геотермы земной коры рассчитывались по уравнению (3).

Как видно по Рисунку 35а-б, тренды изменения температуры на глубине 1000 м в пределах Волго-Уральского субкратона, полученные в результате моделирования, в целом согласуются с таковыми на карте Н.Н. Христофоровой [Христофорова и др., 2004]. Расчетные и измеренные температуры на глубине 1000 м показывают отрицательные аномалии в области Уральских гор и в центре Волго-Уралии. Также в обоих случаях наблюдается положительная температурная аномалия на юге моделируемой области вблизи Прикаспийской впадины, а также локальные положительные аномалии на востоке Волго-Уралии близ системы Камско-Бельских рифтов. В то

же время на Рисунке 35в видно, что между измеренной и расчетной температурами существует расхождение до 15 °С. Чтобы устранить такое расхождение, необходима более высокая вариация теплопроводности земной коры или ее РТО. Однако при проверке температурных профилей в скважинах была обнаружена их относительно хорошая сходимость со средним значением -2.22 °С и среднеквадратичным отклонением 6.15 °С (Рисунок 39). В связи с высокой сходимостью скважинных геотерм и тем фактом, что на карте Н.Н. Христофоровой [Христофорова и др., 2004] множество изолиний температуры имеют предполагаемый характер, было принято решение оставить полученные вариации геотермических параметров в качестве итоговой модели.

Рисунок 33. Вариации тепловых параметров Волго-Уральского субкратона для однослойной модели земной коры: (а) тепловой поток на земной поверхности, (б) теплопроводность земной коры, (в) теплопроводность мантии, (г) РТО земной коры, (д) мантийный тепловой поток. Здесь мантийный тепловой поток рассчитывается по уравнению (4)

Рисунок 34. Коррелограммы структурных (глубины Мохо и LAB) и тепловых параметров (тепловой поток на земной поверхности Q0, РТО ЗК Ас, теплопроводность ЗК кс, мантийная теплопроводность км, мантийный тепловой поток QM) для однослойной модели земной коры. Нижний треугольник показывает двумерные графики ядерной оценки плотности, которые

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

диагонали

Рисунок 35. Разность температур для однослойной модели земной коры. (а) Температура, измеренная на глубине 1000 м в скважинах Поволжья, как показано в [Христофорова и др., 2004] вместе с расположением скважин с существующими температурными профилями, используемыми в настоящем исследовании. (б) Температура на глубине 1000 м, рассчитанная по результатам байесовской инверсии. (в) Разница между рассчитанной и измеренной

температурами на глубине 1000 м

4.2.2. Геотермическая модель Волго-Уралии с многослойной земной корой

В случае трехслойной модели закономерности распределения тепловых параметров практически идентичны таковым у однослойной геотермической модели (Рисунок 36). Например, вариации РТО различных слоев земной коры отражают друг друга и близки к вариациям РТО земной коры из однослойной модели (Рисунок 36д-ж). Единственной очевидной особенностью многослойной модели является то, что среднее значение РТО осадочного чехла и верхней земной коры значительно выше по сравнению с нижней земной корой. Как и в однослойной модели, теплопроводность в трехслойной модели остается практически постоянной для литосферной мантии и демонстрирует большие вариации в слоях земной коры (Рисунок 36а-г). Расчетный мантийный тепловой поток остался практически таким же, как и для однослойной модели (Рисунок 36з). Корреляционные тренды также сохраняются (Рисунок 37). Из Рисунка 37 видно, что в многослойном подходе, как и в однослойном, тепловой поток на земной поверхности положительно коррелирует с РТО, теплопроводностью и мантийным тепловым потоком. Здесь значение коэффициента корреляции выступает в качестве индикатора параметров, оказывающих наибольшее влияние на поверхностный тепловой поток. Таким образом, поверхностный тепловой поток в большей степени определяется теплообразованием слоев земной коры, при этом основную роль играет РТО верхней коры. Значительное влияние на поверхностный тепловой поток оказывает также теплопроводность слоев земной коры. Как теплопроводность

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

При проверке трехслойной модели для расчета температур в земной коре были использованы уравнения (18) и (19). Разница температур на глубине 1000 м между смоделированной и расчетной температурами при многослойном подходе имеет среднеквадратичное отклонение 10.09 °С по сравнению с 5.65 °С для однослойной модели. При этом расчетная карта температур в случае многослойной модели обладает гораздо более высоким уровнем шума, а также локальной высокотемпературной аномалией в районе центральных Уральских гор (Рисунок 38). Проверка многослойной модели путем рассмотрения сходимости температурных профилей в скважинах также не показала каких-либо значительных улучшений по сравнению с однослойным подходом (Рисунок 39). Суммарное среднеквадратическое отклонение и среднее разности температурных профилей в скважинах при многослойном подходе составляют 6.8 и 2.8 °С соответственно, что несколько выше, чем у однослойной модели. Тем не менее, следует отметить, что многослойное решение может отражать изменение градиента температуры при переходе от осадочного чехла к верхнему слою земной коры, чего не может быть достигнуто с использованием единственного слоя земной коры.

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

Рисунок 36. Вариации тепловых параметров Волго-Уральского субкратона для многослойной модели земной коры: теплопроводности (а) осадочного чехла, (б) верхней ЗК, (в) нижней ЗК, (г) верхней мантии; РТО (д) осадочного чехла, (е) верхней ЗК, (ж) нижней ЗК; (з) мантийный тепловой поток. Здесь мантийный тепловой поток рассчитывается по уравнению (21)

Рисунок 37. Коррелограммы структурных (глубины Мохо и LAB) и тепловых параметров

(тепловой поток на земной поверхности и мантийный тепловой поток Q0, QM; теплопроводности и РТО осадочного чехла, верхней и нижней ЗК кг, к2, к3, и Аг, А2, А3; мантийная теплопроводность км) для многослойной модели земной коры

Рисунок 38. Разность температур для многослойной модели земной коры. (а) Температура, измеренная на глубине 1000 м в скважинах Поволжья, как показано в [Христофорова и др., 2004] вместе с расположением скважин с существующими температурными профилями, используемыми в настоящем исследовании. (б) Температура на глубине 1000 м, рассчитанная по результатам байесовской инверсии. (в) Разница между рассчитанной и измеренной

температурами на глубине 1000 м

Скважина 1 (Сергач) Температура (° С) 20 40 60 80 100

Скважина 2 (Кукморская 20010) Температура (° С) 20 40 60 80 100 120

\ - Измеренная температура Однослойная ЗК РМЗ=5.15 ..... Многослойная ЗК ВМЗ=2.78

4 л

Осадка

Верхня я кора

Скважина 3 (Орьебаш 82) Температура(° С) 20 40 60 80 100

\\ - Измеренная температура Однослойная ЗК ПМЭ=2.06 ----- Многослойная ЗК ИМЗ=4.65

\

\

\\

\

Скважина 4 (Минибаевская 20000)

Температура (° С) 0 20 40 60 80 100 120

Измеренная температура Однослойная ЗК ИМ8=5.96 Многослойная ЗК РМ5=2.43

Скважина 5 (Ямашинская 2092) Температура (° С) 0 20 40 60 80 100 120

Скважина 6 (Веденяпинская 217)

Температура (0 С) 0 20 40 60 80 100 120

\ - Измеренная температура Однослойная ЗК ЯМ8=7.02 ----- Многослойная ЗК ЯМЗ=2.64

уЛ

у\ \\

Осад^

Верхня я кора

\ - Измеренная температура Однослойная ЗК Р?МЗ=1.56 ----- Многослойная ЗК !={МЗ=0.98

1

Осадки __

Верхняя кора

Скважина 7 (Приволжская 13)

Температура (° С) 20 40 60 80 100

Скважина 8 (Тепловская 1)

Температура (° С) 20 40 60 80 100