Численное моделирование колебаний жидкости, вызванных мгновенным наклоном основания резервуара

Вычислительные технологии

Том 3, № 1, 1998

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ КОЛЕБАНИЙ ЖИДКОСТИ, ВЫЗВАННЫХ МГНОВЕННЫМ НАКЛОНОМ ОСНОВАНИЯ РЕЗЕРВУАРА *

В. Б. Барахнин, Г. С. Хлкимзянов Институт вычислительных технологий СО РАН Новосибирск, Россия e-mail: bar@net.ict.nsc.ru, khak@adm.ict.nsc.ru

The numerical solutions of the problem on fluid oscillations caused by instantaneous inclination of the base of the reservoir are considered. Calculations are carried out on the basis of potential fluid flows model and also nonlinear and nonlinearly-dispersive models of the shallow water theory. The bounds of applicability of the used models and algorithms are defined from numerical experiments. The qualitative calculation results are given.

Исследования динамики сосудов с жидкостью, имеющей свободную поверхность, интенсивно ведутся в течение последних пятидесяти лет. К числу важнейших работ в этой области следует отнести монографии [3, 4, 10-12]. В указанных работах, а также в статьях [2, 5-7, 9] описаны аналитические и численные методы решения задач о свободных колебаниях жидкости в сосудах различной формы; при этом задачи зачастую рассматриваются в линейном приближении, т.е. в предположении малости колебаний. Нелинейные постановки задач данного класса исследуются весьма редко.

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

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

* Работа выполнена при финансовой поддержке Программы интеграционных фундаментальных исследований СО РАН, проект №43 и Российского фонда фундаментальных исследований, грант №97-01-00819.

© В. Б. Барахнин, Г. С. Хакимзянов, 1998.

или какую-либо нелинейно-дисперсионную модель (НЛД-модель), например, модель Железняка— Пелиновского [8], интересную тем, что при ее получении не учитывается предположение о малости амплитуды волн.

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

1. Математическая постановка

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

Ф** + У2Ф = 0, -К < г < п,

У ФУК + = 0, г = -К, П + УФУп - Ф* = 0, г = п, ф + 1(УФ)2 + 2ф2 + п = 0, г = п,

где Ф(хг,х2,г,£) — потенциал скорости, п(хг,х2, ¿) — возвышение поверхности над плоскостью г = 0, дно описывается функцией г = -К(хг,х2).

Двумерные (плановые) течения жидкости со свободной поверхностью в рамках модели мелкой воды и нелинейно-дисперсионной модели Железняка — Пелиновского описываются следующей системой уравнений для безразмерных переменных:

И + а1у(Я и) = 0,

и + (и  У)и + Уп = в,

где для модели мелкой воды первого приближения

Б = 0,

для модели Железняка — Пелиновского

1 /и3 И 2 \ /И Б = И У ( — Яг + — ДН - УК ( у Яг + Я2

Яг = (&уи)4 + (и  У)а1уи - (а1уи)2, Я2 = и - УК + (и  У)(и  УК),

иа(хг, х2, ¿) — декартовы компоненты вектора скорости и (а = 1, 2), И = п + К — полная глубина.

Все уравнения дополняются начальными и краевыми условиями.

2. Вычислительные алгоритмы

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

При использовании криволинейных сеток в каждый момент времени Ь подвижная область 0(Ь), занятая жидкостью, отображается на неподвижную расчетную область Q, являющуюся в зависимости от размерности задачи единичным отрезком, единичным квадратом или единичным кубом. Такое отображение осуществляется на основе взаимнооднозначного невырожденного преобразования координат. Рассматриваемые уравнения, а также начальные и краевые условия записываются в новой системе координат [1, 13-15].

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

г 100

80

60

40

20

0 10 ж 20

Рис. 1. Траектории узлов подвижной сетки.

Кратко опишем алгоритмы численного решения. Уравнения мелкой воды решались с помощью явной конечно-разностной схемы предиктор — корректор с автоматически настраиваемой аппроксимационной вязкостью [14, 15].

Для уравнений модели Железняка — Пелиновского, содержащих производные высокого порядка искомых функций, в частности, смешанные производные по времени и про-

странственным переменным, нами разработан алгоритм с расщеплением. Он заключается в том, что на каждом временном шаге сначала решается эллиптическое уравнение относительно новой зависимой переменной [1], а затем — неоднородная гиперболическая система с правой частью, зависящей только от производных по пространственным переменным. Эллиптическое уравнение аппроксимируется с помощью схемы типа "косой крест", а для решения гиперболической системы используется упомянутая ранее схема предиктор — корректор, при этом эллиптическое уравнение решается перед обоими шагами указанной схемы.

Наконец, для модели потенциальных течений на п-м шаге вначале определяются значения потенциала Фп+1 в узлах прообраза свободной поверхности, значения Фп+1 в остальных узлах сетки и новое положение свободной поверхности пп+1. Все эти величины рассчитываются на промежуточной сетке, построенной для £п+1. После завершения таких вычислений строится окончательная адаптивная сетка на момент времени 1п+1, и значения Ф и п пере-считываются с промежуточной сетки на окончательную. Цикл вычислений на временном шаге повторяется несколько раз до сходимости величин Ф и п с заданной точностью [13].

3. Результаты расчетов

Перейдем к описанию вычислительных экспериментов. В случае "двумерного" резервуара проводилось сравнение результатов расчетов по всем трем упомянутым моделям. Область течения и начальное положение жидкости показаны на рис. 2. Здесь Ь — длина дна сосуда, ф — угол наклона дна. Заметим, что для упрощения постановки граничных условий в расчетах по приближенным длинноволновым моделям стенки сосуда полагались вертикальными. При малых углах наклона это допущение влияет на результат незначительно.

Рис. 2. Двумерная область течения.

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

В табл. 1 приведены величины заплеска на правую стенку. Отметим, что в расчетах

НЛД-модель

Модель мелкой воды

Модель потенциальных течений

Ь = 1

ф = 5°

Ь = 1

ф = 15°

Ь =10

ф = 1°

Ь =10

ф = 3°

Рис. 3. Хронограммы заплесков на правую стенку.

б

Рис. 4. Профили свободной поверхности жидкости при Ь =10, ф = 3°. а — модель мелкой воды, б — НЛД-модель, в — модель потенциальных течений.

по НЛД-модели и модели потенциальных течений при достаточно большйх Ь и ф максимальный заплеск давала не первая волна (см. рис. 3).

На рис. 4 показаны полученные по использованным моделям профили свободной поверхности в три момента времени. Профиль при I = 9.3 соответствует процессу наката на правую стенку, при I = 21.6 и 29.7 — отражения соответственно от левой и правой стенки. Очевидно, НЛД-модель и модель потенциальных течений дают качественно близкие профили, а модель мелкой воды, как и следовало ожидать, — профили с крутыми передними

а

в

Таблица 1. Величина заплеска на правую стенку.

Модель ф, град

1 2 3 5 10 15

Ь = 1

НЛ 0.017 0.034 0.051 0.085 0.167 0.248

НЛД 0.017 0.035 0.053 0.089 0.182 0.280

ПТ 0.015 0.030 0.046 0.078 0.158 0.251

Ь = 5

НЛ 0.085 0.169 0.250 0.409 — —

НЛД 0.079 0.163 0.256 0.463 — —

ПТ 0.086 0.178 0.280 0.493 — —

Ь = 10

НЛ 0.169 0.331 0.487 — — —

НЛД 0.166 0.387 0.648 — — —

ПТ 0.169 0.392 0.682 — — —

Ь = 20

НЛ 0.331 — — — — —

НЛД 0.341 — — — — —

ПТ 0.329 — — — — —

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

На рис. 5 показана расчетная область для трехмерного резервуара. В этом случае число параметров задачи удваивается, так как теперь имеются длина сосуда Ь1 и ширина Ь2, а наклон плоскости дна описывается двумя углами: ф1 между плоскостью дна и координатной плоскостью ж10ж2 и ф2 между плоскостью дна и плоскостью x10z. Отметим, что точка В во всех расчетах имела аппликату -1.

фронтами.

Рис. 5. Трехмерная область течения.

Как показали численные эксперименты по модели мелкой воды, проведенные при Ь! = Ь2 = 10 и ф! = 3°, отклонение максимального заплеска в одномерном случае от заплеска в плановом случае не превышает 10 %, если боковой наклон резервуара при этом не превосходит 0.3°. Это означает, что при 89.7° < ф2 < 90.3° можно использовать одномерную, а не плановую модель мелкой воды.

Таблица 2. Величина заплеска на ребре В В'.

Модель ф2 — 90°, град

0 0.5 1 1.5 2 2.5 3

НЛ 0.475 0.549 0.609 0.653 0.676 0.666 0.480

ПТ 0.471 0.545 0.606 0.650 0.673 0.660 0.474

а

Рис. 6. Свободная поверхность жидкости, полученная по моделям мелкой воды (а) и потенциальных течений (б).

Из рис. 3 видно, что значения первых заплесков, полученные по всем моделям, близки между собой. При моделировании трехмерных течений этот факт также имеет место, что подтверждается табл. 2, в которой приведены значения первых заплесков на ребре В Б', полученные по модели мелкой воды и модели потенциальных течений при разных углах бокового наклона резервуара. Однако, как и для двумерного сосуда, значения последующих заплесков, рассчитанные по предложенным моделям, весьма отличаются друг от друга. Различия в характере течения, особенно заметные при больших видны из рис. 6, демонстрирующего поведение жидкости в резервуаре с квадратным дном (Ь! = Ь2 = 10), наклоненном симметричным образом (ф! = 3°, ф2 = 92.12°), так что точки А и С лежат на одном уровне (точка В' расположена в верхней части рисунков). При таком наклоне резервуара следует ожидать симметричного относительно плоскости В В' V V развития течения, что и подтверждают картины свободной поверхности жидкости. Заметим, что ввиду разной скорости распространения волн, описываемых в рамках рассматриваемых моделей, сравнительно схожие волновые картины возникают в разные моменты времени.

Резюмируя вышесказанное, можно сделать вывод, что трехмерная задача является

чрезвычайно интересной, хотя и весьма сложной, и требует дальнейших исследований,

характер которых может определяться потребностями практики.

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

[1] Барахнин В. Б., ХАкимзянов Г. С. Об алгоритме численного решения уравнений одной нелинейно-дисперсионной модели мелкой воды. Вычислительные технологии, 1, №3, 1996, 5-20.

[2] Богоряд И. Б. О работах НИИПММ по динамике твердых тел, взаимодействующих с жидкостью. В "Динамика упругих и твердых тел, взаимодействующих с жидкостью. Труды семинара". Изд. Томского гос. ун-та, 1972, 17-25.

[3] Богоряд И. Б. Динамика вязкой жидкости со свободной поверхностью. Изд. Томского гос. ун-та, 1980.

[4] Богоряд И. Б., Дружинин И. А. Лекции по динамике тел с полостями, содержащими жидкость. Изд. Томского гос. ун-та, 1991.

[5] Богоряд И. Б., Дружинина Г. З. О колебаниях вязкой жидкости в прямоугольном сосуде. В "Динамика упругих и твердых тел, взаимодействующих с жидкостью. Труды II семинара". Изд. Томского гос. ун-та, 1975, 32-39.

[6] Богоряд И. Б., Лаврова Н. П., Христенко Г. В. К исследованию колебаний вязкой жидкости, частично заполняющей сосуд. В "Динамика упругих и твердых тел, взаимодействующих с жидкостью. Труды III семинара". Изд. Томского гос. ун-та, 1978, 5-11.

[7] Демин А. В. Численное решение задачи о пространственных колебаниях вязкой жидкости со свободной поверхностью. В "Динамика упругих и твердых тел, взаимодействующих с жидкостью. Труды V семинара". Изд. Томского гос. ун-та, 1984, 31-37.

[8] Железняк М. И., ПЕлиновский Е.Н. Физико-математические модели наката цунами на берег. В "Накат цунами на берег". ИПФ АН СССР, Горький, 1985, 8-33.

[9] Жмурин И. П. Расчетное и экспериментальное исследование колебаний жидкости в прямоугольном наклонном баке. Труды ЦАГИ, вып. 2096, 1981, 3-31.

[10] Луковский И. А. Введение в нелинейную динамику твердого тела с полостями, содержащими жидкость. Наук. думка, Киев, 1990.

[11] МикишЕв Г. Н., Рабинович Б. И. Динамика твердого тела с полостями, частично заполненными жидкостью. Машиностроение, М., 1968.

[12] Моисеев Н.Н., Румянцев В .В. Динамика тела с полостями, содержащими жидкость. Наука, М., 1965.

[13] Шокин Ю. И., ХАкимзянов Г. С. Конечно-разностный метод расчета трехмерных потенциальных течений жидкости со свободной поверхностью. В "Вычислительные технологии". ИВТ СО РАН, Новосибирск, 1, №1, 1992, 154-176.

[14] Barakhnin V. B., Khakimzyanov G. S. On the application of adaptive grids to the numerical solution of one-dimensional problems in the shallow-water theory. Russ. J. Numer. Anal. Math. Modelling, 10, No. 5, 1995, 373-391.

[15] Barakhnin V. B., Khakimzyanov G. S. Numerical modelling of fluid flows in the framework of a two-dimensional shallow-water model with the use of adaptive grids. Ibid., 12, No. 2, 1997, 95-110.

Поступила в редакцию 6 ноября 1997 г.