CC BY
36
24
УДК 519.6
А. В. Фаворская, Д. И. Петров, Н. И. Хохлов, И. Б. Петров
СОВМЕСТНОЕ МОДЕЛИРОВАНИЕ УПРУГИХ И АКУСТИЧЕСКИХ ВОЛН В УСЛОВИЯХ АРКТИКИ
Численно решаются задачи сейсмической разведки в условиях Арктики. Представлены результаты численного моделирования ударного воздействия на айсберг и сейсмического воздействия на нефтехранилище. Решаются полная система уравнений, описывающая состояние линейно-упругого тела и система уравнений, описывающая акустическое поле с помощью сеточно-характеристического метода, на границе ставится контактное условие между жидкостью и твердым телом.
Numerical solutions of seismic prospecting problems in the Arctic is presented. Results of numerical modeling of the explosive impact on the iceberg and earthquake impact on the oil storage are considered as well. The complete system of equations describing the state elastic body and a system of equations describing the acoustic field are solved using grid-characteristic method, the contact condition between the liquid and the solid is used on the border.
Ключевые слова: сеточно-характеристический метод, механика деформируемого твердого тела, акустическое поле, шельфовая сейсморазведка, иерархические сетки, сейсмостойкость.
Key words: grid-characteristic method, the mechanics of deformable solids, acoustic field, shelf seismic prospecting, hierarchical girds, seismic stability.
Введение
В последнее время актуальным стал вопрос освоения запасов углеводородов в Арктических морях [1; 2]. Соответственно, возникает необходимость точного компьютерного моделирования задач сейсмической разведки в шельфовых зонах. Немаловажным является вопрос безопасности судов и нефтедобывающих платформ в ледовые периоды. В связи с этим возникает вопрос раздробления айсбергов и интерес к численному моделированию данных процессов.
В данной работе приведены результаты численного эксперимента по моделированию волновых процессов в системах «лед — вода — грунт» при наличии нефтесодержащего включения сеточно-характе-ристическим методом [3 — 7], также представлены результаты численного моделирования ударного воздействия на айсберг и численного моделирования сейсмостойкости нефтехранилищ. Для моделирования поставленных задач решаются полная система уравнений, описывающих состояние сплошной линейно-упругой среды [8] и полная система уравнений, описывающая акустическое поле [9], при моделировании воздействия землетрясения на нефтехранилище используются иерархические сетки [10], расчет волновых процессов проводится непосредственно от эпицентра землетрясения до нефтехранилища.
© А. В. Фаворская, Д. И. Петров, Н. И. Хохлов, И. Б. Петров, 2015
Вестник Балтийского федерального университета им. И. Канта. 2015. Вып. 10. С. 24 — 31.
1. Постановка задачи
Компоненты скорости движения V и симметричного тензора напряжений Коши ст в линейно-упругой среде описываются следующей системой уравнений [8]:
р^б = (V-о)т, (1)
д т
— о = + ц(У®б +^®б)т). (2)
й
Для численного моделирования морской воды и нефтесодержащих включений использовалось приближение идеальной жидкости [9] и решалась полная система уравнений, описывающая акустическое поле давления р и компонент скорости V :
д
р V = ^р, (3)
дt
д 2
- р = -с 2р( V-V). (4)
дt
В (1), (2) коэффициенты X, ц — параметры Ляме, определяющие свойства упругого материала; а ® Ь — тензорное произведение векторов а и Ь, (а ® ЬУ = а'Ь'; I — единичный тензор второго ранга. В (4) за
с обозначена скорость звука в идеальной жидкости. Скорость продольных волн в линейно-упругой среде можно найти по формуле
Ср =^ , (5)
р
а скорость поперечных волн вычисляется в соответствии с
25
(6)
С, =
2. Численный метод
Для численного решения систем (1), (2) и (3), (4) используется се-точно-характеристический метод на криволинейных структурных сетках, позволяющий строить корректные численные алгоритмы для расчета граничных точек и точек, лежащих на поверхностях раздела двух сред с разными параметрами Ляме и (или) плотностями.
Опишем сеточно-характеристический метод на примере двумерного случая системы, описывающей акустическое поле. Систему (3), (4) в двумерном случае можно представить в следующем виде:
% + А1 ^ + А2 = 0. (7)
дt дх1 дх2
26
В выражении (7) под вектором ц понимается вектор, составленный
из двух компонент скорости и давления: Ц е {с1, с2, р}Т.
Вначале применяется метод расщепления по пространственным координатам, в результате чего имеем две одномерных системы:
^ = Л, . (8)
дt 7 дх] 4 У
Каждая из этих систем — гиперболическая и обладает полным набором собственных векторов с действительными собственными значениями, поэтому каждую из систем можно переписать в виде
дц = а- Ч, О, Л,
дt 1,1 дх,
где О, — матрица, составленная из собственных векторов; Л, — диагональная матрица, элементами которой являются собственные значения. Для всех координат матрица Л имеет вид (индекс , далее опускается, где это возможно): Л = diag{c, -с,0}.
После замены переменных р = ОЦ каждая из систем (8) распадается на три независимых скалярных уравнения переноса:
дР + Лдр = 0.
дt дх
Одномерные уравнения переноса решаются с помощью метода характеристик либо обычными конечно-разностными схемами.
После того как все компоненты перенесены, восстанавливается само решение:
< п+1 = О-1 р"+1.
Решение системы (3), (4) в двумерном случае, а также решение систем (1), (2) и (3), (4) в трехмерном случае проводится по аналогичным алгоритмам.
В программе реализовано применение ТУО-разностных схем 2-го порядка точности [11], 15 различных лимитеров [12], в расчетах в основном используются ограничитель 8ирегЪее [13] и сеточно-характеристи-ческие схемы 2 — 4 порядка точности [7].
Реализованые контактные условия с заданной внешней силой (в том числе частный случай свободной границы), с заданной скоростью границы, неотражающие граничные условия, контактные условия между жидкостью и твердым телом (описание которых можно посмотреть в [14]), а также контактные условия полного слипания и свободного скольжения.
3. Сейсморазведка в условиях Арктического шельфа
Для моделирования сейсмических процессов в условиях северных морей трехмерном случае были проведены расчеты на прямоугольной сетке на 5,4 • 107 узлов. Исследовался отклик от одиночного источника
импульса Рикера в многокомпонентной среде лед — вода — грунт — включение на поверхности льда. Ширина области интегрирования составляла 120 м, глубина — 60 м.
По бокам области интегрирования задавались неотражающие граничные условия.
Было выполнено 15 000 итераций, на каждой проводился шаг по времени, равный 3 -10-5 с. Шаг по координате по был постоянным и равнялся 0,4 м.
Были приняты следующие параметры среды: лед толщиной 4 м с плотностью 917 кг/ м3, скорость распространения продольных волн была взята равной 3940 м/с, поперечных — 2491 м/с. Слой воды имеет толщину 20 м с плотностью 1000 кг/м3 и скоростью звука 1500 м/ с. Для грунта плотности 2500 кг/м3 скорости продольных и поперечных волн равнялись соответственно 6500 м/ с и 3700 м/ с. Плотность включения взята равной 2000 кг/м3. Скорости продольной и поперечной волны суть 4000 и 1250 м/с.
На рисунке 1 изображены сейсмические волны в толще льда, а на рисунке 2 во всех слоях системы «лед — вода — грунт — включение — грунт».
27
Рис. 1. Сейсмические волны в толще льда
Рис. 2. Волновая картина в системе «лед — вода — грунт включение — грунт»
4. Ударное воздействие на айсберг
Система «лед — вода — грунт» решалась в трехмерной постановке.
Рассматривалось ударное воздействие на айсберг размером 100 на 200 м с плотностью, равной 917 кг/м3, скорость продольных волн, равной 3940 м/с, и скоростью поперечных волн 3650 м/с. Вокруг льда моделировался слой воды плотностью 1000 кг/м3 и скоростью звука 1500 м/с. На верхней и боковых границах айсберга, а также на верхней границе воды ставилось условие свободной границы, а между айсбергом и водой ставилось контактное условие, описание которого можно найти в [14]. На внешних границах слоя воды ставились неотражающие граничные условия. На верхнюю поверхность айсберга на площадке шириной 10 м действовала поверхностная плотность сил с амплитудой,
28
равной 0,9 МПа. Шаг по времени брался равным 2,5 -10 5 с, а шаг по пространству — 0,1 м. Расчет проводился в двумерном приближении.
На рисунке 3 приведена волновая картина, сформированная в айсберге и воде в момент времени, равный 0,17 с. Градацией серого показан модуль скорости, черным цветом показана зона локализации трещин в айсберге.
Рис. 3. Волновая картина в воде и айсберге при ударном воздействии на айсберг. Момент времени 0,17 с
5. Воздействие землетрясения на нефтехранилище
Рассматривалась область интегрирования размером 9 на 7 км. На дневной поверхности задавалось условие свободной границы, а по бокам области интегрирования — неотражающее граничное условие. Очаг землетрясения располагался на глубине 3 км по оси у и центрально по оси х. Грунт задавался с упругими характеристиками: скорость продольных волн 4500 м/с, скорость поперечных волн 2250 м/с, плотность 2500 кг/м3. Бетон фундамента задавался с упругими характеристиками: скорость продольных волн 3650 м/с, скорость поперечных — 1825 м/с, плотность — 2300 кг/м3. Шаг по времени брался равным 10-5 с во всех трех постановках, всего было проведено 200 000 шагов по времени. Плотность нефти бралась 900 кг/ м3, а скорость звука — 1470 м/с. Нефтехранилище находилось на 1 км по дневной поверхности правее очага землетрясения. В грунте вокруг нефтехранилища задавалась система из 6 вложенных иерархических сеток, изображенная на рисунке 4. В самой крупной сетке шаг по пространству равнялся 4 м, в самой мелкой — 0,0625 м. В соседних сетках шаги отличались в 2 раза.
Рис. 4. Система вложенных иерархических сеток
На рисунке 5 представлено воздействие продольных сейсмических волн на нефтехранилище в момент времени 0,71 секунды, а на рисунке 6 оно же в увеличенном масштабе.
29
Рис. 5. Воздействие продольных сейсмических волн на нефтехранилище в момент времени 0,71 с
Рис. 6. Воздействие продольных сейсмических волн на нефтехранилище в момент времени 0,71 с. Увеличенный масштаб
На рисунке 7 изображено воздействие поперечных сейсмических волн на нефтехранилище в момент времени 1,43 с, а на рисунке 8 оно же в увеличенном масштабе.
Рис. 7. Воздействие продольных сейсмических волн на нефтехранилище в момент времени 1,43 с
30
Рис. 8. Воздействие продольных сейсмических волн на нефтехранилище в момент времени 1.43 с. Увеличенный масштаб
Можно сделать вывод о том, что использование иерархических сеток позволяет проводить численное моделирование сейсмических воздействий на нефтехранилища, проводя расчет сейсмических волн непосредственно от очага землетрясения.
Заключение
Было рассмотрено моделирование волновых процессов в слоистых средах с совместным решением систем уравнений, описывающих линейно-упругое тело и акустическое поле. Рассматривались задачи сейсмической разведки в условиях Арктического шельфа, задача ударного воздействия на айсберг и моделирования воздействия землетрясения на нефтехранилище с использованием системы иерархических сеток. Во всех расчетах использовался сеточно-характеристический метод.
Результаты получены в рамках работы по гранту РФН №14-11-00263.
Список литературы
1. Etter P. C. Underwater acoustic modelling and simulation. L., 2003.
2. Новиков Ю. Н., Гажула С. В. Особенности оценки месторождений углеводородного сырья арктического шельфа России и их переоценки в соответствии с новой классификацией запасов // Нефтегазовая геология. Теория и практика. 2008(3). С. 1-19.
3. Куликовский А. Г., Погорелов Н. В., Семенов А. Ю. Математические вопросы численного решения гиперболических систем уравнений. М., 2001.
4. Магомедов К. М., Холодов А. С. Сеточно-характеристические численные методы. М., 1988.
5. Иванов В. Д., Кондауров В. И., Петров И. Б., Холодов А. С. Расчет динамического деформирования и разрушения упругопластических тел сеточно-характе-ристическими методами // Матем. моделирование. 1990. Т. 2, № 11. C. 10 — 29.
6. 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.
7. Голубев В. И., Петров И. Б., Хохлов Н. И. Численное моделирование сейсмической активности сеточно-характеристическим методом // Журнал вычисл. матем. и матем. физ. 2013. Т. 53, № 10. С. 1709-1720.
8. Новацкий В. К. Теория упругости. М., 1975.
9. Ландау Л. Д., Лифшиц Е. М. Теоретическая физика. М., 1987. Т. 7.
10. Petrov I. B., Favorskaya A. V., 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, No. 5. P. 409-415.
11. Harten Ami. High resolution schemes for hyperbolic conservation laws // Journal of Computational Physics. 1997. Vol. 135(2). P. 260-278.
12. Петров И. Б., Хохлов Н. И. Сравнение TVD лимитеров для численного решения уравнений динамики деформируемого твердого тела сеточно-характе-ристическим методом // Математические модели и задачи управления : сб. науч. тр. 2011. С. 104-111.
13. Roe P. L. Characteristic-Based Schemes for the Euler Equations // Annual Review of Fluid Mechanics. 1986. № 18. P. 337-365.
14. Петров И. Б., Фаворская А. В., Петров Д. И., Хохлов Н. И. Численное решение арктических задач с помощью сеточно-характеристического метода // Известия ЮФу. Технические науки. 2014. Вып. 12. С. 192 — 200.
Об авторах
Алена Владимировна Фаворская — асп., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Дмитрий Игоревич Петров — асп., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Николай Игоревич Хохлов — канд. физ.-мат. наук, науч. сотр., Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
Игорь Борисович Петров — проф., д-р физ.-мат. наук, чл.-корр. РАН, Московский физико-технический институт, Долгопрудный.
E-mail: [email protected]
About the authors
31
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]
