УДК 536.2 (075)
ПОЛУЧЕНИЕ ТОЧНОГО АНАЛИТИЧЕСКОГО РЕШЕНИЯ ГИПЕРБОЛИЧЕСКОГО УРАВНЕНИЯ ПРИ ГИДРАВЛИЧЕСКОМ УДАРЕ В ТРУБОПРОВОДЕ
И.В. Кудинов, А.В. Еремин, С.В. Колесников, А.Э. Кузнецова
Самарский государственный технический университет 443100, г. Самара, ул. Молодогвардейская, 244 E-mail: [email protected]
Получено точное аналитическое решение гиперболического уравнения, описывающего течение вязкой несжимаемой жидкости в простом трубопроводе в условиях гидравлического удара.
Ключевые слова: вязкая несжимаемая жидкость, гидравлический удар, гиперболическое уравнение, точное аналитическое решение, метод Фурье, ортогональные методы взвешенных невязок.
Фундаментальный вклад в решение проблемы исследования нестационарного течения идеальной упругой жидкости внесен Н.Е. Жуковским [1, 2]. В случае идеальной жидкости решение краевой задачи о распределении давления в трубопроводе по пространственной переменной и времени сводится к интегрированию классического линейного гиперболического уравнения, методы решения которого в настоящее время хорошо разработаны. Трудности решения краевых задач о течении реальных вязких жидкостей связаны с их нелинейностью ввиду зависимости коэффициента гидравлического сопротивления от скорости. Уравнения, описывающие распределение давления и скорости для несжимаемой жидкости, в данном случае имеют вид
е2 82г{х,t) = а2г{х,р + 2аФ(м). (1)
8 л2 8г2 дt '
с 2 8 23(*, t) =8^^ + 2од®М (2)
8 *2 8t2 дt ' В качестве конкретного примера получения аналитического решения уравнений (1), (2) рассмотрим краевую задачу о распространении скачка давления в трубопроводе с неподвижной в исходном состоянии жидкостью. Допустим, что в сечении * = 0 произошло скачкообразное изменение давления, а сечение * = I перекрыто (скорость равна нулю). Требуется найти распространение давления по длине трубопровода во времени [1]. На практике подобная задача встречается при расчете гидравлических регуляторов или передач, когда в сечении * = 0 находится источник давления, а в сечении * = I к трубопроводу присоединен какой-либо прибор (например регулятор расхода или давления), включающийся в работу лишь после того, как давление в данном сечении достигает определенной (заданной) величины. Практический интерес здесь представляет определение запаздывания импульса и его величины, что зависит от длины трубы, вязкости жидкости и коэффициента гидравличе-
Игорь Васильевич Кудинов (к.т.н.), доцент кафедры «Теплотехника и гидравлика». Антон Владимирович Еремин, ассистент кафедры «Теоретические основы теплотехники и гидромеханика».
Сергей Владимирович Колесников (к.т.н.), докторант. Анастасия Эдуардовна Кузнецова, аспирант.
ского сопротивления. К тому же определение величины импульса в сечении x = l представляет решение задачи о гидравлическом ударе.
Математическая постановка задачи по определению давления в данном случае имеет вид
с2 = ЩА + 2ЗМ (t > 0; 0 < x < l) (3)
a x2 at2 dt
p(x, 0) = po = const; (4)
dp(x,0)/dt = 0; (5)
p(0, t ) = pi = const; (6)
dp(l, t) / dx = 0, (7)
где p0 - начальное давление в трубе; pi - давление, приложенное в точке x = 0 и действующее в течение всего времени процесса вплоть до установления стационарного состояния (p1 > p0); l - длина трубопровода.
С целью упрощения математической постановки задачи и процесса получения аналитического решения, а также для нахождения решения наиболее общего вида введем следующие безразмерные переменные и параметры:
2 . 2
гл p" pi С с t x с
© =--; Fo =--; y = -; For =——, (8)
p0 - pi 2al2 l 4a2l2
где © - безразмерное давление; Fo - число Фурье (безразмерное время); y - безразмерная координата; For = const - безразмерный параметр.
С учетом принятых обозначений задача (3) - (7) приводится к виду
d©(y ,Fo) „ d 2©(y ,Fo) d 2©(p,Fo) , n n l)
—^—'- + For-=-—(Fo > 0; 0 < у < 1) (9)
dFo dFo2 d у2
©(y,0) = 1; (10)
d©(y,0) / dFo = 0; (11)
©(0,Fo) = 0; (12)
d©(1,Fo)/dy = 0. (13)
Для удобства получения аналитического решения сделаем замену независимой переменной y по формуле
S = 1 - y. (14)
Относительно переменной задача (9) - (13) примет вид
dMyoLg©^ + For d^2Fo); (Fo>0; 0<^< 1) (15) d^2 dFo r dFo2
©M) = 1; (16)
d©(£, ,0)/ dFo = 0; (17)
©(1,Fo)= 0; (18)
d©(0,Fo) / d^ = 0. (19)
Решение задачи (15) - (19) принимается в виде
©(^,Fo) = 9(Fo)V(^). (20)
Подставляя (20) в (15), получаем
Ford 2ф / dFo2 + dy / dFo + Аф = 0; (21)
d 2у / д^2 + Ху = 0, (22)
где X - некоторая постоянная.
Подставляя (20) в (18), (19), находим
¿у(0) / = 0; (23)
у(1)= 0. (24)
Решение краевой задачи Штурма - Лиувилля (22) - (24) принимается в виде
у(£)= ^(глЕ, /2). (г = 2к -1; к = 1да) (25)
Собственные функции у(^) ввиду однородности уравнения (22) с точностью до постоянного множителя (который в данном случае можно принять равным единице [5]) находятся из (25).
Очевидно, что соотношение (25) удовлетворяет граничным условиям (23), (24). Подставляя (25) в (22), для определения собственных чисел Xк получаем формулу
Хк = г 2л2 /4. (г = 2к -1; к = 1да) (26)
Характеристическое уравнение для уравнения (21) будет
Foгz2 + z + Х к = 0. (27)
Уравнение (27) для каждого собственного числа имеет два корня:
Zik =(-1 ±л/Т-4?01к)/^Ог). (28)
Если дискриминант Б = ^огXк -1)< 0, то из (28) для каждого собственного числа будем иметь два действительных отрицательных корня и (к = 1, да). С учетом найденных значений и 22к решение уравнения (21) для каждого собственного числа будет
Фк (Fo)= С1к ехР(21кFo)+ С2к ехР(22кМ, (29)
где Сд (/ = 1,2; к = 1, да) - неизвестные коэффициенты. Подставляя (25), (29) в (20), находим
©($, Fo) = [С1к ехр(21кFo) + С2к ехр^кFo)]cos^г). (г = 2к -1; к = 1^) (30)
Каждое частное решение (30) точно удовлетворяет уравнению (15) и граничным условиям (18), (19), но ни одно из них не удовлетворяет начальным условиям (16), (17). Для их выполнения составим сумму частных решений
©(^^[Ок ехр(^) + С2к exp(Z2kFo)]cos^г. (г = 2к-1) (31) к=1 V 2 )
Неизвестные коэффициенты С1к и С2к определяются из начальных условий (16), (17). Подставляя (31) в (17), получаем
С1к =-С2к22к /21к. (32) Подставляя (31) в (16), с учетом (32) находим
I С2к (1 - ^г^ 1. (33)
к=1
21к
Соотношение (33) представляет разложение единицы в ряд Фурье по собственным функциям краевой задачи Штурма - Лиувилля на отрезке [0; 1]. Умножим обе части соотношения (33) на ^(/л^/2) и найдем интеграл в пределах от ^ = 0 до й = 1:
J z
о k=1
С
2k
1 -
z2k z1k
cosj^ r ^ j-1 cos^ j ^ jd£ = 0. (j = r = 2k -1; k = 1, да) (34)
Соотношение (34) ввиду ортогональности косинусов принимает вид
С
2k
1-
z2k z1k
\
cos21r §*j-cos^r^
d^ = 0.
Определяя интегралы в (35), находим
C2k = |- I/
т
(
) ■ (r
1 - . (r = 2k -1; k = 1, да
z1k
) ■
(35)
(36)
После определения постоянных интегрирования точное аналитическое решение задачи (15) - (19) находится из (31).
Если в соотношении (28) D = (4ForXk -1)> 0, то будем иметь следующие два комплексных корня:
z1k = У+ip; z2k = У-ip,
где i = V-T ; у = - 0,5 /For; Р = (^4ForXk -1 )/(2For). Частные решения уравнения (21) будут
9tk = exP(У + ip)Fo; Ф2k = ехР(У - ip)Fo ■ (37)
На основе частных решений запишем общий интеграл уравнения (21):
9k(Fo) = C1k exp[(y + i'P)Fo] + C2k exp[(y-ip)Fo], (38)
где Сjk (j = 1,2; k = 1, да) - неизвестные постоянные.
Соотношение (38) можно переписать следующим образом:
9k (Fo) = exp(yFo)[C1k exp(ipFo) + C2k exp(-pFo)]. (39)
При использовании формул Эйлера exp(is) = cos s + i sin s; exp(-is) = cos s - i sin s соотношение (39) приводится к виду
9k (Fo) = exp(yFo)[C1k (cos pFo + i sin pFo) + C2k (cos pFo - i sin pFo)] = = exp(yFo)[(C1k + C2k )cos(PFo) - i (C2k - Cu )sin(pFo)] . (40)
Соотношение (40) с учетом обозначений C1k + C2k = Е1к ; i (C2k - C1k) = B2k будет 9k (Fo) = exp(yFo)[ Bu cos(PFo) - B2k sin(pFo)]. (41)
Подставляя (25), (41) в (20) и составляя сумму частных решений, находим
©fe Fo) = ¿[exp(yFo)[B1kcos(PFo) - B2ksin(PFo)]]cos| r|.
k=1
(42)
(г = 2к -1; к = 1, да)
Для определения постоянных Б^ и Б2к используются начальные условия (16), (17). Подставляя (42) в (17), получаем
Z(YB1k -PB2k)cos| r2^1 = 0.
k=1 V 2 )
Отсюда находим
B2k = yB1k/р.
Подставляя (42) в (16), получаем
(43)
1
0
Соотношение (45) представляет разложение единицы в ряд Фурье по собственным функциям краевой задачи Штурма - Лиувилля на отрезке [0; 1]. Поступая так же, как при определении коэффициента С2к (см. (33)), получим
В1к = 4/(гл). (46)
После определения постоянных интегрирования В^ и В2к точное аналитическое решение задачи (15) - (19) в замкнутом виде находится из (42).
На рис. 1 - 3 даны результаты расчетов по формулам (31), (42) конкретной задачи о распределении давления в находящейся в стальном трубопроводе нефти при следующих условиях:
V = 7,65 • 10-6 м2/с; д = 0,1м; X = 0,01; 30 = 0; 31 = 2 м/с; с = 998 м/с; р = 840кг/м3; к = 1500МПа; Е = 2•105МПа; р0 = 10атм; р1 = 100атм.
Расчеты выполнялись для следующих длин трубопровода:
I = 500 км; I = 100 км; I = 1км.
Соответственно указанным длинам числа Foг для ламинарного режима течения были следующие:
Foг = 6,648•Ю-3; Foг = 0,1662; Foг = 1662,0277.
Аналогично для турбулентного режима
Foг = 0,8964•Ю-3; Foг = 0,02241; Foг = 224,1.
Анализ результатов расчетов позволяет заключить, что изменение давления характеризуется движением гидравлической волны, на фронте которой наблюдается скачок давления от его значения на фронте до величины давления невозмущенного потока. Область, находящаяся за пределами фронта гидравлической волны, оказывается невозмущенной, и давление здесь равно начальному давлению р0 . Отмечается линейная закономерность движения фронта гидравлического возмущения по пространственной переменной во времени = FoFoг 0'5, что подтверждается исследованиями решения уравнения вида (3), выполненными другими авторами [6] применительно к определению температуры в теле с учетом конечной скорости распространения теплоты.
В зависимости от величины числа For наблюдается существенное отличие получаемых результатов. Так, при очень малых его значениях скачок кривых давления имеет место лишь на некоторых начальных участках трубопровода. Например, при
Foг = 6,648•Ю-3 скачок давления наблюдается лишь в диапазоне 0<Fo<0,07 (см. рис. 1). Для всех Fo > 0,07 решение для безразмерного давления полностью совпадает с решением параболического уравнения теплопроводности вида
8©(£, Fo) / 8Fo = 8 2©(^, Fo) / 8^ 2 (47)
при аналогичных начальном (16) и граничных (18), (19) условиях [3, 6, 7].
Из анализа кривых распределения давления следует, что при Fo > 0,07 возрастание давления по длине трубопровода во времени не сопровождается возникновением его скачков - давление с течением времени плавно повышается вплоть до достижения стационарного состояния, при котором давление по всей длине трубопровода устанавливается равным давлению на входер = р\. Отметим, что Fo = 0,07 (при
For = 6,648 -10_3; 1 = 82 км) соответствует времени 1 = 10,46с при турбулентном режиме течения и 1 = 1,28 с при ламинарном режиме.
Для очень малых значений For (длинные трубопроводы или трубопроводы с большим гидравлическим сопротивлением) скачок давления наблюдается лишь на незначительном по длине участке трубопровода вблизи сечения х = 1 и при малых
_7
значениях времени. Например, при For = 10 скачки давления практически заканчиваются при Fo «10_6 . Фронт гидравлического возмущения для данного момента времени перемещается лишь на величину = 0,003, что составляет 0,3 % от всей длины трубопровода.
На всей остальной части трубопровода при Fo > 10_6 повышение давления происходит без скачков вплоть до достижения стационарного состояния, когда давление в сечении х = 1 (^ = 0) становится равным давлению, заданному граничным условием первого рода в сечении х = 0 (^ = 1). Решение задачи (15) - (19) на временном
_6 _7
участке 10 < Fo <да (при For = 10 ) полностью совпадает с классическим точным аналитическим решением параболического уравнения теплопроводности вида (47) при симметричных граничных условиях первого рода [7]. Значения For, близкие к нулевым, могут быть в следующих случаях: большая длина трубопровода; высоковязкая жидкость; большое значение коэффициента гидравлического сопротивления.
1,0 ©
0,8
0,6
0,4
ОД
О ОД 0,4 0,6 0,8 1,0
_3
Р и с. 1. Распределение давления в трубопроводе при Fo г = 6,648 -10
При дальнейшем увеличении времени наблюдается обратная волна гидравлического возмущения со скачком давления в сторону, противоположную скачку давления в прямой волне. По сути это скачок давления, вызывающий гидравлический удар. Причем давление в обратной волне может превышать давление, заданное на входе в трубопровод. С увеличением времени участки трубопровода, на которых на-
блюдается превышение давления по сравнению с давлением р = 100 атм, перемещаются по направлению к точке, где задано это давление, то есть к точке ^ = 1.
С увеличением времени наблюдается периодическое изменение давления в каждой точке трубопровода вплоть до наступления стационарного состояния, при котором давление по всей длине трубопровода становится равным давлению Р1 = 100 атм , заданному на входе.
0,8 § 1,0
Р и с. 2. Распределение давления в трубопроводе при For = 1662, n = 100 (n - число членов ряда решений (31), (42))
-l.ol'-'
500
1000
1500
2000
2500
3000
3500
4000
4500 Fo 5000
Р и с. 3. Распределение давления в точке ^ = 0 во времени Fo при For = 1662 (n = 100)
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Чарный И.А. Неустановившееся движение реальной жидкости в трубах. - М.: Недра, 1975. - 296 с.
2. Гусейнзаде М.А. О характере изменения основных параметров течения жидкости и газа в сложной трубопроводной системе. - М.: ФГУП Изд-во «Нефть и газ» РГУ нефти и газа им. И.М. Губкина, 2005. - 120 с.
3. Кудинов В. А., Кудинов И.В. Получение и анализ точного аналитического решения гиперболического уравнения теплопроводности для плоской стенки // Теплофизика высоких температур. Т. 50. -2012. - № 1. - С. 118-125.
4. Кудинов В.А., Кудинов И.В. Методы решения параболических и гиперболических уравнений теплопроводности. - М.: Книжный дом «Либроком», 2012. - 280 с.
5. Тихонов А.Н., Самарский А.А. Уравнения математической физики. - М.: Наука, 1976. - 724 с.
6. Баумейстер К., Хамит Т. Гиперболическое уравнение теплопроводности. Решение задачи о полубесконечном теле // Теплопередача. - 1969. - № 4. - С. 112-119.
7. Лыков А.В. Теория теплопроводности. - М.: Высшая школа, 1967. - 600 с.
Статья поступила в редакцию 1 апреля 2013 г.
OBTAINING EXACT ANALYTICAL SOLUTIONS OF HYPERBOLIC EQUATION WITHIN SIMPLE PIPELINE IN SURGE
I. V. Kudinov, A. V. Eremin, S. V. Kolesnikov, A.E. Kusnetzova
Samara State Technical University
244, Molodogvardeyskaya st., Samara, 443100
The exact analytical solution of the hyperbolic equation describing viscous fluid flow within a simple pipeline in surge is obtained.
Keywords: viscous incompressible fluid, surge, hyperbolic equation, the Fourier transform, exact analytical solution, orthogonal methods of weighted residuals.
Igor V. Kudinov (Ph.D. (Techn.)), Associate Professor. Anton V. Eremin, Assistant. Sergey V. Kolesnikov, Doctoral Candidate. Anastasiya E. Kusnetzova, Postgraduate Student.