Известия РАН. Энергетика, 2020, № 6, стр. 76-80

В том случае, когда теплофизические коэффициенты переноса материала имеют ярко выраженную зависимость от температуры, и процесс нагрева (или охлаждения) тела осуществляется в широком интервале ее изменения, приходится решать нелинейное дифференциальное уравнение нестационарной теплопроводности [1–3]. Аналитическое интегрирование такого уравнения совместно с заданными краевыми условиями, как правило, сопряжено со значительными математическими трудностями. Получить строгое, в теоретическом смысле, решение таких задач удается только в исключительных случаях. Теоретическому исследованию подобных проблем посвящено сравнительно много работ, систематизированный обзор которых приведен в статье А.В. Лыкова [3]. В ней дан всесторонний анализ приближенных методов, применяющихся в современной инженерной практике. Существующим способам свойственны недостатки, главными из которых являются их математическая сложность и отсутствие достаточных доказательств точности предлагаемых решений поставленных задач.

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

Проиллюстрируем рекомендуемый аналитический способ на примере конкретной теплофизической задачи, записанной для удобства в следующей безразмерной форме

Параметр $a,$ характеризующий динамику изменения коэффициента теплопроводности изделия от температуры может быть как положительным ($a > 0$), так и отрицательным ($a < 0$), что зависит от вида материала [1].

Очевидно, что в частном случае, а именно, когда $a = 0,$ задача (1)–(4) оказывается линейной, решение которой известно [2]

Бесконечный ряд в формуле (5) является быстросходящимся и, начиная с некоторого значения числа Фурье (обычно ${\text{Fo}} \geqslant {\text{0}}{\text{.3}}$), можно ограничиться одним первым слагаемым, т.е. тогда имеем

При малых числах ${\text{Fo}}$ (${\text{Fo}} \leqslant {\text{0}}{\text{.2}}$), т.е. для начальной стадии процесса, целесообразно вместо (5) применить более простое выражение [2]

где

Специальная функция ${\text{erf}}Z$ (интегральная ошибка) достаточно подробно исследована авторами книги [4], в которой приведены ее табличные числовые значения с шагом 0.01 до предельной величины $Z = 3.$ Наличие данной таблицы существенно облегчает выполнение вычислений по соотношению (7). При больших значениях $Z$($Z > 2$) можно использовать асимптотическую формулу [5]

Таким образом, на основе (6) и (7) особенно легко находится безразмерная температура в центре пластины ($X = 0$) соответственно по элементарным выражениям

1) если ${\text{Fo}} \leqslant {\text{0}}{\text{.2}}$

2) если ${\text{Fo}} > {\text{0}}{\text{.2}}$

Нетрудно показать, что расчеты по формулам (6), (7), (10) и (11) будут давать завышенные результаты в случае $a > 0$ и наоборот заниженные при $a < 0.$ Для нахождения противоположной границы для искомого поля температуры $\vartheta \left( {X,{\text{Fo}}} \right)$ целесообразно ввести новую зависимую переменную

Тогда задача (1)–(4) преобразуется к виду

Решая квадратное алгебраическое уравнение (12), получим соотношение

в котором при $a > 0$ в правой части следует брать знак “+”, а при $a < 0$ – знак “–”. Если в левой части дифференциального уравнения (13) заменить сомножитель $\frac{1}{{1 + a\vartheta }}$ на комплекс $\frac{1}{{1 + a}},$ то аналитическое решение системы (13)–(16) можно представить в виде которое при повышенных значениях ${\text{Fo}}$ становится проще

При малых величинах ${\text{Fo}}$ предпочтительнее использовать выражение

Таким образом, границами для действительного температурного поля $\vartheta \left( {X,{\text{Fo}}} \right)$ удовлетворяющего исходной задаче (1)–(4) служат $\vartheta {\text{*}}\left( {X,{\text{Fo}}} \right)$ и $\vartheta {\text{**}}\left( {X,{\text{Fo}}} \right),$ которые рассчитываются по формуле (6) (либо (7)) или по выражению (17), с предварительным определением $\theta \left( {X,{\text{Fo}}} \right)$ на основе (19) (или (20)).

В таблице 1 приведены данные расчетов $\vartheta {\text{*}}\left( {X,{\text{Fo}}} \right)$ и $\vartheta {\text{**}}\left( {X,{\text{Fo}}} \right)$ для центра пластины ($X = 0$) и плоскости с координатой $X = 0.5$ при коэффициенте $a = 0.2.$ Здесь же указаны результаты, полученные численным методом [6]. Таблица 2 содержит подобные вычисления, проведенные для варианта, когда $a = - 0.2.$

Принятые величины параметра $a$ соответствуют характерным случаям весьма значительной зависимости коэффициента теплопроводности от температуры [7]. Из анализа названных таблиц видно, что предлагаемый способ позволяет получать простыми математическими действиями сравнительно узкую зону, в которой располагается искомое поле температуры. С уменьшением параметра $a$ область между верхней и нижней граничными кривыми очевидно будет уменьшаться.

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