УДК 528
Д.М. Вишневский, В.В. Лисица ИНГГ СО РАН, Новосибирск
КОНЕЧНО-РАЗНОСТНОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ РАСПРОСТРАНЕНИЯ ВОЛНОВЫХ ПОЛЕЙ В АНИЗОТРОПНЫХ УПРУГИХ СРЕДАХ
D.M. Vishnevsky, V.V. Lisitsa
Trofimuk Institute of Petroleum Geology and Geophysics SB RAS Acad. Koptyug av., 3, Novosibirsk, 630090, Russian Federation
FINITE DIFFRENCE SIMULATION OF WAVES’ PROPAGATION IN ANISOTROPIC ELASTIC MEDIA
Study of waves’ propagation within complex anisotropic heterogeneous elastic media can be done only via accurate numerical simulation as analytical solution can not be constructed for the considered problems. Numbers of methods exist nowadays to deal with isotropic models while proper treatment of anisotropic ones is under investigation. In the paper we present a second-order conservative finite difference scheme on staggered grids for simulations of waves within arbitrary anisotropic elastic media. We also provide series of the numerical experiments demonstrating efficiency of the proposed approach and typical peculiarities of wavefields caused by the presence of anisotropy.
Введение
Численное моделирование формирования и распространения волновых полей в сложноустроенных существенно неоднородных анизотропных средах является в настоящее время, пожалуй, единственным способом изучения данных процессов. По этой причине разработка эффективных высокоточных алгоритмов такого моделирования является весьма актуальной задачей. В настоящее время существует достаточно большое подходов к численному решению задачи распространения упругих волн в изотропных средах. Прежде всего, следует отметить конечно-разностные методы, в частности, так называемую, схему Вирье [1], которая на данный момент является наиболее распространенной и часто применяемой для решения данного класса задач в силу простоты реализации и, что более важно, консервативности. Последнее фактически означает, что данная схема сохраняет второй порядок сходимости на разрывах (границах слоев), теоретическое обоснование данного факта представлено в [2]. Следует также упомянуть, так называемые, схемы высокого порядка, позволяющие вычислять решение с более высокой точность, но исключительно лишь в случае гладкой скоростной модели (их применение вполне обоснованно, например, для конечно-разностной миграции в обратном времени), в то время как наличие разрывов (контрастных разделов сред) автоматически понижает точность данных методов до первого порядка. Существуют также и другие подходы к численному моделированию, например, спектральные, псевдоспектральные методы, метод Галеркин. Основным преимуществом
данных подходов является возможность применения нерегулярных сеток, что позволяет достаточно точно описывать геометрию модели. Но, в то же самое время, это является и их недостатком, поскольку требует предварительного построения сетки, что само по себе есть достаточно нетривиальная и трудоемкая задача.
Несмотря на огромное количество подходов к моделированию волновых процессов в изотропных средах, их обобщение на случай анизотропии требует существенных, далеко нетривиальных, модификаций. В частности, в настоящее время, имеются некие попытки применения методов типа Г алеркина и спектральных методов для решения данной задачи, например, но в силу своей сложности и недостатка универсальности, данные подходы не находят широкого применения. В отличие от упомянутых подходов, конечноразностные методы оказались более востребованными. Так, в частности, в работе [3] была предложена, так называемая схема на повернутых сетках или RSGS (от английского: Rotated Staggered Grid Scheme) которая является, пожалуй, наиболее применяемой на данный момент. Тем не менее, следует отметить, что данная схема является, во-первых, более сложной в реализации, чем, например, схема Вирье [1], во-вторых, требует хранения существенно большего количества переменных, что снижает ее эффективность. В настоящей работы представлена схема Лебедева, которая является естественным обобщением схемы Вирье [1] на случай анизотропной среды. Как следствие, данная схема обладает всеми свойствами схемы Вирье, а именно простотой применения, консервативностью и пр. В работе также показывается, что схема Лебедева позволяет существенно сократить требования на RAM, в сравнении с RSGS, что, в силу специфики архитектуры современных CPU (доступ к памяти занимает существенно больше времени, чем выполнение арифметических операций), приводит к увеличению эффективности расчетов.
Представленные в работе теоретические результаты подтверждаются серией численных экспериментов моделирующих распространение волновых полей в анизотропных средах.
Теоретическое обоснование
Прежде чем переходить к описанию схемы Ледебева и ее сравнению со схемой на повернутых сетках [3], рассмотрим схему Вирье [1], широко применяемую для моделирования волновых процессов в изотропных средах. В силу специфической формы закона Гука, то есть квази-диагональной структуры тензора упругих модулей, переменные могут быть разнесены как в пространстве, так и по времени. В терминах построения схемы, различные переменные определены разных узлах сетки. Пример сетки используемой для схемы Вирье представлен на рис. 1. Данный способ дискретизации и использование сдвинутых сеток позволяет аппроксимировать уравнения динамической теории упругости со вторым порядком, используя простой двухточечный шаблон и операторы следующего вида:
ги _ ги гп+1/2 _ ги-1/2
/(, = Г, X = х.,^ 2 ^) ~ ;'+1/2]к /-1/2* д ~ ;']к ;']к
дх кх д1 г
Операторы, аппроксимирующие производные по другим пространственным направлениям могут быть получены посредством простой замены индексов.
Рис. 1. Сетки для схемы Вирье (слева), RSGS (в центре) и семы Лебедева
(справа)
В случае общей анизотропии, тензор упругих модулей не обладает какой-либо специфической структурой, по этой причине переменные не могут быть разнесены в пространстве и все компоненты тензора напряжений должны быть определены в одном наборе узлов, в то время как все компоненты скорости - в другом. Попытки применения схем, схожих со схемой Вирье, для решения таких задач, очевидно, приводят к необходимости использования интерполяций, что существенным образом снижает точность вычислений. Иные результаты были получены в работе [3], в которой были предложены схемы на повернутых сетках. Пример такой сетки представлен на рис. 1. При этом пространственные производные аппроксимируются вдоль некоторых специальных направлений (главных диагоналей куба), после чего производные вдоль физических направлений строятся как линейная комбинация вычисленных вспомогательных производных. На первый взгляд, авторам удалось обойти проблему интерполяции, однако более детальный анализ, показывает, что это далеко не так, и операторы, аппроксимирующие пространственные производные представляются в следующем виде:
д/_ад+1,2к-ад г-1/2]к ^|- ^-|п 1 I ^п | /п | £ п | /п ^
^ [ / ]г+1/2 Iк ~~\/г+1/2 /+1 / 2к+1/2 + /г+1/2 /-1/2 к+1/2 + /г+1/2 /+1/2 к-1/2 + /г+1/2 /-1/2 к-1/2/
дх пх 4
В данной работе представлен принципиально иной подход к построению конечно-разностной схемы для анизотропной упругости, называемый схемой Лебедева. Наиболее простое его представление может быть описано как комбинация четырех сеток для схемы Вирье, сдвинутых относительно друг друга, так что все компоненты вектора скорости и тензора напряжений определены в соответствующих узлах. Пример сетки для схемы Лебедева представлен на рис. 1. Несложно понять, что использование данной геометрии сеток позволяет обойти любого рода интерполяцию, дает возможность
использовать простейшие конечно-разностные операторы для аппроксимации пространственных производных и, как следствие, сохранить точность схемы.
Как было показано в работе [4], несмотря на то, что сетка для схемы Лебедева на первый взгляд выглядит существенно плотнее, чем для RSGS, в двумерной постановке рассматриваемые схемы фактически эквивалентны. В трехмерном случае, как это следует из дисперсионного анализа и анализа устойчивости, схема Лебедева позволяет использовать несколько более грубую сетку, чем RSGS, для проведения моделирования с одинаковой точностью. Однако ключевым отличием представленных схем является способ задания модели. Как уже говорилось выше, в силу консервативности обеих схем все параметры среды (тензор упругих модулей и плотность) предполагаются быть заданными в целых точках по пространству и могут быть пересчитаны в прочие точки по мере необходимости, теоретическое обоснование данного подхода представлено в [2]. Несложно понять, что даже применение примерно одинаковых сеток для рассматриваемых схем, приводит к тому, что для схемы Лебедева, число точек, в которых хранятся параметры среды в четыре раза меньше, чем для RSGS. Если быть более точным, применение схемы Лебедева требует всего 30-60% памяти, используемой RSGS для проведения аналогичных расчетов. Если также учесть тот факт, что доступ к памяти может занимать гораздо больше времени, чем непосредственно сами вычислительные операции, то данное преимущество схемы Лебедева становится еще более существенным.
Численные эксперименты
Для иллюстрации эффективности предложенного метода была проведена серия численных экспериментов. Моделирование было проведено для трансверсально-изотропной среды с осью симметрии принадлежащей плоскости Х7, при том в одном случае рассматривалась УТ1 модель, в другом та же модель, но с осью, наклоненной на 45 градусов. Параметры Т1 модели были выбраны следующими:
р = 2500к^/т3, ^ = 4000т /5, ^ = 2055т / 5, е = 0.334, у = 0.575, 8* = 0.93.
Индикатрисы групповых скоростей для УТ1 модели представлены на рис. 2. Для модели с наклоненной осью симметрии индикатрисы, очевидно, могут быть получены посредством поворота на 45 градусов. На рис. 2 представлены рассчитанные волновые поля для УТ1 модели с использованием источника типа центра расширения (слева) и диполя с моментом (справа). Можно видеть, что, как это следует из теории, см. [5], применение источника типа центра расширения не порождает SH волны, в то время как данная волна хорошо видна при использовании диполя с моментом. Результаты аналогичных экспериментов для ТТ1 модели представлены на рис. 3 и 4. Можно видеть, что результаты моделирования (волновые фронты) с высокой точностью совпадают с индикатрисами групповых скоростей.
о
Рис. 2. Индикатрисы групповых скоростей для рассматриваемой УТІ модели (слева). Волновые поля (х-компоненты) для УТІ модели с источником типа центра расширения (в центре) и диполя с моментом (спрева)
о
Рис. 3. Волновые поля х, у и z компоненты соответственно для ТТ1 модели при использовании источника типа центра расширения
Рис. 4. Волновые поля х, у и z компоненты соответственно для ТТ1 модели при использовании источника типа диполя с моментом
Заключение
Представленная в работе консервативная конечно-разностная схема Лебедева для моделирования процессов формирования и распространения волновых полей в анизотропных упругих средах является, по мнению авторов, наиболее универсальным и эффективным подходом к численному моделированию. Прежде всего, в силу простоты в применении - данный подход не требует посторенние сложной сетки для каждой конкретной модели. В сравнение с прочими схемами, в частности RSGS, схема Лебедева является наиболее экономичной и позволяет существенно в 2-3 раза сократить затраты на реализацию.
Данная работа была выполнена при частичной финансовой поддержке РФФИ, гранты № 07-05-00538, 08-05-00265, 09-01-00253 и 09-05-00372.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Virieux J. P-SV wave propagation in heterogeneous media: Velocity-stress finite-difference method // Geophysics. - 1986. - Vol. 51. - P. 889-901.
2. Самарский А.А. Теория разностных схем / А.А. Самарский. - М.: Наука, 1983. -
616 с.
3. Saenger E. H., Gold N., Shapiro S. A. Modeling the propagation of the elastic waves using a modified finite-difference grid // Wave Motion. - 2000. - Vol. 31. - P. 77-92.
4. Lisitsa V.V. Efficient Finite-difference scheme for anisotropic elastic problems // Proceedings of the 8-th International Conference on theoretical and computational acoustics. -2008. - P. 285-297. - ISBN 978-960-89758-4-2
5. Гольдин С.В. Сейсмические волны в анизотропных средах / С.В. Гольдин -Новосибирск: издательство СО РАН. - 2008. -375 с.
© Д.М. Вишневский, В.В. Лисица, 2009