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

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ имени М.В. Ломоносова

На правах рукописи Друца Александр Валерьевич

Конечно-разностные методы решения уравнений мелкой воды на неструктурированных сетках

01.01.07 — Вычислительная математика

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

Москва, 2013

005051779

005051779

Работа выполнена на кафедре вычислительной математики Механико-математического факультета Московского государственного университета имени М.В. Ломоносова

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

профессор Кобельков Георгий Михайлович

Официальные оппоненты: Дианский Николай Ардальенович,

Ведущая организация: ФГБУН Институт прикладной математики

имени М.В. Келдыша РАН.

Защита состоится 6 марта 2013 года в 16 ч. 45 мин. на заседании диссертационного совета Д 501.002.16 при Московском государственном университете имени М.В. Ломоносова по адресу: Российская Федерация, 119991, Москва, ГСП 1, Ленинские горы, д. 1, Механико-математический факультет, аудитория 14-08.

С диссертацией можно ознакомиться в научной библиотеке Московского государственного университета имени М.В. Ломоносова. Автореферат разослан «6» февраля 2013 г.

Ученый секретарь

диссертационного совета Д 501.002.16 доктор физико-математических наук,

профессор Корнев Андрей Алексеевич

доктор физико-математических наук (ФГБУН Институт вычислительной математики РАН, ведущий научный сотрудник)

Кукаркин Алексей Борисович кандидат физико-математических наук, доцент (Научно-исследовательский вычислительный центр МГУ имени М.В. Ломоносова, старший научный сотрудник)

/

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

Актуальность темы. В настоящее время одним из важнейших направлений современной вычислительной математики является моделирование различных природных процессов, связанных с перемещением потоков жидкости в мировом океане. К числу таких процессов можно отнести возникновение приливных волн или цунами. Основной математической моделью, описывающей динамику приливных волн, является система уравнений динамики мелкой воды. Эта система также является частью математической модели динамики океана, она моделирует краевые условия на поверхности океана. Система уравнений динамики приливов является упрощением известной системы уравнений Навье— Стокса при условии, что горизонтальный масштаб много больше вертикального. Вывод уравнений мелкой воды из уравнений крупномасштабной динамики океана в сг-системе координат можно найти в работе В.Б. Залесного1. Теоретическое исследование этих уравнений проводилось в работе Г.И. Марчука и Б.А. Кагана2. Для модели крупномасштабной динамики океана в работах Г.М. Кобелькова3 была доказана теорема существования решения в «целом».

Исследование системы уравнений мелкой воды ведётся не один десяток лет. Для решения данной проблемы использовались такие методы, как метод конечных разностей, метод расщепления, схема Годунова4, поверхностный градиентный метод, метод дробных шагов5 и др. Линеаризованная система уравнений динамики приливов при условии рав-

1 V.D. Zalesny, Mathematical model of яра dynamics in a о -coordinate system, Russ. J. Numer. Anal. Math. Modelling, V. 20, N. 1, pp. 87-113, 2005.

2Марчук Г.И., Каган Б.А. Океанские приливы, Л.: Гидрометеоиздат, 1977.

3Kobelkov G.M.. Exi stence of a solution "in the large" for ocean dynamics equations, J. math, fluid inech., 9, pp. 588-610, 2007.

'Eymard H., Gallouet Т., Herbin R., The finite volume method. Handbook of Numerical Analysis, 2000, Vol. VII, pp. 713-1020

"Marchuk G.I., Gordeev R.G., Rivkind V.Y., Hagan B.A., A numerical method for the solution of tidal dynamics equations and the results of its application, Journal of Computational Physics, 1973, Vol. 13, №1, pp. 15-34

номерности дна тесно связана с задачей движения вязкого слабосжи-маемого газа, поэтому может быть решена методами, применяемыми в газовой динамике6. Все вышеперечисленные подходы к решению поставленной задачи аппроксимируют уравнения на равномерных структурированных сетках. При моделировании приливных волн зачастую необходимо строить триангуляцию, которая учитывает наличие большого количества островов разного размера. Размеры самого моря и островов могут сильно различаться по величине. Например, размеры Охотского моря составляют примерно 2200 х 1500 км, а размеры островов в диаметре составляют от одного километра до десятков километров. Помимо небольших островов, необходимо учитывать наличие узких проливов, ширина которых может составлять несколько километров. При построении триангуляции потребуется, как минимум, два элемента для учёта таких проливов, то есть размер элемента сетки должен быть примерно 2-3 км. Следовательно, равномерная сетка должна быть очень мелкой, а количество элементов в ней для Охотского моря будет превышать полмиллиона. Поэтому обычно решение ищут на неравномерных сетках. При этом такие сетки, как правило, являются неструктурированными. При работе с неструктурированными сетками аппроксимация уравнений строют с помощью методов конечных элементов или методов конечных объёмов. Аппроксимация системы динамики мелкой воды при помощи классических конечных элементов сводит задачу к конечномерной, в которой отсутствует сохранение баланса на ячейке (то есть количество втекающей в ячейку жидкости не соответствует количеству вытекающей). Это приводит к тому, что через некоторое время получаемое в результате расчёта решение существенным образом отличается от решения дифференциальной задачи. Таким образом требование сохранения баланса должно точно выполняться для каждой ячейки сетки. В

вЖуков К.А., Попов A.B., Исследование экономичной разностной схемы для нестационарного движения вязкого слабоежнмасмого газа, Ж. вычисл. матем. и матем. физ., 2005, Том 45, Л"' 4, стр. 701-717

работе Равьяра П.А. и Тома Дж.М.7 для подобной аппроксимации были использованы специальные неконформные конечные элементы (элементы Равьяра-Тома), которые обеспечивают сохранение баланса на ячейке. В работе К.Ю. Богачёва и Г.М. Кобелькова8 такие элементы были применены для аппроксимации задачи приливных волн. Проведённые численные эксперименты показали эффективность этого подхода для решения рассматриваемой системы уравнений. Однако в этом случае элементы являются разрывными функциями, и вычислительные формулы для них довольно сложны. В связи с этим в диссертационной работе был выбран метод построения конечномерной задачи уравнений мелкой воды, основанный на конечно-разностной аппроксимации на неструктурированной сетке.

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

'Haviart P. A., Thomas Л. M., A Mixed Finite Element Method for Second Order Elliptic Problems, Lecture Notes in Mathematics 1977, Vol. 606, pp. 292-315

aBogachev K.Yu., Kobelkov G.M, Numerical solution of a tidal wave problem, ill Proceedings of "Parallel Computational Fluid Dynamics", v.2, 2004, 163-173, J.-Wiley Press

°Popov I.V., Fryazinov I.V., Stanichenko M.Yu., Taimanov A.V., Construction of a difference, scheme for Navies-Stokes equations on unstructured grids, Russ. J. Numer. Anal. Math. Modelling, 2009, v.23, No. 5, 487-503

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

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

Цели работы

1. Построение аппроксимации линеаризованных уравнений динамики мелкой воды на неструктурированных сетках в декартовых и сферических системах координат.

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

3. Исследование сходимости разностных схем в декартовых и сферических системах координат.

4. Построение аппроксимаций нелинейных и вязких уравнений динамики мелкой воды на неструктурированных сетках в декартовой системе координат.

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

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

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

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

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

3. Проведены численные эксперименты на неструктурированных сетках, в том числе для реальных географических объектов, подтверждающие эффективность предложенного подхода. Помимо этого результаты расчётов показали, что наличие порядка 15% тупоугольных треугольников в сетке является допустимым, что позволяет использовать стандартные генераторы сеток.

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

Методы исследований. Была использована методика построения

аппроксимации на неструктурированных сетках, разработанная И.В. Фрязиновым и соавторами10. Для получения априорных оценок применялась методика Г.М. Кобелькова. Кроме того, применялись методы дифференциальной геометрии, общей теории решения систем линейных алгебраических уравнений, а также методы анализа разностных схем. Достоверность теоретических результатов подтверждена детальными численными экспериментами для физически осмысленных задач.

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

• международные конференции молодых ученых «Ломоносов-2009» и «Ломоносов-2011»;

• международная конференция «Современные проблемы математического моделирования» в Абрау-Дюрсо (Россия, 2011 г.);

• международная научная конференция «Современные проблемы вычислительной математики и математической физики», посвященная памяти академика А.А.Самарского (Россия, 2009 г.).

Они неоднократно докладывались и обсуждались на научно-исследовательских семинарах:

• научно-исследовательский семинар кафедры вычислительной математики механико-математического факультета МГУ им. М.В. Ломоносова под руководством д.ф.-м.н. профессора Г.М. Кобелькова (неоднократно в 2009—2012 годах);

• научно-исследовательский семинар ИВМ РАН «Вычислительная математика, математическая физика, управление» под руководством д.ф.-м.н. профессора Г.М. Кобелькова и д.ф.-м.н. профессора A.B. Фурсикова (2010 г.).

10см. примечание 9 на стр. 3

• семинар научно-исследовательского вычислительного центра МГУ имени М.В. Ломоносова под руководством д.ф.-м.н, профессора Я.М. Жилейкина (2012 г.).

• семинар научно-исследовательского вычислительного центра МГУ имени М.В. Ломоносова под руководством чл.-корр. РАН, профессора В.Н. Лыкосова (2011 г.).

Публикации. По теме диссертационной работы опубликовано 4 работы, 3 из которых в журналах из «Перечня ведущих рецензируемых журналов и изданий, в которых должны быть опубликованы основные научные результаты диссертационной работы на соискание степени доктора и кандидата наук» [1-3].

Структура диссертационной работы. Работа состоит из введения, четырёх глав, разбитых на параграфы, заключения и списка литературы. Список литературы включает 33 наименования. Объём диссертационной работы составляет 105 страниц.

Содержание диссертации

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

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

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

центр описанной окружности треугольника, потоковым узлом называется пересечение отрезка, соединяющего центры соседних ячеек, и общей стороны двух ячеек. Поскольку все треугольники остроугольные, центры ячеек лежат внутри каждого треугольника, и потоковый узел — это середина стороны, на которой он лежит. Задаётся нумерация центров ячеек к = 1,2,...,К и нумерация потоковых узлов г = 1,...ЛГ. Центр ячейки с номером к обозначается Ок, а потоковый узел с номером г —

Двойным верхним индексом к, а. обозначаются различные элементы сетки, связанные с потоковым узлом Хг, который лежит на стороне ячейки А; с номером а (то есть Х{ = Хк,а). Двойной нижний индекс г,т обозначает элементы, связанные с ячейкой с номером к, которая содержит потоковый узел Х{ (то есть Оц = Ок). На рисунке 1 изображён узел Ок и основные элементы сетки: длина стороны - 1к,а, длина отрезка, соединяющего центры, - йк'а, площадь ячейки -

вк. На рисунке 2 изображён узел Х^ и отмечены Рис' 11 ^зел ^ и три

соседних узла О а — элементы сетки, связанные с ним: длина ребра _ 1 2 3

и, площадь четырёхугольника - .

Через О обозначено множество центров ячеек, а через X — множество потоковых узлов. Определяются пространства скалярных и векторных сеточных функций на множествах О и X. Сеточные скалярные произведения задаются как

(/,<?) = £ 5*/*-<Д (С,0 = 1>й-й, *=1 ¿=1

где функции /, д заданы в центрах ячеек, а £ — в потоковых узлах.

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

ркЛ

|в*'"| - й1'

заданной в потоковых узлах по формуле

3

о* „=,

где пк,а - внешняя нормаль к стороне 1к'а треугольника Ок, — значение функции £ в точке Хк<а. Компоненты вектора нормали пк'а по а; будем обозначать через пк'а, по у — через пк,а.

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

«Ну* «I = 4 (-Иі,іПі,і - иіі2пі2), (2)

І Лі Оі

где и*,т — значение функции и в точке 0;,т.

В работе И.В. Фрязинова11 показано, что введённые выше операторы сопряжены, то есть (иЛ, УЛСл) = — (сііул ил, Сл) и имеют порядок аппроксимации О(Н); здесь и далее через к обозначен максимальный размер ячейки.

Рис. 2: Потоковый узел

Л", И две соседние ячейки Аппроксимация оператора Лапласа в центрах с центрами ОІЛ и 0(,2. ячеек производится по формуле

з

Ли

ДЛи

- -2_ и а=1

¿к,а

к = 1,К.

(3)

Аппроксимация нелинейного оператора

о*

= ЛГ(и,у)

ок

Лгл(и,у)

о*

(4)

а=1

где к = 1,К, а (у^а)к'а — линейная интерполяция нормальных к 1к-а компонент скорости V в точку Хк'а.

см. примечании 9 на стр. 3

В работе доказана следующая теорема: Теорема 1. Для сеточного оператора дивергенции (Ну'1 верно следу-

ющее соотношение, которое является некоторым аналогом формулы Гаусса-Остроградского:

В работе уделено внимание тому, что для применения вышеописанных аппроксимаций требуется триангуляция, состоящая только из остроугольных треугольников. Существуют методики построения таких сеток12, однако на практике, в силу нетривиальности данной задачи, применяются «универсальные» программные решения, которые строят триангуляции, близкие к требуемой: Тпаг^1е, аш2с1, {^твЬ. В работе введены оценки такой близости сеток: количество тупоугольных треугольников и отношение площади самой большой ячейки сетки к самой маленькой. В конце главы описаны два подхода к построению аппроксимации на «плохих» (тупоугольных треугольниках) ячейках, одним из которых является рассмотрение в качестве центра ячейки тупоугольного треугольника центра масс. В дальнейших главах результаты численных экспериментов показывают, что такая аппроксимация не влияет существенным образом на результат.

Вторая глава диссертационной работы посвящена аппроксимации линейных уравнений мелкой воды на неструктурированной сетке в декартовой системе координат.

Рассматриваются следующие уравнения на двумерной ограниченной области Г2 с границей Г = Гх У Гг:

1ги. в. Попов, С. В, Поляков, Построение адаптивних нерегулярных треугольных сеток для двумерных многосвязных невыпуклых областей, Матем. моделирование, 2002, том 14, К* б, стр.

ди , дС „

дх

(5)

25-35

ft; , 3( „ , ,ч

_ = + + (6)

üí = ~ 4. (74

dt дх ду'

где и, v - горизонтальные компоненты скорости и, а С - высота волны. R,X,g — константы. Граничные и начальные условия записываются в виде:

и • п = ип\ + vn2 = 0 на Гг, (8)

С = 0 наГ2; (9)

и{х, у, 0) = ио(х, у), v(x, у, 0) = щ(х, у), ((ж, у, 0) = Са(г, у). (10)

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

^ = АХ х ц* - Rük + f + gVhC, к = ТГК, (11)

где k — единичный вектор вдоль оси Oz, f = (/ь/г). Аппроксимация уравнения неразрывности (7) во внутренних точках имеет вид

= divhü, г = 1777 и Xi е ПЛ\Г\ (12)

В граничных узлах JV¡ 6 Г^ уравнение неразрывности модифицируется следующим образом:

^^ = div^ü = А (-й4,щи), i = ljy и X¡ 6 rf. (13) т o¡

Для оставшихся узлов X¡ 6 Г§ уравнения имеют вид краевого условия

(9): _

Ci = 0, i = l,N и Xj 6 Г2. (14)

Уравнения (11)—(14) аппроксимируют задачу (5)—(6) и граничные условия (8)-(9).

Путём исключения неизвестных система (11)—(14) сводится к решению системы линейных алгебраических уравнений относительно Матрица этой системы имеет вид 6 2

(ЛаС)І = ~ СО, і = 1

где — множество соседних с Хі потоковых узлов,

Му =

(Д + 1/т)дсгё(фи ¡) + дХ щп{ч{ХіО і,ті Оі,тХ]У) А2 + (Д + 1/г)2

(15)

(16)

где у( а^, — двумерное векторное произведение векторов и ї/. Шаблон для оператора (15) изображён на рисунке 3.

Рис. 3: Шаблон для Л/,. Мі = {Х^.^ Х^, Х^}

В работе доказана следующая теорема для коэффициентов цц системы (15):

Теорема 2. Если все треугольники остроугольные, тогда существует т* > 0 такое, что для всех г, 0 < т < т*, коэффициенты Ціі удовлетворяют условию //у > 0.

С помощью предыдущей теоремы в диссертационной работе доказан следующий результат о матрице системы линейных алгебраических уравнений (15):

Теорема 3. Пусть выполнены условия теоремы 2. Тогда оператор Ад, задаваемый формулой

МІ = - £ ~ б)-

Х3-єП(

положительно определён, и в нём имеет место диагональное преобладание.

Из теоремы 3 следует, что оператор (15) является положительно определённой М-матрицей. Таким образом, разностная задача разрешима. Помимо разрешимости системы, в диссертационной работе доказана устойчивость решения по начальным данным, то есть для всех п ^ О верно

!|иТ + 5Г112^1К||2 + 5||Со112.

где и", £" - значения функций на п-ом слое по времени.

Для построенной разностной схемы доказан следующий результат.

Теорема 4. Пусть существует гладкое решение дифференциальной задачи (5)—(10) и углы всех ячеек меньше п/2 - (1 при некотором положительном ц, тогда решение системы (11)—(Ц) сходится к решению задачи (5)—(10) с порядком 0(т + /г).

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

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

Рассматривается следующая система уравнений с нелинейным членом и вязкостью

и* = 5VС - Яи - \к х и - (и, У)и 4- { + иАи, (17)

С< = ШуЯи, (18)

где R, A, g, и - константы. Н(х,у) € С(П) - глубина, функция координат. Предполагается, что Я » |Ç|. Граничные и начальные условия такие же, как (8)—(10) второй главы .

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

N (u. u)o* = -дъ 2Lj -dM--2-•

a=l

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

= Хк х û* - Rûk + fk + je J2 1к'%к'апк'а

a~\

_ 1 v> П Д ,к.сйк-°-й1< (19)

S1 ¿J ' 7F» 2 Tp^ dk-a >

a=l a-1

к=1,К.

Во внутренних узлах Хг аппроксимация уравнения неразрывности имеет вид

= 4 ("Ямйищл - Яюй^зпу), г = ТТШ и е (20)

т

Для узлов на границе Г^ имеем

^^ = 4 (-Яц^ццд), » = 1737 и € Г?, (21)

т Ь¡ а для узлов на границе :

С» = 0, г = 1777 И X, е Г* (22)

На каждом временном слое система (19)—(22) является системой линейных алгебраических уравнений относительно й и Для полученной системы проведены численные эксперименты, в которых сравнивалось решение нелинейной задачи с решением линейной; моделировалось

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

Далее в этой главе рассмотрена задача с нелинейным уравнением неразрывности:

Се = сЩН + С)и. (23)

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

Задача состоит в аппроксимации уравнений (17), (23) с граничными и начальными условиями (8)—(10) и в проведении численных экспериментов, направленных на исследование влияния нелинейности. Уравнение (17) аппроксимируется, как и в предыдущем случае, с помощью (19). Во внутренних узлах Xi имеем

^^ = (-(Ям + (Скрипи - (Я,2 + (СЬКгпЫ , (24)

Т

где (С)г,ш — линейная интерполяция значений С в центр ячейки ¿,т. А для узлов на границе Г]1:

^ = + < СЫад,,). (25)

В оставшихся точках £ Г£ аппроксимация имеет вид (22). Система (19),(24),(25),(22) также является системой линейных алгебраических уравнений относительно й и

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

Четвёртая глава диссертационной работы посвящена аппроксимации линейных уравнений мелкой воды на неструктурированной сетке в

сферической системе координат при неровном дне. Уравнения мелкой воды в сферической системе координат имеют вид:

^ = Аі; + р^-Ли + /1, (26)

^ = -Хи + дд^-П.у + /2, (27)

SC Jd(Hu) d(*Hv)

где 0 < х < 2тт, 0 < у < тг - координаты на сфере, q — 1/Rq, р — q/s\n(y), Rq — радиус Земли. Коэффициент силы Кориолиса имеет вид А = 2wcos(y). На функцию глубины ставятся те же условия, что и в первой задаче предыдущей главы. В работе считается, что область Г2 не содержит полюса с некоторой окрестностью радиуса 5. Уравнения (26)—(28) дополняются граничными и начальными условиями:

и • п = uni + -vП2 = 0 на Ti ; (29) Р

Ç = 0 наГ2; (30)

и{х,у,0) = и0(х,у), v(x,y,0) = v0{x,y), ф,у,0) = Ço{x,y). (31)

Для этой системы строится чисто неявная по времени разностная схема с шагом т, с использованием сеточных операторов (1) и (2)

= Акк xûk- Rûk + f |о* + gpkWkVhÇ, к = Tj<, (32) где W - диагональная матрица размера 2x2:

W= (1 °

V° î

и W* = ~W(Ok). Уравнение неразрывности (28) во внутренних точках аппроксимируется следующим образом:

^^ = pjdivл (tfWû), i = TJf и Xi 6 Çp\I*. (33)

т

В граничных узлах в Г^ аппроксимация выглядит как

^ = ^ (-Ям\Умиипм), г = 1777 и Х± 6 Г]1, (34)

а для оставшихся узлов 6 Г§ граничные условия имеют вид

6 = 0, г = 1777 и X; <5 Г§. (35)

Путём исключения неизвестных система (32)—(35) сводится к решению системы линейных алгебраических уравнений относительно С- Матрица этой системы имеет вид

СлХс)< = § - ^ £ 25(6 - * = йУ, (36)

где — некоторые коэффициенты, зависящие от углов треугольника, шага по времени т, коэффициента Кориолиса, числа Рейнольдса и ординаты у. Для коэффициентов системы (36) доказана Лемма 1. Пусть все треугольники Г2Л остроугольные и пусть выполнено дополнительное уыовие на углы треугольников сетки:

> \ - ып(*)) , V*. а ф Р. (37)

Тогда существует т* > 0 такое, что для всех т, 0 < т < г", коэффициенты положительны. Из леммы следует

Теорема 5. Пусть выполнены условия лемльы 1. Тогда оператор Л/, из (36) положительно определён, и в нал« имеет место диагональное преобладание.

Таким образом, разностная задача (32)—(35) разрешима. Помимо разрешимости системы, доказана устойчивость решения по начальным данным, то есть для всех тг ^ 0 верно

кнг+гнинионг+м,

17

где || • ||,, || • ||0 — некоторые нормы, эквивалентные стандартным сеточным нормам Ьг-

В итоге, для построенной разностной схемы доказан следующий результат:

Теорема 6. Пусть существует гладкое решение задания (26)—(28), (29)—(31), выполнены условия леммы 1 и углы всех ячеек меньше тх/2— ¡1 при некотором положительном ¡1, тогда решение системы (32)— (35) сходится к решению дифференциальной задачи с порядком 0(т + к).

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

Основные результаты, выносимые на защиту

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

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

2. Построена конечно-разностная аппроксимация уравнений мелкой воды с нелинейными членами и вязкостью на неструктурированной

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

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

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

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

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

Работы автора по теме диссертации

1. Арушанян И.О., Друца A.B., Кобельков Г.М. Метод конечных разностей для решения системы уравнений динамики приливов. Дифференциальные уравнения, 2009, том 45, №7, стр. 965-972. (Дру-це A.B. принадлежит аппроксимация уравнений предложенной схемой; доказательство положительной определённости матрицы системы уравнений для высоты волны.)

2. Drutsa A.V., Kobelkov G.M. Finite difference approximation of tidal wave equations on unstructured grid in spherical coordinates. R.uss. J. Nume.r. Math, and Math. Model, 2010, том 25, №6, стр. 535-544. (Друце A.B. принадлежит аппроксимация уравнений в сферической системе координат и исследование полученной системы уравнений для высоты волны.)

3. Друца A.B. Конечно-разностный метод для решения нелинейной системы уравнений динамики мелкой воды на неструктурированной сетке. Вычислительные методы и программирование, 2012, том 13, №2, стр. 511-516.

4. Арушанян И.О., Друца A.B., Кобельков Г.М. Метод конечных разностей для решения системы уравнений мелкой воды на неструктурированной сетке. Сборник лекций, прочитанных для участников молодёжной школы-конференции, проводимой в рамках Восьмой Всероссийской конференции "Сеточные методы для краевых задач и приложения", посвящённой 80-летию со дня рождения А.Д. Ляшко, Казань, стр. 50-68, 2010. (Друце A.B. принадлежит аппроксимация уравнений в различных системах координат и исследоваг ние полученных систем уравнений.)

Отпечатано в отделе оперативной печати Геологического ф-та МГУ Тираж ¡ОС экз. Заказ № Д

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ

имени М.В. Ломоносова Механико-математический факультет

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

04201357124

Друца Александр Валерьевич

Конечно-разностные методы решения уравнений мелкой воды на неструктурированных сетках

01.01.07 — Вычислительная математика

Диссертация

на соискание учёной степени кандидата физико-математических наук

Научный руководитель:

доктор физико-математических наук профессор Кобельков Георгий Михайлович

Москва, 2012

Оглавление

Введение 4

1 Основные обозначения и определения 10

1.1 Основные обозначения ........................10

1.2 Сеточные функции ..........................13

1.3 Аппроксимация операторов на сетке.................15

1.4 Свойства сеточных операторов....................16

1.4.1 Согласованность сеточных операторов дивергенции и градиента ..............................16

1.4.2 Сеточный аналог формулы Гаусса—Остроградского .... 17

1.5 Вспомогательные утверждения....................19

1.6 Неструктурированные сетки.....................19

1.7 Нерегулярные элементы сетки....................21

1.8 Выводы.................................23

2 Линеаризованные уравнения мелкой воды в декартовой системе координат 24

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

2.2 Разностная аппроксимация уравнений ...............25

2.3 Преобразование сеточной задачи...................26

2.4 Свойства сеточной задачи ......................31

2.5 Итерационный процесс........................34

2.6 Устойчивость решения по начальным данным...........35

2.7 Сходимость...............................36

2.8 Численные эксперименты.......................41

2.8.1 Сравнение с аналитическим решением............41

2.8.2 Численная оценка порядка сходимости............43

2.8.3 Сохранение баланса.......................45

2.8.4 Сила Кориолиса.........................47

2.8.5 Расчёт на реальной географической области.........47

2.9 Выводы.................................54

3 Нелинейные уравнения мелкой воды в декартовой системе координат 55

3.1 Постановка задачи для уравнений мелкой воды в декартовой системе координат с вязкостью и нелинейным членом........55

3.2 Сеточная задача............................56

3.3 Численные эксперименты.......................57

3.3.1 Влияние нелинейности.....................57

3.3.2 Тупоугольные треугольники..................59

3.3.3 Многосвязная область.....................61

3.3.4 Расчёт на реальной географической области.........63

3.4 Постановка задачи с нелинейным уравнением неразрывности . . 68

3.5 Аппроксимация уравнений......................68

3.6 Численные эксперименты.......................70

3.6.1 Неравномерное дно.......................70

3.6.2 Набегание волны на берег...................73

3.7 Выводы.................................76

4 Линеаризованные уравнения мелкой воды в сферической системе координат 77

4.1 Постановка задачи..........................'. 77

4.2 Сеточные функции ..........................78

4.3 Разностная аппроксимация уравнений ...............79

4.4 Преобразование сеточной задачи...................80

4.5 Свойства сеточной задачи ......................85

4.6 Устойчивость решения по начальным данным...........90

4.7 Сходимость...............................91

4.8 Численный эксперимент для реального географического объекта 95

4.9 Выводы.................................99

Заключение 100

Список литературы 102

Введение

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

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

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

масштабной динамики океана в а-системе координат можно найти в работе В.Б. Залесного [6]. Существование решения для крупномасштабной модели динамики океана в «целом» доказана в работах [18, 17].

Систему уравнений мелкой воды рассматривали во многих работах, в которых для нахождения решения использовали различные подходы: метод конечных разностей [23, 22, 26], метод расщепления [25], схема Годунова [27, 23], поверхностный градиентный метод [28], метод дробных шагов [24] и др. Линеаризованная система уравнений динамики приливов при условии ровного дна может быть решена методами, применяемыми в газовой динамике, как это сделано, например, в работах [16, 15]), поскольку данная задача очень похожа на задачу о движении вязкого слабосжимаемого газа. Другое направление в решении поставленной задачи — применение теории оптимального управления. Данное направление было разработано В.И.Агошковым и Е.А.Ботвиновским [4]. Все вышеперечисленные подходы аппроксимируют уравнения на равномерных структурированных сетках. Однако при моделировании приливных волн на реальных географических объектах часто необходимо строить триангуляцию, которая учитывает наличие большого количества островов разного размера. Например, размеры Охотского моря и островов сильно различаются. Размеры Охотского моря составляют примерно 2200 х 1500 км, а размеры островов в диаметре составляют от одного километра до десятков километров. Например, остров Большой Шантар — 65 х 47 км, а остров Утичий — 2, 5 х 0, 5 км. Помимо небольших островов, в Охотском море присутствуют узкие проливы. Например, ширина Татарского пролива, соединяющего Охотское и Японское моря, между островом Сахалин и материком в самом узком месте составляет 7, 3 км. При построении сетки потребуется, как минимум, два элемента для учёта данного пролива, то есть размер элемента сетки должен быть около 2-3 км. В таких условиях равномерная сетка будет очень мелкой, а количество элементов в ней для всего Охотского моря будет превышать полмиллиона. По этой причине такие области часто аппроксимируют неравномерной сеткой. Неравномерные структурированные сетки рассматривали в работе [19]. Однако, как правило, такие сетки являются неструктурированными. Также стоит отметить, что рассматриваемые области существенно многосвязны, то есть количество внутренних «дырок»

(островов) в них достаточно велико (например, для Охотского моря — до 300 островов).

При работе с неструктурированными сетками аппроксимацию уравнений строят с помощью методов конечных элементов [1] или методов конечных объёмов [20]. Аппроксимация системы динамики мелкой воды при помощи классических конечных элементов сводит задачу к конечномерной, в которой отсутствует сохранение баланса на ячейке (то есть количество втекаемой в ячейку жидкости не соответствует количеству вытекаемой). Это приводит к тому, что через некоторое время получаемое в результате расчёта решение начинает существенным образом отличаться от решения дифференциальной задачи. Очевидно, что требование сохранения баланса должно строго выполняться для каждой ячейки сетки. В работе [21] для аппроксимации были использованы специальные неконформные конечные элементы (элементы Равьяра-Тома), которые обеспечивают сохранение баланса на ячейке. В работе [3] такие элементы были применены для аппроксимации задачи приливных волн. Проведённые численные эксперименты показали эффективность этого подхода для решения рассматриваемой системы уравнений. Однако в этом случае элементы являются разрывными функциями, и вычислительные формулы для них довольно сложны. Учитывая вышеизложенное, в диссертационной работе был выбран метод построения конечномерной задачи уравнений мелкой воды, основанный на конечно-разностной аппроксимации на неструктурированной сетке.

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

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

Во второй главе рассмотрена линеаризованная система уравнений динамики мелкой воды с ровным дном на двумерной области в декартовой системе координат:

ди л д(

^ = -Ли + д^ - Яу + /2,

дС, ди ду дх ду'

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

с расчётными данными; подтверждение теоремы о сходимости; моделирование распространения волны цунами на реальном географическом объекте.

Третья глава диссертации посвящена нелинейным уравнениям мелкой воды. В первой части главы рассмотрены уравнения мелкой воды с нелинейной частью (и, \7)и, вязким членом иАи и неровным дном. Для этой задачи строится линеаризованная по временным слоям конечно-разностная схема. Для полученной системы проведены численные эксперименты, в которых сравниваются решения в случае нелинейной и линейной задач, моделируется распространение волны на реальном географическом объекте (Индийский океан). Также проводится численный анализ влияния наличия тупоугольных треугольников в сетке. Предложен метод аппроксимации уравнений на тупоугольных треугольниках. Проведённые тестовые расчёты показали, что наличие даже 15% тупоугольных треугольников в сетке существенно не влияет на результат. Во второй части главы рассмотрена предыдущая задача, в которой уравнение неразрывности имеет вид: = сНу(Я + С)и- Таким образом, новая задача позволяет моделировать динамику мелкой воды, в которой отсутствует требование на глубину |£| <С Н. Проводится численное сравнение данной задачи с предыдущими задачами, также приводятся результаты численных экспериментов, в которых влияние модифицированного уравнения неразрывности играет существенную роль.

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

Результаты, изложенные в диссертации, были представлены на международных конференциях молодых ученых «Ломоносов-2009» и «Ломоносов-2011», на международной конференции «Современные проблемы математического моделирования» в Абрау-Дюрсо (Россия, 2011 г.), на международной

научной конференции «Современные проблемы вычислительной математики и математической физики», посвященной памяти академика А.А.Самарского (Москва, 2009 г.). Помимо этого, результаты работы неоднократно докладывались на научно-исследовательских семинарах механико-математического факультета МГУ имени М.В. Ломоносова, Института вычислительной математики РАН в 2009—2012 годах, Научно-исследовательского вычислительного центра МГУ имени М.В. Ломоносова, Института прикладной математики имени М.В.Келдыша.

Автор выражает глубокую благодарность своему научному руководителю Георгию Михайловичу Кобелькову, за постановку задачи, за помощь и поддержку на протяжении всей научно-исследовательской деятельности, приведшей к получению результатов, изложенных в настоящей диссертации. Автор также выражает глубокую благодарность ведущему научному сотруднику Института прикладной математики имени М.В.Келдыша Игорю Владимировичу Фрязинову за консультации и плодотворные обсуждения аппроксимаций поставленных задач. Кроме того, автор выражает благодарность ведущему научному сотруднику Института прикладной математики имени М.В.Келдыша Игорю Викторовичу Попову за консультации, обсуждения результатов и помощь в построении расчётных сеток. Помимо этого, автор выражает благодарность ведущему научному сотруднику Института вычислительной математики РАН Владимиру Борисовичу Залесному за постановку задачи приливных сил в сферической системе координат, плодотворные консультации и обсуждения результатов. Автор также выражает благодарность сотрудникам кафедры вычислительной математики механико-математического факультета МГУ имени М.В.Ломоносова, в особенности И.О.Арушаняну, О.Б. Арушаняну, К.Ю.Богачёву, В.Д.Валединскому, И.С.Григорьеву, А.А.Корневу, Е.А.Лапшину, М.А.Ольшанскому, А.В.Попову, В.М.Староверову, Е.А.Чижонкову, в прекрасном коллективе которых сформировалось его профессиональное видение проблематики современной науки и методики её преподавания в университете.

Глава 1

Основные обозначения и определения

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

1.1 Основные обозначения

В настоящей работе все уравнения будут рассматриваться в двумерной ограниченной многосвязной области П с границей сЮ = Г. Для построения разностной схемы нам потребуется триангуляция области состоящая только из остроугольных треугольников (ячеек). Это условие необходимо в силу того, что в остроугольном треугольнике центр описанной окружности находится внутри треугольника. Построение такой триангуляции будет рассмотрено в разделе 1.6.

Замечание 1.1.1. Построение аппроксимации на тупоугольных и прямоугольных треугольниках, если их количество в триангуляции невелико,