CC BY
40
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
УДК 519.6
А. В. Фаворская, Д. И. Петров, Н. И. Хохлов, И. Б. Петров
ЧИСЛЕННОЕ РЕШЕНИЕ УПРУГОАКУСТИЧЕСКИХ ЗАДАЧ С ПОМОЩЬЮ СЕТОЧНО-ХАРАКТЕРИСТИЧЕСКОГО МЕТОДА
Сеточно-характеристическим методом моделируются волновые процессы в упругоакустических задачах, в частности в задачах шельфо-вой сейсморазведки. Рассматриваются системы «вода — грунт — угле-водородосодержащий слой» и «лед — вода — грунт — углеводородосо-держащий слой». Разработано условие на контакте между жидкостью и твердым телом.
Wave processes in elastic-acoustic problems are simulated with the help of grid-characteristic method, in particular, for the offshore seismic prospecting. We consider systems of water-soil-hydrocarbon-bearing layer and ice-water-soil-hydrocarbon-bearing layer.
Ключевые слова: сеточно-характеристический метод, механика деформируемого твердого тела, упругоакустические задачи, шельфовая сейсморазведка, модель флюидонасыщенной трещины.
Key words: grid-characteristic method, the mechanics of deformable solids, wave problems in coupled elastic-acoustic media, marine seismology, fluid-filled crack, structural acoustics.
Введение
Проблема поиска перспективных месторождений углеводородов — одна из основных стратегических задач. В последнее время рассматривается возможность добычи в условиях северных морей, что требует совершенствования методов гидросейсмической разведки и их численного моделирования [1].
В данной работе приведены результаты численного эксперимента по моделированию волновых процессов в системах «вода — грунт» се-точно-характеристическим методом [2 — 6]. Была получена форма отклика от резервуара с месторождением углеводородов при различных начальных условиях на больших глубинах. Эксперимент производятся для различных конфигураций геологических отложений и двух типов начальных возмущений. Также исследовались волновые процессы в системе «лед — вода — грунт», причем грунт содержал неоднородное нефтесодержащее включение.
1. Постановка задачи
Решаются полная система уравнений, описывающих состояние сплошной линейно-упругой среды и полная система уравнений, описывающая акустическое поля. Компоненты скорости движения v и
7
© Фаворская А. В., Петров Д. И., Хохлов Н. И., Петров И. Б., 2015
Вестник Балтийского федерального университета им. И. Канта. 2015. Вып. 10. С. 7—12.
симметричного тензора напряжений Коши ст в линейно-упругой среде описываются следующей системой уравнений [7]:
р^б = (V-о)т, (1)
д т
— о = + ц(У®б + ^®б)т). (2)
й
Для численного моделирования морской воды и нефтесодержащих включений использовалось приближение идеальной жидкости [8] и решалась полная система уравнений, описывающая акустическое поле давления р и компонент скорости V :
д
р—V = ^р, (3)
дí
д 2
- р = -с 2р( V-V). (4)
дí
В выражениях (1), (2) X, ц — параметры Ляме, определяющие свойства упругого материала, а ® Ь — тензорное произведение векторов а и Ь , (а ® ЬУ = а'Ь', I — единичный тензор второго ранга. В уравнении
(4) за с обозначена скорость звука в идеальной жидкости. Скорость продольных волн в линейно-упругой среде можно найти по формуле
Ср , <5)
а скорость поперечных волн вычисляется в соответствии с
С, . (6)
2. Численный метод
Для численного решения систем (1), (2) и (3), (4) используется се-точно-характеристический метод на криволинейных структурных сетках, позволяющий строить корректные численные алгоритмы для расчета граничных точек и точек, лежащих на поверхностях раздела двух сред с разными параметрами Ляме и (или) плотностями.
Систему (1), (2) в двумерном случае можно представить в виде
В. + А1 м. + А2 = 0. (7)
дt дх1 дх2
В выражении (7) под вектором ц понимается вектор, составленный из двух компонент скорости и трех компонент симметричного тензора напряжений ц е {с1, б2, ст11, ст22, ст12 }т.
Вначале применяется метод расщепления по пространственным координатам, в результате чего имеем две одномерных системы
^ = А, . (8)
дt ' дх,
8
Каждая из этих систем является гиперболической и обладает полным набором собственных векторов с действительными собственными значениями, поэтому каждую из систем можно переписать в виде
=о- .о. ,
дt 1,1 дх]
где О, — матрица, составленная из собственных векторов, Л, — диагональная матрица, элементами которой служат собственные значения. Для всех координат матрица Л имеет вид (индекс , далее опускается, где это возможно)
Л = ,-ср ,св ,-св ,0}.
После замены переменных р = Оу каждая из систем (8) распадается на пять независимых скалярных уравнений переноса
др + Лдр = 0.
дt дх
Одномерные уравнения переноса решаются с помощью метода характеристик либо обычными конечно-разностными схемами.
После того как все компоненты перенесены, восстанавливается само решение:
у = О р .
Решение системы (3), (4) в двумерном случае, а также решение систем (1), (2) и (3), (4) в трехмерном случае проводится по аналогичным алгоритмам.
В программе реализовано применение ТУО-разностных схем 2-го порядка точности [9], 15 различных лимитеров [10], в расчетах в основном используются ограничитель 8ирегЪее [11] и сеточно-характеристи-ческие схемы 2 — 4 порядка точности [6].
3. Контактное условие между жидкостью и твердым телом
Пусть в части области интегрирования а решается система уравнений (1), (2), а в части области интегрирования Ь решается система уравнений (3), (4). Рассмотрим контактное условие между ними. Система (1), (2) в трехмерном случае имеет три выходящих характеристики, а система (3), (4) в трехмерном случае имеет одну выходящую характеристику. Таким образом, для нахождения всех четырех выходящих характеристик необходимо выполнение следующих контактных условий:
рЬ -"+1 =-(оа '"+1 • р )• р, (8)
оа '"+1 • р -(оа '"+1 • р )• р = 0 (9)
Vа '"+1 • р = VЬ '"+1 • р . (10)
Условие (8) представляет собой равенство нормальной компоненты поверхностной плотности сил со стороны твердого тела давлению в идеальной жидкости, условие (9) обеспечивает равенство нулю тангенциальной компоненты поверхностной плотности сил со стороны твердого тела, а выражение (10) задает равенство нормальных компонент
9
скоростей в идеальной жидкости и твердом теле. В (8) —(10) за р обозначена внешняя нормаль к твердому телу, являющаяся внутренней нормалью к жидкости.
Вычисляется вектор скорости V такой, чтобы условия (8) — (10) были выполнены. Далее для твердого тела используется граничное условие с заданной скоростью, а для жидкости — граничное условие с заданной нормальной компонентой скорости.
4. Результаты численного моделирования
10
В данной работе исследовались волновые картины, возникающие в системах «вода — грунт» и «лед — вода — грунт», причем грунт имел нефтесодержащее включение.
1. Система «вода — грунт». Расчет проводился на трехмерной прямоугольной сетке, количество узлов составляло 10 000, шаг по времени равнялся 10-3 с, шаг по пространству — 10 м.
На рисунке 1 представлено поведение волны с плоским фронтом при прохождении грунта на глубине 1000 м, не содержащего включений. Плотность воды полагалась равной 1000 кг/м3, скорость звука составляет 1500 м/с. Грунт имеет плотность 2500 кг/м3, продольная скорость распространения возмущения 4500 м/с, поперечная — 2250 м/с.
Отклик той же среды на сферическое возмущение представлен на рисунке 2.
Результат моделирования прохождения плоской волны через угле-водородосодержащий слой показан на рисунке 3. Использованы следующие параметры слоя: плотность 2000 кг/м3, продольная скорость распространения звука 4000 м/ с, поперечная — 1250 м/с. Глубина залегания нефтесодержащего пласта составляла 500 м.
Рис. 1. Прохождение плоской волны в системе «вода — грунт»
Рис. 2. Прохождение сферического возмущения в системе «вода — грунт»
Рис. 3. Прохождение плоской волны через углево-дородосодержащий слой
2. Система «лед — вода — грунт». Эта система решалась в двумерном приближении. Шаг по пространству составлял 0,2 м, шаг по времени — 310-5 с, было выполнено 6 тысяч шагов по времени. Рассматривался слой льда толщиной 4 м плотностью 917 кг/ м3. Скорость продольных волн во льду полагалась равной 3940 м/с, а поперечных —
3650 м/ с. Рассматривался слой воды толщиной 100 м, что соответствует Ленинградскому и Русановскому месторождениям Карского моря [12]. Скорость продольных волн в грунте — 6500 м/с, поперечных — 3700 м/ с, плотность задавалась равной 2500 кг/м3. Плотность воды равнялась 1000 кг/ м3, скорость звука — 1500 м/ с.
На рисунке 4 приведена картина сейсмически и гидроакустических волн, возникающих при точечном воздействии на верхнюю кромку льда импульсом Риккера, а на рисунке 5 — формирование отклика от резервуара с глубиной залегания 50 м.
11
Рис. 4. Сейсмические и гидроакустические волны в системе «лед — вода — грунт»
Рис. 5. Отклик включения в системе «лед — вода — грунт» от нефтесодержащего включения
Заключение
В данной статье было рассмотрено решение упругоакустических задач шельфовой сейсморазведки с помощью сеточно-характеристичес-кого метода. Было разработано граничное условие между частью области интегрирования, в которой решается система уравнений, описывающая состояние бесконечно малого объема сплошной линейно-упругой среды и той части области интегрирования, в которой решается система уравнений, описывающая акустическое поле.
Список литературы
1. Etter P. C. Underwater acoustic modelling and simulation. L., 2003.
2. Куликовский А. Г., Погорелов Н. В., Семенов А. Ю. Математические вопросы численного решения гиперболических систем уравнений. М., 2001.
3. Магомедов К. М., Холодов А. С. Сеточно-характеристические численные методы. М., 1988.
4. Иванов В. Д., Кондауров В. И., Петров И. Б., Холодов А. С. Расчет динамического деформирования и разрушения упругопластических тел сеточно-характе-ристическими методами // Матем. моделирование. 1990. Т. 2, № 11. C. 10 — 29.
12
5. Favorskaya A. V., Petrov I. B., Sannikov A. V., Kvasov I. E. Grid Characteristic Method Using High Order Interpolation on Tetrahedral Hierarchical Meshes with a Multiple Time Step / / Mathematical Models and Computer Simulations. 2013. Vol. 5, № 5. P. 409-415.
6. Голубев В. И., Петров И. Б., Хохлов Н. И. Численное моделирование сейсмической активности сеточно-характеристическим методом // Журнал вычисл. матем. и матем. физ. 2013. Т. 53, № 10. С. 1709 — 1720.
7. Новацкий В. К. Теория упругости. М., 1975.
8. Ландау Л. Д., Лифшиц Е. М. Теоретическая физика. М., 1987. Т. 7.
9. Harten Ami. High resolution schemes for hyperbolic conservation laws // Journal of Computational Physics. 1997. Vol. 135(2). P. 260-278.
10. Петров И. Б., Хохлов Н. И. Сравнение TVD лимитеров для численного решения уравнений динамики деформируемого твердого тела сеточно-характе-ристическим методом // Математические модели и задачи управления : сб. науч. тр. 2011. С. 104—111.
11. Roe P. L. Characteristic-Based Schemes for the Euler Equations // Annual Review of Fluid Mechanics. 1986. № 18. P. 337-365.
12. Новиков Ю. Н., Гажула С. В. Особенности оценки месторождений углеводородного сырья арктического шельфа России и их переоценки в соответствии с новой классификацией запасов // Нефтегазовая геология. Теория и практика. 2008. № (3). С. 1—19.
Об авторах
Алена Владимировна Фаворская — асп., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Дмитрий Игоревич Петров — асп., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Николай Игоревич Хохлов — канд. физ.-мат. наук, науч. сотр., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Игорь Борисович Петров — проф., д-р физ.-мат. наук, чл.-корр. РАН, Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
About the authors
Alena Favorskaya — PhD student, ass., Moscow Institute of Physics and Technology (State University), Dolgoprudny.
E-mail: [email protected]
Dmitry Petrov — PhD student, ass., Moscow Institute of Physics and Technology (State University), Dolgoprudny.
E-mail: [email protected]
Dr Nikolay Khokhlov — researcher, Moscow Institute of Physics and Technology (State University), Dolgoprudny.
E-mail: [email protected]
Prof. Igor Petrov — corresponding member of RAS, Moscow Institute of Physics and Technology (State University), Dolgoprudny. E-mail: [email protected]
