王潤(rùn)濤，段明輝，楊 方*，姚 陽(yáng)，張 薇，李友澤

（1.東北農(nóng)業(yè)大學(xué)電氣與信息學(xué)院，哈爾濱 150030；2.天津大學(xué)電氣與自動(dòng)化工程學(xué)院，天津 300072）

隨著計(jì)算機(jī)技術(shù)發(fā)展，數(shù)值模擬水平不斷提高，CFD（Computational fluid dynamics）是可數(shù)值計(jì)算和圖像顯示工具，廣泛應(yīng)用于溫室園藝領(lǐng)域。目前，高效節(jié)能型日光溫室在我國(guó)北方較為普及[1]，但在我國(guó)寒冷地區(qū)（如哈爾濱）冬季，需依靠輔助加熱設(shè)備才能保證室內(nèi)溫度。為保暖冬季溫室中不進(jìn)行通風(fēng)，而溫室中加熱設(shè)備熱量無(wú)空氣流通，單靠熱輻射作用不能有效擴(kuò)散，出現(xiàn)溫度場(chǎng)分布不均勻、多區(qū)域不能達(dá)到植物生長(zhǎng)適宜溫度等問(wèn)題，溫室作物生長(zhǎng)環(huán)境不理想。

利用CFD模擬溫室氣流組織分布，根據(jù)模擬結(jié)果研究分析各影響因素情況，優(yōu)化溫室內(nèi)環(huán)境。本文針對(duì)溫室內(nèi)溫度場(chǎng)問(wèn)題，利用CFD技術(shù)建立溫室數(shù)值模擬模型，與試驗(yàn)測(cè)量值比較驗(yàn)證模型準(zhǔn)確性，利用CFD軟件模擬分析溫室內(nèi)溫度場(chǎng)均勻性優(yōu)化情況，驗(yàn)證優(yōu)化效果良好，可為冬季溫室內(nèi)環(huán)境的優(yōu)化提供參考。

1 模型及計(jì)算條件

1.1 溫室?guī)缀文Ｐ?/h3>

實(shí)驗(yàn)溫室采用東北農(nóng)業(yè)大學(xué)園藝學(xué)院溫室，溫室全長(zhǎng)58 m，中部被pc板隔開(kāi)，形成兩個(gè)溫室環(huán)境，試驗(yàn)選取其中一個(gè)。實(shí)驗(yàn)溫室長(zhǎng)27.5 m，跨度7 m，脊高3.5 m，后墻高2.5 m，后坡水平投影1.5 m，日光溫室?guī)缀文Ｐ腿鐖D1所示。

圖1 溫室?guī)缀文Ｐ虵ig.1 Geometric model of the greenhouse

1.2 溫室CFD模型

1.2.1 湍流模型

湍流流動(dòng)是一種邊界條件保持不變，速度等流動(dòng)特性均隨機(jī)變化的不穩(wěn)定流動(dòng)狀態(tài)，通常用Reynolds數(shù)（用Re表示）作為判斷依據(jù)，當(dāng)Re（雷諾數(shù)）大于臨界值時(shí)為湍流流動(dòng)。加熱條件下溫室內(nèi)氣流可看成湍流流動(dòng)[2-8]。根據(jù)溫室的情況選取標(biāo)準(zhǔn)模型進(jìn)行數(shù)值模擬計(jì)算。標(biāo)準(zhǔn)模型的輸運(yùn)方程為：

式中，Gk為平均速度梯度引起的湍動(dòng)能k的產(chǎn)生項(xiàng)；Gb為浮力引起的湍動(dòng)能k的產(chǎn)生項(xiàng)；YM代表可壓湍流中脈動(dòng)擴(kuò)張的貢獻(xiàn)；C1ε、C2ε、C3ε為經(jīng)驗(yàn)常數(shù)；σk、σε分別是與湍動(dòng)能k和耗散率ε對(duì)應(yīng)的Prandtl數(shù)；Sk、Sε是用戶定義的源項(xiàng)。

高Re數(shù)的湍流模型針對(duì)充分發(fā)展的湍流才有效，對(duì)近壁區(qū)內(nèi)的流動(dòng)，湍流發(fā)展并不充分，目前一般采用壁面函數(shù)法或低Re數(shù)k-ε模型求解。壁面函數(shù)法為軟件默認(rèn)處理辦法，對(duì)各種壁面流動(dòng)有效，計(jì)算效率較高，工程實(shí)用性強(qiáng)，因此采用壁面函數(shù)法處理近壁面區(qū)域[9-10]。

1.2.2 基本控制方程

由于哈爾濱地區(qū)冬季室外氣溫很低，冬季溫室不會(huì)開(kāi)窗通風(fēng)，室內(nèi)主要是由于溫度變化引起的自然對(duì)流，根據(jù)溫室內(nèi)氣流特點(diǎn)，為便于處理溫差帶來(lái)的浮升力項(xiàng)，采用Boussinesq假設(shè)簡(jiǎn)化加熱空氣產(chǎn)生的自然對(duì)流[11]，即假設(shè)流體密度的變化并不明顯改變流體性質(zhì)，即除密度外流體的其他物性不變；密度變化對(duì)慣性力項(xiàng)、壓力差項(xiàng)和粘性力項(xiàng)影響可忽略不計(jì)；對(duì)密度僅考慮動(dòng)量方程中與質(zhì)量力有關(guān)項(xiàng)，其余各項(xiàng)中密度作為常數(shù)。

由于流動(dòng)處于湍流狀態(tài)，系統(tǒng)還要遵守附加的湍流輸運(yùn)方程[12]。控制方程通用形式：

式中，φ為通用變量，可代表u、v、w、T等求解變量；Γ為廣義擴(kuò)散系數(shù)；S為廣義源項(xiàng)。

1.2.3 輻射模型

CFD為用戶提供5種輻射模型，參考文獻(xiàn)[13-16]及試驗(yàn)溫室實(shí)際情況，選取DO模型計(jì)算輻射影響。輻射方程為：

2 網(wǎng)格的劃分及邊界條件選擇

2.1 計(jì)算域網(wǎng)格劃分

網(wǎng)格是CFD模型的幾何表達(dá)形式，也是模擬與分析的載體。網(wǎng)格質(zhì)量對(duì)CFD的計(jì)算精度和計(jì)算效率具有重要影響[17]。根據(jù)溫室情況，選取溫室內(nèi)空氣作為計(jì)算區(qū)域，采用混合網(wǎng)格形式劃分，并對(duì)熱水管道進(jìn)行網(wǎng)格加密。

2.2 邊界條件的選擇

冬季溫室不與外界進(jìn)行通風(fēng)，只存在內(nèi)部流固耦合傳熱，因此邊界條件全部設(shè)置為Wall（壁面邊界條件）；地面給定實(shí)際測(cè)量的固定溫度；由于散熱器為螺旋翅片管型，直接計(jì)算復(fù)雜，參考文獻(xiàn)[18-19]對(duì)散熱器進(jìn)行化簡(jiǎn)；墻體及覆蓋材料為對(duì)流與外部輻射相結(jié)合邊界類型；一般采用空氣有效溫度作為外部輻射溫度，用以反映空氣輻射背景對(duì)覆蓋材料傳熱影響，利用公式Tsky=0.0552T01.5對(duì)外部輻射溫度進(jìn)行計(jì)算[20]；溫室模型中未考慮植物影響。主要邊界條件設(shè)置如表1所示。

表1 邊界類型及取值Table 1 Main boundary types and their values

3 CFD模擬及分析

首先模擬并分析夜間供暖情況下溫室熱環(huán)境分布。室外溫度為-28.4℃（244.75 K），分別取Z=5 m，Z=12 m，Z=20 m三個(gè)剖面，經(jīng)分析發(fā)現(xiàn)沿Z軸方向上，溫度分布基本相同。給出Z=12 m溫度線分布圖2。由圖2可見(jiàn)，在靠近散熱器附近1 m，溫度線分布密集，溫度較高；但是僅靠熱源的熱輻射作用，大部分熱量并不能傳遞到溫室內(nèi)部，靠近溫室中部區(qū)域的溫度較低；遠(yuǎn)離熱源區(qū)域，受外部低溫影響，越靠近溫室上方，溫度越低。溫度分布梯度很大，作物生長(zhǎng)所處大部分區(qū)域的溫度狀況都不理想，溫度場(chǎng)分布均勻性很差。

圖2 Z=12 m截面的溫度線分布Fig.2 Distribution of temperature about the cross section of Z=12 m

溫室沿Y方向（沿溫室高度方向）剖面溫度分布見(jiàn)圖3。在靠近熱源附近，溫度分布梯度較大，遠(yuǎn)離熱源區(qū)域，溫度下降梯度明顯，在0.25 m處大部分溫度處于283 K，在0.5 m處大部分的溫度在280 K及以下。溫度分布在越靠近溫室上方的區(qū)域溫度越低，僅0.25 m已相差3 K（3℃），顯然溫室內(nèi)部分區(qū)域溫度偏低，溫度場(chǎng)分布不均勻，植物生長(zhǎng)環(huán)境不理想。

圖3 Y=0.25 m，Y=0.5 m溫度布云Fig.3 Distribution of temperature in cloud about the two cross sections that Y=0.25 m and Y=0.5 m

對(duì)溫室白天室外溫度為-14.85℃的情況進(jìn)行模擬計(jì)算。選取與上面相同的截面進(jìn)行分析，計(jì)算結(jié)果如下圖4?？梢?jiàn)，由于白天陽(yáng)光的影響，溫室內(nèi)的氣溫比夜間高，但溫度場(chǎng)分布的均勻性仍不理想。白天溫室內(nèi)Y=0.25 m截面層上氣溫可以到達(dá)285 K，Y=0.5 m截面層大部分區(qū)域氣溫在282 K。依舊是靠近溫室頂棚區(qū)域溫度低，靠近溫室熱源附近的溫度分布密集，溫度較高，溫度偏低且分布梯度較大。綜合模擬結(jié)果以及溫室的實(shí)際情況，可以得出，相對(duì)與冬季的外部溫度，溫室內(nèi)較溫暖。室內(nèi)溫度受陽(yáng)光照射的影響很大，白天溫度較高，夜間由于供暖設(shè)備及白天儲(chǔ)存熱量的原因，溫度也不會(huì)很低。但是溫室內(nèi)溫度分布梯度很大，分層明顯，靠近熱源部分的溫度較高，溫室中部溫度偏低，散熱器熱量不能有效充分?jǐn)U散。溫室上部受外界干擾影響較為明顯，尤其是在較冷的時(shí)間，供熱容量沒(méi)有相應(yīng)提升，會(huì)對(duì)溫室內(nèi)溫度影響較大。溫室內(nèi)溫度整體上可以滿足植物生命活動(dòng)的基本需求，但是溫度偏低且溫度場(chǎng)的均勻性并不理想，不能為植物提供良好的溫度環(huán)境，影響溫室作物的質(zhì)量及產(chǎn)量。

圖4 Y=0.25 m，Y=0.5 m截面溫度分布云Fig.4 Distribution of temperature in cloud about the two cross sections that Y=0.25 m and Y=0.5 m

4 模型準(zhǔn)確性驗(yàn)證

為驗(yàn)證模擬的準(zhǔn)確性，將模擬所得溫室溫度與試驗(yàn)測(cè)量溫度進(jìn)行比較，一共選擇9個(gè)測(cè)量點(diǎn)，測(cè)量點(diǎn)的具體選取情況見(jiàn)下表2。其中，沿溫室的寬度方向定義X，沿溫室的高度方向定義Y，沿溫室的長(zhǎng)度方向定義Z（參照?qǐng)D2所示）。選取溫室內(nèi)一個(gè)立方體區(qū)域及其中心點(diǎn)作為測(cè)點(diǎn)進(jìn)行測(cè)量。

表2 溫室內(nèi)測(cè)點(diǎn)分布Table 2 Distribution of measure points of the greenhouse

試驗(yàn)選取較為寒冷的天氣狀況，選取的夜間室外溫度為-28.4℃（244.75 K），與試驗(yàn)測(cè)量室外為-28.4℃時(shí)相對(duì)應(yīng)的溫室內(nèi)部溫度加以比較，如表3所示。由表3可見(jiàn)，由于溫室上部距離熱源較遠(yuǎn)，受外界的影響較大，溫室中靠近上部頂棚的溫度相對(duì)較低，如測(cè)點(diǎn)B、D、F、H，在靠近溫室下部的溫度相對(duì)較高，如測(cè)點(diǎn)A、C、E、G，雖然測(cè)點(diǎn)在高度上相差僅有0.5 m，但相差最多有5℃，溫度分布梯度很大，溫度場(chǎng)均勻性很不理想。經(jīng)比較可以看出，模擬結(jié)果與試驗(yàn)測(cè)量的數(shù)據(jù)之間最大誤差為0.5℃（B點(diǎn)），最大誤差百分比為10%，證明模型建立基本準(zhǔn)確，模擬結(jié)果證明，雖然模擬結(jié)果與溫室的實(shí)際情況存在一定誤差，但是所建CFD模型及邊界條件的設(shè)定可行，能夠較為準(zhǔn)確地反應(yīng)溫室內(nèi)溫度場(chǎng)的變化情況。

表3 溫室實(shí)測(cè)值與預(yù)測(cè)值對(duì)比Table 3 Measured date compared with experimental date

5 氣流組織優(yōu)化及分析

由于溫室內(nèi)的自然微循環(huán)不能產(chǎn)生足夠均勻且適宜的氣溫，為優(yōu)化溫室內(nèi)的環(huán)境，提出在溫室內(nèi)加入內(nèi)循環(huán)風(fēng)扇。內(nèi)循環(huán)風(fēng)扇并不引入外界空氣，并不會(huì)影響溫室整體溫度，并且合理的風(fēng)扇安裝，可以在風(fēng)扇工作時(shí)形成良好的氣流組織（氣流組織是指對(duì)氣流流向和均勻度按一定要求進(jìn)行組織，使工作區(qū)形成比較均勻而穩(wěn)定的溫度、濕度及氣流速度），不僅使溫室的溫度均勻，空氣流通可保持植物葉面干燥，抑制病蟲(chóng)害的發(fā)生。

根據(jù)實(shí)驗(yàn)溫室的實(shí)際情況，經(jīng)計(jì)算溫室所需通風(fēng)量為7500 m3·h-1，初步選擇在溫室的對(duì)角線安置兩臺(tái)內(nèi)循環(huán)風(fēng)扇進(jìn)行試驗(yàn)。對(duì)夜間供暖情況下的溫室進(jìn)行數(shù)值計(jì)算分析，模型選擇與前文夜間供暖相同的室內(nèi)外條件。

選取與前文相同的界面進(jìn)行分析。如圖5所示，加入內(nèi)循環(huán)風(fēng)扇后，受風(fēng)力的影響，雖然三個(gè)截面的溫度線分布各不相同，但比較與圖2溫度分布，溫室內(nèi)溫度場(chǎng)有了很大的改善。溫室主體部分溫度最大浮動(dòng)在285～282 K之間，比較之前溫室內(nèi)溫度分布狀況，溫度分布均勻性有改善，溫度提升。

圖5 3個(gè)Z截面的溫度線分布Fig.5 Distribution of temperature about the three cross sections of Z

溫室沿Y方向剖面溫度分布見(jiàn)圖6，對(duì)比圖3可見(jiàn)，兩個(gè)界面層的溫度均勻性都有很大提升，兩個(gè)截面的大部分區(qū)域都在284 K以上，溫室內(nèi)溫度場(chǎng)的均勻性得到顯著改善。為驗(yàn)證加入內(nèi)循環(huán)風(fēng)扇后，溫室內(nèi)溫度及溫度場(chǎng)均勻性改善狀況，各測(cè)點(diǎn)溫度為A 11.45℃，B 11.15℃，C 11.95℃，D 11.65℃，E 10.95℃，F(xiàn) 10.75℃，G 11.35℃，H 11.15℃，I 11.15℃?？梢?jiàn)各測(cè)點(diǎn)的溫度均有提升，尤其是靠近溫室上部的各點(diǎn)，溫度提升幅度很大，溫室內(nèi)溫度場(chǎng)的均勻性有改善。

對(duì)模型的模擬結(jié)果做顯著性檢驗(yàn)，如表4所示，置信水平為5%?？梢?jiàn)，放置風(fēng)扇后，由于熱量得到有效擴(kuò)散，溫室內(nèi)的溫度較之前平均提升3.35℃，大約為42%左右，有明顯提升。放置風(fēng)扇前，溫度的平均值與期望值存在顯著性差異；放置內(nèi)循環(huán)風(fēng)扇對(duì)溫室的溫度場(chǎng)做優(yōu)化后， 可以證明優(yōu)化后的溫度平均值較期望值無(wú)顯著差異，結(jié)合優(yōu)化后溫室內(nèi)各測(cè)點(diǎn)溫度可見(jiàn)，試驗(yàn)溫室的溫度有顯著提升，溫度場(chǎng)的均勻性得到改善，達(dá)到預(yù)期結(jié)果。

表4 優(yōu)化前后的顯著性檢驗(yàn)Table 4 Measured date compared with experimental date

6 結(jié)論

由于溫室環(huán)境復(fù)雜，模擬過(guò)程中并不能設(shè)置精確的邊界條件，因此數(shù)值模擬結(jié)果與溫室的實(shí)際狀況存在差距，通過(guò)驗(yàn)證可知，溫室模型基本準(zhǔn)確，可反映溫室內(nèi)溫度變化狀況。

從模擬結(jié)果可見(jiàn)，溫室內(nèi)總體溫度分布不均勻，在兩側(cè)靠近熱源附近的溫度分布梯度較大，溫度較高；在中間區(qū)域溫度相對(duì)較低。雖然熱水供熱系統(tǒng)相對(duì)于其他供熱方式優(yōu)勢(shì)大，但存在溫度分布不均勻，不能產(chǎn)生植物生長(zhǎng)適宜的溫度場(chǎng)。尤其是在哈爾濱地區(qū)，不僅冬季溫度很低，而且溫度浮動(dòng)很大，如果不能及時(shí)調(diào)節(jié)鍋爐的供熱狀況，熱量不能有效擴(kuò)散到離供熱設(shè)備較遠(yuǎn)的區(qū)域，這些區(qū)域會(huì)受外界影響較大，容易對(duì)溫室內(nèi)溫度場(chǎng)產(chǎn)生較大影響，不利于植物生長(zhǎng)，內(nèi)循環(huán)風(fēng)扇引入促進(jìn)熱量有效擴(kuò)散，可優(yōu)化溫室溫度場(chǎng)，適當(dāng)減小外部溫度對(duì)溫室的影響。

本試驗(yàn)表明，在溫室內(nèi)部增加內(nèi)循環(huán)風(fēng)扇，優(yōu)化溫室氣流組織形式，有效促使溫室內(nèi)部熱源熱量擴(kuò)散，雖然與期望溫度有差異，但已得到相對(duì)均勻和適宜溫度場(chǎng)，不僅提高溫室適用范圍，產(chǎn)生空氣流通還可在一定程度上優(yōu)化作物生長(zhǎng)環(huán)境，防止由于空氣不流通造成的葉面結(jié)露等問(wèn)題。為在最佳位置安置合適的可調(diào)控風(fēng)源，建立合理的溫室內(nèi)的氣流組織形式，優(yōu)化溫室內(nèi)溫度場(chǎng)及速度場(chǎng)環(huán)境，建立適宜的溫室內(nèi)微氣候環(huán)境提供科學(xué)依據(jù)。

[1] 齊飛.我國(guó)溫室及配套設(shè)備產(chǎn)業(yè)現(xiàn)狀及發(fā)展趨勢(shì)[J].上海農(nóng)業(yè)學(xué)報(bào),2005,21(1):53-57.

[2] 趙雯倩.計(jì)算流體力學(xué)在溫室通風(fēng)中的應(yīng)用研究[J].農(nóng)業(yè)科技與裝備,2011(8):23-27.

[3] Kacira M,Sase S,Okushima L.Optimization of vent configuration by evaluating greenhouse and plant canopy ventilation rates under wind-induced ventilation[J].Trans ASAE,2004,47(6):2059-2067．

[4] 王建,汪小旵,丁為民.風(fēng)壓通風(fēng)的單棟溫室內(nèi)部流場(chǎng)的ANSYSCFD模擬[J]．農(nóng)業(yè)機(jī)械學(xué)報(bào),2007(3):114-116．

[5] Bennisa N,Duplaisb J,Eneab G,et al．Greenhouse climate modeling and robust control[J]．Computers and Electronic in Agriculture,2008,61(2):96-107．

[6] 何國(guó)敏.現(xiàn)代化溫室溫度場(chǎng)數(shù)字化模擬研究[D].南京:南京農(nóng)業(yè)大學(xué),2011．

[7] Mistriotis A．Arcidiacono C,Picuno P,et al．Computational analysis of ventilation in greenhouses at zero-and low-windspeeds[J]．Agricultural and Forest Meteorology,1997,88:121-135.

[8] 楊振超．日光溫室內(nèi)最佳風(fēng)速指標(biāo)與CFD模擬[D]．楊凌:西北農(nóng)林科技大學(xué),2006.

[9] Yakhot V,Orzag S A.Renormalization group analysis of turbulence:basic theory[J].Scient Comput,1986(1):3-11.

[10] Versteeg H K,Malalasekera W.An Introduction to computational fluid dynamics:The finite volume method[M].New York:Wiley,1995.

[11] Molina-Aiz F,Valera D,Alvarez A,et al.A wind tunnel study of airflow through horticultural crops:Determination of the drag coefficient[J].Biosystem Engineering,2006,93(4):447-457.

[12] Takashi Y,Yusuke M,Koki T,et al.A comparative study between flow visualization and comutional fluid dynamic analysis for the sun medical centrifugal blood pump[J].Artif Organs,2004,28(5):25-28.

[13] 李永博.基于CFD模擬模型的溫室溫度場(chǎng)均勻性控制[J].農(nóng)業(yè)機(jī)械學(xué)報(bào),2012,43(4):156-161.

[14] 吳飛青,張立彬,胥芳,等．機(jī)械通風(fēng)條件下玻璃溫室熱環(huán)境數(shù)值模[J]．農(nóng)業(yè)機(jī)械學(xué)報(bào),2010,41(1):153-158．

[15] 程秀花,毛罕平,伍德林,等．玻璃溫室自然通風(fēng)熱環(huán)境時(shí)空分布數(shù)值模擬[J]．農(nóng)業(yè)機(jī)械學(xué)報(bào),2009,40(6):179-183．

[16] Fidaros D K,Baxenanou C A,Bartzanas T,et al.Numerical simulation of thermal behavior of a ventilated arc greenhouse during a solar day[J].Agricultural and Forest Meterology,2009,149(6-7):1050-1062.

[17] 侯翠萍,馬承偉.FLUENT在研究溫室通風(fēng)中的應(yīng)用[J].農(nóng)機(jī)化研究,2007(7):5-9.

[18] 邢蘭.傾斜式日光溫室冬季熱環(huán)境的模擬研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2006.

[19] 陶文銓.數(shù)值傳熱學(xué)[M].西安:西安交通大學(xué)出版社,2001.

[20] Swinbank W C.Long-wave radiation from clear skies[J]．Quarterly Journal of the Royal Meteorological Society,1963,89(381):339-348．