Обратная задача сейсмического зондирования с использованием поверхностных волн тема диссертации и автореферата по ВАК РФ 05.13.18, кандидат наук Ян Цзяньсюнь

  • Ян Цзяньсюнь
  • кандидат науккандидат наук
  • 2019, ФГБОУ ВО «Московский государственный университет имени М.В. Ломоносова»
  • Специальность ВАК РФ05.13.18
  • Количество страниц 131
Ян Цзяньсюнь. Обратная задача сейсмического зондирования с использованием поверхностных волн: дис. кандидат наук: 05.13.18 - Математическое моделирование, численные методы и комплексы программ. ФГБОУ ВО «Московский государственный университет имени М.В. Ломоносова». 2019. 131 с.

Оглавление диссертации кандидат наук Ян Цзяньсюнь

Введение

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

1.1 Уравнение распространения упругих волн в слоистой среде

1.2 Постановка задачи микросейсмического зондирования

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

2.1 Метод расчета характеристики бегущих поверхностных волн в полупространстве

2.2 Метод расчета характеристики бегущих поверхностных волн в слоистой среде

2.3 Метод тензора сейсмического импеданса

2.4 Сравнение между методом МОК-О/П и методом тензора сейсмического импеданса

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

3.1.Предварительная обработка сейсмических данных одиночной станции

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

3.3. Метод измерения дисперсии поверхностных волн

3.3.1 Техника множественных фильтров (Multiple Filter Technique)

3.3.2 Применение взаимнокорреляционной функции в обработке

сейсмических данных двух станций

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

4.1 Постановка обратной задачи сейсмического зондирования

4.1.1 Генетический алгоритм

4.1.2 Регуляризация обратной задачи с учетом контрастности

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

4.1.3. Принцип невязки (метод выбора параметра регуляризации)

4.1.4. Метод «генетический алгоритм, связанный с методом регуляризации Тихонова (ГАМРТ)» для исследования обратной задачи

дисперсии волн Релея

4.2 Численное исследование обратной задачи

4.2.1. Анализ влияния изменения отношения параметров слоистых сред к дисперсии

4.2.2. Экспериментальные результаты инверсии дисперсии поверхностных волн по классическому генетическому алгоритму

4.2.3. Экспериментальные результаты с помощью метода ГАМРТ для

инверсии поверхностных волн в слоистой среде

4.2.4Имитационный эксперимент для определения параметров

регуляризации

4.2.5 Инверсия измеренной дисперсионной кривой с помощью метода

ГАМРТ и мониторинг изменений параметров о слоистой среде во

времени

Заключение

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

Введение

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

Землетрясение - это не только катастрофа для людей, но и яркий свет, который освещает внутреннее строение Земли. Оно является одним из самых важных способов чтобы понять Землю. Фактически, до сих пор многие из основных открытий о внутренней структуре Земли были получены из информации сейсмических волн. История изучения структуры, состава и состояния внутри Земли насчитывает более 100 лет: В 1906 году Олдхэм [59] пытался вывести всю внутреннюю структуру Земли с момента прохождения сейсмической волны через Землю. В 1909 году Андрия Мохоровичич [17] установил существование поверхности раздела на глубине 56 км, основанный на времени движения первого прибытия землетрясения, откуда появилось название границы Мохоровичича. В 1914 году Гутенберг подтвердил существование ядра, основанного на «теневой зоне» сейсмической объемной волны, и измерил глубину границы между земным ядром и мантией до 2900 км. Так как ядро не может распространять поперечные волны, то сейсмологи предположили, что ядро состоит из жидкости. В 1936 году Леман обнаружил, что в жидком ядре имеется твердое ядро Земли [48].

Китай является одной из стран в мире с сильной сейсмической активностью и сильными землетрясениями. Статистика показывает, что около 35% мировых континентальных землетрясений выше 7 баллов происходят в Китае. Среди 1,2 миллиона человек, которые погибли в 20 веке во всем мире из-

4

за землетрясений, на Китай приходилось 590000 человек. Кроме того, почти 1/3 территории Китая расположены в районах с высокой сейсмической активностью выше VII баллов по Шкале Меркалли, а районы с высокой интенсивностью выше 7 баллов в США составляют лишь 12% площади суши. С 1900 года по 2007 год 70 землетрясений магнитудой 7,0-7,9 произошло на материке Китая и 6 землетрясений выше 8 баллов. В результате бедствий, вызванных этими землетрясениями, было охвачено 28 провинций с 590 000 смертей, 760 000 инвалидов и более миллиардов пострадавшего населения.

Рисунок 1. Топографическая карта, включающая материк Китая, и распространение землетрясения в последние годы.

Распространение сильных землетрясений на материке Китая является неравномерным. Самая поразительная особенность заключается в том, что сейсмическая активность на западе Китая сильнее, чем на востоке. Из исторической записи, на 107 градусов к востоку от границы, было 91 сильное землетрясение более 7 баллов на западе и только 27 на востоке. (Исключая глубокие землетрясения в провинции Тайвань и в Северо-Востоке Китая, в котором источник находится в глубине земных недр более 300 км.). Однако из-за густонаселенности и экономического развития в восточном районе страны, бедствия, вызванные землетрясением, намного больше, чем на западе.

Причиной этой сейсмической активности является то, что основная движущая сила тектонической деформации в материке Китая происходит от толчка индийской плиты над Тибетским плато (рис. 1).

Учитывая сейсмическую активность, распределение плотности населения Китая и степень вреда, причиненного землетрясениям для человеческого общества, и другие причины, в конце 90 годы 20-ого века с быстрым развитием сейсмических цифровых технологий китайская сейсмическая администрация начала упорядочивать оборудование сейсмического наблюдения (включают поточные сейсмометры и стационарные сейсмические станции) для мониторинга землетрясения и оценки опасностей в Китае (Чтобы обеспечить научную основу для сейсмического зондирования в Северном Китае, в 2006 году геофизикой институт китайской сейсмической администрации провел крупномасштабное наблюдение за сейсмическими массивами в Северном Китае. В этом районе было установлено 190 комплектов широкополосных сейсмографов, 10 комплектов очень широкополосных сейсмографов и 50 комплектов короткопериодных сейсмографов (рис. 2). Сейсмологический научный национальный массив в севере Китая создает две северо-западные разведывательные (измеряемые) линии: первая линий - от города Тангхая (Tang hai) через города Таншаня (Tangshan), Саньхе (Sanhe), Пекина (Beijing) и Чжанцзякоу (Zhang Jia kou) до Шангду (Shangdu), а другая от города Лингсяна (Lingxian) через города Дечжоу (Dezhou), Хэншуй (Hengshui), Синтан (Xingtang) и Шуочжоу (Shuozhou) до Циншуихэ (Qingshuihe). Расстояние между станциями двух измеряемых линией около 10 км. Кроме того, расстояние остальных станций около 35 км.).

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

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

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

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

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

Краткое содержание диссертации

Ключевыми моментами для построения математической модели сейсмического зондирования является определение местоположения источника

возбуждения и типа возникающих сейсмических волн этих источников. Возбуждение источника происходит при следующих ситуациях (рис. 3):

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

2. Если источник находится далеко от сейсмометра, принимающего сигнал, то сейсмическая волна считается сигналом поверхностной волны после длительного распространения внутри Земли.

3. Если источник расположен внутри Земли (очень медленные землетрясении или микросейсмики), то плоская волна будет возникать и попадать на границу верхних слоев Земли (Сейсмометр близко к источнику).

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

В первой главе мы описали происхождения уравнения упругой волны в слоистой среде. Решая волновое уравнение, мы можем получить параметры, зависимые от слоистой среды (плотности р, коэффициентов Ламе X, л и т. д.). Нахождение этих параметров является целью сейсмической зондирования -''наблюдение изменения параметров слоистой среды и оценка деятельности и опасности землетрясения''. В пункте 1.2 главы 1 были постановлены прямая и

8

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

В главе 2 рассматривается метод обобщенного коэффициента отражения-передачи (МОК-О/П) (The generalized reflection-transmission coefficient method) [24, 26, 39, 60] для расчета характеристики поверхностных бегущих волн слоистых сред и анализируется его недостаточность. В пункте 2.4 главы 2 разработан новый метод тензора сейсмического импеданса, позволяющий быстро вычислить параметры бегущей поверхностной волны для произвольной слоистой среды. В процессе сравнения корней дисперсионного уравнения между методом МОК-О/П и методом тензора сейсмического импеданса (МТСИ) были выявлены некоторые корни, которые не существует при решении методом МТСИ. Чтобы проверить точность корней двух методов, мы сравнили экспериментальные результаты между классическим методом матрицы перехода 'Thomson Haskell' и нашем новым методом МТСИ при заданной модели слоистой среды. Экспериментальные результаты показывают, что корни, полученные методом импеданса, совпадают с классическим методом матрицы перехода. Это также показывает, что некоторые корни, полученных из дисперсионного уравнения в методе МОК-О/П, не имеют физический смысл.

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

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

Далее изложен подробный обзор глав и пунктов диссертации.

В первой главе мы исследуем волновое уравнение и его решение для распространения упругих волн в изотропных слоистых средах. В этом процессе

описываются понятия, связанные с классической теорией упругости (напряжение, деформация и поле смещения), и их отношения, которые описываются законом Гука. Решая уравнение, видно, что в упругой среде могут распространяться две волны с разной скоростью: одна - невращательная волна с большой скоростью. Эта волна, распространяющаяся внутри Земли, часто называется продольной волной. Она также называют первичной волной (P-волны, англ. primary wave), потому что она впервые поступает на записывающую станцию, другая - волна с меньшей скоростью, которая имеет соленоидальное векторное поле. Она называет поперечной волной, также называется вторичной волной (S-волной, англ. secondary wave), так как эта волна является второй записывающейся волной. Затем используя потенциальную функцию для представления уравнений P-волны и S-волны, решение волнового уравнения было упрощено.

В результате из решения волнового уравнения следует, что S-волна может разделиться на две волны (SV волна и SH волна), у которых направления поляризации перпендикулярны для изотропной среды (направление поляризации SV волны считается перпендикулярным направлению распространения сейсмической волны в одной плоскости (рис. 4)). Они обладают одинаковой скоростью. А для анизотропной среды в направления поляризации существуют некоторые различия: скорости SV волны и SH волны разные.

Направление распространения волны _

Поляризация SH волны

^ * Поляризация Б V волны -►

Ограничено направление распространения волны в плоскости

Рисунок 4. направление поляризации различных волн.

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

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

Рисунок 5. Трехкомпонентные сейсмограммы, записанные на станции.

11

В слоистой среде энергия сейсмических волн с одной стороны образует объемные волны и распространяется наружу, а с другой образует поверхностные волны, которые распространяются вдоль параллельной границы. Более того, поверхностная волна представляет собой управляемую волну, энергия которой ограничена волноводом. Границы соответствующих слоев Земли становятся верхней и нижней границами волновода. Между свободной границей и поверхностью Мохо образован волновод. А между свободной поверхностью и границей ядерного тигля образован другой волновод. В отличие от объемных волн скорость поверхностных волн зависит от разной частоты. Это явление представляется дисперсией. Обычно информация о поверхностной волне из сейсмических записей является дисперсия фундаментальной моды.

Сейсмические поверхностные волны делятся на две категории: Рэлея и Лява. Направление поляризации волны Рэлея находится в вертикальном и горизонтальном направлениях (параллельно направлению распространения волны), направление поляризации волны Лява находится в горизонтальной плоскости и перпендикулярно направлению распространения волны (рис. 6).

Surface waves Love wave

ОтвсАол от '.тур лгсрадпсюп

Рисунок 6. Характеристики движения волны Рэлея и Лява.

Сейсмические поверхностные волны обычно обладают следующими свойствами:

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

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

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

4. Затухание поверхностной волны в горизонтальном направлении намного меньше, чем у объемной волны, в то время как в вертикальном направлении энергия быстро убывает (экспоненциально убывает). Глубина проникновения Земли изменяется с длиной волны, и чем больше длина волны, тем больше глубина проникновения. Как правило, поверхностная волна наиболее чувствительна к среде с глубиной около 1/3 длины волны. Чувствительность поверхностной волны с разными периодами к скорости поперечной волны изменяется с глубиной (рис 7).

Для решения распространения поверхностных волн в слоистых средах (Для расчета характеристики поверхностных волн слоистых сред) существует много классических методов, например, метод «Т-матрица», метод «матрица перехода» и другие численные методы. Метод «матрица перехода» обладает физическим смыслом и выведенные результаты легко реализуются с помощью ЭВМ. Этот метод был предложен Томсоном, Hanskell и Thrower [38, 75, 76] и является первой эффективным методом для вычисления фазовых скоростей и собственных значений поверхностных волн. В процессе расчета на ЭВМ для случая высокой частоты очень быстро достигаются вычислительные пределы ЭВМ и переполнение. Для решения этой проблемы исследователи, такие как авторы [45, 46], улучшили форму переходной матрицы. Чтобы изменить форму

13

матрицы переноса и сделать ее более простой, авторы [14, 78] провели дальнейший анализ и обсуждение метода матрицы перехода, изучали взаимосвязи между элементами в матрице и объясняли с точки зрения физики. На основе методов Zhang [19, 84] изменил форму матрицы переноса объяснил некоторые новые физические значения.

Рисунок 7. Чувствительность фазовой скорости волны Рэлея фундаментальной моды к скорости S-волны. Число рядом с кривой представляет период (сек.).

Кеннет и Керри [43, 44] предложили метод коэффициента отражения-передачи (МК-О/П) с другой совершенно другой точки зрения, принципиально отбросив экспоненциальный рост, который вызвал численную неустойчивость в методе матрицы перехода. С помощью метода МК-О/П полностью решил проблему потери точности и численного переполнения. Исходя из этого, исследователи разработали метод обобщенных коэффициентов отражения-передачи (МОК-О/П) Комбинируя преимущества метода МК-О/П, система формул метода МОК-О/П более кратка. Алгоритм обладает лучшей устойчивостью, чем любой другой метод в случае толстого слоя с высокой частотой. Однако, для модели слоистых сред содержащей слой с низкой скоростью трудно находить моды низкого порядка с единичной секулярной функцией (secular function). Для решения этой проблемы авторы [39]

предлагают использовать семейство секулярной функцией (secular function) вместо единичной секулярной функцией. Пей (2008) [26] предложил быстрый алгоритм обобщенного коэффициента отражения-передачи (Fast Generalize Reflection Transmission, FGRT). Этот метод упрощает шаги вычисления и значительно увеличивает скорость вычислений для много-слоистых сред с высокой частотой.

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

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

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

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

u

(x,t) = cos(a.t - kxx) + cos(a21 - k2x) (1)

^ a, + a2 k + k2

Пусть a = ^ 2 , k = , т.е.:

=a + $a, a2 =a-$a, a>

(2)

k1 = k + $k, k2 = k -Sk, k > Sk

Подставив в (1), получим:

u(x,t) = 2cos(at - kx)cos(Sat - dkx) (3)

Первый член в уравнении (3) cos (at - kx) указывает на распространение бегущей волны, т.е. фазовая скорость равна с, а последняя часть cos(Sat - Skx)

представляет собой распространение энергии, огибающей бегущие волны, т.е. групповая скорость равна U, которая приводит к эффектам амплитудной модуляции на первом члене. Поскольку ^и^ намного меньше a и k, по сравнению с групповой скоростью фазовые скорости быстро изменяются как по частоте и по пространству. Из рисунки 8 можно видеть зависимость между фазовой скоростью (соответствующей вершинам и подошвам волны) и групповой скоростью (соответствующей огибающей и энергии). Из уравнения (3) можно получить:

Sa d (ck) dc

U = _®=^—'- = c + k— (4)

Sk dk dk

Если 8а ^ 0,8к ^ 0, то имеем:

а ? к

и _8а _ da _ da _ da _ с _ с _ 8к ~ dk ~ d (ф с)_ ~ 1 _Ф dc_ ~ 1 +

dc С с da с dT

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

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

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

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

(ah cos{шit — /¿q*) ec»i — + Sai, A, — к + &k

cos (tttjf - frjjf) й»г = I» - ¿¡л», k£=k — 8k

Envelope Carrier

Рисунок 8. Понятия фазовой скорости и групповой скорости [74]. (а) две гармонических волн, обладающие аналогичные частоты и волновые чисел. (б). Наложение двух простых гармонических волн при ? = t2, . Толстая линия -огибающая линия волны после наложения. Скорость распространения является групповой скоростью и. Внутренние бегущие волны с амплитудной модуляцией распространяются с фазовой скоростью С.

В литературе [66] введена техника преобразования Фурье в измерение дисперсионных кривых. Этот метод соответствует сигналам шумов (микросейсмики) и сигналам с быстрым изменением. Alexander [16] впервые применил цифровые методы фильтрации к измерениям дисперсии поверхностных волн (в литературе [65] был разработан метод фильтрации переменный по времени (Time-varying filtering method) для измерения фазовая скорость). Применение техники цифровой фильтрации и метода фильтрации переменного по времени имеет эпохальное значение при измерении дисперсии поверхностных волн. Для анализа измеренных сейсмических данных многие методы, были разработаны на основе быстрого преобразования Фурье и численной фильтрации. Landiman et al. (1969) [47] предлагает метод перемещение окон анализ (Moving Windows Analysis). Использование методов фильтрации улучшает отношение сигнал-шум данных, расширяет диапазон периодов измерений дисперсии и улучшает разрешение. Dziewonski [28] предложили технику множественных фильтров (ТМФ) (Multiple Filter Technique). Узкая полосовая фильтрация (Narrow band pass filter) эффективно улучшает отношение сигнал-шум за определенный период и улучшает точность измерения дисперсии групповой скорости. Исходя из этой техники, был создан метод анализа частоты-времени поверхностных волн (The frequency-time analysis method, FTAN [89]). Можно показать, что если функция Гаусса выбрана как оконная функция, метод анализа перемещения окон (Moving Windows Analysis) и метод множественных фильтров (Multiple Filter Technique) являются эквивалентными. Dziewonski [29] доказывает, что метод множественного фильтра выполняется только в приближении первого порядка. В целях повышения точности измерений предлагается метод измерения дисперсии остаточной ошибки (Residua Dispersion Measure). Этот метод не только соответствует измерению групповой скорости поверхностных волн, но и применяется к методу фильтрации временных переменных (Time variable Filtering) и может значительно улучшать вычислительную эффективность. Поскольку разрешение частот-времени ТМФ связано с выбором коэффициентов гаусса, Cara [23] предложил оптимальную схему фильтрации. Существует принцип неопределенности измерения в анализе времени. Разрешение частотной области и временной области является взаимоисключающим. Чтобы поддерживать лучшее временное разрешение в течение более длительных периодов, Nyman & Landisman [58] предложил корректирующий фильтр дисплея (The display equalized filter). Herrin &

18

Goforth[40] предложили технику адаптивного фильтра (ТАФ) для выделения поверхностных волн различных форм колебаний, устраняя эффекты множественных путей и помех более высокого порядка. В настоящее время ТМФ и ТАФ широко используются для обработки сейсмических данных и разработано относительно зрелое программное обеспечение [83].

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

1. Существование ошибки для определения местоположения источника.

2. Существование ошибки во время землетрясения (earthquake origin time).

3. Механизм источника (earthquake mechanism) и начальная фаза неизвестны.

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

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

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

Можно считать, что эти данные, измеряемые из второй станции, являются возбужденными сигналами (нулевой фазовый сдвиг и сигнал 5-импульса) из первой станции. Поэтому фазовая информация взаимно-корреляционной функции может быть использована для расчета фазовой скорости распространения поверхностных волн между двумя станциями. Метод узкополосной фильтрации-взаимной корреляции [20] является широко используемым методом для измерения фазовой скорости между двумя станциями. С помощью методов взаимно-корреляции между двумя станциями значительные результаты были получены при инвертировании стратиграфических структур земной коры [49, 50, 51].

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

1. Предварительно обработка сейсмических данных двух выбранных из двух станций.

2. Отсекания относительно устойчивых соответствующих данных.

3. Нормализация выбранных данных в области времени.

4. Данные, соответствующие времени, обрабатываются двумя интегралами.

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

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

Учитывая прямую задачу для расчета дисперсии поверхностных волн и измеренную дисперсионную кривую (или теоретически рассчитанную дисперсионную кривую), мы можем получить параметры модели слоистой среды (коэффициенты Ламе А, д, и плотность р)

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

Методы инверсии поверхностных волн в основном включают метод наблюдения, метод локальной линейной оптимизации [79] и нелинейный алгоритм глобальной оптимизации [22, 33, 64, 70, 71, 73]. В методе исследования поверхностных волн на основе природного источника измеренный сигнал поверхностной волны может получить дисперсионную кривую, получив скорость распространения на каждой частоте с помощью различных методов разделения поверхностных волн. На основе измеренных

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

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

В соответствии с теорией полуволны этот метод учитывает длину волны 1/2 как глубину обнаружения и вычисляет соотношение между скоростью поперечной волны У'' и фазоывой скоростью Ут в изотропной слоистой среде.

\ Да

Вывод результаты

Рисунок 9. Основные идеи инверсии.

При разработке методов инверсии возникли локальные линейные методы (например, наименьшие квадраты, метод Оккама и т. д.). Однако, если существует большая разница между исходной моделью и измеренной моделью, локальный метод оптимизации не сможет вычислить глобальное оптимальное решение. Это значит, что метод линеаризации не только сильно зависит от выбора исходной модели, но еще требует вычисления частной производной, на которую непосредственно влияет точность матрицы Якоби. Для преодоления этого дефекта был разработан алгоритм нелинейной глобальной оптимизации. В методе нелинейной глобальной оптимизации (например, алгоритм имитации отжига [25, 64], искусственная нейронная сеть [33, 71], метод глобальной оптимизации роя частиц [22, 73], генетический алгоритм (ГА) [37, 55, 56, 63, 70], и т.д.) не нужно вычислять частную производную, этот метод широко применяется в инверсии кривой дисперсии сейсмических поверхностных волн.

В области геофизических исследований, будь то гравитационное поле, магнитное поле, электромагнитное поле или поле сейсмической волны, все можно отнести к операторному уравнению первого рода (operator equation of the first kind). Целью решения обратной задачи является получение оценки параметров модели инверсии путем решения первого типа операторного уравнения. Характеристика обратной задачи является некорректной задачей: решение не обязательно существует; непрерывность или устойчивость решения не могут быть гарантированы. Обратные и некорректные задачи были значительно продвинуты с момента проведения Адамара (Hadamard). Многие ученые работали над теорией некорректных задач и методов их решения, добились плодотворных результатов и предложили ряд методов для решения таких задач (например, метод отбора, метод квазирешения, метод сглаживания и метод регуляризации, и др. [13]). Для решения некорректных задач различные методы и примеры их применения подробно описаны в литературе [11]. Среди них метод регуляризации сыграл важную роль в различных междисциплинарных научных областях.

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

singular value decomposition, TSVD), метод дискретной регуляризации, методы итеративная регуляризация и методы фильтрации и т. д.

Метод регуляризации Тихонова сыграл важную роль в исследовании обратной задачи [2, 11]. При решении практических проблем он может не только преодолеть некорректность, но и значительно уменьшить распространение погрешности и шумов. В пункте 4.2 главы 4 рассматривается новый метод «генетический алгоритм, связанный с мет-дом регуляризации Тихонова (ГАМРТ)» для исследования обратной задачи дисперсии поверхностных волн. Учитывая преимущества двух методов (ГА и Метод регуляризации Тихонова) и сложность обратных задач в теории поверхностных волн в слоистых средах, расчет инверсии проводился с использованием дисперсионной кривой, рассчитанной теоретической моделью слоистой среды для доказательства эффективности метода ГАМРТ. Далее мы добавляем относительную погрешность к теоретической дисперсии для дальнейшей проверки эффективности метода. В конечном эксперименте мы рассчитали структуру скоростей земной коры и верхней мантии в некотором районе Китая, используя дисперсионные кривые, полученные с помощью метода выделения поверхностных волн из сейсмических данных, описанного в главе 3, и достигли хороших результатов.

Цели диссертационной работы:

1. Разработка математических моделей для описания прямой и обратной задач сейсмического зондирования (дисперсия сейсмической поверхностной волны), базирующегося на характеристиках распространения сейсмического уравнения и основном уравнении Ламе для слоистых сред.

2. Проведение численного расчета связи с данными моделями слоистых сред и установленным прямым алгоритмом и сравнение полученных результатов со стандартными методами.

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

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

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

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

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

В диссертации получены следующие основные результаты,

выносимые на защиту:

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

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

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

4. Разработанный программный комплекс был опробован для определения состояния геологических сред в районе Куньмин Китая (станции КМ1) с использованием данных от 13 удаленных землетрясений. Полученные результаты показали эффективность разработанных методов, позволяющих устойчиво определять состояния геологической среды.

Структура и объём диссертационной работы

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

списка литературы. Полный объем диссертации составляет 131 страниц. Список

литературы в диссертации включает 96 наименований.

Рекомендованный список диссертаций по специальности «Математическое моделирование, численные методы и комплексы программ», 05.13.18 шифр ВАК

Введение диссертации (часть автореферата) на тему «Обратная задача сейсмического зондирования с использованием поверхностных волн»

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

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

1. Конференция «Annual Meeting Chinese Geoscience Union», доклад: Research on the mathematical model of microseismic monitoring (in Chinese). (Beijing, China, 2015)

2. Международный российско-китайский учебно-научный семинар МГУ-ППИ по прикладной математике и информатике. Доклад: Методы обработки сейсмических данных. (Москва 2018)

3. Конференция «Annual Meeting Chinese Geoscience Union» (in Chinese). (Beijing, China, 2018)

4. Ежегодная научная конференция «Тихоновские чтения», Москва, 2018

Материалы диссертационной работы опубликованы в 4 статьях в журналах из списка scopus [77, 80, 81, 82] и 1 сборниках [93], 3 тезисах докладов конференций [90-92].

Автор выражает глубокую благодарность своему научному руководителю Владимиру Ивановичу Дмитриеву за помощь в постановке содержательной задачи и постоянное внимание к работе. Автор также выражает благодарность Ингтем Женни Гастоновна за редактирования диссертации.

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

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

1.1 Уравнение распространения упругих волн в слоистой среде

Рассмотрим упругую изотропную среду, которая занимает в момент времени ^ область V трехмерного евклидова пространства с достаточно-

гладкой границей £ ( V = V и £ , V е Е3 )и характеризуется параметрами р=р(г), Л = Л(г), и = и(г), выполняющими условие

р> 0, и> 0, 3Л + 2и> 0, где р - плотность среды, Л, и - коэффициенты Ламе.

Пусть в момент времени t > t0 внешние силы действуют на упругую среду в область V' е Е3. Тогда основная задача классической теории упругости будет заключаться в определении состояния упругой среды в заданных точках области V' в момент времени t, то есть определяется вектор смещения в Декартовых координат и = {их, иу, иг} , тензор деформаций Е, тензор

напряжения Т = {ту }, г, у = 1,2,3 (рис.1.1).

Рисунок 1.1. Напряжения, действующие на грани элементарного куба.

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

3

Ь = Л8Xекк + 2леу, •/ =1,2,3 (111)

к=1

где - символ Кронекера. Соотношение между деформацией и смещением имеет вид:

дих ди дих

_ _ Л | / ; е — л ф

-л XX -л

ду дх дх

Я = 1 {е е = е = -^ ^• е ^• (1 1 2)

^ 2*' еУ2 ^ дх ду ' еу ду ' ( )

ди ди ди„

е = е =—- +---; е

д2 дх д2

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

поверхности (достаточно-гладкая граница объема V) 5, равна .

Общая массовая сила, действующая на данный объем V, равна |Ц fdV .

V

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

д 2и

фь Ш + = 0 (113)

5 V Vй1

(1.1.3) можно переписать следующим образом:

/Ь • ^-Я/А*=0 (114)

5 V Vй1

С помощью интегральной формулы Гауссы поверхностную силу можно представить в виде:

<5^)0-пс® = |Л div0dV (1.1.5)

5 V

переписываем:

Ц^-п^(116)

£ V дХ1

Подставив (1.1.6) в (1.1.4) получим:

Ш

Г доц п д 2м Л

+ Л-р

2

дх дt

V 1

СУ = 0 (1.1.7)

V

Так как точка области V выбирается произвольно, то выражения в скобках должно быть равны нулю, т.е.:

0 + /г-р%=0 (1.1.8)

дху 1 д^

или

д 2м

рд^=0 + л (119)

Уравнение (1.1.9) называется упругим динамическим уравнением. Для решения уравнения (1.1.9) нам нужно знать связь между деформацией и смещением (1.1.2) и соответствующие граничные условия. Для изотропной среды перепишем закон Гука в виде:

01 = Щии * + и(м, ,1 + и1,) (1.1.10)

Подставив (1.1.10) в (1.1.9) получим уравнение в изотропных средах для смещений

д 2«

р~Им=м* *д Л+Ли*м+(щ ,1 + и,, )д и+и(+ и, *)+Л (ил1)

Это уравнение можно представить в векторной форме:

д 2и

р= ^ЛУ-и + (Л + и)УУ-и + Уи (Ум + мУ) + иУ2м + / (1.1.12)

где мУ - транспорт для Ум, если компонент Ум равен и],, то компонент мУ равен и11. Тогда согласно свойству тензора деформации Ум + иУ=2Е, где

Е - тензор деформации. С помощью выражения УхУх и = УУ • и - У и уравнение (1.1.12)можно представить в виде:

д 2и

р= УЛУ • и +(Л + 2и)УУ • и +Уи (Ум + мУ) + иУхУхм + / (1.1.13)

Уравнение (1.1.13) является общей формой упругого динамического уравнения в изотропной среде. Для задачи возбуждения сейсмической волны f -массовая сила, зависящая от сейсмического источника. Если мы рассмотрим задачу свободного колебания Земли, массовая сила f должна включать гравитационное поле Земли. При отсутствии массовой силы уравнение (1.1.13) имеет вид:

д 2и

р^г = +(Я + 2л)УУ-и + Ул-(Уи + uV)-¡^VхVхu (1.1.14)

Уравнение (1.1.14) называется уравнением сейсмической волны вне зоны действия источника. Для решения уравнения (1.1.14) имеется следующие способы:

1. В случае однородной среда существуют аналитическое решения.

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

3. Численные методы.

Задача решения уравнения (1.1.14) является самым главным научным направлением в теоретической сейсмологии, в том числе исследование характеристики сейсмических фаз, вычисление пути сейсмических волн, расчет синтетических сейсмограмм и т.д. Для латеральной неоднородной среды решение (1.1.14) становится очень трудным из-за существования градиента параметра среды, и в настоящее время трудно получить аналитическое решение. Для горизонтально-слоистых одномерных сред градиентный слой среды можно разделить на множество тонких слоев. Каждый слой можно считает однородным. Когда толщина слоя намного меньше длины доминирующей волны, приближенное решение для (1.1.14) может обеспечить достаточную точность. В настоящее время существуют хорошие результаты исследования задачи возбуждения и распространения сейсмических волн в однородных горизонтально-слоистых средах. В однородной изотропной среде уравнение (1.1.14) можно упростить следующим образом д2и (х ?)

р-= (Л + 2л)УУ- и (х, ?) - ¡¡V х V х и (х, ?) (1.1.15)

Для расчета характеристики сейсмических волн в однородной изотропной среде, нужно переписать уравнение (1.1.15) в спектральной форме через преобразования Фурье

оо

и(х,с) = | и (х,г)вшЖ (1.1.16)

—<х>

Подставив (1.1.16) в (1.1.15) получим:

(Л + 2^)УУ-и (х,ю) —^УхУхи (х,а) + рсо2и (х,с) = 0 (1.1.17)

Переобозначив в (1.1.17) через а = \ р= ^ получим:

V Р \Р

а2УУ- и (х,с) — в2УхУхи (х,с) + с2и (х,с) = 0 (1.1.18)

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

и (х,с) = иа(х,о) + ир(х,о) (1.1.19)

и выполняет дополнительное условие:

Ухиа= 0, У-ир= 0 (1.1.20)

Если для иа и ир , выполняются уравнения (1.1.18)-(1.1.20), то они

независимы.

Если иа зависит от ир , то существуют отношение иа = сир , где с —

ненулевая постоянная. Подставляя это соотношение в (1.1.20) получаем:

У - иа= сУ - ир = 0

У хив =1У хиа = 0 с

Тогда Ухиа=Ухир=У-иа=У-ир = 0 , из уравнения (1.1.19) узнаем У-и = 0, Ухи = 0 , подставив выражение в (1.1.18), получим и(х,с) = 0 . Полученный результат противоречит предположению, что и (х,с) не равно нулю, следовательно иа и ир независимы.

Подставив уравнение (1.1.19) в (1.1.18) и согласно условию (1.1.20) получим:

а2УУ-иа — р2УхУхир + с2 (иа+ир) = 0 (1.1.21)

Подставив уравнение (1.1.19) в формуле УхУхи = У У - и — У и , получим:

УхУх( иа+ир) =УУ - (и а + и р ) У2 (иа+ир) (1.1.22)

Согласно условию (1.1.20) уравнение (1.1.22) можно упростить в виде:

1р+У2ир=УУ- иа—У\

У х У х и В+У2и в=УУ - иа—У2иа (1.1.23)

Так как иа и ир независимы, то левая и правая часть уравнения (1.1.23) должны быть равны нулю, т.е.

УхУхи в = -У2и^

УУ- и а = У2иа

(1.1.24)

Подставив (1.1.24) в (1.1.22), получим:

(а2У 2иа + (о2иа) + (Р2У 2ир + (о2ир) = 0 (1.1.25)

Согласно выражения (1.1.25) в скобках должны быть равны нулю

а2У 2иа + о)2иа = 0

(1.1.26)

в У ив + —ив = 0

Пусть ка= — , кв= — , уравнение (1.1.26) имеет вид: а р

У 2 иа+ к2а иа= 0 (1.1.27)

У2ир + крир = 0 (1.1.28)

Где а - скорость продольной волны, в - скорость поперечной волны, ка -волновое число продольной волны, кр - волновое число поперечной волны.

Таким образом, для однородной изотропной среды задачу решения уравнения сейсмической волны (1.1.18) можно свести к задаче решения уравнения Гельмгольца векторной формы (1.1.27) - (1.1.28). Уравнения (1.1.27) и (1.1.28) являются волновыми уравнениями по частоте. В однородной изотропной среде поперечные и продольные волны существуют независимо и

распространяются с различными скоростями в и а , причем а . В

жидкости не существует поперечной волны, так как модуль сдвига равен нулю О = 0 (Скорости продольной волны а и поперечной в зависят от модуля

сдвига: а = ^2О(1 - у)/р(1 - 2v) , в = ^О/р , где О - модуль сдвига, V -коэффициент Пуассона, р - плотность среды). Так как У х иа = 0 , то продольные волны распространяются без вращающегося поля, а поперечные волны создают соленоидальное векторное поле, т.е. У-ир= 0 . При этом

характер поля смещения продольной и поперечной волны можно представить в виде:

иа = Уф? ив = Ух^, У- ^ = 0 (1.1.29)

где р - скалярный потенциал для продольной волны, ^ - векторный потенциал для поперечной волны. Поле смещения и имеет вид:

и = Ур + Ух^, У- ¥ = 0 (1.1.30)

Введем скалярные потенциалы фиф, удовлетворяющие скалярному уравнению Гельмгольца, т.е.:

У 2р + к2ар = 0 (1.1.31)

У2у + к2ру = 0 (1.1.32)

Определяются поля смещений продольной и поперечной волн в виде:

иа= ~Т УР

ка

(1.1.33)

ивн =Ух( ау)

(1.1.34)

и

= — У хУ х( ау)

к

(1.1.35)

где вектор с отличается в разных системах координат. В декартовых и цилиндрических координатах а = 2, в сферической системе координат а = г. Поперечные волны определяют два решения (1.1.34) и (1.1.35). Правые части уравнения (1.1.33)-(1.1,35) можно представить в виде:

Ь = —Ур, М =Ух(ау), N = —УхУх(ау)

ка кв

(1.1.36)

Векторы Ь, М, N называется векторами Хансена. Когда скалярные потенциалы ф и ф удовлетворяют скалярному уравнению Гельмгольца, вектор Хансена удовлетворяет векторному уравнению Гельмгольца. Тогда имеем следующие отношения:

1

N = —УхМ, У-Ь = -кр( к

(1.1.37)

в

Здесь векторы М, N, Ь независимы друг от друга. Векторы ХансенаМ, N, Ь являются линейно-независимыми векторами, удовлетворяющие векторному уравнению Гельмгольца, поэтому решение уравнения сейсмических волн можно представить в виде:

и = АЬ + ВМ + CN (1.1.38)

где А, В, С - постоянные зависимых от частоты.

Если волновой фронт сейсмической волны является плоскостью, то волна называется плоской волной. Гармоническая плоская волна во времени можно представить в виде:

и (г ) = а (—) е ~—г+'гх (1.1.39)

где у = — - характеристика плоской волны, с - фазовая скорость

с

распространения плоской волны.

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

В декартовой системе координат решение скалярного уравнения Гельмгольца (1.1.31) (1.1.32) может быть получено с помощью метода переменных разделения. Пусть (1.1.32) имеет решение по форме:

(р = X (х)У (у) 7 (I) (1.1.40)

Подставив (1.1.40) в (1.1.32), с помощью метода переменных разделения получим:

X" + а X = 0

у"+к2у = о

7" + к2аз7 = 0

(1.1.41)

где ка1, ка2, ка3 - постоянные независимых от пространственных координат и выполняют условие:

ка=а+ка2+ка (1.1.42)

Решение системы(1.1.41) имеет вид:

X « е±ках, У « е±к*2У, 7 « е±1а, (1.1.43)

Скалярный потенциал ф можно представить в виде:

ф = Ае'( ка1х1 + ка 2 х2 +ка3х3 ) = Де'ка'х (1 1 44)

где ка=—п = ( ка1, ка2, ка3)- волновое число продольной волны, А -а

постоянный независимого от пространственных координат. Подставив (1.1.44) в (1.1.33) получим:

= А—) каеКх (1.1.45)

иа

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

В декартовых координатах согласно (1.1.41) -(1.1.46) с помощью метода разделения переменных решение скалярного уравнения Гельмгольца поперечных волн имеет:

Н = Ве(кв1Х1 +кв2Х2+кв3Х3) (1146)

или

Н = ВеквХ (1.1.47)

где кв=вП = ( кв1, кв2, кв3)- волновое число поперечной волны, В -

постоянный вектор, независимый от пространственных координат. Подставив (1.1.47) в (1.1.34), получим:

ивн = В (ю)(кв2,-кв,0) е&в (1.1.48)

где В(ю)-комплексные амплитуды независимые от координат. Известно, что

направление поляризации волны SH находится в горизонтальном направлении и перпендикулярно направлению распространения волны [96]. Подставив (1.1.47) в (1.1.36), получим другую плоскую поперечную волну, называемая БУ-волной

ив = С (ю)(-кв1кв3, - кв2кв3, к2в! + к2в2) е^ (1.1.49)

где С(ю)- комплексные амплитуды независимые от координаты. Известно,

что направление поляризации БУ-волны перпендикулярно направлению распространения сейсмической волны и поляризации БН-волны.

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

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

Если направление распространения плоской волны выбрано в плоскости х1 - х3, то в этом случае выражение трех плоских волн представляет в виде:

= А—)(каХ, 0, ка3)е'к"х (1.1.50)

иа

и8рн = В—)( 0,1,0 ) еквх (1.1.51)

ив = С—)(-кв3,0, кр1) еЛвХ (1.1.52)

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

где q - коэффициент пропорциональности. Легко понять, что поляризация БИ-волны находится в направлении х2, направление поляризации волны SV находится в плоскости х1 - х3 и перпендикулярно направлению поляризации

волны Р, а направление поляризации продольной волны параллельно направлению распространения сейсмических волн [95].

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

ка2 = ка3 = кР2 = кР3 = 0 П = х1, т.е.:

иа

ив

ив

= А—) хек«х (1.1.53)

•5н = В—) х2е'квх (1.1.54)

■8Г = С—) х3екрх (1.1.55)

Рассматривается изотропная идеально-упругая плоско-слоистая среда, в которой характеристики Ламе Я, д, и плотность р являются произвольными кусочно-непрерывными функциями координаты 2. Слои не ограничены по простиранию, число слоев конечно и равно N, снизу модель ограничена полупространством г > гп . Полупространство г > гп является упругой средой, в которой характеристики Я, д, р положим постоянными. Построенную модель среды будем называть кусочно-градиентной (рис. 1.2).

Для j-ого слоя волновое уравнение имеет вид:

Р

(= »+ ^уу^и(^ х,г )-М( ;)УхУхм(')(х,г) + / (х,г) (1.1.56)

Данное уравнение удовлетворяет следующим граничным условиям:

1. На земной поверхность тензор напряжения равен нулю

(х>' )|г=0 = 0

2. На внутренней границе смещение и напряженности непрерывны

и

(х, г) (х, г)

и

(/+1)

г=г/

=т+1)

(х, г) (х, г)

г=г(/)

г=г (/)

(1.1.57)

(1.1.58)

(1.1.59)

3. В нижнем полупространстве не имеет отраженную (уходящую) волну.

Свободная поверхность

„ (Ч .

z(1) Г«

лю

;(1) //(1) п® ^ , М ■ Р

Ли\ци\ри)

1 источник

— —

X

Рисунок 1.2. Модель слоистой среды. Волновое уравнение по частоту имеет вид:

-р(/ )а2и{])(х,т) = (Я(/) + 2//(/))УУ• и(/ )(х, т) -^)УхVхи(/ )(х,т) + /(х,т) (1.1.60) Тогда граничные условия перепишем в виде:

1. На земной поверхности тензор напряжения равен нулю

Т (х—) ,=0 = 0

(1.1.61)

2. На внутренней границе смещение и напряженности непрерывны

и (у)( х,— = и (7+1)( х, — х, — ( = Т+1)( х, —

г=гу '

(1.1.62)

(1.1.63)

3. В нижнем полупространстве не имеет отраженную (уходящую) волну.

1.2 Постановка задачи микросейсмического зондирования

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

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

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

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

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

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

Рисунок 1.3. Возможные источники измеряемых данных.

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

Пусть имеется изотропная слоистая среда при г е[0, Н ] (г = 0 - земная

поверхность), упругие параметры который зависят только от глубины z . При г > Н имеем однородную среду из которой на границу г = Н падает плоская волна. Уравнения распространения упругих колебаний в однородном случае имеют вид:

Р- =1 [«- >£

'( г ^ = 4 ы г > + г

(1.2.1) (1.2.2)

где p( z)- плотность среды, /л(z) и Л( z) -коэффициенты Ламе для упругой среды. Обычно рассматривают преобразование Фурье:

<х>

ul(z,a>) = J ul(z,t)eiatdt, l = (x,z) (1-2.3)

—<X>

Для которого уравнения (1.2.1) и (1.2.2) можно записать в общем виде:

—|ф(z)du^\ + Юp(z)u(ю) = 0, z е[0,H] (1.2.4)

dz ^ dz ^

где для поперечных волн ф(z) = ju(z), u(ю) = ux (ю), для продольных волн ф(z) = Л(z) + 2^(z) , u(ю) = uz (ю) . на земной поверхности выполняется условие:

dU (z = 0) = 0 (1.2.5)

dz

На нижней границе ( z = H ) выполняется условие возбуждения поля плоской волной, а в точках z = zn разрыва коэффициента ф(z) должны выполняться условия сопряжения:

U ( zn + 0 ) = U ( zn — 0 ) ф n + 0)du(zn + 0) = ф(zn — 0)du(zn — 0) (1.2.6)

dz dz

Выведем условия сопряжения поля плоской волной. При z > H имеем однородную среду ф = фн = const, p = pH = const . Тогда уравнение (1.2.4) имеем вид:

^ + kHu = 0, кн = alp (1.2.7)

Тогда при z > H имеем решение в виде:

u (z) = Ae~ikH (z—H) + BeikH (z—H) (1.2.8)

где первый часть падающая волна, а второй часть отраженная волна. Вычислив производную u(z) и исключив отраженную волну, получим:

u'( z) — ikHu (z ) = —2 ikHAe~ikH (z—H), при z > H (1.2.9)

Тогда при z = H + 0 имеем вид:

u'( H + 0) — ikHu (H + 0) = —2 ikHA (1.2.10)

Учитывая условия сопряжения (1.2.6), получим окончательное условие возбуждения при г = Н - 0:

ф(г = Н - 0> и'(Н - 0> - ¡кНфНи (Н + 0> = -2¡кНфНЛ (1.2.11)

Таким образом мы получим формулировку прямой задачи (1.2.4) -(1.2.11), которая позволяет определить и (г> при г е [0, Н ]. Следовательно, мы можем вычислить частотные характеристики поля смещения на земной поверхности:

и (г = 0> = и0 (с,р( г> ,ф( г>> (1.2.12)

и0 -нелинейный оператор зависящий от частоты с и действующий на р(г>,ф(г> .Обратная задача сейсмического зондирования заключается в

определении р(г>,ф(г> по измеренному на земной поверхности и(т^(с) в зависимости от частоты. Рассмотрим метод расчета и (г = 0 >. Заметим, что и (г >, в силу линейности задачи, с помощью метода переменных разделения можно представить в виде:

и (г> = и0 (с>у(г>, где ио (с) = и (г = 0) (1.2.13)

Функция V(г>, согласна (1.2.4) -(1.2.6), является решением следующей задачи Коши:

(ф( г > v'( г ))> +с2р( г > V ( г > = 0

V(г = 0> = 1; V (г = 0> = 0 (1.2.14)

(г> и ф(г>V '(г>-непрерывны

Полученная задача легко решается численно, а для кусочно-постоянных коэффициентов р(г>,ф(г> из системы (1.2.14) можно получить аналитическое

решение. Пусть коэффициенты уравнения заданы следующим образом:

р(г) = рп; ф(г)=фп при ге[гп-1>гп],пе[1,Щ,г0 = 0

/ ч / ч (1.2.15)

р(г = Н) = рн; ф(г = Н) =Фн при х>хм = Н

Тогда при г е [гп-1, гп ] решение задачи (1.2.14) имеем вид:

v(г) = V! СОБкп (г - гп-1) + кп (г - гп-1) (1.2Л6)

кп

гдеVn-! = V(гп-1), = гп_{ + 0), кп = . Продифференцировав (1.2.16)

и положив г = гп = гп-1 + кп, получим рекуррентную формулу перечета:

vn = V1cos knhn + knhn

К

V' ( ^n - 0) = - Vl sin knhn + V'-l cos knhn

Согласно условия (1.2.4) -(1.2.6) имеем вид:

ФгУ' ( ^n - 0 ) = Фп+/ ( ^n )

v ( Zn - 0 ) = ф v ( zn ) = ф v'n Kn ) фп Kn) фп n

(1.2.17)

(1.2.18)

(1.2.19)

Подставив (1.2.19) в (1.2.18), получим:

кф

v = - v ,

n n-1

sin ±

v' , cos knhn

n-1 n n

(1.2.20)

/ п п 1 п—1 Н Н V

ФН+1 ФН+1

Зная у0 = 1, у'0 = 0 , по формулам (1.2.18) и (1.2.20) находим последовательно у1, у{ до ,. Зная отношения = V(г = Н + 0) и = V '(г = Н + 0) ,

согласно (1.2.13), можно определить,

и (г = Н + 0) = и0 (ю)уу, и (г = Н + 0) = и0(т)^ (1.2.21)

Подставив полученные выражения в (1.2.20), найдем:

(1.2.22)

U0 (ю) = -

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

Прежде чем рассмотреть общий случай решение обратной задачи. Рассмотрим поведение и0 (т) для простейших моделей. Пусть данной слой имеет параметры фх, р1,^ , который находится на полупространстве с параметрами фН, рН , .Тогда согласно формулы (1.2.18) и (1.2.20) имеем:

v1 = cos k1h1; v| = - sin k1h1

фн

Подставив эти выражения в (1.2.22), найдем:

2ikHфHA

Uо М = -

или

Uо М =

ikнфн cos k1h1 ± k1ф1 sin k1h1

2 A(cos k1h1 ± i61H sin k1h1) cos2 k1h1 ± 61H sin2 k1h1

(1.2.23)

(1.2.24)

(1.2.25)

где вш = 1ж

\РнФи

Из уравнения (1.2.25) следует:

, ч ImU0(о) Л 7 , Q(o) = тт°) = вш tan kh ReU0 (о)

(1.2.26)

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

вш Р hi

ф <чФ 1 (1.2.27)

рНфН \ф

Рассмотрим два слоя с параметрами ф1, р1, к1 и ф2, р2, , на однородном

полупространстве имеет параметры фН,рН. Тогда, зная v!,V' согласно (1.2.23) вычислим по формулам (1.2.18) и (1.2.20):

где в12 =

i

Р1Ф1

Р2Ф2

v2 = cos k1h1 cos k2h2 - 6l2 sin k1h1 sin k2h2 . Продифференцировав выражение, получим:

(1.2.28)

v

к2ф2

Ф

(sin k2h2 cos k1h1 + 6l2 cos k2h2 sin k1h1)

н

Подставив (1.2.28) и (1.2.29) в (1.2.22), найдем:

U (о)= 2A(-12 + Í@2H^ 12 )

U 0 V»)- -2 + a2 2

12 2н 12

(1.2.29)

(1.2.30)

где

-12 = cos k1h1 cosk2h2 - 012 sin k1h1 sin k2h2 gX2 = cosk1h1 sin k2h2 + 012 sin k1h1 cosk2h2

a

Р2Ф2

2 H

Рнф

H

легко видно, что

Q (о) =

Im U0 (°)_a2н^ 12

ReU0 (о) -12 Тогда для n-ого слоя, Q(о) можно приставить в виде:

(1.2.31)

QM = = (1.2.32)

где

с = £ sin kh + 0 с cos kh

' mn = pm n n mn^pm n n

£ = £ cos kh -0 с sin kh

mn pm n n mn pm n n

m = n -1,p = n - 2,n = 1..N -1 Это означает, что зная Q(ю), можно определить следующие величины: knhn, k

0 H ,0 , — или следующие относительные величины:

К

hn Рп Фп

К' Рт

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

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

проверим осуществимость подхода, хотя наш метод установлен на теоретическом выводе.

Похожие диссертационные работы по специальности «Математическое моделирование, численные методы и комплексы программ», 05.13.18 шифр ВАК

Список литературы диссертационного исследования кандидат наук Ян Цзяньсюнь, 2019 год

1 источник

— —

X

Рисунок 1.2. Модель слоистой среды. Волновое уравнение по частоту имеет вид:

-р(/ )а2и{])(х,т) = (Я(/) + 2//(/))УУ• и(/ )(х, т) -^)УхVхи(/ )(х,т) + /(х,т) (1.1.60) Тогда граничные условия перепишем в виде:

1. На земной поверхности тензор напряжения равен нулю

Т (х—) ,=0 = 0

(1.1.61)

2. На внутренней границе смещение и напряженности непрерывны

и (у)( х,— = и (7+1)( х, — х, — ( = Т+1)( х, —

г=гу '

(1.1.62)

(1.1.63)

3. В нижнем полупространстве не имеет отраженную (уходящую) волну.

1.2 Постановка задачи микросейсмического зондирования

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

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

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

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

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

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

Рисунок 1.3. Возможные источники измеряемых данных.

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

Пусть имеется изотропная слоистая среда при г е[0, Н ] (г = 0 - земная

поверхность), упругие параметры который зависят только от глубины z . При г > Н имеем однородную среду из которой на границу г = Н падает плоская волна. Уравнения распространения упругих колебаний в однородном случае имеют вид:

Р- =1 [«- >£

'( г ^ = 4 ы г > + г

(1.2.1) (1.2.2)

где p( z)- плотность среды, /л(z) и Л( z) -коэффициенты Ламе для упругой среды. Обычно рассматривают преобразование Фурье:

<х>

ul(z,a>) = J ul(z,t)eiatdt, l = (x,z) (1-2.3)

—<X>

Для которого уравнения (1.2.1) и (1.2.2) можно записать в общем виде:

—|ф(z)du^\ + Юp(z)u(ю) = 0, z е[0,H] (1.2.4)

dz ^ dz ^

где для поперечных волн ф(z) = ju(z), u(ю) = ux (ю), для продольных волн ф(z) = Л(z) + 2^(z) , u(ю) = uz (ю) . на земной поверхности выполняется условие:

dU (z = 0) = 0 (1.2.5)

dz

На нижней границе ( z = H ) выполняется условие возбуждения поля плоской волной, а в точках z = zn разрыва коэффициента ф(z) должны выполняться условия сопряжения:

U ( zn + 0 ) = U ( zn — 0 ) ф n + 0)du(zn + 0) = ф(zn — 0)du(zn — 0) (1.2.6)

dz dz

Выведем условия сопряжения поля плоской волной. При z > H имеем однородную среду ф = фн = const, p = pH = const . Тогда уравнение (1.2.4) имеем вид:

^ + kHu = 0, кн = alp (1.2.7)

Тогда при z > H имеем решение в виде:

u (z) = Ae~ikH (z—H) + BeikH (z—H) (1.2.8)

где первый часть падающая волна, а второй часть отраженная волна. Вычислив производную u(z) и исключив отраженную волну, получим:

u'( z) — ikHu (z ) = —2 ikHAe~ikH (z—H), при z > H (1.2.9)

Тогда при z = H + 0 имеем вид:

u'( H + 0) — ikHu (H + 0) = —2 ikHA (1.2.10)

Учитывая условия сопряжения (1.2.6), получим окончательное условие возбуждения при г = Н - 0:

ф(г = Н - 0> и'(Н - 0> - ¡кНфНи (Н + 0> = -2¡кНфНЛ (1.2.11)

Таким образом мы получим формулировку прямой задачи (1.2.4) -(1.2.11), которая позволяет определить и (г> при г е [0, Н ]. Следовательно, мы можем вычислить частотные характеристики поля смещения на земной поверхности:

и (г = 0> = и0 (с,р( г> ,ф( г>> (1.2.12)

и0 -нелинейный оператор зависящий от частоты с и действующий на р(г>,ф(г> .Обратная задача сейсмического зондирования заключается в

определении р(г>,ф(г> по измеренному на земной поверхности и(т^(с) в зависимости от частоты. Рассмотрим метод расчета и (г = 0 >. Заметим, что и (г >, в силу линейности задачи, с помощью метода переменных разделения можно представить в виде:

и (г> = и0 (с>у(г>, где ио (с) = и (г = 0) (1.2.13)

Функция V(г>, согласна (1.2.4) -(1.2.6), является решением следующей задачи Коши:

(ф( г > v'( г ))> +с2р( г > V ( г > = 0

V(г = 0> = 1; V (г = 0> = 0 (1.2.14)

(г> и ф(г>V '(г>-непрерывны

Полученная задача легко решается численно, а для кусочно-постоянных коэффициентов р(г>,ф(г> из системы (1.2.14) можно получить аналитическое

решение. Пусть коэффициенты уравнения заданы следующим образом:

р(г) = рп; ф(г)=фп при ге[гп-1>гп],пе[1,Щ,г0 = 0

/ ч / ч (1.2.15)

р(г = Н) = рн; ф(г = Н) =Фн при х>хм = Н

Тогда при г е [гп-1, гп ] решение задачи (1.2.14) имеем вид:

v(г) = V! СОБкп (г - гп-1) + кп (г - гп-1) (1.2Л6)

кп

гдеVn-! = V(гп-1), = гп_{ + 0), кп = . Продифференцировав (1.2.16)

и положив г = гп = гп-1 + кп, получим рекуррентную формулу перечета:

vn = V1cos knhn + knhn

К

V' ( ^n - 0) = - Vl sin knhn + V'-l cos knhn

Согласно условия (1.2.4) -(1.2.6) имеем вид:

ФгУ' ( ^n - 0 ) = Фп+/ ( ^n )

v ( Zn - 0 ) = ф v ( zn ) = ф v'n Kn ) фп Kn) фп n

(1.2.17)

(1.2.18)

(1.2.19)

Подставив (1.2.19) в (1.2.18), получим:

кф

v = - v ,

n n-1

sin ±

v' , cos knhn

n-1 n n

(1.2.20)

/ п п 1 п—1 Н Н V

ФН+1 ФН+1

Зная у0 = 1, у'0 = 0 , по формулам (1.2.18) и (1.2.20) находим последовательно у1, у{ до ,. Зная отношения = V(г = Н + 0) и = V '(г = Н + 0) ,

согласно (1.2.13), можно определить,

и (г = Н + 0) = и0 (ю)уу, и (г = Н + 0) = и0(т)^ (1.2.21)

Подставив полученные выражения в (1.2.20), найдем:

(1.2.22)

U0 (ю) = -

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

Прежде чем рассмотреть общий случай решение обратной задачи. Рассмотрим поведение и0 (т) для простейших моделей. Пусть данной слой имеет параметры фх, р1,^ , который находится на полупространстве с параметрами фН, рН , .Тогда согласно формулы (1.2.18) и (1.2.20) имеем:

v1 = cos k1h1; v| = - sin k1h1

фн

Подставив эти выражения в (1.2.22), найдем:

2ikHфHA

Uо М = -

или

Uо М =

ikнфн cos k1h1 ± k1ф1 sin k1h1

2 A(cos k1h1 ± i61H sin k1h1) cos2 k1h1 ± 61H sin2 k1h1

(1.2.23)

(1.2.24)

(1.2.25)

где вш = 1ж

\РнФи

Из уравнения (1.2.25) следует:

, ч ImU0(о) Л 7 , Q(o) = тт°) = вш tan kh ReU0 (о)

(1.2.26)

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

вш Р hi

ф <чФ 1 (1.2.27)

рНфН \ф

Рассмотрим два слоя с параметрами ф1, р1, к1 и ф2, р2, , на однородном

полупространстве имеет параметры фН,рН. Тогда, зная v!,V' согласно (1.2.23) вычислим по формулам (1.2.18) и (1.2.20):

где в12 =

i

Р1Ф1

Р2Ф2

v2 = cos k1h1 cos k2h2 - 6l2 sin k1h1 sin k2h2 . Продифференцировав выражение, получим:

(1.2.28)

v

к2ф2

Ф

(sin k2h2 cos k1h1 + 6l2 cos k2h2 sin k1h1)

н

Подставив (1.2.28) и (1.2.29) в (1.2.22), найдем:

U (о)= 2A(-12 + Í@2H^ 12 )

U 0 V»)- -2 + a2 2

12 2н 12

(1.2.29)

(1.2.30)

где

-12 = cos k1h1 cosk2h2 - 012 sin k1h1 sin k2h2 gX2 = cosk1h1 sin k2h2 + 012 sin k1h1 cosk2h2

a

Р2Ф2

2 H

Рнф

H

легко видно, что

Q (о) =

Im U0 (°)_a2н^ 12

ReU0 (о) -12 Тогда для n-ого слоя, Q(о) можно приставить в виде:

(1.2.31)

QM = = (1.2.32)

где

с = £ sin kh + 0 с cos kh

' mn = pm n n mn^pm n n

£ = £ cos kh -0 с sin kh

mn pm n n mn pm n n

m = n -1,p = n - 2,n = 1..N -1 Это означает, что зная Q(ю), можно определить следующие величины: knhn, k

0 H ,0 , — или следующие относительные величины:

К

hn Рп Фп

К' Рт

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

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

проверим осуществимость подхода, хотя наш метод установлен на теоретическом выводе.

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

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

В данной главе подробно рассмотрен метод обобщенного коэффициента отражения-передачи (МОК-О/П)) (The generalized reflectiontransmission coefficient method) решения прямой задачи для поверхностных бегущих волн в слоистой среде. предложен новый метод тензор сейсмического импеданса, который позволяет быстро и устойчиво получить характеристики поверхностных бегущих волн.

2.1 Метод расчета характеристики бегущих поверхностных волн в полупространстве

Рассматривается изотропная упругая среда, в которой характеристики Ламе Л, ju и плотность р являются произвольными функциями координаты z. Ось z направлена вниз, полупространство z > 0 , является упругой средой, в которой характеристики Л, ju, р положим постоянными. Тогда уравнение упругой волны для изотропной среды имеем вид:

где и - смещение.

Если бегущая поверхностная волна бежит по оси ОХ как в'7Х то имеем:

(Л + 2j)grad divU - jurot rotU + со2pU = 0

(2.1.1)

U = u (z )eirx

(2.1.2)

где и (z) их (г) ,0,иг (г))-вектор зависимый от глубины г, у распространения бегущей поверхностной волны. Тогда имеем:

постоянная

divU = U + dUz = (iyux (z) + u'z (z)) eiYX dx dz

grad div U = (ly( lyux (z) + u'z (z)); 0; lyu'x (z) + u'' (z)) e'

dU.

elYx

(2.1.3)

(2.1.4)

л

ix ly iz

(2.1.5)

iyx

rot rotU = (~u"x + iyu'z; - u"y + y2uy; iy(u'x - z>uz))< Подставим (2.1.3) - (2.1.6) в уравнении (2.11) получим:

(Л + 2ц»( iYux (z) + uZ (z))- м(-и'Х + iYuZ) + aPux = 0

-ц(К -Y2uy) + ®2Puy = 0

(Л + 2м)(iYuX(z) + <(z))- ц(ir(uX - iYuz)) + aPuz = 0

(2.1.6)

(2.1.7)

Волна распространяет по оси ОХ , поэтому компонента ОУ может отсутствовать. Тогда система уравнений (2.1.7) может быть переписана для двух компонент в виде:

иХ + аХ + =0

где

uz + auz + b2u' = 0

z z z 2 x

2 pa)2 - (Л + 2ц) Y , (Л + ц)iy 2 pa2 - НУ2 и (Л + ц)iy

(2.1.8)

ax =

bi =

> az =

b2 =

ц ц Л + 2 ц Л + 2 ц

Приz ^ да, ux, uz ^ 0, надо, чтобы Imax ф 0 и Imaz ф 0, т.е.:

(2.1.9)

ax = 1

. 1(Л + 2ц)у -paa_ a = . \цу--ра

ц

Л + 2ц

(2.1.10)

Из (2.1.10) легко видно отношение между характеристикой бегущей волны и волновым числом поперечной и продольной волны:

2 2 р ap 2 „2 „2

Y > a — > -—J— ^ у > Ks > Kp

(2.1.11)

Л Л + 2ц

где К, Кр - волновые числа поперечных и продольных волн.

Пусть V = (у1 = их, у2 = и2, у3 = их , у4 = иг, тогда систему уравнений (2.1.8)

можно представить в виде:

dz dz

2

= V,

dz dz

= - аЛ - Ъ^4

'аг~^2 - b2V3

(2.1.12)

Перепишем в векторную форму:

су

где А =

0 0 1 0

0 0 0 1

- ах2 0 0 - Ъ]

0 - ах2 -Ъ2 0

с?

л

= АV, V = (VI,

V , V2, Vз, V4

(2.1.13)

. Пусть V = Увп, подставив (2.1.13) получим:

(А -пЕ)г = 0

Если полученная система (2.1.14)существует решение, то имеем:

А - пЕ) = ёй

г -п 0 1 0 Л 0 -п 0 1 -а2х 0 -п

0 - а 2

п

0

п4 + (а2х + а"; - Ъ^ )п2 + а2ха] = 0 Корни данного уравнения получаются:

п = ±

ЪЪ2 - (аХ + аг2) а_2 + а\ - ЪЪ )2 - 4«

2 2 :2а

х z

2

(2.1.14)

(2.1.15)

(2.1.16)

(2.1.17)

Теперь нам надо определить полученные корни. Тогда подставив (2.1.9) в (2.1.17), получим:

Ъ1Ъ2 - (а] + а^= 2 Ъ1Ъ2 - (а2 + а:2 ) + ^( а2х + а] - ЪХЪ2 )2 - 4а2ха] =2

г2 Ср

V

У

м

/

2

с р (Л + 2м)

(2.1.18)

(2.1.19)

Согласно условию (2.1.11) легко видно, что (2.1.18) и (2.1.19) являются положительными числами при постоянных Л, м, Р> У,® . Тогда корни п имеют действительные значения:

2 С Р

\У---, п2

2 С Р

У

м

м

(2.1.20)

пз = АУ

С р

п = - У

С р

(Л + 2м)' V (Л + 2м) Следовательно, общее решение уравнение (2.1.8) имеет вид:

их = Схеп + С2е~п + Съеп + САе~п

= -пз

(2.1.19)

и г = Охеп + П2еп + П3епз + ПАе~пз

Для полупространства бегущая волна не имеет отраженную волну, т.е.:

их = Схе~т + С2е

(2.1.20)

и = Це~п + Де"п

z 1 2

Подставив (2.1.20) в систему уравнений (2.1.8), параметр Д, П2 можно представить в виде:

-2 , „2

П = ах +п С 1 _ и 1

2 , „2

(2.1.21)

П2=ах+п С Ь2п4

Тогда в решении (2.1.21) имеется только 2 неизвестных параметра Сг.и С2. Для полупространства имеем граничное условие: на свободной поверхности нормальной и касательной компоненты тензора напряжения

Т = {<Гу } = ЛЗуШуи + /

ди ди. —L + —-

д. д/

=0, т.е.:

а

Х2

I = 0

л

ди ди '

_Х +_I

д1 дх

V у

л( />и1 (1) + и'х (1))е

I = 0

¡ух

= 0

I = 0

а

22

I = 0

Л

(ди ди Л

х + I

дх дI

ди + 2/— дI

(2.1.22)

I = 0

((Л + 2/) иI (I) + Л/уих (I) е/^х

= 0

I = 0

Подставив (2.1.20) - (2.1.21) в граничном условии (2.1.22) получим:

Л ■ 2 Л /

V 2 2 Л ( 2 2 Л

п + ь4 С1 п + п4 + ь4

V

ЬП

У

2 , „2 Л

ь3 + 3 ¿1

С^ +

ь2П4 2 2 Л

Ь3 + 02+32

Ь2 У

С2 еп1 = 0

у

С2 еп = 0

(2.1.23)

и Л/у

где Ь3 =1,4 , Ь4 = Л + 2/

Зная существование коэффициенты С1, С2, детерминант системы уравнения (2.1.23) равен нулю:

А) =

,2 Л

ЬП

av + п , п г от + п ,

П2 + ———Ь4 еп п4 + ——— Ь

2

ь2П4

ь3 + 0.1П2 3 ¿1

2

ь.+

2

72 У

= 0

(2.1.24)

Тогда мы можем ввести функцию:

¡(Л,/,р,у,а) = ёе1( А) = 0 Упростив детерминант (2.1.24) получим:

/(Л р, ^ со) =(г2 + П22 )2 - 4ППГ2 = 0 49

(2.1.25)

(2.1.26)

Таким образом, мы получили дисперсионное уравнение для расчета у-постоянная распространения бегущих поверхностных волн. В полупространстве коэффициенты Ламе Л, и и плотность р являются постоянными, у всегда зависит от частоты о . Для вычисления характеристики у введем примерные коэффициенты Л, и, р, о и с помощью ЭВМ получим следующую таблицу (Таб.2.1):

= 10584 VS =л\м!р= 2.1 км / сек, VP = у/л + 2^¡р = 3.5 км / сек

л = 8232 KS = аЦ^р = 0.043, KP = аЦл + 2^р= 0.026

р = 2400 кг / м3

Y = 0.05796

а = 0.09 рад / сек

Таблица. 2.1

Из таблицы. 1 легко видно, что полученный результат у соответствует условию (2.1.11). Для проверки устойчивости задачи (2.1.25) задаем ряд параметров о и вычисляем у(Рис.2.\). В рис.2.1 показано сравнение между полученными характеристиками бегущих волн (у) и волновыми числами поперечных и продольных волн (К8 и Кр) по каждому узлу и доказано что решение задачи (2.1.25) устойчиво.

з 0 25

5

„ 0.2

к

^ 0.1В

"3

ч \\ 1 1 1 1 - * - характеристики бегущих волн (у) —волновые чисел поперечных волн (Ks) --■о-волновые чисел продольных волн (Кр)

м м

u о W - \ 1 S

\ N \ N 4 Л*

Ъ,

...........А ..........<■ 1 >.......

...........<9..........< >

Период(Т/сек.)

Рисунок 2.1. Сравнение между полученными характеристиками бегущих волн (у) и волновыми числами поперечных и продольных волн (К8 и Кр ) для полупространства.

Предположив у = — , где с - фазовая скорость бегущей волны, уравнение

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