CN103455676A - 一种利用流体力学模拟室内热环境的方法 - Google Patents

Abstract

本发明公开了一种利用流体力学模拟室内热环境的方法，包括以下步骤：根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导空调工程施工。本发明通过流体力学的计算机仿真，具有模拟真实条件的能力，完善模型以及边界条件保证模拟的准确性，能够较全面的了解不同边界条件下室内环境的情况，流体力学的计算可以为建筑设计提供大量的附加信息。此外，本发明符合日益重视绿色、环保和人群健康，强调以人为本的要求，在建筑设计中可以得到大力的发展。

Description

一种利用流体力学模拟室内热环境的方法

技术领域

本发明属于建筑室内热环境研究技术领域，尤其涉及一种利用流体力学模拟室内热环境的方法。

背景技术

中庭类建筑室内热环境问题一直是国内和国际设计师所关注研究的重点，受到“烟囱效应”和“温室效应”的双重作用，有大面积玻璃幕墙或透明屋顶的中庭室内热环境主要存在两方面的问题：一是夏季中庭上部温度很容易过高，从而影响与中庭相邻房间的热舒适性；二是当中庭发生火灾时，在热压作用下，火灾会迅速猛烈的向上蔓延。中庭的特殊性使得传统的设计方法很难对此类问题提出有根据的解决方案，因此越来越多的中庭实际问题需要借助数值模拟的方法来解决。

目前对于建筑群风环境的预测以及室内热、湿、风环境、设备性能的预测都是通过实验来完成的，如风洞试验、热工实验等，通过实物模型化导致不真实、实际对象的复杂性导致实验困难、无法对各个参数进行研究，并且费用高、周期长。

发明内容

本发明实施例的目的在于提供一种利用流体力学模拟室内热环境的方法，旨在解决现有对建筑室内热环境的测试通过实物模型化存在的模拟不真实、实际对象的复杂性导致实验困难、无法对各个参数进行研究，实验费用高和周期长的问题。

本发明实施例是这样实现的，一种利用流体力学模拟室内热环境的方法，所述利用流体力学模拟室内热环境的方法包括以下步骤：

根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导空调工程施工。

进一步、所述利用流体力学模拟室内热环境的方法具体步骤为：

步骤一、根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

数值模拟:购物中心中庭周边的大多区域送风方式采用上送上回，部分区域吊顶内布置侧送风口和球形喷口向中庭送风，空调方案为全空气系统；

步骤二、经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

步骤三、分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导工程施工，针对施工图设计当中造成的室内温度分布不均匀，风量偏大的问题，提出重新分配各层风口风量的修改方案；

边界条件：中庭夏季热源，根据负荷计算表，把人员负荷赋值在各层地面，灯光负荷赋值在各层吊顶上，计算供冷量；

模拟结果：模型网格数为100×100×100，计算迭代次数为3000次，判定是否收敛。

进一步、所述步骤一数据模拟还包括：

数学模型：

控制方程，为便于编程和计算机识别与运算，PHOENICS将流体流动的控制方程写为通用形式；

湍流模型：高大空间室内气流的流动状态为湍流，湍流模型方程选用k-εCHEN湍流模型；

数学模型的离散方法：使用的PHOENICS软件采用的离散方法为有限体积法；

离散方程组的求解方法：采用的计算方法为改进的压力校正法；

物理模型：购物中心地下三层，地上六层，建筑高度37.25米，建筑面积约28.8万平米，地下二三层为汽车库、设备用房及人防区域，地下一层及地上均为大型商业与餐饮、娱乐，购物中心中庭包括地下一层到地上六层公共区域。

进一步、所述流体流动的控制方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ\Φ}\right)+\mathrm{div}\left(\mathrm{\ρU\Φ}\right)=\mathrm{div}\left(\mathrm{\Γ}\mathrm{grad}\mathrm{\Φ}\right)+S---(2-1)$

其中：ρ为密度，Γ为扩散系数，S为源项，U为速度矢量，对应于Φ的特定意义，Γ和S应具有特定形式；

通用微分方程包含变化率项、对流项、扩散项和源项，当Φ代表不同的物理量，并对扩散系数Γ和源项S作对应的调整变化可导出下述各项方程：

连续性方程：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\mathrm{div}\left(\mathrm{\ρU}\right)=0$

式2-1中Φ＝1，Γ＝0，S＝0

即连续性方程为：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\frac{{\∂\mathrm{\ρu}}_x}{\∂x}+\frac{{\∂\mathrm{\ρu}}_y}{\∂y}+\frac{{\∂\mathrm{\ρu}}_y}{\∂z}=0$

动量方程：

$\frac{{\∂U}_i}{\∂t}+U_i\frac{{\∂U}_i}{{\∂x}_i}=-\frac{1}{\mathrm{\ρ}}\frac{\∂P}{{\∂x}_i}+v\frac{{\∂}^2{U}_i}{{\∂x}_j{\∂x}_j}+B_i$

式2-1中Φ＝Ui，Γ＝ρν，

其中：ν为运动粘性系数，P为压力，Bi为体积力项。

即动量方程为：

$\left\{\begin{array}{c}\mathrm{\ρ}\frac{\mathrm{Du}}{\mathrm{d\τ}}=\mathrm{Fx}-\frac{\∂p}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dv}}{\mathrm{d\τ}}=\mathrm{Fy}-\frac{\∂p}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dw}}{\mathrm{d\τ}}=\mathrm{Fz}-\frac{\∂p}{\∂z}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zz}}}{\∂z}\end{array}\right.$

浓度输运方程：

$\frac{\∂C}{\∂t}+U_j\frac{\∂C}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_c}{\mathrm{\ρ}}\frac{{\∂}^{2}C}{{\∂x}_j{\∂x}_j}+S_c$

式2-1中Φ＝C，Γ＝Γc，S＝Sc

其中：C为质量浓度，Γc为C的扩散系数，Sc为源项，

能量方程：

$\frac{\∂T}{\∂t}+U_j\frac{\∂T}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_T}{\mathrm{\ρ}}\frac{{\∂}^{2}T}{{\∂x}_j{\∂x}_j}+S_T$

式2-1中Φ＝T，Γ＝ΓT，S＝ST

其中：T为温度，ΓT为T的扩散系数，ST为源项

即能量方程为：

${\mathrm{\ρc}}_p\frac{\mathrm{DT}}{\mathrm{d\τ}}=\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)+\frac{\∂}{\∂y}\left(\mathrm{\λ}\frac{\∂T}{\∂y}\right)+\frac{\∂}{\∂z}\left(\mathrm{\λ}\frac{\∂T}{\∂z}\right)+{\mathrm{\α}}_{V}T\frac{\mathrm{Dp}}{\mathrm{d\τ}}+\mathrm{\μ\Φ}.$

进一步、所述k-εCHEN湍流模型为：

标准k-ε模型的单时间尺度的缺点而进行修正，通过在ε方程中增加额外的时间尺度项来提高ε方程的动态响应，如下：

$s_{\mathrm{\ϵ}}=\frac{{-C}_{3e}{P}_k^2}{k}$

其中C3e＝0.25。

进一步、所述有限体积法为：

将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制体积，将待解的微分方程对每一个控制体积积分，便得出一组离散方程，其中的未知数是网格点上的因变量Φ值的数值，为了求出控制体积的积分，必须假定Φ值在网格点之间的变化规律，即设定Φ值分段的分布剖面。

进一步、所述压力校正法采用的离散格式，缺省格式是上风与中心相结合的一阶精度杂交型格式，但备有多种一阶、二阶精度格式供选择，一阶上风格式守恒型低耗散格式,原理是将每个四方形单元用对角线拆分为4个三角形单元，称为X单元，再按上风格式计算对流通量。

进一步、所述步骤三的边界条件的供冷量计算：

根据送风量及送风温差，可计算供冷量为：

Q＝CρVΔt/3600

式中，C为空气比热[1400J/(kg·℃)]，Ρ为空气密度(1.2kg/m³)，V为送风量(m³/h)，Δt为送风温差(℃)。

进一步、所述步骤三的模拟结果通过监测控制方程的迭代残余量以及各表面压力变化来判定是否收敛。

本发明的利用流体力学模拟室内热环境的方法，通过流体力学的计算机仿真，具有模拟真实条件的能力，完善模型以及边界条件保证模拟的准确性，能够较全面的了解不同边界条件下室内环境的情况，流体力学的计算可以为建筑设计提供大量的附加信息。此外，本发明符合日益重视绿色、环保和人群健康，强调以人为本的要求，在建筑设计中可以得到大力的发展。

附图说明

图1是本发明实施例提供的利用流体力学模拟室内热环境的方法的流程图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

图1示出了本发明提供的利用流体力学模拟室内热环境的方法。为了便于说明，仅仅示出了与本发明相关的部分。

本发明的利用流体力学模拟室内热环境的方法，该利用流体力学模拟室内热环境的方法包括以下步骤：

根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导空调工程施工。

作为本发明实施例的一优化方案，利用流体力学模拟室内热环境的方法具体步骤为：

步骤一、根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

数值模拟:购物中心中庭周边的大多区域送风方式采用上送上回，部分区域吊顶内布置侧送风口和球形喷口向中庭送风，空调方案为全空气系统；

步骤二、经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

步骤三、分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导工程施工，针对施工图设计当中造成的室内温度分布不均匀，风量偏大的问题，提出重新分配各层风口风量的修改方案；

边界条件：中庭夏季热源，根据负荷计算表，把人员负荷赋值在各层地面，灯光负荷赋值在各层吊顶上，计算供冷量；

模拟结果：模型网格数为100×100×100，计算迭代次数为3000次，判定是否收敛。

作为本发明实施例的一优化方案，步骤一数据模拟还包括：

数学模型：

控制方程，为便于编程和计算机识别与运算，PHOENICS将流体流动的控制方程写为通用形式；

湍流模型：高大空间室内气流的流动状态为湍流，湍流模型方程选用k-εCHEN湍流模型；

数学模型的离散方法：使用的PHOENICS软件采用的离散方法为有限体积法；

离散方程组的求解方法：采用的计算方法为改进的压力校正法；

物理模型：购物中心地下三层，地上六层，建筑高度37.25米，建筑面积约28.8万平米，地下二三层为汽车库、设备用房及人防区域，地下一层及地上均为大型商业与餐饮、娱乐，购物中心中庭包括地下一层到地上六层公共区域。

作为本发明实施例的一优化方案，流体流动的控制方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ\Φ}\right)+\mathrm{div}\left(\mathrm{\ρU\Φ}\right)=\mathrm{div}\left(\mathrm{\Γ}\mathrm{grad}\mathrm{\Φ}\right)+S---(2-1)$

其中：ρ为密度，Γ为扩散系数，S为源项，U为速度矢量，对应于Φ的特定意义，Γ和S应具有特定形式；

通用微分方程包含变化率项、对流项、扩散项和源项，当Φ代表不同的物理量，并对扩散系数Γ和源项S作对应的调整变化可导出下述各项方程：

连续性方程：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\mathrm{div}\left(\mathrm{\ρU}\right)=0$

式2-1中Φ＝1，Γ＝0，S＝0

即连续性方程为：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\frac{{\∂\mathrm{\ρu}}_x}{\∂x}+\frac{{\∂\mathrm{\ρu}}_y}{\∂y}+\frac{{\∂\mathrm{\ρu}}_y}{\∂z}=0$

动量方程：

$\frac{{\∂U}_i}{\∂t}+U_i\frac{{\∂U}_i}{{\∂x}_i}=-\frac{1}{\mathrm{\ρ}}\frac{\∂P}{{\∂x}_i}+v\frac{{\∂}^2{U}_i}{{\∂x}_j{\∂x}_j}+B_i$

式2-1中Φ＝Ui，Γ＝ρν，

其中：ν为运动粘性系数，P为压力，Bi为体积力项。

即动量方程为：

$\left\{\begin{array}{c}\mathrm{\ρ}\frac{\mathrm{Du}}{\mathrm{d\τ}}=\mathrm{Fx}-\frac{\∂p}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dv}}{\mathrm{d\τ}}=\mathrm{Fy}-\frac{\∂p}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dw}}{\mathrm{d\τ}}=\mathrm{Fz}-\frac{\∂p}{\∂z}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zz}}}{\∂z}\end{array}\right.$

浓度输运方程：

$\frac{\∂C}{\∂t}+U_j\frac{\∂C}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_c}{\mathrm{\ρ}}\frac{{\∂}^{2}C}{{\∂x}_j{\∂x}_j}+S_c$

式2-1中Φ＝C，Γ＝Γc，S＝Sc

其中：C为质量浓度，Γc为C的扩散系数，Sc为源项，

能量方程：

$\frac{\∂T}{\∂t}+U_j\frac{\∂T}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_T}{\mathrm{\ρ}}\frac{{\∂}^{2}T}{{\∂x}_j{\∂x}_j}+S_T$

式2-1中Φ＝T，Γ＝ΓT，S＝ST

其中：T为温度，ΓT为T的扩散系数，ST为源项

即能量方程为：

${\mathrm{\ρc}}_p\frac{\mathrm{DT}}{\mathrm{d\τ}}=\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)+\frac{\∂}{\∂y}\left(\mathrm{\λ}\frac{\∂T}{\∂y}\right)+\frac{\∂}{\∂z}\left(\mathrm{\λ}\frac{\∂T}{\∂z}\right)+{\mathrm{\α}}_{V}T\frac{\mathrm{Dp}}{\mathrm{d\τ}}+\mathrm{\μ\Φ}.$

作为本发明实施例的一优化方案，k-εCHEN湍流模型为：

标准k-ε模型的单时间尺度的缺点而进行修正，通过在ε方程中增加额外的时间尺度项来提高ε方程的动态响应，如下：

$s_{\mathrm{\ϵ}}=\frac{{-C}_{3e}{P}_k^2}{k}$

其中C3e＝0.25。

作为本发明实施例的一优化方案，有限体积法为：

将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制体积，将待解的微分方程对每一个控制体积积分，便得出一组离散方程，其中的未知数是网格点上的因变量Φ值的数值，为了求出控制体积的积分，必须假定Φ值在网格点之间的变化规律，即设定Φ值分段的分布剖面。

作为本发明实施例的一优化方案，压力校正法采用的离散格式，缺省格式是上风与中心相结合的一阶精度杂交型格式，但备有多种一阶、二阶精度格式供选择，一阶上风格式守恒型低耗散格式,原理是将每个四方形单元用对角线拆分为4个三角形单元，称为X单元，再按上风格式计算对流通量。

作为本发明实施例的一优化方案，步骤三的边界条件的供冷量计算：

根据送风量及送风温差，可计算供冷量为：

Q＝CρVΔt/3600

式中，C为空气比热[1400J/(kg·℃)]，Ρ为空气密度(1.2kg/m³)，V为送风量(m³/h)，Δt为送风温差(℃)。

作为本发明实施例的一优化方案，步骤三的模拟结果通过监测控制方程的迭代残余量以及各表面压力变化来判定是否收敛。

下面结合附图及具体实施例对本发明的应用原理作进一步描述。

如图1所示，本发明实施例的利用流体力学模拟室内热环境的方法包括以下步骤：

S101：根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

S102：经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

S103：分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导空调工程施工。

本发明的具体步骤如下：

步骤一、根据设计方案，对全空气系统即各层顶部布置风口（双层百叶风口或球形喷口）、气流场（送风、排风气流组织情况）、空气龄的气流组织方式进行室内CFD模拟；

数值模拟:购物中心中庭周边的大多区域送风方式采用上送上回，部分区域吊顶内布置侧送风口和球形喷口向中庭送风，空调方案为全空气系统，

1、数学模型：

①控制方程：为便于编程和计算机识别与运算，PHOENICS将流体流动的控制方程写为通用形式；

②湍流模型：高大空间室内气流的流动状态为湍流，本报告的湍流模型方程选用k-εCHEN湍流模型；

③数学模型的离散方法：本发明使用的PHOENICS软件采用的离散方法为有限体积法（FVM）；

④离散方程组的求解方法：本发明采用的计算方法为改进的压力校正法（SIMPLEST）；

2、物理模型：购物中心地下三层，地上六层，建筑高度37.25米，建筑面积约28.8万平米，地下二三层为汽车库、设备用房及人防区域，地下一层及地上均为大型商业与餐饮、娱乐，购物中心中庭包括地下一层到地上六层公共区域；

由于地下一层到地上六层各层平面不同，且公共区域的面积和形状都不相同，所以需要根据建筑图分别建筑各层的几何模型，鉴于模型较复杂，在不影响计算结果的情况下对整个几何模型进行简化，建立各层连通区域，不建立各层的上下自动扶梯模型。

步骤二、经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值：

①购物中心中庭地下一层

a.温度场

温度场分布较为均匀，右边中段靠近自动扶梯B107和B108区域和左下中庭区域温度得到改善，在25.6℃左右，其他区域温度分布均匀，且在25℃左右，此温度场分布满足设计要求，

b.速度场，室内空气速度场分布较均匀，平均风速均不超过0.3m/s，满足人的舒适性要求，属于正常设计范围内，

c.PMV，室内热舒适预测指标PMV分布与温度场分布相似，右边中段靠近自动扶梯B107和B108区域和左下中庭区域温度得到改善，大部分区域在-1到1之间，满足人的舒适性要求，属于设计范围之内，

d.PPD室内大部分区域人员不满意率PPD在20%以下，满足人的舒适性要求，属于设计范围之内，

e.空气龄AGE负一层空气龄分布良好，绝大部分区域均不超过15分钟；

步骤三、分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导工程施工，

针对施工图设计当中造成的室内温度分布不均匀，风量偏大等问题，提出重新分配各层风口风量的修改方案，

1）边界条件

地下一层：增加两个黑框区域内送风口风量，使得K-B1-4框中的五组风口（共20个）单个送风量由1000m³/s增加到1050m³/s，K-B1-4框中的五组风口（共20个）单个送风量由625m3/s增加到700m³/s；

表1 地下一层风口数量及送风量

一层：增加黑框区域内送风口风量，使得K-R-01.02框中的三组风口（共12个）单个送风量由1000m³/s增加到1250m³/s，

表2 一层风口数量及送风量

中庭夏季热源，根据负荷计算表，我们把人员负荷赋值在各层地面，灯光负荷赋值在各层吊顶上，

根据送风量及送风温差，可计算供冷量为：

Q＝CρVΔt/3600

式中，C为空气比热[1400J/(kg·℃)]；Ρ为空气密度(1.2kg/m³)；V为送风量(m³/h)，Δt为送风温差(℃)；

2）模拟结果

模型网格数为100×100×100，计算迭代次数为3000次，通过监测控制方程的迭代残余量以及各表面压力变化来判定是否收敛（当相对迭代残余量均小于5×10-3，且各表面风压系数基本不发生变化时，认为所得流场进入了稳态），

为了得到热舒适性指标PMV和PPD，辐射温度设为25℃，服装系数为0.5clo，新陈代谢率为轻微活动值，外部做功为0；

本发明的工作原理：

PHOENICS的模型基础：

1）流体流动的控制方程：

流体流动的控制方程的通用形式为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ\Φ}\right)+\mathrm{div}\left(\mathrm{\ρU\Φ}\right)=\mathrm{div}\left(\mathrm{\Γ}\mathrm{grad}\mathrm{\Φ}\right)+S---(2-1)$

其中：ρ为密度；Γ为扩散系数；S为源项；U为速度矢量，对应于Φ的特定意义，Γ和S应具有特定形式，

通用微分方程包含变化率项、对流项、扩散项和源项，当Φ代表不同的物理量，并对扩散系数Γ和源项S作对应的调整变化可导出下述各项方程：

连续性方程：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\mathrm{div}\left(\mathrm{\ρU}\right)=0---(2-3)$

式2-1中Φ＝1，Γ＝0，S＝0

即连续性方程为：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\frac{{\∂\mathrm{\ρu}}_x}{\∂x}+\frac{{\∂\mathrm{\ρu}}_y}{\∂y}+\frac{{\∂\mathrm{\ρu}}_y}{\∂z}=0---(2-4)$

②动量方程：

$\frac{{\∂U}_i}{\∂t}+U_i\frac{{\∂U}_i}{{\∂x}_i}=-\frac{1}{\mathrm{\ρ}}\frac{\∂P}{{\∂x}_i}+v\frac{{\∂}^2{U}_i}{{\∂x}_j{\∂x}_j}+B_i---(2-5)$

式2-1中Φ＝Ui，Γ＝ρν，

其中：ν为运动粘性系数，P为压力，Bi为体积力项，

即动量方程为：

$\left\{\begin{array}{c}\mathrm{\ρ}\frac{\mathrm{Du}}{\mathrm{d\τ}}=\mathrm{Fx}-\frac{\∂p}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dv}}{\mathrm{d\τ}}=\mathrm{Fy}-\frac{\∂p}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dw}}{\mathrm{d\τ}}=\mathrm{Fz}-\frac{\∂p}{\∂z}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zz}}}{\∂z}\end{array}\right.---(2-6)$

③浓度输运方程：

$\frac{\∂C}{\∂t}+U_j\frac{\∂C}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_c}{\mathrm{\ρ}}\frac{{\∂}^{2}C}{{\∂x}_j{\∂x}_j}+S_c---(2-7)$

式2-1中Φ＝C，Γ＝Γc，S＝Sc

其中：C为质量浓度，Γc为C的扩散系数，Sc为源项，

④能量方程：

$\frac{\∂T}{\∂t}+U_j\frac{\∂T}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_T}{\mathrm{\ρ}}\frac{{\∂}^{2}T}{{\∂x}_j{\∂x}_j}+S_T---(2-8)$

式2-1中Φ＝T，Γ＝ΓT，S＝ST

其中：T为温度，ΓT为T的扩散系数，ST为源项

即能量方程为：

${\mathrm{\ρc}}_p\frac{\mathrm{DT}}{\mathrm{d\τ}}=\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)+\frac{\∂}{\∂y}\left(\mathrm{\λ}\frac{\∂T}{\∂y}\right)+\frac{\∂}{\∂z}\left(\mathrm{\λ}\frac{\∂T}{\∂z}\right)+{\mathrm{\α}}_{V}T\frac{\mathrm{Dp}}{\mathrm{d\τ}}+\mathrm{\μ\Φ}---(2-9)$

2）湍流模型

PHOENICS内置了21种适合于各种雷诺数场合的湍流模型，包括雷诺应力模型、多流体湍流模型和通量模型及k-ε模型的各种变异，其中LVEL模型和MF（多项流）模型是PHOENICS独有的，下面介绍几种PHOENICS中常用于模拟室外热环境的湍流模型：

①k-ε双方程模型：

k-ε模型是目前工程问题中应用最为广泛的湍流模型，该模型由Launder和Spalding提出，经许多学者和工程技术人员的研究与应用证明k-ε：模型优于除雷诺应力模型以外的其他各种模型，

②RNGk-ε湍流模型：

Yakhot与Orszag在重正化群方法的基础上对标准k-ε模型进行了改进，该模型用RNG方法发展一个适用于大尺度的理论，而小尺度的影响效果则通过改进传输系数来重新描述，RNG k-ε对于标准k-ε中的几个常数的值进行了轻微的改动，并在ε传输方程中加了一项：

$s_{\mathrm{\ϵ}}=\frac{{-\mathrm{a\ϵ}}^2}{k}---(2-15)$

其中：

$a={\mathrm{C\μ\η}}^3\frac{1-\frac{\mathrm{\η}}{{\mathrm{\η}}_0}}{1+{\mathrm{\β\η}}^3}---(2-16)$

③Chen-Kim k-ε湍流模型：

该模型是Chen与Kim在1987年对k-ε模型修正后的模型，其中针对标准k-ε模型的单时间尺度的缺点而进行修正，通过在ε方程中增加额外的时间尺度项来提高ε方程的动态响应，该项如下：

$s_{\mathrm{\ϵ}}=\frac{{-C}_{3e}{P}_k^2}{k}---(2-18)$

其中C3e＝0.25，

Chen-Kim k-ε模型既保持了RNG k-ε模型的优点，而且在模拟射流与羽流时与标准k-ε模型效果一样好，大空间建筑的模拟多采用此模型，

④LVEL湍流模型：

LVEL模型属于零方程湍流模型，是PHOENICS所特有的模型，特别适合暖通空调领域的模拟，标准k-ε模型不太适用于空间中存在大量固体物体的情况，而LVEL模型则特别适合，其采用Spalding墙体法则来计算νt，对于三维问题也有很好的经济性，一般适用于雷诺数较低的情况，

3）PHOENICS中的壁面函数

虽然标准k-ε两方程模型是应用最为广泛的湍流模型，它忽略了分子粘性，对于离开壁面一定距离以上的湍流区域具有相当高的适用性，但在与壁面相邻的粘性底层，由于Re数很低，必须充分考虑分子粘性的影响，所以必须对标准k-ε两方程模型进行改进，对于与壁面接触的这部分区域的处理一般有两种方法，一种是采用低Re数模型，另一种则是用壁面函数法，由于粘性底层中的速度梯度极为陡峻，为了得到理想的解答，必须在粘性底层中布置极为细密的网格，需要花费较多的计算机资源，因此一般更为常用的是采用壁面函数法即在粘性底层的外边缘规定相应的边界条件，PHOENICS中提供了以下4种壁面函数：①布拉修斯公式（Blasius）；②对数律（Log-law）；③通用对数律（General-log-law）；④完全粗糙对数律（Fully-rough），

（5）PHOENICS的离散方法——有限体积法

PHOENICS的离散方法采用有限体积法（FVM，即Finite VolumeMethod），也称为控制体积法，其基本思路是：将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制体积，将待解的微分方程对每一个控制体积积分，便得出一组离散方程，其中的未知数是网格点上的因变量Φ值的数值，为了求出控制体积的积分，必须假定Φ值在网格点之间的变化规律，即设定Φ值分段的分布剖面，从积分区域的选取方法看来，有限体积法属于加权余量法中的子区域法；从未知解的近似方法看，有限体积法属于采用局部近似的离散方法，简言之，子区域法加离散，就是有限体积法的基本方法，

（6）PHOENICS的求解方法——SIMPLEST算法

SIMPLEST算法是由SIMPLE（即Semi-Implicit Method forPressure-Linked Equations，压力校正法）算法改进得来的，SIMPLE算法自1972年问世以来，在计算流体力学及计算传热学中得到广泛的应用，同时也得到不断的改进与发展，成为目前CFD通用软件求解低速问题的主要方法，该方法所采用的离散格式，缺省格式是上风与中心相结合的一阶精度杂交型格式（Hybrid Scheme），但备有多种一阶、二阶精度格式供选择，此外，PHOENICS还提供了一种改进的一阶上风格式CLDS(守恒型低耗散格式),其原理是将每个四方形单元用对角线（面）拆分为4个三角形（三棱柱）单元，称为X单元，再按上风格式计算对流通量，可有效降低流动方向与网格偏斜导致的数值耗散。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种利用流体力学模拟室内热环境的方法，其特征在于，所述利用流体力学模拟室内热环境的方法包括以下步骤：

根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导空调工程施工。

2.如权利要求1所述的利用流体力学模拟室内热环境的方法，其特征在于，所述利用流体力学模拟室内热环境的方法具体步骤为：

步骤一、根据设计方案，对全空气系统即各层顶部布置风口、气流场、空气龄的气流组织方式进行室内CFD模拟；

数值模拟:购物中心中庭周边的大多区域送风方式采用上送上回，部分区域吊顶内布置侧送风口和球形喷口向中庭送风，空调方案为全空气系统；

步骤二、经过模拟计算后，得出室内温度场、湿度场、PMV值及PPD值；

步骤三、分别进行分析、比较和评价，选出最优方式，对设计不当的地方提出修改意见，指导工程施工，针对施工图设计当中造成的室内温度分布不均匀，风量偏大的问题，提出重新分配各层风口风量的修改方案；

边界条件：中庭夏季热源，根据负荷计算表，把人员负荷赋值在各层地面，灯光负荷赋值在各层吊顶上，计算供冷量；

模拟结果：模型网格数为100×100×100，计算迭代次数为3000次，判定是否收敛。

3.如权利要求2所述的利用流体力学模拟室内热环境的方法，其特征在于，所述步骤一数据模拟还包括：

数学模型：

控制方程，为便于编程和计算机识别与运算，PHOENICS将流体流动的控制方程写为通用形式；

湍流模型：高大空间室内气流的流动状态为湍流，湍流模型方程选用k-εCHEN湍流模型；

数学模型的离散方法：使用的PHOENICS软件采用的离散方法为有限体积法；

离散方程组的求解方法：采用的计算方法为改进的压力校正法；

物理模型：购物中心地下三层，地上六层，建筑高度37.25米，建筑面积约28.8万平米，地下二三层为汽车库、设备用房及人防区域，地下一层及地上均为大型商业与餐饮、娱乐，购物中心中庭包括地下一层到地上六层公共区域。

4.如权利要求3所述的利用流体力学模拟室内热环境的方法，其特征在于，所述流体流动的控制方程为：

$\frac{\∂}{\∂t}\left(\mathrm{\ρ\Φ}\right)+\mathrm{div}\left(\mathrm{\ρU\Φ}\right)=\mathrm{div}\left(\mathrm{\Γ}\mathrm{grad}\mathrm{\Φ}\right)+S---(2-1)$

其中：ρ为密度，Γ为扩散系数，S为源项，U为速度矢量，对应于Φ的特定意义，Γ和S应具有特定形式；

通用微分方程包含变化率项、对流项、扩散项和源项，当Φ代表不同的物理量，并对扩散系数Γ和源项S作对应的调整变化可导出下述各项方程：

连续性方程：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\mathrm{div}\left(\mathrm{\ρU}\right)=0$

式2-1中Φ＝1，Γ＝0，S＝0

即连续性方程为：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\frac{{\∂\mathrm{\ρu}}_x}{\∂x}+\frac{{\∂\mathrm{\ρu}}_y}{\∂y}+\frac{{\∂\mathrm{\ρu}}_y}{\∂z}=0$

动量方程：

$\frac{{\∂U}_i}{\∂t}+U_i\frac{{\∂U}_i}{{\∂x}_i}=-\frac{1}{\mathrm{\ρ}}\frac{\∂P}{{\∂x}_i}+v\frac{{\∂}^2{U}_i}{{\∂x}_j{\∂x}_j}+B_i$

式2-1中Φ＝Ui，Γ＝ρν，

其中：ν为运动粘性系数，P为压力，Bi为体积力项，

即动量方程为：

$\left\{\begin{array}{c}\mathrm{\ρ}\frac{\mathrm{Du}}{\mathrm{d\τ}}=\mathrm{Fx}-\frac{\∂p}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{xz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dv}}{\mathrm{d\τ}}=\mathrm{Fy}-\frac{\∂p}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{yz}}}{\∂z}\\ \mathrm{\ρ}\frac{\mathrm{Dw}}{\mathrm{d\τ}}=\mathrm{Fz}-\frac{\∂p}{\∂z}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zx}}}{\∂x}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zy}}}{\∂y}+\frac{{\∂\mathrm{\τ}}_{\mathrm{zz}}}{\∂z}\end{array}\right.$

浓度输运方程：

$\frac{\∂C}{\∂t}+U_j\frac{\∂C}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_c}{\mathrm{\ρ}}\frac{{\∂}^{2}C}{{\∂x}_j{\∂x}_j}+S_c$

式2-1中Φ＝C，Γ＝Γc，S＝Sc

其中：C为质量浓度，Γc为C的扩散系数，Sc为源项，

能量方程：

$\frac{\∂T}{\∂t}+U_j\frac{\∂T}{{\∂x}_j}=\frac{{\mathrm{\Γ}}_T}{\mathrm{\ρ}}\frac{{\∂}^{2}T}{{\∂x}_j{\∂x}_j}+S_T$

式2-1中Φ＝T，Γ＝ΓT，S＝ST

其中：T为温度，ΓT为T的扩散系数，ST为源项

即能量方程为：

${\mathrm{\ρc}}_p\frac{\mathrm{DT}}{\mathrm{d\τ}}=\frac{\∂}{\∂x}\left(\mathrm{\λ}\frac{\∂T}{\∂x}\right)+\frac{\∂}{\∂y}\left(\mathrm{\λ}\frac{\∂T}{\∂y}\right)+\frac{\∂}{\∂z}\left(\mathrm{\λ}\frac{\∂T}{\∂z}\right)+{\mathrm{\α}}_{V}T\frac{\mathrm{Dp}}{\mathrm{d\τ}}+\mathrm{\μ\Φ}.$

5.如权利要求3所述的利用流体力学模拟室内热环境的方法，其特征在于，所述k-εCHEN湍流模型为：

标准k-ε模型的单时间尺度的缺点而进行修正，通过在ε方程中增加额外的时间尺度项来提高ε方程的动态响应，如下：

$s_{\mathrm{\ϵ}}=\frac{{-C}_{3e}{P}_k^2}{k}$

其中C3e＝0.25。

6.如权利要求3所述的利用流体力学模拟室内热环境的方法，其特征在于，所述有限体积法为：

将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制体积，将待解的微分方程对每一个控制体积积分，便得出一组离散方程，其中的未知数是网格点上的因变量Φ值的数值，为了求出控制体积的积分，必须假定Φ值在网格点之间的变化规律，即设定Φ值分段的分布剖面。

7.如权利要求3所述的利用流体力学模拟室内热环境的方法，其特征在于，所述压力校正法采用的离散格式，缺省格式是上风与中心相结合的一阶精度杂交型格式，但备有多种一阶、二阶精度格式供选择，一阶上风格式守恒型低耗散格式,原理是将每个四方形单元用对角线拆分为4个三角形单元，称为X单元，再按上风格式计算对流通量。

8.如权利要求2所述的利用流体力学模拟室内热环境的方法，其特征在于，所述步骤三的边界条件的供冷量计算：

根据送风量及送风温差，可计算供冷量为：

Q＝CρVΔt/3600

式中，C为空气比热[1400J/(kg·℃)]，Ρ为空气密度(1.2kg/m³)，V为送风量(m³/h)，Δt为送风温差(℃)。

9.如权利要求2所述的利用流体力学模拟室内热环境的方法，其特征在于，所述步骤三的模拟结果通过监测控制方程的迭代残余量以及各表面压力变化来判定是否收敛。

Images