ЧИСЛЕННЫЙ АНАЛИЗ ЛИНЕАРИЗОВАННОЙ ОДНОМЕРНОЙ МОДЕЛИ ГИДРАВЛИЧЕСКОГО СКАЧКА Текст научной статьи по специальности «Физика»

УДК 532.5

DOI: 10.14529/ mmp210406

ЧИСЛЕННЫЙ АНАЛИЗ ЛИНЕАРИЗОВАННОЙ ОДНОМЕРНОЙ МОДЕЛИ ГИДРАВЛИЧЕСКОГО СКАЧКА

Е.В. Семенко1, Т.И. Семенко1

1 Новосибирский государственный технический университет, г. Новосибирск, Российская Федерация

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

Ключевые слова: гидравлический скачок; теория мелкой воды; сдвиговая теория мелкой воды; ударная волна; фронт ударной волны.

Введение

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

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

канале с запрудой на сливе

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

Работа организована следующим образом: в первой части кратко приводится необходимые теоретические сведения из [6], во второй части описываются методы численных расчетов, третья часть посвящена результатам численных расчетов и их интерпретации.

1. Постановка задачи, аналитические формулы для решения

Рассматривается одномерное (осесимметрическое) течение в канале с запрудой, искомыми являются функции Н(х, Ь) - уровень воды над дном канала, и(х, Ь) - усредненная по вертикали горизонтальная скорость, Ф(х, Ь) - энстрофия, /(Ь) - координата точки гидравлического скачка. Здесь Ь > 0 - время, х - пространственная переменная, х Е [—Ь\, Ь2], где х = -Ь\ соответствует входу в канал, х = Ь2 - точка расположения запруды.

Одномерная нелинейная модель гидравлического скачка [1-5] основана на системе дифференциальных уравнений

дк д(Ни) --Ь ——- = О,

дЬ дх

^ V 2 ) дх\ V 2 2

= овНи - ( Сг + Сг——— ) и3,

V ^ +

параметрами здесь являются следующие константы: д - ускорение свободного падения, в - тангенс угла наклона дна канала, Cf - коэффициент трения о дно, Сг - коэффициент сопротивления [2, 3], (р - пристеночная энстрофия. Величина с = л/дН + 3(</? + Ф)Н2 является аналогом скорости звука для системы (1), при и > с течение называется сверхкритическим, при и < с - докритическим.

При х = -Ь\ задаются постоянные значения всех искомых величин, соответствующие постоянному сверхкритическому решению системы (1):

Н = НХ, и = щ = у^Нъ гдеи = Р/С,>1, Ф = 0, (2)

при х = Ь2 задается довольно распространенное граничное условие перелива [3-5]:

( 0, к<й, _ , ,

иН = д2{К)={ 2 ,__Са =-+ 0.08——, (3)

\ И > с1, тг + 2 ^

где й есть высота запруды, рис. 1.

(1)

При х = f (¿) (в точке разрыва решения) ставятся стандартные условия сопряжения [2-6] (аналог условий Ренкина - Гюгонио):

[hu] = f '(t)[h],

2

f '(t)[hu],

(4)

u

3(p + Ф)^

hu I у + gh + 2

f '(t)

h

и2 + gh + 0 + Ф)/г2

2

где квадратные скобки означают скачок функции на разрыве. Положение точки скачка х = f (¿) называют фронтом ударной волны, область х < f (¿) - областью перед фронтом, а х > f (¿) - областью за фронтом. Условия (4) связывают сверхкритическое течение перед фронтом и докритическое за фронтом.

Кроме того, задаются начальные условия: предполагается, что в начальный момент времени точка разрыва расположена в нуле, т.е. f (0) = 0, а начальные значения перед фронтом - те же константы, что и на входе, см. условие (2). Начальные значения за фронтом, как показано в [6], можно не задавать, так как их влияние на решение исчезает за конечное время.

Как отмечено в [6], при таких начальных и граничных условиях решение в области перед фронтом будет постоянным сверхкритическим решением системы (1) вида (2), и остается лишь найти решение за фронтом. Переход к линейной задаче за фронтом состоит в следующем [6]. Во-первых, замена переменной х := ¿2 сводит ее к задаче в постоянной области х £ (0,Ь2). Далее вводятся безразмерные переменные, функции и параметры (с сохранением исходных обозначений): все величины, имеющие размерность длины, заменяются на отношение их к в качестве единицы скорости для величин и, с, V = берется л/дТц] единицей энстрофии становится д/Н 1, а единицей времени - \Jh\jд. Безразмерное решение задачи записывается в виде вектора ф(х,£) = (Л,(х,^) и(х, ¿) Ф(х,£))т. Безразмерные параметры

v> 1, ^о = ^hi/g > 0, Cf = в/v, Cr > 0, L = L2/hi > 0, d0 = d/h (5)

называются управляющими параметрами, именно они определяют режим течения. Наконец, задача линеаризуется на неподвижном фронте f0 (t) = 0, при этом система дифференциальных уравнений и условие при x = L линеаризуются на специальном постоянном докритическом решении

ф2 = (h2 u2 0)T = const, (6)

где U2 = Vuh,2, а значение задается из физических соображений, о которых речь пойдет ниже. В свою очередь, система граничных условий на разрыве (при x = 0) линеаризуется на векторе ф* = (h* u* Ф*), являющемся ее решением при неподвижном фронте. В результате получается линейная задача для неизвестных функций f (t) и G(x,t) = -0(x,t) — ф2. К полученной линейной задаче применяется преобразование Фурье, которое приводит к линейной системе алгебраических уравнений для функций

те ь

С(£,ш)

(2п)2

йЬ С(х,Ь)ег йх,

0 0

СгМ = ^ I С(0,¿1, б2{и) = С(Ь, ¿1, ¡(и) = ^ / 0 0 0 Для решения алгебраической задачи в [6] получены формулы:

, ц{ш){ф*-ф2)-д2т

/И = - 777 \-^оИ,

гшУ (ш)

(7)

А 1 - ег(?;'

о;) = /оМ ^ ^ , - ф2)~

з=1

2пг(£ - &(ш))

(8)

"/оМ уй и

где /0(ш) - это одностороннее (при Ь > 0) преобразование Фурье от единицы,

/ 3 \ / и2 Н2 0

и(ш) = ^ ^ б*'(ш)ьЦ- (ш) А, А = I с2/к2 и2 Н2 | , и2 = (и2 - ^2 к2 0),

з=1 0 0 и2

= (ш+гС)/м2, б.гМ = ЩШ + 2г^1 ± 70, 7оМ = \!+

гвш (и2 + 2с2)

22

9в2

и2 - с2

и2

^2

(к^ йк

Ь,=Ь,2

с2 = + 3</?о/г|, С

2Сгг^ _ 2С>3/2

^2 = (к2) - и2к2,

^0) = (А*)-1^, ^ = (к* - 1)1 0

к*

т

А*

и

* 2 * 2

4

к* 0

(и )2 + (с* )2 2к*и* (к*)3

и*((и*)2 + 3(с*)2 + дЬ*) Ь*{?> (и*)2 + (с*)2 + дЬ*) З(к*)3и*

~2 2 2 / а матрицы (ш), П2 (ш), П3 (ш) вычисляются по формулам

V

(

1

^(ш)

\

270 и2С1(ш) ~ ш 2/г27о6М 0

® (ш - и2^1 (ш) + 2гв/и2 к2^1(ш)

Н2е?

гС + ш - и2^1(ш)

1

4

^2 (ш)

^з(ш) =

( J_ \

27о

21г2ъЫш) 0

( НЦш + гС)2 \ иР3(ш) ¿Ch2(w + iC)

0 (ш — U26M + 2гв/и2 h2&(ш)

h2C2

iC + ш — U2^2 (ш)

Рз(ш)

V

з 1

u2

0(0 0 1), P3(ш) = с2ш2+(2^с2+^и2)ш+C2(u2—c2) —3^u2C.

/

Наконец, функция Y(ш) в (7), (8) вычисляется по формуле

Y (ш) = U (ш )F.

(9)

Как видим, в формулы (7), (8) для решения входит неизвестная константа Л,2 из (6), ее значение определяется из следующих соображений. Вообще говоря, функция Д (ш) вида (7) имеет в нуле полюс второго порядка: Д(ш) = а/ш2 + 0(1/ш), что в физических переменных даст Д(£) = + 0(1), и, следовательно, точка гидравлического скачка за конечное время выйдет за пределы физической области. Чтобы этого избежать, необходимо, чтобы полюс Д (ш) был не выше первого порядка, что означает и(0)(ф* - ^2) = $2, или

^e-CL/u2-3fiL/(c2i-«!) _ ^

Ф* h2 Cv2

C(c2 _ u2) + 3ви2

Q^b/icl-ul) + ^ Л 2 + _

3 u2

+

(10)

+ (ъь* - к* -

V «2 / 3

где С, и2, с2, , есть функции от Л,2, см. выше. Итак, (10) - это уравнение для определения Л,2. Найдя Л,2 из (10), далее по формулам (7), (8) можно полностью построить решение линеаризованной задачи.

Качественные свойства решения задачи, как видно из (7), зависят от наличия и местоположения корней функции У (ш). Если У (ш) не имеет корней в верхней полуплоскости и на вещественной оси, то решение задачи имеет асимптоту по времени: при £ ^ то точка разрыва (фронт) х = Д (£) стремится к постоянному значению хо, где

и'(0)(ф* - ф2)

xo

iY (0)

(11)

а решение ф(х,£) за фронтом - к стационарному решению линеаризованной задачи. Такой случай называется случаем устойчивости. Если функция У(ш) имеет вещественные корни и/или корни в верхней полуплоскости, то к «устойчивому:», решению задачи добавляются осциллирующие слагаемые. При этом, если хотя бы один корень

шк лежит в верхней полуплоскости (1т ш^ > 0), то решение ф(х,£) экспоненциально растет по времени. Этот случай называется неустойчивостью. Если же все корни вещественны (1т шк = 0), то экспоненциального роста нет, а есть асимптотические колебания фронта Д (£) по времени с частотами шк и амплитудами

А

и (шк )(ф* - - $2

ш2У '(шк )

(12)

и асимптотические колебания решения ф(х,£) по времени и по пространству [6]. Последний случай называется нейтральной устойчивостью.

Случаи устойчивости/неустойчивости/нейтральной устойчивости разделяются в зависимости от знака функции

0о(^, ^о,С7,СГ= Ао - |Л| - |А2|,

(13)

где Л = ехр ( Ау = ехр Г Щ^-щ) ^^

А

ехр

2«2С2(«2 - С2)

-ЬС

«2

С2 - «2

1 1 ~2С2 2^2

2«2^2 («2 + С2) Т

0

с2 И2 —

С2

1

1 1 ^

Т

О - ) 0(0 о 1).

С2«2 «2

«2 + С2 V 2С2 2^2 У V С2

В [6] показано, что условие д0 > 0 определяет область устойчивости, д0 < 0 -область неустойчивости, и эти области разделяет область нейтральной устойчивости -так называемая переходная поверхность

^о(V, ^0,^7,СГ,Ь,^о) = 0.

(14)

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

ео8(а1ш) = -31^(^4), еов(а2ш) = -sign(A2), где а1

2ЬС2

Со По

ЬС2

, «2

«2(С2 - «2)

. (15)

1

2. Схема численных расчетов

Формулы (7), (8), (10) позволяют выразить решение линейной задачи о гидравлическом скачке через набор управляющих параметров V, ^о,С*7, СгОпишем методы и последовательность расчетов по этим формулам. Остановимся на следующих четырех задачах:

1. Решение уравнения (10). Итак, требуется найти значение Л,2 из нелинейного уравнения (10), имеющего вид

Т (Я, Л) = 0, (16)

где Я = Л,2, Л = (V, <^о, С7, Сг, Ь, ^о) - набор заданных параметров (ниже - параметр). Для численного решения уравнения мы используем метод Ньютона. Проблема выбора начального приближения в методе Ньютона решается с применением метода

продолжения по параметру, суть которого кратко состоит в следующем. Предположим, что известно решение (Н0, Л0) уравнения (16). Далее мы меняем параметр с постоянным шагом АЛ от значения Л0 до заданного в (16) значения Л, то есть рассматриваем значения Лк = Л0 + к АЛ, к = 1, 2,..., п, где Лп = Л, и последовательно для к = 1, 2,..., п находим методом Ньютона решение Нк нелинейного уравнения

Т (Нк, Лк) = 0,

вычисляя итерационные приближения для Нк по формулам:

40) = Нк_ъ Я« = Я- , г =1,2,...

Таким образом, в качестве нулевого приближения для решения Нк, к= 1 ,п, используется вычисленное на предыдущем шаге решение Нк-1, соответствующее параметру Лк-1. Величина шага АЛ должна быть достаточно мала для того, чтобы обеспечить на каждом шаге продолжения сходимость метода Ньютона. В результате на п-ном шаге будет вычислено решение Н уравнения (16), соответствующее заданному параметру Л.

Определив значение к2 из (10), мы можем далее по формулам (7), (8) полностью найти решение линейной задачи в плоскости преобразования Фурье.

2. Вычисление точек на переходной поверхности (14). Пусть имеются два набора параметров: Л1 = (V1, С^, С^, Ь1 ,й0) и аналогичный набор Л2, причем д0(Л1) > 0, а д0(Л2) < 0. Соединив точки Л1 и Л2 отрезком, методом деления отрезка пополам находим точку Л, удовлетворяющую уравнению д0(Л) = 0 (точку на переходной поверхности) .

3. Построение профилей переходной поверхности. Мы можем построить сечения переходной поверхности (14) любыми двумерными плоскостями, то есть исследовать зависимость на переходной поверхности между какими-либо двумя управляющими параметрами Л1,Л2 при фиксированных значениях остальных четырех параметров. Это также делается методом продолжения по параметру: мы находим на переходной поверхности начальную точку, как это описано выше, затем маленькими шагами меняем параметр Л1 и методом Ньютона вычисляем соответствующее значение параметра Л2 как решение уравнения (14).

4. Вычисление частот и амплитуд колебаний в случае нейтральной устойчивости. Для вычисления шк мы используем формулы (15). Из формул следует, что отношение

а = ^ = (17)

сх\ 2 2 V У К 1

должно быть рациональным числом, поэтому мы используем следующий алгоритм: приближаем число а вида (17) с заданной степенью точности рациональным числом, т.е. полагаем а ~ п2/п1, причем четность/нечетность чисел п1>2 должна соответствовать условиям (15): пз- - четное (нечетное) ^^ Аз- < 0 (Аз- > 0), ] = 1,2. Тогда

искомый корень ш может быть вычислен по формуле ш = —-. Очевидно, число

а1

решений ш для (15) будет, вообще говоря, бесконечно, причем ш ^ то. Мы будем вычислять только первую (самую маленькую положительную) частоту ш, при этом

особый интерес для нас будут представлять не конкретные значения частот, а их зависимость от управляющих параметров. Аналогично описанному выше способу движения вдоль переходной поверхности мы строим зависимости частоты от любого из управляющих параметров при фиксированных значениях четырех из них; последний, шестой, параметр определяется из условия нейтральной устойчивости (14).

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

Как отмечалось выше, решение линейной задачи о гидравлическом скачке [6] однозначно определяется шестью безразмерными управляющими параметрами: V, <0, Cf, Сг, Ь и й0. Для перехода к размерным величинам как в решении, так и в исходных данных, необходимо задать еще Н1 - уровень воды на входе. Для удобства сопоставления наших расчетов с имеющимися в [3,4] мы будем в число исходных данных включать размерные величины: = Н1и1 - поток на входе (м2/с), Ь2 = ЬН1 -точка расположения запруды (м), й = й0Н1 - высота запруды (м). Итак, в конечном счете решение линейной задачи будет определяться семью параметрами: V, <0, С/, Сг, Ь2, ^1, й, и, в частности, уравнение переходной поверхности (14) примет вид

д0(V, <0, С/,СГ,Ь2,дьй) = 0, (18)

где д0 - это функция д0 вида (13) после указанной замены параметров.

В работе [3] была численно решена нелинейная задача для трех наборов управляющих параметров. Во всех вариантах выбирались одинаковые значения четырех параметров: <0 = 0, 005, С/ = 0, 00177, Ь2 = 10 м, = 0, 0835 м2/с, и варьировались значения V (характеризует свойства дна канала - уклон, трение), Сг - коэффициент сопротивления, й - высота плотины:

01: V =1, 69, Сг = 0, 0977, й = 0, 0103 м. Обнаружена устойчивость.

02: V = 4, Сг = 0,174, й = 0, 026 м. Обнаружены колебания фронта.

03: V = 126, 563, Сг = 1, 74, й = 0,135 м. Обнаружены колебания фронта.

На основании численных расчетов в [3] были выдвинуты две гипотезы о качественном поведении решения. Во-первых, при возрастании числа Фруда Г г = л/и (т.е. с ростом V) происходит переход от устойчивости к нейтральной устойчивости (колебаниям), причем переходное число Фруда, по мнению авторов, находится между 1, 3 и 2 (считается примерно равным 1, 5). Во-вторых, с ростом числа Фруда частота колебаний фронта убывает, а амплитуда растет.

Наличие аналитических формул для решения линейной задачи позволяет нам исследовать закономерности более детально, в частности - при любых заданных шести управляющих параметрах из (18) методом продолжения по параметру можно вычислить переходное значение последнего, седьмого, параметра. Численно были решены следующие задачи:

1. При фиксированных значениях параметров <0, С/, Ь2, , Сг, й, соответствующих вариантам 01, 02, были вычислены переходные (от устойчивости к неустойчивости) значения параметра V как решения уравнения (18). Расчеты демонстрируют следующее: переходные значения параметра V в вариантах 01, 02 не очень отличаются друг от друга, для 01 получаем V ^ 1, 475, что соответствует числу Фруда

^г ~ 1, 214, для варианта 02 имеем V « 1, 557, что соответствует ^г ~ 1, 248. Как видим, переход к неустойчивости в линейной задаче происходит при несколько меньших, чем предполагалось в [3], значениях, но при этом качественная картина сохраняется: с ростом числа Фруда при неизменных значениях остальных параметров происходит переход от устойчивости к неустойчивости.

2. Численно исследована зависимость между V и д1, а также между V и d на переходной поверхности, а именно: проделаны расчеты профилей переходной поверхности в плоскости переменных (V, д1) при разных фиксированных значениях d и, соответственно, в плоскости переменных (V, d) при разных фиксированных для следующего варианта значений остальных параметров:

Ei :

^о = 0, 005, Cf = 0, 00177, Cr = 0,15, L2 = 10 м.

Результаты изображены на рис. 2: в плоскости (V, ) области неустойчивости расположены выше соответствующих графиков, а в плоскости (V, d) - ниже графиков.

V

ё=\8мм ё=26мм й=34Амм

V

4(м)

0.10 0.08 0.06 0.04

0.02 f

V

д=0,03286м 2/c

д=0,0835м2/с

д=0Д445м2/с

Рис. 2. Переходная поверхность в плоскости (V, д1) при разных d ив плоскости (V, d) при разных , вар. Е1

2

0.08

0.06

0.04

Графики отражают влияние на характер течения параметров V, д1, d при фиксированных значениях прочих параметров в ситуации нейтральной устойчивости. Из рис. 2 мы видим, что с ростом V происходит переход из области устойчивости в область неустойчивости, причем переходное значение V увеличивается с ростом высоты запруды d и уменьшается с ростом напора

3. Проводился поиск переходного значения V при остальных заданных параметрах из варианта 03 (см. выше). Расчеты показали, что при всех V > 1 в этом случае выполнено неравенство (о < 0, т.е. неустойчивость в линейной задаче не достигается. Аналогичная ситуация возникает при расчетах с рассмотренными в [4] параметрами

Б: V = 30, 9136, ^о = 0, 005, С7 = 0, 00236, Сг = 0, 682, ¿2 = 10 м, 51 = 0, 02286 м2/с, d = 0, 0344 м;

в [4] в этом случае обнаружены колебания фронта. Отсутствие переходного значения V в линейной задаче означает, как и для случая 03, что при данных наборах параметров прямая, параллельная оси V, не пересекает переходную поверхность (18). Чтобы поймать переходную поверхность, нужно сдвинуть эту прямую, изменив значения некоторых управляющих параметров. На рис. 3 приведены профили переходной поверхности в плоскости переменных (V, 51) при разных значениях d и, соответственно, в плоскости переменных (V, d) при разных ^ для следующего варианта значений остальных четырех параметров:

Е2: ^о = 0, 005, Cf = 0, 00236, Сг = 0, 85, ¿2 = 10 м.

0.038 1 0.036 1 0.034 1 0.032 1 0.030 1 0.028 1 0.026 1

2.0 2.5 3.0 3.5

д=0,03103м2/с ^1=0,0367м2/с 41=0,0437м2/с

Рис. 3. Переходная поверхность в плоскости (V, 51) при разных d ив плоскости (V, d) при разных , вар. Е2

Как видим, для варианта Е2 сечения переходной поверхности плоскостью (V, д1) имеют минимум в точке V ^ 2, 45 (т.е. число Фруда ^У ~ 1, 57). Другими словами, если задать значение меньше некоторого минимального, то прямая, параллельная оси V, не пересекает переходную поверхность, при минимальном значении она касается переходной поверхности, а при дальнейшем увеличении пересекает ее в двух точках. То же самое мы наблюдаем и в плоскости (V, d), только здесь функция d(v) имеет максимум. В частности, это означает, что гипотеза [3] об обязательном переходе от устойчивости к неустойчивости при возрастании V не вполне справедлива: переходная поверхность устроена несколько более сложно - по крайней мере, для линейной задачи.

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

4. Проведены вычисления характеристик осцилляций в случае нейтральной устойчивости. На рис. 4 и 5 для вариантов Е1 и Е2 соответственно изображены зависимости первых частот ш от V при разных фиксированных значениях d (графики слева) и при разных фиксированных (графики справа). Недостающий седьмой параметр в обоих случаях определялся из (18). Из рис. 4 слева мы видим, что при достаточно малых

V

ё(м)

V

0.040

0.035

0.030

0.025

d=18мм ё=26мм d=34,4мм

d зависимость частоты от V не является монотонной: сначала частоты убывают с ростом V, затем - растут. Во всех остальных случаях, как следует из рис. 4, 5, наблюдается рост частот с ростом V, и это не соответствует упомянутой выше гипотезе [3] об уменьшении частот с ростом V. Что касается влияния на частоты высоты запруды d и величины входного потока 51, то при фиксированном V увеличение как d, так и 51 приводит к росту первой частоты ш.

а

а

V

V

с/=18лш с/=26лш

^=0,03286Л/2'С ^=0,083 5,м2'С г/, 0.1445лг с

Рис. 4. Зависимость первой частоты от V при разных d и при разных д1, вар. Е1

а

а

V

V

с/=18лш <3=26мм £/=34,4лш

^=0,03103Л(2'С ^=0,0367л/2/С г/, 0.0437л/"г

Рис. 5. Зависимость первой частоты от V при разных d и при разных д1, вар. Е2

5. Исследовано поведение точки гидравлического скачка (асимптотики хо фронта) в случае устойчивости. На рис. 6 и 7 для вариантов Е1 и Е2 соответственно изображена зависимость хо от V (при заданном и различных заданных d) и зависимость хо от d (заданном и различных заданных V). Графики в плоскости (V, хо) заканчивается при переходном (из области устойчивости) значении V, графики в плоскости ^,хо) начинаются при переходном (в область устойчивости) значении d. Как видно из рис. 6 слева, при малых значениях V точка гидравлического скачка смещена влево от начального положения, с увеличением V точка смещается вправо до некоторого максимального значения, зависящего от d, и затем опять влево - до тех пор,

25

20

20

10

5

5

пока V не переходит в область неустойчивости. Зависимость от высоты d плотины здесь следующая: чем выше плотина, тем меньше размах положений точки гидравлического скачка при разных V. Из рис. 6 справа можно заключить, что при малых фиксированных V (т.е. при малом уклоне канала либо при очень большом трении о дно) положение точки скачка мало зависит от d: точка остается немного смещенной вправо от начального положения. При больших фиксированных значениях V с увеличением высоты d запруды точка гидравлического скачка смещается вправо, т.е. в сторону запруды.

х0(м)

Х0(м)

д1=0,1445м7с

й=18мм й=26мм с= 34,4мм

#1=0,1445м2/с

- у=1,6 ■ у=1,65 . у=2,0 . у=2,4

Рис. 6. Зависимость асимптотики фронта от V и d в случае устойчивости, вар. Е1

Что касается варианта Е2, то здесь мы видим монотонную зависимость х0 от V и ту же качественную картину зависимости х0 от d, что и для варианта Е1, рис. 6.

х0(м)

х0(м)

д1=0,0437м2/с

<?1=0,0437м2/с

¿=18мм с1=26мм с!=34,4мм

1.14 1.16

-у=1,6

---------у=2,0

-т=2,5

-т=3,0

=- с1(м) 1.01 1 7

Рис. 7. Зависимость асимптотики фронта от V и d в случае устойчивости, вар. Е2

V

Заключение

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

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

Литература

1. Тешуков, В.М. Газодинамическая аналогия для вихревых течений со свободной границей / В.М. Тешуков // Прикладная механика и техническая физика. - 2007. - Т. 48, № 3. -C. 8-15.

2. Richard, G.L. A New Model of Roll Waves: Comparison with Brock's Experiments / G.L. Richard, S.L. Gavrilyuk // Journal of Fluid Mechanics. - 2012. -V. 698. - P. 374-405.

3. Richard, G.L. The Classical Hydraulic Jump in a Model of Shear Shallow-Water Flows / G.L. Richard, S.L. Gavrilyuk // Journal of Fluid Mechanics. - 2013. -V. 725. - P. 492-521.

4. Delis, A.L. Numerical Simulation of Hydraulic Jumps with Shear Shallow Water Model / A.L. Delis, H. Guillard, Tai Yih-Chin // Journal of Computational Mathematics. - 2018. -V. 4. - P. 319-344.

5. Ivanova, K.A. Structure of the Hydraulic Jump in Convergent Radial Flows / K.A. Ivanova, S.L. Gavrilyuk // Journal of Fluid Mechanics. - 2019. - V. 860. - P. 441-464.

6. Semenko, E.V. The Analytic Description of Hydraulic Jump in the Linear Theory of the Shear Shallow-Water Flows / E.V. Semenko // Physics of Fluids. - 2019. - V. 31. - Article ID: 016101.

Евгений Вениаминович Семенко, доктор физико-математических наук, профессор, кафедра инженерной математики, Новосибирский государственный технический университет (г. Новосибирск, Российская Федерация), semenko54@gmail.com.

Татьяна Ивановна Семенко, кандидат физико-математических наук, доцент, кафедра алгебры и математической логики, Новосибирский государственный технический университет (г. Новосибирск, Российская Федерация), setativ3@gmail.com.

Поступила в редакцию 14 мая 2021 г.

MSC 76L05 DOI: 10.14529/mmp210406

NUMERICAL ANALYSIS OF A LINEARIZED ONE-DIMENSIONAL HYDRAULIC JUMP MODEL

E.V. Semenko1, T.I. Semenko1

Novosibirsk State Technical University, Novosibirsk, Russian Federation E-mails: semenko54@gmail.com, setativ3@gmail.com

In this paper, a one-dimensional linear hydraulic jump problem is numerically analyzed. The stability of the solution is investigated, the investigation is based on the known formulas

for the linear problem's solution, which express it in the terms of the model parameters. A system of calculations, that allows to separate the regions of stability and instability and calculate their common boundary - the surface of neutral stability, is constructed. For the cases of stability and neutral stability, the principal asymptotic characteristics of the hydraulic jump are calculated: the front asymptotics and, accordingly, the first oscillation frequency. The results are compared with the available numerical results for the solution of the general (nonlinear) problem.

Keywords: hydraulic jump; shallow water theory; shear shallow water theory; shock wave; shock wave front.

References

1. Teshukov V.M. Gas-Dynamic Analogy for Vortex Flows with a Free Boundary. Applied Mechanics and Technical Physics, 2007, vol. 48, no. 3, pp. 8-15. (in Russian)

2. Richard G.L., Gavrilyuk S.L. A New Model of Roll Waves: Comparison with Brock's Experiments. Journal of Fluid Mechanics, 2012, vol. 698, pp. 374-405.

3. Richard G.L., Gavrilyuk S.L. The Classical Hydraulic Jump in a Model of Shear Shallow-Water Flows. Journal of Fluid Mechanics, 2013, vol. 725, pp. 492-521.

4. Delis A.L., Guillard H., Tai Yih-Chin. Numerical Simulation of Hydraulic Jumps with Shear Shallow Water Model. Journal of Computational Mathematics, 2018, vol. 4, pp. 319-344.

5. Ivanova K.A., Gavrilyuk S.L. Structure of the Hydraulic Jump in Convergent Radial Flows. Journal of Fluid Mechanics, 2019, vol. 860, pp. 441-464.

6. Semenko E.V. The Analytic Description of Hydraulic Jump in the Linear Theory of the Shear Shallow-Water Flows. Physics of Fluids, 2019, vol. 31, article ID: 016101. DOI: 10.1063/1.5072772

Received May 14, 2021