Численное решение пространственной динамической задачи теории упругости на многопроцессорных вычислительных системах тема диссертации и автореферата по ВАК РФ 05.13.18, кандидат физико-математических наук Кучунова, Елена Владимировна

• Кучунова, Елена Владимировна
• кандидат физико-математических наук
• 2008

Оглавление диссертации кандидат физико-математических наук Кучунова, Елена Владимировна

Введение

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

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

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

1.3 Обзор методов моделирования неотражающих граничных условий

1.4 Обзор основных тенденций в области параллельных вычислений

II Численное решение задачи о распространении упругих волн в блочной среде

2.1 Постановка пространственной задачи.

2.2 Построение сетки в области.

2.2.1 Разбиение области на блоки.

2.2.2 Нахождение направляющих векторов.

2.2.3 Построение сетки в блоке.

2.3 Численный метод решения задачи.

2.3.1 Метод двуциклического расщепления.

2.3.2 Решение одномерных систем

2.3.3 Предельная реконструкция инвариантов решения

2.3.4 Склейка решения на внутренней границе двух блоков

III Технологии реализации параллельного алгоритма на многопроцессорных вычислительных системах

3.1 Проблема эффективной реализации алгоритмов па многопроцессорных вычислительных системах.

3.2 Сравнение MPI- и DVM-технологий распараллеливания на примере решения одномерной задачи теории упругости.

3.3 Технология распределения данных по процессорам

3.4 Обмен граничными значениями.'

3.4.1 Обмен значениями внутри одного блока.

3.4.2 Обмен значениями между блоками.

3.5 Программный комплекс.

3.5.1 Основные идеи программного комплекса

3.5.2 Параллельная реализация алгоритма построения расчетной сетки.

3.5.3 Параллельная реализация алгоритма численного решения задачи

3.5.4 Методы визуализации решения.

3.6 Эффективность параллельной реализации алгоритмов

IV Результаты исследования на модельных задачах

4.1 Вычислительный эксперимент

4.1.1 Модельная задача 1. Действие сосредоточенной нагрузки на упругое полупространство.

4.2 Прохождение волн через различные границы раздела сред

4.2.1 Модельная задача 2. Прохождение волн через наклонную границу раздела двух сред.

4.2.2 Модельная задача 3. Прохождение волн через горизонтальную границу раздела двух сред.

4.2.3 Модельная задача 5. Прохождение волн через границу раздела двух сред со сбросом.

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

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

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

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

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

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

• ускорение расчетов при увеличении числа процессоров;

• необходимость синхронизации при остановке (завершении) итераций;

• снижение затрат времени на синхронизацию;

• периодический сбор данных для записи и анализа;

• обмен данными между процессорами для продолжения расчета;

• снижение объема передаваемой информации при обмене;

• снижение числа актов приема-передачи (то есть числа обменов между процессорами);

• уменьшение влияния размера задачи на ускорение;

• обеспечение масштабируемости и переносимости программы.

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

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

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

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

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

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

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

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

• математическая формулировка задачи,

• построение приближенного (численного) метода решения задачи, написание вычислительного алгоритма;

• программирование на ЭВМ вычислительного алгоритма,

• проведение расчетов на ЭВМ,

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

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

Заключение диссертации по теме «Математическое моделирование, численные методы и комплексы программ», Кучунова, Елена Владимировна

Результаты исследования на модельных задачах

4.1 Вычислительный эксперимент

Вычислительный алгоритм, описанный в главе 2, реализован комплексом параллельных программ на алгоритмическом языке Fortran90 с использованием библиотеки MPI (message passing interface). Реализация параллельного алгоритма была основана на основных идеях, представленных в главе 3. В данной главе представлены результаты численных тестов для различных моделей геофизических сред.

Во всех представленных моделях исходная нагрузка действует нормально на верхнюю грань области в некоторой точке Хо. Функция источника (си-пусообразный импульс) имеет вид: P(t) = Pq sin (аЛ), если t < Т, где Pq -амплитуда и ш - частота генерируемых волн в среде. Верхняя грань области, за исключением зоны приложения импульсной нагрузки, считается свободной от напряжений. На остальных гранях поставлены неотражающие условия, которые соответствуют бесконечной протяженности массива и формулируются для одномерных систем с помощью уравнений на характеристиках.

Расчеты производились на суперкомпьютере МВС-1000/16 и МВС-1000 Института Вычислительного Моделирования СО РАН. Суперкомпьютер МВС-1000/16 введен в эксплуатацию в начале 2002 года. Общая производительность МВС-1000/16 составляет 14 млрд. операций в секунду, скорость пересылок - ЮОМЬ в секунду. Суперкомпьютер МВС-1000, введенный в эксплуатацию в 2005 году, имеет пиковую производительность 220.4 млрд. операций в секунду и производительность по тесту LinPack 134 ГФлоп.

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

4.1.1 Модельная задача 1. Действие сосредоточенной нагрузки на упругое полупространство

Изотропная упругая среда занимает прямоугольную область ft = {х G R3\ 0 < xi < Я, -Н < х2 < Н, -Н <х3< Н} о в декартовой системе координат. Параметры среды: р = 7850 кг/м , ср — 6000 и cs = 3210 м/с. В точке хо = (0,0,0) приложена точечная импульсная нагрузка, действующая перпендикулярно к плоскости xi = 0. Размеры расчетной области1: 200 X 400 X 400 м. Расчеты проводились на сетке из 200 X 400 х 400 узлов. Шаг по времени выбирался в соответствии с услови

1 Здесь и далее размеры области перечисляются в следующем порядке: первая цифра - по xi, вторая -но Х2, третья - по х3. ем Куранта-Фридрихса - Л еви и составил 0.0016 с. Время действия нагрузки Т — 0.005 с, то есть около четырех шагов по времени. Величина амплитуды колебаний Pq — 250 МПа. Среднее время расчета одного шага по времени на 32-х процессорах кластера MB С-1000 - 2 минуты. Количество расчетных ячеек, приходящихся на один процессор, равно 106.

Рассмотрим полученную волновую картину в вертикальном сечении области, проходящем через точку приложения нагрузки хц. На рис. 3 приложения Б приведена сейсмограмма вертикальной компоненты вектора перемещения вдоль профиля, проходящего через точку Х(). На рис. 4 приложения Б представлена визуализация вертикальной компоненты вектора перемещения вдоль вертикального сечения области в фиксированный момент времени 0.016 с. Здесь введены следующие цифровые обозначения: 1 - продольная волна, 2 - поперечная волна, 3 - головная поперечная волна, 4 - поверхностная волна Релея.

Результаты соответствуют точному решению задачи для случая плоского деформированного состояния: видны продольная и поперечная волны, менее интенсивные головные волны и волны Релея, которым на рисунках соответствуют движущиеся точки на поверхности х\ = 0. На рис. 4.1 приведено волновое поле нормального напряжения.

Рис. 4.1: Поверхности уровня напряжения в моменты времени 0.02, 0.03 и 0.04 с

4.2 Прохождение волн черёз различные границы раздела сред

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

4.2.1 Модельная задача 2. Прохождение волн через наклонную границу раздела двух сред

Рассмотрим модель двухслойной среды с однородным изотропным верхним слоем и скачком акустической жесткости на границе с подстилающим полупространством. Под однородной покрывающей средой с параметрами р = 500кг/м , ср=500 и cs=300 м/с распложена вторая среда с параметрами ро = 750 кг/м3, Ср—750 и с°=450 м/с (рис. 4.2.а). Горизонт, разделяющий две среды, представляет собой плоскость, наклоненную по линии АС под углом ipi = 5°, а по линии DB - под углом ip2 = 7°. Заметим, что на границе раздела сред выполняются условия рср ф роср и Pcs ф Pocs: гарантирующие образование отраженных волн.

Размеры расчетной области - АВ — AD = 600, АА2 = 400 м. Глубина залегания наклонного горизонта: АА\ — 175, ВВ\ = 87,5, СС\ = 100, DD\ = 187, 5 м. Внешний импульс на среду действует в точке О. Глубина залегания границы раздела в точке О составляет 00\ = Н = 156,25 м.

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

Величина приложенной нагрузки - Pq=250 МП а, время действия - 0,01 с, радиус действия - 3.5 м. Число узлов сетки - 200 х 200 х 200, временной шаг - т — 0.002 с. Среднее время расчета одного шага на 8-ми процессорах составило 2.63 минуты. Количество расчетных ячеек, приходящихся на один процессор - 106.

На рис. 5-7 приложения Б представлены сейсмограммы вертикальной составляющей вектора перемещения вдоль профилей ОС, ОМ, DB, EF и OQ. Здесь введены следующие обозначения: 1 - падающая Р-волна; 2 - падающая S-волна; 3 - отраженная РР-волна; 4 - отраженная PS-волна; 5 - преломленная РР-волна; 6 - преломленная PS-волна; 7 - головная волна, образованная РР-волиой; 8 - отраженная SS-волна; 9 - отраженная SP-волна.

Исходя из законов геометрической оптики, годографы прямой Р-волны вдоль профилей AC, MN и PQ: проходящих через источник импульса, представляют собой прямые линии: ti(г) = г/ср, где г - расстоянии от точки на профиле до источника волн. Годографы Р-волны вдоль профилей DB и EF, не проходящих через источник импульса, описываются уравнением параболы, симметричной относительно центра сейсмограммы: t(r) = \Д2 + г2/ср, где г - расстояние от точки профиля до центра, I - расстояние от центра профиля до источника импульса (I — OL для линии EF, I = О К для линии DB). Картина распространения падающей поперечной волны (S-волпы) полностью эквивалентна картине продольной волны и отличается только скоростью. Годографы S-волны вдоль профилей АС, NM и PQ имеют вид: t(r) = r/cs, вдоль профилей DB и EF: ^г(^) — у/12 + r2/cs.

Годографы отраженной РР-волиы вдоль профилей АС, NM и PQ, характеризуются уравнением: t(r) = — л/г2 - 2Hr sin 2<р + 4Я2 cos2 ip, (4.1)

Ср где ip - угол наклона горизонта (для профиля АС угол ср равен 5°, MN -ср = 8,56°, MN - ip = —1,19°). Годографы (4.1) представляют собой параболу, сдвинутую относительно центра сейсмограммы в сторону подъема горизонта на расстояние г\ = iJ sin время прихода отраженной волны ti = 2Н cos (p/ср. Эта точка отмечена на сейсмограммах как Р\. Для профиля АС: г\ = 27 м, t\ — 0.623 с, для профиля MN: ri=50 м, t\ = 0.618 с, для PQ: п -7, 3 м, ti = 0.625 с.

Вдоль профилей DB и EF уравнения годографов отраженной волны усложняется тем, что угол горизонта зависит от г, расстояния до центра сейсмограммы: ч tan + (г/Л tancpo ip(r) = arctan- , w ; — (4.2) y/1 + r2/l2 K J

Годографы отраженной волны вдоль профилей DB и EF: t(r) = -\fr2 + l2- 2Hy/l2 + г2 sin 2<p(r) + 4Я2 cos2 <p(r).

Ср

На рис. 4.3 представлены годографы отраженной РР-волны для профилей EF и DB.

Рис. 4.3: Годограф отраженной Р Р-волны вдоль профилей EF (а) и DB (б)

Годографы отраженной SS-волны вдоль профилей AC, NM и PQ описываются уравнениями t(r) = — \Jr2 - 2Яг sin 2tp + 4Я2 cos2 (р, cs где (р - угол наклона горизонта. Отраженная волна смещена относительно центра на расстояние r2 = Я sin 2с/?, время прихода отраженной волны t2 = 2Нcos(p/cs. Для профиля АС: г2=27 м, t2 = 1, 04 с; для MN: Г2=50 м, ^2=1,'03 с; для PQ: г2=-7.3 м, t2 = 1,04 с. Годографы отраженной SS-волны вдоль профилей DB и EF: t(r) = -^r2 + l2- 2Ну/12 + г2 sin 2ip(r) + 4Я2 cos2 сp{r), cs где угол <р(г) определяется по формуле (4.2).

При переходе продольной Р-волны через границу раздела двух сред кроме отраженной продольной волны образуется вторая отраженная волна - поперечная PS-волна. Аналогично при прохождении поперечной S-волны кроме отраженной поперечной SS-волны образуется отраженная продольная SP-волна. Эти волны обозначены на сейсмограммах цифрами 4 и 9 соответственно.

Угол отражения поперечной PS-волны (/3) связан с углом падения Р-волны (а) по закону Снеллиуса: sin в = — sin а. скорость поперечных волн в среде (cs) меньше скорости продольных (ср), поэтому (3 < а. Следовательно, при любом угле падения а (0 < а < отраженная волна выйдет на верхнюю границу среды на расстоянии: ч Я (tan а + tan/З) ,, . г(а) = -—, (4.3)

1 + tan <р tan а где ip - угол наклона горизонта. Время прихода отраженной Р5-волны определяется: t(a) -- Н

Cg tan а + Ср tan /3 cpc2s sin а(1 + tan <р tan а)' где ip - угол наклона горизонта, для профилей DB и EF определяется соотношением (4.2). Годограф отраженной Р5-волны представляет собой параболу не симметричную относительно центра сейсмограммы.

Угол отражения продольной 5Р-волпы (5) связан с углом падения S-волны (7) соотношением Снеллиуса: г ср • sm д = — sin 7.

Cs

В этом случае существует критический угол 7* = arcsin — = 36, 87° при котоср ром угол отражения становится максимальным (90°). При 7 < 7* расстояние выхода отраженной 5Р-волны определяется аналогично (4.3). Годограф отраженный 5Р-волны 7 < 7* описывается уравнением: cl tan 7 + tan 5 t(j) = Я£--, (4.4)

CpCs sin 7(1 ~Ь tan (р tan 7)' где (р - угол наклона горизонта, для профилей DB и EF определяется соотношением (4.2).

Критический угол г, при падении под которым луч продольной Р-волны преломляется под углом в 90° к нормали, определяется по закону Снеллиуса: sin г = Ср/Ср. При этом преломленная продольная волна начинает распространяться вдоль границы раздела сред. Фронт проходящей продольной

РР-волны, скользя вдоль границы раздела, возбуждает в верхнем слое колебания, которые вызывают появление головной преломленной волны. Одной стороной фронт головной волны касается фронта отраженной из критической точки волны, другой примыкает к фронту скользящей преломленной волны. На рис. 6 приложения Б представлены фрагменты визуализации вертикальной компоненты вектора перемещения для случаев падения продольной Р-волны под углом, меньшим критического угла i (а) и большим критического (б). При этом на рис. 6 (б) прослеживается образование головной продольной волны. Ее годографы для профилей АС, MN и PQ имеет вид: t(r) = — (г sm(i <р) + 2Н cos <р cos г), (4.5)

Ср где знак "—" берется для поднимающегося горизонта, а знак "+" - для опускающегося. Точка выхода головной волны на поверхность (точка Р2 на сейсмограмме) расположена от источника на расстоянии: гз = 2Н cos <рCos(i-p); время выхода головной волны на поверхность составляет: ts = с^оф-р) • Для профиля АС: гз=261 м, £3=0.776 с; для профиля MN: г3=246 м, £3=0.726 с; для профиля PQ: гз=287 м, ^з=0.866 с. Для профилей DB и EF годографы продольной головной волны аналогичны (4.5) с заменой г на л/12 + г2 и ip -на <р(г) по формуле (4.2): t(r) = — (\//2 + г2 sin(z =f ip) + 2Н cos ip cos г).

Ср

Из этого уравнения можно получить, что головная волна вдоль линии DB (.EF) выходит на поверхность в двух точках: на расстоянии 84 м (225 м) слева от центра сейсмограммы и на расстоянии 45 м (183 м) справа. Расстояние выхода головной продольной волны: 2Н cos (р(ф) sin г

-Г--ГмГ"' cos (г — фуф)} где ф - угол подъема горизонта: р(ф) = arctan tan (f \ + tan ф ■ tan y/l + tan2 ф

Картина распространения головной поперечной SS-волны аналогична. Уравнение годографа поперечной головной волны имеет вид: t = — (г sin(i ip) + 2Н cos <р cos г). cs

Точка выхода (точка на сейсмограмме) на поверхность расположена от источника на расстоянии Г4 = 2iiTcos(^i sin г/соз(г—ipi), время выхода составляет = 2Н cos2 ipi/(cs cos(i — ip)). Для профиля AC: 7-4 = 261 м, = 1.475 с. Головная поперечная волна выходит на поверхность в тех же точках, что и головная продольная волна, отличается только время.

На рис. 4.4 кривой линией представлена зависимость расстояния от точки приложения импульса (т. О) до точки выхода головной продольной волны на поверхность в зависимости от угла подъема горизонта (ф).

X,

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

4.2.2 Модельная задача 3. Прохождение волн через горизонтальную границу раздела двух сред

Рассмотрим модель двухслойной среды (аналогичную модельной задачи 2) с горизонтальной границей раздела сред (рис. 4.5.а). Параметры верхнего слоя: р=500 кг/м3, ср—500 и cs=300 м/с, параметры нижнего слоя: ро—750 кг/м3, Ср=750 и с®—450 м/с. Поверхность раздела двух сред представляет собой горизонтальную плоскость па некоторой глубине Н от верхней границы. Выполняются условия неравенства акустической жесткости двух сред: рср ф р0Ср и pcs ф р0с°3.

Рис. 4.5: Модель области с горизонтальной поверхностью раздела двух сред: а - геометрия среды, б - расположение вертикальных профилей в горизонтальной плоскости Х2Х3

Размеры расчетной области — АВ = AD = 400 м, АА2 = 200 м. Глубина залегания горизонта: АА\ = Н — 100 м. Внешний импульс действует в точке О. Величина приложенной нагрузки — Pq=250 МПА, время действия — 0.01 с, радиус действия — 2,5 м. Число узлов сетки — 200 х 200 х 200, временной шаг — т=0,002 с. Среднее время расчета одного шага на 8-ми процессорах — 1.26 минуты. Количество расчетных ячеек, приходящихся на один процессор

На рис. 8-9 приложения Б представлены сейсмограммы вертикальной сох. х

106 ставляющей вектора перемещения вдоль профилей ОС, ON, EF, DB и QN. На рис. 10-11 приложения Б представлены визуализации вертикальной компонент вектора перемещения вдоль профилей AC, MN и DB в момент времени 0.4 с. Задача является симметричной относительно оси АС.

Здесь введены следующие числовые обозначения: 1 - падающая Р-волна; 2 - падающая S-волна; 3 - поверхностная волна Релея; 4 - отраженная РР-волна; 5 - отраженная Рб'-волна; 6 - преломленная РР-волна; 7 - преломленная Рб'-волна; 8 - головная волна, образованная РР-волной; 9 - вторичные отраженные волны от верхней границы области; 10 - отраженная б'б'-волна; 11 - отраженная б'Р-волна; 12 - преломленная 5Р-волна; 13 - преломленная 55-волна; 14 - головная волна, образованная SP-волной; 15 - головная волна, образованная 55-волной.

Годографы прямых Р- и S-волн описываются аналогичными уравнениями, представленными в описании задачи 2. Годографы отраженных и головных волн упрощаются за счет того, что поверхность раздела сред является горизонтальной. Годографы отраженной РР-волны вдоль профилей АС и MN описываются уравнением параболы, симметричной относительно положения источника импульса (точка О): t — \/г2 + 4Н2/ср. Вдоль профилей EF, DB и QN годографы отраженной РР-волны характеризуются уравнением: t = у/г2 +12 + 4Я2/ср, где I - расстояние от центра профиля до источника импульса (/ = OL для профиля EF, I — О К для профиля DB, I = OG для профиля QN). Уравнение годографа отраженной 55-волны вдоль профилей ОС, ON: t = у/г2 + 4H2/cs, и вдоль линий EF, DB и QN: t = y/x2 + l2 + 4H2/cs.

Годограф продольной головной волны для профилей AC, MN и PQ описывается уравнением: t = (г sin г + 2Н cos i)/cp, где i - критический угол угол под которым продольная волна преломляется под углом 90° к нормали). Согласно формуле Снеллиуса величина критического угла составит г — arcsin ^ = 42°. Головная волна выйдет на поверхность (точка на сейсмограммах помечена Р) на расстоянии от источника волны, равным хя = 2Н sin г/ cosi= 180 м, вместе с отраженной волной. Однако далее головная волна начинает опережать отраженную волну.

Модельная задача 4. Прохождение волн через вертикальную границу раздела двух сред

Рассмотрим модель среды, состоящей из двух однородных изотропных блоков с вертикальной границей между ними и скачком акустической жесткости на вертикальной границе (рис. 4.6). Параметры одного блока: ^=500 кг/м3, ср=500 и ся=300 м/с, параметры второго блока: ро=750 кг/м3, c^^SO и Сд=450 м/с.

Рис. 4.6: Модель области с вертикальной поверхностью раздела двух сред

Размеры расчетной области - АВ = AD — 300 м, АА\ = 200, АР = BQ = I = 200 м. Внешний импульс действует в точке О. Величина приложенной нагрузки - Pq-250 МПА, время действия нагрузки - 0.01 с, радиус действия нагрузки - 1,5 м. Число узлов сетки - 100 х 300 х 300, временной шаг - т=0.01 с. Среднее время расчета одного шага по времени на 9-ти процессорах -1.13 минуты. Количество расчетных узлов на один процессор - 106.

На рис. 13 приложения Б представлены сейсмограммы вертикальной составляющей вектора перемещения вдоль профилей ОС и О2С2. На рис.12 приложения Б представлена визуализация компонентов вектора перемещения вдоль профиля АС в момент времени 0.195 с. Здесь введены следующие числовые обозначения: 1 - падающая Р-волна; 2 - падающая S-волна; 3 -поверхностная волна Релея; 4 - отраженная продольная РР-волна; 5 - отраженная поперечная PS-волна; 6 - преломленная продольная РР-волна; 7 - преломленная поперечная Р5-волна; 8 - отраженная поперечная волна 55-волна; 9 - отраженная продольная волна 5 Р-волна; 10 - преломленная поперечная б'З'-волна; 11 - преломленная продольная волна 5 Р-волна.

Вдоль профиля АС годографы всех волн являются прямыми линиями (рис. 13.а прил. Б): t = г/ср (Р-волна), t = r/cs (S-волна), t = 1/ср + г/сРр (отраженная РР-волна), t = l/cp + r/c® (отраженная Р5-волна), t = l/cs г/с® (отраженная S Р-волна), t = l/cs + r/c® (отраженная б'б'-волна), где I - расстояние от центра импульса до границы раздела сред (I = OR).

Отраженные РР— и SP—волны распространяются с одинаковой скоростью (с®). Аналогично PS— и SS—волны соответственно со скоростью с®. Поэтому линии годографов для этих пар воли представляют собой параллельные линии. На сейсмограмме хорошо просматриваются преломленные РР- и З^-волпы, в то время как преломленные PS- и 5Р-волны выражены крайне слабо. Большая часть сейсмической энергии при пересечении границы переходи в.волну той же поляризации, какая у падающей волны. Отраженных волн на сейсмограмме практически не видно.

Вдоль профиля О2С2, расположенного на глубине 75м, отраженные волны прослеживаются более явно (6,7,10,11). Кроме того, более четко прослеживается преломленная поперечная Р5-волна.

4.2.3 Модельная задача 5. Прохождение волн через границу раздела двух сред со сбросом

Рассмотрим модель трехмерного сброса (рис. 4.7, а) с амплитудой L = 15 м. Параметры материала верхнего слоя: р=500 кг /м3, ср=500 и cs=200 м/с, подстилающего слоя: pq—750 кг/м3, с^=750 и с^=300 м/с. Размеры расчетной области - АВ = AD = 750, АА2 = 600 м. Глубина залегания поверхности раздела материалов: АА\ — DD\ — 200, ВВ\ = СС\ = 185 м. Сброс проходит по линии ran. Внешний импульс на среду действует в точке О. Величина приложенной нагрузки - Ро—250 МПа, время действия - 0,01 с, радиус действия - 3.5 м. профилей в горизонтальной плоскости Х2Х3

Число узлов сетки - 300 х 300 х 300, временной шаг - т = 0.002 с. Число используемых процессоров - 27, количество расчетных ячеек, приходящихся на один процессор - 106. Число шагов по времени - 1000, среднее время расчета одного шага - 2.1 минуты.

Полученное в ходе расчетов волновое поле представлено на серии вертикальных плоскостей, секущих сброс под различными углами: 0° (EF), 11° (GH), 26.5° (KL), 37° (PQ), 45° (АС) и проходящими через источник импульса.

На рис. 14, 15 приложения Б представлены сейсмограммы вертикальной составляющей вектора перемещения вдоль профилей OF (а), ОН (б), OL (в), OQ (г) и ОС (д). На рис. 16 приложения Б представлены фрагменты визуализации вертикальной компоненты вектора перемещения для тех же профилей. Здесь введены следующие обозначения: 1 - падающая Р-волна; 2 - падающая S-волна; 3 - отраженная РР-волна; 4 - отраженная PS-волна; 5 - преломленная РР-волна; 6 - преломленная PS-волна; 7, 8 - дифрагированные волны.

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

Заключение

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

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

2. Алгоритм реализован в виде программного комплекса на языке программирования Fortran-90 с использованием функций библиотеки передачи сообщений MPI. Технология распараллеливания вычислений основана на разделении исходной области по процессорам, исходя из требований равномерной загрузки. Программный комплекс позволяет получать результаты, которые могут быть использованы в прикладных геофизических программах, таких как Promax или SeisView для анализа синтетических сейсмограмм.

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

Список литературы диссертационного исследования кандидат физико-математических наук Кучунова, Елена Владимировна, 2008 год

1. Якобский, М.В. Распределенные системы и сети / М.В. Якобский. - М. : МГТУ "Станкин", 2000. - 118с.

2. Самарский, А. А. Теория разностных схем: учебное пособие / А.А. Самарский. М. : Наука : Главная редакция физико-математической литературы, 1983. - 616с.

3. Domain decomposition method for partial differential equations. Proc. of the 1st Int. Symp. (Eds.R.Glowinsti et al.), SIAM, Philadelphia, 1988.

4. Алексеев, А.С. О лучевом методе вычисления интенсивности волновых факторов / А.С. Алексеев, Б.Я. Гельчинский // Вопросы динамической теории распространения сейсмических волн. JL, 1961. - N 5. - С.3-24.

5. Бабич, В.М. Асимптотические методы в задачах дифракции коротких волн. / В.М. Бабич, B.C. Булдырев. М. : Наука, 1982. - 248с.

6. Бабич, В.М. Пространственно-временной лучевой метод /В.М. Бабич, B.C. Булдырев, И.А. Молотков. JI. : Изд. Леи. Универ., 1985.

7. Клем-Mycamoe, К.Д. Теория краевых волн и ее применение в сейсмике / К.Д. Клем-Мусатов. Новосибирск : Наука, 1980.

8. Keller, J.B. Geometrical theory of diffraction /J.B. Keller //J. Opt. Soc. Am. 1981. - Vol.52. - N 4. - P. 175-188.

9. Киселев, А.П. Высшие приближения лучевого метода и "нелучевые явления" в неоднородных вязкоупругих средах /А.П. Киселев // Изв. АН СССР. Физика Земли. 1992. - N 11. - С.35-38.

10. Popov, М.М. A new method of computation of wave fields in the high-frequency approximation /М.М. Popov. Leningrad, 1981. - 20p. -(Preprint AN SSSR)

11. Алексеев, А.С. О лучевом методе вычисления полей волн в случае неоднородных сред с криволинейными границами раздела / А.С. Алексеев, Б.Я. Гельчинский // Вопросы динамической теории распространения сейсмических волн. JL, 1959. - N 3. - С. 16-47.

12. Cerveny, V. Ray method in seismology / V. Cerveny, I.A. Molotkov, I. Psencik. Prague : Varlovar. Univ.,1977.

13. Keller, J.B. Surface waves on waves on water of non-uniform depth /

14. J.B. Keller //J. Appl. Phys. 1980. - Vol.31. - N 6. - PP.1039-1046.

15. Tygel, M. Transient waves in layered media / M. Tygel, P. Hubral //Elsevier, Amsterdam. 1987.

16. Hron, F. Numerical modeling of nongeometrical effects by the Alekseev-Mikhailenko method / F. Hron, B.G. Mikhailenko //Bulletion of the seismological of America. 1981. - Vol.71. - N 4. - P.1011-1029.

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

18. Numerical methods used un atmospherical models // GARP Publication Series. 1979. - Vol.11. - N17.

19. Fornberg, В. On a Fouriese mehod for the integration of hyperbolic equation / B. Fornberg // Soc. Industr. Appl. Math., J.Numer. Anal. 1975. - N 12. - P.509-528.

20. Kreiss, H. Comparison of accurate methods for the integration of hyperbolic equation / H. Kreiss, J. Oliger // Tellus. 1972. - N 24. - P.199-215.

21. Orszag, S.A. Comparison of pseudospectral and spectral approximation /S.A. Orszag // Stud. Appl. Math. N 51. - P.253-259.

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

23. Kosloff, D. Elastic wave calculation by the Fourier method /D.Kosloff, M. Reshef, D. Loewenthal // Bull. Seis. Am. 1984. - N 74. - P.875-891.

24. Михайленко, Б. Г. Численное моделирование сейсмических полей в двумерно-неоднородных средах / Б.Г. Михайленко, В.И. Корнев // Моделирование волновых полей. Новосибирск, 1983 - С.41-70.

25. Mikhailenko, B.G. Synthetic seismograms for complex threedimensional geometries using an analytical-numerical algorithm /B.G. Mikhailenko // Geopys. J. R. Astr. Soc. 1984. - N 79,3. - P. 963-986.

26. Mikhailenko, B.G. Numerical experiments in seismic investigations /B.G. Mikhailenko // J. Geophys. 1985. - N 58. - P.101-124.

27. Алексеев, А. С. О задаче Лэмба для неоднородного полупространства / А.С. Алексеев, Б.Г. Михайленко// Докл. АН СССР. 1974. - Т.214. -N 11. - С.84-86.

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

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

30. Алексеев, А. С Метод вычисления теоретических сейсмограмм для сложно построенных моделей сред/А.С. Алексеев, Б.Г. Михайленко // ДАН СССР. 1978. - т.240. - N 5. - С.1062-1065.

31. Михайленко, Б.Г. Сейсмические поля в сложнопостроенных средах (Атлас численных снимков и теоретических сейсмограмм) / Б.Г. Михайленко; под ред. А.С. Алексеева. Новосибирск, 1988. - 312с.

32. Петпрашенъ, Г.И. Распространение упругих волн в слоисто-изотропных средах, разделенными параллельными плоскостями / Г.И. Петрашень //Уч. Зал. ЛГУ, 1952. N 162. - С.3-189.

33. Петрашень, Г.И. Общая количественная теория отраженных и головных волн, возбуждающихся в слоистых средах с плоско-параллельными границами раздела /Г.И. Петрашень //Вопросы динамической теории распространения сейсмических волн. JL, 1957. - С. 70-164.

34. Огурцов, К. И. Динамические задачи для упругого полупространства в случае осевой симметрии /К.И. Огурцов, Г.И. Петрашень. Учен. зап. ЛГУ, 1951. - N 149. - С.3-117.

35. Петрашень, Г.И. О распространении волн в сложно-изотропных средах / Г.И. Петрашень, И.Н. Успенский. Учен. зап. ЛГУ, 1956. - N 208.1. C.58-141.

36. Петрашень, Г.И. Элементы динамической теории распространения сейсмических волн / Г.И. Петрашень //Вопросы динамической теории распространения сейсмических волн. Л., 1959. - Сб.З. - С.11-106.

37. Thomson, W. Т. Transmission of elastic waves through a stratified solid medium /W.T. Thomson // J. Appl. Phys. 1959. - N 21. - P.89-93.

38. Haskell, N.A. The dispersion of surface waves a multilayered media / N.A. Haskell // Bull. Seismol. Soc. Amer. 1953. - Vol.43. - N 1. - P. 17-34.

39. Harkrider D.G. Surface waves in multilayered elastic media. 1. Rayleigh and Lowe waves from buried sources in a multilayered elastic half-space. /

40. D.G. Harkrider // Bul.Seimol.Soc.Amer. 1964. Vol.54. - N 2. - P.627-673.

41. Молотков, Л. А. Матричный метод в теории распространения волн в слоистых упругих и жидких средах / J1.A. Молотков. JI. : Наука, 1984. - 366с.

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

43. Alterman, Z. Propagation of elastic waves layered media by finite-difference method / Z. Alterman, F.C. Karal // Bull. Seism.Soc. Am. 1968. - Vol. 58. - P.367-398.

44. Boore, D.M. Finite-difference methods for seismic wave propagation in heterogeneous materials / D.M. Boore // Methods in Computational Physics. 1972. - N 11. - P.l-37.

45. Smith, W.D. The application of finite element analysis to body wave propagation problems / W.D. Smith // Geophys. J.R.Astr.Soc. 1975. - N 42. - P.747-768.

46. Zahradnik, J. Finite-difference solution to certain diffraction problem / J. Zahradnik // Stud, and Gead. 1975. - Vol.19. - P.233-244.

47. Кошур, В.Д. Континуальные и дискретные модели динамического деформирования элементов конструкций /В.Д. Кошур, Ю.В. Немиров-ский. Новосибирск: Наука. Сиб. отд.-ние, 1990. - 198с.

48. Магомедов, К.М. Сеточно-характеристические численные методы /К.М. Магомедов, А.С. Холодов. М.: Наука, 1988.

49. Магомедов, К.М. О построении разностных схем для уравнений гиперболического типа на основе характеристических соотношений /К.М. Магомедов, А.С. Холодов // Журнал вычислительной математики и математической физики. 1969. - Т.9. - N 2. - С. 373-386.

50. Петров, И. Б. Численное решение некоторых динамических задач деформируемого твердого тела сеточно-характеристическим методом / И.Б. Петров, А.С. Холодов // Журнал вычислительной математики и математической физики. 1984. - N.5. - С. 722-739.

51. Рахматулин, Х.А. Применение метода пространственных характеристик к решению задач по распространению упругопластических волн / Х.А. Рахматулин, Т.Д. Каримбаев, Т. Байтелиев // Изв. Кав. СССР. Сер. физ.-мат. 1973. - N. 1. - С. 141-152.

52. Сабодаш, П. О. Применение методы пространственных характеристик к решению осесимметричных задач по распространению упругих волн / П.О. Сабодаш, Р.А. Чередниченко // ПМТФ. 1971. N. 4. - С. 101-109.

53. Годунов, С.К. Численное решение многомерных задач газовой динамики /С.К. Годунов, А.В. Забродин, А.В. Иванов и др. М.: Наука, 1976. -400с.

54. Lax, P.D. Systems of conservation laws / P.D. Lax, Wendroff // Communs Pure and Apppl. Math. 1960. - V.13. - P.217-237.

55. Lax, P.D. Weak solutions of nonlinear hyperbolic equations and their numerical computation /P.D. Lax // Communs Pure and Apppl. Math.- 1954. V.7. - P.159-193.

56. Русанов, В. В. Разностная схема третьего порядка точности для сквозного счета разрывных решений / В.В. Русанов // ДАН СССР. 1968.- N 6. С.1303-1305.

57. Courant, R. On the solution of nonlinear hyperbolic differential equations by finite differences / R. Courant, E. Isaacon, M. Rees // Communs Pure and Apppl. Math. 1952. - vol.5.

58. Roe, P.L. Approximate Riemann solvers, parameters vectors and difference schemes / P.L. Roe // Journal of Computational Physics. 1981. - vol. 43.

59. Roe, P.L. Characteristic based schemes for Euler equations / P.L. Roe // And. Rev. Fluid Mechanics. 1986. - Vol. 18.

60. Roe, P.L. Efficient construction and utilization of approximate Riemann solutions / P.L. Roe, J. Pike // Computing Method in Applied Sciences and Engineering. Amsterdam, North - Holland, 1984. - Vol.6

61. Osher, S. Numerical solution of singular perturbation problems and hyperbolic systems of conservation laws / S. Osher // North Holland Mathematical Studies. 1981. - Vol. 47.

62. Steger, J.L. Flux vector splitting of the inviscrid gas dynamic equations with application to finite difference methods / J.L. Steger, R.F. Warming // Journal of Computational Physics. 1981. - Vol. 40. - No 2.

63. Harten, A. High resolution schemes for hyperbolic conservation laws / A. Harten // Journal of Computational Physics. 1983. - Vol.49. - P.367-393.

64. Swety, P.K. High resolution schemes using flux limiters for hyperbolic conservation laws / P.K. Swety // SIAM J. Numer. Analys. 1984. - Vol.21. - P.995-1011.

65. Yee, H.C. Application of TVD-schemes for the Euler equations of gas dynamics / H.C. Yee, R.F. Warming, A. Harten // Lecture of Applied Mathematics. 1985. - Vol. 22.

66. Harten, A. On a class of high resolution total-variation-stable finite-difference schemes / A. Harten // NYU Report New York, NYU, 1982.

67. Harten, A. The artificial compression method for computation of shocks and contact discontinuities: III. Self-adjusting hybrid schemes / A. Harten // Math. Comput. 1978. - Vol. 32. - P.363-389.

68. Chakravarthy, S.R. Computing with high-resolution upwind schemes for hyperbolic equation / S.R. Chakravarthy, S. Osher // Lectures in Applied Mathematics. 1985. - vol. 22.

69. Годунов, С.К. Разностный метод численного расчета разрывных решений гидродинамики / С.К. Годунов // Математический сборник. 1959. - Т.47. - N 3. - С.271-306.

70. Van Leer В. Toward the ultimate conservation difference scheme. V. A second-order segcul to Godunov's method / Van Leer B. // Journal of Computational Physics. 1979. - V.32. - P. 101-136.

71. Harten, A. Uniformly high-order accurate nonoscillatory schemes / A. Harten, S. Osher // SIAM J. Numer. Analys. 1987. - V.24. - P.279-309.

72. Richtmycr, R.D. A survey of difference method for nonsteady flyid dynamics / R.D. Richtmycr // NCAR Technical Note 63-2-Colorado, Boulder, 1963.

73. Mac Cormack R. W. The effect of viscosity in hypervelocity impact cratering / Mac Cormack R.W. // AIAA Paper N 69 354, 1969.

74. Колган, В. П. Применение принципа максимальных значений производной к построению конечно-разностпых схем для численного анализа разрывных течений гидродинамики / В.П. Колган // Ученые записки ЦАРИ. 1972. - Т 3.

75. Van Leer В. Toward the ultimate conservation difference scheme. I. The quest of monotonicity / Van Leer В. // Lecture Notes in Physics. 1973. -vol.18.

76. Van Leer B. Towards the ultimate conservative difference scheme II. Monotonicity and conservation combined in a second-order scheme / Van Leer B. // Journal of Computational Physics. 1974. - Vol.14.

77. Van Leer B. Towards the ultimate conservative difference scheme III. Upstream-centered finite-difference schemes for ideal compressible flow /

78. Van Leer В. // Journal of Computational Physics. 1977. - Vol. 23. - P. 263-275.

79. Van Leer B. Towards the ultimate conservative difference scheme IV. A new approach to numerical convection / Van Leer B. // Journal of Computational Physics. 1977. - Vol. 23. - P. 276-299.

80. Harten, A. The method of artificial compression / A. Harten // CIMS Report COO 30077-50- New York, Courant Institute, NYU, 1974.

81. Harten, A. Self-adjusting hybrid schemes for shock computation / A. Harten, G. Zwas // Journal of Computation Physics. 1972. - Vol.6.

82. Beam, R. An implicit finite-difference algorithm for hyperbolic systems in conservation-law-form / R. Beam, R.F. Warming // Journal of Computational Physics. 1976. - Vol.22.

83. Boris, J.P. Flux-corrected transport. I. SHASTA, a fluid transport algorithm that works / J.P. Boris, D.L. Book // Journal of Computational Physics. 1973. - Vol. 11.

84. Boris, J.P. Flux-corrected transport. II Generalization of the method / J.P. Boris, D.L. Book, K. Hain // Journal of Computational Physics. -1975. Vol. 18.

85. Boris, J.P. Flux-corrected transport. II Minimal-error FGT algorithms / J.P. Boris, D.L. Book // Journal of Computational Physics. 1976. - Vol. 20.

86. Коновалов, A.H. Численные методы в динамических задачах теории упругости /А.Н. Коновалов. Сиб. матем. журн. - Новосибирск, 1997. - 38:3. - С. 551 - 568.

87. Белоцерковский, О.М. Метод расщепления для исследования течений стратифицированной жидкости со свободной поверхностью / О.М. Белоцерковский, В.А. Гущин, В.Н. Копынин // Журнал вычислительной математики и математической физики. 1987. - Т. 27.

88. Коновалов, А.Н. Вариационная оптимизация итерационных методов расщепления / А. Н. Коновалов. Сиб. матем. эюурн. - Новосибирск, 1997. - 38:2. - С. 312-325.

89. Софронов, И.Л. Искусственные граничные условия адекватные волновому уравнению вне сферы / И.Л. Софронов. Препринт Института Прикладной Математики им. М.В.Келдыша Академии Наук СССР. -Москва, 1992. - No. 42.

90. Софронов, И.Л. Условия полной прозрачности на сфере для трехмерного волнового уравнения / И.Л. Софронов // Доклады Академии Наук.- 1992. Т. 326. -С. 953-957.

91. Keller, J.B. Exact nonreflecting boundary conditions / Joseph B. Keller, Dan Givoli // J.Comput. Phys. 1989. - Vol. 82. - P.172-192.

92. Givoli, D. nonreflecting boundary conditions based on Kirchhoff-type formulae / Dan Givoli, D.Cohen. // J. Comput. Phys. 1995. - Vol. 117.- P.102-113.

93. Grote,M.J. Exact nonreflection boundary conditions for the time-dependent wave equation / M.J. Grote, J.B. Keller. // SIAM J. Appl. Math. 1985.- Vol 55. P. 280-297.

94. Grote M.J. Nonreflecting boundary conditions for time-dependent scattering /M.J. Grote, J.B. Keller. // J. Comput. Phys. 1996. - Vol 127. - P.52-65.

95. Engquist, B. Far field boundary conditions for computation over long time / B. Engquist, L. Halpern. // Appl. Numer. Math. 1988. - Vol. 4• -P.21-45.

96. Engquist, B. Absorbing boundary conditions for the numerical simulation of waves / B. Engquist, A. Majda // Math. Comput. 1977. - Vol.31. -P. 629-651.

97. Higdon, R.L. Absorbing boundary conditions for difference approximations to the multidimensional wave equation / Robert L. Higdon // Math. Сотр. 1986. - Vol. 47(176). - P.437-459.

98. Higdon, R.L. Absorbing boundary conditions for accoustic and elastic waves in stratified media / Robert L. Higdon // J. Comput. Phys. 1992. - Vol. 101(2). - P.386-418.

99. Bayliss, A. Radiation boundary conditions for wave-like equations / Alvin Bayliss, Eli Turkel // Comm. Pure Appl. Math. 1980. - Vol. 33(6). -P.707-725.

100. Givoli, D. Nonreflection boundary conditions for elastic waves / Dan Givoli, Joseph B. Keller // Wave Motin. 1990. - Vol. 12(3), P.261-279.

101. Hedstrom, G.W. Nonreflecting voundary conditions for nonlinear hyperbolic systems / G.W. Hedstrom // J. Comput. Phys. 1979. - Vol. 30(2). -P.222-237.

102. Atkins, H. Nonreflective boundary conditions for high-order methods / H. Atkins, J. Casper // AIAA Journal. 1994. - Vol. 32. - P.512-518.

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

104. Крюков, В.А. Разработка параллельных программ для вычислительных кластеров и сетей / В.А. Крюков //Информационные технологии и вычислительные системы. 2003. - N 12.

105. Марчук, Г.И. Методы растепления / Г.И. Марчук. М. : Наука, 1988.

106. Садовская, О. В. Метод сквозного счета для исследования упругопла-стических волн в сыпучей среде / О.В. Садовская // Ж. вычисл. математики и матем. физики. 2004• - N 10. - С. 1909-1920.

107. Куликовский, А.Г. Математические вопросы численного решения гиперболических систем уравнений / А.Г. Куликовский, Н.В. Погорелое, А.Ю. Семенов. М. : Физматлит, 2001.

108. Магометов, К.М. Сеточно-характеристические численные методы / К.М. Магометов, А.С. Холодов. М. : Наука, 1988.

109. Воеводин, В.В. Параллельные вычисления / В.В. Воеводин, Вл.В. Воеводин Спб. : ВХВ-Перетбург, 2002.

110. Ортега, Дою. Введение в параллельные и векторные методы решения линейных систем : пер. с англ.] / Дж. Ортега. М. : Мир, 1991.

111. Кучунова, Е.В. Параллельный алгоритм расчета упругих волн в блочной среде / Е.В. Кучунова // Материалы VI межрегиональной школы-семинара "Распределенные и кластерные вычисления". Красноярск, 2006г. - С. 79-89.

112. Кучунова, Е.В. Численное исследование распространения сейсмических волн в блочных средах на многопроцессорных вычислительных системах / Е.В. Кучунова, В.М. Садовский // Вычислительные методы и программирование. 2008. - Т.9. - N 1. - С.70-80.