Вычислительные технологии
Том 9, № 6, 2004
ОБ ОДНОЙ СИСТЕМЕ НЕЛИНЕЙНО-ДИСПЕРСИОННЫХ УРАВНЕНИЙ ГИДРОДИНАМИКИ С ПОЛЕЗНЫМ СВОЙСТВОМ*
З. И. Федотова
Институт вычислительных технологий СО РАН, Новосибирск, Россия
e-mail: [email protected]
The paper is devoted to derivation and investigation of nonlinear dispersive equations of hydrodynamics, which admit the analogues of Riemann invariants in the case of a flat bottom.
Введение
В работе рассмотрена система уравнений из класса приближенных слабо нелинейных моделей гидродинамики, описывающих распространение длинных поверхностных волн на воде с учетом дисперсионных эффектов. Такие модели успешно применяются для решения задач гидродинамики длинных волн, когда применение полных трехмерных моделей сталкивается с необоснованно большими затратами вычислительных ресурсов.
В научной литературе о волнах на воде описано несколько разновидностей нелинейно-дисперсионных (НЛД) моделей. Большинство из них выведено из "базовой" модели, приведенной к нормированным переменным и содержащей два параметра: параметр нелинейности а и параметр дисперсии ß (см., например, [1, 2]). В работах [1, 3, 4] выполнен сравнительный анализ и дана классификация этих моделей по свойствам, важным для моделирования задач волновой гидродинамики. Одной из причин разнообразия модельных уравнений является возможность выбора зависимых переменных. Например, под "скоростью" в НЛД-уравнениях можно понимать не только скорость воды на дне канала (как в "базовой" модели), но и скорость на свободной поверхности, или усредненную по глубине скорость, или скорость на некоторой высоте Y над основанием канала (что использовано в данной работе) и т.д. В результате можно добиться тех или иных свойств от системы НЛД-уравнений. Ранее [4, 5] мы обратили внимание на полезное свойство одной системы НЛД-уравнений, состоящее в возможности привести ее к переменным, аналогичным инвариантам Римана [6]. В работах [5, 7] изучена соответствующая одномерная модель и продемонстрирован ряд ее свойств. Ниже мы покажем ряд преимуществ этой модели перед другими из класса НЛД-моделей и приведем ее вывод для произвольного дна в двумерном случае.
* Работа выполнена при поддержке Российского фонда фундаментальных исследований (грант № 0305-64108), Программы интеграционных исследований СО РАН (грант № 2003-5) и Программы поддержки ведущих научных школ РФ (грант № НШ-2314.2003.1).
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2004.
1. Вывод модели
Для вывода НЛД-модели с указанным свойством рассмотрим уравнения безвихревого течения несжимаемой жидкости со свободной поверхностью и жестким непроницаемым дном. Чтобы учесть длинноволновой характер течения и проследить за нелинейностью, используем, как обычно, разномасштабные характерные размеры в горизонтальном и вертикальном напрвлениях: I — длина, Но — глубина, по — амплитуда возмущения поверхности воды и параметры в = (Но/1)2 (параметр дисперсии) и а = По/Но (параметр нелинейности). Далее нормируем переменные с помощью следующей замены:
х'г = 1х (г = 1, 2), у' = Ноу, £ = -, у' = ^у, Н' = НоН, п' = поП-
Со Со
Здесь х^х2 — пространственные координаты на горизонтальной плоскости; у — вертикальная координата; £ — время; п — возвышение свободной поверхности над ее равновесным положением; Н — глубина невозмущенной жидкости; у — потенциал скорости жидкости; со = (дНо)1/2, д — ускорение свободного падения. Штрихами помечены размерные переменные. Соответствующие уравнения гидродинамики в указанных безразмерных переменных можно записать в виде
А у + 1 Ууу = 0, -Н (XI ,Х2) < у < ап(х1,Х2,Ь); (1)
УуУН + 1 уу = 0, у = -Н (Х1,Х2); (2)
П + аУуУп - 1 Уу = 0, у = ап(х1,Х2,£); (3)
а /-г—, \2 а 2 / \ /л\
уг + 2(уУ) +2вУу + у = 0, У = ап(х1, х2, - (4)
^ ( д д \ А „2 д2 д2 Здесь V = ——, —— и А = V2 = ——2 + — определенные в горизонтальной плоскости
\дх1 дх2 / дх1 дх2
операторный вектор-градиент и оператор Лапласа соответственно.
Для получения НЛД-модели примем разложение потенциала по вертикальной коорди-
^ _
нате [1, 8]: у = ф(х1,х2,£) + ^ фк(х1,х2,£)Ук, где ф, фк (к = 1, то) — функции, зависящие
к=1
только от х1,х2 и У = ау + Н(х1,х2). Подстановка суммы в уравнение (1) и граничные условия (2)-(4) приводят к рекуррентным соотношениям между этими функциями, что в результате позволяет перейти от бесконечного числа уравнений к системе двух уравнений бесконечного порядка для функций ф и п.
В целом указанная процедура весьма громоздка, поэтому мы приведем окончательные уравнения, в которых члены порядка О(в2) отброшены:
Пг + V (^ф) =
= в (V (й ^^ф) VH) + 1V (й2 (V (VHVф) + АфVH)) + 1V ; (5)
1 2 1 /1 2 фг + 2а|Vф|2 + - (й - Н) = в -а(VHVф)2 + й ^Н ^ф)4 + аVфV ^^ф)) + 2 а 2
+у ((Дф)4 + аVфV (Дф) - а (Дф)2) ). (6)
Эти уравнения называют "базовыми" НЛД-уравнениями мелкой воды (это уравнения мелкой воды второго приближения). По аналогии с классическими уравнениями мелкой воды первого приближения первое из этих уравнений будем называть уравнением неразрывности, второе — уравнением движения.
Уравнения (5), (6) являются основой для получения многих других описанных в литературе нелинейных уравнений с дисперсией. В частности, они содержат подкласс моделей Буссинеска, у которых параметры нелинейности и дисперсии имеют одинаковый порядок малости. В этом специальном случае исходная модель упрощается ввиду пренебрежения членами порядка О(ав, а2) и принимает вид
п + V (^ф) = 3вV3 (Я^ф) ; (7)
ф + 2а^ф|2 + а (Ь - я) = 2вV (я^ф)4. (8)
2 а 2 4 ъ
Здесь Ь = ап + Я — полная глубина жидкости. Преобразуя эти уравнения, удается получить ряд НЛД-систем с улучшенными свойствами [1, 3, 4]. Например, применив к уравнению (6) оператор V и заменив в системе (5)-(6) производные фХ1 и фХ2 переменными и>2, которые можно трактовать как компоненты "скорости" (на самом деле, это и есть скорость течения на дне канала), получим отличную от (5)-(6) систему НЛД-уравнений. В свою очередь, скорость w = Vф может быть заменена некоторой другой "скоростью", в частности скоростью на высоте У над основанием канала. Если рассмотреть в качестве "скорости течения" в приближенной модели истинную скорость течения и на глубине У = лД/ЗЯ:
и = Vф - вYV (ЯVЯVф) - в у V (Я2Дф) + О (в2) ,
где 72 = 2/3, то после надлежащих преобразований членов порядка О (в) и приведения полученных уравнений по mod (в2) получим уравнения
|п + V{Ы + 3 ДИ) = 61
д в д ( )
— и + а(и ■ V)u + Vn = 6—VФ[u],
где ( )
Я[и] = -кЯ2 ((и ■ V)VЯ + (VЯ х V)u±) - ЯVЯ(VЯ ■ и),
Ф[п] = V (Я2п) - (к + (Я2) п,
Ф[и] = V (Я2и) + ^ (Я2) ■ и, к = 2 - л/б.
В случае ровного дна Я(ж1,ж2) = Я0 модель (9) принимает компактный вид:
^п + V (Ьи) = вЯо2дДп,
6 (ю)
д в д ( —и + а(и ■ V)u + Vn = 6Яо2—V ^и).
Если воспользоваться предположением о потенциальности исходного течения, то на основе равенства их2 — vХl = О (в2) (и, V — компоненты скорости) в уравнениях движения можно сделать замену правой части на выражение в—— что делает уравнения более
6 от
привлекательными с точки зрения конструирования вычислительных алгоритмов. Стандартная процедура линейного дисперсионного анализа исследуемой модели приводит к дисперсионному соотношению и выражениям для мод (к = (к\, к2) - волновой вектор):
2 —0|к|2 —01/2|к| Ш = ^-„„О.. . , & =
(1+ в—02|к|/6)2' 1 + в—02|к|/6"
Таким образом, система уравнений (10) имеет две моды, которым соответствуют уравнения конечного порядка со смешанной производной третьего порядка. Нетрудно видеть, что фазовая и групповая скорости при |к| ^ то стремятся к нулю, откуда следует несущественность высокочастотных гармоник (свойство, полезное при численной реализации алгоритмов). Другие свойства будут рассмотрены ниже.
2. Численные эксперименты
Как следует из вида дисперсионного соотношения, полученные уравнения замечательны тем, что в линейном одномерном случае возможна замена переменных, при которой система уравнений распадается на два независимых уравнения, каждое из которых содержит только одну независимую функцию. В определенном смысле эту замену можно трактовать как переход к инвариантам Римана [6]. Из общего вида уравнений в случае ровного дна следует, что если такую систему уравнений записать в матричном виде, то матрица с дисперсионными членами будет диагональной. Это позволяет в полной мере использовать методику исследования целого класса разностных схем, разработанную для гиперболической системы уравнений мелкой воды, технически основанную на возможности изучения одного скалярного уравнения.
В одномерном случае уравнения (9) имеют вид
щ + ((к + в/Щи)х = в/6(Ф[п])х*,
иг + аиих + Пх = в/6 (Ф[и])хг,
Ф[п] = (—2п)х — (к + 1) (—2)х П,
Ф[и]= (—2и)х + к(—2(х u,
Е[и] = Яи, Я = —— (к——хх + (—х)2) , к = 2 — у/6.
Отметим, что Я зависит только от формы дна, т.е. это известная функция. Правая часть уравнения неразрывности не зависит от и, а правая часть уравнения движения не зависит от п. Чтобы введенная скорость имела физический смысл, следует потребовать выполнения неравенства л/2/3— < к. По существу, это есть ограничение на величину амплитуды волн понижения. Напомним, что модель (9) получена при условии, что а и в имеют один и тот же порядок малости, т. е. уже при выводе модели заложено ограничение на размах амплитуды. Указанная модель рассмотрена в [4], где, в частности, приведена вариационная формулировка и выписаны законы сохранения импульса и энергии.
Рассмотрим некоторые свойства построенной модели в сравнении с аналогичными моделями. Для построения численного алгоритма вводится замена, такая, что исходную систему уравнений можно представить в виде
п - 6(ф[п])х = 0, и - §№])* = и,
©4 + й + рщ и
0, и + Пх + 1 а(и2)х = 0.
(11)
Пусть Д* и Дж — шаги сетки в направлении осей по времени * и пространству соответственно. Введем ряд обозначений. Пусть Е — тождественный оператор, Т — оператор сдвига: Т±30(ж,*) = 0(ж ± Дж,*); Т±га0(ж,*) = 0(ж, * ± пД*). Используя комбинации этих операторов, введем следующие разностные операторы: Д1 = Т — Е, Д-1 = Е — Т_1, Д1 = Т1 — Е, ^ = (Т1 + Т_1)/2. Сеточную функцию с п-го слоя по времени будем обозначать ^^ = Дж). Уравнения системы (11) аппроксимируем следующим образом:
пГ1 — вДх2 [Д_1(н1пГ1) — (к + 1)^пП+1Д_1Я|] = ©Г1,
и
га+1
— вД2 [Д_1(Я|и"+1) + кК+1Д_1Я|] = иГ+1
Д1©" , Д1 + Д_1
д*
+
2Дж Д1^
(«(пГ1 + п" )/2 + я,- + в/3Д,- )ип
Д*
+
Д2ДД_ 3 + а</2] =0,
—--
я,
3 Дж2
[АН,Д1Д_1 я, + ((Д1 + Д_1)Яз)2/4 .
Эта схема подробно исследована в работе автора [5]. Она имеет второй порядок аппроксимации и условно устойчива, однако ограничение на шаг по времени менее строго, чем для соответствующей разностной схемы без дисперсионных членов.
Как правило, тестирование НЛД-моделей начинается с проверки, насколько точно они передают движение по ровному дну Я (ж) = Но уединенной волны
п = по 8ееЬ2[(3по)1/2/(2Яо(Яо + по)1/2)(ж — жо — С*)], С = [д(Яо + по)]1/2. (12)
Предложенная модель, как и все модели типа Буссинеска, в отличие от полных НЛД-моделей [8], не допускает точного решения (12). Если взять начальное возмущение жидкости в виде (12), то в численном решении при * > 0 происходит сброс осциллирующего хвоста и формируется солитонообразная волна, соответствующая конкретной модели.
Существенным преимуществом предложенной модели является возможность простой конструкции граничных условий, обеспечивающих выход волны за границу области без отражения, — свойство, присущее гиперболическим уравнениям.
Рассмотрим задачу: при * = 0 на отрезке [0, 30] задается возвышение свободной поверхности в виде (12) для по = 0.2, Но = 1, жо = 15, д =1 с нулевой скоростью. Границы полагаются открытыми. В точках ж = 0 (левая граница) и ж = 5 установлены виртуальные мареографы, фиксирующие высоту проходящей волны. При * > 0 волна распадается на две, идущие в противоположных направлениях. Вершина волны достигает границ за
х
0
Рис. 1. Выход волны из области (x = 0).
Рис. 2. Поведение вблизи границы (x = 5).
t =
15. На границах заданы "мягкие условия" простейшего вида: rqon = nin, Uon
П\
Пм-in, uNn
uN-1n; 0,... , N — номера узлов сетки.
На рис. 1, 2 показаны мареограммы решений, рассчитанных по однотипным разностным схемам для ряда моделей: "треугольниками" отмечено решение гиперболических уравнений мелкой воды (по описанной выше разностной схеме с в = 0), кривые с "кружками" и "ромбиками" получены соответственно по предложенной модели и известной модели Перегрина [8, 9]. На рис. 1 показана форма выходящей за пределы области волны. Видно, что в последнем случае возникает существенная по амплитуде отраженная волна, тогда как для первых двух моделей граница оказывается "прозрачной". На рис. 2 видно, что до взаимодействия с границей решения, полученные по НЛД-моделям, близки.
Мы рассмотрели и другие известные НЛД-модели типа моделей Буссинеска. Все они показали примерно тот же результат, что и модель Перегрина. Используя простейшие граничные условия (выпуск волны за пределы области и пересчет линии уреза с помощью экстраполяции), мы попытались применить НЛД-модели для расчета наката волны на берег. Работающий алгоритм удалось получить только в случае описанной здесь модели именно благодаря свойству правильно выпускать волну за пределы области.
Список литературы
[1] Уизем Дж. Линейные и нелинейные волны. М.: Мир, 1977.
[2] Fedotova Z.I., Karepova Е.Д. Variational principle for approximate models of wave hydrodynamics // Russ. J. Numer. Anal. Math. Modelling. 1996. Vol. 11, N 3. P. 183-204.
[3] Марчук А. Г., Чубаров Л.Б., Шокин Ю.И. Численное моделирование волн цунами. Новосибирск: Наука, 1983.
[4] Fedotova Z.I., Pashkova V.Yu. On the numerical modelling of the dynamics of weakly nonlinear waves with dispersion // Russ. J. Numer. Anal. Math. Modelling. 1995. Vol. 10, N 5. P. 407-424.
[5] Fedotova Z.I., Pashkova V.Yu. Methods of construction and the analysis of difference schemes for nonlinear dispersive models of wave hydrodynamics // Russ. J. Numer. Anal. Math. Modelling. 1997. Vol. 12, N 2. P. 127-149.
[6] Рождественский Б.Л., Яненко Н.Н. Системы квазилинейных уравнений и их приложения к газовой динамике. Второе доп. издание. М.: Наука, 1968.
[7] Chubaroy L.B., Fedotoya Z.I., Shokin Yu.I., Einarsson B.G. Comparative analysis of nonlinear dispersive models of shallow water // Intern. J. of Computational Fluid Dynamics. 2000. Vol. 14, N 1. P. 55-73.
[8] ВольцингЕР Н.Е., Клеванный К.А., Пелиновский Е.Н. Длинноволновая динамика прибрежной зоны. Л.: Гидрометеоиздат, 1989.
[9] Peregrine D. Long Waves on a Beach // J. Fluid Mech. 1967. Vol. 27. P. 815-827.
Поступила в редакцию 10 сентября 2004 г.