CC BY
9ВЕСТН. САМАР. ГОС. ТЕХН. УН-ТА. СЕР. ТЕХНИЧЕСКИЕ НАУКИ. 2017. № 4 (56)
УДК 536.2(075)46
АНАЛИТИЧЕСКОЕ РЕШЕНИЕ НЕСТАЦИОНАРНОЙ ЗАДАЧИ ТЕПЛООБМЕНА ПРИ ЛАМИНАРНОМ ТЕЧЕНИИ ЖИДКОСТИ В ЦИЛИНДРИЧЕСКОМ КАНАЛЕ*
В.А. Кудинов, А.В. Еремин, Е.В. Стефанюк, В.В. Жуков, Т.Б. Тарабрина
Самарский государственный технический университет Россия, 443100, г. Самара, ул. Молодогвардейская, 244
На основе использования двукратного интегрального преобразования Лапласа - Карсона задача нестационарного теплообмена для движущейся в цилиндрическом канале жидкости разделяется на две независимые задачи - нестационарную и стационарную, каждая из которых определена лишь в определенных диапазонах временной и пространственной координат. Нестационарная задача совпадает с задачей теплопроводности для бесконечного цилиндра (теплообмен протекает как бы в неподвижной жидкости), точные аналитические решения которой известны. Для решения стационарной задачи (задача Гретца - Нуссельта) применяется интегральный метод теплового баланса при использовании дополнительных искомых функций и дополнительных граничных условий. Благодаря использованию дополнительных искомых функций задача Гретца - Нуссельта разделяется на две задачи, в каждой из которых решение уравнения в частных производных сводится к интегрированию обыкновенного дифференциального уравнения. Использование дополнительных граничных условий позволяет, выполняя искомым решением уравнение лишь на границах, получить его выполнение и внутри области с точностью, зависящей от числа приближений.
Ключевые слова: течение жидкости в цилиндрическом канале, нестационарная задача теплообмена, интегральные преобразования Лапласа - Карсона, ортогональный метод Бубнова - Галеркина, интегральный метод теплового баланса, дополнительные искомые функции, дополнительные граничные условия.
Введение
В процессе решения нестационарной задачи теплообмена при ламинарном течении жидкости в цилиндрическом канале происходит ее физически обоснованное разделение на нестационарную и стационарную задачи, границы применимости которых определяются соответствующими диапазонами временной и продольной пространственной переменных. Последующее решение каждой из этих задач выполняется раздельно. Нестационарная задача теплообмена, по сути, представляет собой задачу теплопроводности для сплошного цилиндра (теплообмен протекает как бы в неподвижной жидкости), точные аналитические реше-
*Работа выполнена при финансовой поддержке Минобразования и науки РФ в рамках базовой части государственного задания ФГБОУ ВО «СамГТУ» (код проекта: 1.5551.2017/8.9).
Василий Александрович Кудинов (д.ф-м.н., проф.), заведующий кафедрой «Теоретические основы теплотехники и гидромеханики».
Антон Владимирович Еремин (к.т.н.), доцент кафедры «Теоретические основы теплотехники и гидромеханики».
Екатерина Васильевна Стефанюк (д.т.н., доц.), профессор кафедры «Теоретические основы теплотехники и гидромеханики».
Виталий Владимирович Жуков, аспирант.
Тамара Борисовна Тарабрина (к.п.н.), доцент кафедры «Теоретические основы теплотехники и гидромеханики».
ния которой известны. Наибольшую трудность представляет решение стационарной задачи с учетом движения жидкости (задача Гретца - Нуссельта) [1-5]. Эта задача впервые была решена Гретцем в 1885 г. [4] и независимо от Гретца Нуссельтом в 1910 г. [5]. Уточненное решение задачи Гретца - Нуссельта, данное в [1], представляет бесконечный ряд, плохо сходящийся при малых значениях продольной координаты. К тому же в зависимости от величины поперечной координаты решение содержит функции Бесселя различного (в том числе и дробного) порядка. Такое решение малопригодно для инженерных приложений. В настоящей работе приводится метод решения задачи Гретца - Нуссельта, основанный на определении дополнительных искомых функций и дополнительных граничных условий в интегральном методе теплового баланса [7, 8].
Для решения нестационарной задачи теплообмена при течении жидкости в цилиндрическом канале введем следующие допущения: течение жидкости и процесс теплообмена стационарны; жидкость несжимаема, ее физические свойства постоянны; профиль скорости не изменяется по длине трубы; температура жидкости на входе в трубу неизменна по сечению и равна to; температура внутренней поверхности стенки трубы постоянна и равна tс, причем tс Ф ^; начальная температура жидкости равна ^; внутренние источники тепла и диссипация энергии не учитываются; переносом теплоты теплопроводностью вдоль оси трубы пренебрегается (рис. 1).
Рис. 1. Схема стабилизированного ламинарного течения жидкости в круглой трубе
Математическая постановка задачи имеет вид [1-4]
л, т) Э^Х, Л, т)_ а Э^Х, Л, т) + Э 2^Х, Л, т)
Эт х Эл (т > 0; п > 0; 0 < ^ < г)
Х эх
+ а-
ЭХ2
t (Х, л,0)=^.
?
t (Х ,0, т) = ^
t (г, л, т) = ^,
?
Эt (0, л, т)/ЭХ = 0,
(1)
(2)
(3)
(4)
(5)
где 122
t - температура;
X - продольная и поперечная координаты; г - радиус трубы;
to - температура жидкости на входе в трубу (при ^ = 0); ^ — температура стенки; ^ - начальная температура; a - коэффициент температуропроводности;
юх = 2юср (1 - X2 / г 2) - распределение скорости по координате X (0 < X £ г ); Юср = 0,5 т
ср ' тах — средняя скорость; ®тах - максимальная скорость ламинарного течения (см. рис. 1). Введем следующие безразмерные переменные и параметры:
©=& х=Ро=01
^ - ^ ; у =Х ; 2ре Г ; ре = V ; Г2, (6)
г а
где 0 - относительная избыточная температура;
у, х - безразмерные поперечная и продольная координаты; Ре - число Пекле; Ро - число Фурье. С учетом принятых обозначений задача (1) - (5) будет
Э0(у, х,Ро) + - 2 Э0(у, х,Ро) = 1 Э0(у, х,Ро) + Э 20(у, х,Ро) ЭРо Э х у Э у Э у2
(Бо > 0; х > 0; 0 < у < 1)
0( у , х,0) = Тн.
?
0( у,0, Ро) = 1,
?
0(1, х,Бо) = 0;
?
Э0(0, х,Бо)/Эу = 0
где Т. =(*н - ¿с )/(^ - tc).
(7)
(8) (9)
(10) (11)
1. Математическая постановка задачи
Для решения задачи (7) - (11) применим двукратное преобразование Лапласа -Карсона по безразмерному времени Ро и относительной координате х [2, 3]:
Т *( у, 5,Бо) = 5| 0( у, х, Бо) ехр(-5х)ох; (12)
Т (у, 5, Р) = Р | Т *( у, 5 ,Ро)ехр(-РРо)^ Ро, (13)
0
где Р, 5 - параметры двукратного преобразования Лапласа - Карсона.
0
Соотношение (13) с учетом (12) примет вид
Т (у, Б, Р) = РБД 0( у, х,Бо)ехр(-Бх - РБо)й№о.
0 0
Задача (7) - (11) в области изображений будет
0
2
Р[Т(у,Б,Р) - Тн] + Б(1 - у2)[Т(у,Б,Р) - Т]-
1 ЭТ(у, Б, Р) Э Т(у, Б, Р) = 0,
у ЭУ ЭУ2 (15)
Т(1, Б,Р) = 0; (16) ЭТ (0, Б, Р)/Эу = 0. (17)
Решение задачи (15) - (17), следуя ортогональному методу Бубнова - Галер-кина, в первом приближении примем в виде
Т (у, Б, Р) = Ц(Б, Р)ф(у), (18)
где ^(Б, Р) - неизвестный коэффициент - изображение; Ф1 (у) - координатная функция, принимаемая в виде
фМ =1 - у2. (19)
Очевидно, что соотношение (18) при использовании координатной функции вида (19) точно удовлетворяет граничным условиям (16), (17). Для определения неизвестного коэффициента П^Б,Р) составим невязку уравнения (15) и потребуем ортогональности невязки и координатной функции ф1 (у), то есть
| {уР[Т (у, Б, Р) - Тн ] + уБ (1 - у2)[Т (у, Б, Р) - Т ]-
0
- уу* . (20)
эу Эу2 ]
Подставляя (18) в (20), находим 1
| {уР^(Б, Р)(1 - у 2) - Тн ]+ уБ (1 - у 2)[у1(Б , Р)-Т ] +
0
+ 2У1 (Б, Р)у + 2У1 (Б, Р)у }(1- у2)йу = 0. (21)
Определяя интегралы в (21), относительно неизвестного коэффициента \^(Б, Р) получаем формулу
, Р) = РР1 +Б1'2 (22)
Р^3 + £Р4 +
1
где ^ = Тн|у(1 - у2)йу ; ¥2 = |у(1 - у2)2йу; Я) = |у(1 - у2)2йу ;
0
0
0
Р _ | У(1 - У2 )3 dy ; Р = 4| у(1 - у2 )йу .
(23)
0 0
В пространстве оригиналов соотношение (18) принимает вид
( * \ Р
Ф1(у) при X »-¡^Бо;
0( у, х,Бо) =
Ррхр
- Р5Ро
рз .
Р2
(
ехр
- х
Р Х
Л
Рз р
ф1(у) при х <—Бо.
(24)
3
Верхняя строка формулы (24) совпадает с решением в первом приближении задачи (7) - (11) при равенстве нулю конвективного члена (второе слагаемое в левой части) уравнения (7), то есть с решением нестационарной задачи теплопроводности для бесконечного цилиндра. Нижняя строка этой формулы совпадает с решением в первом приближении стационарной задачи конвективного теплообмена, то есть при равенстве нулю первого члена в левой части уравнения (7). Таким образом, для областей теплообмена, которых не достигло возмущение, вызванное начальной температурой жидкости на входе в канал (при х = 0), нестационарный теплообмен происходит как бы в неподвижной жидкости, то есть перенос теплоты происходит только теплопроводностью. Данная задача представляет нестационарную задачу теплопроводности для бесконечного цилиндра, классические точные аналитические решения которой известны.
Для областей, подверженных влиянию температурных условий на входе в трубу (жидкость при течении достигла этих областей), теплообмен не зависит от времени (не зависит от начальной температуры и полностью определяется течением среды, то есть задача становится стационарной с учетом конвективного переноса теплоты по оси х.
Подобное разделение решения нестационарной задачи конвективного теплообмена при течении жидкостей в трубах и плоских каналах на стационарную и нестационарную составляющие физически обосновано и в других работах [1, 2].
2. Решение задачи Гретца - Нуссельта
Математическая постановка стационарной задачи теплообмена (задача Грет-ца - Нуссельта) имеет вид
(1 - У2)
Э0(у, X) _Э2 0(у, х) +1 Э0(у, х).
Э х Э у2 у Э у (х > 0; 0< у < 1)
0( у ,0) = 0.
?
Э0(0, х)/Эу _ 0.
?
0(1, х)_ 1
(25)
(26)
(27)
(28)
Процесс нагрева среды разделим на две стадии по координате х: 0 < х < Х1 и х1 < х < ¥. Для этого введем движущуюся по координате у границу (фронт температурного возмущения), разделяющую исходную область на две подобла-
сти: прогретую ql(x) < у < 1 и непрогретую 0 < у < q1(х) . Здесь х) - функция, определяющая продвижение границы раздела по координате у в зависимости от продольной координаты х (фронт температурного возмущения), являющаяся дополнительной искомой функцией (рис. 2 а). При этом в области, не затронутой нагревом, сохраняется температура (0, заданная на входе в канал. Первая стадия процесса заканчивается при достижении движущейся границей центра трубы (у = 0). Во второй стадии изменение температуры происходит по всему объему
среды 0 < у < 1. Здесь в рассмотрение вводится новая функция q2 (х), характеризующая изменение температуры во времени в центре трубы q2(х) = ©(0, х).
2.1. Первая стадия процесса
Для упрощения получения решения сделаем замену переменной у, отсчитываемой от центра трубы, новой переменной г = 1 - у , отсчитываемой от ее стенки. Тогда задачу нагрева жидкости для первой стадии процесса можно сформулировать в виде (рис. 2 б)
(2г - г2 )Э©М = Э!©М - Л. Э©(г, х) ; (29)
Эх Эг2 1 - г Эг
(0 < х < х1; 0 < г < ^(х))
©(0, х) =1; (30)
©(^, х) = 0; (31)
Э©(^, х)/Эг = 0, (32)
где соотношения (31), (32) представляют условия тепловой изоляции подвижной границы (условия сопряжения прогретой и непрогретой зон). Первое из них означает равенство температуры жидкости на движущейся границе ее температуре на входе в канал ©(у,0) = 0. Второе - отсутствие теплового потока за пределы фронта температурного возмущения (за пределами подвижной границы). Отметим, что в задаче (29) — (32) отсутствует граничное условие вида (26), так как за пределами фронта температурного возмущения она вообще не определена. В задаче (29) — (32) отсутствует также граничное условие вида (27), так как оно не влияет на процесс теплообмена в первой его стадии.
Потребуем, чтобы искомое решение задачи (29) — (32) удовлетворяло не исходному уравнению (29), а осредненному по толщине термического слоя 0 < г < ql( х) уравнению. Определяя интеграл от уравнения (29) по переменной 2 в пределах от г = 0 до г = ql(x) , получаем следующее интегральное уравнение (интеграл теплового баланса):
} (2 г - г 2 ^ йг = ]
Э 2©(г, х)__1_ Э©(г, х)
Эг2 1 - г Эг
йг
. (33)
0
0
Рис. 2. Расчетная схема теплообмена
Решение задачи (29) - (32) принимается в виде
п
&{г, х) = £ ак (3 )Х, (34)
к=0
где ак(31) - неизвестные коэффициенты.
Для определения коэффициентов ак(31) используются граничные условия (30) - (32). Подставляя (34) (ограничиваясь тремя членами ряда) в (30) - (32), относительно ак (#1), (к = 0,1,2) будем иметь систему трех алгебраических линейных уравнений. После определения из решения этой системы коэффициентов ак (#1) соотношение (34) принимает вид
0^, х) = (1 - х / ^)2. (35)
Подставляя (35) в (33), получаем
31 , У _ Л О 31
Ь \ (2 х - х2)
0
1
V
, 2 1 2 ах =
Г .
31 0 х -1
31
Определяя интегралы в (36), находим
32 (2з2 - 9дх +10)031 / ах = 60.
- —
V 312 31
0х
(36)
(37)
Уравнение (37) - нелинейное. Для получения его приближенного аналитического решения введем новую, независимую переменную п = х / х1, где Х1 - расстояние, на которое перемещается бывшая на входе в канал жидкость по оси х после достижения фронтом температурного возмущения 31(х) координаты х = 1 , то есть это граница раздела возмущенного граничным условием на входе в канал участка трубы по оси х (0 £ х £ Х1) и невозмущенного (х1 £ х < ¥, где теплообмен протекает как бы в неподвижной жидкости. Напомним, что для невозмущенного по оси х участка (х > х1) температурное состояние жидкости определяется из решения нестационарной задачи теплопроводности для бесконечного сплошного цилиндра, коэффициент теплопроводности материала которого равен теплопроводности жидкости при ламинарном течении (см. выше решение, полученное с использованием преобразования Лапласа - Карсона и метода Бубнова -Галеркина).
С учетом новой независимой переменной уравнение (37) будет
q12 (2д22 - 9q1 +10)йц1 = 60 х2йц ,
?
(0 < п < 1; 0 < д,^ 1). (39)
Решение уравнения (38) при граничном условии ql (0) = 0 принимается в виде
^ (П) = Цп", (40)
где ц , " - некоторые постоянные.
Решение (40) должно также удовлетворять условию ql (1) = 1, следующему из неравенств (39). Отсюда находим ц = 1.
ql (п)=п".
Потребуем, чтобы соотношение (40) удовлетворяло осредненному по переменной п уравнению (38):
q12(2q12 - 9q2 +10) ^ - 60 х2 йц
йп = 0 .
(41)
Подставляя (40) в (41), находим 1
2
йц1
п2"(2п2" -9п" +10)^--60х,)
йц = 0.
(42)
Вычисляя интегралы в (42), получаем
89
60х2--= 0.
60 (43)
Из соотношения (43) находим х2 = 0,0247.
Для определения 1 потребуем, чтобы соотношение (40) (при ц = 1) удовлетворяло уравнению (38) в точке п = 1. Подставляя (40) в (38), применительно к точке п = 1 будем иметь
^(1)(2^(1) - 9^(1) +10)
йд2
йп
-60 х = 0.
п=2 (44)
Последнее соотношение приводится к следующему уравнению относительно " :
31- 60х2 = 0. (45)
Из решения уравнения (45) находим " = 0,4944.
Соотношение (40) с учетом найденных значений х, и " для независимой переменной х приводится к виду
ql( х) =
V х2
= 6,2305 х
0,4944
0
0
Уточняя полученные значения х, и ", в соответствии с методом, изложенным в [6], получаем
q1(x) = 7,1464х °'54. (46)
Результаты расчетов по формуле (46) практически совпадают с точным решением уравнения (38). Следовательно, формула (46), являясь приближенным аналитическим решением уравнения (38), с достаточной точностью может быть использована в соотношении (35). Подставляя (46) в (35), находим
©( г, х) = (1 - 0.13993гх "°'54)2
(47)
Формула (47) представляет решение задачи (29) - (32) в первом приближении. Результаты расчетов безразмерных температур в первом (по формуле (47)) и
последующих приближениях при х = 0,005; 10 ; 10-5 в сравнении с точным
решением [1, 2] представлены на рис. 3-5. Расхождение результатов при х = 10-4 составляет более 15 % и увеличивается при уменьшении х .
0,24
0
0,18
0,12
0,06
\
5, 10, 15 N \ 1,3
\ V 1
3 1 ^^ 5, 10, 15
0,20
0,27
0,34
0,41
0,48
Рис. 3. Распределение температуры при х = 0,005 :
1, 2, 3, 5, 10, 15 - соответственно первое, второе, третье, пятое, десятое,
пятнадцатое приближения;----численное решение;
--приближенное аналитическое решение
1,0
0
0,8
0,6
0,4
0,2
\ 3 \ 2 \Л \ 1 \ XV /15 10
/ 5
0
0,025
0,05
0,075
* 0,1
Л-4
Рис. 4. Распределение температуры при х = 10- : 1, 2, 3, 5, 10, 15 - соответственно первое, второе, третье, пятое, десятое,
пятнадцатое приближения;----численное решение;
--приближенное аналитическое решение
1,0 0 0,8
0,6
0,4
0,2
V \\\ \
\ 2 VN \ 1 \ т 10 /
\ \ \ 5
0
0,01
0,02
0,03
0,04
Рис. 5. Распределение температуры при х = 10 -5: 1, 2, 3, 5, 10, 15 - соответственно первое, второе, третье, пятое, десятое,
пятнадцатое приближения;----численное решение;
--- - приближенное аналитическое решение
Для повышения точности найдем решение задачи (29) - (32) во втором приближении. При этом будем использовать шесть членов ряда (34), для нахождения неизвестных коэффициентов которого основных граничных условий (30) - (32) будет недостаточно. В связи с этим необходимо ввести три дополнительных граничных условия, определяемых из основных граничных условий (30) - (32) и уравнения (29) и имеющих вид [6, 7]
^©0^-^©^ = 0; Э 2©(^, х)/Эг2 = 0; Э3©(^, х)/Эг3 = 0. (48) Э г 2 Эг
Подставляя (34) во все основные и дополнительные граничные условия, для определения коэффициентов ак , (к = 0,5) получаем цепочную систему шести алгебраических линейных уравнений. После нахождения из решения этой системы неизвестных коэффициентов ак, (к = 0,5) соотношение (34) принимает вид
©(г х) = 2 -20г -20г2 + 20(Ц + 2)г3 - 5(3ц + 8)г4 + 4^1 + 3)г5
' qlr qlr q13r q14r q15r
(49)
г = q2 + 8 где .
Подставляя (49) в интеграл теплового баланса (33), получаем
g22 (-5^ - 41д23 + 150д22 + 552д2 ) ^ + 2520д2 + 20160 = 0
йх . (50)
Используя метод, изложенный выше, при начальном условии ql(0) = 0 находим следующее решение уравнения (50):
q1( х) = 5,8676х0,43302 (51)
Положив в (51) gl(х) = 1, находим х = х, = 0,0168.
Анализ результатов расчетов по формуле (51) позволяет заключить, что с увеличением числа приближений п величина продольной координаты х, при которой фронт температурного возмущения достигает координаты г = 1, уменьшается. И в пределе при п ® ¥ х ® 0 . Этот результат находится в полном соответствии с гипотезой о бесконечной скорости распространения теплоты, лежащей в основе вывода параболического уравнения (29).
Соотношения (49), (51) определяют решение задачи (29) - (32) во втором приближении. Результаты расчетов по формуле (49) в сравнении с точным решением даны на рис. 3-5. Их анализ позволяет заключить, что во втором приближении происходит существенное уточнение решения. И, в частности, при
х = 10-4 расхождение результатов расчетов составляет менее 10 %.
Для дальнейшего повышения точности ограничимся в соотношении (34) девятью членами ряда, следовательно, для нахождения девяти неизвестных коэффициентов ак к трем основным (30) - (32) и трем дополнительным (48) необходимо использование трех дополнительных граничных условий, которые в данном случае имеют вид
Э3©(0, х) Э2©(0, х) Э©(0,х) 0 Э4©( )/э 4 0 ——3---—2---^-= 0; Э ©(^,х) Эг = 0;
Эг Эг2 Эг
Э5©(^, х)/Эг5 = 0. (52)
Для нахождения решения задачи (29) - (32) в третьем приближении используются основные (30) - (32) и дополнительные (48), (52) граничные условия. Подставляя (34) во все перечисленные условия, относительно неизвестных коэффициентов ак ^1), (к = 0,8) получаем цепочную систему девяти алгебраических уравнений. После подстановки (34) с учетом найденных значений коэффициентов ак ^1) в интеграл теплового баланса (33) находим
10g18 + 105q17 + 80q16 - 1803g15 - 12078g14 - 42840g13 +
+138600^ -498960?, = 3492720. (53)
Используя изложенный выше метод, при начальном условии (0) = 0 получаем
ql( х) = 5,56428х0,39974. (54)
Положив в (54) ?2(х) = 1, находим х = х2 = 0,013654.
Результаты расчетов по формуле (34) в третьем приближении (п = 8) даны на
рис. 3-5. Их анализ позволяет заключить, что при х = 10-4 расхождение с точным аналитическим решением [1, 2] не превышает 5 %. С увеличением числа приближений точность существенно возрастает. Так, в пятнадцатом приближении точное и приближенное решения практически совпадают.
2.2. Вторая стадия процесса
Вторая стадия теплового процесса, соответствующая продольной координате х — х,, характеризуется изменением температуры уже по всему сечению канала вплоть до состояния, когда вся жидкость прогреется до температуры стенки. Для этой стадии понятие термического слоя теряет смысл, и в качестве дополнительной искомой функции принимается температура в центре канала ©(1, х) = g2(х) (рис. 2 б). Математическая постановка задачи в данном случае имеет вид
(2г - г2 )Э©(г, х) = Э2©(г, х)__1_ Э©(г, х)
Эх Эг2 1 - г Эг ; (55)
(х - х1; 0 < г < 1)
©(0,х) = 1; (56)
©(1, х)=?2(х); (57)
Э©(1, х) / Эг = 0 (58)
Начальным условием задачи (55) - (58) будет формула (35) при ?1(х1) = 1:
©(г, х1 ) = (1 - г)2 . (59)
Осредняя дифференциальное уравнение (55) по всему объему среды (0 £ г £ 1), получаем интеграл теплового баланса вида
}(1 - г)(2г - г2 ^ = }
„ чЭ2©(г, х) Э©(г, х) (1 - г)-Ч—------
Эх2 Эг
(60)
0
0
(62) (63)
Решение задачи принимается в виде
п
9(г, х) = £ Ък (^) гк
к=0 , (61)
где Ък (#2) _ неизвестные коэффициенты.
После определения коэффициентов Ък, (к = 0, 1, 2) из граничных условий (56) - (58) соотношение (61) будет
9(г, х) = 1 + (#2-1)(2 - г)г
Подставляя (62) в интеграл теплового баланса (60), находим (д2(х)/ (х+12[д2(х)-1] = 0
Интегрируя уравнение (63), получаем
#2 (х) = -~ехр[-12(С - х)] +1
12 , (64)
где С - постоянная интегрирования, определяемая из начального условия хО = 0. После определения постоянной интегрирования соотношение (64) принимает вид
(х) = 1 - ехр[12(х1 - х)], (65)
где = 0,02731 получено в первом приближении первой стадии процесса.
Подставляя (65) в (62), находим
9(г, х) = - 2)ехр[12(хх - х)]+1 (66)
Непосредственной подстановкой можно убедиться, что соотношение (66) точно удовлетворяет граничным условиям (56) - (58), начальному условию (59) и интегралу теплового баланса (60). Уравнение (55), как это следует из (60), в данном случае удовлетворяется лишь в среднем.
Результаты расчетов по формуле (66) в сравнении с точным решением [1, 2], даны на рис. 6. Из анализа следует, что их расхождение составляет 16-18 %.
Для повышения точности найдем решение задачи (55) - (58) во втором приближении с привлечением дополнительных граничных условий. Первое дополнительное условие в данном случае совпадает с первым соотношением из (48), а второе и третье имеют вид
Э 20(1, х) / Эг2 = 0,5(д2 (х) / (х; Э3©(1, х) / Эгъ = 0 (67)
Подставляя (61) в основные и дополнительные граничные условия и ограничиваясь шестью членами ряда, относительно неизвестных коэффициентов
Ък (#2), (к = 0,5) получаем цепочную систему шести алгебраических уравнений.
После определения неизвестных Ък (#2) соотношение (61) будет
9( г, х) = 1 + (16г4 + 55г3 + 60г2 + 10г + 20) +
г (а2 ( 4 13 3 2 1 1 ^
+--^I - г4 + — г3 - г2 + — г + — I, (68)
3 (х { 4 4 2)
где г = -1.
1,0
0
0,8
0,6
0,4
0,2
1 1
/ 3 7 / 4
X 1
/2
/ 3 4
0
0,1
0,2
0,3
0,4
0,5
0,6
0,7
0,8
1,0
Рис. 6. Распределение температуры при х = 0,05; х = 0,15. 1, 2, 3 - соответственно первое, второе, третье приближения; 4 - точное решение [1];----численное решение
Подставляя (68) в интеграл теплового баланса (60), находим
129 <*2х) + 6484 ^^ + 40320q2 (х) = 40320 Граничные условия к уравнению (69) имеют вид
q2(х1) = 0. Э^2(х1)/Эх = 0
?
Интегрируя (69), находим
(69)
(70)
д2 (х) = 1 + С ехр(-7,26986х) + С2 ехр(-42,99371х), (71)
где С1, С2 - постоянные интегрирования, определяемые из граничных усло -вий (70).
Подставляя (71) в (70), получаем
С1 = ехр(7,26986х1); С2 = ехр(42,99371х1)
Подставляя (72) в (71), находим
(72)
(73)
q2 = 1 + ехр[- 7,26986л -х)]+ехр[- 42,99371(х1 - х)]
где х1 = 0,0168 - получено во втором приближении первой стадии процесса.
Соотношения (68), (73) представляют решение задачи (55) - (58) во втором приближении.
Отметим, что константы, расположенные под знаком экспонент, незначительно отличаются от первых двух собственных чисел, получающихся при решении задачи (25) - (28) точными аналитическими методами [1, 2]. Точные значения этих чисел 1 = 7,3136; 1,2 = 44,6095.
Результаты расчетов по формуле (68) даны на рис. 6. Их анализ позволяет заключить, что максимальное расхождение результатов второго приближения по сравнению с точным решением не превышает 3 %.
Для повышения точности найдем решение задачи (55) - (58) в третьем приближении. Для этого будем использовать еще два дополнительных условия вида
4 Э49(1,г) = 1 (#2(х) - 2 (#2(х). Э59(1,2) = 0 3 Эг4 2 (х2 (х ' Эг5 .
Подставляя (61) во все перечисленные условия и ограничиваясь девятью членами ряда, относительно неизвестных коэффициентов Ък, (к = 0, 8) получаем цепочную систему девяти алгебраических уравнений. Подставив (61) (с учетом найденных Ък, (к = 0,8)) в интеграл теплового баланса (60), получаем
15649 (За? 382469 (2а7 72767 96 96 --11 +--11 +--11 + — #2 = —
683020800 (х3 56918400 (х2 254100 (х 55 2 55 . Решение уравнения (75) имеет вид
(х) = С1 ехр(у1х) + С2 ехр(у2 х) + С3 ехр(у3 х) +1
(75)
(76)
где С1, С2, С3 - константы интегрирования;
=-7,321467. =-42,789754.
п3 =-243,174485
Константы интегрирования С1, С2, С3 находятся из граничных условий
#2(х1) = 0. (#2(х1)/(х = 0. (2#2(х1)/ёх2 = 0 (77)
После их определения соотношение (76) принимает вид
х) = ехр^ (х - х1)] + А2 ехр^2 (х - х^] + А3 ехр[у3 (х - х1)] +1
(78)
где А =-1,243873; А2 = 0,250502; А3 =-0,006629.
Соотношения (61) и (78) с учетом найденных значений коэффициентов Ък, (к = 0,8) представляют решение задачи (55) - (58) в третьем приближении. Результаты расчетов по формуле (61) в третьем приближении (п = 8) даны на рис. 6. Их анализ позволяет сделать заключение о существенном повышении точности решения. И, в частности, в диапазоне х1 = 0,013654£ х<¥ расхождение с точным решением не превышает 1,5 %.
Выводы
1. С использованием двукратного интегрального преобразования Лапласа -Карсона, ортогонального метода Бубнова - Галеркина и интегрального метода теплового баланса получено аналитическое решение нестационарной задачи теплообмена при течении жидкости в цилиндрическом канале. Математически показано, что решение задачи разделяется на два решения, одно из которых представляет решение нестационарной задачи, а второе - стационарной. Пространственно-временные границы раздела этих двух задач определяются на основе использования интегрального преобразования Лапласа - Карсона и ортогонального метода Бубнова - Галеркина. Дальнейшее решение выполняется раздельно для нестационарной и стационарной задач.
2. Для решения стационарной задачи (задача Гретца - Нуссельта) используется интегральный метод теплового баланса, согласно которому процесс теплообмена разделяется на две стадии по продольной пространственной переменной. На каждой из этих стадий вводятся дополнительные искомые функции и дополнительные граничные условия. Введение дополнительных искомых функций позволяет в обеих стадиях сводить решение уравнения в частных производных к интегрированию обыкновенного дифференциального уравнения.
3. Дополнительные граничные условия находятся в таком виде, чтобы их выполнение искомым решением в граничных точках и на фронте температурного возмущения было эквивалентно выполнению исходного дифференциального уравнения (в частных производных) и производных от него различного порядка по переменной z в граничных точках. Показано, что выполнение уравнения в граничных точках приводит к его выполнению и внутри области с точностью, зависящей от числа приближений (числа дополнительных граничных условий) [7].
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Петухов Б.С. Теплообмен и сопротивление при ламинарном течении жидкости в трубах. - М.: Энергия, 1967. - 412 с.
2. Цой П.В. Системные методы расчета краевых задач тепломассопереноса. - М.: Изд-во МЭИ, 2005. - 568 с.
3. Еремин А.В., Стефанюк Е.В., Рассыпнов А.Ю., Кузнецова А.Э. Нестационарный теплообмен в цилиндрическом канале при ламинарном течении жидкости // Вестник Самарского государственного технического университета. Сер. Физико-математические науки. - 2013. - № 4 (33). - С. 122-130.
4. Graetz L. Über die Wärmeleitungsfähigkeit von Flüssigkeiten // Ann. Phys. Chem. 1885. No. 25. S. 337-357.
5. Nusselt W. Die Abhängigkeit der Wärmeübergängszahl von der Rohrlange // Z. Ver. Deutsch. Ing. 1910. Bd. 54. S. 1154-1158.
6. Кудинов В.А., Стефанюк Е.В., Антимонов М.С. Аналитические решения задач теплообмена при течении жидкости в плоскопараллельных каналах на основе определения фронта температурного возмущения // Инженерно-физический журнал. - 2007. - Т. 80. - № 5. - С. 176-186.
7. Стефанюк Е.В., Кудинов В.А. Дополнительные граничные условия в нестационарных задачах теплопроводности // Теплофизика высоких температур. - 2009. - Т. 47. - № 2. - С. 269-282.
8. Федоров Ф.М. Граничный метод решения прикладных задач математической физики. - Новосибирск: Наука, 2000. - 220 с.
Статья поступила в редакцию 5 сентября 2017 г.
ANALYTICAL SOLUTION OF THE NON-STATIONARY HEAT TRANSFER PROBLEM FOR A LAMINAR FLUID FLOW IN A CYLINDRICAL CHANNEL
V.A. Kudinov, A. V. Eremin, E. V. Stefaniuk, V. V. Zhukov, T.B. Tarabardina
Samara State Technical University
244, Molodogvardeyskaya st., Samara, 443100, Russian Federation
On the basis of the use of the twofold Laplace-Carson integral transformation, the problem of non-stationary heat transfer for a fluid moving in a cylindrical channel is divided into two independent problems - non-stationary and stationary, each of which is defined only in certain ranges of time and space coordinates. The non-stationary problem coincides with the problem of thermal conductivity for an infinite cylinder (heat transfer proceeds as if in a fixed liquid), the exact analytical solutions of which are known. To solve the stationary problem (the Gretz-Nusseltt problem), an integral heat balance method is used when using additional unknown functions and additional boundary conditions. Owing to the use of additional unknown functions, the Gretz-Nusselt problem is divided into two problems, in each of which the solution of the partial differential equation is reduced to the integration of an ordinary differential equation. The use of additional boundary conditions makes it possible, by carrying out the desired solution, the equation only at the boundaries, and obtaining its fulfilment within the region with an accuracy that depends on the number of approximations.
Keywords: fluid flow in a cylindrical channel, non-stationary heat transfer problem, Laplace -Carson integral transforms, Bubnov - Galerkin orthogonal method, integral heat balance method, additional unknown functions, additional boundary conditions.
Vasily A. Kudinov (Dr. Sci. (Phys.& Math.)), Professor. Anton V. Eremin (Ph.D. (Techn.)), Associate Professor. Ekaterina V. Stefaniuk (Dr. Sci. (Techn.)), Professor. Vitaly V. Zhukov, Postgraduate Student.
Tamara B. Tarabardina (Ph.D. (Pedag.)), Associate Professor.
