Дульнев Г.Н., Парфенов В.Г., Сигалов А.В. Применение ЭВМ для решения задач теплообмена (Дульнев Г.Н., Парфенов В.Г., Сигалов А.В. Применение ЭВМ для решения задач теплообмена.djvu), страница 6

На основе принятой нумерации для системы (1.26) — (1.28) в подпрограмме ЬУВТТ производится формирование массивов А и Т, соответствующих матрице и столбцу свободных членов, а затем с помощью стандартной подпрограммы ОЕЕО (см. 2 1.3) решается система уравнений н находится столбец неизвестных. Отметим, что после решения системы найденные температуры размещаются в том же массиве Т, в котором находились свободные члены.

В уравнениях баланса в общем случае могут быть связаны температуры 7'и 7ь У, с любыми номерами 1, 1', 1, поэтому матрица А является матрицей общего вида и для решения системы используется подпрограмма ОЕЕО. Рассмотрим процедуру формирования матрицы А и столбца Т. Сначала двумерный массив А н одномерный массив Т обнуляются, а затем производится расчет их ненулевых элементов путем последовательного суммирования отдельных членов, входящих в формулы (1.26) — (1.28).

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

С целью экономии машинной памяти и сокращения объема исходных данных целесообразно приводить информацию только об отличных от нуля оц В приведенной программе эта информация задается с помощью двух одномерных массивов 1,) и $13. В массиве ц длиной 2 ь И1), где И1Л вЂ” число связей между телами, парами записаны номера взаимодействующих тел: элементы 1Л (2 ч М вЂ” 1), 1Л (2 ч М) указывают номера тел, участвующих в М-м взаимодействии. Массив 5!3 длиной И13 содержит значения тепловых проводимостей о ь соответствующих этим тепловым связям. Совершенно аналогичным образом представлена информация о тепловых связях от~ между телами и теплоносителями с помощью массивов П. и $11., а также о связях оД между телами и средами с помощью массивов 1К и Я К.

Движение потоков теплоносителей описывается с помощью массива М1- длиной 2 ь ИО и массива ОС длиной ИО, где ИΠ— число по- 19, 8 Т(ЬЬе()-т(ЬЬ«1)+СО ТСОНа) 198 19 СЯ)(!НЕЕ (ев с 2. Рхмци(е сис)хим УРмнщи)Н ие и САьь Оцл(Т,А,И,(,!.3-7,1н() и1 1Р((хн.не.е)72121 12,1и( И2 12 УОИМАТ(ОХ,'Ш( ',12) из с з. Оечдть РееудьтАТОЙ ИА РН1ИТ 13 И5 13 МНМАТ(/ОХ, 'ТЙИМРАТНЭ) ТЕА'! ИВ ВО Ы 1-(,Й! И7 14 РН1ИТ 15,1,Т(1) ИВ 15 УОИМАТ(22,13,812.3) ИЭ РН1ИТ 18 129 18 уснмАт(/Ох, 'твниРАт»РН тиюнсситеААЙ'/ !2! «1Х, 'НОМЕР О.ВХОДА О.ВЫХОДА') 122 ОО 17 )~1,ИЬ 123 Ы~Н«Ь«И! 124 17 ПЦИТ 18лиТ(ЬЬ),Т(Ы 1) 125 18 УОНМАТ(15,«о!2.3) 128 НЕП)НМ 127 ж(О Рис.

1.4 Продолжение токов теплоносителей из одного объема в другой. Массив МЬ содержит номера связанных объемов: из объема номер МЬ (2 * К вЂ” 1) теплоноситель втекает в объем с номером М1 (2:» К), К = 1, ..., МС«. В массиве С!С записаны произведения удельных теплоемкостей на массовые расходы с 6 ь Поток может вытекать либо из среды с заданной температурой 8» [6», в уравнении (1.28)1, либо из объема с неизвестной температурой Ь)' " (6,).

Чтобы в программе различить эти две разные ситуации, использован следующий прием: номера сред с заданной температурой в массивах МЬ и 1К начинаются не с единицы, а с 1)Ь + 1, где МЬ вЂ” число объемов с теплоносителями. Для формирования матрицы А и столбца свободных членов Т организуются пять циклов: по телам; по связям между телами (а(/); по связям между телами и теплоносителями (а() ); по связям между телами и средами (о(А) и по потокам теплоносителя между объемами (6 ь 6»,).

В каждом цикле последовательно выбирается связь, определяются номера взаимодействующих тел и объемов, находятся номера строк и столбцов, определяющие коэффициенты матрицы и столбца свободных членов, в которые вносит «вклад» данная тепловая связь, и, наконец, к ранее вычисленному значению коэффициента матрицы прибавляется «вклад» от рассматриваемого взаимодействия. Например, в цикле по связям между телами и объемами с тепло- носителями на основе массивов 1Ь и 81Ь проводится учет вкладов в матрицу А проводимостей ай, которые входят в ее первые У, рок [уравнения (1.26)! и в строки с номерами Ж, + 1, У, + 3, ...

А!,+ 2 э Ф вЂ” 1 !уравнения (1.27)). Номер тела 1 совпадает с номером температуры Т, в столбце неизвестных, а по номеру объема 1 определяются номера температур теплоносителя: (! !! + 2 ч 1. — 1)— для температуры У~"" и (741 + 2 э 1.) — для температуры 01". далее в коэффициенты матрицы А с указанными номерами добавляются слагаемые пп, оГ~ь пй (1 — 7',) согласно уравнениям (1.26), (1.27).

После формирования матрицы А и столбца Т выполняется обращение к стандартной подпрограмме СЕЕС, выходным параметром коюрой является'массив Т. В этом массиве первые М, элементов соответствуют температурам тел Т„а затем парами расположены температуры теплоносителей У;"", У,"". Эти температуры выводятся на печать. Если в подпрограмме СЕ1С параметр ошибки 1ЕВ принимает значение, отличное от нуля, то печатается соответствующее сообщение об ошибке.

$ КЗ. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ ЗАДАЧ КОШИ ДЛЯ ОВЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ РРАВНЕНИЯ Расчет нестационарного теплового режима по моделям с сосредоточеннымн параметрами сводится к решению систем уравнений теплового баланса вида (1.2), (1.3) с начальными условиями (1.6), т. е. к решению задачи Коши для систем обыкновенных дифференциальных уравнений первого порядка. В случае линейных уравнений решение удается представить в аналитическом виде при числе уравнений А! ( 4. Для нелинейных задач н в случае М ) 4 точное решение в аналитическом виде получить не удается, за исключением некоторых частных случаев.

Поэтому при расчетах нестационарных тепловых режимов систем тел широко применяют численные методы, которые мы сначала рассмотрим применительно и одному уравнению вида (1.29) — =1(т, Т), 0<т < т,„ дт с начальным условием Т! т=а — — Тм В дальнейшем методы, изложенные на примере задачи (!.29), (1.36), будут обобщены для системы А! уравнений. Основные понятия теории численных методов решения дифференциальных уравнений будут достаточно подробно рассмотрены в главе 3 на примере дифференциального уравнения теплопроводности. Сейчас лишь кратко сформулируем ряд понятий, которые понадо- 27 бятся для изложения практических вопросов численного решения обыкновенных дифференциальных уравнений.

При численном решении задачи непрерывная область изменения независимой переменной (О, т,„) заменяется множеством значении (тт) . з „которые будем называть узлами сетки. В случае равномерной сетки тз = 1Лт, 1 = 1, ..., /; Лт=т,„/l — шаг по времени. Вместо задачи отыскания непрерывной функции Т (т) ставится задача определения дискретного множества значений функции в узлах сетки: Тз = Т(тт). Величина Тг называется сеточной функцией точного решения. Как мы увидим дальше, точные значения Тз найти не удается, а вместо них в результате численного решения задачи получаются приближенные значения искомой функции в узлах сетки, которые будем обозначать и' и называть сеточной функцией разностного решения или просто разностным (численным) решением. Погрешность численного решения определим как разность сеточных функций точного и разностного решений: е~ = Т~ — и~.

Для нахождения разностного решения вместо задачи (1.29), (1.30) рассматривают задачу решения системы алгебраических уравнений относительно искомых значений (ит),.~,. Такую систему алгебраических уравнений называют разностной схемой. Прн измельчении сетки (при Лт- 0) погрешность численного решения должна в случае удачной схемы стремиться к нулю, т, е. 1пп )е/(=О. вч о Это требование к разностной схеме называют условием сходимости. Для сходимости разностной схемы необходимо н достаточно выполнения двух других условий — аппроксимации и устойчивости, которые будут пояснены ниже на примере схем Эйлера.

Если при достаточно малых Лт начинает выполняться неравенство (е('(=) Т'' — и!) к АЛт", А =сопз1, (1.31) то говорят, что разностная схема сходится со скоростью 0 (Лтг) или численный метод имеет порядок точности р, т. е. порядок точности характеризует скорость убывания погрешности при измельчении сетки. Теперь перейдем к рассмотрению наиболее распространенных методов численного решения задачи (1.29), (1.30).

К ним относятся: методы, основанные на разложении функции Т (т) в ряд Тейлора (наиболее распространенная схема этого вида — схема Эйлера); методы Рунге — Кутта; линейные многошаговые методы. Схемы Эйлера. Используя разложение решения Т (т) в ряд Тейлора в точке т,, можно записать следующее выражение для значения решения в точке тры — — тз + Лт: Т!+ ' = Тз+ Т'(ту) Лт+ 7" (т~) Лтд(2 — , '.... (1.32) 28 Пренебрегая в (1.32) членами порядка О (Лт') и выше и учитывая, что согласно исходному уравнению (1.29) производная Т' (т;) равна ~(тп Т!), получим следующую схему: ит'ь ~ =-ш+ Лт! (т~, и!), (1.33) называемую схемой Эйлера. Обратим внимание, что вместо ТГ мы записали и~, так как отбросив члены порядка О (Лт») в (! .32), мы уже не получим точное решение Т! = Т (т!). Из (!.ЗЗ) вытекает, чтозначениеи~-' вычисляется на основе значения иг для предыдущего момента времени ло лвной формуле.

Подобные схемы называются явными. Зная из начального условия (1.3О) значение и» =- Т, в момент времени т == О, легко по формуле (!.33) последовательно вычислить значения и', и-", ..., из во все последующие моменты времени. Кроме явных существуют неявные схемы, в которых значение искомой функции на новом временном слое из+' находится в результате решения уравнения, включающего это значение и значения для предыдущих моментов времени.