Теплообмен при ламинарном течении жидкости в круглой трубе

Полный текст

ВВЕДЕНИЕ

Изучение процессов теплообмена при ламинарном режиме течения среды в каналах посвящено значительное число научных работ, например [1, 2, 9–14]. В этой области достигнуты существенные теоретические результаты, имеющие важное значение для практики конструирования современных высокоэффективных теплообменных устройств.

В указанных публикациях сравнительно подробно на математической основе рассмотрены случаи ламинарного движения жидкости в трубах различной геометрической конфигурации при различных краевых условиях. В монографиях [2, 4] для решения сформулированных задач были использованы специальные вырожденные гипергеометрические функции [3, 5, 15]. На основе этих функций были выведены аналитические зависимости для определения собственных функций и собственных значений поставленных задач.

Однако полученные аналитически строгие расчетные зависимости в математическом отношении являются сравнительно громоздкими при проведении по ним конкретных инженерных вычислительных операций. В связи с этим целесообразно к методам расчетов, рекомендуемым в [2–4], разработать более упрощенные способы проведения технических вычислений характеристик изучаемого вида теплообмена.

ПОСТАНОВКА ЗАДАЧИ

С учетом принятых выше допущений математическая постановка рассматриваемой задачи в безразмерной форме имеет вид

$\begin{array}{l}\frac{{\partial }^2\vartheta (X,R)}{\partial R^2}+\frac{1}{R}\frac{\partial \vartheta (X,R)}{\partial R}=(1-R^2)\frac{\partial \vartheta (X,R)}{\partial X},\\ 0\le X\le \infty ;0\le R\le 1\end{array}$, (1)

$\frac{\partial \vartheta \left(X\mathrm{,0}\right)}{\partial R}=0$, (2)

$\frac{\partial \vartheta \left(X\mathrm{,1}\right)}{\partial \gamma }=-Bi\vartheta \left(X\mathrm{,1}\right)$, (3)

$\vartheta \left(\mathrm{0,}R\right)=1$, (4)

где ϑ(X, R) – безразмерная температура жидкости (искомое температурное поле),

$\vartheta (X,R)\Rightarrow 0$ при $x\Rightarrow \infty$,

где X, R – соответственно осевая и радиальная координаты канала.

Bi – число Био, характеризующее интенсивность теплообмена на наружной поверхности канала.

РЕШЕНИЕ

В монографии [1] изложены математические методы расчета теплообмена при ламинарном течении несжимаемой жидкости в каналах разной геометрической формы при граничных условиях первого, второго и третьего рода. В работах [2–4] представлены результаты исследований, дополняющие материалы, рассмотренные в книге [1].

Так, в частности, опубликованы в [2] аналитические решения задач теплообмена в ламинарных потоках жидкости при граничных условиях третьего рода на внешних поверхностях плоских и круглых труб. Для расчета температурного поля в движущейся жидкости использована безразмерная зависимость вида

$\vartheta (X,R)=\sum _{n=1}^{\infty }{A}_n{\psi }_n\left(R\right)\exp (-{\mu }_n^{2}X)$, (5)

где ψ_n(R) и μ_n – собственные функции и числа рассматриваемой задачи. Для их определения в случае круглого канала необходимо провести исследование следующей задачи

$\psi ´´+\frac{\psi ´}{R}+{\mu }^2(1-R^2)\psi =0$, (6)

$\psi ´=0$ при R = 0, (7)

$\psi ´=-Bi\psi$ при R = 1. (8)

Здесь использована безразмерная форма записи, как более удобная с математической точки зрения.

Представить решение данной задачи через элементарные функции в общем случае не удается, поэтому необходимо использовать специальные функции. Тогда можно решение задачи (6)–(8) записать в виде [5]

$\psi =\exp (-\mu \frac{R^2}{2})F_a(\alpha ,\gamma ,\mu R^2)$, (9)

где F_a(a, g, mR²) – конфлюэнтная гипергеометрическая функция, определяемая как бесконечная сумма [6]

$F_a(\alpha ,\gamma ,\mu R^2)=1+\frac{\alpha }{\gamma }\mu R^2+\frac{\alpha (\alpha +1)}{\gamma (\gamma +1)}\frac{{\mu }^2{R}^4}{2!}+\frac{\alpha (\alpha +1)(\alpha +2)}{\gamma (\gamma +1)(\gamma +2)}\frac{{\mu }^3{R}^6}{3!}+\mathrm{...}$ (10)

Для круглой трубы [5]

$\alpha =\frac{1}{2}-\frac{\mu }{4}$ и $\gamma =1$.

Тогда окончательный вид соотношения (9) будет следующим

$\psi =\exp \left(-\mu \frac{R^2}{2}\right)\left[\begin{array}{l}1+\frac{(2-\mu )\mu }{4}{R}^2+\frac{(2-\mu )(6-\mu ){\mu }^2}{16}\frac{R^4}{{(2!)}^2}+\frac{(2-\mu )(6-\mu )(10-\mu ){\mu }^3}{64}\frac{R^6}{{(3!)}^2}+\mathrm{...}\\ +\frac{(2-\mu )(6-\mu )(10-\mu )\mathrm{...}(4m-2-\mu ){\mu }^m}{4^m}\frac{R^{2m}}{{(m!)}^2}+\mathrm{...}\end{array}\right]$ (11)

где m = 1, 2, 3, …

Произведя подстановку (11) в граничное условие (8), получим уравнение для нахождения характеристических чисел μ_n

$\begin{array}{l}1-\left[\begin{array}{l}\frac{2-\mu }{2}+\frac{(2-\mu )(6-\mu )\mu }{4{(2!)}^2}+\frac{(2-\mu )(6-\mu )(10-\mu ){\mu }^2}{4^3{(3!)}^2}+\mathrm{...}\\ +\frac{2m(2-\mu )(6-\mu )(10-\mu )\mathrm{...}(4m-2-\mu ){\mu }^{m-1}}{4^m{(m!)}^2}+\mathrm{...}\end{array}\right]\times \\ {\left[\begin{array}{l}1+\frac{(2-\mu )\mu }{4}+\frac{(2-\mu )(6-\mu ){\mu }^2}{16{(2!)}^2}+\frac{(2-\mu )(6-\mu )(10-\mu ){\mu }^3}{64{(3!)}^2}+\mathrm{...}\\ +\frac{(2-\mu )(6-\mu )(10-\mu )\mathrm{...}(4m-2-\mu ){\mu }^m}{4^m{(m!)}^2}+\mathrm{...}\end{array}\right]}^{-1}=Bi\end{array}$ (12)

Формулы (11) и (12) в общем случае являются сложными для выполнения по ним инженерных расчетов. Однако они позволяют сравнительно просто определить некоторые важные параметры процесса в промежуточных состояниях. Так, из выражения (11) вытекает, что при Bi = 2 первый корень характеристического уравнения (11) равен μ₁ = 2. Далее следует, что μ₂ = 6 при Bi = 3.6 и μ₃ = 10 при Bi = 4.839. Эти значения можно рассматривать как эталонные величины и на их основе установить по две зоны, в которых находятся искомые корни уравнения (11). Очевидно, что

0 ≤ μ₁ ≤2, если 0 ≤ Bi ≤ 2,

2 ≤ μ₁ ≤ 2.7044, если 2 ≤ Bi < ∞,

5.0675 ≤ μ₂ ≤ 6, если 0 ≤ Bi ≤ 3.6,

6 ≤ μ₂ ≤ 6.6790, если 3.6 ≤ Bi ≤ ∞,

9.1576 ≤ μ₃ ≤ 10, если 0 ≤ Bi ≤ 4.839,

10≤ μ₃ ≤ 10.6734, если 4.839 ≤ Bi < ∞.

Подобные границы могут быть достаточно легко найдены и для последующих корней μ_n.

Также несложно установить некоторые собственные функции на основе уравнения (11). В частности,

${\psi }_1\left(R\right)=\exp (-R^2)$ для Bi = 2, (13)

${\psi }_2\left(R\right)=(1-6R^2)\exp (-3R^2)$ для Bi = 3.6, (14)

${\psi }_3\left(R\right)=(1-20R^2+50R^4)\exp (-5R^2)$, если Bi = 4.839. (15)

Графическое представление собственных функций (13)–(15) представлено ниже.

Коэффициент A_n решения (5) определяется на основе выражения [2]

 

Рис. 1. График собственной функции для уравнения y₁(R) = exp(–R²), если Bi = 2

 

Рис. 2. График собственной функции для уравнения y₂(R) = (1 – 6R²)exp(–3R²), если Bi = 3.6

 

Рис. 3. График собственной функции для уравнения y₃(R) = (1 – 20R² + 50R⁴)exp(–5R²), если Bi = 4.839

 

$A_n=\frac{{\displaystyle {\int }_0^{1}R(1-R^2){\psi }_n\left(R\right)dR}}{{\displaystyle {\int }_0^{1}R(1-R^2){\psi }_n^2\left(R\right)dR}}$. (16)

Используя зависимости (13)–(16) можно вычислить сравнительно просто некоторые контрольные A_n. Например, если Bi = 2, то коэффициент A₁ будет равен

$A_n=\frac{{\displaystyle {\int }_0^{1}R(1-R^2){\psi }_n\left(R\right)dR}}{{\displaystyle {\int }_0^{1}R(1-R^2){\psi }_n^2\left(R\right)dR}}=\frac{2R^2\exp (-R^2)\left|{}_0^1\right.}{\left[R^2\exp (-2R^2)-\frac{1}{2}\exp (-2R^2)\right]\left|{}_0^1\right.}=\frac{4\exp (-1)}{\exp (-2)+1}=1.2961$. (17)

Это значение полностью соответствует величине, приведенной в работе [2], в которой использовался численный метод расчета первых трех корней μ_n уравнения (12), собственных функций ψ_n и коэффициентов A_n.

Таким образом можно указать границы для A₁, а конкретно при 0 ≤ Bi ≤ 2, 1 ≤A₁ ≤ 1.2961, а если 2 ≤ Bi ≤ ∞, то будет 1.2961 ≤ A₁ ≤ 1.4764.

Итак, строгое аналитическое решение типа (11) позволяет достаточно просто определить важные локальные характеристики изучаемого процесса.

Однако для инженерной практики целесообразно в дополнение к расчетным зависимостям (11)–(12), (16) получить упрощенные математические соотношения. Это особенно актуально при изучении процесса теплообмена в зоне канала, следующей за начальным термическим участком, т. е. когда осевая координата X превышает значение 0.1. Тогда в выражении (5) можно ограничиться первым слагаемым, так как остальные становятся пренебрежимо малыми. Исходя из этого, авторами разработана следующая приближенная математическая формула для нахождения упрощенным способом собственной функции ψ₁(R)

${\psi }_1\left(R\right)=1-\frac{{\mu }^2}{4}{R}^2(1-\frac{R^2}{4})+\frac{{\mu }^4}{64}{R}^4(1-\frac{5}{9}{R}^2+\frac{R^4}{16})-\frac{{\mu }^6}{64\cdot 36}{R}^6(1-\frac{7}{8}{R}^2+\frac{89}{400}{R}^4-\frac{R^6}{64})$. (18)

Рекомендуемое выражение (18) удовлетворяет с высокой степенью точности исходному дифференциальному уравнению (6) и полностью соответствует условию симметрии (7). Подставляя (18) в граничное условие (8), удается получить характеристическое уравнение для определения первого корня μ₁. Необходимо отметить, что учет первых трех членов в правой части предлагаемого решения (18) позволяет определить функцию ψ₁(R) с некоторым незначительным превышением, а расчет по полному выражению (18) характеризуется весьма малым занижением ψ₁(R).

Для подтверждения такой особенности решения (18) приведем расчет для варианта μ₁ = 2, что соответствует величине Bi = 2.

В первом случае согласно (18) получим

${\psi }_{1\max }\left(1\right)=1-\frac{1}{4}(1-\frac{1}{2})+\frac{16}{64}(1-\frac{5}{9}+\frac{1}{16})=0.376736$,

а во втором

${\psi }_{1\min }\left(1\right)=\mathrm{0,376736}-\frac{1}{36}(1-\frac{7}{8}+\frac{89}{400}-\frac{1}{64})=0.367517$.

В работе [2] для рассматриваемого варианта приводится табличное значение ψ₁(1) = 0.3679, рассчитанное численным методом. Сравнение этой величины с ψ_1max(1) и ψ_2min(1) свидетельствует, что невязка невелика.

Нужно также сказать, что с уменьшением параметра Bi и радиальной координаты R точность вычислений по формуле (18) повышается.

При проведении инженерных расчетов наибольший интерес представляет температура поверхности канала. Из зависимости (18) следует, что при R = 1 функция ψ₁(1) может быть записана в виде

${\psi }_1\left(1\right)=1-\frac{3}{16}{\mu }_1^2+\frac{73}{144}\frac{{\mu }_1^4}{64}-\frac{59}{6400}\frac{{\mu }_1^6}{64}$. (19)

Определим с помощью (19) ψ₁(1) для варианта Bi = 1 (μ₁ = 1.6413)

${\psi }_1\left(1\right)=1-\frac{3}{16}{1.6413}^2+\frac{73}{144}\frac{{1.6413}^2}{64}-\frac{59}{6400}\frac{{1.6413}^6}{64}=1-0.5051+0.057482-0.002816=0.549566$.

Табличное значение [2]

${\psi }_1\left(1\right)=0.5497$.

Отсюда можно сделать вывод, что в области 0 ≤ Bi ≤ 1 рекомендуемое простое аналитическое решение (18) обладает достаточной точностью.

Подставляя в граничное условие (8) рекомендуемое решение (18), получим, учитывая соотношение (19), характеристическое уравнение для определения первого собственного значения μ₁

$\frac{{\mu }_1^2}{4}-\frac{7{\mu }_1^4}{6\cdot 64}+\frac{83{\mu }_1^6}{63\cdot 16\cdot 180}=Bi(1-\frac{3}{16}{\mu }_1^2+\frac{73}{144\cdot 64}{\mu }_1^4-\frac{59}{6400\cdot 64}{\mu }_1^6)$. (20)

Нетрудно показать, что зависимость (20) сводится к алгебраическому уравнению третьей степени, решение которого может быть найдено, например, с помощью применения формулы Кардано [7].

При умеренных величинах параметра Bi (Bi ≤ 1) в выражении (20) допустимо не учитывать комплексы, содержащие μ₁⁶. Тогда из (20) несложно получить модификацию алгебраического уравнения второй степени, решение которого общеизвестно [7].

ЗАКЛЮЧЕНИЕ

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

Об авторах

Ю. В. Видин

Федеральное государственное автономное образовательное учреждение высшего образования «Сибирский федеральный университет»

Email: roman.v.kazakov@gmail.com
Россия, Красноярск

Р. В. Казаков

Федеральное государственное автономное образовательное учреждение высшего образования «Сибирский федеральный университет»

Автор, ответственный за переписку.
Email: roman.v.kazakov@gmail.com
Россия, Красноярск

Список литературы

Дополнительные файлы