基于CFD的玻璃溫室環境數值模擬及優化分析

賈鶴鳴，韓駿騁，張森，孫康健，李瑤

東北林業大學 機電工程學院，黑龍江 哈爾濱 150040

溫室植物的產量取決于溫室內環境的溫度、空氣濕度等小氣候要素的適宜度[1]，尤其是較好的溫度和氣流分布，對于植物生長有著極其重要的作用。在所有環境因子中，適宜的溫度是植物生長的首要條件，過低的溫度甚至會造成溫室內植物的死亡。因此，溫室通風口的設置和風速都會對溫室內氣流的分布情況有較大的影響。程秀花等[1]提出室內氣流的速度和分布情況較易受到植物的影響；汪鄭邦[2]認為植物溫室的氣流分布對作物生長有較大的影響，不均勻的氣流分布容易在作物周圍形成換流，影響作物的蒸騰和光合作用；劉文合等[3]的研究表明通過對溫室進行升溫，可以使得溫室內氣溫顯著升高，但對溫室內土壤深處的溫度影響不大；周偉等[4]通過試驗證明了通風口的設置和風速都會對溫室內氣流的分布情況有較大的影響，同時溫室的氣流分布對植物蒸騰和光合作用有較大的影響。因此，如何合理優化溫室結構，得到適宜的風速和溫度并使兩者均勻性提高是當前急需解決的重要問題。

1 材料和方法

1.1 溫室模型

試驗在東北林業大學進行，模擬溫室建立在室內，受外界環境影響因素較小。溫室的長×寬×高為2 m×0.83 m×1.8 m，由于溫室是密閉的，因此忽略溫室的門，并對溫室四周及頂部進行簡化，即將溫室的頂部和四周均視為以浮法玻璃為材料。考慮到植物的冠層阻力，將植物模型簡化為長1.50 m、寬0.6 m、高0.5 m的多孔介質，置于溫室正中央。暖氣扇規格為0.4 m×0.4 m，置于溫室西側，離地面1.3 m高處。采用Gambit進行以溫室西北角為零點的溫室模型的建立和網格化，得到的玻璃溫室模型如圖1所示。

圖1 玻璃溫室模型

1.2 監測點布置

在溫室中選取Z=0.415 m處的截面，均勻分布6個監測點，高度分別為Y=0.6 m和1.2 m，將監測點依次記為P1,P2,···,P6。在平面X=1 m處，共建立4個監測點，高度分別為0.6 m和1.2 m，記為P7、P8、P9、P10。通過將在監測點處得到的實際值與模擬值進行對比驗證。采用風速儀和溫度計對溫室內風速和溫度進行測量，溫室環境模擬值由Fluent軟件進行模擬，監測點的布置如圖2、3所示。

圖2 Z=0.415 m截面監測點布置（Y=0.6 m和Y=1.2 m）

圖3 Z=0.415 m截面監測點布置（X=1 m）

1.3 數學模型

1.3.1 CFD模型

當Fluent模擬溫室環境時，使用標準k-ε模型[5]將溫室內的氣流視為遵循質量守恒方程、動量守恒方程和能量守恒方程的湍流。

動量守恒方程指出，所有的流動在實際情況下都應該遵循動量守恒定律，定律用公式表示為：

式中的u、v、w分別為在x、y、z方向上溫室內氣流流速的分量

動量守恒用公式表現為：

式中：μeあ代表有效黏度；β代表空氣熱膨脹系數；g代 表重力加速度；T代表空氣溫度；Tref代表空氣的參考溫度。

能量守恒用公式表示為

式中：q代 表熱源；代表比熱容； λeあ代表有效導熱系數[6]。

在處理封閉室內因溫差引起的浮升力項時采用Boussinesq假設[7]。Boussinesq假設是對流體流動現象的一種簡化，假設認為：流體中的黏性耗散忽略不計；除密度外其他物性為常數；對密度僅考慮動量方程中與體積力有關的項，其余各項中的密度亦作為常數。重力項中的密度如若采用作參考溫度，則在不同溫度T下溫室內部的空氣密度ρ可以用數學方程表示為：

式中：β是溫室內氣流的熱膨脹系數；T為溫室內的溫度；T0為溫室內部空氣的參考溫度。

流經多孔介質的空氣可通過Darcy-Forchheimer理論描述為

式中：Sφ代表動量源項；Kp代表多孔介質的滲透性系數；CF代表非線性動量損失系數；μ代表空氣的動力黏度；ρ代表空氣密度；u代表空氣流速[8]。

將植物視為各向同性的多孔介質[9]，并建立其與溫室內部風速的數學關系模型：

式中：ILAV代 表葉面積指數；CD是植物冠層的阻力系數。

對于多孔介質動量源項，使用冪律模型，自定義經驗系數C0為0.395，C1為0.2。植物通過呼吸作用產生的熱量與溫室中氣流進行熱交換，將植物通過呼吸作用產生的熱量視為常值，在進行數值模擬時將植物區域看作“體積熱源”。

1.3.2 計算域與網格劃分

在使用Gambit劃分網格時將整個玻璃溫室作為計算域，采用非結構化四面體網格(混合TGrid網格類型以及Tet/Hybrid混合單元)將溫室整體劃分為四面體混合網格單元[10]，并局部加密暖風扇出風口和植物區域，以適應氣流和溫度可能出現的較大變化[11]。

1.3.3 邊界條件

因為溫室密閉性較好且置于室內，不考慮太陽輻射和漏風的影響，對模型進行簡化。將暖風扇側的邊界條件設為速度入口條件，吹出水平方向25 ℃暖風，平均速度為1 m/s；覆蓋層全部簡化為玻璃，設為壁面條件，其與整個溫室進行能量交換的方式主要為對流傳熱和傳導[10]，初始溫度設為20 ℃；為簡化能量交換過程，在進行模擬時不考慮其他的傳熱過程。將擁有大的熱容量的土壤的邊界條件設為壁面[12]，初始溫度設為20 ℃。

1.3.4 求解運算

在進行模擬時，為修正速度，選用SIMPLEC算法：梯度項選用最小二乘法；壓力項采用標準算法；能量和動量采用二階迎風；為加快收斂速度，黏性和湍流選用一階迎風。

2 FLUENT仿真結果與驗證

2.1 FLUENT仿真結果

用Fluent進行模擬的結果如圖4、5所示。在吹入熱風后，暖風口周圍區域的氣體流速變化明顯，且存在較大梯度。氣流從暖風扇水平送入溫室后，略微向下方運動，在溫室內部的東面、暖風扇正下方和溫室底部形成渦流，且溫度與風速有著直接的關系。

圖5 Y=0.5 m和Y=0.1 m處溫度云圖

選取溫室內植物冠層平面Y=0.5 m和植物內部平面Y=0.1 m的速度云圖（圖4）進行分析。當風速為0.30～1.00 m/s時，植物生長較好[13]。在使用Fluent進行模擬計算后，在距離地面0.5 m處的平均空氣流速為0.152 m/s，最低速度為0.009 m/s，低于最佳風速；而植物冠層部分沒有處于最適風速的區域出現。在平面Y=0.1 m處氣流平均流速為0.195 m/s，最小流速為0.013 m/s，低于適宜風速的最小值0.3 m/s，最大風速為0.552 m/s，雖略大于0.3 m/s，但植物內部達到適宜植物生長的風速的部分僅占整體的20%。綜上所述，此時的氣體流動不適于植物的生長。

選取溫室內植物表面Y=0.5 m和植物內部Y=0.1 m的溫度云圖（圖5）進行分析。在暖風口附近溫度最高，溫度與其他區域相比有很大的梯度。在2個平面內最高溫度可達24.45 ℃，平均溫度23.17 ℃，溫差的最大值為1.28 ℃，此溫度較適于大部分植物的生長。

2.2 對比驗證

在試驗期間比較10個監測點的實際測量值和使用Fluent模擬得出的計算值，其結果如圖6、7所示。

圖6 風速在監測點模擬值與實測值的比較

圖7 溫度在監測點模擬值與實測值比較

根據圖6，對比不同監測點的風速，所有觀測點的絕對誤差的平均值為0.008 m/s，相對誤差的平均值為5.54%，絕對誤差和相對誤差的最大值出現在監測點P6處，分別為0.027 m/s和24.228%。雖然溫室內風速的模擬值小于測量值，但大部分情況下，模擬結果與實際情況一致。因此可以認為模型的速度計算值可以較好地體現溫室內部氣體流速的情況。

根據圖7，對比不同監測點的溫度，10個監測點絕對誤差均值為0.44 ℃，相對誤差均值為4%，P10處出現最大絕對誤差1.66 ℃和最大相對誤差2.78%，模擬值與從監測點測出的值差別不大，出現誤差的原因可能是溫室模型密封不夠好而造成漏風，導致室溫略低。結合儀器位置和人員干擾對實驗結果造成的影響[14-15]，從整體來看，溫度的實際值和計算值誤差很小，因此可以將溫度的模擬值視作溫室內溫度的實際值。

3 溫室內環境的優化設計

3.1 優化方案設計

溫室內的速度場、溫度場的均勻性受暖風扇的影響較大，因此通過在溫室中改變進風口的速度和改變暖風口的位置來設計出3種改進方案，以滿足植物生長的需要。方案一，升高暖風口速度，即暖風口速度從1 m/s升高到2 m/s；方案二，下移暖風口的位置，將暖風口從1.3 m處下移到1 m處，如圖8所示；方案三，下移暖風口的位置至0.7 m處，如圖9所示。將原溫室作為方案零來進行對比，溫室、植物和暖風扇的尺寸視為定量不做改變，且方案一、方案二和方案三的網格劃分與方案零相同，其余邊界條件的設置也相同。通過對不同方案使用Fluent進行模擬計算，比較植物冠層表面和內部的風速及溫度分布情況，依此來選擇擁有適宜風速溫度的最適宜植物生長的方案。

圖8 暖風扇離地1 m時的玻璃溫室模型

圖9 暖風扇離地0.7 m時的玻璃溫室模型

3.2 優化模擬分析

3.2.1 速度場模擬分析

經過改良后的方案氣流分布模式基本相同：在溫室上部和底部都有渦流的產生，暖風扇口附近的氣流流速最大，并隨著流動方向減弱，直到在東側壁面附近開始沿溫室頂部返回，溫室西側的植物表面部分和溫室壁面附近風速較小，且會在植物內部的渦流處出現整個植物區域中風速的最大值。

在方案一中，從圖10速度云圖上看，選取的2個平面空氣流速平均值為0.695 m/s，均高于方案零的0.152 m/s，且滿足植物生長的最佳風速。風速在南北兩側墻處達到最大，為1.562 m/s，遠高于最適宜植物生長的風速1 m/s，這樣的速度容易對植物造成損傷，此時植物表面風速位于0.3～1.0 m/s的部分約占整個植物冠層的63%。在植物內部Y=0.1 m處最大風度為1.581 m/s，平均風速0.928 m/s，處在最佳氣流流速范圍內的區域約占27.9%，剩下的大部分區域風速均大于1 m/s，不利于植物生長，由此可知不能僅增大風速來優化溫室內部的氣流分布情況。

圖10 方案一速度場模擬

從圖11方案二的云圖可知，y=0.5 m處氣流的流速均值為0.284 m/s，最大為0.632 m/s；y=0.1 m處平均風速為0.423 m/s，最小風速為0 m/s，最大風速為0.754 m/s，均小于方案一。處于適宜風速的占25%，雖然相對于方案零適宜風速的占比有所增加，氣流的均勻性與方案一相比較好，但適宜風速的占比仍小于方案一。

圖11 方案二速度場模擬

分析圖12方案三的云圖，選取分析的平面平均氣體流速為0.383 m/s，最高流速可達0.834 m/s，均符合最佳風速。在植物部分的最低風速為0.07 m/s，大于方案零和方案二，但仍低于0.3 m/s。這時的植物冠層風速達到適宜風速的占70%，高于其他方案，且均勻性最好。

圖12 方案三速度場模擬

3.2.2 溫度場模擬分析

分析方案一、方案二和方案三的溫度情況（圖13～15），并將其與速度云圖進行比對，發現氣流速度較大的區域溫度上升較多。

圖13 方案一溫度場模擬

圖14 方案二溫度場模擬

圖15 方案三溫度場模擬

3個方案的最高溫度分別為24.49、24.45和24.72 ℃，平均溫度為 23.93、23.15 和 23.49 ℃，均高于方案零，且植物區域的溫度分布均勻性優于方案零。就溫度而言，方案一的均值和最大值都為最高，其次為方案三。

通過模擬分析實驗中設計的3種溫室的優化方案，可以發現：1）方案一的氣流均勻性稍差，且有較多的區域風速高于最佳風速，不利于植物生長；2）方案二雖然略優于方案零，但對于氣流場和溫度場并不能提供有效的優化；3）方案三的氣流均勻性最好，不存在流速大于1 m/s的區域，從而不會對植物的正常生長造成威脅，且風速大于方案零和方案二?？傮w來講，方案三的暖風扇位置更適宜植物的生長。

4 結論

1）溫度數值模擬結果與實驗實測結果吻合良好，表明建立的CFD模型有效，模擬玻璃溫室內環境分布是可行的；

2）設計的3種結構優化方案的風速值及均勻性均優于方案零。綜合考慮，方案三氣流循環模式最佳，氣流均勻性最好，適宜植物生長。

本文僅基于風速和暖風扇位置進行溫室的優化，在之后的研究中，可考慮將添加排氣扇和增加暖風扇的數量作為溫室的優化方式，考慮植物在陽光下的光合作用產生的熱量和外部輻射等會影響溫室內部氣流和溫度的條件，探究溫室內部氣流與溫度的分布情況。