Метод расчета тепловых режимов многослойных полупрозрачных материалов Текст научной статьи по специальности «Физика»

УЧЁНЫЕ ЗАПИСКИ Ц А Г И Том XVII 1986

М2

УДК 536.33

МЕТОД РАСЧЕТА ТЕПЛОВЫХ РЕЖИМОВ МНОГОСЛОЙНЫХ ПОЛУПРОЗРАЧНЫХ МАТЕРИАЛОВ

В. П. Тимошенко, М. Г. Тренев

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

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

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

ничных поверхностей слоев. В работе [3], где рассматривалась задача прохождения внешнего излучения сквозь многослойные материалы, предлагается сращивать решения уравнений переноса излучения в каждом слое с использованием рекуррентных соотношений. Процессы теплопроводности и собственное излучение материалов авторами не рассматривались, однако есть основания ожидать, что при некоторых условиях использование рекуррентных соотношений не обеспечит устойчивости вычислительного процесса при решении сопряженной задачи радиационно-кондуктивного теплопереноса.

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

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

Рассмотрим вначале процессы переноса тепла теплопроводностью и излучением в одном полупрозрачном слое с номером т (1 </псМ). Распределение температуры по толщине слоя описывается уравнением

где С — объемная теплоемкость, qR — суммарный радиационный тепловой поток внутри слоя, А, Та, <3 — коэффициенты, определяющие интенсивность объемных источников тепла в линеаризованном виде относительно температуры Т.

Граничные условия для левой и правой границ слоя Х°т, Х1т запишем в виде

дТ _ д дt дх

^е. + А(та-Т) + (1,

(1)

где <7т, Яи> — потоки, обусловленные конвективным теплообменом, теплопроводностью и поверхностными источниками тепла, причем потоки <7К и q^í одновременно не реализуются.

Граничные условия первого рода (Т = Та) специально не выделяются, так как могут быть реализованы при достаточно больших значениях а. Тепловая емкость тонких покрытий считается пренебрежимо малой величиной, в противном случае эти покрытия могут быть рассмотрены в качестве самостоятельных слоев.

Для построения алгоритмов расчета введем в рассматриваемом слое две неравномерные разностные сетки (рис. 1)

Хп^={Х°т = Хт,0, Хт,1, ..., Хт,н = Ххт),

Зт = {5т = 5т, о, • • • , Зт, р, . . . , 5т, Л? +1 == 5т}»

координаты узлов которых Хт, ¡, 5т,р и шаги Нт,и Нт,р удовлетворяют соотношениям

Ьт, I = Хщ,» Xт, I—1, Нот, р — 5т, р " $т, р—1, , р == Хт< I (/2,

Нт, р+1 = (Ат,г+ Лт, «+1)/2, (/ = /?; г,/?=. 1, УУ), (3)

Лт,0 == Лт—1, Л?, 5т,0==-‘^пг, 0 — Лт—1,Ту/2, 5т, Л/ + 1 — ЛГ “Ь ^т+1, 1 /2.

Сетку АГт, в узлах которой А'т,,- вычисляются значения температур, будем считать основной, а сетку 5т, в узлах которой 5т, р (сдвинутых влево относительно узлов Л"т,,) вычисляются значения потоков, — вспомогательной. При этом всегда ¿=р, если г и р используются одновременно.

В соответствии с интегро-интерполяционным методом построения разностных схем [4] запишем уравнения (1) в виде баланса энергии для интервала (¿у 5р+1), содержащего узел сетки Х1 (1 = р, 1 <(/, —1). Обозначим /ч-, р1+ средние значения

7П -/

хк-1

ТП,0

т.

,1-1

N

¿-1

р-1

н,

т.,р

рЧ

И,

т,р*1

¿+/ N-1

N

N

т + 1 X 0

лт + !

N+1

Нт,Ы+1

"т., р-1

г,р*1

°т,Х+1

Рис. 1

произвольной функции F(x) на интервалах (Sp, Xt) и (Xly Sp+\), тогда, проинтегрировав (1) на интервале (5р, Sp+i), получим

±- (С, _л + с,+ А,+1) = (x|^+i -{>■■£)-

— р+1 + qR<p -f -g-(A- ТА, /-Л/ + Ai+T А> i+hi+i) —

---2~ (Ai— hi + Ai+ hi^ i) Tt -]—(Qi_ Лг + Q/+(4)

Уравнение (4) остается справедливым для зоны стыка двух слоев, теплообмен между которыми определяется теплопроводностью и потоками qw. Действительно, проинтегрировав (1) на интервалах (Sm-l, N, Х]п_х), (Хт, Sm, 1) И СЛОЖИВ результаты, МОЖНО

получить с учетом (2) уравнение, аналогичное (4), с дополнительным членом в правой части:

q1 , + о0 (5)

*w, га—1 1 “w,m v '

Если слой т представлен жидкой или газообразной фазой и его теплообмен со слоями т—1 и т+1 осуществляется конвекцией, то уравнение (4) здесь также можно использовать, если задать в (1) и (2) следующие значения коэффициентов:

С = А = Q = О, А = 1, TA=.fA(x, t), {Xl<x<Xlm),

1) — q^i-Km), qR {Sm, n) = qR (Xm),

qR(Sm,p) = 0, (p — 2, A/—1), m =f i'i,,m = 0,

где /а (x) — заданный профиль температуры в слое т.

Для границ раздела слоев (т—1, т) и {т, т+1) в уравнение (4) вместо (5) войдут следующие дополнительные члены:

Таким образом, показано, что уравнение переноса тепла вида (4) с дополнительными членами в правой части (5) или (6) справедливо ДЛЯ всех узлов СКВОЗНЫХ и разностных сеток Хг, (1 = 0, Ы) и

М \

р = 0, Л/+ 1, N = 2 Мт )> составленных из локальных сеток Хт

т — \ /

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

Перепишем (4) для простоты в виде

с, 4^ - Х,+, (4^, - X, (+ о, Т, + Л/,+1 Е,. (7)

Выражения для коэффициентов Си Ей б* не приводятся, так как легко могут быть получены из (4) с учетом (3), (5) и (6).

Обозначим

№р+. = Хр+, т‘+' ~ П , И, р = бГло, л/+1

(6)

тогда разностная аппроксимация (7) с порядком 0(А2'+т,) имеет вид

где

\Рр - \Ур+1 -Р1Т1 = - £>„ (/, р = О, ЛО, (8)

Т1 - Т-1+1 = ЧГр+%/Ьр+и (/, Р = О, М-\), (9)

Р1 = (С,1ъ-<ЗдНр+1, £г==(ГгСА4-£,)Яр+ь Ьр+1 = Ьр+11Н1+1,

V

т,— шаг счета повремени, Т1 — значение Т1 в предыдущий момент,

V

£ = / — Т, .

Форма записи уравнения (4) и соответственно (7) позволяет учесть граничные условия на внешних поверхностях рассматриваемого пакета с помощью коэффициентов объемных источников тепла. Поэтому недостающие разностные уравнения для системы (8), (9) получим, используя уравнение (8) ^для р = 0, р = Ы\ и соотношения"

1^0 + 1^ п ( ■1дт\ _и7л? + ^+1

дх )х- 2 —и’ I кдх )хы— 2

(П)

Таким образом, на границах пакета имеем Ш0-Т0Г0/2 = -а012, )

№7г — Ты /^/2 == — Ду/2, / (1и)

По аналогии с [5] можно получить расчетные формулы для решения системы (8), (9), (10).

При 1 <&р+1<;оо (I, р=О, N—1) прогоночные коэффициенты ар, фр вычисляются следующим образом:

ао = 2//70> Ро — О01Р0,

ар+1 — (ар + уьр+1)1г,

Фр+1=(ар £—\ + /"¿+1 (ар 4- 1/^р+1).

Если 0<6р+1<1, то вместо (11) используются формулы 'Л

ар+1 = (14- V1 ар)/г, Рр+1 = Ьр+1 (ар Д + рр)/£,

Z= Ьр+1 + 1+1 (^р+1 !)•

Формулы обратной прогонки имеют вид

Н^+1 = [— -Олг/2 + /^/2 фдг + “я Оя)]1(Рц/2 — 1),

^=(1Гр+1-А)(1-/=■;«,)+.^рл> (¿,

7\ = — (№¡,+1 — Д) *р + Рр.

Расчет радиационного теплового потока <7д является сопряженной задачей по отношению к решению уравнения теплопроводности. Предполагая поле температур в текущий момент времени известным, по-

строим алгоритм сквозного расчета потока <7Н в рассматриваемой многослойной системе.

Разобьем весь интервал длин волн излучения (0—оо) на спектральных интервалов

д/-, = ^.с/- 1, М),

в каждом из которых оптические характеристики слоев можно считать постоянными. Интегральный по длинам волн Я-поток результирующего излучения qmR в слое т записывается через направленные потоки селективно серого излучения следующим образом:

ЛЬ 1

Я ту)<

]=1

д+. = 2іг /+} !^, Яті = } ІтічЛіЬ [а = СОЭ Є,

О о

(12)

где /то/ — интенсивности излучения в положительном и отрицательном направлении оси х, 6 — угол между направлением луча Я осью х. Ось х проходит вдоль направления нормали к поверхности слоя.

Уравнение переноса селективного серого излучения без учета рассеяния имеет вид

/ в„

йх

[1Л^-==Хт/((13)

где хті — коэффициент поглощения, Вті — поток полусферического излучения абсолютно черного тела в интервале ДЯ3-, определяемый согласно закону излучения Планка по формулам [1, 2]:

ВШ) = [? («ту — У Г)1-

пХТ

Г сг/9 .,

<Р (лХУ*) = J (лА.Г)5 [ехр (с2Цп\Т)) - 1] din.IT).

Введем оптическую координату Тпц- слоя т по формуле

*

хт/ (■*•) = *т)(Хт) = Х/п/' — О,

Х°т

т) — хт/.

Будем считать, ЧТО граничные поверхности Л"т И Л'т вместе с имеющимися на них покрытиями пропускают, отражают и испускают излучение внутри слоя т полусферически изотропно и их свойства со стороны, обращенной внутрь слоя, определяются значениями спектральных полусферических характеристик г®, е®, ггт,

(пропускная способность, отражательная способность, степень черноты; й+Г + Е=\).

Тогда с учетом (12) решение уравнения (13) для произвольного спектрального диапазона Дк] можно записать через интегральные экспоненциальные функции Еп(т) следующим образом:

41

Я+т Ы = 2?+ (0) Ег (ТИ) + 2 / Вт (6) Е2 (*т - 5) Л,

о

*

Я~т (*«) = 2^- (^) Е, рт - хт) + 2 | Вт (?) Ег (6 - *т) Л,

(14)

Еп (*) = / ^П-2 ехр (- */ц) ¿1»,

о

1

£„<0) =

п — 1

» ^я(°°) = 0-

(15)

Соотношения между потоками излучения на границах смежных слоев (/я—1, т) и (т., т+ 1) имеют вид

Я+т (0) = г°т (0) + Вт (0) е0т + ¿1т_1 д+_г (х1т_1), Ю = Г« <7т Ю + Ю + ¿°т+1 ?т+1 (°)'

(16)

Для того, чтобы выражения (16) можно было использовать и на внешних границах X?, Хм рассматриваемого пакета слоев, введем два дополнительных слоя с номерами от = 0 и т— М + \, полагая в каждом из них потоки излучения q\ = ^ (т£), д°м+1 =

= Ям+\Ф) и оптические характеристики ¿¿, го, 4, с1%+и А+ъ ®ж-н— заданными в соответствии с физической постановкой задачи.

Как видно из (14), для вычисления в любой точке х значений потоков Я%[1т{х)] необходимо знать граничные потоки <7+(0), д-(г]п), направленные внутрь соответствующего слоя т. Обозначим для краткости

пт= <7т(0), лт = Я^т) (РИС- 2),

Рис. 2

тогда выражения (14) для = ^ и т„ = 0, а также соотношения (16) можно представить в виде

?«<XL) = /?mnr» + Gm.

Ят (°) — Лт + Нт, пт = г°т д~ (0) + Фт_, д+_, (Х^_х) + Вт (0) £°т. лт = r'm q+ (Х^) + df> q~+1 (°) + Вт (z'J е'т>

(17)

(18)

где

Rm=2E3(,]n),

J "т

Gm = 2/BmEt[?m-%)d\, Нт = 2 j BmE2{t~^m)dl

*

(19)

Исключив ИЗ (18) ПОТОКИ д%_ит, т+Р подходящие к внутренним границам соответствующих слоев, с использованием (17) и аналогичных выражений для слоя /га — 1, можно получить

где

Пт= rfii-l 1 Пт_1 +Гт /?т Лт -}- Кт, ТП— 2, Ж, •^т — Гт d0mjf.\RmjriJ\m^.\-\- Nm, т— 1, М—1,

Кт= Гт Нт + (0) гт -j- Gm-l,

Nm = Гт Gm -f- Вт (Тт) sm Ч~ dm + l •

(20)

(21)

Выражения для потоков и Лм получаются аналогичным образом из (17), (18) с учетом заданных характеристик для дополнительных слоев т = 0 и т = М~\-Аш.

(22)

где

Л, + к».

Лм = ггм Ям Пм + Мм,

к, = г? Нх + ВХ (0) 8? + ¿1 т - 1,

Мм — г]ц Ом + Вм (х’и) + еРм+1 <рм+1, т = М.

Соотношения (20), (22) представляют замкнутую систему 2М линейных уравнений относительно потоков Пот, Лт, (т= 1, М). Для дальнейшего анализа запишем ее в матричном виде

Вд = /

или

1 1 0 А\ В% С2 Qi /,

X _ /2

Ат Вт Ст Qm /ж

0 .... . . . . . .

Ам Вм Ям fм

где

dm-1 Rm-1 0 1 r° f? ' т ‘\т

II гт Rm 1

О О

0 0 Пт Кт

0 dm-^-i > Чт , /« =

Лт Nm

Для существования решения системы (23) необходимо, чтобы ранг блочной трехдиагональной матрицы В совпадал с рангом расширенной матрицы (В, /), [6], т. е.

rang(,B,/) = rang(Z?) = я. (25)

При этом, если п = 2-М, (йе1(В)ф 0), решение системы (23) единственно.

Из (24) следует,, что поэтому

м

м

det (В) = П det (Вт) = П (1 - Г°т гЖ) ,

т =1 /га = 1

М

п = rang (В) = £ rang (Вт), [1 < rang (Вт) < 2].

/72 = 1

Таким образом, п<2-М, если имеется по крайней мере одна вырожденная матрица Вт [det (Вт) = 0], что в, соответствии с (15), (19) возможно только при

г° __ г1 _ Е> -- 1

• т • т Д/п—*>

(в® «el, = d° = .

Физически это эквивалентно оптически прозрачному слою (т^ = 0) с идеально отражающими внутренними поверхностями.

При наличии в (23) вырожденных матриц Вт условие (25) для соответствующих значений т принимает вид rang (Вт, f т) = 1, что с учетом (24), (26) может быть выполнено только при Km=Nm. Из (19), (21), (26) следует, что при произвольных значениях оптических тол-

щин смежных слоев “tm-i, Tm+1 ЭТО ВОЗМОЖНО ПрИ УСЛОВИИ, ЧТО

dm-1 = d°m+1 = 0 (fm = 0),

(27)

т. е. излучение не должно проникать извне в слои т.

Непосредственной проверкой с использованием (26), (27) легко убедиться, что матрицы Ст_ь Ат, Ст, Ат+\. окружающие вырожденную матрицу Вт, являются нулевыми, следовательно, система (23) распадается на ряд независимых подсистем линейных уравнений, в каждой из которых матрица является либо блочной трехдиагональной, либо вырожденной матрицей Вт. В последнем случае множество решений определяется соотношением Пт=Лт.

Поскольку слой т является диатермическим и не обменивается излучением со смежными слоями, удобно принять что Пт = Лт=0.

В подсистемах линейных уравнений с невырожденными матрицами решение может быть найдено одним из прямых численных методов, например методом исключения Гаусса. Вычисленные таким образом значения потоков Лт, Пт позволяют рассчитать потоки (т:т) в любой точке рассматриваемого пакета слоев.

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

те профилей температуры и радиационных тепловых потоков на каждом шаге счета по времени. В связи с нелинейностью задачи целесообразно применение итераций, число которых обычно не превышает 2—3 для достижения сходимости с погрешностью меньше 1% .

В качестве иллюстрации рассмотрим результаты решения модельной задачи прогрева пакета из пяти слоев. Средний из них с номером III является диатермическим зазором, а оба примыкающих к нему слоя II я IV я внешние слои пакета I я V имеют попарно одинаковые свойства, которые в расчетах принимались близкими к свойствам кварца. Различия между слоями состоят в том, что слой V считался непрозрачным, а коэффициенты теплопроводности слоев / и У заданы в пять раз большими, чем для слоев II и IV. Радиационный перенос энергии рассчитывался в трех спектральных интервалах с границами раздела 2,7 и 4,5 мм.

Начальная температура слоев пакета задавалась равной 300 К-Нагрев внешних поверхностей пакета осуществлялся постоянным тепловым потоком 50 кВт/м2 в течение 200 с, затем на этих поверхностях задавалось условие полной теплоизоляции и расчет проводился до полного выравнивания поля температур в пакете.

Определяющие параметры задачи были подобраны так, чтобы точное конечное значение температуры составляло 800 К. Целью расчета являлся анализ характера перераспределения поля температур в пакете и проверка точности выполнения баланса энергии.

На рис. 3 показаны профили температур по толщине пакета в моменты 200, 400 с и асимптотическое точное значение 7=800 К при ¿->оо. Видно, что в оптически полупрозрачном слое / перенос энергии осуществляется значительно интенсивнее, чем в непрозрачном слое V,

Рис. 3. Распределение температур по толщине пакета в различные моменты времени 1—У—номера слов; 0—5—номера граничных точек

Рис. 4. Изменение по времени температур в граничных точках слоев

приводя к более быстрому выравниванию температуры в слоях / и II по сравнению со слоями IV и V. Это подтверждается и характером стремления температур в граничных точках слоев к асимптотическому значению (рис. 4). Наиболее медленно поле температур выравнивается в непрозрачном слое V.

При использовании шага счета по времени, равного 1 с, точность выполнения баланса энергии в любой момент и отклонение асимптотического значения температуры от точного значения 800 К не превышает 0,5%. При шаге счета 5 с соответствующие значения погрешностей возрастают до 0,7 %, что свидетельствует о достаточно высокой точности предложенного метода.

ЛИТЕРАТУРА

1. Зигель Р., Хауэлл Дж. Теплообмен излучением. — М.: Мир,

1975.

2. Спэрроу Э. М., Сесс Р. Д. Теплообмен излучением. — Л.: Энергия, 1971.

3. Ильясов С. Г., Красников В. В. Перенос энергии излучения в многослойных рассеивающих и поглощающих материалах.— ИФЖ т. XXVII, № 6.

4. Самарский А. А. Введение в теорию разностных схем. —

М.: Наука, 1971.

5. Дегтярев Л. М., Фаворский А. П. Потоковый вариант метода прогонки. — Ж. вычисл. мат. и матем. физ., 1983, т. 8, № 3.

6. Бронштейн И. Н., Семендяев К. А. Справочник по математике.— М.: Наука, 1980.

Рукопись поступила 15/У 1984 г.