CC BY
83
МАТЕМАТИЧЕСКАЯ МОДЕЛЬ СВОБОДНОЙ КОНВЕКЦИИ ВОЗДУХА В КОМНАТЕ
В.А. КУЗНЕЦОВ, В.П. КОЖЕВНИКОВ
Белгородский государственный технологический университет им. В.Г. Шухова
Дифференциальные уравнения Навье-Стокса представлены в форме, учитывающей влияние силы тяжести неизотермической среды на каждую компоненту скорости. Предложен способ расчета слабой турбулентности при естественной конвекции. Составлена трехмерная математическая модель, разработана компьютерная программа и выполнено численное моделирование свободной конвекции воздуха в обогреваемой комнате, обладающей диагональной симметрией. Представлены распределения температуры и скорости воздуха в помещении.
Для решения проблемы повышения энергетической эффективности зданий необходимо развивать методы математического моделирования теплообмена при свободной конвекции воздуха. Математическую модель составляют на основе дифференциальных уравнений движения (неразрывности и Навье-Стокса) и конвективного переноса теплоты, полагая, что с помощью их численного решения возможно найти распределение температуры *, давления р и компонент скорости и, V, м> в декартовых координатах х, у, z. Но обычно задачу упрощают, переходя к преобразованным переменным - функции тока и завихренности, что приводит к ограничению применения математических моделей областью двухмерных задач, делает их малоинформативными.
Одной из причин, затрудняющих решение трехмерных задач свободной конвекции, является особенность структуры уравнений Навье-Стокса. Сила тяжести входит только в уравнение для вертикальной компоненты скорости V, а на горизонтальные компоненты скорости и и сила тяжести оказывает влияние лишь косвенно посредством уравнения неразрывности, связывающего между собой все три компоненты скорости. Такая отдаленная взаимосвязь функций препятствует нормальной сходимости итераций при численном решении задачи.
Другой причиной может быть малая вязкость воздуха при комнатной температуре, приводящая к появлению турбулентности [1]. Встает естественный вопрос, почему известные методы расчета турбулентности не используются в математических моделях свободной конвекции воздуха в помещениях? Возможно, ответ заключается в том, что пока еще не вполне ясно, как рассчитывать слабую турбулентность вблизи стенки.
В численных исследованиях процессов вынужденной конвекции применяется двухпараметрическая диссипативная модель турбулентности, позволяющая вычислять турбулентный аналог кинематического коэффициента вязкости ут в соответствии с полуэмпирической гипотезой Колмогорова:
Введение
© В.А. Кузнецов, В.П. Кожевников Проблемы энергетики, 2008, № 7-8
где кт - кинетическая энергия турбулентных пульсаций; гт - скорость
диссипации энергии турбулентности; с ^ - эмпирический коэффициент (с ^ » 0,09);
fv^ - поправочная функция для слаборазвитой турбулентности вблизи стенки [2].
Для расчета турбулентности вблизи стенки используют также гипотезу Прандтля, из которой следует безразмерное квадратное уравнение:
Ґ \2
ди
ди+
дУ+
+ —--1 = 0,
ду+
где и + = и / и* - безразмерная осредненная скорость, параллельная поверхности стенки; у+ = уи* /V - безразмерное расстояние от стенки до расчетной точки;
l + = Ы* /V - безразмерная длина пути смешения; V - кинематический коэффициент вязкости; и* - динамическая скорость, характеризующая величину касательного напряжения на стенке.
Его решение дает безразмерные выражения для производной от скорости
(2)
Зу + ^ л
и для турбулентной вязкости вблизи стенки
— = -12. (3)
Формулы (2) и (3) становятся пригодными для расчета турбулентного движения только в том случае, если известно изменение длины пути смешения ^ по нормали к стенке. При расчете параметров развитой турбулентности обычно учитывают, что вблизи поверхности стенки в турбулентном потоке образуется равновесный слой, в котором, в соответствии с гипотезой Прандтля, длина пути смешения прямо пропорциональна расстоянию до стенки у + [2]:
I + = КУ + ,
где к - универсальная постоянная, имеющая значение 0,41 [2].
В случае слабой турбулентности целесообразно использовать эмпирическую формулу, определяющую безразмерную длину пути смешения в вязком подслое и переходной области (при у+< 60), где турбулентные пульсации заглушаются, при параметре A+ = 30 [3]:
l + = ку + { 1 - {(- у + / (+ - 0,25 у + / А.+ ))}. (4)
В итоге появляются предпосылки для учета слабой турбулентности. Преобразование дифференциальных уравнений движения
Известно [4], что при постоянной плотности р слагаемое, представляющее силу тяжести в уравнениях Навье-Стокса, можно заменить на производную от
функции ряк, где к - координата по вертикальной оси у, а затем эту функцию объединить с давлением р, введя для этой суммы, например, обозначение р = Р + ряк.
Чтобы в неизотермических условиях при переменной плотности воздуха производные от давления р в уравнениях Навье-Стокса заменить проекциями градиента функции давления р, следует производные от ряк записать по правилам дифференцирования произведения функций. В результате в дифференциальные уравнения Навье-Стокса вместе с производными от функции ~ войдут также и производные от переменной плотности воздуха [5]:
дрии дpvu дрми д
+ . +
дх ду дг
др^ дрvv др^
+ +
дх ду дг
дрим др^ дрмм
- + - +
дх ду дг
( ди' д ( Л ди д ( ди І
Д— Д— Д—
V дх Ч ду V дУ V дг V дг >
др др
=---------+ — Як,
дх дх
( дv" д ( \ дv д Ґ \ дv Р
д— д— д— = + — як
хд ду V дУ у дг V дг у ду ду
г \ дм д ( N дм д " дм'
д— д— Д—
хд ду V ду > дг V дг >
др Эр
=-----+ — Як ,
дг дг
где ц - динамический коэффициент вязкости воздуха; я - ускорение силы тяжести.
В этом варианте записи системы дифференциальных уравнений Навье-Стокса для стационарной задачи влияние силы тяжести на движение воздуха учитывается по каждой из трех координат. Влияние переменной плотности непосредственно передается каждой компоненте скорости, что значительно улучшает сходимость итераций.
Математическое моделирование слабой турбулентности
Для нахождения турбулентной вязкости по формуле (1) во внутренних точках расчетной области стационарной задачи решают численно дифференциальные уравнения переноса кинетической энергии турбулентности к т :
дрик т дрvkт дрмк т -------+---------+--------
дх
ду
дг
д
дх
ґ
р(V+vт)-
дк т
дх
+
д ( д і д ( дк т N
+ — р Г—\ < + V т )— + р( V + -'Ч ) р +
ду V ду J дг V дг
и скорости диссипации энергии турбулентности гт
( ( л Л
дриг т дрvг т дрмгт д 'т дгт
— + + — = р V + — +
дх ду дг дх V V У дх ,
+ -
ду
V + -
д£т
д + —
дг
Р
V + -
д£т
дг
+ Р (£ 1 О т с£ 2 £ т )5
к т
где эмпирические постоянные имеют значения: аЕ » 1,3; се1 » 1,44 и се2 » 1,92 [2] и учитывается скорость генерации энергии турбулентности От.
Граничные условия к этим уравнениям целесообразно формулировать в узлах сетки, ближайших к стенам, с помощью полуэмпирической гипотезы Прандтля. Численно интегрируя уравнение (2) совместно с экспоненциальным выражением длины пути смешения (4), можно установить уточненное распределение безразмерной скорости и + в пристенной области (кривая 2 на рис. 1), которая при у+< 40 лежит значительно ниже логарифмической кривой 1, соответствующей линейной зависимости длины пути смешения I + в гипотезе Прандтля.
Рис. 1. Безразмерная зависимость скорости и+ от расстояния у+ по закону:
1 - логарифмическому; 2 - уточненному; 3 - хорда; 4 - касательная
Динамическую скорость и* предлагается вычислять в процессе итераций как величину, заменяющую касательное напряжение на стенке:
и*2 = V (ди/ду )у=0.
Производная в этом выражении равна тангенсу угла наклона касательной в начале координат к графической зависимости скорости и от расстояния у. Обычно при дискретном представлении производной касательная приближенно заменяется хордой, что вносит в расчет заметную погрешность. Так, например, в безразмерных координатах на рис. 1 угол наклона хорды 3 заметно меньше, чем касательной 4. Для устранения расчетной погрешности достаточно дискретное значение производной и/у поделить на соотношение безразмерных величин и + / у+, найденное по зависимостям рис. 1, так как тангенс угла наклона касательной 4 определяется соотношением и + / у+, равным единице:
д
V
V
т
т
Р
о
о
£ У
£ У
2 иУ +
и* = V----------
У и+
Уравнение (5) пригодно для расчета любого режима течения среды вблизи стенки: и ламинарного, и турбулентного, и переходного. Единственным
ограничением его применимости является предположение о том, что в расчетных условиях сохраняется классическая структура турбулентного пограничного слоя. По найденному значению динамической скорости и* уточняется
безразмерное расстояние у + до ближайшего к стенке узла сетки, вычисляется турбулентная вязкость vт по уравнению (3), а затем и кинетическая энергия турбулентности к т на границе :
В математической модели принято, что в равновесном слое поправочный множитель /к = 1 и, в соответствии с экспериментальными данными [2], по мере приближения к границе вязкого подслоя повышается до /к * 1,2, а затем уменьшается до нуля на поверхности стенки.
Граничное значение скорости диссипации энергии турбулентности можно рассчитать по формуле, следующей из зависимости (1):
Как показал расчет, граничные значения вт оказывают заметное влияние на вычисляемые по формуле (1) значения коэффициента турбулентной вязкости вдали от стен, задавая таким образом уровень турбулентности в центральной части расчетной области. В связи с этим возрастает роль поправочной функции fv^ в формуле (7).
Из анализа выражения (1) следует, что поправочная функция (кривая 1 на рис. 2) вблизи стенки равна соотношению значений турбулентной вязкости, рассчитанных по уравнению (3) при длине пути смешения, определяемой разными зависимостями - экспоненциальной (4) и линейной:
Чтобы поднять расчетный коэффициент турбулентной вязкости во внутренней области до уровня, обеспечивающего сходимость итераций, значения поправочной функции при численном решении задачи пришлось несколько уменьшить (см. кривую 2 на рис. 2). Для сравнения представлена также эмпирическая поправочная функция [2] (кривая 3 на рис. 2).
(6)
(7)
(8)
Рис. 2. Поправочная функция /ц: 1 - уточненная; 2 - расчетная; 3 - Лаундер-Шарма [2]
Математическая модель свободной конвекции воздуха
Математическая модель стационарного движения воздуха состоит из преобразованных дифференциальных уравнений Навье-Стокса, в коэффициентах которых к молекулярной вязкости ц добавлена турбулентная вязкость у,,., и дифференциального уравнения неразрывности
dpu dpv dpw
---+-----+-----= 0.
dx dy dz
Турбулентная вязкость в комнате вычисляется по формуле (1) везде, кроме ближайших к стенкам узлов сетки, для чего численно решаются дифференциальные уравнения переноса кинетической энергии турбулентности кт
и скорости диссипации энергии турбулентности £т. Для постановки граничных условий к этим уравнениям в ближайших к стенам узлах сетки использованы выражения (3) - (8).
При переходе к дискретным алгебраическим аналогам дифференциальное уравнение неразрывности преобразуется в уравнение для расчета поля давления по схеме метода SIMPLER [4]. Но, в отличие от этого метода, в математической модели поправка давления заменена корректирующими выражениями с расчетными источниками массы. После их подстановки в дискретные уравнения движения и небольших преобразований, разностная схема становится полностью неявной [6] и приобретает свойство безусловной устойчивости. В итоге сходимость итераций улучшается.
Построение разнесенной (шахматной) сетки производится по рекомендациям метода SIMPLER так, чтобы избежать вычисления граничных значений поля давления. Для однозначности решения задачи необходимо в одном из узлов сетки задать значение преобразованной функции давления, например p =0. Скорость воздуха на поверхности стен принимается равной нулю, а граничные коэффициенты переноса для компонент скорости, параллельных стенам, определяются по вычисленным значениям динамической скорости и*.
Математическая модель учитывает, что побудителем движения в комнате является конвективный нагрев воздуха обогревающими устройствами при
одновременном отводе теплоты через стены и окна. Стационарное температурное поле во внутренних узлах сетки рассчитывается с помощью дифференциального уравнения переноса теплоты конвекцией и теплопроводностью без учета теплового излучения:
дрисі дргсі дрп>а д
-------+--------+--------= —
дх
ду
ді дх
ГГ > \
ді
к + +
V Рг V X 1 т ) дх,
д + Ґ Ґ \ к + ді д + — ґ Ґ \ , СДт к + \ ді
ду V Рг V ГІ т ) ду, ді V Рг V X 1 т ) ді ,
где * - температура; с - удельная теплоемкость воздуха при постоянном давлении; к - коэффициент теплопроводности воздуха; ц - турбулентный аналог динамического коэффициента вязкости воздуха; Ргт - турбулентное число Прандтля для воздуха.
Для вычисления плотности теплового потока теплоотдачи к внутренней поверхности стен, окон и нагревателей используется аналогия интенсивностей переноса количества движения и теплоты:
*
Ч гр = к (д*1 % )у=о = Рс—№/&
Рг
где Рг - число Прандтля воздуха.
Аналогично формуле (5), производная в этом выражении вычисляется по тангенсу угла наклона хорды в графической зависимости разности температуры воздуха в пограничном слое и температуры на поверхности нагревателя, стены или окна (* — *гр) от расстояния у и корректируется делением на соотношение
Т+ /(у + Рг), равное в начале координат единице:
і — і
ч гр = Рсу-
гр
у+
(9)
У Т+
где Т+ = ри* (* — *гр)/Чгр - безразмерная температура, зависимость которой от безразмерного расстояния у + определяется формулой [7]
дТ+
ду+
( 2 '\—1
1 1+ ди+
— +------------
Рг Ргт ду + )
от
Этой формуле придадим следующий вид, показывающий, что производную безразмерной температуры можно, вообще говоря, заменить
скорректированной производной от безразмерной скорости:
дТ+
= Ргт
д + і Ґ Ргт т 1 1+ / Ргт ..1і т ( \ Ргт т 1 >+
.ду+ г Р / ] / 1 Рг г Р / \
ду+
Действительно, как следует из экспериментальных данных [7], для турбулентного движения воздуха вблизи стенки можно принять, что за пределами вязкого подслоя безразмерные производные от температуры и скорости пропорциональны (с погрешностью порядка 5 % при Ргт и 0,75) и, следовательно,
-» Ргт
(10)
у+
у+
Подставив соотношения (10) и (5) в формулу (9), получим выражение для граничных условий теплоотдачи в математической модели:
.2
( - *гр ),
реи
9 гр
(11)
Ргт и
где * - температура в ближайшем к стенке узле сетки; и - скорость в этом же узле сетки; * гр - температура поверхности стены, окна или нагревателя.
Температура поверхности окон и наружных стен комнаты находится из условия равенства плотности локального теплового потока в выражении (11) и в известной формуле теплопередачи через плоскую многослойную стенку. Температура поверхности внутренних стен, пола и потолка комнаты принимается равной температуре воздуха в ближайших к ним узлах сетки, так что теплоотдача к этим поверхностям не рассчитывается.
Обсуждение численных результатов
Цель вычислений состояла в оценке адекватности математической модели свободной конвекции воздуха, а не в исследовании особенностей обогрева реальных объектов. Это объясняет упрощенную геометрическую постановку задачи и исключение радиационного теплообмена из математической модели.
Численный эксперимент выполнен применительно к квадратной комнате длиной и шириной 5 м, высотой 2,4 м. Комната занимает угловое положение в здании и на каждой из двух ее взаимно перпендикулярных наружных стен расположен нагреватель высотой 0,6 м и над ним - окно высотой 1 м. Ширина окна и нагревателя 3 м. Нагреватель и окно считались плоскими и составляли часть внутренней поверхности стены. Таким образом, задаче придана симметрия относительно вертикальной плоскости, пересекающей комнату по диагонали.
Наружные стены, выполненные из глиняного кирпича толщиной 250 мм, утеплены изнутри слоем пенобетона толщиной 120 мм и создают термическое сопротивление 0,88 м2К/Вт. Для окон задавалось термическое сопротивление
0,464 м2К/Вт. Температура нагревателя принята 60 °С, температура атмосферного воздуха -20 °С при отсутствии ветра.
Сетка по оси х, перпендикулярной передней стене, и по оси г, перпендикулярной боковой стене, выполнена с переменным шагом, который, по мере удаления от наружных стен, возрастает от 10 мм до 100 мм. По вертикальной оси координат у шаг сетки принимается постоянным, равным 50 мм. В итоге весь объем комнаты разделен по осям х и г на 55 слоев, по оси у - на 48 слоев. Пересечение слоев
и
+
+
создает контрольные объемы, в центре которых находятся узлы сетки. В них определяются значения искомых функций (температура, давление) и физические свойства воздуха (теплоемкость, вязкость, плотность). В узлах сетки, лежащих на гранях контрольных объемов, определяются компоненты скорости.
Численное решение систем дискретных уравнений математической модели выполнялось полинейной прогонкой с итерациями. Сходимость итераций контролировалась по абсолютной величине расчетных источников массы. Расчет завершался, когда суммарный расчетный источник массы снижался до 10 мг/м3.
Полученные результаты рассматриваются здесь, прежде всего, с точки зрения их достоверности. Предполагается, что отклонения полученных результатов от диагональной симметрии будет указывать на существование определенных погрешностей в решении задачи.
На рис. 3 представлены кривые, показывающие изменение температуры воздуха на разных уровнях в средней части комнаты, по оси х - перпендикулярно передней стене и по оси і - перпендикулярно боковой стене. Вдали от наружных стен комнаты температура воздуха повышается в направлении от пола к потолку, причем на каждом горизонтальном уровне она изменяется незначительно. Резкое изменение температуры воздуха происходит в довольно тонком пристенном слое, толщиной порядка 10 мм, прилегающем к наружным стенам комнаты.
и
О
Я
6
я
I
г
с
г
и
о
я
І.
■
и
с
ш
л
й
5
£■
к
г
н
20 40 60 80 100
Расстояние от стены, мм
а)
1000 2000 3000 4000 5000
Расстояние от стены, мм
б)
Рис. 3. Температура в вертикальной плоскости посередине комнаты по осям х (сплошные линии) и г (пунктирные линии) вблизи наружных стен и вдали от них: 1 - на 25 мм выше пола;
2, 3 - на уровне середины нагревателя и окна; 4 - на 25 мм ниже потолка
Температура воздуха над полом (кривая 1) является наиболее низкой. На внутренней поверхности стены она составляет 6,7 °С, а вдали от стен возрастает до 16,2 °С. Температура на уровне нагревателя (кривая 2) уменьшается от 60 °С на его поверхности до 28,4 °С на расстоянии 5 мм от него и затем до 16,8 °С вдали от стен.
На уровне окна в восходящем потоке нагретого воздуха температура сначала повышается от 10,2 °С на его внутренней поверхности до 19,9 °С на расстоянии 17...33 мм от окна, затем снижается до 18...18,4 °С вдали от стен (кривая 5). Температура под потолком (кривая 4) сначала ниже, чем на уровне окна, а вдали от стен становится выше, достигая 18,5 °С. Характерно, что температурные кривые по осям х и г почти везде совпадают, что подтверждает неплохую диагональную симметрию температурных полей в комнате.
Распределение вертикальной скорости воздуха характеризуется более заметными локальными отклонениями от диагональной симметрии вблизи наружных стен. Так, согласно рис. 4, на расстоянии 5 мм от передней стены скорость составляет 0,44 м/с, в то время как на таком же расстоянии от боковой стены только
0,37 м/с. Следовательно, локальная погрешность вычисления скорости воздуха достигает здесь 20 %. Вдали от наружных стен вертикальная скорость, как и следовало ожидать, близка к нулю.
0,5
0,4
0,3
0,2
0,1
0,1
V. /
2
1000 2000 3000 4000 5000
Расстояние от стены, им Рис. 4. Изменение вертикальной компоненты скорости V посередине комнаты (удаление от соседних стен 2,5 м, от пола 1,2 м) вдоль осей: 1 - х; 2 - г
Более наглядно локальное отклонение от диагональной симметрии показано на рис. 5, где представлена вертикальная компонента скорости воздуха на уровне середины нагревателя в плоскостях, отстоящих на 5 мм от боковой и передней стен. Рис. 5 дает представление о том, что практически по всей ширине нагревателя существует подъемный поток нагретого воздуха, который способствует повышению температуры вблизи расположенного над ним окна. Левее и правее нагревателя возникает опускное движение воздуха, причем максимальная скорость направленного вниз потока составляет 0,6 м/с в углу комнаты, образованном наружными стенами. Более слабое опускное движение образуется под нагревателем.
е.
с
2
X
в
н
а.
«
0,4
0,3
0.2
0,1
о
0,1
0,2
0,3
0,4
0,5
0.6
1
1 2 - "-А
■'
/ 1/
1
1000 2000 3000 4000 5000
Расстояние вдоль стены, мм Рис. 5. Вертикальное движение воздуха на уровне нагревателя на расстоянии 5 мм: 1 - от боковой стены; 2 - от передней стены
Этот нисходящий поток воздуха встречает на своем пути пол комнаты, поворачивает и продолжает свое движение в горизонтальной плоскости непосредственно над полом. Кривые на рис. 6 показывают, что горизонтальная скорость воздуха над полом может достигать 0,2 м/с. В то же время распределение горизонтальных компонент скорости в комнате в целом подчиняется диагональной симметрии.
Рис. 6. Изменение горизонтальной скорости в плоскости на 25 мм выше пола:
1 - и по оси г; 2 - по оси х;
На рис. 7 отложены расчетные значения турбулентного аналога кинематического коэффициента вязкости. Хотя кривые 1 и 2, представляющие его значения по осям х и г, в своей средней части не вполне совпадают друг с другом, все же можно утверждать, что и здесь приближенно воссоздается диагональная симметрия комнаты.
0,003
Гч|
Е
А
В
О
2
ва
9
X
н
X
ф
£
о
£
0,002
0,001
[
'
0 1000 2000 3000 4000 5000
Расстояние от стены, ми
Рис. 7. Изменение турбулентного аналога кинематического коэффициента вязкости ут посередине комнаты: 1 - по оси х; 2 - по оси г
Как отмечалось ранее, расчетные значения коэффициента турбулентной вязкости вдали от стен несколько завышены из-за введения в математическую модель скорректированной поправочной функции f. Заметные расхождения
кривых на рис. 7 также возникают вдали от стен комнаты. Однако отмеченные неточности не приводят к сколько-нибудь заметному искажению распределения температуры и скорости воздуха. Таким образом, в какой-то мере подтверждается положение, принятое в работе [8], что в основном объеме помещения влиянием турбулентной вязкости можно пренебречь.
Выводы
Решены две основные проблемы, препятствующие в настоящее время численному исследованию трехмерной свободной конвекции воздуха при обогреве помещений: улучшена форма представления силы тяжести в дифференциальных уравнениях Навье-Стокса и корректно учтены слабые проявления турбулентности.
Получены формулы для вычисления динамической скорости на стенке и уточнены граничные условия в математических моделях турбулентности и теплоотдачи к воздуху. В итоге построена трехмерная математическая модель свободной конвекции воздуха в обогреваемой комнате, разработана компьютерная программа и обеспечена сходимость итераций при ее численной реализации.
Представленные распределения температуры и скорости воздуха в комнате, обладающей диагональной симметрией, позволяют оценить точность численного эксперимента и свидетельствуют о приемлемой достоверности полученных результатов.
Summary
Navier-Stokes differential equations are presented in a form that takes into account an influence of the nonisothermal medium gravity force directly upon each velocity component. A method is suggested to compute any small turbulence while there is some natural convection. New computing equations for skin-friction velocity have been derived, and turbulence and heat transfer boundary conditions are defined more precisely. A threedimensional mathematical model and a computer program have been elaborated, and a numerical simulation has been realized for the natural convection in a heated room with some diagonal symmetry. Temperature and velocity distributions in the room are adduced.
Литература
1. Математическое моделирование конвективного тепломассообмена на основе уравнений Навье-Стокса / В.И. Полежаев, А.В. Бунэ, Н.А. Верезуб и др. -М.: Наука, 1987. - 270 с.
1. Пейтел В.К., Роди В., Шойерер Г. Модели турбулентности для течений в пристеночной области с малыми числами Рейнольдса // Аэрокосмическая техника. - 1986. - № 2. - С. 184-197.
2. Кузнецов В.А. Уточнение гипотез пристенной турбулентности // ИФЖ. -1986. - Т. 50. - № 6. - С. 917-922.
3. Ferziger J.H., Peric’M. Computational methods for fluid dynamics. - Berlin, Heidelberg, New York: Springer, 2003. - 423 p.
4. Кузнецов В.А., Наренкова Я.В. Математическая модель свободной конвекции воздуха в помещениях // Научные исследования, наносистемы и ресурсосберегающие технологии в стройиндустрии. - БГТУ им. В.Г. Шухова. -2007. - Ч. 8. - С. 72-74.
5. Кузнецов В.А. Основы математического моделирования
теплотехнологических процессов. - Белгород: Изд-во БГТУ им. В.Г. Шухова, 2004. - 95 с.
7. Жукаускас А. А. Конвективный перенос в теплообменниках. - М.: Наука, 1982. - 472 с.
8. Табунщиков Ю.А., Бродач М.М. Математическое моделирование и оптимизация тепловой эффективности зданий. - М.: АВОК-ПРЕСС, 2002. - 194 с.
Поступила 01.04.2008
