Динамика кондуктивно-ламинарного свободно конвективного движения жидкости в квадратной области Текст научной статьи по специальности «Физика»

УДК 536.25

ДИНАМИКА КОНДУКТИВНО-ЛАМИНАРНОГО СВОБОДНО КОНВЕКТИВНОГО ДВИЖЕНИЯ ЖИДКОСТИ В КВАДРАТНОЙ ОБЛАСТИ

М.И. Слюсарев, Е.Д. Чертов, В.И. Ряжских, В.Г. Стогней, А.А. Богер

Использованием конечного интегрального преобразования Фурье и одностороннего интегрального преобразования Лапласа получено аналитическое решение линеаризованной системы уравнений Обербека-Буссинеска в переменных Гельмгольца для квадратной области с разными значениями постоянных температур на вертикальных стенках и отсутствием теплового потока через верхнее и нижнее основания

Ключевые слова: нестационая свободная конвекция, аналитическое решение, квадратная область

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

Такие модельные представления позволили уже получить ряд точных решений задач о свободной конвекции у бесконечной вертикальной стенки [6], в плоском вертикальном канале [7], в прямо -угольной каверне с отношением высоты к ширине намного больше единицы [8]. В связи с этим нет

Слюсарев Михаил Иванович - ВГТА, канд. техн. наук, доцент, тел. (4732) 55-35-54

Чертов Евгений Дмитриевич - ВГТА, д-р техн. наук, профессор, тел. (4732) 55-35-54

Ряжских Виктор Иванович - ВГТА, д-р техн. наук, профессор, тел. (4732) 55-35-54

Стогней Владимир Григорьевич - ВГТУ, канд. техн. наук, профессор, тел. (4732) 55-53-54

Богер Андрей Александрович - ВГТА, канд. техн. наук, доцент, тел. (4732) 55-35-54

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

Постановка задачи. Пусть имеется квадратная плоская область (рис. 1), полностью заполненная вязкой несжимаемой жидкостью, находящейся в начальный момент в состоянии покоя с температурой tc. Далее, например, температура левой стороны области скачком изменяется до температуры 4, при этом верхнее и нижнее основания теплоизолированы, т.е. тепловой поток отсутствует д=0, а темпера-

Рис. 1. Расчетная схема

тура правой стенки остаётся равной первоначельной температуре жидкости (4>4)-

В соответствии с выбранной системой координат уравнения Обербека-Буссинеска в переменных Гельмгольца примут следующий вид:

да дТ да дТ дП _ д2П + д2П + 5Т ;

50 + дУ дХ дХ дУ ~ дХ2 + дУ2 + Г дХ ’

дТ дТ дТ дТ дТ _ _1_( д2Г д2Т ^ ;

50+ дУ дХ ~дХ~дУ ~ Рг[дХ2 + дУГ) ’

п(Х,У,0) _Т(Х,У,0) _ Т(Х,У,0) _ 0 ; (4)

Т(0,У,0) _ Т(1, У,0) _ Т(Х,0,0) _ Т(Х, 1, 0) _ 0; (5)

5T(G,Y,0) _дТ(1, Y,0)

CX CX

d¥(X,G,0) d¥(X,1,0)

CY CY

CT (X,G,0) CT (X,1,0)

_ G;

_ G;

(6)

(7)

д7 д7

Т (0,7,0) = 1, Т (1,7,0) = 0, (8)

где по определению функции вихря и тока есть

_ ЗУ ди т 5Т у

0 =----------, и =-------, V =--------; и = фф;

дХ д7 д7 дХ ^

X = х/Н , 7 = у/к ; и = и/*, V = и/ * ;

Т = (t-tc)/(-tc); *=ЧН; т=N*;

вг = PgН3 (tН - tc )у2 - критерий Грасгофа; Рг = V а

- число Прандтля; т, х, у - текущее время и декартовы координаты; и , и, t - локальные компоненты скорости и температура; к - длина стороны каверны; а, V, в - коэффициенты температуропроводности, кинематической вязкости и температурного расширения жидкости; g - ускорение свободного падения.

В случае кондуктивно-ламинарного режима система (1) - (8) декомпозируется на две последовательно решаемые задачи: тепловую

дТ 1 д Т

(9)

д0 Pr CX 2 ’

Т(X,G)_ G, Т(G,0)_ 1, Т(1,0) _ G; (1G)

и гидродинамическую

д Г д2Ф д2Ф^ д4Ф

30І3Х2 CY2 J dX4 dX2dY2

I д4Ф дТ (X,0)

+~37T дХ ’

Ф(X,Y,G )_ G, Ф^^^) _ Ф(1^,0) _ _Ф(Х, ,0)_Ф( X,1,0) _G ; дФ(G,Y,0) _ дФ(1Х0) _

(11)

(12)

(13)

CX CX

дФ(ХД0) _ дФ(Х,1,0)

CY

CY

_ G.

(14)

где Ф _ Т/ Gr .

Решение. Решение тепловой задачи (9), (1G) известно [9]

2 v(-l)P sln

T (X, 0)_ 1 - X

-sin

[(1 - X )np ]>

x exp

2 2 П p , Pr

(15)

и поэтому сосредоточим основное внимание на системе (11) - (14), принимая во внимание (15).

Применив одностороннее интегральное преобразование Лапласа [10] по переменной 0 с символьным обозначением оператора 2 к (11) - (14), получим

3Ф, 3Ф,

S I S --

3X2 CY2

3 4Ф

3Y4

IZ

CX4 3X23Y2

дТ (X, 0)"

дХ

Ф- (G, Y, s )_Ф- (1, Y, s )_ _Ф- (X,G,s) _ Ф- (X,1,s)_ G ; дФL (G,Y,s) _ дФL (l,Y,s) _

дХ _ дХ _

дФL (X,G,s) _ дФL ^JX,1, s)

(16)

(17)

CY

CY

_ G.

(1В)

где 8, Фь - изображения 0 и Ф.

Конечное интегральное синус-преобразование [11] (16) - (18) по переменной Х с символьным оператором РХ имеет вид:

—sX Ф -X I S

d 2Ф— dY2

32ФL (G,Y,s)

X

d 4Ф

-_-X

32ФL (1,Y, s)

CX2

cos X-

IX4Ф„ -2X

dY2

dY

f-1F X Z

3T (X, 0)

CX

Фlx (М_Ф-х (1,s) _ G;

dФLX (G,s) _ dФLX (1,s)

dY

dY

_G,

(19)

(2G)

(21)

где ФLX - изображение ФL, к - характеристическое число, определяемое из уравнения sin X _ G.

Конечное интегральное синус-преобразование (19) - (21) по переменной Y с символьным оператором FY таково

-sX2Ф,^ -sp2Ф,^ _ -XE(s)(cosX - 1)i

IX4ФLXY 12X2ц2ФLXY -рiF(s)(cosц- 1)i

3T (X, 0)"

+l4фLXY I FYFXZ

CX

(22)

где ФLXY - изображение Ф1Х, ц - характеристическое число, определяемое из уравнения sin ц = 0 ;

Е (5 ) = Fy

[32 Ф L (1, Y, s ) _ F ~3 2Ф(G, Y, s )

CX2 Y CX2

F (s )_

32ФLX (1,s) _ 32ФLX (G,s)

CY CY

Из (22) следует решение в изображениях

ФLXY _{XE(s)(cosX-l) + pF(s)(cosц-l)-

+Fr F X Z

~ дТ (X, 0) 1/

CX і

(2 +i2 )2 + s (X2 + ц2)

. (23)

Для того, что бы получаемое решение удовлетворяло начальному условию (12), положим

Z

Z

1 [E[1 -exp(-a0)]; (24) F(s)]_xFі0^[1 -exp(-a0)] , (25)

X + ц

где Z 1 - обратный оператор Лапласа; Е(9), F (9)

- ограниченные функции в области решения; a -величина, подлежащая определению.

Оригинал выражения (23), имея в виду (24) и (25), есть

(2os Х2)) Е №- exp (-a9)]+

(Х2+ц2)

Ц 2 ) F (9)[! - exP (-a9)] +

(Х2 + ц2)

~dT (X, Л)]

ФXY _-

дХ

Т.к.

xexp[-(X2 +ц2)0-Л)]dЛ . ^ = -! - 2Z (-!)' cos [(1 — X )„p]!

(26)

p_I

П2p2 Л

<exp| -^тЛ

fX {cos [(! - X )V )] ,

FX [-1] = Х(C0sХ -1) , Fy [-1] = -1- (C0s Ц - 1)

x ц

то решение в изображениях будет

Ф XY =X(^C20S Х2)2) ) (9)[1 - eXP (-a9)] +

(Х2 +ц2)

+ц(ц^ F (9)[1 - eXP (-a9)]+

(Х2 +ц2)

+2Х(-1)

p _i

(cos X - і) (cos ц- і)

(X2 + ц2) Xi x

<{і - exp [-(X2 +ц2) 0]} +

X [cos X - cos (np)] (cos ц-1)

n2p2 -X2 (X2 +ц2)

-n2 p2/Pr + X2 +ц2 exp(-^Pp-) - exp[-(X2 +ц2 ЦЛ . (27)

Применяя формулу обращения синус-

преобразования [11], имеем

да да

Ф(Х,Y,0)_ 4 £ £ (x. (cosXn - і) (0)

n_I,3,5... m_I,3,5... x [1 - exp ( anm0)] + И"m (C0S И"m - Fm (0) x

x[i—exp (—*„ 0Д—(cos Xn "X^005 ц-—1) x

x{I -exp[-(x. +ц. )0]} +

Xn [C0S Xn - C0S (p)] (cos цm - 1)

+2£(-1)

p _I

П2 p2 -X n2

(x n +ц. )

-n2 pVPr +x. +ц.

exp

2 2 П p

Pr

-ехР[-{Х2+*)0]^!!п<^(Ж), (28)

(Xп + ут )

где Хп =пп ; цт =пт , а в силу одинаковой скорости убывания слагаемых в (28) из физических соображений

апт = Х П + Ут ■

Нетрудно проверить, что решение (28) удовлетворяет условиям (12) и (13). В силу несимметричности решения (28) относительно характеристических чисел лп и мп коэффициенты Ёп и Рт не равны между собой и могут быть найдены из систем уравнений, полученных подстановкой (28) в граничные условия (14):

X ([лп (сов л п - 1)Ёп (и) + мт (0в мт - 1)^т (и) -

(cos X„ - 1)(cos ц. - 1)

Xn ц.

+2£(-1)

p _i

Iі - exp [-( +ц. )0]} +

Xn [C0S Xn - C0S (np)] (cos цm - 1) .

П2 p2 -Xn2 ц.

(x n +ц. )

-П2 p2/Pr+Xn +ц

-exp [-(X n +ц. )0]})

exp

2 2 П p

Pr

(x n+ц m)

_ G.

(29)

m

и

где m = 1,3,5,...;

Е \[л» (C0S л» - 1)E%n (и)+ Mm (C0S Mm - 0) (и) -

(C0s Х„ - - !)

Х (Цш

+2Е(-1)

Р=1

I1 - exp [-(х (+ц m )9]}+

Хя [C0s Хя - C0s (ПР)] (c0S цm - 1)

П Р2 -Х(2 Цm

(х(+цm)

-л2 р2/Рг+Х( +ц

- exp [-(х(+цm )9])^

exp

2 2 П Р

Pr

Ц»

• = 0.

(30)

(х (+ц m)

где ( = 1,3,5,__

Для того, что бы решение (28) удовлетворяло исходному уравнению (11) необходимо после его подстановки выполнение (11) тождественно, что предполагает наличие Е( (9) и Fm (9) с принадлежностью к С1 [0, да). Такому свойству отвечают, например, полиномиальные представления.

Анализ решения. Заметим, что в (28) чётные слагаемые рядов равны нулю, что очевидным образом распространяется и на систему (29)-(30), решение которой осуществлено итерационным методом неполной релаксации [12]. Кроме того, при проведении вычислений с использованием (28)-(30) необходимо учитывать, что

C0S (nn )-C0S (пр) lim----^= 0 .

р^( п2 (р2 - (2 )

В качестве примера результаты расчётов приведены в табл. для Рг=7 при 0 = 0,15.

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

Эволюция функции тока для срединного сече-

Рис. 2 Поле температур для Рг=0,7 в различные моменты времени и: 1 - 0,001; 2 - 0,01; 3 - 0,1; 4 - 0,5

ния каверны при Рг=0,7 и Рг=7 представлена на рис. 4 и 5. Анализ полученных результатов показывает, что в связи с очень большой величиной начального

теплового импульса вблизи стенки (рис. 6, 7) возникает волновое движение жидкости с инверсией (обращением на противоположное направление) первоначально возникшего течения с нисходящим движением у горячей стенки. При этом величина функции тока неинверсированного течения тем больше, чем больше число Рг.

n Е( (0,15) m Fm 20,15)

1 5,8399-10-3 1 -0,045335

3 -3,8923-10-4 3 -6,1799-10-3

5 2,0071-10-3 5 -2,7301-10-3

7 1,713310-3 7 -1,7924-10-3

9 1,3718-10-3 9 -1,3635-10-3

11 1,1286-10-3 11 -1,1099-10-3

13 9,5524-Ю-4 13 -9,3882-10-4

15 8,2698-Ю-4 15 -8Д440-10-4

17 7,2878-Ю-4 17 -7,1942-10-4

19 6,5134-Ю-4 19 -6,4437-10-4

21 5,8876-Ю-4 21 -5,8352-10-4

23 5,3717-Ю-4 23 -5,3316-10-4

25 4,9390-Ю-4 25 -4,9079-10-4

27 4,5711-Ю-4 27 -4,5464-10-4

29 4,2542-Ю-4 29 -4,2344-10-4

31 3,9786-Ю-4 31 -3,9623-10-4

33 3,7366-Ю-4 33 -3,7230-10--4

35 3,5225-Ю-4 35 -3,5109-10-4

37 3,3316-Ю-4 37 -3,3216-10-4

39 3Д604-10-4 39 -3,1516-10-4

41 3,0059-Ю-4 41 -2,9981-10-4

43 2,8659-Ю-4 43 -2,8589-10-4

45 2,7384-Ю-4 45 -2,7320-10-4

47 2,6218-Ю-4 47 -2,6158-10-4

49 2,5147-Ю-4 49 -2,5091-10-4

51 2,4160-Ю-4 51 -2,4108-10-4

53 2,3248-Ю-4 53 -2,3198-10-4

55 2,2403-Ю-4 55 -2,2355-10-4

57 2Д616-10-4 57 -2,1570-10-4

59 2,0884-Ю-4 59 -2,0839-10-4

Рис. 3 Поле температур для Рг=7 в различные моменты времени и: 1 - 0,01; 2 - 0,1; 3 - 1; 4 - 2

С течением времени по мере прогревания жидкости вблизи вертикальных стенок формируются циркуляционные контуры с восходящим движением жидкости у горячей стенки и нисходящим - у хо-

4

Рис. 4 Эволюция функции тока при Рг=0,7: 1- и=1 • 10-4; 2- и=0,001;3- и=0,005; 4- и=0,01; 5- и=0,02;6- и=0,05; 7- и=0,1; 8- и=0,5

Рис. 5 Эволюция функции тока при Рг=7: 1- и=0,001; 2- и=0,01; 3- и=0,1;4- и=0,15; 5- и=0,2; 6- и=0,5; 7- и=2

Рис. 6 Поле градиента температур для Рг=0,7 в различные моменты времени и: 1 - 0,001; 2 - 0,01;

3 - 0,1; 4 - 0,5

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

0.5 Х 1

Рис. 7 Поле градиента температур для Рг=7 в различные моменты времени и:1 - 0,01; 2 - 0,1; 3 - 1; 4 - 2

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

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

Рис. 8 Динамика изменения нормированных среднеинтегральных значений функциитока (------) и

температур (..) при Рг=0,7 (1);Рг=1 (2); Рг=7 (3)

и температуры не эквивалентны при Рг=1, в отличие от уже известных решений аналогичных задач. Более детальный анализ динамики изменения нормированных среднеинтегральных значений функций тока и температуры при Рг=7 в области малых значений безразмерного времени 0 показывает (рис. 8), как и ожидалось, преимущественный рост температуры по сравнению с функцией тока, при 0» 0,73 их значения совпадают, а затем уже преимущественный рост над температурой имеет функция тока. Это означает, что первопричиной конвективного

движения жидкости является тепловая неравновес-ность системы. Такой результат находится в согласии с общепринятыми физическими представлениями о естественной конвекции [3].

Важным заключительным моментом анализа является вопрос удовлетворения найденного решения (28) исходному уравнению (11), что было проверено прямой подстановкой. Так для центральной точки с координатами X=Y=0,5 для различных 0 уравнение (11) выполнялось при (=m=150 с относительной погрешностью менее чем 0,04, что вполне приемлемо при выполнении практических оценок эффективности и точности различных вычислительных процедур.

Литература

1. O0sthuizen P.H. An Intr0duCti0n to C0nveCtive Heat Transfer Analysis - Singap0re: WCB/MCGrow-Hill,

1999. - 620 pp.

2. Latif M. Jiji. Heat C0nveCti0n. - Springer, 2006. - 443 pp.

3. Гебхарт Б., Джалурия Й., Махаджан Р., Самакия Б. Свободноконвективные течения, тепло- и массооб-мен. В 2-х кн., кн. 1. - М.: Мир, 1991. - 678 с.

4. Ландау Л.Д., Лифшиц Е.М. Теоретическая физика: В 10 т. T.VI. Гидродинамика. - М.: Наука 1988.-736 с.

Воронежская государственная технологическая академия Воронежский государственный технический университет

5. Коган В.Б. Теоретические основы типовых процессов химических технологий. - Л.: Химия, 1977. - 589 с.

6. SChetz J. A., Eighh0rn R. Unsteady natural C0nveCti0n in the viCinity 0f a d0ubly infinite vertiCal plate // TransaC-ti0ns 0f the ASME J. 0f Heat Transfer. - 1962. - № П.-p. 334-338.

7. Ryazhskikh V.I., Slyusarev M.I., B0ger A.A., Ryab0v

S.V. DynamiCs 0f laminar free-C0nveCtive fl0w 0f a new-t0nian fluid between vertiCal plane is0thermal walls // J0urnal 0f Engineering PhysiCs and Therm0physiCs.-2009.-V. 82.-№ 6.-Pp. 1163-1170.

8. Богер А.А., Рябов С.В., Ряжских В.И., Слюсарев М.И. Расчёт кондуктивно-ламинарного движения термоконвекции ньютоновской среды в прямоугольной каверне с вертикальными изотермическими стенками // Изв. РАН Механика жидкости и газа.- 2010.- № 3. С. 17-21.

9. Лыков А.В. Тепломассообмен: (Справочник). - М.: Энергия, 1978. - 480 с.

10. Дёч Г. Руководство к практическому применению преобразования Лапласа и Z - преобразования. - М.: Физматгиз, 1971. - 288 с.

11. Корн Г., Корн Т. Справочник по математике для научных работников и инженеров. - М.: Наука, 1973. -832 с.

12. Фаддеев Д. К., Фаддеева В. Н. Вычислительные методы линейной алгебры. - М. - Л.: Гл. изд-во физ.-мат. лит-ры, 1963. -734 с.

CONDUCTIVE LAMINAR FREE CONVECTION LIQUID FLOW DYNAMICS IN SQUARE REGION M.I. Slyusarev, E.D. Chertov, V.I. Ryazhskih, V.G. Stogney, A.A. Boger

Using the finite integral Fourier transform and one-sided integral Laplace transform, an analytical solution of the linearized Oberbeck-Boussinesq equations in Helmholtz’s variables for square region with different values of constant temperatures vertical walls and no heat flow through the upper and lower base is received

Key words: non-staionary free convection, analytical solution, cquare redion