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

004614925

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

А

РЕШЕТОВА ГАЛИНА ВИТАЛЬЕВНА

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

05.13.18 - математическое моделирование, численные методы и комплексы программ

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

НОВОСИБИРСК-2010

- 2 ЛЕК 2010

004614925

Работа выполнена в Учреждении Российской академии наук Институте вычислительной математики и математической геофизики Сибирского отделения РАН (ИВМиМГ СО РАН)

Научный консультант: доктор физико-математических наук, академик РАН Михайленко Борис Григорьевич

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

доктор физико-математических наук, профессор Лаевский Юрий Миронович

доктор физико-математических наук, профессор Роменский Евгений Игоревич

доктор физико-математических наук, профессор Крукиер Лев Абрамович

Ведущая организация: Учреждение Российской академии наук Институт вычислительного моделирования Сибирского отделения РАН (ИВМ СО РАН, г. Красноярск)

Защита состоится 15 февраля 2011 г. в 15 час. на заседании диссертационного совета Д 003.061.02 при Учреждении Российской академии наук Институте вычислительной математики и математической геофизики СО РАН.

Адрес: пр-т Ак. Лаврентьева, 6, Новосибирск, 630090. Факс: (383) 3308783, e-mail: kosova@rav.sscc.ru

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

Автореферат разослан _ ноября 2010 г.

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

С.Б. Сорокин

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

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

Актуальность

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

К ним относится, например, изучение сейсмоакустических волновых полей при акустическом каротаже скважин. Действительно, даже в простейшем случае открытой (необсаженной) скважины приходится иметь дело с двумя характерными размерами - диаметром скважины 0.1-0.2 м и максимальным расстоянием источник - приемник 10-12 м. При наличии обсадки приходится также учитывать влияние стальной трубы, отделяющей заполненную жидкостью скважину от вмещающей среды. Толщина стенок трубы составляет первые сантиметры и является еще одним масштабом задачи.

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

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

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

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

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

Научные задачи

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

2. Разработать численный метод расчёта сейсмических волновых полей в трёхмерных неоднородных разномасштабных упругих средах (с учетом системы трещин, каверн) с использованием конечно-разностных схем с локальным измельчением сеток по пространству и времени. Создать на этой основе научно-исследовательский вариант параллельного программного обеспечения.

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

Теория и методы исследования

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

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

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

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

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

разностной сетки с сохранением порядка аппроксимации конечно-разностной схемы;

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

- теория параллельных вычислений, в том числе методы пространственной декомпозиции областей для разработки параллельных версий алгоритмов, а также использование последних эффективных разработок программного интерфейса обмена сообщениями в стандарте MPI-2, в частности, новые возможности параллельных коллективных обменов и ввода/вывода в среде MPI-2 I/O.

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

Окончательная верификация программного обеспечения проводилась в нефтяных и сервисных компаниях (Schlumberger, Total, ООО «РН-КрасНИПИНефть») и Югорском НИИ Информационных технологий путём сопоставления результатов численного моделирования и лабораторных и натурных наблюдений.

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

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

2. Теоретически и экспериментально обоснованный численный метод расчёта сейсмических волновых полей в трёхмерных неоднородных разномасштабных средах (с учетом системы трещин, каверн), разработанный на основе конечно-разностных схем с локальным пространственно-временным измельчением сеток. Научно-исследовательский вариант параллельного программного обеспечения.

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

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

Научная новизна и личный вклад

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

1.1. Интегро-интерполяционным методом (методом баланса) построена консервативная конечно-разностная схема для численного моделирования обобщенной модели стандартного линейного тела GSLS для нескольких релаксационных механизмов, проведен анализ устойчивости схемы и исследованы ее дисперсионные свойства для минимизации численной дисперсии скорости.

1.2. В ходе реализации метода разработаны оригинальные подходы для эффективного решения ряда ключевых задач:

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

- ограничение расчетной области слабо отражающими слоями в цилиндрической системе координат с помощью модифицированного метода PML (аббревиатура английского Perfectly Matched Layer) без расщепления переменных по азимуту;

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

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

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

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

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

Создан научно-исследовательский вариант параллельного программного обеспечения для численного расчёта сейсмических волновых полей в трёхмерных неоднородных разномасштабных упругих средах, основывающийся на новых возможностях параллельных коллективных обменов и параллельного ввода/вывода больших массивов данных с использованием библиотек MPI-2 I/O (Input/Output):

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

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

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

3. Разработан, теоретически и экспериментально обоснован и протестирован новый спектрально-разностный метод численного моделирования сейсмических волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела, основывающийся на отделении времени с использованием интегрального преобразования Лагерра (совместно с Б.Г.Михайленко).

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

но реализован новый оригинальный способ конструирования PML без расщепления, основанный на преобразовании Лагерра по времени.

Создан научно-исследовательский вариант параллельного программного обеспечения для численного моделирования волновых полей в трёхмерных неоднородных средах с резкоконтрастными границами раздела. Численными экспериментами установлено существование предсказанной ранее теоретически поверхностной волны Стоун-ли-Шолтэ и модифицированной волны Лэмба.

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

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

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

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

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

- для оценки качества выполненного гидроразрыва пласта;

- для исследования скользящих волн (creeping waves) на стенке скважины для определения анизотропии прискважинной зоны.

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

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

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

- при моделировании волновых полей в средах с детальным описанием неоднородных включений (конструкций).

3. Разработанный новый спектрально-разностный метод численного моделирования волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела на основе интегрального преобразования JIareppa расширяет рамки использования численных методов при решении ряда важнейших задач прикладного характера, таких как:

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

- раннее обнаружение цунамигенного землетрясения;

- контроль подземных ядерных взрывов.

Реализация результатов

Результаты диссертации составили основу проекта «Численное моделирование взаимодействия сейсмических волн с кавернозно-трещиноватыми коллекторами в трёхмерно-неоднородных средах реалистичного строения», занявшего первое место в Конкурсе проектов в сфере высокопроизводительных вычислений под девизом «Невозможное стало возможным», организованного корпорацией Intel совместно с Российской корпорацией нанотехнологий («РОСНАНО») в 2010 году.

Представленный в диссертации спектрально-разностный метод на основе интегрального преобразования JIareppa получил развитие в серии работ зарубежных и отечественных исследователей (Läcik, Chung, Sarkar, Jung, A.A. Михайлов, А.Ф. Мастрюков, В.Н.Мартынов). Эффективность этого метода для численного моделирования электромагнитных волновых полей продемонстрирована в работе Chung et al. (2003), в которой убедительно показано, что применение преобразования JIareppa дает ускорение времени вычислений в 100 раз по сравнению с применением обычных конечно-разностных схем при численном моделировании двумерной системы уравнений Максвелла.

Результаты, представленные в диссертации, были включены в Основные результаты исследований Сибирского отделения РАН по приоритетным направлениям развития науки и техники и в Список достижений СО РАН в 2003, 2005, 2006 годах.

Основным направлением научно-исследовательских работ, выполняемых в ИВМиМГ СО РАН с участием диссертанта, являлось развитие данного подхода к численному моделированию подземных и наземных ядерных взрывов в рамках проекта по программе СО РАН 1.4.1.1. «Математическое моделирование природных и техногенных геофизических полей в средах сложной геометрии и реологии».

Все исследования, проводимые по теме диссертации, являются составной частью планов НИР Института, а их выполнение постоянно поддерживалось Российским фондом фундаментальных исследований в рамках проектов № 00-05-65323, № 04-05-64177, № 06-0564149, № 07-05-00538 № 10-05-00233.

Внедрение научных результатов

Созданный диссертантом программный продукт по моделированию сейсмических волн в двумерно-неоднородных сложнопостро-енных упругих средах был внедрен в мобильный программно-аппаратный комплекс, созданный в рамках комплексного проекта «Разработка комплексной технологии поиска и разведки углеводородов в сложнопостроенных, глубокозалегающих месторождениях» по Госзаказу 2005-РИ-00.0/009/202, шифр - RU.IBMMG.00103-01 34 01 (Акт о внедрении от 25.11.2006).

Разработанное автором программное обеспечение для моделирования трёхмерных сейсмоакустических волновых полей в средах с затуханием (акустический каротаж) в настоящее время используется инженерными технологическими центрами компании БсЫитЬе^ег (в рамках контракта с ИНГГ СО РАН, Акт о внедрении от 01.08.2009).

Разработанное диссертантом программное обеспечение для моделирования сейсмических волновых полей в трехмерных неоднородных средах с субсейсмическими неоднородностями (кавернозно-трещиноватые коллекторы) используется ООО «Геола» с целью изучения проявлений ориентации коридоров трещиноватости по полю рассеянных волн (Акт о внедрении от 14.09.2010).

Апробация, публикации, объем и структура диссертации

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

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

- Генеральные Ассамблеи Европейского геофизического общества (EGU) - 2000, 2002, 2004;

- Коференции Европейской ассоциации геофизиков (EAGE) - 2000, 2003,2007, 2008,2010;

- Конференции сообщества геофизиков-исследователей (SEG) - 2006, 2008, 2010;

- Международные конференции "Математические и численные аспекты теории распространения волн" - 2007, 2009;

- 7-я Международная европейская конференция по вычислительной математике, Австрия, Грац, 2007;

- 8-я Международная конференция "Математические и численные аспекты теории распространения волн", Великобритания, Рединг, 2007.

Диссертация состоит из введения, трех глав и заключения. Содержит 298 страниц, в том числе 98 рисунков. Библиография содержит 171 наименование.

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

Успешному проведению исследования способствовала поддержка академиков РАН А.Н. Коновалова и Б.Г. Михайленко, оказавших большое влияние на формирование научных взглядов соискателя.

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

Автор ценит всестороннюю поддержку, постоянное внимание к работе и благодарит всех сотрудников Лаборатории численного моделирования сейсмических полей Института вычислительной математики и математической геофизики Сибирского отделения РАН, а также сотрудников Лаборатории вычислительных методов геофизики Института нефтегазовой геологии и геофизики им. A.A. Трофимука Сибирского отделения РАН, особенно В.В.Лисицу и Е.В. Лыся.

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

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

Глава 1. Численное моделирование сейсмоакустических волновых полей при акустическом каротаже в трёхмерных неоднородных средах с поглощением

Первый раздел посвящен истории вопроса и описанию современных методов численного моделирования сейсмоакустических волновых полей при акустическом каротаже. Сформулированы особенности постановки, связанные с наличием в упругом пространстве скважины, заполненной жидкостью (Biot, 1952; Крауклис П.В. и Крауклис JI.A., 1976). Введено понятие поглощения сейсмической энергии и прокомментированы математические модели, позволяющие описывать волновые процессы с затуханием (Chistensen, 1982; Robertsson et al., 1994).

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

Во втором разделе формулируется математическая постановка задачи и вводятся основные модели среды с поглощением: модель Максвелла, модель Фойгта, стандартная линейная модель твёрдого тела (Curtin and Sternberg, 1962; Christensen, 1982). Наиболее подробно исследуется последняя модель, для которой получены частотные зависимости фазовой скорости и добротности. Далее вводится модель обобщённого стандартного линейного твёрдого тела и доказывается, что для разумных значений добротности (Q> 10) для её описания достаточно двух релаксационных механизмов.

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

диг _ дагг + 1 даГ9 ^ даГ2 | сг„ _ 8t дг г дер & г

dt дг г dtp д: г 1 да™ За

СГ

dt дг г dtp д:

Закон Гука:

= [л/(1 + Ьт") - 2//(1 + Ьт" )]© + 24 + Ь* Ъ- + £ггг1 + Ф; от дг /=1

^ = [л/(1 + 1тр) - 2/л{\ + Ьт* )]© + 2 ц [1 + Ы* + X ^ + Ф;

ОГ СЕ (=1

ГА/(И-^гр>-2//(1 + /,гл')1©+2//[1 + ^г-л'| I +Ф;

д( 1 1 \г дф г ) % "

(2)

да,* д(

■ Ц

\ + " г

1 8иг 8нр ^ г д<р &

(тЬ

да

а/ 1 \ 5г г дер г

где введены обозначения 0 = ¿/п>я,Ф =

<5(Г-Гй,<р~(р0,2-2й)

2 ЛГ

Уравнения на переменные памяти:

дг, д1

п' - - \(Мтр -2цт*)© + 2//гл ^ + г„,

Зг

^ = _2Ит«)0 + ^ + ,,, у.

дг„

—\(Мтр -2цт*)<д + 2цт* & ты

1 д119 , «,

г д<р г

+ г„

«у*

дгн 1 л-( диг ди

(3)

дгк I 1 I .V

(\_ди^+ди1л г д(р дг

+ К

дг,

д1

1 \цтБ

ди? | 1 диг _ % " дг г д<р г

+ г.

Ц1

Система уравнений (1)-(3) решается с нулевыми начальными данными

иг\/=0~ ир|г=0 = |/=0 = 11=0 = °гг|<=0 = = = °«г|г=0 = = 0

и условиями контакта на границе скважины г = /?0, заключающимися

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

["г]'г=я„ = 0;[г^]|г=до = [г„]|г=яо ;[тг9]\г=к0 = 0. (5)

Для ограничения расчётной области автором реализована оригинальная версия PML, разработанная при её непосредственном участии (Костин, Решетова, Чеверда, 2006).

Конечно-разностная аппроксимация задачи (1)-(5) строится на сдвинутых сетках (Virieux, 1986) методом баланса (Самарский, 1977). Полученная при этом явная конечно-разностная схема является условно устойчивой (доказан спектральный признак устойчивости) и имеет второй порядок аппроксимации всюду вне границ разрыва упругих параметров среды. На границах разрыва она аппроксимирует условия согласования типа (5) с первым порядком аппроксимации, однако её решение сходится к решению исходной задачи со вторым порядком (Самарский, 1977; Lisitsa, Podgornova and Tcheverda, 2010).

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

Разбухание ячейки по азимуту. Для обеспечения примерно одинаковых размеров ячеек сетки по радиусу и по азимуту предлагается уменьшать в два раза дискретизацию по азимуту каждый раз, как радиус увеличивается в два раза, но при этом необходимо согласовывать сетки, заданные при разном азимутальном шаге. На рис.1 приведено взаимное расположение узлов при переходе с одной сетки на другую. Принимая во внимание 2п -периодичность всех участвующих в рассмотрении переменных, восполнение неизвестных значений (в красных узлах) выполняется с помощью интерполяционной процедуры, основанной на быстром преобразовании Фурье (БПФ), что обеспечивает экспоненциальную точность и тем самым позволяет избежать сколько-нибудь существенных артефактов при переходе с сетки на сетку. Кроме того, использование для этих целей процедуры БПФ обеспечивает максимальное быстродействие при выполнении такой интерполяции.

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

гш измельчения азимутального шага

лк#нмя измельчения азимутального шага

кону: /г, = (рис.2). Конечно, при этом аппроксимация становится уже первого порядка. Однако за счёт выбора оптимального для заданного частотного диапазона коэффициента изменения шагов удаётся добиться приемлемого уровня погрешности._

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

002 0 04 0 06 0 08 0.1 012 0.14 0.16

Рис.2. Сетка с периодическим измельчением по азимуту и квазиравномерным изменение шагов по радиусу в окрестности границы скважины.

Устранение математической особенности при г = О. Используемая сетка начинается с окружности г = 0.5ИГ и, следовательно, не содержит ось г = 0. Однако для вычисления некоторых компонент вектора скорости смещений и тензора напряжений необходимо знание

14

компонент иг,сгг9,<гп при г = 0. Так как скважина заполнена жидкостью, внедиагональные компоненты тензора напряжений внутри неё равны нулю, поэтому в нуле надо определить только радиальную компоненту вектора скорости смещений. Для этого используется интерполяция полиномами Лагранжа по шести соседним узлам.

Организация параллельных вычислений производится декомпозицией расчётной области на непересекающиеся дискообразные элементы, каждый из которых приписывается своему процессорному элементу. На каждом временном шаге соседние процессорные элементы обмениваются текущими значениями векторов скорости смещений и тензора напряжений, рассчитанными на соприкасающихся гранях элементарных дисков (рис.3). Высокая эффективность и масштабируемость параллельного программного обеспечения была достигнута благодаря организации асинхронных вычислений. Для этого использовалось сочетание неблокирующих МР1-процедур межпроцессорного обмена ЬепсИгесу и специальной организации вычислений внутри элементарной области - из её самой внутренней точки по направлению к соседним процессорным элементам. Оценка точности вычислений как для идеально-упругой среды, так и среды с поглощением проводилась путём сравнения численного решения с точным, а также выделением таких характеристик волнового поля, как времена вступления отдельных волн (продольных и поперечных) и уменьшение их амплитуды.

В заключение представлены результаты расчётов для нескольких реалистичных сейсмоакустических моделей, одна из которых показана на рис. 4. Расчётная область - цилиндр радиусом 1 м и высотой 6 м, радиус скважины (1) - 0.1 м, толщина охватывающей её стальной обсадной колонны (2) - 0.01 м, цементного кольца (3) - 0.04 м, толщина вертикальной трещины (7) - 0.02 м. Параметры модели:

1. буровой раствор в скважине: р = \Шкг/мъ,Ур = 1500л</с,£> = 65;

2. обсадная колонна: р = 7830 кг/м\Ур = 5600 м/с, V, = 327О.«/с,0 = 100;

3. цементное кольцо: р = 2400 кг/м\\'р = 4200м/с, V„ = 2425м/с,<2 = 80;

4. упругая среда №4: р = 2400 кг/м3,¥р = 4989м/с, V= 2605м/с,д = 100;

5. упругая среда №5: р = 2400 кг/м\Ур = 3208м/с, У,. = 1604и</с,0 = 60;

6. упругая среда №6: р = 2400 кг/м\Ур = 2650м/с, У, = 1219л</с,б = 15;

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

к процессорному элементу к-1

к процессорному элементу к+1

Рис.3. Декомпозиция расчётной области.

Рис. 4. Модель строения среды и конструкции скважины для моделирования проявления вертикальной трещины в сейсмоакустических полях: а) вертикальное сечение; б) горизонтальное сечение; в) увеличение в окрестности скважины.

1

£ ш

1 ¡7

■ -Т

0 1

ПИ

: !

Рис. 5. Серия моментальных снимков разности компонент а1Г тензора напряжения, рассчитанных для модели с вертикальной трещиной.

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

Глава 2. Конечно-разностные методы численного моделирования распространения сейсмических волн в трёхмерных неоднородных разномасштабных средах.

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

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

ванный на локальном пространственно-временном измельчении шагов сетки.

Конечно-разностное моделирование с использованием сеток с локальным измельчением по пространству и по времени. Типичная задача трёхмерного численного моделирования сейсмических волновых полей требует десятки и сотни гигабайт оперативной памяти, и следовательно, её решение невозможно без использования вычислительных систем с параллельной архитектурой. На сегодня наиболее распространённой и, пожалуй, самой эффективной является организация таких вычислений с использованием явных конечно-разностных схем, аппроксимирующих систему уравнений динамической теории упругости в скоростях-напряжениях. В работе используется стандартная схема на сдвинутых сетках (У теих, 1986), адаптированная для неоднородных сред с границами раздела с использованием метода баланса (Самарский, 1977). Переход с одного временного слоя на другой осуществляется в два шага: от скоростей для момента времени / к напряжениям в момент времени / + 0.5# и от напряжений в момент времени /+0.5Я к скоростям в момент времени 1+8.

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

Г7ГГГ

ЁН ггу:::: Р

— „-п....

1гф„ ТТГТ V

г 5 \

\

а)

юммиемм по простраспу

.Н—Д.—*

I

б)

хи

-Л-Х-*

> ►

4

В)

Рис. 6. а) Взаимное расположение поверхности с измельчением шага по времени (красный прямоугольник) и шага по пространству, б) Согласование ша-_гов по времени, в) Согласование шагов по пространству._

Измельчение шага по времени. Для одномерного случая измельчение шага по времени при фиксированном шаге по пространству представлено на рис. 66 и достаточно очевидно. Его модификации в двумерной и трёхмерной постановках отличаются только использованием соответствующего пространственного шаблона, оставляя неизменной структуру вычислений (Ие8Ье1:оуа е1 а1., 2009).

Измельчение шага по пространству. Для согласования сеток с различными шагами по пространству используется интерполяция на основе БПФ. Чтобы упростить объяснение, рассмотрим двумерную реализацию этого подхода. На рис. 6в изображено взаимное расположение крупной и мелкой сеток. Треугольниками обозначены узлы, в которых задаются компоненты векторов скорости смещения (ориентация треугольника совпадает с направлением компоненты вектора), а прямоугольниками и косыми крестиками - узлы, в которых задаются соответственно диагональные (а^ и а..) и внедиагональные (о^ ) компоненты тензора напряжений. Рассмотрим шаг, на котором по заданным компонентам вектора скоростей смещения вычисляются компоненты тензора напряжений. Как видно, для того чтобы вычислить компоненту а^ на границах мелкой сетки, необходимо знать векторы скоростей смещения в узлах, отмеченных красным треугольником. Но эти узлы не принадлежат крупной сетке, поэтому искомые значения в них должны быть получены интерполированием, для чего используется преобразование Фурье. Определяющими факторами этого выбора явились экспоненциальная точность такой интерполяции и высокое быстродействие БПФ. Единственное отличие в трёхмерном случае -необходимость выполнения двумерной интерполяции.

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

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

¿Ш

А

/7Ш

ю

2D

200 SOO 400 500 800 700 ВОО

число процессорных элементов

б)

Рис. 7. а) Сверху вниз: одномерная, двумерная и трёхмерная декомпозиции области. б) Ускорение при использовании одномерной (чёрная), двумерной (синяя) _и трёхмерной (красная) декомпозиции расчётной области._

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

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

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

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

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

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

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

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

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

Рис. 8. а) Две группы процессоров, б) Взаимодействие между ними

группа 1

Процессоры мелкой сетке

Процессоры на крупной сетке

Резервуар на мелкой сетке

Вмещающая среда на крупной сетке

месторождений нефти и газа в Восточной Сибири, в частности в Юрубчено-Тохомской зоне. В основу построенной численной модели резервуара положена схематическая модель, представляющая погруженный в однородную среду однородный слой, содержащий в себе два слоя с коридорами трещиноватости. Параметры вмещающей среды выбраны соответствующими типичной модели карбонатных пород: Ур = 4500 л//с, Vs = 2500м/с, р = 2500кг/мъ.

Резервуар содержит область трещиноватости коридорного типа с изменчивостью концентрации трещин от 0 в цельных фациях до максимального значения 0.3 (рис.9). Для дискретизации резервуара использовался пространственный шаг 0.5 м. В итоге получается массив из 4000x4000x400 точек. Декомпозиция области была выполнена таким образом, чтобы на каждый из процессорных элементов приходился куб 200x200x200 точек (примерно 1Gb RAM), что потребовало 800 процессорных элементов. Вмещающая среда заполняет куб 2000x2000x2000 м и дискретизуется с шагом 2.5 м, что соответствует массиву 800x800x800 точек. Выбирая такую же загруженность процессорных элементов, получаем необходимость использования 64 процессорных элементов. Таким образом, для моделирования взаимодействия сейсмических волн с тонкой структурой данного резервуара необходимо 864 процессорных элемента.

Vp в резервуаре Vp в резервуаре

Ось Z ах О

и( 41 осьХ и

у штшяш 1 IX*

Ось Z - 13« иве

а) ось Y б) Ось Y

Рис.9. Распределение коридоров трещиноватости в резервуаре, а) вид сбоку, б) вид сверху.

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

Рис. 10. Трёхкомпонентные сейсмограммы, вычисленные на профилях, параллельных (а) и ортогональных (б) коридорам трещиноватости.

Глава 3. Спектрально-разностный метод численного моделирования сейсмических волновых полей в резкоконтрастных средах на основе преобразования Лагерра по времени

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

|верх

^резервуара

основание [резервуара

[обменные

есть доминирующая частота сигнала, а к - некоторый множитель, зависящий от выбранной схемы (так, для схемы на сдвинутых сетках второго порядка к=20). Из приведённых соотношений видно, что в средах с высокими контрастами скоростей применение явных конечно-разностных схем, особенно при необходимости проведения моделирования волновых полей на большие расстояния и, следовательно, для больших времён, связано со значительными затратами вычислительных ресурсов.

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

Указанного недостатка лишён подход с применением преобразования Лагерра по времени. В этом случае получается строго отрицательно определённый эллиптический оператор. Более того, этот оператор не зависит от параметра разделения, что существенно уменьшает затраты на вычисление и хранение соответствующей матрицы. Этот подход предложен и обоснован при непосредственном участии автора (Конюх (Решетова), Михайленко, Михайлов; 2001; М1к1ш1епко, М1к-ИаПоу, КевЬе^уа, 2003). Одно из приложений подхода, опирающегося на преобразование Лагерра, подробно рассматриваемое в диссертации, связано с моделированием взаимодействия сейсмических и акусто-гравитационных волновых полей для модели неоднородной среды литосфера - атмосфера.

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

тальные работы, в которых приводятся экспериментальные данные о существовании взаимозависимости волновых процессов в литосфере и атмосфере. В частности, A.C. Алексеевым с соавторами описан эффект акусто-сейсмической индукции (Alekseev et al., 1996), заключающийся в возбуждении акустической волной от мощного вибратора интенсивных сейсмических поверхностных волн, распространяющихся на расстояние в десятки километров. Теоретическим исследованиям волновых процессов на границе упругого однородного полупространства с изотермической моделью однородной атмосферы посвящены работы A.B. Разина (1993), JI.A. Гасиловой и Ю.В. Петухова (1999). В них с использованием метода контурного интегрирования и результатов исследования дисперсионного уравнения для собственных решений построенной модели установлено существование поверхностной волны Стоунли-Шолтэ и модифицированной волны Лэмба.

Постановка задачи. Система линеаризованных уравнений, описывающая распространение акусто-гравитационных волн для модели неоднородной атмосферы в цилиндрической системе координат (г, z), имеет вид:

8uz 8Р 8иг дР дР 2

dt dz dt 8r ■ dt

dt dz "z

8P0

--и.

ÖP n

dur 1 du. —- + -ur +—-dr r dz

Uz+F(r,z,t).

dz

(6)

Здесь иг, и- - компоненты скорости движения частиц воздуха в волне, р, р - соответственно, возмущения давления и плотности под действием источника массы /г(г,г,о = ^(г-го)ло, - скорость звука в атмосфере, константа % - ускорение силы тяжести. Нулевыми индексами обозначены величины, относящиеся к среднему невозмущенному состоянию атмосферы. Зависимость значений гидростатического давления /'о и плотности атмосферы р0 от высоты г определяется соотношением Р0, р0 ~ехр(-2/Я), где Н - высота изотермической однородной атмосферы. Для покоящейся атмосферы в однородном поле силы тяжести имеют место формула: Р>0=-р^ и ра(?) = р^ ехр(-г///), где р -

плотность при 2=0 (граница литосфера - атмосфера).

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

du. da,, der„ er Po = -.+ + /г (r>Zj/);

dt dr dz r

dur d(j„ derrT &rr-<jm „ ,

Po -1Г = -Г-+-Г-+—-e± + PaFr(r,z,t)-,

dt dr dz r

г ;: со

der.. _ \dur J du. ur Л

--^НЬт)

dt v г l dz dr

der„ Г Эиг

—— = d —-+—i

3/ Ч & й-

В этих уравнениях и,., и. - компоненты скорости смещения упругой среды, а::.аrr,аqq - компоненты тензора напряжений, pQ(z) - плотность упругой среды, л(:) и m{z) - коэффициенты Ламе. При этом Fz,Fr - составляющие силы F(;.r,t) = F-(:,r,t)e: + Fr(:,rj)er, описывающие распределение локализованного в пространстве источника типа «центра давления», «вертикальной силы» или «диполя без момента» (если источник находится в упругой среде). Например, в случае источника типа «центра давления» компоненты силы F(z,r,t) имеют вид:

Fz{z,rJ)=S}rl~S(z-z0)f(t)- Fr(z,r,t) = ±^ls{z-ZQ)f(t).

2 яг dz аг 2лг

Функция f(t) задаёт изменчивость источника по времени. Для формального описания источников используются выражения, содержащие обобщенные функции и производные от них. При численном решении для обобщенных функций вводятся соответствующие гладкие аналоги. Задачи (6) и (7) решаются при нулевых начальных данных в атмосфере и упругом полупространстве, а также и при условии согласования на границе литосфера - атмосфера: да„\ (дР . all

1 о ^

dt ! I dt ^ "z Ii '

i /а im v

-I— + I , cr„\ = 0. (8)

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

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

Спецификой поставленной задачи является необходимость моделирования совместных волновых процессов в атмосфере и упругом неоднородном полупространстве, скорости распространения волн в которых отличаются в 10-15 раз. Поэтому для изучения волновых явлений необходимо рассматривать волновые поля для больших промежутков времени, а следовательно, и для больших пространственных областей. Для ограничения расчетной области в атмосфере и упругом полупространстве вводятся специальным образом сконструированные достаточно узкие (меньше длины волны) поглощающие слои, учитывающие специфику задачи, связанную с использованием преобразования Лагерра (Решетова, Чеверда, 2006). С их помощью удаётся обеспечить устойчивость алгоритма вне зависимости от длительности расчетов по времени.

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

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

V, = 1760.«/с, плотность р0 =2.3г/о13. Скорость распространения акустических волн в атмосфере со = 340м/с, плотность ро = 0001225г/ш3

(при г = 0). Источник типа «центра давления» расположен в упругом полупространстве.

На рис. 11а приведен численный снимок в момент времени г = 20с для горизонтальной компоненты скорости смещения в упругом полупространстве и давления в атмосфере. Источник типа «центра давления» расположен в полупространстве на расстоянии 3000м («Л, где я - доминирующая длина волны Р) от границы г=0. В первых вступлениях в упругом полупространстве регистрируется продольная р - волна от источника и отраженная рр - волна от границы с атмосферой и далее обменная ръ- волна. В атмосфере регистрируется коническая волна.

Волновая картина усложняется, если источник типа «центра давления» находится вблизи свободной поверхности на расстоянии 750м (1/4Я ). На рис. 116 приведен численный снимок в тот лее момент для той же самой горизонтальной компоненты скорости смещения в упругом полупространстве и давления в атмосфере. В упругом полупространстве, кроме классических волн, возникает «нелучевая» 8*-волна. Как видно из рисунка, она возникает и в атмосфере. Кроме того, по границе литосфера-атмосфера распространяется поверхностная волна Стоунли-Шолтэ, а в атмосфере - серия конических волн._

а)

б)

Puc.ll. Численный снимок горизонтальной компоненты волнового поля в момент времени 20 с для двух положений источника: 3000 м (а) и 750 м (б).

ЗАКЛЮЧЕНИЕ

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

1. На основе конечно-разностного подхода разработан и теоретически и экспериментально обоснован новый метод, создан ориги-

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

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

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

- использование для согласования сеток пространственной интерполяции на основе преобразования Фурье, что обеспечило её экспоненциальную точность и низкий уровень возникающих при этом артефактов (менее 0.1% от амплитуды падающей волны);

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

3. Разработан, теоретически и экспериментально обоснован и протестирован новый устойчивый эффективный спектрально-разностный метод численного моделирования волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела на основе интегрального преобразования Лагерра. Применение интегрального преобразования Лагерра по времени позволило:

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

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

- обеспечить возможность ограничения расчетной области устойчивыми PML без расщепления.

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

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

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

др.;

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

Основные публикации по теме диссертации

1. Кошох (Решетова) Г.В., Михайленко Б.Г., Михайлов A.A. Численное моделирование сейсмических полей в вязкоупругих средах на основе спектрального метода Jlareppa // Математическое моделирование. -2001.-Т. 13 (2). -С.61-70.

2. Ельцов И.Н., Кашеваров A.A., Решетова Г.В., Чеверда В.А. Проявление неоднородностей зоны проникновения в геофизических полях вдоль ствола скважины // Геофизика. -2004. -Т.6. -С. 17-21.

3. Ельцов И.Н., Решетова Г.В., Чеверда В.А. Численное моделирование процессов распространения сейсмоакустических полей с учетом неоднородности зоны проникновения // Физическая мезомехани-ка. -2005.-Т.8. -С.99-105.

4. Михайленко Б.Г., Решетова Г.В. Численно-аналитический метод решения задачи о распространении сейсмических и акусто-гравитационных волн для неоднородной модели Земля-Атмосфера // Сибирский журнал вычислительной математики. -2006. -Т.9. -С.37-46.

5. Михайленко Б.Г., Рсшетова Г.В. Математическое моделирование распространения сейсмических и акусто-гравитационных волн для неоднородной модели Земля-Атмосфера // Геология и геофизика, 2006, т.47(5), с. 547-556.

6. Решетова Г.В., Чеверда В.А. Использование преобразования Jla-герра для построения идеально подходящих поглощающих слоев без расщепления // Математическое моделирование. -2006. -Т. 18 (10). -С.91-101.

7. Костин В.И., Решетова Г.В., Чеверда В.А. Численное моделирование трехмерного акустического каротажа с использованием многопроцессорных вычислительных систем // Математическое моделирование. -2008. -Т.20(9). -С.51-66.

8. Гилбо Ж., Ланда Е., Решетова Г.В., Хайдуков В.Г., Чеверда В.А. Численное моделирование сейсмических волновых полей в двумерно-неоднородных упругих разномасштабных средах (карстовые включения)//Технологии сейсморазведки. -2008. -Т.З. -С.19-28.

9. Konyukh (Решетова) G.V., Krivtsov Y.I., Mikhailenko B.G. Numerical-Analytical Modeling of Seismic Wave Propagation in Vertically Inhomogeneous Media//Applied Mathematics Letters. -1998. -V.ll(l). -P.99-104.

10. Konyukh (Решетова) G.V., Mikhailenko B.G., Mikhailov A.A. Application of the integral Laguerre transforms for forward seismic modeling//Journal of Computational Acoustics. -2001. -V.9(4). -P.1523-1541.

11. Mikhailenko B.G., Mikhailov A.A., Reshetova G.V.Numerical vis-coelastic modeling by the spectral Laguerre method // Geophysical Prospecting. -2003. -V.51. -P.37-48.

12. Mikhailenko B.G., Mikhailov A.A., Reshetova G.V.Numerical modeling of transient seismic fields in viscoelastic media based on the Laguerre spectral method // Pure and Applied Geophysics. -2003. -V.160. -P. 1207-1224.

13. Pissarenko D.V., Reshetova G.V., Tcheverda V.A. 3D finite-difference synthetic acoustic logging in cylindrical coordinates // Geophysical Prospecting. -2009. -V.57(3). -P. 367-377.

14. Mikhailenko B.G., Reshetova G.V. Mathematical modeling of seismic and acousto-gravitational waves in a Heterogeneous earth-atmosphere model // Journal of Computational and Applied Mathematics. -2010.-V.234.-P. 1678-1684.

15. Pissarenko D.V., Reshetova G.V., Tcheverda V.A. 3D finite-difference synthetic acoustic log in cylindrical coordinates // Journal of Computational and Applied Mathematics. -2010. -V. 234. -P. 1766-1772.

Подписано к печати 03.11.2010 Формат 60x84/16. Объем 2 печ.л. Тираж 100 экз. Зак. 172. ООО «Омега Принт», пр-т. Ак. Лаврентьева, 6, Новосибирск 630090

ВВЕДЕНИЕ.

Глава 1. ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ СЕЙСМОАКУСТИЧЕСКИХ ВОЛНОВЫХ ПОЛЕЙ ПРИ АКУСТИЧЕСКОМ КАРОТАЖЕ В

ТРЁХМЕРНО-НЕОДНОРОДНЫХ СРЕДАХ С ПОГЛОЩЕНИЕМ

1.1. Изученность решения задачи.

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

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

1.2.2. Математические модели для сред с затуханием.

1.2.3. Переменные памяти для моделей БЬБ и ОЗЬ8.

1.2.4. Дифференциальные уравнения в цилиндрической системе координат.

1.2.5. Построение конечно-разностной схемы интегро-интерпо-ляционным методом.

1.2.6. Конечно-разностные уравнения.

1.2.7. Ограничение расчетной области для уравнений вязко-упруго сти.

1.2.8. УстохЧчивость и дисперсионный анализ.

1.2.9. Радиальное и азимутальное измельчение сетки.

1.2.10. Согласование сеток.

1.2.11. Аппроксимация в нуле.

1.2.12. Реализация параллельных вычислений.

1.2.13. Оценка скорости сходимости конечно-разностного решения.

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

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

1.3.1. Решение задачи с учетом особенности строения при-скважинной зоны.

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

1.3.3. Численные эксперименты: описание и анализ.

Глава 2. КОНЕЧНО-РАЗНОСТНЫЕ МЕТОДЫ ЧИСЛЕННОГО

МОДЕЛИРОВАНИЯ РАСПРОСТРАНЕНИЯ СЕЙСМИЧЕСКИХ ВОЛН В ТРЕХМЕРНО-НЕОДНОРОДНЫХ РАЗНОМАСШТАБНЫХ СРЕДАХ.

2.1. Изученность решения задачи.

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

2.3. Известные подходы к разработке конечно-разностных методов с использованием сеток с локальным пространственно-временным измельчением.

2.4. Конечно-разностное моделирование с использованием сеток с локальным измельчением по пространству и времени.

2.4.1. Измельчение сетки по времени.

2.4.2. Теоретический коэффициент отражения при измельчении сетки по времени.

2.4.3. Измельчение сетки по пространству.

2.4.4. Теоретический коэффициент отражения при измельчении сетки по пространству.

2.4.5. Экспериментальная оценка коэффициента отражения для пространственного измельчения сеток.

2.4.6. Обоснование устойчивости метода.

2.5. Параллельное программное обеспечение.

2.6. Численное моделирование сейсмических волновых полей с учетом мелкомасштабных неоднородностей.

2.6.1. Современное численное моделирование полей рассеянных волн.

2.6.2. Описание численных экспериментов.

Глава 3. СПЕКТРАЛЬНО-РАЗНОСТНЫЙ МЕТОД ЧИСЛЕННОГО

МОДЕЛИРОВАНИЯ СЕЙСМИЧЕСКИХ ВОЛНОВЫХ ПОЛЕЙ В РЕЗКОКОНТРАСТНЫХ СРЕДАХ НА ОСНОВЕ

ПРЕОБРАЗОВАНИЯ ЛАГЕРРА ПО ВРЕМЕНИ.

3.1. Изученность решения задачи.

3.2. Теоретическое обоснование.

3.2.1. Обобщенные многочлены и функции Лагерра.

3.2.2. Интегральное преобразование Лагерра.

3.3. Преобразование Лагерра по времени для спектрально-разностного метода. Постановка задачи.

3.3.1. Одномерная система уравнений акустики.

3.3.2. Одномерное волновое уравнение.

3.3.3. Трехмерная система динамической теории упругости.

3.3.4. Сравнительный анализ интегральных преобразований

Фурье и Лагерра.

3.3.5. Оценка точности приближенного решения.

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

3.3.7. Оптимизация вычислений.

3.3.8. Программная реализация метода. Описание численных экспериментов.

3.4. Ограничение расчетной области при использовании преобразования Лагерра.

3.4.1. Способы ограничения расчетной области.

3.4.2. Реализация PML во временной области.

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

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

3.5.1. Изученность решения задачи.

3.5.2. Анализ известных методов и подходов.

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

3.5.4. Спектрально-разностный метод решения задачи.

3.5.5. Ограничение расчетной области.

3.5.6. Численные расчеты.

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

Актуальность

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

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

Наиболее интересным является изучение сейсмоакустических волновых полей при акустическом каротаже скважин. Учёт цементного кольца и металлической обсадной колонны, наличия глиняной корки и системы трещин гпдроразрыва, а также сложного геологического строения прискважинной зоны с необходимостью приводит к решению задачи в резкоконтрастной разномасштабной среде. Действительно, даже в простейшем случае открытой (необсаженной) скважины приходится иметь дело с двумя характерными размерами - диаметром скважины 0.1-0.2 м и максимальным расстоянием источник - приемник 10-12 м. При наличии обсадки приходится также учитывать влияние стальной трубы, отделяющей заполненную жидкостью скважину от вмещающей среды. Толщина стенок трубы составляет первые сантиметры и является еще одним масштабом задачи.

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

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

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

Научные задачи

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

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

3. На основе интегрального преобразования Лагерра разработать спектрально-разностный метод численного моделирования сейсмических волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела.

Теория и методы исследования

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

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

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

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

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

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

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

- теория параллельных вычислений, в том числе методы пространственной декомпозиции областей для разработки параллельных версий алгоритмов, а также использование последних эффективных разработок программного интерфейса обмена сообщениями в стандарте MPI-2, в частности, новые возможности параллельных коллективных обменов и ввода/вывода в среде MPI-2 I/O

Тестирование разработанных алгоритмов и программ проводилось в несколько этапов. Сначала оно выполнялось на простейших моделях, допускающих точные аналитические решения. Для более сложных моделей проверялась корректность описания ряда ключевых характеристик изучаемых физических процессов: времён первых вступлений, связь добротности с амплитудой волн в вязкоупругой среде, наличие в результатах численных расчетов теоретически предсказанных волн, таких как волна Стоунли-Шолтэ и модифицированной волны Лэмба при моделировании взаимодействия литосфера-атмосфера и семейства трубных волн в задаче моделирования акустического каротажа. По возможности проводилось сравнение с численными расчетами, полученными другими исследователями по программам, основанным на других алгоритмах (Вишневский, 2002; Lys, 2008). Окончательная верификация программного обеспечения проводилась в нефтяных и сервисных компаниях (Schlumberger, Total, ООО "РН-КрасНИПИНефть") и Югорском НИИ Информационных технологий путём сопоставления результатов моделирования и лабораторных и натурных наблюдений.

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

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

2. Численный метод расчёта сейсмических волновых полей в трёхмерно-неоднородных разномасштабных средах (с учетом системы трещин, каверн) разработанный на основе конечно-разностных схем с локальным пространственно-временным измельчением сеток. Научно-исследовательский вариант параллельного программного обеспечения для расчета и изучения волновых полей, рассеянных от скоплений микрон еоднородно-стей (трещин, каверн).

3. Спектрально-разностный метод численного моделирования сейсмических волновых полей для трёхмерных неоднородных сред с резкокон-трастными границами, разработанный на основе интегрального преобразования Лагерра по времени.

Научная новизна и личный вклад

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

1.1. По результатам изучения известных математических моделей была выбрана и обоснована модель стандартного линейного тела GSLS (аббревиатура английского Generalized Standard Linear Solid) для корректного описания сейсмоакустических волновых полей в средах с поглощением сейсмической энергии (вязкоупругие модели среды).

1.2. Интегро-интерполяционным методом (методом баланса) построена консервативная конечно-разностная схема для численного моделирования обобщенной модели стандартного линейного тела GSLS для нескольких релаксационных механизмов, проведен анализ устойчивости схемы и исследованы ее дисперсионные свойства для минимизации численной дисперсии скорости.

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

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

- ограничение расчетной области слабо отражающими слоями в цилиндрической системе координат с помощью модифицированного метода PML (аббревиатура английского Perfectly Matched Layer) без расщепления переменных по азимуту;

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

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

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

1.4. Диссертантом создан научно-исследовательский вариант параллельного программного обеспечения для численного моделирования сей-смоакустических волновых полей при акустическом каротаже для трёхмерно-неоднородных вязкоупругих сред, оснащённый дружественным графическим интерфейсом пользователя. Программный продукт был применен для изучения проявления неоднородностей призабойной зоны скважины в сейсмоакустических волновых полях по данным акустических методов ГИС. На основании выполненной серии численных экспериментов была сформулирована гипотеза о зоне проникновения как "динамическом волноводе".

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

Лично Решетовой Г.В. создан научно-исследовательский вариант параллельного программного обеспечения для численного расчёта сейсмических волновых полей в трёхмерно-неоднородных разномасштабных упругих средах, основывающийся на новых возможностях параллельных коллективных обменов и параллельного ввода/вывода больших массивов данных с использованием библиотек MPI-2 I/O (Input/Output):

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

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

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

3. Разработан, теоретически и экспериментально обоснован и протестирован новый спектрально-разностный метод численного моделирования сейсмических волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела, основывающийся на отделении времени с использованием интегрального преобразования Лагерра (совместно с Б.Г.Михайленко). Метод опробован и показал свою эффективность для:

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

Лично Г.В.Решетовой создан научно-исследовательский вариант параллельного программного обеспечения для численного моделирования волновых полей в трёхмерных неоднородных средах с резкоконтрастными границами раздела. С использованием программного обеспечения проведены численные расчеты сейсмических и акусто-гравитационных волновых полей для неоднородной резкоконтрастной среды литосфера-атмосфера. Численно установлено существование теоретически предсказанной поверхностной волны Стоунли-Шолтэ и модифицированной волны Лэмба.

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

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

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

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

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

- для оценки качества выполненного гидроразрыва пласта;

- для исследования скользящих волн (creeping waves) на стенке скважины и определения анизотропии околоскважинного пространства.

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

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

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

- при моделировании волновых полей в средах с детальным описанием неоднородных включений (конструкций).

3. Разработанный новый спектрально-разностный метод численного моделирования волновых полей для трёхмерных неоднородных сред с резкоконтрастными границами раздела на основе интегрального преобразования Лагерра расширяет рамки использования численных методов при решении ряда важнейших задач прикладного характера, таких как:

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

- раннее обнаружение цунамигенного землетрясения с удаленным от берегов очагом;

- контроль подземных ядерных взрывов.

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

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

- поверхностных волновых полей в зонах малых скоростей, в рыхлых приповерхностных слоях и зонах выветривания;

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

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

Реализация результатов

Представленный в диссертации спектрально-разностный метод на основе интегрального преобразования Лагерра получил развитие в серии работ зарубежных и отечественных исследователей (Läcik, Chung, Sarkar, Jung, A.A. Михайлов, А.Ф. Мастрюков, В.В.Мартынов). Эффективность этого метода для численного моделирования электромагнитных волновых полей продемонстрирована в работе Chung et al. (2003), в которой убедительно показано, что применение преобразования Лагерра дает ускорение времени вычислений в сто раз по сравнению с применением обычных конечно-разностных схем при численном моделировании двумерной системы уравнений Максвела.

Эффективность применения преобразования Лагерра для численного моделирования сейсмических волновых полей в неоднородных вяз-коупругих средах с произвольными функциями последействия в интегральных соотношениях Больцмана обоснована в кандидатской диссертации Михайлова A.A. "Численное моделирование нестационарных сейсмических полей в неоднородных упругих и вязкоупругих средах", защищенной в Институте вычислительной математики и математической геофизики СО РАН в 2006 году с участием автора в качестве научного консультанта, соавтора публикаций и научно-исследовательских работ по этой теме.

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

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

- "Численно-аналитический метод моделирования распространения и взаимодействия сейсмических и акусто-гравитационных волн для неоднородной модели Земля-Атмосфера", 2005г.;

- "Математическое моделирование вибросейсмического мониторинга сложнопостроенных сред", 2006 г.

Основным направлением научно-исследовательских работ, выполняемых в ИВМиМГ СО РАН с участием диссертанта, являлось развитие данного подхода к численному моделированию подземных и наземных ядерных взрывов в рамках проекта по программе СО РАН 1.4.1.1. "Математическое моделирование природных и техногенных геофизических полей в средах сложной геометрии и реологии".

Все исследования, проводимые по теме диссертации, являются составной частью планов НИР Института.

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

- 2000-2002, проект РФФИ № 00-05-65323 "Теоретические и экспериментальные исследования возможности фокусирования сейсмической энергии распределенных поверхностных виброисточников во внутренних областях неоднородной среды" (руководитель Г.В. Решетова);

- 2004-2006, проект РФФИ № 04-05-64177 "Теоретические и экспериментальные исследования возможности обнаружения и локализации областей изменчивости геомеханических параметров по вызванным ими вариациям стационарных вибросейсмических полей" (руководитель Г.В. Решетова);

- 2006-2008, проект РФФИ № 06-05-64149 "Численное моделирование распространения сейсмических и акусто-гравитационных волн для неоднородной модели Земля-Океан-Атмосфера" (руководитель Б.Г.Михай-ленко);

- 2007-2009, проект РФФИ № 07-05-00538 "Изучение проявления скоплений неоднородностей субсейсмического масштаба в волновых сейсмических полях" (руководитель Г.В. Решетова);

- 2010-2012, проект РФФИ № 10-05-00233,"0т микроструктуры к макромодели трещиновато-кавернозных резервуаров: гомогенизация, проявление флюидонасыщенности, описание рассеянных волн" (руководитель Г.В. Решетова).

Внедрение научных результатов

Созданный диссертантом программный продукт по моделированию двумерно-неоднородных сложнопостроенных упругих сред был внедрен в мобильный программно-аппаратный комплекс в рамках комплексного проекта "Разработка комплексной технологии поиска и разведки углеводородов в сложнопостроенных, глубокозалегающих месторождениях" по Госзаказу 2005-РИ-00.0/009/202, шифр - 1Ш.1ВММС.00103-01 34 01 (Акт о внедрении от 25.11.2006).

Разработанное автором программное обеспечение для моделирования трёхмерных сейсмоакустических волновых полей (при акустическом каротаже) в средах с затуханием в настоящее время используется инженерными технологическими центрами компании БсЫитЬе^ег (в рамках контракта с ИНГГ СО РАН, Акт о внедрении от 01.08.2009).

Разработанное диссертантом программное обеспечение для моделирования сейсмических волновых полей в трехмерно-неоднородных средах с субсейсмическими неоднородностями (кавернозно-трещиноватые коллекторы) используется ООО "Геола" с целью изучения проявлений ориентации коридоров трегциноватости рассеянных волнах (Акт о внедрении от 14.09.2010).

Апробация работы и публикации

Результаты диссертационной работы известны научной общественности. Всего по теме диссертации автором лично и в соавторстве опубликовано более 100 работ, в том числе 38 статей, из которых 15 - в ведущих рецензируемых научных журналах из перечня ВАК. Результаты диссертационной работы докладывались и обсуждались на Международных научных конференциях в России и за рубежом. Основные:

- Генеральные Ассамблеи Европейского геофизического общества EGU, Франция, Ницца, 2000, 2002, Австрия, Вена, 2004;

- Международный конгресс математиков, Китай, Пекин, 2002;

- 62-я Коференция EAGE, Великобритания, Глазго, 2000;

- Генеральная Ассамблея ассоциации сейсмологии и физики Земли IASPEI, Вьетнам, Ханой, 2001;

- 65-я Конференция EAGE, Норвегия, Ставангер, 2003;

- Международные конференции "Мониторинг ядерных испытаний и их последствий Казахстан, Боровое, 2004, 2008;

- Международная конференция по вычислительной математике ICCM, Россия, Новосибирск, 2004;

- 2-й Международный симпозиум "Активный геофизический мониторинг литосферы Земли", Россия, Новосибирск,. 2005;

- Ежегодная международная конференции сообщества геофизиков-исследователей (SEG), США, Новый Орлеан, 2006;

- Международная конференция "Науки о Земле - открыть и разработать Россия, Санкт-Петербург, 2006;

- Всероссийская конференция по вычислительной математике КВМ, Россия, Новосибирск, 2007;

- 69-я Международная конференция EAGE, Великобритания, Лондон, 2007;

- 8-я Международная конференция "Математические и численные аспекты теории распространения волн Великобритания, Рединг, 2007;

- 7-я Международная европейская конференция по вычислительной математике, Австрия, Грац, 2007;

- Десятая международная научно-практическая конференция "Геомодель-2008", Россия, Геленджик, 2008;

- Международная конференция геофизиков и геологов "К эффективности через сотрудничество Россия, Тюмень, 2007;

- Международная конференция "ГеоСибирь - 2008 Россия, Новосибирск, 2008;

- 70-я Международная конференция EAGE, Италия, Рим, 2008;

- Международная конференция "Математическое моделирование в геофизике Россия, Новосибирск, 2008;

- 9-я Международная конференция "Математические и численные аспекты теории распространения волн Франция, По, 2009;

- 72-я Международная конференция EAGE, Испания, Барселона, 2010;

- Международная конференция EAGE, Россия, Санкт-Петербург, 2010.

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

Успешному проведению исследования способствовала поддержка академиков РАН А.Н. Коновалова и Б.Г. Михайленко, оказавших большое влияние на формирование научных взглядов соискателя.

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

Автор ценит всестороннюю поддержку, постоянное внимание к работе и благодарит всех сотрудников Лаборатории численного моделирования сейсмических полей Института вычислительной математики и математической геофизики Сибирского отделения РАН, а также сотрудников Лаборатории вычислительных методов геофизики Института нефтегазовой геологии и геофизики им. A.A. Трофимука Сибирского отделения РАН, особенно В.В.Лисицу и Е.В. Лыся.

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

ЗАКЛЮЧЕНИЕ

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

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

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

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

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

2. С использованием конечно-разностных схем с локальным пространственно-временным измельчением сеток разработан, теоретически и экспериментально обоснован и программно реализован новый метод численного моделирования волновых полей в трёхмерно-неоднородных разномасштабных упругих средах. Выгодным отличием пространственно-временного измельчения сеток от известных ранее способов является:

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

- использование для согласования сеток пространственной интерполяции на основе преобразования Фурье, что обеспечило её экспоненциальную точность и низкий уровень возникающих при этом артефактов (менее 0.1% от амплитуды падающей волны);

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

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

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

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

- возможность ограничения расчетной области устойчивыми РМЬ без расщепления;

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

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

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

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

1. Алексеев A.C., Бабич В.М., Гельчииский Б.Я. Лучевой метод вычисления интенсивности волновых фронтов // Вопросы динамической теории распространения сейсмических волн. Л.: Наука -1961. -№5. -С.3-24.

2. Алексеев A.C., Гельчинский Б.Я. О лучевом методе вычисления полей волн в случае неоднородных сред с криволинейными границамиураздела // Вопросы динамической теории распространения сейсмических волн. Л.: Наука-1959. -№3. -С.16-47.

3. Алексеев A.C., Глинский Б.М. и др. Эффект акусто-сейсмической индукции при вибросейсмическом зондировании // Докл. РАН. -1996. -Т.346. -№5. -С.664-667.

4. Алексеев A.C., Михайленко Б.Г. Решение задач Лэмба для вертикально-неоднородного полупространства // Изв. АН СССР, сер. Физика Земли. -1976. -№12. -С. 11-25.

5. Алексеев A.C., Михайленко Б.Г. Численное моделирование распространения сейсмических волн в радиально-неоднородной модели Земли // Докл. АН СССР. -1977. -Т.235. 1. -С.46-49.

6. Алексеев A.C., Михайленко Б.Г. "Нелучевые" эффекты в теории распространения сейсмических волн // Докл. АН СССР. -1982. -Т.267. -С.1079-1083.

7. Бабич В.М., Алексеев A.C. Лучевой метод вычисления интенсивности волновых фронтов // Изв. АН СССР. Сер. Геофизическая. -1958. -т. -С.9-15.

8. Бахвалов Н.С. Численные методы. -М.: Наука, 1973. -Т.1. -632с.

9. Бахвалов H.С., Жидков H.П., Кобельков Г.М. Численные методы. -М.: Наука, -1987. -598с.

10. Бреховских JI.M., Годин O.A. Акустика слоистых сред. -М.: Наука, 1989. -214с.

11. Вишневский Д.М., Костин В.И., Чеверда В.А. Возбуждение сейсмических волн источником, расположенным в скважине, заполненной жидкостью // Физическая мезомеханика. -2002. -Т.5. -К2 5. -С.85-92.

12. Гасилова JI.A., Петухов Ю.В. К теории поверхностных волн, распространяющихся вдоль разных границ раздела в атмосфере // Известия РАН. Физика атмосферы и океана. -1999. -Т.35. -№1. -С.14-23.

13. Гик Л.Д. Определение коэффициента пористости в методе акустического каротажа // Наука и технология углеводородов. -2002. -№3. -С.3-11.

14. Гилбо Ж., Ланда Е., Решетова Г.В., Хайдуков В.Г., Чеверда В.А. Численное моделирование сейсмических волновых полей в двумерно-неоднородных упругих разномасштабных средах (карстовые включения) // Технологии сейсморазведки. -2008. -Т.З. -С. 19-28.

15. Гольдин C.B., Поздняков В.А., Смирнов М.Ю., Чеверда В.А. Построение сейсмических изображений в рассеянных волнах как средство детализации сейсмического разреза // Геофизика, Специальный выпуск к 40-летию ОАО "Тюменнефтегеофизика". -2004. -С.23-29.

16. Ельцов И.Н., Кашеваров A.A., Решетова Г.В., Чеверда В.А. Проявление неоднородностей зоны проникновения в геофизических полях вдоль ствола скважины // Геофизика. -2004. -Т.б. -С.17-21.

17. Ельцов И.Н., Решетова Г.В., Чеверда В.А. Численное моделирование процессов распространения сейсмоакустических полей с учетомнеоднородности зоны проникновения // Физическая мезомеханика. -2005. -Т.8. -С.99-105.

18. Зиатдинов С.Р., Каштан Б.М., Бакулин A.B., Троян В.Н. Влияние конечной трещины на распространение трубной волны // Геофизика. -2006. -Т.5. -С.22-27.

19. Исимару А. Распространение и рассеяние волн в случайно-неоднородных средах. М.: Мир, 1981. -Т.1. -280с.

20. Кашеваров A.A., Ельцов И.Н., Эпов М.И. Гидродинамическая модель формирования зоны проникновения при бурении скважин // ПМТФ. -2003. -Т.44. -т. -С. 148-157.

21. Коновалов А.Н. Решение задач теории упругости в напряжениях: Учеб. пособие./ Новосибирский гос. ун-т, 1979. -92с.

22. Коновалов А.Н. Введение в вычислительные методы линейной алгебры. -Новосибирск: Наука. Сибирское отделение, 1993. -158с.

23. Коновалов А.Н. Динамическая задача теории упругости в постановке "скорости-напряжения" // Дифференциальные уравнения. -1999. -Т.35. -т. -С.238-248.

24. Конюх (Решетова) Г.В., Михайленко Б.Г. Применение интегрального преобразования Лагерра при решении динамических задач сей-смики // Труды ИВМиМГ СО РАН. Сер. Математическое моделирование в геофизике. -1998. -№5. -С.106-123.

25. Конюх (Решетова) Г.В., Мнхайленко Б.Г., Михайлов A.A. Численное моделирование сейсмических полей в вязкоупругих средах на основе спектрального метода Лагерра // Математическое моделирование. -2001. -Т.13. -№2. -С.61-70.

26. Костин В.И., Решетова Г.В., Чеверда В.А. Численное моделирование трехмерного акустического каротажа с использованием многопроцессорных вычислительных систем // Математическое моделирование. -2008. -Т.20. -т. -С.51-66.

27. Крауклис П.В., Крауклис Л.А. Волновое поле точечного источника в скважине // Вопросы динамической теории распространения волн. -Л.: Наука, 1976. -Вып. XVI. -С.41-53.

28. Крауклис П. В., Крауклис Л. А. О спектре продольной волны в скважине с зацементированной обсадной колонной // Вопросы динамической теории распространения сейсмических волн. -Л.: Наука, 1976. -Вып. 17. -С.156-164.

29. Крукиер Л.А., Лапшина O.A. Численное сравнение вариационных методов решения СЛАУ, получаемых при копечно-разностной аппроксимации уравнения конвекции-диффузии // Математическое моделирование. -2004. -Т. 16. -С.23-32.

30. Кузнецов О.Л., Курьянов Ю.А., Чиркин H.A., Шленкин С.И. Сейсмический локатор бокового обзора / / Геофизика. Спецвыпуск к 40-летию ОАО "Тюменьнефтегеофизика". -2004. -С.17-22.

31. Курант Р., Гильберт Д. Методы математической физики. М.-Л.: ГИТТЛ, 1951. -Т.1. -476с.

32. Левянт В.Б., Тропов Ю.А., Шустер В.Л. Использование рассеянной компоненты сейсмического поля для дифференциации кристаллического фундамента на коллекторские и монолитные зоны // Геофизика. -2003. -№3. -С.17-26.

33. Левянт В.Б., Антоненко М.Н., Антонова И.Ю. Исследование методами численного моделирования сейсмического поля, обусловленного рассеиванием на зонах диффузной кавернозности и трещинова-тости // Геофизика. -2004. -№2. -С.8-20.

34. Линьков Е.М., Петрова А.Н., Осипов К.С. Сейсмогравитационные пульсации Земли и возмущения атмосферы как возможные предвестники сильных землетрясений // Докл. АН СССР. -1990. -Т.313. -№5. -С.1095-1098.

35. Макаров А.Н., Багринцева К.И., Садыбеков Д.Т. Особенности строения рифейских карбонатных коллекторов Юрубченского газонефтяного месторождения // Геология нефти и газа. -1998. -№4. -С. 1830.

36. Михайленко Б.Г. Численное решение задачи Лэмба для неоднородного полупространства // Математические проблемы геофизики. ВЦ СО РАН, Новосибирск. -1973. -С.273-297.

37. Михайленко Б.Г. Расчет теоретических сейсмограмм для многомерно-неоднородных моделей сред // Условно-корректные задачи математической физики в интерпретации геофизических наблюдений. ВЦ СО РАН, Новосибирск. -1978. -С.75-88.

38. Михайленко Б.Г. Метод решения динамических задач сейсмики длядвумерно-неоднородных моделей сред // Докл. АН СССР. -1979. -Т.246. -М. -С.47-51.

39. Михайленко Б.Г. Сейсмические поля в сложнопостроенных средах // Новосибирск: Изд-во СО РАН. -1988. -311с.

40. Михайленко Б.Г., Решетова Г.В. Численно-аналитический метод решения задачи о распространении сейсмических и акусто-гравитаци-онных волн для неоднородной модели Земля-Атмосфера // Сибирский журнал вычислительной математики. -2006а. -Т.9. -С.37-46.

41. Михайленко Б.Г., Решетова Г.В. Математическое моделирование распространения сейсмических и акусто-гравитационных волн для неоднородной модели Земля-Атмосфера // Геология и геофизика. -20066. -Т.47. -№5. -С.547-556.

42. Молотков J1.A. Матричный метод в теории распространения волн в слоистых упругих и жидких средах. -JL: Наука, -1984. -270с.

43. Никифоров А.Ф., Уваров В.Б. Основы теории специальных функций. -М.: Наука, 1974. -303с.

44. Поздняков В.А., Чеверда В.А. Фокусирующие преобразования сейсмических данных для площадных стационарных систем / / Геология и геофизика. -2005. -№3. -С.328-338.

45. Прудников А.П., Брычков Ю.А., Марпчев О.И. Интегралы и рады: Специальные функции. М.: Наука, 1983. - 749 с.

46. Разин A.B. Распространение сферичного акустического дельта-импульса вдоль границы газ-твердое тело // Известия АН. Физика Земли. -1993. -№2. -С.73-77.

47. Решетова Г.В., Чеверда В.А. Использование преобразования Лагер-ра для построения идеально подходящих поглощающих слоев без расщепления // Математическое моделирование. -2006. -Т.18. -№10. -С.91-101.

48. Рытов С.М., Кравцов Ю.А., Татарский В.И. Введение в статистическую радиофизику. Ч. 2. Случайные поля. -М.: Наука, 1978. -464с.

49. Самарский A.A. Теория разностных схем. -М.: Наука, 1983, -616с.

50. Суетин П.К. Классические ортогональные многочлены. -М.: Наука, 1974. -327с.

51. Технология исследования нефтегазовых скважин на основе ВИКИЗ // Методическое руководство. Ред. Эпов М.И., Антонов Ю.Н. Составители: Ельцов И.Н., Жмаев С.С., Петров А.Н. и др. -Новосибирск: НИЦ ОИГГМ СО РАН, Изд-во СО РАН, -2000, -121с.

52. Фатьянов А.Г., Михайленко Б.Г. Численное решение задачи Лэм-ба для неоднородной среды Больцмана с упругим последействием // Математические методы интерпретации геофизических наблюдений. -Новосибирск. -1979. -С. 115-160.

53. Фатьянов А.Г., Михайленко Б.Г. Метод расчета нестационарных волновых полей в неупругих слоисто-неоднородных средах // Докл. АН СССР. -1988. -Т.301. -Ш. -С.1024-1027.

54. Шемякин Е.И. Динамические задачи теории упругости и пластичности. -Новосибирск: НГУ, 1968. -337с.

55. Шемякин Е.И. Введение в теорию упругости.-М.: МГУ, 1993. -96с.

56. Abarbanel S., Gottlieb S., Hesthaven J.S. Long time behavior of the perfectly matched layer equations in computational electromagnetics // Journal of Scientific Computing. -2002. -V.17. -P.405-422.

57. Alford R.M., Kelly K.R., Boore D.AL Accuracy of finite-difference modeling of the acoustic wave equation // Geophysics. -1974. -V.39. -P.834-842.

58. Alekseev A.S., Mikhailenko B.G. The solution of dynamic problems of elastic wave propagation in inhomogeneous media by a combination of partial separation of variables and finite-difference method // Geophysics. -1980. -V.48. -P.161-172.

59. Alterman Z., Karal F.G. Propagation of elastic waves layered media by finite-difference methods// Bull. Seism. Soc. Am. -1968. -V.58. -P.367-398.

60. Angerer E., Lanfranchi P., Rogers S.F. Fractured reservoir modeling from seismic to simulator: A reality? // The Leading Edge. -2003. -P.684-689.

61. Artru J. Ducic V. and others. Ionospheric detection of gravity waves induced by tsunamis // Geophys. J. Int. -2005. -V.160. -P.840-848.

62. Asvadurov S., Knizhnerman L. and Pabon, J. Finite-difference modeling of viscoelastic materials with quality factors of arbitrary magnitude // Geophysics. -2004. -V.69. -№.3. -P.817-824.

63. Asvadurov S., Knizhnerman L. and Pabon J. Finitedifference modeling of viscoelastic materials with quality factors Q of arbitrary magnitude // Schlumberger-Doll Research preprint OFSR/ RN/ 2001 /151 /RTMI /. -2001. -151p.

64. Bamberger A., Glowinski R., Tran Q.H. A domain decomposition methodfor the acoustic wave equation with discontinuous coefficients and grid change // SIAM J. Num. Anal. -1997. -V.34. -№.2. -P.603-639.

65. Bayliss A. et.al. A fourth order accurate finite-difference scheme for the computation of clastic waves // Bull. Seis. Soc. Am. -1986. -V.76. -P.1115-1132.

66. Becache E., Petropoulos P.G., Gedney S.G. On the long-time behavior of unsplit perfectly matched layers // INRIA Rapport de Recherchc, -2002. -№4538. -18p.

67. Berenger J.P. A perfectly matched layer for the absorption of electromagnetic waves // J. Comput. Phys. -1994. -V.114. -P.185-200.

68. Biot M.A. Propagation of elastic waves in a cylindrical bore containing a fluid // Journal of Applied Physics. -1952. -V.23. -P.997-1005.

69. Bland D.R. The theory of linear viscoelasticity. -Pergamon press, -1960. -125p.

70. Blanch J., Robertsson 0., and Symes W. Modeling of a constant Q: Methodology and algorithm for an efficient and optimally inexpensive viscoelastic technique // Geophysics. -1995. -V.60. -P.176-184.

71. Bohlen T. Parallel 3-D viscoelastic finite-difference seismic modeling // Computers & Geoscients. -2002. -№.28. -P.887-899.

72. Boore D.M. Finite-difference methods for seismic wave propagation in heterogeneous materials // Methods in Computation Physics. -1972 -V.ll. -P. 1-37.

73. Carcione J.M., Wang J.P. A Cheb^yshev collocation method for the elastodynamic equation in generalized coordinates // Comput. Fluid Dynamics. -1993. -V.2. -P.269-290.

74. Carcione J.M., KoslofF D. and Kosloff R. Wave propagation simulation in a linear viscoelastic medium // Geophys. J. R. Astr. Soc. -1988. -V.95. -P.597-611.

75. Carcione J., Kosloff D. and Kosloff R. Wave propagation simulation in a linear viscoacoustic medium // Geophysical Journal of the Royal Astronomical Society. -1988. -V.93. -P.393-407.

76. Cerveny V., Ravindra R. Theory of seismic head waves // Toronto Univ. Press, Toronto, 1971. -250p.

77. Cerveny V., Molotkov I.A., Psencik I. Ray method in seismology // Prague, Varlovar. Univ., 1977. -281p.

78. Chew W. C., Liu Q. H. Perfectly matched layers for elastodynamics: A new absorbing boundary condition // Journal of Computational Acoustics. -1996. -V.4. -№.4. -P.341-359.

79. Christensen R.M. Theory of viscoelasticity An introduction // Academic Press, New-York, seconde edition, 1982. -114p.

80. Collino Francis. Perfectly matched absorbing layers for the paraxial equations // Journal of Computational Physics. -1997. -V.131. -№.1. -P.164-180.

81. Collino F., Tsogka C. Application of PML absorbing layer model to the linear elastodynamic problem in anisotropic heterogeneous media // Geophysics. -2001. -V.66. -P.294-307.

82. Collino F., Fouquet T., Joly P. Analyse numerique d'une methode de raffinemeiit de maillage espace-temps pour l'equation des ondes // Technical Report 3474, INRIA, Aout, 1998.

83. Collino F., Fouquet T., Joly P. A conservative space-time mesh refinement method for the 1-D wave equation. Part I: Construction // Numer.

84. Math. -2003a -V.95. -P. 197-221.

85. Collino F., Fouquet Т., Joly P. A conservative space-time mesh refinement method for the 1-D wave equation. Part II: Analysis. Numer. Math. -20036 -V.95. -P.223-251.

86. Dablain M.A. The application of high-order differencing to the scale wave equation // Geophysics. -1986. -V.51. -P.54-56.

87. Day S. and Minster J. Numerical simulation of attenuated wavefields using a Padi"e approximant method // Geophysical Journal of the Royal Astrological Society. -1984. -V.78. -P.105-118.

88. Diaz J., Grote M.J. Energy conserving explicit local time-stepping for second-order wave equation // INRIA Research Report, 2007. -№.6377.

89. Dong Z., McMechan G.A. 3-D viscoelastic anisotropic modeling of data from a multi-component, multi-azimuth experiment in northeast Texas // Geophysics. -1995 -V.60. -P.1128-1138.

90. Emmerich H., and M. Korn. Incorporation of attenuation into timedomain computations of seismic wave fields // Geophysics. -1987. -V.52. -P.1252-1264.

91. Engquist В., Majda A. Absorbing Boundary Conditions for Numerical Simulation of Waves // Proceedings of the National Academy of Sciences of the United States of America. -1977. -V.74. -№.5. -P. 1765-1766.

92. Francis S.H. Global propagation of atmospheric gravity waves: A review // J. Atmos. Terr. Phys. -1975. -V.37. № 6-7. -P.1011-1054.

93. Fuchs K., Muller G. Computation of synthetic seismograms with the reflectivity method and comparison with observations // Geophys. J. R. Astr. Soc. -1971. -V.23. -P.417-433.

94. Futterman W.I. Dispersive body waves //J. Geophys. Res. -1962. -V.67. -P.5279-5291.

95. Gossaro E. Hooken. Waves in the Atmosphere // Elsevier Scientific Publishing Company. -1975. 532p.

96. Graves R.W. Simulating seismic wave propagation in 3d elastic media using staggered grid finite difference // Bulletin of seismic Society of America. -1996. -V.86. -P.1091-1106.

97. Gurtin M.E. and Sternberg E. On the linear theory of visscoelasticity // Arch. Ration. Mech. Anal. -1962. -V.ll. -P.291-356.

98. Haskell N.A. The dispersion of surface waves in multi-layered media // Bull. Seism. Soc. Am. -1953. -V.43. -P.17-43.

99. Heiner I. Wave propagation in three-dimensional spherical sections by the Chebyshev spectral method // Geophys. J. Int. -1999. -V.136. -P.559-566.

100. Hestholm S. et al. Quick and accurate Q parameterization in viscoelastic wave modeling // Geophysics. -2006 -V.71. -doi:10.1190/1.2329864.

101. Holdberg O. Computational aspects of the choice of operator and sampling interval for numerical differentiation in large-scale simulation of wave phenomena // Geophys. Prosp. -1987. -V.35. -P.629-655.

102. Hron F. Introduction to the ray theory in a broader sense: application to seismology // Textbook of Laboratorie de Physique de l'Écola Normale Supérieure, Université de Paris, Paris, 1968. -180p.

103. Hron F., Mikhailenko B.G. Numerical modeling of nongeometrical effects by the Alekseev-Mikhailenko method // Bull. Soc. Am. -1981. -V.71. №.4. -P.10011-1099.

104. Hu F.Q. A stable, perfectly matched layer for linearized Euler equations in unsplit physical variables// J. Comp. Physics. -2001. -V.173. -P.455-480.

105. Hustedt B., Operto S., Virieux J. Mixed-grid and staggered-grid finite-difference methods for frequency-domain acoustic wave modelling // Geophysical Journal International. -2004. -V.157. -P. 1269-1296.

106. Jo C.H., Shin C.S., Suh J.H. An optimal 9-point., finite-difference, frequency-space, 2D scalar wave extrapolator // Geophysics. -1996. -V.61. -P.529-537.

107. Karal F.C., Keller J.B. Elastic wave propagation in homogeneous and inhomogeneous media //J. Acoust. Soc. Am. -1959. -V.31. -P.694-705.

108. Kelly K.R. et. all. Synthetic seismograms: a finite-difference approach // Geophysics. -1976. -V.41. -P.2-27.

109. Kim S.I., Hoefer W.J.R. A local mesh refinement algorithm for the timedomain finite-difference method to solve Maxvell's equations // IEEE Trans. Microwave Theory Tech. -1990. -V.38. -P.812-815.

110. Konyukh G.V. and Mikhailenko B.G. Application of the Laguerre integral transform for solving dynamic seismic problems //Bull, of the Novosibirsk Computing Center, Series: Mathematical Modeling in Geophysics, Novosibirsk. -1998c. -V.4. -P.79-91.

111. Konyukh G.V., Krivtsov Y.I., Mikhailenko B.G. Numerical-Analytical Modeling of Seismic Wave Propagation in Vertically Inhomogeneous Media // Applied Mathematics Letters. -1998. -V.ll. -№1. -P.99-104.

112. Konyukh G.V., Mikhailenko B.G., Mikhailov A. A. Integral Laguerre transform as applied to forward seismic modeling // Bull, of the Novosibirsk Computing Center, series Mathematical Modeling in Geophysics, Novosibirsk. -1999. -V.5. -P.71-91.

113. Konyukh G.V., Mikhailenko B.G. and Mikhailov A.A. Application of the integral Laguerre transforms for forward seismic modeling // Journal of Computational Acoustics. -2001. -V.9. №.4. -P.1523-1541.

114. Kosloff D., Baysal E. Forward modeling by a Fourier method // Geophysics. -1982. -V.47. -P. 1402-1412.

115. Kosloff D., et. al. Solution of equations of dynamic elasity by a Chebyshev spectral method // Geophysics. -1990. -V.55. -P.734-748.

116. KoslofT D., Reshef M., Loewenthal D. Elastic wave calculations by the Fourier method // Bull. Seis. Soc. Am. -1984. -V.74. -P.875-891.

117. Kostin V.I., Pissarenko D.V., Reshetova G.V., Tcheverda V.A. Numerical Simulation of 3D Acoustic Logging // Applied Parallel Computing. State of the Art in Scientific Computing. Lecture Notes in Computer Sciences. -2007. -V.4699. -P. 1045-1054.

118. Kovalevsky V.V., Reshetova G.V. Simulation of superpower shaft hydro-resonance vibrational source // Bull, of the Novosibirsk Computing Center, series Mathematical Modeling in Geophysics, Novosibirsk. -2002. -V.7. -P.45-57.

119. Levander A.R. Fourth-order finite-difference P-SV seismograms // Geophysics. -1988. -V.53. -P. 1425-1436.

120. Liu H.P., Anderson D.L., Kanamori H. Velocity dispersion due to anelas-ticity; implications for seismology and mantle composition // Geophys. J. R. astr. Soc. -1976. -V.47. -P.41-58.

121. Liu E., Maultzsch S., Chapman M., Li X.Y., Queen J.H., Zhang Z. Frequency-dependent seismic anisotropy and its implication for estimating fracture size in low porosity reservoirs // The Leading Edge. -2003. -P.662-665.

122. Levander A.R. Fourth order velocity-stress finite-difference scheme // Proc. 57th SEG Annual Meeting, New Orleans. -1987. -P.234-245.

123. Levin F.K. Seismic velocities in transversely isotropic media // Geophysics. -1979. -V.44. -P.918-936.

124. Leviant V.B., Petrov I.B., Chelnokov F.B., Antonova I.Y. Nature of the scattered seismic response from zones of random clusters of cavities and fractures in a massive rock // Geophysical Prospecting. -2007. -V.55. -P.507-524.

125. Lysmer B., Drake N. A finite-element method for seismology // In Bolt B.A. Eds., Methods in computational physics, Seismology: Surface waves and Earth oscillations, Academic Press Inc., 1972. -V.ll. -P.181-216.

126. Marfurt K.J. Accuracy of finite-difference and finite-element modeling of the scalar and elastic wave equations // Geophysics. -1984. -V.49. -P.533-549.

127. Marfurt K.J., Shin C.S. The future of iterative modeling in geophysical exploration //In Eisner E. Ed., Handbook of geophysical exploration: Seismic exploration, supercomputers in seismic exploration, Pergamon Press. -1989. -V.21. -P.203-228.

128. McDonal F.J., Angona F.A., Mills R.L., Sengbush R.L., van Nostrand R.G., White J.E. Attenuation of shear and compressional waves in Pierre shale // Geophysics. -1958. -V.23. -P.421-439.

129. Mikhailenko B.G. Spectral Laguerre method for the approximation solution of time-dependent problems // Appl. Math. Lett. -1999. -V.12. -P.105-110.

130. Mikhailenko B.G. Synthetic seismograms for complex 3D geometries using an analytical-numerical algorithm // Geophys. J.R. Astr. Soc. -1984. -V.79. -P.963-986.

131. Mikhailenko B.G. Numerical experiment in seismic investigations // Geophys. -1985. -V.58. -P.101-124.

132. Mikhailenko B.G. Elastic wave propagation in complex 3D media // The 61st Annual Internal. Mtg. Soc. Expl. Geophys., Expanded Abstracts. -1991 -P.1569-1572.

133. Mikhailenko B.G. Seismic modeling by the spectral-finite difference method // Physics of the Earth and Planetary Interiors. -2000. -V.119. -P. 133-147.

134. Mikhailenko B.G. Spectral Laguerre method for the approximate solution of time dependent problems // Appl. Math. Lett. -1999. -V.12. -P. 105110.

135. Mikhailenko B.G., Mikhailov A.A. and Reshetova G.V. Numerical Modeling of transient seismic fields in viscoelastic media based on the Laguerre spectral method // Pure Appl. Geophys. -2003. -V.160. -P. 1207-1224.

136. Mikhailenko В., Reshetova G. Application of the perfectly matched absorbing layer method for the spectral Laguerre modeling // Труды международной конференции по вычислительной математике, МКВМ-2004, Новосибирск. -2004а. -Р.901-906.

137. Mikhailenko B.G., Mikhailov А.А., Reshetova G.V. Numerical viscoelastic modeling by the spectral Laguerre method // Geophysical Prospecting -2003. -V.51. -P.37-48.

138. Mikhailenko B.G., Mikhailov A.A., Reshetova G.V.Numerical modeling of transient seismic fields in viscoelastic media based on the Laguerre spectral method // Pure and Applied Geophysics. -2003. -V.160. -P. 12071224.

139. Mikhailenko В., Reshetova G. The perfectly matched layers technique in 3D seismic modeling // Proceedings of International Symposium on Mathematical Modeling of Dynamic Processes in Atmosphere, Ocean, and Solid Earth, Novosibirsk. -20046. -P.91-97.

140. Mikhailenko B.G., Reshetova G.V. Mathematical modeling of seismic and acousto-gravitational waves in a Heterogeneous earth-atmosphere model // Journal of Computational and Applied Mathematics. -2009. -V.234. -№6. -P.1678-1684.

141. Mikhailenko B.G., Reshetova G.V. Mathematical modeling of seismic and acousto-gravitational waves in a Heterogeneous earth-atmospheremodel // Journal of Computational and Applied Mathematics. -2010. -V.234. -P. 1678-1684.

142. Muller G. Rheological properties and velocity dispersion of a medium with power-law dependence of Q on frequency // Geophysics. -1983. -V.54. -P.20-29.

143. Muller G. The reflectivity method; a tutorial // Geophysics. -1985. -V.58. -P. 153-174.

144. Pissarenko D.V., Reshetova G.V., Tcheverda V.A. 3D finite-difference synthetic acoustic logging in cylindrical coordinates // Geophysical Prospecting. -2009. -V.57. -№3. -P.367-377.

145. Pissarenko D.V., Reshetova G.V., Tcheverda V.A. 3D finite-difference synthetic acoustic log in cylindrical coordinates // Journal of Computational and Applied Mathematics. -2010. -V.234. -P.1766-1772.

146. Pratt R.G., Shin C., Hicks G.J. Gauss-Newton and full Newton methods in frequency-space seismic waveform inversion // Geophys. J. Int. -1998. -V.133. -P.341-362.

147. Robertsson J., Blanch J., Symes W. Viscoelastic finite-difference modeling // Geophysics. -1999. -V.64. -P.1444-1456.

148. Ruger A. Variation of P-wave reflectivity with offset and azimuth in anisotropic media // Geophysics. -1997. -V.63. -№3. -P.935-947.

149. Ruger A., Tsvankin I. Using AVO for fracture detection: Analytic basisand practical solution // The Leading Edge. -1997 October. -P. 14291434.

150. Saad Y. Iterative methods for sparse linear systems // PWS Publishing Company. -1995. -lOOp.

151. Saad Y., Van der Vorst H.A. Iterative solution of linear systems in the 20th century // J. of Computational and Applied Mathematics. -2000. -№123. -P.l-33.

152. Schoenberg M., Sayers C.M. Seismic anisotropy of fractured rock // Geophysics. -1995 -V.60. -№1. -P.204-211.

153. Siddiqui S.A. Dispersion analysis of seismic data // M.S. thesis, Univ. of Tulsa. -1971.

154. Singh S.K., Abu-Habbiel H., Khan B., Akbar M., Etchecopar A., Monatron B. Mapping fracture corridors in naturally fractured reservoirs: an example from Middle East carbonates // First Break. -2008 May. -P. 109113.

155. Sonneveld P. CGS, a fast Lanczos-type solver for nonsymmetric linear system // SIAM Journal of Scientific and Statistical Computing. -1989. -V.10. -P.36-52.

156. Stekl I., Pratt R.G. Accurate viscoelastic modeling by frequency-domain finite-difference using rotated operators // Geophysics. -1998. -V.63. -P. 1779-1794.

157. Tal-Ezer H., Kosloff D., Koren Z. An accurate scheme for seismic forward modeling // Geophys. Prosp. -1987. -V.35. -P.479-490.

158. Tal-Ezer H., Carcione J.M, Kosloff D. An accurate and efficient scheme for wave propagation in linear viscoelastic media // Geophysics. -1990. -V.55. -P.1366-1379.

159. Thomsen S. Weak elastic anisotropy // Geophysics. -1986. -V.51. -№10. -P. 1954-1966.

160. Thomson W.T. Transmission of elastic waves through a stratified solid // J. Appl. Phys. -1950. -V.21. -P.89-93.

161. Trefethen L. and HalpernL. Well-posedness of one-way wave equations and absorbing boundary conditions // Math. Comput. -1986 -V.47. -P.421-435 .

162. Virieux J. P-, SV- wave propagation in heterogeneous media: velocity-stress finite-difference method // Geophysics. -1986. -V.51. -P.889-901.

163. Wuenschel P.C. Dispersive body waves an experimental study // Geophysics. -1965. -V.30. -P.539-551.

164. Yi W., Nakagawa S., Nihei K.T., Rector J.W., Myer L.R., Cook N.G.W. Numerical investigation of fracture-induced seismic anisotropy // SEG Expanded Abstracts. -1997. -P.960-963.