О применении и решении задачи Стефана в строительной теплофизике Текст научной статьи по специальности «Физика»

4./2011 ВЕСТНИК _7/202J_МГСУ

О ПРИМЕНЕНИИ И РЕШЕНИИ ЗАДАЧИ СТЕФАНА В СТРОИТЕЛЬНОЙ ТЕПЛОФИЗИКЕ

THE APPLICATION AND STEFAN PROBLEM IN BUILDING

THERMAL PHYSICS

H.A. Парфентьева, О.Д. Самарин, В.Л. Кашинцева N. Parfentyeva, O. Samarin, V. Kashinceva

ГОУ ВПО МГСУ

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

The paper offers analytical and numerical solutions of the Stefan problem, except this simple formulas for estimating the depth and time of freezing, which is of great practical interest.

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

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

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

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

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

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

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

Как мы говорили, одной из основных целей решения является определение скорости перемещения фазовой границы.

В простейшем случае одномерной задачи требуется решить систему уравнений

[1]:

!=о,х^

(1)

- = а2-г 8(.)< X<8^

при краевых условиях:

© = 0; ^(8с™,т)=; г1(8,т)=^(8,т)=(2)

Здесь х - промежуток времени с момента начала промерзания; ^ и ^ - температура среды соответственно в оттаявшей и мерзлой зоне (ось ОХ направлена вдоль теплового потока, т.е. в сторону падения температуры); а1 и а2, м2/с - коэффициенты температуропроводности материала в этих зонах, принимаемые с учетом данных [2]. В частности, а1 = (Уср)ь где с1, Дж/(кг-К) - удельная теплоемкость материала в оттаявшей зоне, а Вт/(м-К), и р1, кг/м3 - соответственно его теплопроводность и плотность. Величина ^ соответствует температуре наружного воздуха (пренебрегаем ее незначительным отличием от температуры на поверхности); 1пл - температуре фазового перехода. В дальнейших рассуждениях принято, что 1пл равна нулю. Величина 5СТ представляет собой толщину материала в целом, а параметр 5 - толщину оттаявшего слоя. Ясно, что при полном промерзании 5 = 0. Кроме того, должно быть задано начальное

4./2011 ВЕСТНИК _4/20|Т_МГСУ

распределение температуры по сечению материала, соответствующее стационарному режиму.

Обозначим массовое содержание замерзающей влаги в материале и ее удельную теплоту плавления соответственно как кг/кг, и Гщ,, Дж/кг. Тогда удельная теплота фазового перехода, отнесенная к единице объема материала, составит '№р1гпл, Дж/м3. На границе раздела фаз (на фронте промерзания или оттаивания) справедливо условие Стефана [1]:

= -(&) (3)

дх I дх )х _*_„ I дх

Знак плюс берется при затвердевании, знак плюс при оттаивании. В общем случае систему (1)-(3) можно решать, используя численное моделирование задачи с использованием ЭВМ. При этом необходимо использовать конечно-разностную аппроксимацию дифференциального уравнения теплопроводности и условия Стефана (1) - (3). Поскольку решается задача с подвижной границей, шаг по времени должен быть очень малым, а в этом случае для аппроксимации можно использовать только явную схему для обеспечения устойчивости численного решения [1]. Помимо существенного упрощения вычислений, такая схема позволяет за счет надлежащего подбора соотношения шагов по времени Ах и по пространственной координате h повысить порядок аппроксимации до уровня Ах2 + h4. Как показано в [3], для этого достаточно выполнения равенства aAx/h2 = 1/6, где а = 1/cpo, м/с - коэффициент температуропроводности материала стенки.

Тогда температура i-ro слоя стенки в j+1-й момент времени определится через температуры i-ro, i-1-го и i+1-го слоев для j-ro момента по выражению [1], [3]:

ki+1 = 6-у + 4tj +1,+u )/6. (4)

Перемещение фронта промерзания при этом можно вычислить по соотношению, являющемуся конечно-разностной аппроксимацией условия (3):

ХАх hwporm

Первоначально соответствующая программа для ЭВМ тестировалась для случая одностороннего промерзания стенки при ее нулевой начальной температуре и отличном от нуля значении tH. На рис. 1 пунктирными линиями приведены графики зависимости положения границы фаз от времени для слоя утеплителя из пенополистирола при следующих значениях теплофизических параметров [2]: с = 1340 Дж/кг-К, X = 0.052 Вт/(м-К), w = 0.1, р1 = 100 кг/м3, Вт/(м2-К) и темперетуре на границе tH = -5°С (это сопоставимо с суточной амплитудой для Москвы) и для кирпичной стенки ( с = 880 Дж/кг'К, X = 0.81 Вт/(м-К), w = 0.02, р1 = 1800 кг/м3). Сплошные линии соответствуют аналитическому решению, согласно которому фронт оттаивания определяется по формуле [4]:

8= ¡-^^(м), где k =-2- (6)

V wpr» 1 - ctul 2wrm

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

AS = - 2t. i + ti+1J) .(5)

5% в относительном выражении. Тем не менее, окончательное суждение о достоверности полученных результатов требует сопоставления с экспериментом.

0,250.20,150,1-O.OS-

00 5 10 15 20 25 30

Однако, например, процесс сезонного промерзания грунта, для которого существует достаточное количество экспериментальных данных, не вполне эквивалентен случаю, для которого получена формула (6), т.е. режиму tH = const, поскольку на самом деле наружная температура вначале понижается до некоторого минимального значения, а затем снова растет. Разумеется, существует ее средняя величина tx, °С, за период с отрицательными температурами, продолжительность которого равна z*, сут. Значения tx и Zx для различных населенных пунктов приведены в Таблице 1 [6]. Не вызывает сомнения, что иной характер изменения t„ будет сказываться только на числовом коэффициенте в выражении для 5.

Многочисленные расчеты с использованием конечно-разностной схемы (4)-(5) показали [7], что для такого режима формула (6) по-прежнему справедлива, но числовой коэффициент в выражении для k вместо 2 следует принимать равным л/2 = 1.57. Следует также заметить, что в рассматриваемом случае, в отличие от режима tH = const, поправка на влияние теплоемкости грунта, предложенная в [5] и замедляющая промерзание, не требуется. Во всяком случае, расхождения между результатами, полученными при разных теплоемкостях и прочих равных условиях, не выходят за пределы погрешности аппроксимации.

Сопоставление с фактическими данными по глубине промерзания грунта для ряда населенных пунктов России и СНГ [6] показано на Рис.2. Сплошной линией показана аппроксимация по формуле (7):

8 = 0.0457^,м (7)

Легко видеть, что разброс точек не превосходит 5-10% -обычной погрешности задания теплофизических параметров. При этом числовой коэффициент с точностью порядка 1% совпадает с вычисляемым из выражения (6) для глинистого грунта с учетом замены 2 на л/2 и очевидного соотношения z = 86400т. Таким образом, получен-

Рис. 1. Зависимость толщины промерзшего слоя (метры)

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

^ / 2

*

1 1

4/2011 ВЕСТНИК _4/2011_МГСУ

ная авторами здесь, а также в [1], [4] и других работах приближенная методика решения задачи Стефана полностью подтверждается экспериментальными данными.

Рис.2. Связь глубины промерзания грунта и произведения -1х*1х.

♦ По данным [6] По формуле (7)

Для практики представляет интерес получить возможность рассчитать время промерзания труб в условиях аварии в зимний период года.

В этом случае мы решаем осесимметричную задачу и для замерзшей зоны уравнение имеет вид:

а д ( дг \ дг

--1 г— ,- —

г дг ^ дг) дх

х = 0: г = 0, гр = Го; г = Го: г = гн; (8)

СгЕ ( дг Л

Г=Гр: рг се +о

Здесь г0, м - радиус трубопровода; гр, м - радиус фронта промерзания. Считая в первом приближении, как и при выводе (6), температурное поле в промерзшем слое подчиняющимся стационарному распределению, для безразмерного радиуса г р = гР/го получаем дифференциальное уравнение:

~Гр 1п[гЕ^сСгр = сСЕв (9).

Параметр Бо'' представляет собой модифицированное число Фурье (безразмерное время), определяемое в данном случае по выражению (10):

(10)

Рг„Л

Интегрируя (10) в пределах соответственно от 1 (поскольку при х = 0 гр = г0) до г р и от 0 до х, окончательно находим:

ВЕСТНИК МГСУ

4/2011

2Fo =- + r, 2 F

In lrF

(11)

Если теперь положить гБ = 0, для времени полного промерзания получаем Бо'' = %. Учет теплоемкости промерзшего слоя выполняется аналогично (6) умножением

Бо'' на поправочный коэффициент 1 — &н/-\Г2м>гт (в данном случае числовой коэффициент 2 стоит под знаком корня из-за цилиндрической симметрии) [8]. Полученные результаты дают возможность располагать заметным резервом времени при ремонте аварийных участков тепловой сети в холодный период года. Соответствующие зависимости имеют достаточно простой вид и удобны в инженерной практике.

Литература

1. Прусаков Г.М. Математические модели и методы в расчетах на ЭВМ. М.: Наука, 1993, 144 с.

2. СП 23-101-2004 «Проектирование тепловой защиты зданий». М.: ГУП ЦПП, 2004.

3. Теория тепломассообмена / Под ред. А.И.Леонтьева. М.: Изд-во МГТУ, 1997, 684 с.

4. Богословский В.Н. Строительная теплофизика. М.: Высшая школа, 1982, 415 с.

5. Самарин О.Д. Теплофизика. Энергосбережение. Энергоэффективность. М.: Изд-во АСВ, 2009, 296 с.

6. СНиП 23-01-99* Строительная климатология. М.: ГУП ЦПП, 2004.

7. Парфентьева Н.А., Самарин О.Д. Об экспериментальном подтверждении решения задачи Стефана // Кровля и изоляция. 2005, №3, с. 51 - 52.

8. Парфентьева Н.А., Самарин., О.Д. Решение задачи Стефана при промерзании трубопроводов // Вестник МГСУ. 2007, № 1, с. 67 - 70.

References

1. G. Prusakov. Mathematical models and methods in calculations on EVM. M.: Science, 1993, 144 pp.

2. SP 23-101-2004 "Design of thermal protection of buildings". M.: GUP ZPP, 2004.

3. Theory of Heat and Mass Transfer. Ed. A. Leontiev. Moscow: Publishing House of Moscow State Technical University, 1997, 684 pp.

4. V. Bogoslovsky. Building Thermal Physics. M.: Higher school, 1982, 415 pp.

5. O. Samarin. Thermo physics. Energy conservation. Energy efficiency. Moscow: Publishing House ASV, 2009, 296 pp.

6. SNiP 23-01-99* "Building Climatology. M.: GUP ZPP, 2004.

7. N. Parfentyeva, O. Samarin. An experimental confirmation of the decision of the Stefan problem // Roofing and Insulation, 2005, №3, pp. 51-52.

8. N. Parfentyeva, O. Samarin. The solution of the Stefan problem with freezing pipes // Vestnik MGSU. 2007, № 1, pp.67-70.

Ключевые слова: Задача Стефана, глубина промерзания, время, аппроксимация точное решение, численные методы

Key words: Stefan problem, the depth of freezing time, an approximation of the exact solution, numerical methods

Парфентьева H. A.: +7-495-4704244; +7-925-0544560 n.grigorieva.13@gmail.com

Рецензент: Либерман А. А., д.т.н., зав. лабораторией, ФГУП «ВНИИОФИ»