Бесплатный автореферат и диссертация по наукам о земле на тему

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

ВАК РФ 25.00.10, Геофизика, геофизические методы поисков полезных ископаемых

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

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

КАВАННИК Артем Валерьевич

ОПРЕДЕЛЕНИЕ СКОРОСТНОЙ СТРУКТУРЫ СРЕДЫ С ИЗВЕСТНЫМИ И НЕИЗВЕСТНЫМИ ИСТОЧНИКАМИ МЕТОДОМ СЕЙСМИЧЕСКОЙ ТОМОГРАФИИ БЕЗ ТРАССИРОВКИ ЛУЧЕЙ

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

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

НОВОСИБИРСК 2005

Работа выполнена в Геофизической службе Сибирского отделения Российской академии наук

Научный руководитель:

доктор геолого-минералогических наук Селезнев Виктор Сергеевич

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

доктор технических наук, профессор Зеркаль Сергей Михайлович доктор физико-математических наук Санина Ирина Альфатовна

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

Международный институт теории прогноза землетрясений и математической геофизики РАН (МНТП РАН, г. Москва)

Защита состоится 6 июля 2005 г. в 12 час. на заседании диссертационного совета Д 003.050.05 при Объединенном институте геологии, геофизики и минералогии им. A.A. Трофимука СО РАН, в конференц-зале.

Адрес: 630090, Новосибирск-90, пр-т Ак. Коптюга, 3

Факс: (383) 333-27-92

С диссертацией можно ознакомиться в библиотеке ОИГГМ СО РАН

Автореферат разослан мая 2005 г.

Ученый секретарь

диссертационного совета

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

Ю.А. Дашевский

и^-ч зашгг

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

Актуальность темы

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

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

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

^ОС НАЦИОНАЛЬНАЯ I

БИБЛИОТЕКА 1

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

Научная задача исследований — определение скоростного строения среды по временам прихода волн для случаев известных и неизвестных источников. Поставленная задача предполагает выполнение следующих этапов:

1 Изучить область определения оператора линейной сейсмической томографии и выяснить возможность его матричного представления

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

3 Адаптировать итерационную методику оценки точности скоростных моделей ЬЭС^ПА к разработанному методу сейсмической томографии без трассировки лучей.

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

Фактический материал и методы исследования

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

С помощью разработанного метода была проведена обработка

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

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

- времен вступлений Р- и S-волн, зарегистрированных временной сетью цифровых станций Акташского полигона от афтершоков Чуй-ского землетрясения;

Основной метод исследования — математическое моделирование. Исследование включало в себя: построение разностных схем решения систем дифференциальных уравнений, решение плохо обусловленных систем линейных уравнений с разреженными несимметричными матрица^ ми с помощью LSQR-алгоритма, определение параметров гипоцентров землетрясений методом Гейгера, применение пакета программ VELEST для определения одномерной оптимальной скоростной модели, нелинейная минимизация методом наименьших квадратов с изменяющимися весами, решение прямой кинематической задачи конечно-разностным способом, вычисление сингулярного разложения матриц и расчет матриц ковариации, проектирование баз данных сейсмологической информации и разработку пользовательского интерфейса к ним.

Защищаемые научные результаты

1. Доказана ограниченность и компактность линейного томографического оператора в двухмерном случае, в предположении о регулярности поля лучей.

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

3. Адаптирована к разработанному алгоритму и программно реализована итерационная методика оценки точности скоростных моделей LSQRA (Zhang, McMechan, 1995), основанная на вычислении матриц ковариации и разрешающей способности.

4. Алгоритм опробован на данных от афтершоков Чуйского землетрясения: выделены аномалии скорости продольных (до 0.15-0.25

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

Новизна работы. Личный вклад

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

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

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

- используя тот факт, что на каждой итерации по методу Ъ8<311 достаточно знания действия прямого и сопряженного операторов на конкретный вектор, развит подход, не требующий знания матричного представления томографического оператора;

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

- на основе метода разделения задачи на дискретную и непрерывную части метод сейсмической томографии без трассировки лучей адаптирован к задаче определения скоростного строения с неизвестными источниками;

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

- с применением гауссовского сглаживания разработатаны ре-гуляризующие процедуры, максимально использующие априорную информацию;

- с использованием модификации метода ЬЗС^И, реализован алгоритм расчета матриц ковариации и разрешающей способности.

2. Разработанный подход реализован в виде комплекса прикладных программ и библиотек.

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

Теоретическая и практическая значимость результатов

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

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

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

Апробация работы

Основные результаты работы докладывались на различных российских и международных конференциях, наиболее значимыми из которых являются: Генеральная ассамблея IUGG (Саппоро, Япония, 2003 г.), Международная конференция "Проблемы сейсмологии Ш-ого тысячелетия" (Новосибирск, 2003 г.), Международная конференция "Математические методы в геофизике" (Новосибирск, 2003 г.), Международная научная конференция, посвященная 90-летию академика H.H. Пу-зырева (Новосибирск, 2004 г.), а также на геофизическом семинаре ИГФ ОИГГМ СО РАН, 29 сентября 2004 г.

Объем и структура работы

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

Автор выражает глубокую благодарность и особую признательность своим учителям: доктору геолого-минералогических наук B.C. Селезневу и кандидату физико-математических наук В.А. Чеверде. Автор благодарен доктору физико-математических наук, академику РАН C.B. Голь-дину за постоянное внимание, ценные советы и замечания, а также кандидату геолого-минералогических наук В.М. Соловьеву, доктору технических наук А.Ф. Еманову и доктору геолого-минералогических наук В.Д. Суворову за многочисленные конструктивные замечания и помощь в работе над диссертацией.

Автор выражает признательность А.Г. Филиной и М.С. Резепиной за моральную поддержку, а также В.И. Самойловой за неоценимую помощь в работе над рукописью.

ОСНОВНОЕ СОДЕРЖАНИЕ РАБОТЫ

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

Глава 1. Локальная сейсмическая томография

^ В главе проводится анализ работ, посвященных задачам одновремен-

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

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

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

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

• Отбор и предварительная обработка входных данных

• Построение начальной референтной скоростной модели

• Параметризация модели

• Решение прямой кинематической задачи

• Разделение системы на дискретную и непрерывную части

• Обращение систем линейных уравнений

• Оценка разрешающей способности метода и ошибок в решении

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

Глава 2. Определение скоростной структуры среды с известными и неизвестными источниками

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

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

-о

с,

Здесь со - распределение скорости в невозмущенной среде, известная величина. Таким образом, обращение невязок времен пробега в рамках линейной сейсмической томографии сводится к необходимости решения линейного операторного уравнения первого рода (1). Стандартным подходом к переходу от континуальной к дискретной постановке является покрытие целевой области равномерной прямоугольной сеткой, в каждой ячейке которой неизвестная скорость предполагается постоянной. Однако для возможности его применения и обеспечения сходимости решения получаемой на этом пути системы линейных алгебраических уравнений к решению исходного континуального уравнения томографический оператор А должен быть компактным из £2 в ¿2- В разделе доказывается непрерывность и компактность оператора в пространстве ¿2 в предположении о регулярности поля лучей.

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

Действие прямого оператора п = А{5с), ищется в виде решения следующего дифференциального уравнения, полученного из уравнения эйконала:

а(х,у, г)^ + Ь(х,у, г)^ + с(х, у, г)= п0(х,у, г)щ(х,у, г)

пи = о

Здесь, то(х,у, г) - функция поля времен в невозмущенной среде с распределением медленности По = 1/со, а(х,у,г) = Ь(х, у, г) =

Q%j-,c(x,y,z) = ^ - частные производные поля времен по пространственным координатам, а щ = — Ц - возмущения медленности. Анало-

со

гично, можно выписать дифференциальное уравнение для сопряженного томографического оператора v = А*{т{)-.

Ф, У, z) ^ + Ь(х, у, z) ^ 4- с(х, у, z) ^ + vAt0 = О

v\t = ti

Дифференциальные уравнения (2) и (3) решаются конечно-разностными методами. Это позволяет обойтись без обычной в таких случаях трассировки лучей и не требует хранения в оперативной памяти элементов матрицы А, аппроксимирующей оператор А. Более того, это не требует вычисления элементов матрицы.

Важно заметить, что в дифференциальные уравнения (2) и (3) входят времена пробега волн, рассчитанные для невозмущенной среды. Для расчета времен пробега волн используется конечно-разностный метод решения уравнения эйконала (Podvin, Lecompte, 1993), с помощью которого определяются времена прихода волн в узлах равномерной сетки, покрывающей исследуемую область. При этом для обеспечения хорошей точности численного дифференцирования необходимо достаточно детальное разбиение области, а это ведет к значительным затратам оперативной памяти компьютера. Чтобы преодолеть это ограничение исследуемая область разбивается на равные блоки. Времена прихода волн рассчитывается последовательно от блока к блоку, а в оперативной памяти сохраняются только времена в узлах, принадлежащих границам блоков, которые считаются вторичными источниками волн. Наиболее часто используемые блоки хранятся в кэше — области памяти фиксированного размера, что позволяет оптимизировать процедуру решения приямой задачи.

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

(3)

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

В разделе "Смешанная задача: определение скоростной модели среды с одновременным уточнением параметров гипоцентров землетрясений" изложенный подход к определению скоростного строения среды обобщается на случай неизвестных источников. Задачи подобного типа называются смешанными задачами, так как содержат два типа неизвестных параметров: дискретные координаты гипоцентров землетрясений и непрерывную функцию распределения скорости в среде (Pavlis, Booker, 1980). Задача решается методом Ньютона: она сводится к решению системы линейных уравнений относительно поправок к гипоцентральным параметрам Sh и скоростным парамерам 8с на каждом шаге итерационного процесса:

7 = H8h + А(5с),

(4)

где 7 - вектор невязок времен прихода волн. Для разделения системы (4) на дисекретную и непрерывную части в работе используется идея Павлиса и Букера, состоящая в подборе такого преобразования системы линейных уравнений, которое позволило бы получить две независимых системы относительно гипоцентральных и скоростных параметров. Таким преобразованием служит умножение системы (4) справа и слева па матрицу левых сингулярных векторов и матрицы Н:

0

иЬ

XJjUSh О

+

и 1А

U %А

{6с).

(5)

Матрица Ч^ содержит первые г столбцов и, а матрица Иц — оставшиеся 771 — г столбцов (г - ранг матрицы Н). Таким образом, смешанная задача решается с помощью итерационного метода Ньютона, состоящего в поочередном уточнении параметров источников и скоростной модели на каждом шаге:

Ск+1 = Ск 4- 8ск+1 = ск + (и0Л)+ (й - д(кк, ск)) кк+1 =Ьк- Н¡А(ёск+1) (6)

= 1,2...)

Глава 3. Обращение синтетических и реальных данных

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

Рис. 1: Обращение данных физического моделирования: а) исходная телесейсмическая модель, б) решение, полученное за 20 итераций ЬвС^Я

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

Первая модель отображает телесейсмическую ситуацию - в однородной среде находится объект, выполненный в виде плюма, с пониженной скоростью распространения волн (см. Рис. 1 (а)). Параметры системы наблюдения следующие: пункты возбуждения в количестве 21 расположены на глубине 800 км с шагом 80 км. Шаг между пунктами приема, расположенными на дневной поверхности составлял 40 км. Референтной моделью являлся однородный слой, имеющий скорость продольной волны во вмещающей среде.

Анализ результатов применения метода ЬЭСЩ показывает, что достаточно уверенно выделяются контуры объекта, размещенные на фоне неоднородной яркости освещения целевой области (см. Рис. 1 (б)). При этом кроме восстановления геометрии объекта, восстанавливается и значение скорости распространения волн внутри него.

а) б)

Рис. 2: Обращение данных физического моделирования: а) исходная межскважинная модель, б) решение, полученное за 7 итераций ЬЭСЗН

Вторая модель, соответствующая межскважинной томографии изображена на Рис. 2 (а). Модель включает в себя субгоризонтальные слои различной формы с повышенной (V = 2.48 км/с) и пониженной (V = 2.29 км/с) скоростью продольной волны относительно вмещающей среды (V = 2.38 км/с) и имеет нарушение в виде сброса под углом 45°. Мощности слоев лежат в интервале от 4 до 60 м. Источники и приемники расположены в двух скважинах, разнесенных на 200 м. Параметры системы наблюдения следующие: пункты взрыва в количестве 45 расположены с шагом 10 м; пункты приема расположены с шагом 5 м.

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

Раздел "Сейсмическая томография эпицентральной зоны Чуйского землетрясения" посвящен опробованию разработанного алгоритма на синтетических и реальных данных. Записи афтершоков Чуйского землетрясения, произошедшего 27 сентября 2003 года, зарегистрированные сетью станций Акташского сейсмологического полигона АСФ СО РАН (Еманов А.Ф., Еманов A.A., Филина А.Г. и др., 2003), несут в себе уникальную информацию о глубинном строении земных недр в эпицентральной зоне. Однако то обстоятельство, что эпицентральная область Чуйского землетрясения оказалась довольно слабо окружена станциями привело к необходимости проведения численных экспериментов с синтетическими данными, чтобы получить ответы на следующие вопросы:

• Какова разрешающая способность данных и как она связана с системой наблюдений?

• Как зависит решение от выбора стартовой скоростной модели?

Истинная скоростная модель, Vp

0 10 20 30 40 50 60 70 80 км Рис. 3: Численные эксперименты: исходная модель

Входные данные для обращения насчитывались следующим образом: глубинам всех очагов землетрясений из каталога присваивалось некоторое фиксированное значение, а затем времена прихода волн из каталога заменялись на времена пробега от "обновленных" гипоцентров, рассчитанные для тестовой однородной модели среды с включением, имеющим положительную скоростную аномалию (см. Рис. 3). Скорость Р-волн во вмещающей среде составляла 6.1 км/с (3.5 км/с для Б-волн), а в аномалии - 6.3 км/с (соответственно, 3.61 км/с). Значения скоростей, полученные на втором шаге свидетельствуют о расходимости ньютоновского процесса, хотя границы аномалии восстанавливаются достаточно уверенно. Непосредственно над востановленной неоднородностью появляется артефакт с отрицательной скоростной аномалией, связанный с недостаточной освещенностью области в горизонтальном направлении. В общем случае обе эти аномалии являются равнозначными и отдать предпочтение какой-либо из них отсутствии дополнительной априорной информации нельзя. Чтобы ответить на вопрос, какая из аномалий является "истинной", £1 Кд>КйЯ — артефактом, из решения после первой итерации поочередно исключаются выделенные аномалии, а полученные модели используются в качестве начальных приближений. Сравнение полученных за четыре итерации решений, для обоих стартовых моделей позволяет сделать вывод о том, что отрицательная аномалия в верхней части является артефактом: положительная аномалия сохраняется в обоих моделях, а отрицательная проявляется только в первом случае. Использование этого результата позволяет ввести дополнительные ограничения на решение в виде более сильного сглаживания в районе возникновения артефакта. Решение, полученное на десятом шаге показывает, что итерационный процесс сходится и даже позволяет получить близкие к исходным значения скоростей (см. Рис. 4).

Для изучения влияния выбора стартовой скоростной модели на решение, были заданы две начальные однородные модели со скоростями распространения волн выше (6.3 км/с) и ниже (6.0 км/с) чем во вмещающей среде в "истинной" модели. В Таблице 1 приведены сравнительные значения средневзвешенной квадратичной невязки для различных начальных скоростных моделей полученные на десятой итерации. Минимальное значение среднеквадратичной невязки, соответствует скоростной модели с "истинной" вмещающей средой. А максимальное, - случаю с завышенной начальной скоростной моделью, что вероятно, объясняет-

Vp, итерация 10

0 10 20 30 40 50 60 70 80 км

Рис. 4: Численные эксперименты, решение со сглаживанием

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

Начальная скоростная модель (км/с) Квадратичная невязка (с)

6.1 ("истинная") 0.0344

6.0 ("заниженная") 0.0814

6.3 ("завышенная") 0.1487

Таблица 1: Значения среднеквадратичной невязки RMS, полученные для различных референтных моделей на десятой итерации

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

Эпицентральная зона Чуйского землетрясения расположена в приграничной зоне на юге республики Алтай, между Северо-Чуйским хребтом и Курайской и Чуйской межгорными впадинами, разделенными Чаган-Узунским блоком (см Рис. 5).

Для построения скоростных моделей из каталога землятрясений были отобраны времена первых вступлений волн Pg и Sg, от 291 события, зарегистрированного временной сетью из 32 станций Алтайского сейсмологического полигона в 2003 году.

Рис. 5: Эпицентральная зона Чуйского землетрясения, сеть станций расположение профиля, а также карта разломов по данным ГИН РАН (Ема-нов А.Ф., Еманов A.A., Филина А.Г. и др., 2003).

Решение задачи определения скоростного строения среды по невязкам времен вступлений имеет сильную зависимость от выбора начальной скоростной модели. Мы используем концепцию так называемой "минимальной одномерной модели" (Kissling, et al., 1994), получаемой в процессе минимизации суммы квадратов невязок времен вступлений, реализованном в пакете программ VELEST. Процесс состоит в пошаговом уточнении одномерной скоростной модели методом проб и ошибок с учетом априорной информации и наших представлений о среде. Скоростная модель, соответствующая минимуму суммы квадратов невязок, и рассчитанные для нее гипоцентры землетрясений используются как стартовая модель для последующих уточнений.

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

Разрез модели Ур, итерация 3

Рис. 6: Скоростная модель Р-волн эпицентральной зоны, полученная за три итерации. Через N обозначено количество землетрясений, для которых погрешность определения глубины сгг попадает в соответствующий интервал.

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

• обратно пропорционально эпицентральным расстояниям расстоянию;

• обратно пропорционально абсолютной величине невязки.

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

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

вышающую 5 км и поэтому используется для уточнения скоростной модели.

Разрез распределения разрешенное™ модели Ур итерация 3

Рис. 7: Разрешенность модели Ур.

Сопоставив модель продольных волн, полученную за три итерации (Рис. 6), с диагональными элементами матрицы разрешающей способности (Рис. 7), можно сделать вывод о том, что отрицательная аномалия на юге области скорее всего является артефактом, т.к. разрешенность в этой области практически нулевая, в то время как положительные аномалии Р-волн на севере соответствуют освещенным областям. Отсутствие в южной части отрицательной аномалии поперечных волн, вероятно, объясняется только недостатком данных, зарегистрированных на станциях ЕЬСТ ("Елангаш") и БЫ^ ("Солонешенская"). Для этих станций снято в два раза меньше поперечных волн, чем продольных.

ЗАКЛЮЧЕНИЕ

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

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

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

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

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

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

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

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

в области вычислительных технологий.

ПУБЛИКАЦИИ ПО ТЕМЕ ДИССЕРТАЦИИ

1. Кабанник A.B., Орлов Ю.А., Чеверда В.А. Численное решение задачи линейной сейсмической томографии на проходящих волнах: случай неполных данных // Сибирский журнал индустриальной математики. - 2004.- т.7.- № 2. С. 54-67.

2. Кабанник A.B. Трехмерная нелинейная сейсмическая томография без трассирования лучей: одновременное восстановление скоростной структуры и гипоцентров землетрясений в случае неполного набора данных / ГС СО РАН - Новосибирск. 2004. - 9 с. - Деп. В ВИНИТИ 20.10.2004 г., N1640-B2004.

3. Kabannik А. V., Tcheverda V. А. 3D linear seismic tomography without ray tracing // IUGG 2003 Abstracts, Sapporo, Japan. P. 405.

4. Кабанник A.B., Чеверда В.А. Трехмерная нелинейная сейсмическая томография без трассирования лучей: одновременное восстаг новление скоростной структуры и гипоцентров землетрясений // Математические методы в геофизике: Труды междунар. конф. -Новосибирск: Изд. ИВМиМГ СО РАН, 2003. Т. 1. С. 219-224.

5. Кабанник A.B., Чеверда В.А. Трехмерная нелинейная сейсмическая томография без трассирования лучей: одновременное восстановление скоростной структуры и гипоцентров землетрясений // Проблемы сейсмологии Ш-ого тысячелетия: Мат. междунар. конф. - Новосибирск: Изд. СО РАН, 2003. С. 382-389.

6. Фатеев A.B., Филина А.Г., Кабанник A.B. Представительность и точность определения параметров эпицентров в АСОМСЭ // Проблемы сейсмологии Ш-ого тысячелетия: Мат. междунар. конф. -Новосибирск: Изд. СО РАН, 2003. С. 145-153.

7. Кабанник A.B. Сейсмическая томография афтершоковой зоны Чуйского землетрясения // Междунар. науч. конф. посвящ. 90-летию акад. Н. Н. Пузырева: Сб. докл. - Новосибирск: Изд. СО РАН, 2004. С. 450-456.

8 Кабанник А В. БУЛ-анализ обратной кинематической задачи в линейном приближении // XXXV МНСК. Геология: Сб. докл. - Новосибирск: НГУ, 1997. С. 48-50.

9. Кабанник А.В. SVD-анализ обратной кинематической задачи в линейном приближении // Геофизикаг97: Тез. докл. Всеросс. научно-практ. конф. молод, уч. и спец. - СПб.: ВИРГ-Рудгеофизика, 1997. С. 12-13.

10 Селезнев В С , Кабанник А.В. Томографический алгоритм определения трехмерного скоростного строения земной коры по данным локальных землетрясений с одновременным уточнением параметров гипоцентров //50 лет ГСЗ: прошлое, настоящее и будущее: Тез. Междунар. совещ. - М.: Изд. ОИФЗ РАН, 1999. С. 22.

11 Seleznev V , Kabannik A. An algorithm for reconstruction of three-dimensional velocity structure form local earthquake data with simultaneous determination of earthquake coordinates and origin time / / PAPHOS99: An International Workshop on Tomographic Imaging of 3-D Velocity Structures and Accurate Earthquake Location, Paphos, Cyprus, 1999. P. 501.

«

«

г

к

lr

_Технический редактор P.M. Вараксина_

Подписано к печати 13.05.2005 « Бумага 60x84/16. Бумага офсет № 1. Гарнитура 'Тайме". Печать

офсетная

_Печ. л. 1,4. Тираж 100. Заказ № 142_

Издательство СО РАН. 630090, Новосибирск, Морской пр , 2 Филиал "Гео". 630090, Новосибирск, пр-т. Ак Коптюга, 3

1 3 ^^ Русский фонд

2006-4 10618

Содержание диссертации, кандидата физико-математических наук, Кабанник, Артем Валерьевич

ВВЕДЕНИЕ

Глава 1 ЛОКАЛЬНАЯ СЕЙСМИЧЕСКАЯ ТОМОГРАФИЯ

1.1 Обзор исследований в области обратных кинематических задач сейсмики.

2.1 Определение глубинного строения но данным от близких землетрясений

2.1.1 Отбор н предварительная обработка входных данных

2.1.2 Построение начальной референтной скоростной модели

2.1.3 Параметризация модели.

2.1.4 Решение прямой кинематической задачи.

2.1.5 Разделение системы на дискретную и непрерывную части

2.1.6 Обращение систем линейных уравнений, регуляризация

2.1.7 Оценка разрешающей способности метода и ошибок в решении

Глава 2 ОПРЕДЕЛЕНИЕ СКОРОСТНОЙ СТРУКТУРЫ СРЕДЫ С

ИЗВЕСТНЫМИ И НЕИЗВЕСТНЫМИ ИСТОЧНИКАМИ

1.2 Численное решение задачи линейной сейсмической томографии

1.2.1 Исследование прямого томографического оператора в двухмерном случае.

1.2.2 Компактность линейного томографического оператора

1.2.3 Выбор численного алгоритма для решения задачи линейной сейсмической томографии.

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

1.2.5 Метод ЬЭС^КА, получение оценок матриц разрешающей способности, информационной плотности и ковариации.

2.2 Смешанная задача: определение скоростной модели среды с одновременным уточнением параметров гипоцентров землетрясений

2.2.1 Постановка задачи.

2.2.2 Линеаризация системы.

2.2.3 Разделение системы.

Глава 3 ОБРАЩЕНИЕ СИНТЕТИЧЕСКИХ И РЕАЛЬНЫХ ДАННЫХ

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

1.3.1 Физическое моделирование.

1.3.2 Аппроксимация линейного томографического оператора. Его сингулярный спектр.

1.3.3 Результаты обработки данных физического моделирования . 51 2.3 Сейсмическая томография эиицентральной зоны Чуйского землетрясения

2.3.1 Численные эксперименты.

2.3.2 Обращение данных афтершоков Чуйского землетрясения

Введение Диссертация по наукам о земле, на тему "Определение скоростной структуры среды с известными и неизвестными источниками методом сейсмической томографии без трассировки лучей"

Актуальность темы

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

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

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

Научная задача исследований — определение скоростного строения среды но временам прихода волн для случаев известных и неизвестных источников. Поставленная задача предполагает выполнение следующих этапов:

1. Изучить область определения оператора линейной сейсмической томографии и выяснить возможность его матричного представления.

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

3. Адаптировать итерационную методику оценки точности скоростных моделей ЬБСЗНА к разработанному методу сейсмической томографии без трассировки лучей.

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

Фактический материал и методы исследования

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

С помощью разработанного метода была проведена обработка

- данных, полученых в результате физического моделирования;

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

- времен вступлений Р- и Б-волн, зарегистрированных временной сетью цифровых станций Акташского полигона от афтершоков Чуйского землетрясения;

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

Защищаемые научные результаты

1. Доказана ограниченность и компактность линейного томографического оператора в двухмерном случае, в предположении о регулярности поля лучей.

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

3. Адаптирована к разработанному алгоритму и программно реализована итерационная методика оценки точности скоростных моделей LSQRA (Zhang, McMcchan, 1995), основанная на вычислении матриц ковариацин и разрешающей способности.

4. Алгоритм опробован на данных от афтершоков Чуйского землетрясения: выделены аиомалии скорости продольных (до 0.15-0.25 км/с) и поперечных (до 0.1 км/с) волн, со средним горизонтальным размером около 30 км при разрешающей способности обусловленной неравномерной плотностью расположения станций не более 0.16.

Новизна работы. Личный вклад

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

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

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

- используя тот факт, что на каждой итерации по методу ЬБРИ достаточно знання действия прямого и сопряженного операторов на конкретный вектор, развит подход, не требующий знания матричного представления томографического оператора;

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

- на основе метода разделения задачи на дискретную и непрерывную части метод сейсмической томографии без трассировки лучей адаптирован к задаче определения скоростного строения с неизвестными источниками;

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

- с применением гауссовского сглаживания разработатаны регуляризую-щие процедуры, максимально использующие априорную информацию;

- с использованием модификации метода иЗС^Ы реализован алгоритм расчета матриц ковариации и разрешающей способности.

2. Разработанный подход реализован в виде комплекса прикладных программ и библиотек.

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

Теоретическая и практическая значимость результатов

Предложенный томографический алгоритм является инструментом, позволяющим определять глубинное строение но кинематическим данным в сложных трехмерных средах применительно к широкому диапазону задач геофизики: от построения инженерных сейсмических разрезов до получения трехмерных скоростных моделей по сейсмологическим данным. Применение итерационного метода ЬБС^И. позволяет производить обработку больших объемов данных без значительных затрат оперативной памяти и ресурсов процессора. Адаптированная итерационная методика вычисления матриц ковариации и разрешающей способности делает возможной оценку точности трехмерных скоростных моделей для задач сейсмической томографии с большим количеством наблюдений и неизвестных параметров.

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

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

Апробация работы

Основные результаты работы докладывались на различных российских и международных конференциях, наиболее значимыми из которых являются: Генеральная ассамблея IUGG (Саппоро, Япония, 2003 г.), Международная конференция "Проблемы сейсмологии 111-ого тысячелетия" (Новосибирск, 2003 г.), Международная конференция "Математические методы в геофизике" (Новосибирск, 2003 г.), Международная научная конференция, посвященная 90-летию академика H.H. Пузырена (Новосибирск, 2004 г.), а также на геофизическом семинаре ИГФ ОИГГМ СО РАН, 29 сентября 2004 г.

Объем и структура работы

Данная диссертационная работа состоит из трех глав, заключения и приложения.

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

Глава 2 посвящена задаче определения скоростного строения среды по кинематическим данным для случаев известных и неизвестных источников сейсмических волн. В первой части рассматривается задача сейсмической томографии в линеаризованном виде. Описываются свойства линейного томографического оператора для случая двух пространственных переменных; показывается возможность его конечномерной аппроксимации. Однако основное внимание в данной главе уделено изложению конечно-разностного метода сейсмической томографии без трассировки лучей, основанного на применении итерационного метода LSQR: описывается способ вычисления прямого и сопряженного томографических операторов, а затем способ получения их конечномерной аппроксимации. Приводится адаптированный к разработанному методу способ оценки точности скоростных моделей, основанный па итерационном методе расчета матриц разрешающей способности и ковнри-ации [81]. Во второй части главы изложенный подход к определению скоростного строения среды обобщается на случай неизвестных источников: рассматривается задача одновременного определения скоростного строения среды и координат гипоцентров землетрясений. Следуя [63], эта задача называется смешанной, и решается путем разделения системы линейных уравнений на дискретную и непрерывную части.

Основными целями третьей главы являются: тестирование разработанного комплекса программ на данных физического моделирования; опробывание метода на синтетических моделях, реальных данных. Первая часть главы посвящена обращению данных физического моделирования [18], позволивших дать дополни-телиую оценку качества работы метода. Во второй части главы основное внимание уделено обращению времен пробега от афтершоков Чуйского землетрясения, зарегистрированных временной сетью станций Алтайского полигона АСФ ГС СО РАН. Проведенные численные эксперименты с использванием синтетических моделей, позволили определить границы применимости алгоритма для обращения данных, зарегистрированных на сети полигона [14]. В конце главы приводятся трехмерные скоростные модели Р- и S-волн и параметры гипоцентров землетрясений для эпицентральной зоны, а также оценки их точности.

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

В приложение вынесены: краткое описание сингулярного разложения матриц; последовательность шагов метода LSQR; описание итерационного метода Гейгера.

Автор выражает глубокую благодарность и особую признательность своим учителям: доктору геолого-мипералогических паук B.C. Селезневу и кандидату физико-математических наук В.А. Чеверде. Автор благодарен доктору физико-математических наук, академику РАН C.B. Гольдину за постоянное внимание, цепные советы и замечания, а также кандидату геолого-мипералогических паук В.М. Соловьеву, доктору технических наук А.Ф. Еманову и доктору геолого-мипералогических наук В.Д. Суворову за многочисленные конструктивные замсчания и помощь в работе над диссертацией.

Автор выражает признательность А.Г. Филиной и М.С. Резепиной за моральную поддержку, а также В.И. Самойловой за неоценимую помощь в работе над рукописью.

Локальная сейсмическая томография

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