CC BY
29
УДК 502/504: 627.13:519.852.6
М. А. ВОЛЫНОВ, И. В. ГУГУШВИЛИ
Государственное научное учреждение
«Всероссийский научно-исследовательский институт гидротехники и мелиорации имени А. Н. Костякова»
Н. М. ЕВСТИГНЕЕВ
Институт системного анализа РАН
ПРИМЕНЕНИЕ ЧИСЛЕННЫХ МЕТОДОВ ИНТЕГРИРОВАНИЯ ТРЕХМЕРНЫХ НЕСТАЦИОНАРНЫХ УРАВНЕНИЙ ГИДРОДИНАМИКИ ПРИ РАСЧЕТЕ РАСПРОСТРАНЕНИЯ ВОЛНЫ ПРОРЫВА
Изложен метод численного моделирования волн прорыва в трехмерной по пространству постановке. Приведены результаты апробации метода на решении тестовой задачи и представлены некоторые результаты моделирования.
Метод численного моделирования, волны прорыва, параметры потока, интенсивность размыва русла, сила динамического воздействия, взаимодействие руслового и пойменного потоков.
There is stated a method of numerical simulation of breakthrough waves in three-dimensional on space equations. There are given results of the method approbation on the test task solution and presented some results of simulation.
Method of numerical simulation, waves of breakthrough, parameters of flow, intensity of the channel washing out, force of dynamic impact, interaction of the channel and flood-plain flows.
Проблема расчетов режимов распространения прорывных волн по речной долине, возникающих при разрушении напорных фронтов гидроузлов, приобрела актуальность в связи с необходимостью проведения расчетов для составления деклараций безопасности гидротехнических сооружений. Число гидроузлов, нуждающихся в декларировании, достигает десятков тысяч, и расчеты волн прорыва приобрели массовый характер.
Существует несколько апробированных методик расчетов режимов распространения прорывных волн: программные комплексы семейства MIKE 11 и MIKE 21, Каскад, программы, разработанные специалистами Министерства чрезвычайных ситуаций Российской Федерации, Научно-исследовательского института энергетических сооружений, Всероссийского научно-исследовательского института гидротехники имени Веденеева и ряда других учреждений.
Как правило, эти методики основаны на численном моделировании волн
прорыва в одномерной по пространству постановке и позволяют с удовлетворительной точностью решать задачи в рамках декларирования безопасности гидротехнических сооружений IV класса. Речь идет об определении времени добегания, скорости волны, глубины и продолжительности затопления в заданных точках речной долины.
Определение геометрии расчетной области и интерпретация результатов численного моделирования - такие действия представляли ранее значительные трудности для исследователей. В настоящее время с появлением программного обеспечения, позволяющего создавать геометрию области решения, в частности поперечников речной долины, с использованием электронных карт местности или спутниковой информации в автоматическом режиме, а также визуализировать результаты расчетов с использованием геоинформационных технологий, эти трудности уменьшились. В тех случаях, когда
необходимы более точные расчеты, применяются методики численного моделирования волн прорыва в двумерной постановке. Однако таких методик известно значительно меньше, чем для одномерного моделирования. В ранге программного продукта известны зарубежные комплексы MIKE 21, HEC RAS, ANSYS. Отечественные методики могут использоваться только в сопровождении разработчиков.
Одномерные и двумерные методы бессильны дать ответы именно на те вопросы, которые находятся в рамках стилизации, с помощью которых они были получены из трехмерных эволюционных уравнений Навье-Стокса.
Одномерная система уравнений Сен-Венана (система мелкой воды), полученная из исходной системы путем осреднения потока по глубине и ширине не может определить его поведение, например на резком повороте русла или при резком изменении геометрии поперечного сечения.
С помощью численного моделирования в двумерной по пространству постановке нельзя рассчитать параметры потока по его глубине, определить интенсивность размыва русла, силу динамического воздействия на находящиеся на пойме гидротехнические и другие сооружения и препятствия, взаимодействие руслового и пойменного потоков. Перечисленные вопросы могут быть решены, когда для описания трехмерных потоков применяется трехмерное моделирование.
В настоящей статье изложен метод численного моделирования волн прорыва в трехмерной по пространству постановке. Приведены результаты апробации метода на решении тестовой задачи и даны некоторые результаты моделирования.
Описание численного метода. Численное моделирование волны прорыва в трехмерной по пространству постановке требует применения математической модели, основанной на трехмерных эволюционных уравнениях Навье-Стокса с учетом интенсивно изменяющейся свободной поверхности в
режиме турбулентного течения:
jqdW + ф fyids - §7wnds - Fsg = 0; (1)
an
an
g=[0; и; v; w~J; ¿л = [(0); (u0 + /ixP); (u0 + rayP); (u>0 + nzP)]T;
tn = [О; (пхт„ + путху + nxzz + /nxu + mxz);
(nx +n T +nx +m +m к
xkyx У УУ zkyz yx yz p
+ ™zx + Wzy
f.
(п Т +П Т +ПХ
\ Х1Х у гу г г
где вектор-функция скорости V - W х [0, Т] ^ И3, скалярная функция давления Р - W х [0, Т] ^ И; «х» - топологическое декартовое произведение; t - время; W - объем; 8 - площадь; п - вектор нормали |гах; пу; пг | = 7Т.
Поиск решения в произвольно ограниченной области П е £3.
где х. - направляющие в Е3 (х, у, 2 - декартова система координат); © = п^и + геуи + п^ш; g -
ускорение силы тяжести; и; ы)}=у - декартовы составляющие вектор-функции скоро-
сти; Тц = и
дч1 öv,
- тензор вязких напря-
. дх, dxi I
жений ньютОновсКой жидкости, где п -кинематическая вязкость; i = 1...3, j = 1...3; i = {x; y; 2}; j = {x; y; 2};
= VjVj -v,v, -O.ö^y^Vj -VjSjVj)- тензор
турбулентных напряжений; v. —i-я составляющая v; v. - j-я составляющая v.
Интегральная форма (1) подразумевает наличие слабо обусловленных (турбулентных, разрывных) решений в R4 [1, 2].
Свободная поверхность описывается модифицированным VOF (Volum of Fluid) методом [3]; в котором Fs = okVF - поверхностное напряжение; а - коэффициент смачивания; k - кривизна свободной поверхности; F - цветовая функция, определяющая уровень свободной поверхности.
Система уравнений (1) дополняется кинетическими и динамическими условиями для свободной поверхности (непрерывность нормальных напряжений свободной поверхности), описываемой уравнениями для уровня свободной поверхности конвективного типа и градиентом:
76
№ 5' 2009
ÔF ÔF ÔF ÔF л — + u — + v— + w— = 0; F g [0; 1];
ot ÔX oy 02
k=V
'VF^
.M.
(2)
Предполагается, что цветовая функция уровня F равна 0 для воздуха и 1 для жидкости. Этот предложенный метод отличен от метода VOF, где F > 0 для жидкости и F < 0 для воздуха [3]. Данный метод может быть назван методом установления уровня с процедурой расчета векторов напряжений на свободной поверхности (аналогичен методу VOF).
Система уравнений (1) должна быть дополнена начально-краевыми условиями, чтобы получить корректную задачу в Q. Для уравнений Навье-Стокса такими условиями являются граничные условия типа Дирихле и Неймана [1]. В физической области можно определить следующие граничные условия: на входе задаются значения всех переменных; на выходе градиент каждой из переменных по нормали к границе задается равным 0; на твердых стенках задаются условия прилипания и непротекания.
Для получения решений для турбулентных течений необходимо разрешать все масштабы движения, таким образом явно выражая тензор турбулентных напряжений m...
Это может быть сделано следующими путями:
разрешением всех масштабов течений из нестационарных уравнений Навье-Стокса (прямым численным моделированием - DNS) и, как следствие, тензор m.. равен нулю;
введением дополнительной сложной модели турбулентности через явные уравнения для m.. с осреднением по О. Рей-нольдсу, т.е. разрешая только средние значения гидродинамических функций;
рассмотрением больших структур турбулентного потока (где турбулентное течение анизотропно) через уравнения Навье-Стокса (1) и созданием простой модели турбулентности (под-сеточной модели) для небольших
(изотропных) масштабов турбулентности (динамика больших вихрей LES);
объединением двух выше описанных методов (DES).
LES-метод для турбулентных течений. Полное описание LES-метода и его применение представлено в [4]. Вкратце: метод LES пользуется осред-ненными уравнениями (1) по времени и пространству с произвольным масштабом L, связанным с размером ячейки области дискретизации.
Введем некоммутативный оператор осреднения в R4 общего вида:
v = L(v), (3)
для которого имеет место разложение v = L_1(v) с коммутативными свойствами для дифференциального оператора:
dtYl=L-1(dtyl) = dt(L~1(yl)). (4)
я д
В (4) обозначим ot = — и т. д. для всех х.
С/ь»
В общем смысле m.. является некоммутативным тензором, который появляется из конвективной части уравнений количества движения (импульса) путем подстановки (3) в (1):
щ5 =L(vjL"1(vi))-vjvi.
Численное интегрирование (2) проводится с помощью двух методов. Один из методов используется в аэродинамике, и для уравнения неразрывности в (1) записывается в следующем виде [4]:
—§PdW + j p0ds = 0.
Для классических уравнений Навье-Стокса используется метод «предиктор-корректор», где в - параметр сжимаемости, равный 80...120 [5]. Данный алгоритм более надежный и быстрый благодаря отсутствию корректирующих уравнений Пуассона для функций давления, однако алгоритм «предиктор-корректор» более точный. Таким образом, объединение этих двух методов обеспечивает быстрый и точный расчет течений для сложных условий.
Интегрирование системы (1) состоит из двух этапов: нахождение
пространственных интегралов (интегрирование дифференциального оператора в И3) и применение процедуры продвижения решения по времени (решение ОДУ в И1 с найденными интегралами в И3).
Учитывая, что течение обладает свободной поверхностью, поиск решения (1) производится только в той области, где ¥ > 0. Следовательно, цветовая функция должна быть интегрирована достаточно точно. Более детально метод интегрирования (1), а также особенности численной реализации показаны в [5, 7].
Расчет тестовой задачи о мгновенном обрушении столба жидкости. Данная тестовая задача также является классической. С ее помощью проверяется корректность моделирования свободной поверхности.
В расчетах прорывных волн часто используется предположение о мгновенном образовании прорана в теле плотины. Таким образом, данная тестовая задача наиболее близко подходит к моделированию волны прорыва.
0,437
0,874
0,229 б
0,458
0,305
0,653
Рис. 1. Задача о мгновенном обрушении столба воды: а, б - геометрия, векторы скорости в разные моменты времени соотвествен-но; в - цветовая функция
На рис. 1 показана геометрия области расчета, положение столба жидкости в начальный момент времени, момент обрушения столба жидкости, выраженный в виде векторов скоростей, цветовая функция.
Результаты численного расчета были сопоставлены с экспериментальными данными из [8]. Рис. 2 наглядно иллюстрирует хорошее совпадение результатов.
Представленный численный метод был применен при моделировании волны прорыва по пойме с расположенным
а
в
Рис. 2. Сопоставление полученных численных результатов и физического эксперимента [8]: а - высота начального положения столба жидкости по оси Z (направление действия гравитации); Ь - текущая высота в зависимости от времени; £ - время, с; 4 -эксперимент [8]; - численный эксперимент
на ней препятствием - автодорожной насыпью с мостовым строением в русле реки. Площадь поперечного сечения мостового пролета была принята существенно меньше той, которая необходима для пропуска тела волны. Задача расчета - определение динамического воздействия волны прорыва на препятствие.
На рис. 3. изображен вид с верхнего бьефа на автодорожную насыпь в момент подхода к ней фронта волны прорыва.
По эпюрам давления видно, что динамическое давление испытывают устои моста. При высоте фронта волны
до 10 м и принятой схеме мостового строения это давление достигает 7,698-105 Па. Полученный результат говорит о том, что устои моста испытывают не только гидростатическую нагрузку от разности уровней перед и за насыпью, но и существенную динамическую нагрузку от воздействия волны прорыва.
На рис. 4. изображено обтекание мостового пролета при прохождении волны прорыва.
Видно, что часть струи, обтекающая мостовой пролет над мостовым строением, падает на дно русла реки и оказывает на него динамическое воздействие, могущее вызвать существенные размывы.
■ ' v' ••Vi.
■ .-ц, _ ^ _ • ■ . ...
-!—----——"-' -= - - tv'-^ —-Jf . 1 Ч-*т_ иг Т —- JK^gr1-щ. - - i
0,75
1,5
Рис. 3. Автодорожная насыпь в момент подхода фронта волны. Вид с верхнего бьефа. Эпюры динамического давления
Рис. 4. Обтекание мостового пролета при прохождении волны прорыва. Продольный разрез. Векторы скорости движения воды
Трехмерное моделирование частного случая - обтекания автодорожного моста волной прорыва - показало, что береговые устои, мостовое строение, а также дно русла реки за мостом испытывают значительные динамические нагрузки. Такие нагрузки, как правило, не являются расчетными при проектировании и могут привести к разрушению исследуемых объектов.
Список литературы
1. Temam, R. Navier-Stokes equations -theory and numerical analysis [Text] / R. Temam. - North Holland Publishing Company. - Oxford, 1981.
2. Yanenko, N. N. On the approximation of the Navier-Stokes equations for an incompressible fluid by evolutionary-type equations [Text] / N. N.Yanenko, B. G. Kuznetsov, Sh. Smagulov // Numerical Methods in Fluid Dynamics. - Moscow, 1984. - P. 290-313.
3. Sussman, M. Level Set Approach for Computing Solutions to Incopressible Two-Phase Flow [Text] / M. Sussman, P. Smerka, S. A. Osher // J. Comput. Phys. 114. -1994. - P. 146-159.
4. Chunlei, L. A study of kinetic energy conserving scheme using finite volume collocated grid for LES of a channel flow [Text] / L. Chunlei, N. Evstigneev // Proceedings of the international conference on numerical methods in fluid dynamics. King's Colleg London, - WC2R. -2006. - P. 46-59.
5. Evstigneev, N. M. Solution of 3D nonviscous compressible gas equations on unstructured meshes using the distributed computing approach [Text] / N. M. Evstigneev // J. of Comp. Math. And Math. Physics. - 2007. - V.8. - P. 252-264.
6. Chorin A. A Numerical Method for solving Incompressible Viscous Flow Problems [Text] / A. Chorin // J. Comput. Phys. 2. - 1967. - P. 12-26.
7. Evstigneev, N. M. Integration of 3D incompressible free surface Navier-Stokes equations on unstructured tetrahedral grid
using distributed computation on TCP/IP networks [Text] / N. M. Evstigneev // Proc. Of the VII International conf. «Advances in Fluid Mechanics». - Oxford, 2008. - P. 134-208.
8. Минаков, А. В. Моделирование многофазных течений несжимаемой жидкости со свободной поверхностью [Текст] / А. В. Минаков, А. А. Гаврилов, А. А. Дектерев // Гидродинамика больших скоростей и численное моделирование : Третья международная научная летняя школа. - Кемерово, 2006.
Материал поступил в редакцию 02.11.09. Волынов Михаил Анатольевич, кандидат технических наук, доцент, заместитель директора Тел. 8 (495) 153-21-33
Гугушвили Иракли Викторович, младший научный сотрудник Тел. 8 (495) 153-21-33 E-mail: Gugushvili_I@mail.ru
Евстигнеев Николай Михайлович, доктор физико-математических наук, старший научный сотрудник ИСА РАН E-mail: EvstigneevNM@yandex.ru
80
№ 5' 2009
