ВЕСТН. САМАР. ГОС. ТЕХН. УН-ТА. СЕР. ТЕХНИЧЕСКИЕ НАУКИ. 2017. № 3 (55)
УДК 536. 2. 075
ПОЛУЧЕНИЕ ТОЧНОГО АНАЛИТИЧЕСКОГО РЕШЕНИЯ ГИПЕРБОЛИЧЕСКОГО УРАВНЕНИЯ ПРИ ГИДРАВЛИЧЕСКОМ УДАРЕ В ТРУБОПРОВОДЕ*
Э.М. Карташов1, Е.В. Стефанюк2, С.В. Колесников2, А.Н. Бранфилева2, Г.Н. Максименко2
1 Московский технологический университет 119454, г. Москва, Проспект Вернадского, д. 78
2 Самарский государственный технический университет 443100, г. Самара, ул. Молодогвардейская, 244
E-mail: totig@yandex.ru
При использовании интегрального косинус-преобразования по ограниченной области определения эллиптической координаты и интегрального преобразования Лапласа по неограниченной области изменения односторонней параболической координаты получено точное аналитическое решение гиперболического уравнения, описывающего распределение давления по длине трубопровода во времени применительно к движению реальной вязкой жидкости в условиях гидравлического удара. Показано, что процесс колебаний давления сопровождается движением гидравлической волны, на фронте которой происходит скачок от величины давления в текущей точке до давления невозмущенного потока. После достижения фронтом гидравлической волны конца трубопровода, перекрытого задвижкой, происходит смена направления течения жидкости на противоположное со скачком давления противоположного знака на фронте гидравлической волны. Колебательный процесс изменения давления в любой точке трубопровода во времени происходит с затухающей по экспоненциальной зависимости амплитудой.
Ключевые слова: вязкая жидкость, гидравлический удар, гиперболическое уравнение, интегральное косинус-преобразование, преобразование Лапласа, точное аналитическое решение, волновые колебания давления и скорости.
Краевая задача определения давления в трубопроводе для идеальной жидкости сводится к интегрированию линейного волнового уравнения, методы решения которого хорошо разработаны. Трудности решения краевых задач для реальных жидкостей связаны с их нелинейностью. Уравнения для давления и скорости в данном случае имеют вид [1]
* Исследование выполнено при финансовой поддержке РФФИ в рамках научного проекта № 16-38-00059 мол а.
Эдуард Михайлович Карташов (д. ф.-м. н.), профессор кафедры «Высшая и прикладная математика».
Екатерина Васильевна Стефанюк (д.т.н.), профессор кафедры «Теоретические основы теплотехники и гидромеханика».
Сергей Владимирович Колесников (к.т.н.), докторант кафедры «Теоретические основы теплотехники и гидромеханика».
Анастасия Николаевна Бранфилева (к.т.н.), ассистент кафедры «Теоретические основы теплотехники и гидромеханика».
Галина Николаевна Максименко, аспирант кафедры «Теоретические основы теплотехники и гидромеханика».
ф(x, t) _
dx
S3( x, t) А(3)3
dt 2d
dp(x, t)_ , аз(x, t)
= 6 P - 5 (2)
5 5х
где р - давление;
х - продольная координата; О - скорость; / — время; р — плотность;
А(3) — коэффициент гидравлического сопротивления;
1 / P P d „ с = i/J— н--- скорость звука в капельной упругой жидкости, текущей
V k 5Е
в трубе с упругими стенками (справедливо при 3/с << 1 и ((p -p0)/Е << 1)); k — модуль упругости жидкости; 5 — толщина стенки трубы; d — диаметр трубы;
Е — модуль упругости материала трубы.
Уравнения неустановившегося движения (1) и (2) определяют зависимости между значениями скорости, плотности и давления потока жидкости, осреднен-ными по его сечению. Кроме того, они содержат среднюю величину касательного напряжения на стенке трубопровода. Система уравнений для средних величин оказывается замкнутой в случае, если указана зависимость между касательным напряжением и скоростью. Для ее замыкания предполагается, что при неустановившемся движении касательное напряжение является такой же функцией скорости и свойств жидкости, как и при установившемся движении, т. е. используется гипотеза квазистационарности, в соответствии в которой считается, что закон Ньютона для касательного напряжения выполняется и при нестационарном процессе изменения скорости. Точно так же предполагается, что соответствующие соотношения сохраняются и для коэффициента гидравлического сопротивления. Строгое обоснование этого допущения отсутствует, поэтому его справедливость может быть подтверждена лишь путем сравнения с экспериментальными данными [1]. При малых изменениях давления и скорости во времени теория квазистационарности удовлетворительно подтверждается результатами экспериментов, однако при значительных их изменениях обнаруживаются существенные расхождения. Следовательно, теория квазистационарности может быть использована лишь для некоторого ограниченного круга задач, когда распределение скорости при нестационарном течении жидкости незначительно отличается от такового при стационарном течении.
Ввиду нелинейности в правой части уравнения (1) решение системы уравнений (1), (2) возможно лишь путем численного интегрирования. Способ линеаризации уравнения (1) разработан И.А. Чарным [1], который принял множитель A3 /(2d) в виде постоянной величины, равной его среднему значению по длине
трубы и времени A3 /(2d) = 2a = const. Величина 2a для ламинарного режима, учитывая формулу Пуазейля A = 64/Re = 64v/(3d) , будет 2a = 32v / d , где
Re = &d / v — число Рейнольдса; v — кинематическая вязкость.
Для турбулентного режима течения второе слагаемое в правой части уравнения (1) осредняется в некотором интервале скорости <&<&j при замене
кривой изменения функции y = /(2d) соответствующим отрезком прямой.
2 2
Отсюда формула для величины 2a будет 2a = Л( + ЗД — 2^2)/
/[3d (&! — &о)1, где принимается Х = const.
С учетом линеаризации уравнение (1) принимает вид
ф( х, t) -- = р
дх
^ + 2a&( х, t) dt
(3)
С целью упрощения получения аналитического решения уравнений (2), (3) они сводятся к одному гиперболическому уравнению (для давления или скорости) [1]. Дифференцируя уравнение (2) по переменной t, а уравнение (3) - по переменной х, и сравнивая полученные соотношения, получаем
2 д2p(х, t)_ д2p(х, t) ^ ^д&(х, t) дх2 dt2 дх
Уравнение (2) с учетом (4) будет
^2 д2p(х, t) = д2p(х, t) 2ддР(х, t) (5)
дх2 дt2 ^ '
Аналогично получается и уравнение для скорости
с2 д2&(х, t) = д2а(х, t) д&(х, t)
дх2 дt2 дt • W
В качестве примера нахождения решений уравнений (5), (6) рассмотрим задачу о движении скачка давления в трубопроводе с жидкостью, неподвижной в исходном состоянии. Предположим, что в сечении х = 0 происходит скачок давления, а сечение х = l закрыто (скорость равна нулю). Следует найти распределение давления по длине трубы во времени. Такого рода задача встречается в расчетах гидравлических регуляторов, когда в сечении х = 0 расположен источник давления, а в сечении х = l находится регулятор расхода или давления, включающийся лишь после достижения в этом сечении заданного давления. Требуется найти запаздывание импульса и его величину, что зависит от длины трубы, вязкости жидкости и коэффициента трения. К тому же нахождение значения давления в сечении х = l сводится к задаче о гидравлическом ударе. Математическая постановка задачи имеет вид
c2 aVx^ = aVxIo + 2a dpM (t > 0; 0 < х <I) ; (7)
д х2 дt2 дt
p(х, 0) = p0 = const; (8)
дp(х,0)/ дt = 0; (9)
p(0, t) = pj = const; (10)
дp(l,t)/дх = 0 , (11)
где Po — начальное давление в трубе;
Pi — давление в точке x = 0 (pj > p0); l - длина трубопровода. Решение задачи (7)—(11) позволяет выполнить анализ изменения давления по длине трубопровода от времени. Скорости при найденном давлении находится путем интегрирования уравнения (2).
Методы решения задачи (7)—(11), основанные на применении метода разделения переменных [2, 3], а также метода Бернулли — Фурье [1], технически сложны. К тому же решение в [1] получено в размерных переменных, что не позволяет сделать выводы наиболее общего характера. Ниже излагается более простой способ получения точного аналитического решения, основанный на использовании интегральных методов [4, 5].
Для упрощения постановки задачи введем безразмерные переменные:
2 2 ~ P - Pi . „ c2t x c2
0 =-; Fo=—2; y = T ; Fo r =—2ГГ, (12)
Po - Pi 2al2 l 4a l2
где 0 - безразмерное давление;
Fo - число гомохронности (безразмерное время); y - безразмерная координата; For = const - безразмерный параметр. С учетом безразмерных переменных задача (7)—( 11) будет
50(У,Fo) + For 520(У-F0) = а20(*Fo) (Fo > 0; 0 < y < 1) ; (13)
5Fo r 5Fo2 5y2
0(y,0) = 1; (14)
50( y, 0) / 5Fo = 0; (15)
0(0, Fo) = 0; (16)
50(1, Fo)/ cy = 0 . (17)
Выполним замену независимой переменной y по формуле
4 = 1 - y. (18)
Относительно переменной 4 задача (13)—(17) примет вид
52QYo) = 50(4,Fo) + For 520(4-Fo) (Fo > 0; 0 < 4 < 1) ; (19) 542 5Fo r 5Fo2
0(4,0) = 1; (20)
50(4,0)/ 5Fo = 0 ; (21)
0(1, Fo) = 0 ; (22)
50(0,Fo)/54=0. (23)
Для приведения начального условия (20) к однородному обозначим: W(4 ,Fo) = 0(4, Fo) -1. Задача (19)—(23) примет вид
5W(4,Fo) = 5W(4,Fo) + Fo д2W(4,Fo) (Fo > 0; 0 < 4 < 1) ; (24)
542 5Fo r 5Fo2
W(4,0) = 0 ; (25)
^(4,0)/^ = 0 ; (26)
^(0,Ро)/Э4 = 0; (27)
W(1,Fo) = -1. (28)
Для исключения переменной 4 введем интегральное преобразование вида [4, 5]
W(n, Fo) = ]Щ4, Fo)cos (2П - 1)П 4 й4 . (29)
0
При этом
- РгЛ Пп _1\,г Пп-Л\г V (2п - 1)П ]р ^д
2^4 ^ cos(2n -1)14й4=(-1)п(2п -1)1 ¿542 2 ^ ' ' 2
Формула обращения имеет вид
то
W(4= 2^ W(n,Fo)cos(2n 1)п 4 . (31)
п = 1 2
Применяя (29)-(30) к задаче (24)-(28), находим
Fo г---+-—-+ у п W(n, Fo) = 0; (32)
dFo
W(n,0) = 0 ; (33)
й W(n,0)/dFo = 0. ^
(2п - 1)п
2
У п =
(35)
2
К задаче (32)-(35) применим преобразование Лапласа:
то
W(n,Р) = | W(n,Fo)exp(-PFo)dFo . (36)
0
Решение задачи (32)-(35) в пространстве изображений имеет вид
W(n,Р) = (-1)п<2п -'>(п/2> . (37)
^..Р2 + Р + у п ) Р
Рассмотрим случай, когда квадратный трехчлен
Fo ГР2 + Р + у2„ (38)
не имеет действительных корней, то есть дискриминант О = 1 - 4A,иFo < 0 , что
возможно при 4ynFo > 1.
Решение задачи (32)-(35) в области изображений в этом случае будет
W(n, P) =
(—1)n (2n — 1)(n /(2Fo))
P +
J_ 2Fo,
+
4y2For — 1 4Fo2
(39)
P
Переходя к оригиналу, находим
W(n,P) =
(—1)n (2n — 1)n
л/4г
Fo — 1 0
Fo
J sin
V4Y
Fo r — 1
2Fo
"T/(2Fo -) dT
(40)
Решение задачи (19)-(23) при
4у2For>1
принимает вид
0(4, Fo) = 1 +
(—1)n (2n — 1) (2n — 1)п
1 V4r
-cos-
n F°r — 1
2
Fo
4 J sin
Fo- —1
2Fo
e—T(2Fo- )dx . (41)
Рассмотрим случай, когда квадратный трехчлен (38) имеет действительные
корни
Pu = (—1+ V1 — 4уnFor /(2Fo r )
что возможно, если 4уиРог < 1.
В этом случае решение в изображениях будет
W(n,P) =
(—1)n (2n — 1)(п /(2Fo r)) (—1)n (2n — 1) / 2
(P — P1)(P — P2)P
a/1 — 4y *FO-
1
1
(P — P1)P (P — P2 )P
Переходя к оригиналу по переменной Fo, находим
W(n,P) =
(—1)n (2n —1)/2
а/1 — 4у nFo -
11
Л (
V P2
P
+
1У
1 eP1Fo — ^ eP2Fo
VP1
P
(42)
. (43)
(44)
Решение задачи (19)-(23) при 4у*Foг < 1 приводится к виду 0(4,Fo) = 1 + ]Г (—1)*(2n —1 cos(2n -1)п4 х
n=1y/1 — 4у
n Fo -
1 1
V P2
P
+
1У
V P1
-eP1Fo — 1 eP2Fo
P
(45)
На рис. 1-9 даны результаты исследований по формулам (41), (45) задачи
0 распределении давления нефти в стальном трубопроводе при следующих условиях [3]:
V = 7,65-10~6 м2/с; й = 0,1 м; х = 0,01; э0 = ом/с; э1 = 2м/с; с = 998м/с ;
р = 840кг/ м3 ; к = 1500МПа ; Е = 2-105 МПа; />0 = 10атм; />1 = 100 атм.
Расчеты выполнялись для длин трубопровода I = 1000км ; I = 100 км;
1 = 100 м . Соответственно указанным длинам были получены следующие значе-
2
T
e
T
0
v
У
1
х
ния числа Бог: Бог = 6,648-10"3; Бог = 0,1662; Бог = 1662,0277 - для ламинарного и Бо г = 224,1 - для турбулентного режима течения.
Из анализа результатов следует, что изменение давления происходит в виде движения гидравлической волны, на фронте которой наблюдается скачок давления от его величины на фронте до давления невозмущенного потока. Область за пределами фронта остается невозмущенной с начальным давлением Ро . Движение фронта возмущения по координате 4 во времени Бо подчиняется линейной
закономерности 4 = БоРог0'5 (см. рис. 1), что подтверждается исследованиями, выполненными другими авторами [6].
Величина числа Бог оказывает существенное влияние на получаемые результаты. Так, при малых его значениях скачок скорости наблюдается лишь на
_-2
начальном участке трубы. Например, при Бог = 6,648-10 скачок давления отмечается лишь в диапазоне 0 < Бо < 0,07 (см. рис. 2). Для очень малых значений Бо г скачки давления наблюдаются лишь на некотором малом участке трубы
вблизи сечения х = I и при малых значениях времени. Например, при Бог = 10_7
скачки давления практически заканчиваются при Бо «10_6 . Фронт возмущения за это время перемещается лишь на величину 4 = 0,003, что составляет 0,3 % от
всей длины трубопровода. На остальной части трубы при Бо > 6,648-10 изменение давления происходит без скачков вплоть до полного прекращения колебаний.
При числах Бог > 0,07 скачок давления происходит вплоть до того момента времени, когда фронт возмущения достигает координаты 4 = 0 (см. рис. 4-9). При дальнейшем увеличении времени наблюдается обратная волна с противоположным скачком давления. После достижения обратной волной сечения 4 = 1 вновь наблюдается прямая волна, но уже с уменьшенным давлением. Такое чередование прямых и обратных волн происходит до полного затухания процесса колебаний.
1,0
4 0,8
0,6
0,4
0,2
0 0,01 0,02 0,03 0,04 0,05 0,06 Бо 0,08
Рис. 1. Перемещение фронта гидравлической волны по координате 4 во времени (Бо г = 6,648 -103)
1,0 e 0,8
0,6
0,4
0,2
Рис. 4. Распределение давления в трубопроводе при Fo r = 0,1662
1,0 ®
0,8 0,6 0,4 0,2
- 0,2 -0,4 - 0,6 -0,8 -1,0
40,7 о т О <N Fo=1(
о 1Т1 0 \D 0
0,2
0,4
0,6
0,8
1,0
0
0
Рис. 5. Распределение давления в трубопроводе при Fo r = 0,1662 (n = 100)
Выполним анализ изменения скоростей при гидравлическом ударе. Выше было указано, что уравнения (1), (2) сводятся к одному волновому (гиперболическому) уравнению относительно давления или скорости (уравнение (6)). Решение уравнения для давления было получено выше. Найдем решение задачи об изменении скорости жидкости, движущейся в трубе длиной l, при начальной скорости
$0 = const. В момент времени t = 0 задвижка в сечении x = 1 мгновенно закрывается и, следовательно, скорость здесь будет равной нулю. Найти изменение скорости по длине трубы вплоть до полного прекращения колебаний.
1,0 ®
0,8 0,6 0,4 0,2
- 0,2
0,4
0,6
7000
5000
3000
2000
0,2
0,4
0,6
0,8 \ 1,0
Рис. 6. Распределение давления в трубопроводе при Бо = 1662 (п = 100)
1,0 ®
0,8
0,6 0,4 0,2 0 - 0,2
■ 0,4 - 0,6
■ 0,8 1,0
500
1000
1500
2000 Бс 2500
Рис. 7. Распределение давления для точки £ = 0 Бс при Бо г = 1662 (п = 100)
0
0
1,0
0 0,8
0,6 0,4 0,2 0 0,2 0,4 0,6 0,8 1,0
0 0,2 0,4 0,6 0,8 % 1,0
Рис. 8. Распределение давления в трубопроводе при Бо г = 224,1 (п = 100)
1,01—-----
0 __
0,8---------
0,6------------¡=1
0,4-------------
0,2-------------
0-------------
- 0,2-------------
- 0,4-------------
- 0,6----------и--
- 0,8----J---
- 1,0 _____
0 50 100 150 200 Го 250
Рис. 9. Распределение давления для точки ^ = 0 при Бо г = 224,1 (п = 100)
Математическая постановка в данном случае включает уравнение (6) с краевыми условиями
о Т) II о Рч
Го=15
1000
о о Т) м II
Рч
3(х,0) = Зо; (46)
д3( х, 0)/ дг = 0; (47)
<33(0, г)/ дх = 0; (48)
д3(/, г) = 0, (49)
где $0 - начальная скорость.
После определения решения задачи (6), (46)-(49) давление в сечении X = I,
где происходит его максимальное увеличение, определяется интегрированием уравнения (1), в результате которого находим
1 гд3(х,0)
p(l,t) = p(0, t) - p j
■ + 2a3o( x, t) dt
dx, (50)
o
где p(0, t) = pj = const - давление в сечении x = 0 (предполагается известным). Обозначим:
Э c 2t x c2
0Fo=^7; ^ = т; Fo, =
V 2а/2' " / г 4а2/2 ' где © - безразмерная скорость; Бо - число гомохронности; £ - безразмерная координата; Бог - безразмерный параметр.
С учетом принятых обозначений задача (6), (46) - (49) будет
-= ^ + *ог -Г 2 ; (ро > 0; 0 < С < 1) (51)
д£ дБо дБо
0(4,0) = 1 ; (52)
д0(С, 0) / дБо = 0 ; (53)
30(0, Бо)/ дС = 0 ; (54)
0(1, Бо) = 0 . (55)
Задача (51)-(55) полностью совпадает с задачей (19)-(23). Поэтому результаты исследований, представленные на рис. 1-9, могут быть применены и к задаче (51)-(55).
ВЫВОДЫ
1. На основе интегрального косинус-преобразования и преобразования Лапласа получено точное аналитическое решение гиперболического уравнения, описывающего распределение давления реальной вязкой жидкости по длине трубопровода во времени в условиях гидравлического удара. Уравнение получено из условия квазистационарности течения жидкости, согласно которому считается, что закон Ньютона для касательного напряжения выполняется и при нестационарном процессе изменения скорости.
2. Показано, что процесс колебания давления сопровождается движением гидравлической волны, на фронте которой происходит скачкообразное изменение давления от его величины в текущей точке до давления невозмущенного потока. Процесс скачкообразных колебаний давления в любой точке трубопровода во времени происходит с затухающей по экспоненциальной зависимости амплитудой вплоть до полного прекращения колебаний.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Чарный И.А. Неустановившееся движение реальной жидкости в трубах. - М.: Недра, 1975. - 296 с.
2. Гусейнзаде М.А. О характере изменения основных параметров течения жидкости и газа в сложной трубопроводной системе. - М.: Изд-во «Нефть и газ» РГУ нефти и газа им. И.М. Губкина, 2005. - 120 с.
3. Кудинов В.А., Кудинов И.В., Ерёмин А.В. Исследование распределения скорости течения вязкой жидкости в трубопроводе при гидравлическом ударе // Инженерно-физический журнал. - 2013. - Т. 86. - № 2. - С. 387-393.
4. Карташов Э.М., Кудинов В.А. Математические модели теплопроводности и термоупругости: 2-е изд., перераб. и доп. - Самара: Самар. госуд. техн. ун-т, 2013. - 877 с.
5. Карташов Э.М. Аналитические методы в теории теплопроводности твердых тел. - М.: Высш. шк., 2001. - 550 с.
6. Баумейстер К., Хамим Т. Гиперболическое уравнение теплопроводности. Решение задачи о полубесконечном теле // Теплопередача. - 1969. - № 4. - С. 112-119.
Статья поступила в редакцию 14 июня 2017 г.
OBTAINING EXACT ANALYTICAL SOLUTIONS OF HYPERBOLIC EQUATIONS IN THE HYDRAULIC SHOCK IN THE PIPELINE
E.M. Kartashov1, E.V. Stefanyuk2, S.V. Kolesnikov2, A.N. Branfileva2, G.N. Maksimenko2
1 Moscow Technological University 78, Vernadskogo pr., Moscow, 119454, Russian Federation
2Samara State Technical University
244, Molodogvardeyskaya st., Samara, 443100, Russian Federation
Using integral cosine - transform for the restricted definition area of elliptic coordinate and the Laplace integral transform for the unlimited area changes unilateral parabolic coordinate, it was obtained the exact analytical solution of the hyperbolic equation describing the pressure distribution along the length of the pipeline in time with respect to the motion of a real viscous fluid in the hydraulic shock. It was shown that the process of pressure fluctuation is accompanied by the movement of a hydraulic wave, front in which the pressure surge of the pressure value at the current point up to the pressure of the disturbed flow. After reaching the front hydraulic wave end of the pipeline, closed valve, changing direction offluid flow with a pressure jump of opposite sign on the at the front hydraulic wave. The oscillatory process of pressure change at any point of the pipeline in time occurs with decreasing amplitude according to the exponential dependence.
Keywords: viscous, hydraulic shock, hyperbolic equation, integral cosine - transform, the Laplace integral transform, exact analytical solution, wave oscillations of pressure and velocity.
Eduard M. Kartashov (Dr. Sci. (Phys.& Math.)), Professor. Ekaterina V. Stefanyuk (Dr. Sci. (Techn.)), Professor. Sergey V. Kolesnikov (Ph.D. (Techn.)), Associate Professor. Anastasiya N. Branfileva (Ph.D. (Techn.)), Associate Professor. Galina N. Maksimenko, Postgraduate Student.