А.А. Самарский, Ю.П. Попов - Разностные методы решения задач газовой динамики (1161630), страница 30
Текст из файла (страница 30)
Разностные схемы для уравнения теплопроводности 1. Консервативная схема для уравнения теплопроводности. Выше были рассмотрены некоторые способы построения разностных схем, аппроксимирующих систему одномерных нестационарных уравнений газодинамики без учета реальных диссипативных процессов. Обратимся тепарь к случаю, когда в задаче присутствуют процессы теплопроводности. Изменится лишь уравнение энергии, которое для одномерного плоского случая имеет вид дз ди дГà — — 77 — — —, дс дг дз' (6.1) где е(Т) — удельная внутренняя (тепловая) энергия, Иг = = — крдТ/дз — поток тепла. Введем обозначение )г = рк(р, Т). С учетом сказанного имеем (6.2) Иг = — гг —.
дт (6.3) дг Здесь 7' = — р ди7дг — источник энергии газодинамического происхождения (з этом параграфе его мощность мы будем считать известной). Кзазилинейное (коэффициент й, вообще говоря, есть Функция температуры) уравнение параболического типа (6.2) и является уравнением теплопрозодности. В интегральной форме зто уравнение выглядит следующим образом '~' г г!д — Игг77 -. ( 7(д, Г) сЬпд (6.4) где С вЂ” заикядтый кусо ~по-гладкийг контур, ограничипающяй на плоскосттг (г. г) область С.
В соотзетстзпп с идеей иктогро-нптерполяцнонпого метода (гм. 1 2, и. 4) используем интегральное уравнение (6.4) для построения консервативной схемы. Введем в рассматриваемой части плоскости (з, г) сетку е7 = ((г,, 1,). (з~„~ в Г,), гы1 = з, + йе гыпд = з;+ О 57Ьь г О, 1, ..., У вЂ” 1, гн~ =г,+ть 7'=О, 1, ...) 143 Сеточную функцию потока тепла И'в ~будем относить к узлам сетки (гь й), а к полуцелым точкам (зв+нм 1,) — температуру Т,'+ем = Т';, внутреннюю энергию а„'+пв = с'„коэффициент теплопроводности йв+нэ — — Хр который определяется значениями температуры Т,', и плотности рв, и источник газодинамической энергии )вв+ма= Я (см.
рис. 2.22). Далее черточки вад функциями вс, ~, с, Т и т. д. мы онусим. Применим интегральное уравнение (6.4) к контуру С, охватывающему элементарную ячейку разностной сетки (рис, 2.22). При этом используем следующий способ вычисления пптегралов, входящих в (6.4): в; .в с(г, Г;)овт звв+М,Ь» =- г',!во Вв О+1 И'(г„1)вЬ ж [оИ'",+'+ (1 — о)Ивет; — И'; т;, увг Г-l в-в" Г, ~ взт ввх Г . а,'. Р ', Г' Раа.
2.22 параграфах. Учитывая все сказанное, приходим к следующему разностному уравнению, аппроксимирующему дифференциальное уравнение (6.2): е(+~ — з' Г яв(+1 вав+1 ив( иве ~ в+Ма в+1/2 в+1 в вв в в+1 в,в + Д+ню в 1 ~в или в безындексной форме зв = И'в +в' (а) (6.5) (6.5') '+в 6+1 ) ) (э, Е) вЬ Ж = ),'Л;т,. в. вв Здесь 0 ~ о - 1 — свободный параметр, позволяющий варьировать интерполяцию сеточной функции И'. Способ вычисления газодивамического источника )( подробно описан в предыдущих ! Построим теперь развостное уравнение для теплового потока (6.3). Проинтегрируем для этого на временвбм слое 6 по отрезку 21 112 ~ г < г,+1,2 уравнение дТ/дг = — '41///4, вытекающее из (6.3).
Получим 41 + 1/З И/(4, 1,) Т,''„„Т', „,= ,) 2 (4, 1.) 44 — 1/3 и далее (6.0) Напомним, что беэывдексвое обозпачепне у; соответствует раз- востной левой производвой, вычисленной на неравномерной сет- ке. Здесь '1 ';1/З -1 1 ! 64 О 5(Ь/т Ьз,) .) 2(4, 1;) 41,/., (6,7) '1+1/З Величиву ) //г,'/4 называют тепловь м сопротивлением отрезка 4; 1/2 !21- 1~2, 2<+1/2! 2. Аппроксимация кОзффициснта теплОпроводпости. Интеграл в выражевии (6.7), аппрокснмирующем коэффициент теплопроводности, можно заменить различвыми разноствыми соотношевиями, например: 44.1.1/2 Вг — — й1 = 0,6 (/11 + /11 1), /с' ° 1+1/З 4; 6+1/з чь ! 2 ! 4 4 44 — 1/2 41 гДе /4/зп/2 = /2(Рем/2, Т, 1/2). 10 Л, Л С221Г4222, Ю П.
П4244 145 Считая, что функция И/ в промежутке г; пзбг(г,+1,2 имеет постоянное звачевие, равное И/1 имеем 44.1.1/2 / / / ( гз Т/+1/, — Т4 !/2 — — — И71 ,! 2 (4, 1 )' '1-1/з Соответственно коэффициент а"; будет иметь вид / и/ й/ /г/ г!/ю ь/ й/ $ — 1/и 1+1/2 (-//с г+//г( / ч -1) а' (6.8) Ьс где /г'„= ', — нн герполяция коэффициевта те- / /-1 «+1/$ 1 1 — 1/с Ьг-,-/г, плопронодностк в целу/о точку на кгранкомерпой сетке.
Заметим, по приведенные формулы (6.8) для вычисления козффициевта теплопроводпости а; в некоторых условиях могут оказаться неэффективными. Действительно, пусть / = О, а коэффициент теплопроеодпости й зависит от температуры так, что при Т = 0 обращается в нуль (й = 0). Рассмотрим задачу о распространении тсмксратурпой волны по холодному фону (Т=О). Она имеет в этом случае конечный фронт (см. подробнее в п. 4). Предположим, что в некоторый момент времени фрОнт тепловой (температурной) волны, распространяющейся направо, совпадает с гчм узлом сет'.
т „, ар г к, -// ки, так что Тгы/и = О, а Т/-//з чь 0 (рис. 2.23). В этом случае разностный каэффициент теплопроводности Ь;, вы- 2 г '-', /;+у численный по второй и третьей фор- мулам (6.8), тождественно равен нулю, Рис. 2,23 ибо й;+ьа = й(Т,+ьа) = О. Поэтому по- ток тепла И/( в узле 4 также будет нулевым. Это означает, что тепло никогда не проникнет в следугощий массовыи интервал Ь| сетки и температура Т,„~/и здесь так и останется нулевой. Такая картина, полученная в рамках дискретвой модели, противоречит физическому смыслу.
Очевидно, первая и третья формулы (6.8) приведут к аналогичному результату при расчете тепловой волны, распростравяющейся налево. Избежать подобпых недостатков можво, воспользовавшись для вычисления //( например, формулой, представляющей собой полугумму первых двух соотношений на равпомервой сетко н (6.8): //; = (!'; г/и + 1'/+г/с)/2. (6.9) В силу того, что коэффициент теплопронодности зависит от тем- пературы и плотпостп й,.г~/з = й(Т,.г1/и рчю/з) = Й(То Р ) коэф фнцкспт а'; разпостпого уранпепия и /чм уэлс и (6.9) будет определяться значепиями температуры и плотности в соседних полуцелых точках и/.—.— а(Т(, Т,',, р,'.
р,',), $46 Сводя вместе разностные уравнения (6.5) и (6.6), имеем е/ = — И/, + ), И' = — аТ, (101О) или,после исключения функции И', е/ = (аТ-) + 1. (6.11) или, в раэностной форме, Т14+1/э 44 (84+1/2) 4 0 1» / 1 (6.12) Краевые условия в зависимости от конкретных физических особенностей задачи могут иметь различный вид.
Весьма распространенными являются задачи, где на граннцах рассматриваемой области заданы режимы изменения со временем температуры (1 краевая аадача): Т(0, 1) = 0" (1), Т(м, 1) = 6""(г), (6.13) или теплового потока (П краевая задача): И'(О, г)= "(/), И'(М, г)= ш* (1). (0.1/4) Встречаются и более сложные случаи, когда на границе за- дается связь между тепловым потоком и температурой: И'(О, 1)= шэ[Т(0, 1), Е), И/(М, 1)= шээП[Т(м, Г), г). (615) Функции шэ и ш**, вообще говоря, являются нелинейными. Примером может служить задача, в которой масса газа излучает с границы как абсолютно черное тело: И (О, 1)=.—,Т4(О, г), И/(М, г)= а,т (М, г), где ао — постоянная Стефана — Больцмана.
104 147 Использование интегро-интерполяционного метода позволило построить консервативную схему для уравнения теплопроводиости (6.2), (6.3), имеющую второй порядок аппроксимации по пространству 0(Ь~) на равномерной сетке. Значению а = 0 в (6.10), (6.11) отвечает явная схема, а значению а=1 — чисто неявная схема; обе они имеют порядок аппроксимации 0(т).
При а = 0,5 получаем схему второго порядка аппроксимации по времени 0(тт); ее называют симметричвой схемой или схемой /4ра/4ка — Николасе/4а. Заметим, что разностные уравнения (6.10) и (6.11) записываются одинаково во всех точках сетки, даже если коэффициент теплопроводности разрывен. Поэтому сформулированная схема является однородной. 3. Формулировка краевых условий. Постановка задач для уравнения теплопроводности в ограниченной пространственной области включает формулировку краевых и начальных условий. В качестве начальных данных обычно задается распределение температуры среды Т(г, 0) = 4с(г), 0 < г < М, Возможны также комбинации условий (6.13) — (6Л5) на правой и левой границах. Обратимся к способам записи краевых условий в разностной задаче.
Уравнение (6.10), или (6Л1), записанное, как и исходное уравнение энергии, в полуцелой точке, рассматривается на расширенной сетке, включающей фиктивные крайние интервалы Ь ~ = 0 и йл = О. Значения температуры в этих интервалах фактически совпадают со значениями температуры в граничных узлах сетки (рис. 2.24): Т' т = Т» Ть. = Ти Поэтому краевые Ф-1 Ф Ф+Т Рис.