Численное моделирование и анализ расчета теплообмена в наружных ограждающих конструкциях с особенностями

ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ

УДК 519.6:697.1

Н. Н. Павлов, В. Ю. Шадрин

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ И АНАЛИЗ РАСЧЕТА ТЕПЛООБМЕНА В НАРУЖНЫХ ОГРАЖДАЮЩИХ КОНСТРУКЦИЯХ С ОСОБЕННОСТЯМИ

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

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

N. N. Pavlov, V Yu. Shadrin

Numerical modeling and analysis of heat exchange calculation in outer bounding constructions with peculiarities

The article is about working out and analysis of digital method of finite-difference equation solution approximating the third boundary value for two-dimensional Laplace equation. The effectiveness of the conjugate gradient method for computation of stationary temperature fields of outer bounding buildings’ constructions and constructions with angular elements consisting of materials with thousand-fold differences of heat passage values is showed.

Key words: bounding constructions, math modeling, differential scheme, the conjugate gradient method, approximation error.

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

Как показал опыт работы с программным продуктом «SHADDAN» [1], разработанным авторами, возникают некоторые проблемы при расчете температурных полей наружных ограждающих конструкций с угловыми элементами и наличия в ограждениях соседних материалов с несоизмеримыми размерами и

ПАВЛОВ Никифор Никитич - к. ф.-м. н., доцент кафедры информационных технологий Института математики и информатики СВФУ им. М.К. Аммосова.

E-mail: pnn10@mail.ru

ШАДРИН Василий Юрьевич - к. ф.-м. н., профессор кафедры высшей математики Института математики и информатики СВФУ им. М.К. Аммосова.

E-mail: vshadr@mail.ru

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

Для решения систем уравнений, возникающих при дискретизации третьей краевой задачи для двумерного уравнения Лапласа с переменными коэффициентами, в программе «8НАЭБА№> применялся итерационный метод Зейделя, который неплохо зарекомендовал себя для достаточно широкого класса наружных ограждающих конструкций и был выбран из-за простоты реализации и удобства для организации интерфейса.

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

тями.

1. Постановка задачи

Рассмотрим третью краевую задачу для двумерного уравнения Лапласа с кусочно-постоянными коэффициентами теплопроводности:

дТ

дТ

= 0 , (х1, х2) еО,

лдТ+аср (Г - Тср )=0, (X! , х2) єда,

где Т (х1, х2) - температура точки конструкции с координатами ( Х1,Х2 ) , Л = Л(х1,Х2) - коэффициент теплопроводности, ЭО - граница области О, а - ко-

эффициент теплообмена со средой, Тср - постоянная

дт

дп

температура внешней среды,

производная по

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

Для решения задачи (1)-(2) применим метод сеток.

Введем неравномерную сетку

Ък = {(Xа , х2а), іа= 0,1,..., Иа, а = 1,2}, щ =ак + ук,

Г\

.... н П И н г- / \уЬ.

1 .

Пг

ТІ '

— — — и «—

Рис. 1. Сеточная область

где ($ - множество внутренних узлов, у - мно-

Н Н

жество граничных узлов, = х1а — х1а 1,

К = %= 0.5( Н + к1/1).

На рис. 1 знаком о обозначены внутренние узлы, а знаками х, р, Ц ^, ^, ^, ^, ^, ^ обозначены граничные узлы. Для множества граничных узлов, обозначенных знаком X, введем множества левых

и правых граничных узлов у^ 1 и у1аг), а = 1,2.

Обозначим ^<+) = уг + Уп + Уп + У г,

7<-) = У# + У" + у# + у".

В области а рассмотрим сеточные функции

у== у(\, у^ =мх(±1з-аЬ

и(±) = и( Х^3-а )).

Дискретизацию задачи проведем с помощью интегро-интерполяционного метода, который приводит к консервативной разностной схеме со вторым порядком аппроксимации по обеим пространственным переменным [3]:

Аи = / ,

где оператор А = А + А2, / = / + /2:

12((^<а+ % ) Х + (аа У г ) Х )> х ’

2 ха ха ха ха

2(Г(а«+)Уха -«<+)У) + ^-(а?Уха -у)), х

Ь1 (-«<+> у - «<+> ух ) + К-(-а<-> у - а? Ух )), 2 Па К

1Кг(а<‘+)Уха -а<+у')’ хе г#>

а

2т(а°}Уха)-а<-у)’х’

2К-(-«<+)У-а(а+ Ух )’ х

2 П а Х а

\Ь(-а( )У-№Ух )’ хе ^’

2 К а ха

2(К-(а*}Уха -У) + К+ Ух )х.), хе 7#,

2 К а а Ла ха

2(-Г-(аа+ Ух - «<+) У) + (аа} Ух )„)’ х еУ"’

2 Ка а Ха ха

2((а«+)Ух ) х + К~(-а< )У - а<а-Ух ))’ хе У!,

2 Ха Ха Па Ха

2((аа} УХ ) Х + Т-(-а<Л> У - ^ УХ ))’ Х еУг ’

2 Ха ха К Ха

АаУ =

/а =

О, х є сок,

2 Т {-а^ + а^ Т. х єУа, 2 к-( «£> + С%> )Тср , х єуа,

——а

2 к а ср Тср ’ х ' ’

ІГ ^ . Х ЄГ-

Систему (3) будем решать методом сопряженных градиентов [2]:

Вуы =ак+1(В-тк+1А)ук + (1 -а^Ву1-1 + ак+1тк+1/, к = 1,2,...,

Ву1 = (В -Ті А) у0 +т/,

где ак+1 и Тк+1 - итерационные параметры. Вычисления по (4) основаны на представлении

у*+1 = ак+і ук +(1 -ак+1) у к-1 + ак+1тк+1ак.

Здесь ($к = В ~ггк.

Итерационные параметры рассчитываются по формулам:

(ок, гк) , П1

Тк+і = Л к , к = 0Д,...,

(АО,О)!

ак+1 = (1

(°, гк) 1_)-^ к = 0,1,2,...,

тк (О 1, гк 1) ак

а1 =1.

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

2. Численные эксперименты

Проведены численные эксперименты для четырех видов наиболее распространенных типов наружных ограждающих конструкций:

2.1. Однородная деревянная стена (точное решение известно).

2.2. Неоднородная стена с большими перепадами коэффициентов теплопроводности (точное решение известно).

2.3. Неоднородная стена с теплопроводным включением.

2.4. Угловое соединение.

Введем обозначения:

е =| Тк-Тк-11 - разность между соседними итерациями Тк и Тк-1 в равномерной норме, к - номер итераций, ДР=||Р„ар-Рвнутр ||- разность между

тепловыми потоками через наружную поверхность

ограждения, 1 - время расчета, Тнар - вычисленная

температура на наружной границе ограждения, Твнутр

- вычисленная температура на внутренней границе ограждения.

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

2.1. Рассматривалась прямая вертикальная деревянная стена толщиной 0,2 метра, длиной 1 метр и коэффициентом теплопроводности Х=0,14 Вт/(м^°С). Точное решение легко вычисляется, температуры на наружной и внутренней границах соответственно

равны Т =-51,945 °С, Т =15,568 °С. Результаты

А нар 7 7 внутр 7 ^

расчетов приведем в табл. 1-5.

Таблица 1

Результаты расчета методом Зейделя для е=10-3

сетка к дд і Т нар Т внутр

10Ч50 163 7.20^10-2 <1 с -51.943 15.572

20Ч100 533 3.5810-1 <1 с -51.937 15.585

40Ч200 1579 1.5110+0 3 с -51.911 15.648

100Ч500 5313 9.66^10+0 1 мин. -51.730 16.104

200Ч1000 8795 3.44^10+1 6 мин.30 с -51.078 17.222

Т

Таблица 2

Результаты расчета методом Зейделя для е=10-4

сетка к дд 1 Т нар Т внутр

10Ч50 227 7.58-10'3 <1 с -51.945 15.568

20Ч100 765 3.50-10-2 <1 с -51.944 15.570

40Ч200 2502 1.49 10-1 4 с -51.942 15.576

100Ч500 11071 9.60-10-1 2 мин. -51.924 15.621

200Ч1000 30471 3.87-10+° 23 мин. -51.860 15.786

Таблица 3

Результаты расчета методом Зейделя для е=10-5

сетка к дд 1 Т нар Т внутр

10Ч50 285 7.41-10"4 <1 с -51.945 15.568

20Ч100 996 3.45-10-3 <1 с -51.945 15.568

40Ч200 3425 1.4810-2 6 с -51.945 15.569

100Ч500 16825 9.56-10'2 3 мин. -51.943 15.573

200Ч1000 53471 3.87-104 40 мин. -51.937 15.590

Таблица 4

Результаты расчета методом Зейделя для е=10-6

Как видно из таблиц, для метода Зейделя е = 10-3 недостаточно для получения точных результатов на любой сетке. Тем не менее, он достаточно хорошо работает для е=10-6, хотя не следует увлекаться слишком мелким разбиением области расчета. Результаты табл. 5 свидетельствуют о большой скорости сходимости метода сопряженных градиентов, когда уже при е =10-3 достигается точное решение и число итераций линейно зависит от числа узлов сетки.

2.2. Рассмотрим температурное поле ограждения, которое состоит из плоского пенополистирола (Х2=0,041 Вт/(м^°С)) толщиной 10 сантиметров,

облицованного с обеих сторон алюминиевой пластиной (Х1=221 Вт/(м^°С)) толщиной 2 миллиметра.

Перепад значений коэффициентов теплопроводности составляет 1/^=221/0,041=5390,24. Точное решение также легко вычисляется, температуры на наружной и внутренней границах соответственно равны Т =-52,745 °С, Т =17,681°С. Сетка фиксирована:

нар 7 7 внутр 7 'гг

104^10, то есть шаг сетки по оси Ох И =0,001, шаг

7 х 7 7

сетки по оси Оу Ь=0,01.

Результаты расчетов приведены в табл. 6-7.

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

Таблица 6

Результаты расчета методом Зейделя

сетка к дд 1 Т нар Т внутр

10Ч50 343 7.25-10-5 <1 с -51.945 15.568

20Ч100 1227 3.4110-4 1 с -51.945 15.568

40Ч200 4347 1.46 10-3 8 с -51.945 15.568

100Ч500 22577 9.5310-3 4 мин. 10 с -51.945 15.568

200Ч1000 76464 3.86-10'2 56 мин. -51.944 15.570

е к дд 1 Т нар Т внутр

10-3 19451 3.36-10+1 5 с -33.822 6.288

10-5 171502 4.46-10-1 46 с -52.752 17.189

10-7 398143 4.62-10-3 1 мин. 45 с -52.745 17.676

10-8 511490 4.62-10-4 2 мин. 15 с -52.745 17.681

Таблица 5

Результаты расчета методом сопряженных градиентов для е=10-3, 10-4, 10-5, 10-6

Таблица 7

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

сетка к дд 1 Т нар Т внутр

10Ч50 12 8.86-10-12 <1 с -51.945 15.568

20Ч100 23 1,26'10-12 <1 с -51.945 15.568

40Ч200 43 1.06-10-10 <1 с -51.945 15.568

100Ч500 103 1.59-10-10 2 с -51.945 15.568

200Ч1000 203 1.83 10-9 15 с -51.945 15.568

е к дд 1 Т нар Т внутр

10-3 106 5.5110-6 <1 с -52.745 17.681

10-5 106 5.5110-6 <1 с -52.745 17.681

10-7 106 5.5110-6 <1 с -52.745 17.681

10-8 106 5.5110-6 <1 с -52.745 17.681

*1 ^ ^1 0,002 1 ! 0,20

^3 0,19 1

/ ! 0,01

К ТИ(=21°С

Рис. 2. Ограждение с теплопроводным включением

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

В горизонтальном слое из пенополистирола (Х2=0,041 Вт/(м^°С)) имеется тонкий вертикальный металлический стержень (Х1=58 Вт/(м^°С)). Точное решение неизвестно, о результатах расчета можно судить лишь по качеству температурного поля. Размеры ограждения указаны на рис. 2. Сетка в расчетах бралась постоянной: 42*40.

Введем обозначения:

Т - максимальная температура наружной по-

макс.нар. 1 ^ 1 1 ^

верхности, которая, очевидно, находится посередине наружной поверхности, Т - минимальная темпе-

1 ^ 1 у мин.вн.

ратура внутренней поверхности, которая, очевидно, находится посередине внутренней поверхности ограждения.

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

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

Приведем результаты расчета для прямоугольного двумерного угла из сосны (X = 0,14 Вт/(м^°С)) толщиной 0,2 м и длиной 1 м по наружной поверхности. Данный пример интересен тем, что расчетная область имеет внутренний угол больший, чем 180 градусов.

Введем обозначения:

Т - вычисленная температура внутреннего угла ограждения. Сетка равномерная и состоит из NxN узлов.

Результаты расчетов приведены в табл. 10-12.

По результатам экспериментов видно, что расчет угловых соединений представляет определенные трудности. Это объясняется большими вторыми производными вблизи внутренних угловых узлов [4], которые в данном случае являются величинами порядка 3ДТ0+5. Поэтому, если требуется большая точность температур внутренних угловых точек, то необходимо внимательно следить за размерами сетки, точностью соседних итераций, разностью тепловых потоков через внутреннюю и наружную поверхности ограждающих конструкций. Естественно, что время расчета, погрешность приближенного решения будут возрастать с усложнением геометрии

Таблица 8

Результаты расчета методом Зейделя

є к Ґ Т мин.вн. Т макс.нар.

10-3 3608 8.37'10+0 2 с 18.447 -7.305

10-5 122598 2.28-10-1 49 с 17.642 -20.845

10-7 272234 2.27-10-3 1 мин. 50 с 17.615 -21.144

10-8 347052 2.27-10-4 2 мин. 19 с 17.615 -21.147

Таблица 9

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

є к Ґ Т мин.вн. Т макс.нар.

10-3 502 2.42-10-2 <1 с 17.613 -21.147

10-5 560 4.1410-6 1 с 17.615 -21.147

10-7 608 3.72-10-8 1 с 17.615 -21.147

10-8 631 1.77-10-8 1 с 17.615 -21.147

Таблица 10

Сравнение результатов расчета для е=10-3

Метод Зейделя Метод сопряженных градиентов

N к дд X Т угл к дд X Т угл

50 192 1.610-1 <1 с 8.306 63 2.3^10-3 <1 с 8.298

100 624 7.110-1 1 с 8.002 107 3.310-3 1 с 7.967

500 7204 9.810+0 6 мин. 8.721 428 1.7^10-2 24 с 7.770

1000 14654 3.9-10+1 48 мин. 11.330 791 3.5^10-2 3 мин. 7.765

2000 19013 1.110+2 4 часа 20 мин. 13.912 1397 6.7^10-2 20 мин. 7.779

4000 - - - - 2279 9.310-2 2 часа 8 мин. 7.772

Таблица 11

Сравнение результатов расчета для е=10-6

Метод Зейделя Метод сопряженных градиентов

N к дд X Т угл к дд X Т угл

50 375 1.810-4 1 с 8.298 101 4.2^10-7 <1 с 8.298

100 1358 6.4^10-4 3 с 7.966 173 8.110"7 1 с 7.966

500 25516 8.610-3 21 мин. 7.767 702 1.6 10-5 39 с 7.766

1000 87404 3.5^10-2 4 часа 45 мин. 7.759 1369 1.7^10-5 5 мин. 7.755

2000 290847 1.510-1 75 часов 7.767 2661 2.3^10-5 35 мин. 7.751

4000 - - - - 5114 4.110-5 4 часа 32 мин. 7.750

Таблица 12

Сравнение результатов расчета для е=10-8

Метод Зейделя Метод сопряженных градиентов

N к дд X Т угл К дд X Т угл

50 496 6.810-7 1 с 8.298 116 7.4^10-9 <1 с 8.298

100 1845 2.9^10-6 4 с 7.966 197 1.610-8 1 с 7.966

500 37715 8.0^10-5 31 мин. 7.766 853 1.5^10-8 48 с 7.766

1000 136227 3.310-4 7 часов 30 мин. 7.755 1665 6.510-8 6 мин. 7.755

2000 - - - - 3276 3.110-7 47 мин. 7.751

4000 - - - - 6412 5.410-7 5 час 43 мин. 7.750

области расчета, наличием неоднородных материалов.

Данная работа нашла применение в совершенствовании программы «8НАВБА^2Б». Новая версия программы существенно улучшена не только в отношении совершенствования математического аппарата и методов решения, но и в отношении удобства пользования программы с учетом пожеланий и замечаний строителей-проектировщиков.

Л и т е р а т у р а

1. Данилов Н. Д., Шадрин В. Ю., Павлов Н. Н.

Прогнозирование температурного режима угловых соединений наружных ограждающих конструкций // Промышленное и гражданское строительство. - 2010. - № 4. С. 20-22.

2. Самарский А. А., Вабищевич П. Н. Численные методы решения задач конвекции - диффузии. - М.: Эдито-риал УРСС, 1999. - 248 с.

3. Самарский А. А. Теория разностных схем. - М.: Наука. Гл. ред. физ.-мат. лит., 1977. - 656 с.

4. Шайдуров В. В. Многосеточные методы конечных элементов. - М.: Наука. Гл. ред. физ.-мат. лит., 1989. -288 с.