Решение линейного дифференциального уравнения теплопро-

водности для одномерного нестационарного процесса

Рассмотрим численный метод на примере расчета одномерного темпера-турного поля в стержне. Начальное распределение температуры T ₀ (x) задано численной таблицей или графиком. Боковая поверхность стержня теплоизоли-рована Требуется определить температурное поле в процессе выравнивания температуры в стержне. Математическое описание задачи сформулируем так:

1) ; 0 < x< l; t> 0 (24)

2) T (x, 0 ) = T ₀ (x); (25)

3) ; (26)

4) T (l, t) = T ₀; (27)

Рис. 6. К пояснению численного метода конечных разностей

При решении задачи численным методом дифференциальное уравнение теплопроводности и уравнения краевых условий заменяют уравнениями в ко-нечных разностях – приближенными выражениями, связывающими не беско-нечно малые, а небольшие конечные приращения температуры Δ Т, времени Δ t и координаты Δ x. Весь стержень разбивают на n слоев одинаковой толщины Δ x=l/ n. Время процесса также делят на равные промежутки длительностью Δ t. Слои нумеруются индексами i, а промежутки времени – индексами j. При этом делают следующие допущения. Температуру в пределах одного слоя считают постоянной и равной температуре середины слоя, т. е. непрерывное распреде-ление температуры по стержню заменяют ступенчатым. Температуру каждого слоя в пределах одного промежутка времени считают постоянной, т. е. непре-рывное изменение температуры во времени заменяют скачкообразным изме-нением при переходе от одного промежутка времени к другому.

Обозначим температуру слоя i в промежутке времени j через T _{i, j}. Тогда температура соседних с ним слоев в том же промежутке времени будет T _{i-1, j}и T _{i+1, j}, а температура слоя i в последующем промежутке времени будет T _{i, j+1}.

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

. (28)

Вторая производная температуры также заменяется соответствую-щим ей выражением в конечных разностях. Сперва выражается в конечных раз-ностях градиент в двух смежных сечениях слева и справа от расчетного слоя i:

;

.

Затем вторая производная заменяется средней скоростью изменения гра-диента температуры на отрезке Δ х:

(29)

Подставляя приближенные выражения (28) и (29) для производных тем-пературы в уравнение (24), получаем систему конечно–разностных уравнений:

(30)

Следовательно, расчет выравнивания температуры в стержне сводится к решению системы уравнений.

Однако, из этой системы уравнений можно в явном виде выразить неизвестную температуру для последующего интервала времени T _{i, j+1}:

. (31)

Эта формула позволяет рассчитать температуру любого слоя, расположенного внутри стержня (i=2, 3 … n - 1), в последующий интервал времени (j+1) по зна-чениям температуры этого слоя и соседних с ним слоев в предшествующий интервал времени j.

Для расчета температуры граничных слоев стержня (i = 1 и i = n) необхо-димо учесть граничные условия. Так, при расчете температуры граничного слоя j = 1 в соответствии с адиабатическим условием (26) следует считать темпера-туру слева от него всегда равной температуре этого слоя. Тогда уравнение (31) применительно к слою j = 1 будет иметь вид:

. (32)

При расчете температуры граничного слоя i = n в соответствии с изотер-мическим условием (27) следует считать температуру справа от него всегда равной Т ₀. Тогда уравнение (31) применительно к слою i = n будет иметь вид:

. (33)

Таким образом, интегрирование дифференциального уравнения (24) сво-дится к последовательному расчету по формуле приближенных значений тем-пературы T _{i, j+1} слоев стержня в последующий интервал времени по распреде-лению температуры T _{i, j} в предшествующий интервал времени. Такой вариант метода конечных разностей, при котором расчет ведется по формулам, называют методом конечных разностей по явной схеме.

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

Условие устойчивости для уравнения (31):

. (34)

Если оно не соблюдается, то разность между решениями уравнений (24) и (31) возрастает по мере увеличения t.

Наиболее простая расчетная формула получается при выполнении ус-ловия:

. (35)

Она запишется так: , т. е. температура слоя i в последующий промежуток времени (j+1) есть средняя арифметическая температур двух соседних слоев (i-1) и (i+1) в предыдущий промежуток времениj.

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

Разработаны другие конечно-разностные уравнения, которые лишены этого недостатка.

Рис. 7. Узлы сеточной модели

На рис. 7 показаны узлы сеточной модели, которые используются в ко-нечно-разностных уравнениях. Рассмотренная явная схема метода конечных разностей построена по шаблону 1.

Для шаблона 2 уравнение в конечных разностях запишется:

(36)

Для шаблона 3: (37) Расчетную схему по шаблону 3 называют схемой Кранка-Николсона. Она име-ет широкое применение благодаря большей точности расчетов.

Расчетные схемы по шаблонам 2 и 3 называют неявными схемами, поскольку расчет ведется путем решения систем линейных алгебраических уравне-ний, а не по формулам. Схемы 2 и 3 безусловно устойчивы. Для них можно ис-пользовать любые шаги интегрирования по времени. Однако, при этом возни-кают проблемы точности.

При неявном методе требуется на каждом шаге по времени решать систему линейных алгебраических уравнений, что по сравнению с явным методом резко увеличивает число арифметических операций, а, следовательно, растут ошибки округления. Если Dt _неявн./ Dt _явн< 3, то количество алгебраических операций в неявном методе будет больше, чем в явном. Таким образом, неявный метод выгодно использовать, если шаг по времени более чем в 3 раза больше по сравнению с явным методом.

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

более выгодно использовать явную схему расчета.

Рассмотрим решение уравнения в конечных разностях для неявной схемы по шаблону 2.

После преобразований получим уравнение:

(38)

По этой формуле составим n-2 уравнений для внутренних узлов, температура которых в интервале времени j+1 неизвестна. В результате получим систему линейных алгебраических уравнений, которую удобно записать в матричной форме: [A]× {T}_j+1={B}, где [A] – матрица коэффициентов, {T}_j+1 – вектор-столбец температуры в узлах, {B} – вектор-столбец свободных членов.

Матрица [ A ] симметричная, ленточная, трехдиагональная. Только на трех диагоналях элементы матрицы отличны от нуля. Такая система уравнений легко решается итерационными методами.

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

(39)

Условие устойчивости в этом случае имеет вид:

(40)