Авторефераты по всем темам  >>  Авторефераты по земле

На правах рукописи

Тихоцкий Сергей Андреевич

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

25.00.10 - Геофизика, геофизические методы поисков полезных ископаемых

АВТОРЕФЕРАТ

диссертации на соискание ученой степени доктора физико-математических наук

Москва - 2011

Работа выполнена в Учреждении Российской академии наук Институте физики Земли им. О.Ю.Шмидта РАН.

Официальные оппоненты:

д.ф.-м.н., профессор Винник Лев Павлович Учреждение Российской академии наук Институт физики Земли им. О.Ю.Шмидта РАН д.ф.-м.н., профессор Глазнев Виктор Николаевич Воронежский государственный университет, Геологический факультет д.ф.-м.н., профессор Яновская Татьяна Борисовна Санкт-Петербургский государственный университет, Физический факультет

Ведущая организация:

Московский государственный университет имени М.В.Ломоносова, Геологический факультет

Защита состоится 29 сентября 2011 г. в 14 часов на заседании диссертационного совета Д 002.001.01 при Учреждении Российской академии наук Институте физики Земли им. О.Ю. Шмидта РАН, расположенном по адресу: 123995, Москва, Д-242, ГСП-5, ул. Грузинская Большая, д. 10, стр. 1. ИФЗ РАН.

С диссертацией можно ознакомиться в библиотеке Учреждения Российской академии наук Института физики Земли им. О.Ю.Шмидта РАН.

Автореферат разослан л 2011 г.

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

Ученый секретарь диссертационного совета, к.ф.-м.н.

Пилипенко О.В.

Общая характеристика работы

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

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

2. Появлением качественно новых видов экспериментальной информации, включая данные, получаемые космическими летательными аппаратами. В частности, спутниковая система GRACE (запущена в 2002 г.) предоставляет информацию об аномалиях поля силы тяжести 1-2 раза в месяц, что позволяет изучать не только сами аномалии, но и их временные вариации, отражающие динамику развития планеты. Спутниковые навигационные системы (ГЛОНАСС, GPS) позволяют с точностью в миллиметры наблюдать современные движения земной поверхности, а данные радарной спутниковой интерферометрии (InSAR) - определять смещения больших участков земной поверхности, в том числе - сопровождающие сильные землетрясения.

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

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

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

1. Наличие предположений о точности и природе измеряемых величин, основанных на приближениях и аппроксимациях, которые были вполне верными 1-десятилетия назад, но более не выполняются на современном уровне точности.

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

3. Упрощённые модельные представления, лежащие в основе используемых методов интерпретации.

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

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

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

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

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

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

2. Разработан и опробован метод восстановления гармонического аномального магнитного поля по результатам измерений аномалий модуля магнитного поля T.

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

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

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

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

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

8. Алгоритм локальной сейсмической томографии применён к интерпретации данных, полученных в очаговой зоне Рачинского землетрясения 1991 г. Построены скоростная модель очаговой зоны и распределение гипоцентров афтершоков. Дана сейсмотектоническая интерпретация полученных результатов.

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

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

2. Показана возможность, предложен метод и разработан алгоритм восстановления гармонического аномального магнитного поля по результатам измерений аномалий T.

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

4. Предложено использовать разреженное разложение по вэйвлет-функциям Хаара для построения адаптивной параметризации среды при решении обратной задачи лучевой сейсмической томографии.

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

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

7. В результате обработки данных, полученных в ходе проекта TOMOVES, получена более детальная (чем представленные в предшествующих исследованиях) скоростная модель строения среды под вулканом Везувий в интервале глубин 2-4 км. В модели выделяется низкоскоростная аномалия непосредственно под вулканической постройкой, которая может быть связана со слоями, нарушенным во время подъёма магматического материала и вулканических землетрясений.

8. Построена новая скоростная модель очаговой зоны Рачинского землетрясения 1991 г. Её анализ позволил уточнить представления о механизме очага землетрясения.

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

В частности:

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

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

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

4. Метод нелинейной минимизации функционала может применяться при изучении механизма землетрясений по данным о деформациях земной поверхности, получаемым методами радарной спутниковой интерферометрии (InSAR) и систем глобального позиционирования (ГЛОНАСС, GPS).

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

6. Алгоритм решения задачи локальной сейсмической томографии может применяться при изучении очаговых зон и механизмов сильных землетрясений, а также при изучении геодинамически активных регионов. Алгоритм также может найти применение при изучении техногенной сейсмичности, в частности - при мониторинге гидроразрыва пласта.

На защиту выносятся следующие основные положения:

1. Современный уровень точности данных гравиразведки и магниторазведки и требования к достоверности и детальности результатов интерпретации требуют пересмотра традиционно используемых приближённых трактовок наблюдаемых аномалий и создания новых математических методов их обработки и решения обратных задач. При высокоточных гравиметрических наблюдениях необходимо более точно, чем это предписывает действующая инструкция по гравиметрической разведке, подходить к вопросу выбора систем высот и учитывать смешанный характер аномалий в свободном воздухе, образованных с использованием нормальной либо ортометрической систем высот. При интерпретации аномалий модуля магнитного поля T необходимо принимать во внимание их негармонический характер. Гармоническое аномальное поле амплитудой до 15 000 нТл (т.е. до 1/3 средней амплитуды нормального поля) может быть восстановлено по наблюдённым значениям аномалий T, после чего возможно использование широко распространённых методов интерпретации, в своей основе опирающихся на представление о гармоничности изучаемого поля.

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

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

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

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

При изучении вулканических построек возможно выделение магмаподводящих структур. При изучении очаговых зон сильных землетрясений Ч уточнение их механизма. В частности, обнаружено, что под конусом вулкана Везувий имеются области повышенной скорости, связанные с магмоподводящими каналами, а глубже - область пониженной скорости, которая может быть связана с нарушенными (во время подъёма магматического материала и вулканических землетрясений) слоями. Также показано, что разрыв в очаге Рачинского землетрясения 1991 г. развивался вдоль сложной поверхности, причем основная часть энергии выделилась в ходе подвижки по субгоризонтальной плоскости, залегающей в основании Рача-Лехумского прогиба, а меньшая - в ходе подвижки по двум более круто падающим поверхностям, обрамляющим основной очаг с севера и юга.

Апробация работы. Основные результаты диссертации докладывались на следующих конференциях: V (2004 г.), VII (2008 г.) и VIII (2010 г.) международных конференциях УПроблемы геокосмосаФ (Problems of Geocosmos); Сессиях генеральной ассамблеи Европейского союза наук о Земле (European Geoscience Union) - в 20и 2009 годах; 30 (2003 г.) и 31 (2004 г.) сессиях международного семинара УВопросы теории и практики интерпретации гравитационных, магнитных и электрических полейФ им. Д.Г.Успенского; X геофизических чтениях им. В.В. Федынского в 2008 г.

Материалы диссертации опубликованы в 28 печатных работах, из них 15 статей в рецензируемых журналах, рекомендуемых ВАК РФ для публикации материалов докторских и кандидатских диссертаций.

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

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

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

Структура и объем диссертации. Диссертация состоит из введения, 3 глав, заключения и списка литературы. Общий объем диссертации 241 страница, из них 2страниц текста, включая 45 рисунков. Список литературы включает 245 наименований на 25 страницах.

Благодарности. Автор с глубокой благодарностью вспоминает своего учителя В.М.Гордина, в беседах и совместных трудах с которым родились многие идеи и подходы, воплощённые в данной работе. Автор искренне благодарен своим ближайшим коллегам - И.В.Фокину и Д.Ю.Шур, принимавшим участие в обсуждениях, проведении расчётов, подготовке данных и написании программного кода. Автор глубоко признателен У.Ахауеру (Франция) за многочисленные консультации, неоднократную возможность работать в университете г. Страсбурга и предоставленные данные по проекту TOMOVES. Автор благодарит В.О.Михайлова, Е.А.Киселёву и Е.И.Смольянинову за продолжительный период плодотворного сотрудничества и поддержку, а также - за обработку данных InSAR по очаговой зоне Алтайского землетрясения, результаты которой были использованы в настоящей работе. Автор благодарен С.С. Арефьеву за важные консультации по вопросам очаговой сейсмологии и локальной сейсмотомографии и предоставленные данные по Рачинскому землетрясению. Автор весьма признателен Ю.О. Кузьмину за плодотворные и стимулирующие обсуждения различных аспектов геофизики и за важные замечания по данной работе. Автор признателен С.С.Полякову за написание программного кода, реализующего предложенную в работе модификацию метода Монте-Карло.

Автор выражает глубокую благодарность директору ИФЗ РАН академику А.О. Глико за поддержку и предоставленные широкие возможности для ведения научной работы. Автор также искренне признателен заместителям директора ИФЗ РАН О.Н. Галаганову, В.Н. Конешову, А.В. Пономарёву, Е.А. Рогожину, помощнику учёного секретаря Т.Н. Филатовой, помощникам директора Г.Н. Михайловой и С.В. Ляпуновой за прекрасную, доброжелательную и творческую атмосферу в коллективе, понимание и помощь. Автор весьма благодарен за поддержку и внимание директору ГЦ РАН, чл.-корр. РАН А.Д. Гвишиани. Автор искренне признателен за поддержку своей научной работы академику Д.В.Рундквисту.

Автор также благодарит за многочисленные плодотворные обсуждения различных аспектов научной работы, важные советы и рекомендации своих учителей и коллег: С.М. Агаяна, П.С. Бабаянца, Ю.И. Блоха, Ш.Р. Богоутдинова, Б.Г. Букчина, А.А. Булычёва, М.Л. Владова, А.С. Долгаля, А.Д. Завьялова, М.К. Кабана, А.В.

Калинина, В.В. Калинина, Н.К. Капустян, А.С. Кобрунова, С.С. Красовского, И.Ю.

Кулакова, М.В. Минца, А.А. Никитина, В.А. Рашидова, Ю.Л. Ребецкого, А.В. Старовойтова, В.И. Старостенко, М.Ю. Токарева и многих других.

Автор благодарит свою супругу и всю свою семью за постоянную поддержку и понимание.

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

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

Последние десятилетия характеризовались качественным изменением методов и технологий измерения потенциальных полей, в результате чего стало возможно выполнять съёмки с весьма высокой (сотые доли мГал - для аномалий силы тяжести и первые нТл - для аномалий магнитного поля) точностью. Появление таких данных стимулировало интерес к выделению и интерпретации низкоамплитудных аномалий, в частности - при поисковых работах на нефть и газ [Бычков, 2010]. Вместе с тем большинство используемых методов обработки и интерпретации опираются на приближённые трактовки измеряемых величин, которые были справедливы в эпоху создания этих методов. Не изменились и нормативные и справочные документы. Цели исследований состояли в анализе адекватности традиционно используемых трактовок наблюдаемых аномалий потенциальных полей и в развитии методов и алгоритмов обработки и интерпретации, адекватных современному уровню точности и потребностям геофизической практики.

Часть первая главы 1 посвящена анализу влияния системы высот, используемой при обработке результатов гравиметрических наблюдений, на результаты дальнейшей интерпретации аномалий силы тяжести. Стимулом к проведению соответствующих исследований стали экспериментально обнаруживаемые в гравиметрической практике различия между высотами точек наблюдения, определяемыми посредством нивелировочных геодезических работ и при помощи спутниковых систем глобального позиционирования. Задачи поиска структурных и неструктурных нефтегазовых ловушек требуют выделения аномалий силы тяжести амплитудой в Рис. 1. Экспериментально обнаруженная зависимость между аномалией в свободном воздухе (1) и аномалией высоты (2). Оценка аномалии высоты получена как разница высоты пунктов по результатам нивелировки и данным GPS. Аномалия в свободном воздухе образована с использованием высот, полученных в результате нивелировки сотые-десятые мГал. В практике современных крупномасштабных исследований достигается среднеквадратическая погрешность в (0,02 - 0,04) мГал [Бычков, 2010].

Экспериментально обнаруживаемые различия высот по данным нивелировок и спутниковых систем достигают в равнинных районах 40-50 см на расстоянии в 50-60 км (рис. 1), что соответствует невязкам аномалий в свободном воздухе в 0,12Ц0,15 мГал, т.е. в 3-7 раз превышает погрешность современных съёмок.

Обнаруживаемые различия аномалий силы тяжести, образованных с использованием различных систем высот, с высокой точностью объясняются влиянием так называемого косвенного эффекта, возникающего из-за наличия аномалии высоты = H - Hn, где H - геодезическая высота (определяемая GPS), Hn - нормальная высота (определяемая нивелировкой). Указанные особенности образования аномалий силы тяжести хорошо известны в физической геодезии, они составляют основу всей теории определения формы поверхности Земли по данным геодезических и гравиметрических измерений. С другой стороны, при крупномасштабных поисковых геофизических исследованиях косвенным эффектом практически всегда пренебрегают. Эта практика закреплена, в частности, в основных нормативных и справочных изданиях по гравиразведке: инструкции по гравиметрической съемке, государственном стандарте, справочнике геофизика, которые вообще не вводят различия систем высот, молчаливо предполагая, что все высоты представляют собой превышения относительно эллипсоида, т.е. геодезические высоты H. Исторически это связано с отсутствием детальных карт аномалий высоты (геоида) на территорию суши, а методологически обосновывается представлением, что все аномалии геоида носят длинноволновый характер, соответственно косвенный эффект искажает лишь фоновые значения аномального поля, которые не принимаются во внимание при интерпретации. Как видно из рис. 1 и приведённых выше параметров точности современных съёмок, в настоящее время указанное приближение более не справедливо.

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

Важно, что смешанная аномалия силы тяжести не равна вертикальной производной аномального потенциала u, как это всегда предполагается при решении обратных задач гравиметрии в практике разведочной геофизики. В действительности, смешанная аномалия силы тяжести удовлетворяет (в приближении сферической Земли) уравнению Стокса: g = - (u/r + 2u/R) [Шимбирёв, 1975]. Таким образом, оказывается, что применение традиционного набора математических методов интерпретации аномалий силы тяжести к данным современных высокоточных гравиметрических съёмок некорректно, если аномалии вычислены с использованием высот, определённых в результате нивелировок.

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

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

На необходимость принимать во внимание тот факт, что аномалии модуля магнитного поля T не являются гармонической функцией, неоднократно указывал В.Н.Страхов, который, в частности, предлагал использовать при решении обратных задач приближение второго порядка. В работах Ю.И.Блоха с соавторами, воплощённых в широко используемом программном комплексе СИГМА-3D, используются нелинейные алгоритмы решения обратных задач, адекватно трактующие аномалию T. Таким образом, при решении обратной задачи методом подбора указанная особенность аномалий T может быть адекватно учтена.

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

Обозначим u - скалярный потенциал, соответствующий аномальному магнитному полю. За нулевое приближение гармонического компонента T0 = u t, где t = TN/|TN|, примем измеренное в эксперименте поле T. Найдется вспомогательное распределение магнитных масс, порождающее скалярный магнитный потенциал v(0), такой, что на совокупности точек наблюдения соответствующий гармонический компонент поля совпадает с T :

N v(0) t = T = u + TN - T, (1) Негармоническая часть поля T, соответствующего потенциалу v(0), будет равна ( (0) N N d(0) = T - T00) = v(0) + TN - T - v(0) t = v(0) + TN - T - T. (2) Тогда в качестве следующего (первого) приближения гармонического компонента естественно принять ( N T01) = T - d(0) = 2T - v(0) + TN + T. (3) ( Полученное приближение гармонического компонента T01)вновь аппроксимируем полем вспомогательного распределения источников, т.е. найдём распределение магнитных масс, порождающее скалярный магнитный потенциал v(1) такой, что на со( вокупности точек наблюдения v(1) t = T01). Повторяя операции (2-3) получим следующее (второе) приближение гармонического компонента.

Продолжая описанный процесс по индукции, получим следующее общее выражение для k-го приближения гармонического компонента:

k-1 k- ( N T0k) = T + k T + T - v(i) + TN = T + k u + TN - v(i) + TN, (4) i=0 i=где u - потенциал аномального поля Ta, а v(i) - потенциалы, порождаемые вспомогательным распределением масс и определяемые на каждом шаге итерационного процесса из условия:

( v(i) t = T0i). (5) В работе доказано следующее утверждение: если итерационный процесс, опре( деляемый соотношениями (4-5) сходится, то последовательность T0k) сходится к T0. Иначе говоря, оценка гармонического компонента, получаемого в ходе итерационного процесса (4-5), является состоятельной.

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

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

Наглядно демонстрируется необходимость предварительного восстановления гармонического аномального поля и возможность такого восстановления с использованием предложенного алгоритма (рис. 2).

Рис. 2. Применение алгоритма восстановления гармонического компонента поля T в задаче определения положения особых точек аномалообразующего объекта методом полного нормированного градиента. а - графики значений аномального поля вдоль профиля вкрест простирания модельного объекта: сплошная линия - восстановленное поле T0, прерывистая - поле T ; б - поле полного нормированного градиента в сечении XZ, рассчитанное в результате продолжения вниз исходного негармонического поля T (сплошная линия - контур аномалообразующего объекта); в - поле полного нормированного градиента в сечении XZ, рассчитанное в результате продолжения вниз исходного гармонического поля T0 ; г - поле полного нормированного градиента в сечении XZ, рассчитанное в результате продолжения вниз восстановленного гармонического компо нента T0.

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

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

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

Начиная с 2003 г. система GRACE каждые 30 дней генерирует новую модель геопотенциала. Основной целью проекта является изучение глобальных изменений климата Земли и циркуляции мирового океана. Согласно проекту, кумулятивная точность определения высот геоида по данным GRACE в диапазоне сферических гармоник до порядка l = 90 (что соответствует длине волны около 300 км) должна была составлять около 3,5 миллиметров. Однако фактическая точность 30-дневных моделей геопотенциала почти на два порядка ниже, чем было заложено в проекте.

Спутник GOCE после 20 месяцев работы на орбите привёл к созданию одной модели геопотенциала. На настоящий момент не известно, будут ли созданы модели, отвечающие последующим временным интервалам.

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

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

Конвергенция плит на этих участках возобновляется после накопления напряжения, достаточного для преодоления сил трения. Если запертый участок зоны субдукции можно аппроксимировать одной плоскостью, то для расчета смещений поверхности [Okada, 1985] необходимо задать девять параметров, характеризующих её положение и относительное смещение. Большую часть этих параметров можно определить априори: по данным об ориентировке погружающейся плиты или механизме очага произошедшего землетрясения. В случае землетрясения сейсмологические данные дают ограничения на размеры области разрыва и амплитуду сдвига на ней. Амплитуда смещения в единицу времени a, в случае запертой зоны субдукции, пропорциональна скорости субдукции, которая обычно известна приближённо.

В исследованиях, описываемых в настоящей работе, предполагалось, что гравитационный эффект геодинамических процессов совпадает с возмущениями силы тяжести, возникающими из-за деформации земной поверхности в ходе этих процессов. Данное приближение обосновывалось тем, что контрасты плотности на внутренних границах как минимум на порядок меньше контраста плотности на поверхности твердой Земли. Кроме того, известно, что при землетрясении на Аляске в 1964 г. с магнитудой Mw=9,2 наблюдённые изменения гравитационного поля практически полностью совпали с эффектом от изменения высоты пунктов наблюдений [Barnes, 1997]. Вместе с тем анализ временных вариаций поля силы тяжести по данным GRACE, сопровождавших катастрофическое землетрясение на Суматре (Mw=9,Ц 9,3) 2004 г., показал, что они не могут быть объяснены исключительно деформацией земной поверхности: значительную роль сыграла объёмное разуплотнение пород земной коры [Panet at al., 2007]. Тем не менее, амплитуда и пространственный масштаб наблюдённых вариаций хорошо совпали с оценкой, полученной на основе предположения о доминирующем влиянии деформации земной поверхности. Поэтому приводимые в настоящей работе результаты, касающиеся принципиальной возможности выделения косейсмического сигнала в данных GRACE, опубликованные до момента Суматранского землетрясения, остаются справедливыми и подтверждаются данными, полученными при этом землетрясении.

Рассмотрим следующую модель экспериментального материала. Пусть измеренные значения аномалий высот геоида - u (, , t) слагаются из постоянной во времени составляющей этих аномалий u (, ), вариации высот геоида геодинамического происхождения uG (, , t) и помехи u (, , t), т.е.

u (, , ti) = u (, ) + uG (, , ti) + u (, , ti), i = 1,..., M, (6) где: - коширота, - долгота, ti - моменты времени, в которые производятся наблюдения, M - число серий наблюдений. От равенства (6) можно перейти к аналогичной системе равенств для коэффициентов разложений геопотенциала по сферическим функциям.

Будем предполагать, что данные спутниковых измерений были предварительно исправлены за атмосферные и гидросферные эффекты, а также за эффект гляциоизостатического выравнивания [Wahr et al., 2001]. Тогда помеха u (, , ti) слагается из погрешности, обусловленной конечной точностью самих наблюдений, вариаций высот геоида, связанных с ошибками при введении поправок за атмосферу и гидросферу и влияния других, неизвестных нам источников временных вариаций.

С достаточной степенью общности геодинамический сигнал можно представить в виде произведения:

uG (, , t) = f (, ) g (t). (7) Это справедливо, в частности, для большинства тектонических процессов, характерное время развития которых велико по сравнению с длительностью спутниковых наблюдений. Для запертых участков зон субдукции можно положить g (t) = v (t - t0), где величина v равна модулю вектора скорости конвергенции плит в зоне субдукции. Для сейсмотектонических деформаций, сопровождающих крупные землетрясения зависимость от времени может быть описана функцией скачка: g (t) = a (t - t0).

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

Выполнены оценки вероятности обнаружения косейсмического сигнала, соответствующего землетрясениям с параметрами, аналогичными землетрясениям на Аляске 1964 г. (Mw=9,2) и на Хоккайдо 2003 г. (Mw=8,1). Показано, что такое крупное событие, как землетрясение 1964 года на Аляске, вполне может быть обнаружено по данным GRACE, при фактической точности этой системы. Для обнаружения землетрясения на о. Хоккайдо нынешней точности измерений GRACE недостаточно, но при выходе на проектную точность обнаружение сигнала вполне реально.

Для землетрясения в Чили 21-22 мая 1960 г. предложено три модели плоскости разрыва, параметры которых сильно отличаются [Plafker, Savage, 1970]. Это связано с тем, что геодезические данные характеризуют области, далекие от областей главных деформаций. Результаты расчётов, приводимые в диссертации, показывают, что при фактической точности наблюдений данные GRACE позволят различить две наиболее сильно отличающиеся модели, различение других пар моделей проблематично. Однако если точность наблюдений увеличить хотя бы на порядок, то различение всех моделей возможно с вероятностью ошибок не хуже первых процентов.

При оценке возможности обнаружения накопления деформаций на запертом участке зоны субдукции использовались параметры запертого участка зоны субдукции на Аляске из работы [Zweck et al., 2002]. Предполагалось, что оценивание производится по 60 сериям наблюдений для спутника GRACE, разделенных интервалами времени длиной в один месяц. В результате расчётов показано, что фактической точности системы GRACE оказывается недостаточно для уверенного обнаружения запертой зоны субдукции, но уменьшение погрешности на один порядок позволит уверенно выделить искомый сигнал.

Результаты оценок возможности обнаружения геодинамического сигнала, связанного с крупнейшими землетрясениями, приведённые в настоящей работе и опубликованные в 2004 г. [8], затем были подтверждены при катастрофическом землетрясении на Суматре в декабре 2004 г. В работе [Panet at al., 2007] показано, что косейсмический сигнал был не только обнаружен в данных GRACE, но и позволил прояснить механизмы деформации земной коры, сопутствующие сильнейшим землетрясениям.

Другим важнейшим источником данных о косейсмических деформациях является метод радарной спутниковой интерферометрии InSAR, который рассматривается во второй части главы 2.

Пусть для области землетрясения имеются интерферометрические данные о смещениях в направлении на спутник (WLOS (i, i). Пусть также имеются данные GPS по системе пунктов повторных наблюдений, характеризующие горизонтальные косейсмические смещения на север и восток {VN (i, i), VE (i, i)}. Аппроксимируем поверхность разрыва набором из M плоскостей и используем для описания поля смещений на дневной поверхности решение задачи для плоской дислокации в упругом полупространстве [Okada, 1985].

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

Поверхность разрыва аппроксимируется совокупностью M плоскостей (M задано). Поиск решения осуществляется путём минимизации функционала:

obs mod obs mod obs mod (p) = VE - VE (p) 2 + VN - VN (p) 2 + WLOS - WLOS (p) 2 + L2 L2 L 2 2 2 2 M (8) i 2 i 2Un 2Un 2Uss Uss i i , + + + x2 y2 x2 y2 i=mod mod где расчетные смещения в пунктах наблюдений GPS (VN, VE ) и в узлах сети, mod где заданы смещения по данным интерферометрии (WLOS ), зависят от вектора параметров p = (f, g). Здесь f - составной вектор неизвестных смещений; g - вектор параметров, определяющих геометрию поверхности разрыва (угловые координаты плоскостей, углы падения и простирания, глубины верхней и нижней кромок). Приведённое разбиение вектора параметров p обусловлено тем, что смещения поверхности линейно зависят от компонент вектора f и нелинейно - от компонент g.

Зависимость смещения поверхности от совокупности компонент вектора g носит весьма сложный характер. Функционал (8) может иметь несколько локальных минимумов, а результат его минимизации с использованием методов типа градиентного спуска будет в значительной степени определяться выбором модели начального приближения. Более полную картину поведения подобных функционалов даёт использование методов типа Монте-Карло, которые позволяют опробовать более широкий диапазон пространства параметров.

Для решения поставленной задачи была разработана следующая модификация метода Монте-Карло:

1. Для каждого компонента gj вектора g задаётся тип распределения wj (нормальное либо равномерное на ограниченном интервале, в зависимости от физического смысла параметра) и начальные значения математического ожидания aj и дисперсии D(0).

j 2. На основе заданных распределений случайным образом генерируются K пробных векторов k (т.е. наборов нелинейных параметров модели).

g 3. Для каждого пробного вектора k на основе линейной минимизации определяg ются вектора f, доставляющие минимум функционала (8) при заданном фик k).

сированном векторе g = k (обозначим эти найденные вектора f g 4. Отбираются Q < K пробных векторов q, соответствующих Q наименьшим g q, q значениям функционала (8) при q = f g.

p 5. Значения математических ожиданий компонент gj вектора g приравниваются последовательно к значениям соответствующих компонент отобранных пробных векторов 1,..., Q и в каждом случае случайным образом генерируется g g K/Q новых пробных векторов k. Тип распределений компонент остаётся неизg менным, а значения дисперсий определяются по формуле: D(1) = kD D(0), где j j kD 1 - параметр алгоритма.

6. Осуществляется переход к пункту 3 и процесс повторяется.

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

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

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

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

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

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

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

В обзоре приводится историческая справка и даётся подробный обзор физических и математических основ метода лучевой сейсмической томографии. Далее рассматриваются основные проблемы метода лучевой сейсмической томографии, на решение которых направлена настоящая работа. Как показано в обзорных работах Т.Б.Яновской (1997) и Г.Ноле (2008), одной из наиболее важных проблем является оценка разрешающей способности данных и построение адекватных этой разрешающей способности моделей среды.

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

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

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

Наиболее распространённым подходом к решению таких задач является использование регуляризации [Тихонов, Арсенин, 1986]. В качестве регуляризирующего функционала в задачах сейсмической томографии обычно используется та или иная мера гладкости искомого распределения скорости [Дитмар, 1993; Hobro et al., 2003; Tondi, de Franco, 2005 и др.]. Такой подход подавляет возникновение осцилляций скорости, не подкреплённых данными, и улучшает обусловленность систем уравнений, однако требование гладкости одновременно применяется ко всем частям модели, независимо от локальной разрешающей способности, а следовательно, сглаживается (хотя и в меньшей степени) и решение в хорошо освещённых лучами областях, что не позволяет получить здесь модели с максимально возможной степенью детальности.

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

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

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

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

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

Основная идея алгоритма сводится к следующему. Рассмотрим в качестве примера двумерное разложение аномалии медленности ds(i, j) в некотором интервале глубин (z, z + z), где i, j - индексы ячеек равномерной сетки в направлениях осей X и Y соответственно. Разложение по системе вэйвлет-функций Хаара имеет вид:

J 2J-1-1 ds (i, j) = c0 + ck k (i, j), (9) l,m,n l,m,n l=1 m,n=0 k=где k (i, j) = 2-l(k)((i - 1) /2l-1 - 2m, ( j - 1) /2l-1 - 2n), i, j = 1,.., N = 2J, (10) l,m,n где x обозначает целую часть x, (k) (p, q) - "материнскиеУ функции Хаара, представленные на рис. 3. Индексы m и n в формуле (10) определяют сдвиг базисной функции в направлении, задаваемом индексами i и j соответственно. Индекс l соответствует масштабу вэйвлета: чем больше l, тем больше размер носителя базисной функции.

Аналогичным образом разложение (9) используется и для представления вариаций глубин сейсмических границ dz(x, y).

Рис. 3. Материнские вэйвлет-функции Хаара (k) (p, q). а - k = 1; б - k = 2; в - k = 3.

Выбор базиса Хаара связан с финитностью их носителя - каждая функция q (i, j) вносит вклад в величину медленности в ограниченном объеме среды, за l,m,n пределами которого q 0. При этом чем больше величина l, тем меньше размер l,m,n носителя и тем больше детальность модели в области, определяемой положением носителя, т.е. индексами m, n. Это свойство вэйвлет-функций Хаара может быть эффективно использовано для адаптивной параметризации среды, согласованной с разрешающей способностью: исключая из разложения (9) те члены, коэффициенты при которых не могут быть устойчиво определены из-за недостаточной разрешающей способности в соответствующем объеме, получим аппроксимацию, эквивалентную неоднородной параметризации, детальность которой адаптирована к локальной разрешающей способности.

Актуальной остаётся проблема построения эффективного правила оценки локальной разрешающей способности в задачах сейсмической томографии. Наиболее распространённым методом подобного рода является оценка разрешающей способности на основании плотности сейсмических лучей в конкретном объёме среды [Kulakov et al., 1995; Bijwaard et al, 1998 и др.]. Однако увеличение числа лучей с близкими параллельными траекториями не приводит к увеличению разрешающей способности: соответствующая аномалия скорости может быть локализована в любом конечном объёме, который пересекают все эти лучи. Поэтому другим важнейшим параметром, характеризующим разрешающую способность, является Уугловое покрытиеФ сейсмических лучей, т.е. степень разброса их направлений в объёме среды. В качестве меры углового покрытия будем использовать представление о среднеквадратическом разбросе азимута лучей относительно их среднего направления, применяемое в статистике Фишера.

Предложенные меры разрешающей способности являются эмпирическими. Для исследования вопроса об их адекватности была проведена серия численных экспериментов на синтетических примерах, в ходе которых одновременно со значениями числа лучей Nr(l, m, n) и углового покрытия (l, m, n) вычислялась и матрица разрешения R [Яновская, 1997; Жданов, 2007]. В результате экспериментов установлено, что для получения единой для вэйвлетов всех порядков l корреляционной связи необходимо нормировать значения числа сейсмических лучей на величину, пропорциональную длине стороны носителя соответствующей вэйвлет-функции, т.е.

рассматривать значения Nr = Nr/2(J-l). Во всех случаях имеется статистически зна чимая зависимость между эмпирическими мерами Nr и и значениями диагональных элементов R. Показано, что для достижения наилучших результатов необходимо одновременно рассматривать обе меры разрешающей способности: нормированное число лучей Nr и их угловое покрытие в пределах носителя. Предложено соответствующее решающее правило, позволяющее отбраковывать плохо обусловленные неизвестные по данным о Nr и .

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

Далее в работе делается анализ возможности увеличения быстродействия алгоритма при помощи современных технологий параллельных вычислений. Одной из основных проблем лучевой сейсмической томографии является необходимость привлекать чрезвычайно большой объём данных [Яновская, 1997]. На практике интерпретация данных сейсмических экспериментов всегда требует проведения не одной серии расчётов, поэтому быстродействие алгоритма лучевой сейсмической томографии имеет принципиальное значение. Показано, что возможности его оптимизации с использованием параллельных вычислений весьма широки. Реализация параллельных вычислений была проведена для многопроцессорных (многоядерных) вычислительных систем с общей памятью, с использованием языка cilk++ - расширения языка С++.

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

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

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

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

Результаты трассирования лучей (рис. 4) демонстрируют высокую неоднородность лучевого покрытия как в плане, так и по глубине. Лучше всего освещена сейсмическими лучами область непосредственно под конусом вулкана: здесь максимально не только число лучей, но и их угловое покрытие. Кроме того, из-за рельефа вулканического конуса длина сейсмических лучей в этой области не превышает 5-км, что, вкупе с низкими скоростями сейсмических волн, приводит к малому размеру зоны Френеля. На флангах модели плотность сейсмических лучей уменьшается. Поэтому результаты исследований по проекту TOMOVES - весьма интересный объект с точки зрения опробования алгоритма с адаптивной параметризацией среды.

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

Рис. 4. Система наблюдений в ходе проекта TOMOVES и траектории сейсмических лучей, рассчитанные в модели начального приближения, построенной по априорным данным. а - горизонтальное сечение; б - вертикальное сечение В результате инверсии (рис. 5) обнаружено, что скорости сейсмических волн в пределах вулканического конуса повышены и достигают 4,25 км/с. Эти результаты находятся в хорошем согласии с результатами, полученными в работах [di Stefano, Chiarabba, 2002; Tondi, de Franco, 2003, 2006]. По имеющимся представлениям [Tondi, de Franco, 2003, 2006] эта аномалия маркирует лаву, застывшую в подводящих каналах во время предшествующего извержения.

Под северо-восточными склонами вулкана проявлена низкоскоростная аномалия (Vp = 1,5 - 2,0 км/с), охватывающая вулкан полукольцом. Низкоскоростные аномалии меньшей интенсивности (Vp = 2,0 - 2,25 км/с) также обнаруживаются под южным и западным склонами. Аналогичная структура низкоскоростных аномалий обнаружена в работах [di Stefano, Chiarabba, 2002; Tondi, de Franco, 2003, 2006]. В работе [Lomax et al., 2001] низкоскоростная аномалия обнаружена исключительно под северным склоном, что объясняется двумерной интерпретацией данных и низкой детальностью этой модели. Эти аномалии хорошо соответствуют рыхлым пеплово-туфовым толщам.

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

Наиболее интересная черта, обнаруженная по результатам данной работы, - низкоскоростной ФкореньУ вулкана в пределах известнякового слоя на глубинах 1,5 - 4,км. В предшествующих работах детальность модели на соответствующих глубинах не позволила выделить какие-либо структуры. Можно предположить, что эта аномалия соответствует слоям, нарушенным во время подъёма магматического материала и вулканических землетрясений. Это предположение подтверждается тем, что положение данной аномалии хорошо согласуется с положениями очагов вулканических землетрясений, определёнными в работе [Lomax et al., 2001].

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

Четвёртая часть главы 3 посвящена описанию и применению алгоритма пассивной сейсмической томографии, основанного на принципе адаптивной параметризации среды.

Задача пассивной сейсмической томографии распадается на три этапа: первичная локация очагов в рамках горизонтально-однородной модели v(0) (z); уточнение модели v (z) с одновременной перелокацией очагов; построение трёхмерно-неоднородной модели среды с дальнейшим уточнением параметров событий. Разработанный алгоритм позволяет решать весь комплекс перечисленных задач. Переход от этапа к этапу может осуществляться как на основе интерактивного взаимодействия с сейсмологом, по результатам визуального и статистического контроля моделей и Рис. 5. Строение вулкана Везувий и прилегающей области по данным лучевой сейсмической томографии. a - горизонтальное сечение модели на абсолютной отметке Z = 0; б, в - вертикальные сечения модели.

невязок, так и автоматически.

Алгоритм пассивной томографии применён к изучению строения очаговой зоны Рачинского землетрясения 1991 г. (Mw = 7.0). В работе приводятся сведения о сейсмотектонической позиции очага землетрясения и даётся обзор имеющихся представлений о его механизме. Вопрос о точном положении поверхности разрыва (в частности об угле её наклона) остаётся дискуссионным до настоящего времени (см.

монографию Т.П.Белоусова (2009)).

Временная сейсмологическая сеть, развернутая в эпицентральной зоне Рачинского землетрясения в 1991 г., состояла из 52 станций и зарегистрировала 3792 события [Арефьев и др., 1993]. Для томографической обработки были отобраны 24события, для которых удается устойчиво определить первые вступления по данным не менее чем шести сейсмических станций. Число независимых уравнений (число сейсмических лучей) составило около 27 000.

Размер базовой ячейки модели среды был выбран исходя из частотного состава записей сейсмических волн от афтершоков. Спектр колебаний P-волн, регистрируемых в первых вступлениях, достигает 10 Гц, скорость P-волн в первых 5Ц10 км разреза может составлять 3 км/с, соответственно длина волны может достигать 0,км и более. Поэтому размер базовой ячейки параметризации скорости принимался равным 0,3х0,3х0,3 км. Лучевое покрытие в объеме модели весьма неоднородно.

Максимальное лучевое покрытие приходится на глубины 8Ц10 км; в верхней части модели лучевое покрытие фрагментарно и наследует геометрию расположения сейсмических станций; глубже 12Ц13 км число очагов уменьшается и соответственно уменьшается лучевое покрытие.

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

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

На рис. 6 представлены результаты инверсии. Наиболее яркая черта полученного распределения скорости - наличие на глубине около 8 км зоны пониженных скоростей, которая протягивается с запада-северо-запада на восток-юго-восток вдоль оси Рача-Лехумского прогиба. По имеющимся представлениям именно в этой области происходило вспарывание во время основного события.

Гипоцентры афтершоков разбиваются на два вытянутых в том же направлении кластера, которые оконтуривают зону пониженных скоростей соответственно с юга и севера. При этом афтершоки, составляющие южный кластер, группируются преимущественно на глубинах менее 8Ц10 км; в вертикальных сечениях вдоль юго-запад - северо-восточного направления они образуют вытянутое облако. Ось этого облака наклонена на северо-восток под углом около 45; его нижняя часть достигает зоны пониженных скоростей. Афтершоки, составляющие северный кластер, группируются на глубинах более 8-10 км, образуя облако, также наклоненное на северо-восток под углом около 45.

Вероятно, три указанные структуры: зона пониженных скоростей и два кластера афтершоков - маркируют сложную структуру поверхности разрыва в очаге Рачинского землетрясения 1991 г. В сечениях по направлению юго-западЦсеверо-восток она выглядит как ступенчатая линия, отвечающая поверхности надвигания структур Рис. 6. Скоростная модель очаговой зоны Рачинского землетрясения. а,б - сечения в направлении ЮЗ-СВ (на рис. а дана сейсмотектоническая интерпретация полученной модели); в - сечение в направлении CЗ-ЮВ (положение сечений см. на рис. г.); г - сейсмотектоническая схема очаговой зоны и срез скоростной модели на глубине 9,5 км; цифрами 1 - 4 показаны механизмы и положение субочагов основного толчка, по данным из работы [Fuenzalida et al.

,1997] южного склона Большого Кавказа на Грузинскую глыбу Дзирульского массива. Эта гипотеза подтверждается и более высокими относительными скоростями в южной части модели.

Обнаруженные особенности строения очаговой зоны и распределения афтершоков подтверждают гипотезу, что основное смещение в очаге землетрясения апреля 1991 г. происходило по субгоризонтальной поверхности разрыва, которая, повидимому, сопряжена с поверхностью раздела осадочных пород, выполняющих РачаЛехумский прогиб, и нижележащего кристаллического основания. Эта поверхность, вероятно, представляла собой зону пониженной (относительно смежных участков) прочности, которая способствовала распространению разрыва. Нарушенная приразломная зона в томографическом изображении маркируется пониженными скоростями. Разрядка тектонических напряжений на этом участке, вероятно, произошла полностью, что объясняет отсутствие здесь афтершоковой активности. К югуЦюго-западу и северуЦсеверо-востоку от зоны основного разрыва смещение продолжилось по более круто падающим поверхностям в направлении вверх на южном фланге разломной зоны и вниз на северном. На этих участках энергии, выделившейся во время основного события, было недостаточно для полной разрядки напряжений, поэтому вдоль соответствующих поверхностей сейсмический процесс продолжился афтершоковой активностью.

Полученные результаты достаточно хорошо согласуются с представлениями о механизме очага, развитыми в [Fuenzalida et al.,1997] на основе моделирования записей объемных волн, а также подтверждают основные особенности ранее построенных томографических моделей [Арефьев и др., 2006].

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

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

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

В результате проведённых исследований сделаны следующие выводы.

По первой главе:

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

2. При определении высот пунктов наблюдений на основе нивелировок вычисляемые аномалии силы тяжести в свободном воздухе являются смешанными, причём величина косвенного эффекта сопоставима с величинами аномалий, интерпретируемых при поисковых работах на нефть и газ. Это требует модификации методов их количественной интерпретации. Определение высот пунктов при помощи систем глобального позиционирования (ГЛОНАСС, GPS) снимает проблему косвенного эффекта.

3. Ныне действующие нормативные и справочные документы по гравиметрической съёмке не соответствуют актуальному уровню точности и применяемым технологиям, необходимость их модификации широко признаётся в гравиметрическом сообществе [Бычков, 2005, 2010; Конешов и др., 2010].

4. Различие между величинами аномалии модуля магнитного поля T и гармонического поля T0 может достигать сотен нТл при обычных для платформенных районов аномалиях T в первые тысячи нТл, что может превосходить не только точность современных съёмок, но и амплитуду полезного сигнала, например аномалий, связанных с намагниченными горизонтами осадочного чехла и др.

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

6. Предложенный в настоящей работе алгоритм восстановления гармонического компонента T0 по полю T позволяет с высокой точностью решать эту задачу в диапазоне значений амплитуд аномального поля до 15 000 нТл (до 1/величины нормального поля), что покрывает значительную часть потребностей практики.

Выводы 1 - 6 обосновывают первое защищаемое положение.

По второй главе:

7. В результате статистических оценок, основанных на данных о точности спутниковых систем GRACE и GOCE, показано, что поставляемые этими спутниковыми системами данные о временных вариациях поля силы тяжести в диапазоне периодов от 30 суток и более могут быть использованы для выделения косейсмического эффекта при землетрясениях магнитудой 9 и уточнения их механизма.

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

9. Задача определения геометрических параметров сложной поверхности разрыва в очаге землетрясения по данным о косейсмических деформациях земной поверхности (которые могут быть получены при помощи радарной спутниковой интерферометрии - InSAR и систем глобального позиционирования - ГЛОНАСС, GPS) характеризуется значительным множеством Цэквивалентных решений и сложной формой поверхности минимума функционала невязки.

10. Разработанный алгоритм определения параметров поверхности разрыва, основанный на методе нелинейной минимизации Монте-Карло составного функционала невязки, позволяет устойчиво лоцировать положение глобального минимума при уровне сложности модели поверхности разрыва, адекватной современным представлениям. Алгоритм эффективен при использовании на современных ЭВМ с параллельной архитектурой.

Выводы 7 - 10 обосновывают второе защищаемое положение.

По третьей главе:

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

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

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

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

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

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

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

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

Выводы 11 - 17 обосновывают третье защищаемое положение.

18. Скорости сейсмических волн в пределах вулканического конуса Везувия повышены и достигают 4,25 км/с. Под склонами вулкана проявлена низкоскоростная аномалия (Vp = 1,5 - 2,25 км/с), охватывающая вулкан полукольцом. Эти особенности строения хорошо согласуются с результатами предшествующих интерпретаций [Lomax 2001; Di Stefano, Chiarabba, 2003; Tondi, de Franco, 2003, 2006], однако полученная в настоящей работе модель отличается несколько большей детальностью.

19. Разработанный алгоритм позволил достичь большей, сравнительно с предшествующими исследованиями, детальности в модели строения Везувия на глубинах более 1,5 км. В результате обнаружена низкоскоростная аномалия, расположенная под вулканической постройкой в диапазоне абсолютных глубин 1,5 - 4,0 км. Её положение соответствует области, где группируются очаги вулканических землетрясений.

20. Основная черта полученного распределения скорости в очаговой зоне Рачинского землетрясения - наличие на глубине около 8 км зоны пониженных скоростей, которая протягивается с запада-северо-запада на восток-юго-восток вдоль оси Рача-Лехумского прогиба.

21. Гипоцентры афтершоков Рачинского землетрясения разбиваются на два вытянутых в том же направлении кластера, которые оконтуривают зону пониженных скоростей соответственно с юга и севера. Афтершоки, составляющие южный кластер, группируются преимущественно на глубинах менее 8Ц10 км;

в вертикальных сечениях вдоль юго-запад - северо-восточного направления они образуют вытянутое облако, с осью, наклонённой на северо-восток под углом около 45. Афтершоки, составляющие северный кластер, группируются на глубинах более 8Ц10 км, образуя аналогичное облако.

22. Обнаруженные особенности строения очаговой зоны хорошо согласуются с геологическими данными и ранее полученными результатами [Арефьев и др., 2006а, б], но отличаются большей детальностью.

23. Полученные результаты подтверждают гипотезу, согласно которой разрыв в очаге Рачинского землетрясения 1991 г. развивался вдоль сложной поверхности, причем основная часть энергии выделилась в ходе подвижки по субгоризонтальной плоскости, залегающей в основании Рача-Лехумского прогиба, а меньшая - в ходе подвижки по двум более круто падающим ( 45) поверхностям, обрамляющим основной очаг с севера и юга. Этот вывод хорошо согласуется с моделью очага, предложенной ранее в [Fuenzalida et al., 1997] на основе моделирования записей объемных волн, и позволяет объяснить имеющуюся неоднозначность в результатах определений механизма очага в рамках модели точечного (дипольного) источника.

Выводы 18 - 23 обосновывают четвёртое защищаемое положение.

Список публикаций по теме диссертации Статьи в рецензируемых журналах, рекомендуемых ВАК РФ для публикации материалов докторских и кандидатских диссертаций:

1. Тихоцкий С. А., Фокин И. В., Шур Д. Ю. Активная лучевая сейсмическая томография с использованием адаптивной параметризации среды системой вэйвлетфункций // Физика Земли. 2011. № 4. С. 67Ц86.

2. Tikhotsky S., Achauer U. Inversion of controlled-source seismic tomography and gravity data with the self-adaptive wavelet parametrization of velocities and interfaces. // Geophys. J. Int. 2008. Vol. 172. Pp. 619Ц630.

3. Тихоцкий С. А., Ахауер У. Строение вулкана Везувий по данным активной сейсмической томографии - новые результаты интерпретации данных, полученных в ходе проекта TOMOVES // Вестник Камчатской региональной ассоциации Учебно-научный центр. Серия: Науки о Земле. 2011. № 1. С. 34Ц44.

4. Тихоцкий С. А., Фокин И. В., Шур Д. Ю., Арефьев С. С. Строение очаговой зоны Рачинского (1991 г.) землетрясения по данным локальной сейсмической томографии с адаптивной параметризацией среды // Геофизические исследования. 2011.

Т. 12, № 1. С. 5Ц33.

5. Левин Г. С., Тихоцкий С. А. О влиянии выбора системы высот на результаты высокоточных гравиметрических съёмок // Геофизика. 2003. № 5. С. 55Ц59.

6. Гордин В. М., Тихоцкий С. А., Шур Д. Ю. О восстановлении гармонического компонента аномалий модуля магнитного поля // Физика Земли. 2006. № 4.

С. 69Ц79.

7. Михайлов В. О., Тихоцкий С. А., Диаман М., Пане И. Исследование возможности обнаружения и изучения вариаций силы тяжести геодинамического происхождения по современным спутниковым гравиметрическим данным // Физика Земли. 2005. № 3. С. 18Ц32.

8. Mikhailov V., Tikhotsky S., Diament M., Pannet I., Ballu V. Can tectonic processes be recovered from new gravity satellite data? // Earth Plan. Sci. Letter. 2004. Vol.

228. Pp. 281Ц297.

9. Гордин В. М., Тихоцкий С. А., Курихина О. А., Платонова С. А. Применение пуассоновой модели источников аномального гравитационного поля для изучения океанской литосферы // Физика Земли. 2000. № 7. С. 76Ц88.

10. Kaban M. K., Schwintzer P., Tikhotsky S. A. A global isostatic gravity model of the Earth // Geoph.J.Int. 1999. Vol. 136. Pp. 519Ц563.

11. Фокин И. В., Басакина И. М., Капустян Н. К., Тихоцкий С. А., Шур Д. Ю. Опыт применения сейсмической томографии для археологических исследований оснований и фундаментов зданий // Вопросы инженерной сейсмологии. 2011. № 2.

С. 21Ц34.

12. Михайлов В. О., Назарян А. Н., Смирнов В. Б., Диаман М., Шапиро Н. М., Киселева Е. А., Тихоцкий С. А., Поляков С. А., Смольянинова Е. И., Тимошкина Е. П.

Совместная интерпретация данных дифференциальной спутниковой интерферометрии и GPS на примере Алтайского (Чуйского) землетрясения 27.09.2003 г. // Физика Земли. 2010. № 2. С. 3Ц16.

13. Tiberi C., Diament M., Daverchure J., Petit Mariani C., Mikhailov V., Tikhotski S., Achauer U. Deep structure of the Baikal rift zone revealed by joint inversion of gravity and seismology data // J.Geoph.Res. 2003. Vol. 108, no. B3. Pp. ETG1, 1Ц15.

14. Соловьев А. А., Шур Д. Ю., Гвишиани А. Д., Михайлов В. О., Тихоцкий С. А.

Определение вектора магнитного момента при помощи кластерного анализа результатов локальной линейной псевдоинверсии аномалий T // Докл. РАН. 2005.

Т. 404, № 1. С. 1Ц4.

15. Widiwijayanti C., Mikhailov V., Diament M., Deplus C., Louat R., Tikhotsky S., Gvishiani A. Structure and evolution of the Molucca Sea area: Constraints based on interpretation of a combined sea-surface and satellite gravity dataset // Earth Plan.Sci.Lett. 2003. Vol. 215. Pp. 1365Ц150.

Прочие публикации:

16. Tikhotsky S. Determination of the Sublithospheric Component in the EarthТs Anomalous Gravity Field // Cahiers of ECGS. 2003. Vol. 20. Pp. 79Ц85.

17. Тихоцкий С. А., Ашауер У. Комбинированная инверсия данных сейсмологии и гравиметрии в задаче определения положения геологической границы в трёхмерном случае // Геоинформатика. 2006. № 3. С. 25Ц28.

18. Тихоцкий С. А. Решение обратной кинематической задачи активной сейсмической томографии с использованием адаптивной параметризации среды системой вэйвлетов Хаара // Тез. докл. X геофиз. чтений им. В.В.Федынского. 2008. С. 75.

19. Tikhotsky S. A., Shur D. Y. Modern parallel computing technologies for the traveltime seismic tomography inversion // 8th international conference Problems of geocosmos.

2010. Pp. 183Ц184.

20. Tikhotsky S., Achauer U. Active seismic tomography inversion with the self-adaptive wavelet parameterization: algorithm and its application for the Vesuvius volcano structure // 7th international conference Problems of geocosmos. 2008. Pp. 251Ц252.

21. Tikhotsky S., Achauer U., Fokin I. Controlled-Source Seismic Tomography with Wavelets: Inversion Algorith and its Application to Vesuvius Volcano // Geophysical Research Abstracts / EGU General Assembly. Vol. 11, EGU2009-8626. 2009.

22. Tikhotsky S. A., Mikhailov V. O. Can the geodynamically induced gravity variations be detected by modern sattelite missions? // 6th international conference Problems of geocosmos. 2004. P. 257.

23. Гордин В. М., Тихоцкий С. А., Курихина О. А. Модели случайных источников в задачах интерпретации морских гравимагнитных данных // Материалы межд.

школы-семинара Вопросы теории и практики компл. геол. инт-ции грав., магн.

и электр. полей. М.: ОИФЗ РАН, 2001. С. 242Ц252.

24. Tikhotsky S. A., Fokin I. V., Shur D. Y., Arefiev S. S. Local traveltime tomography with the adaptive wavelet perameterization: algorithm and its application for the 19Racha (M=7,0) earthquake source area study // 8th international conference Problems of geocosmos. 2010. P. 166.

25. Fokin I. V., Tikhotsky S. A., Basakina I. M., Kapustyan N. K., Shur D. Y. The traveltime seismic tomography for the archaeology and engineering geophysics: methodological considerations and applikation for the Solovki island archaeological site // 8th international conference Problems of geocosmos. 2010. Pp. 184Ц185.

26. Гордин В. М., Тихоцкий С. А., Шур Д. Ю. О восстановлении гармонической компоненты поля скалярных магнитных аномалий // Вопросы теории и практики интерпретации гравитационных, магнитных и электрических полей: мат. 31-й сессии межд.сем. им. Д.Г.Успенского. М.: ОИФЗ РАН, 2004. С. 21Ц22.

27. Mikhailov V., Tikhotsky S., Pannet I., M. D. On the recovery of geodynamic signals from data of temporal variations of the global gravity field // Geophys.Res.Abstr.

Vol. 6. 2006.

28. Михайлов В. О., Тихоцкий С. А. Исследование возможности обнаружения и изучения вариаций силы тяжести геодинамического происхождения по современным спутниковым гравиметрическим данным // Вопросы теории и практики интерпретации гравитационных, магнитных и электрических полей: мат. 31-й сессии межд.сем. им. Д.Г.Успенского. М.: ОИФЗ РАН, 2004. С. 45Ц46.

Тихоцкий Сергей Андреевич Разработка математических методов и алгоритмов решения обратных задач геофизики и обработки геофизических данных. Автореф. дисс. на соискание учёной степени доктора физ.-мат. наук. Подписано в печать 17.06.2011 г. Формат 60х90/42. Усл.печ.л. 2. Тираж 200 экз. Изд. ИФЗ РАН.

Авторефераты по всем темам  >>  Авторефераты по земле