CC BY
41
Вычислительные технологии
Том 17, № 5, 2012
Математическое моделирование движения свободной границы волны, обрушающейся под воздействием
u >к
ветровой нагрузки*
С. П. Бдутин
Уральский госудаpственный университет путей сообщения, Екатеринбург, Россия
e-mail: SBautin@math .usurt.ru
Проведено приближённое математическое моделирование обрушения волны под воздействием ветровой нагрузки. При этом рассматривается только движение свободной границы волны без учёта движения среды внутри волны. Для моделирования используются строго установленные математические факты о движении свободной границы газ — вакуум. В данном случае для описания движения поверхности обрушающейся волны необходимо строить решения соответствующей системы дифференциальных уравнений с частными производными, которые выписываются в явном виде. Приведены два факта, обязательно сопровождающие процесс обрушения волны: 1 — возникновение "языка" волны, когда свободная граница становится неоднозначной функцией от горизонтальной переменной, 2 — заострение конца "языка" волны. Моделирование проведено до момента полного обрушения волны на подстилающий горизонтальный слой воды.
Ключевые слова: система уравнений газовой динамики, свободная граница, разрыв сплошности, обрушение волны.
Введение
В работе рассматривается случай двух пространственных переменных х, г, когда в направлении убывания последней действует сила тяжести с постоянным значением ускорения свободного падения д.
Численное моделирование течений жидкости с поверхностными волнами является актуальной проблемой [1], имеющей ряд приложений, в том числе по смягчению последствий катастрофического воздействия волн цунами. В настоящее время для аналитического и численного исследования соответствующих задач часто используется модель мелкой воды. В этом случае внешним давлением воздуха пренебрегается, а движение свободной границы описывается с учётом законов сохранения и действия силы тяжести вдоль оси Ог. Одним из основных элементов модели мелкой воды является постоянство скорости среды при фиксированном значении х, т. е. независимость от г. Поэтому при возникновении у волны вертикального участка свободной границы в данном случае возможны только две ситуации: 1 — движение гидродинамического скачка (так называемый бор или бора) до момента выхода этого профиля на берег [2], 2 — распад соответствующего разрыва в момент выхода волны на берег, причём у течения, обра-
* Работа выполнена при финансовой поддержке РФФИ (грант № 11-01-00198).
зовавшегося в результате распада разрыва, свободная граница задается однозначной функцией от х [3].
В связи с указанными обстоятельствами в модели мелкой воды невозможно описать наблюдающееся в природе явление обрушения волн, ещё не дошедших до берега (рис. 1 [4]). Как видно из рисунка и других многочисленных наблюдений обрушения волны, этот процесс сопровождается двумя обязательными элементами: 1 — образование "языка" волны, когда свободная граница становится неоднозначной функцией от горизонтальной переменной, 2 — заострение "языка" волны и, как следствие, возникновение пенного потока на этом участке обрушающейся волны.
Проблеме движения волн на поверхности воды посвящено достаточно много работ [4-20]. Подробную библиографию можно найти в [4, 20]. Условно все эти работы можно разделить на две группы. В первой (см. [4]) исследуются различные методы описания процесса обрушения волн. Однако причины обрушения волн ни в монографии [4], ни в работах, указанных в ее библиографическом списке, не раскрыты. В работе [20], относящейся ко второй группе исследований, производится математическое моделиро-ваниие движения волн по поверхности воды (в том числе собраны известные факты о поведении волн, включая факты о волнах-убийцах), приведены основные теоретические модели возникновения волн: дисперсионное сжатие волновых пакетов, геометрическое фокусирование, нелинейная модуляционная неустойчивость. Системами дифференциальных уравнений и отдельными дифференциальными уравнениями, решения которых используются в [5-20], являются кинематическое уравнение для частоты волнового вектора, нелинейные неодномерные уравнения Шрёдингера и Кортевега де Вриза (уравнение Кадомцева — Петвиашвили), нелинейные уравнения мелкой воды для случая наклонного берега, а также системы уравнений для потенциальных течений несжимаемой жидкости. Очевидно, что в рамках большинства перечисленных моделей невозможно описать процесс непосредственного обрушения волны и, следовательно, указать причины, его вызывающие.
Среди большого числа работ, в которых обсуждается необходимость описания процесса обрушения волн, имеется единственная [7], где предложена методика реального описания начала процесса обрушения волны для модели несжимаемой жидкости. В основе этой методики лежит переход к новым независимым переменным, фактически позволивший задать свободную границу волны в параметрическом виде. В новом функциональном пространстве с использованием паде-аппроксимаций численно решается задача Дирихле. При восстановлении свободной границы волны в физическом пространстве данная методика описывает начало образования "языка" волны, когда
указанная граница начинает задаваться неоднозначной функцией от горизонтальной переменной. Однако данная методика всегда даёт закруглённый конец "языка" волны (т. е. по результатам расчётов заострение "языка" не наблюдается) и непосредственно процесс обрушения волны не моделирует. К тому же используемый в явном виде закон изменения давления на поверхности волны, благодаря которому и образуется "язык" волны, не обоснован никакими физическими причинами.
1. Математическая модель
Приведённый библиографический обзор указывает на то, что для описания процесса обрушения поверхностных волн необходимо привлекать более сложные модели. В настоящей работе в качестве математической модели используется система уравнений газовой динамики (СУГД) для идеального политропного газа в условиях действия внешних сил [21]
pt + V ■ Vp + p divV = 0,
p [Vt + (V ■ V)V] + Vp = pf, (1)
St + V ■ VS = 0,
являющаяся дифференциальной формой законов сохранения массы, импульса и энергии, где
p =1 A2(S )pY. (2)
Y
Здесь p — плотность газа; V — вектор скорости газа; p — давление; S — энтропия; A 2 (S) — известная энтропийная функция; pf — вектор внешних сил; y = const > 1 — показатель политропы (адиабаты) газа; t — время; V и div — операторы градиента и дивергенции по декартовым переменным x, y, z; точка — скалярное произведение.
Для описания движения воды можно использовать [22, 23] уравнение состояния (2) политропной среды с большим значением показателя политропы: y = 7.02. Естественно, в этом случае требуется достаточно большой перепад давления для незначительного изменения значений плотности, больших единицы. Поэтому традиционно считается, что модель несжимаемой среды достаточно адекватно описывает течения воды. Однако необходимо учитывать два принципиальных момента: 1 — при использовании условия p = const решения соответствующих систем уравнений с частными производными не удовлетворяют либо закону сохранения энергии, либо основному термодинамическому тождеству [24], 2 — имеются примеры [25, 26] очень важных в прикладном аспекте течений газа с незначительными изменениями плотности, которые описываются решениями СУГД (1) и которые невозможно описать в модели несжимаемой среды.
В настоящей работе математически смоделировано движение только свободной границы волны, обрушающейся под воздействием ветровой нагрузки, без описания движение среды внутри волны.
Для течений, описываемых решениями СУГД, строго установлены [27] следующие математические факты о движении свободной границы газ — вакуум. В случае гладкого примыкания газа к вакууму возникает характеристическая задача Коши (ХЗК) стандартного вида [24] с заданными значениями плотности и компонент вектора скорости газа (рассматривается случай изоэнтропических течений) на свободной границе Г0 в момент времени t = 0:
p|r0,t=o = 0, u|r0,t=o = uo (С), w|ro,t=o = Wo(C), (3)
а сама поверхность Г0 в момент времени £ = 0 задается в параметрическом виде:
Го : х = Хо(С), £ = ¿о(£),
где £ — параметр. Тогда необходимые условия разрешимости этой ХЗК имеют вид двух дифференциальных уравнений с частными производными [27]
щ + иих + ши2 = /1(х, £), и + и^х + = /з(х, г),
(4)
рассматриваемых на свободной границе Г0 в моменты времени £ > 0. В системе (4) и, и — первая и третья компоненты вектора скорости газа на Г0; /1(х, г), /3(х, г) — координаты вектора внешней силы ^ = (/1, 0, /3), действующей на рассматриваемый объём среды. От (4) с начальными условиями (3) делается стандартный переход к системе обыкновенных дифференциальных уравнений [27]
д£ (5)
= /з(х,г)
с начальными данными вида
х(0,£) = хо(£), х4(0,£)= ио(£),
г(0,£ ) = *,(£), **(0£) = щ,(£).
(6)
Из формул (5), (6) следует строго установленный математический факт о движении свободной границы газ — вакуум [27]: движение каждой отдельной точки свободной границы происходит так же, как движение отдельной материальной точки в условиях действия заданных внешних сил при задании в начальный момент времени как положения точки на свободной границе, так и скорости движения точки свободной границы. Данный закон движения свободной границы и конфигурация возникшего течения сохраняются до момента появления бесконечного градиента в какой-либо точке свободной границы, в том числе вызванного занулением хотя бы одного из радиусов главных кривизн в отдельной точке свободной поверхности [27]. Эти моменты времени определяются при решении систем транспортных уравнений, описывающих выводящие с Г0 производные [27].
Будем полагать, что движение свободной границы, разделяющей воздух и воду, также определяется сформулированными выше законами. Естественно, что данная гипотеза может использоваться только для приближённого описания движения волны, имеющей границу раздела воздух — вода.
2. Описание процесса обрушения волны
Простой анализ свойств решений задачи (5), (6) указывает на то, что если для различных точек на свободной границе начальные скорости разные, то даже при /1 = /2 = 0 движение различных точек будет осуществляться с разными скоростями и это приведёт к изменению вида свободной границы. Естественно, что помимо непосредственного описания движения свободной границы желательно обосновать причину различия начальных скоростей для различных точек на свободной границе и вид функций /1, /3,
отвечающих за внешнее воздействие на волну. Предполагаемых причин, вызывающих обрушение волны, достаточно много: ветровое воздействие, неровности дна, выход волны на мелководье, подводные течения и т. д. В данной работе в качестве причины разных скоростей различных точек на свободной границе принято ветровое воздействие и при моделировании движения свободной границы рассматривается процесс обрушения волны, вызванный воздействием ветра на одну сторону волну.
Положим, что в начальный момент времени £ = 0 свободная граница (кривая АВСД С1В1А1 на рис. 2) сначала строго повышается от точки А до самой верхней точки Д, а затем строго снижается от точки Д до точки А1. Предполагается, что ниже линии г = 0 находится горизонтальный слой воды, называемый далее подстилающим. Пусть при £ = 0 каждая точка свободной границы имеет некоторую начальную положительную скорость в горизонтальном направлении и (ж).
Для обоснования вида функций /1, /3 из задачи (5), (6) приведём следующие факторы. Пусть с момента £ = 0 на волну слева начинает действовать ветровое давление, изображённое на рис. 2 горизонтальными стрелками. Можно положить, что ветровое давление одинаково по высоте, и тогда горизонтальная составляющая силы, действующая на горизонтальные слои сплошной среды, например на слои ВВ1, СС1 (см. рис. 2), пропорциональна ширине этих слоев ДЛ. Также можно положить, что т — масса данных слоев пропорциональна их длине. Поскольку с ростом г длина слоев убывает, то т(г) в случае роста г тоже убывает. Поэтому с увеличением г растет а = а(г) — ускорение слоев, вызванное ветровой нагрузкой
^ р рДЛ ^
а(г) = —тт = —тт. (7)
т(г) т(г)
Таким образом, предполагается, что в начальный момент времени £ = 0 каждая точка свободной границы помимо начальной положительной скорости и (ж) вдоль оси Ож имеет в том же направлении и начальное ускорение а(г) > 0. А при £ > 0 в (5) можно положить
/1(ж,г) = ф).
Исходя из общего вида обрушающейся волны (рис. 3) введём в качестве гипотезы ещё одно условие движения свободной границы: для каждой точки на свободной границе, движущейся вдоль оси Ож, действие силы тяжести учитывается только в том случае, если вниз по вертикали от рассматриваемой точки имеется разрыв сплошности среды, т.е. фактически исчезает опора (см., например, точки С и Д на рис. 3). Следовательно, /3(ж,г) = —д для точек, имеющих под собой разрыв сплошности, и /3(ж, г) = 0 для точек без такого свойства.
Рис. 2. Движение свободной границы при об- Рис. 3. Общий вид обрушающейся волны рушении волны
Для численного моделирования движения свободной границы рассматриваются координаты (xi(t), Zj(t)), i = 0,1, 2, ...,M, точек, составляющих эту границу. Значения (xi(0), ^¿(0)) определяются её исходным положением в момент времени t = 0. Для каждой точки с номером i, 0 < i < M, необходимо задать начальную скорость v > 0 и ускорение aj > 0 при заданной ветровой нагрузке. В частности, можно полагать V = const > 0. В этом случае значения aj должны возрастать при увеличении номера точки от нулевого до номера самой высокой точки на начальной свободной границе. Затем с ростом номера точки значения а^ должны убывать так, чтобы для точек, имеющих одинаковую начальную высоту, значения а^ были также одинаковыми.
С учётом сформулированных предположений о виде функций /i, / и начальных условиях (5) координаты Xj(t) с течением времени изменяются по закону
xi(t) = жД0)+ Vit + 0.5ait2, 0 < i < M.
Координаты Zj(t), 0 < i < M, изменяются только с тех моментов t = to, начиная с которых при каждом конкретном i значение Xi(tO) становится больше хотя бы одного значения Xj(to), i < j < M. Следовательно, закон изменения координат Zj(t) имеет следующий вид:
z (t) = i t < tO, 0< i<M
Zi(t) = \ zi(0) - 0.5g(t - tO)2, t > tj, 0 < i < M
Замечание 1. Предложенная схема движения свободной границы предполагает наличие давления на волну только с одной стороны — наветренной — и отсутствие давления с другой стороны — подветренной. При необходимости наличие ненулевого давления с подветренной стороны можно учесть при выборе значений aj, например, положив в формуле (5) в качестве давления p разность давлений с наветренной и подветренной сторон.
Замечание 2. Если при учитываемых закономерностях движения свободной границы значение Zj(t) начинает уменьшаться, то следует изменять значения aj с учётом новой длины горизонтального слоя, в который попала точка с номером i. При этом, если начиная с какого-то момента времени на горизонтале z = Zj(t) имеется разрыв плотности левее точки с номером i, то с данного момента значение aj надо полагать равным нулю, поскольку в таком случае на указанный горизонтальный слой не действует ветровая нагрузка, приложенная к волне слева. Однако для простоты вычислений значения aj с течением времени не менялись.
Замечание 3. Если в какой-то момент времени t = t* какая-либо координата Zj(t) станет равной нулю: Zj(t*) = 0, то и далее её значение в приведённых ниже расчётах полагалось равным нулю. Последнее означает, что при попадании какой-либо точки с номером i со свободной границы на свободную поверхность подстилающего горизонтального слоя воды, вдоль которого движется волна (см. на рис. 2, 3 область, расположенную при z < 0), данная точка там и остается: Zj(t) = 0 при t > t*. При этом в работе не рассматривается результат взаимодействия горизонтального слоя воды и упавшей на него части волны.
На рис. 4 приведены результаты численного моделирования движения свободной границы для следующего случая: число точек на свободной границе равно 201, начальное положение свободной границы есть часть косинусоиды, поднятой на единицу,
0.1, значения ai вычислялись в виде
а0(г/100)2, а0(2 - г/100)2,
0 < г < 100, 101 < г< 200,
полагалось а0 = 0.5д = 0.5, т.е. безразмерное значение д равнялось единице. Следовательно, при введении безразмерных переменных за масштабные значения расстояния и времени приняты соответственно 9.8 м и 1 с. На рис. 5 для сопоставления показаны рис. 1 и вид обрушающейся части свободной границы в момент £ = 1.6 согласно проведённому расчёту.
Анализ полученных результатов приводит к следующим выводам.
1. Качественно одинаковое обрушение волны наблюдается как в описанном случае, так и при нулевых значениях а^ Однако при этом обязательным условием должен быть рост значений vi при увеличении номера точки от нулевого до номера самой высокой точки на начальной свободной границе с последующим монотонным убыванием vi при росте номеров точек до последнего. Таким образом, эффект обрушения волны, ещё не дошедшей до берега, вызывается двумя причинами: а — более высокие точки на свободной границе движутся в горизонтальном направлении быстрее точек, расположенных на свободной границе ниже, б — точка на свободной границе начинает движение вниз только тогда, когда под ней исчезает опора, т. е. ниже вдоль вертикали имеет место разрыв сплошности.
Рис. 4. Результаты численного моделирования
V
Рис. 5. Сравнение данных [4] и результатов настоящей работы по обрушению волны
2. В процессе обрушения волны в верхней части её свободной поверхности всегда образуется заострённый "язык". При этом чем больше горизонтальное ветровое воздействие, тем дольше существует возникшее заострение, вплоть до падения соответствующей части свободной границы на подстилающий горизонтальный слой воды.
Исходя из закономерностей движения сжимаемой сплошной среды, моделируемого решениями СУГД в окрестности свободной границы газ — вакуум, конфигурация течения в области заострения, возникающего на свободной границе, должна измениться [27]. Это происходит в силу того, что масса движущейся среды попадает во всё уменьшающийся объём. Представляется, что наблюдаемый в природе эффект возникновения "барашков" (вспененных бурунов) на гребне волны — переход части непрерывного потока в пенный поток (капельно-пузырчатый) — является следствием подобного заострения, возникающего в верхней части свободной границы. Расчёты показали полное качественное совпадение поведения модельной свободной границы с наблюдаемым в природе эффектом появления на гребне волны "барашков", а именно: если ветровое воздействие на волну незначительно (а не велики), то заострение свободной границы имеет место только в самой верхней части волны — при малом ветре "барашки" на волнах малы, и обрушение этого небольшого заострённого "языка" волны происходит не на подстилающий горизонтальный слой воды, а на ближайшую к "языку" часть волны. Если же в расчётах ветровое воздействие (значения а^) положить большим, то область свободной границы с заострением обгоняет волну и падает на подстилающий горизонтальный слой воды — в природе подобное обрушение волн, имеющих вплоть до момента обрушения в голове волны пенную область, наблюдается при достаточно сильном ветре.
3. Учёт свойств решений СУГД в окрестности свободной границы газ — вакуум даже с достаточно простыми зависимостями fl, f3 позволяет достаточно адекватно математически смоделировать движение свободной границы в процессе обрушения поверхностной волны, вызванном ветровой нагрузкой. Подобное моделирование возможно и в трёхмерном пространстве, поскольку приведённые строго доказанные математические факты о движении свободной границы газ — вакуум имеют место в самых общих случаях [27].
Список литературы
[1] ХАкимзянов Г.С., Шокин Ю.И., Барахнин В.Б., ШокинА Н.Ю. Численное моделирование течений жидкости с поверхностными волнами. Новосибирск: Изд-во СО РАН, 2001.
[2] Шокин Ю.И., ХАкимзянов Г.С. Схема предиктор-корректор, сохраняющая гидравлический скачок // Вычисл. технологии. 2006. Т. 11, Спец. выпуск. Часть 2. С. 92-99.
[3] Баутин С.П., Дерябин С.Л., Соммер А.Ф., ХАкимзянов Г.С. Исследование решений уравнений мелкой воды в окрестности подвижной линии уреза // Там же. 2010. Т. 15, № 6. С. 19-41.
[4] Шарков Е.А. Обрушающиеся морские волны: Структура, геометрия, электродинамика. М.: Научный мир, 2009.
[5] Стокер Дж. Дж. Волны на воде. Математическая теория и приложения. М.: Иностр. литература, 1959.
[6] Banner M.Z., Phillips O.M. On the incipient breaking of small scale waves //J. Fluid Mech. 1974. Vol. 65. P. 647-657.
[7] Longuet-Higgins M.S., Gokelet E.D. The deformation of steep surface waves on water.
I. A numerical method of computation // Proc. Royal Society. London. 1976. A 350. P. 1-26.
[8] Тернер Дж. Эффекты плавучести в жидкостях. M.: Мир, 1977.
[9] Сретенский Л.Н. Теория волновых движений жидкости. М.: Наука, 1977.
[10] Longuet-Higgins M.S., Gokelet E.D. The deformation of steep surface waves on water.
II. Growth of normal-mode instabilities // Proc. Royal Society. London. 1978. A 364. P. 1-28.
[11] Филлипс О.М. Динамика верхнего слоя океана. Л.: Гидрометеоиздат, 1980.
[12] ЛАйтхил Дж. Волны в жидкостях. М.: Мир, 1981.
[13] Koga M. Direct production of droplets from breaking wind-waves its observation by a multy-colored overlapping // Tellus. 1981. Vol. 33. Р. 552-563.
[14] Forbes L.K., Schwartz L.W. Free-surface flow over a semicircular obstruction // J. Fluid Mech. 1982. Vol. 114. P. 299-314.
[15] Melville W.K. The instability and breaking of deep-water waves // Ibid. 1982. Vol. 115. P. 165-185.
[16] Nakagawa T. On characteristics of the water-particle velocity in a plunging breaking // Ibid. 1983. Vol. 126. P. 251-268.
[17] Madsen P.A., Svendsen J.A. Turbulent bores and hydraulic jumps // Ibid. 1983. Vol. 129. P. 1-25.
[18] Белоцерковский О.М., Бурынин Е.Е., Гущин В.А. и др. Поле скоростей в окрестности зоны опрокидывания // Проблемы гидромеханики в освоении океана. Часть II. А. Киев: АН УССР, 1984. С. 156-157.
[19] Гущин В.А., Коньшин В.Н. Численное исследование опрокидывания поверхностных волн // Проблемы гидромеханики в освоении океана. Часть II. A. Киев: АН УССР, 1984. С. 164-165.
[20] Куркин А.А., ПЕлиновский Е.Н. Волны-убийцы: Факты, теория и моделирование. Н. Новгород: Нижегородский гос. техн. ун-т, 2004.
[21] Овсянников Л.В. Лекции по основам газовой динамики. М.; Ижевск: Ин-т компьютерных исследований, 2003.
[22] Нигматуллин Р.И., Болотнова Р.Х. Широкодиапазонное уравнение состояния воды и пара. Метод построения // Теплофизика высоких температур. 2008. Т. 46, № 2. С. 206-218.
[23] НигмАтуллин Р.И., БолотновА Р.Х. Широкодиапазонное уравнение состояния воды и пара. Результаты расчётов // Там же. 2008. Т. 46, № 3. С. 362-373.
[24] Баутин С.П. Характеристическая задача Коши и её приложения в газовой динамике. Новосибирск: Наука, 2009.
[25] Баутин С.П. Торнадо и сила Кориолиса. Новосибирск: Наука, 2008.
[26] Баутин С.П., Обухов А.Г. Математическое моделирование разрушительных атмосферных вихрей. Новосибирск: Наука, 2012.
[27] Баутин С.П., Дерябин С.Л. Математическое моделирование истечения идеального газа в вакуум. Новосибирск: Наука, 2005.
Поступила в 'редакцию 21 февраля 2012 г., с доработки — 16 июля 2012 г.
