基于CFD的植物工廠管道通風模擬及優化

吳晨溶 程瑞鋒 方 慧* 楊其長 張 晨

(1.中國農業科學院 農業環境與可持續發展研究所，北京 100081;2.農業農村部設施農業節能與廢棄物處理重點實驗室，北京 100081)

植物工廠是指通過設施內高精度環境控制實現農作物周年連續生產的高效農業系統，利用計算機對植物發育過程的溫度、濕度、光照、CO2濃度以及營養液等環境條件進行自動控制，使設施內植物生長發育不受或很少受自然條件制約的省力型生產方式[1]。由于植物工廠充分運用了現代工程、生物、環境和信息等手段，技術高度密集，近年來一直被國際上公認為是現階段設施農業的最高發展階段，是衡量一個國家農業技術水平的重要標志，受到世界各國的高度重視[2]。

植物工廠內部氣流場的分布情況對植物生長有重要影響[3]。研究表明[4-6]植物工廠內部氣流速度在0.01～1.30 m/s范圍內，葉片光合作用及蒸騰作用隨著氣流速度的增長逐漸增強，更適宜作物生長。但在傳統通風模式的植物工廠中，內部氣流往往分布不均，尤其是栽培架對風速的阻擋，作物光合作用和蒸騰作用受到抑制，導致作物產生疾病，其中以燒尖現象最為普遍，嚴重影響作物的產量和品質，降低了植物工廠經濟效益[7]。大量研究表明，增加作物冠層氣流可在一定程度上改善作物冠層微環境，促進作物光合作用和蒸騰作用。Goto等[8]設計了一種在作物冠層上方的垂直通風模式，采用單根主管連接6根支管的方式，將氣流從主管通向支管氣孔自上而下吹入作物冠層，該技術在一定程度上改善了作物冠層氣流環境。Shibata等[9]同樣設計了垂直氣流系統，該系統主要由1臺風機、1根集氣管和6根通風管組成，每根通風管分布9個氣孔，風機吹出的空氣經過集氣管吹入通風管，最后通過氣孔吹出，該方法提高了葉片的蒸騰速率并在一定程度上抑制了燒尖現象的發生。Nishikawa等[10]將作物放置于獨立的旋轉平臺上，以2 r/min的速度勻速旋轉，一定程度上打破了冠層邊界層阻力，作物產量比靜置時提高了20%。Hesham等[11]提出多風扇系統，將7個風扇分別安裝在栽培層的前后側，產生來自相反水平方向的氣流，改善了氣流環境，促進了植物工廠內作物生長，緩解作物燒尖現象。李琨等[12]設計了一套根際通風系統，空氣經風機導入栽培板與營養液面中間的空氣層中，經栽培板定植孔周圍預留的通風孔自下而上排出，該系統在解決傳統環控方式通風溫控不均勻的同時，對地上部及地下部多種微環境參數調控起到了積極作用。盡管上述通風方式對氣流均勻性有一定積極作用，但試驗周期較長，且有的需借助復雜的機械裝置來實現通風。

近年來，計算流體力學軟件(Computational fluid dynamics，CFD)以周期短、成本低、準確性高的特點被廣泛應用于溫室和植物工廠的氣流場與溫度場的模擬[13]。Lim T等[14]利用CFD模擬了植物工廠內中間過道4種不同進出口位置下氣流分布情況。Baek M S等[15]借助CFD在人工光植物工廠中將空調、排風扇、內部風扇的開啟設置為不同案例，模擬出植物工廠內部溫度場及氣流場，并選出最適宜植物生長的通風案例。張芳等[16]通過CFD找到了自然通風條件下大跨度溫室內溫度場和氣流場的分布規律，以及不同室外風速條件下通風口開度對大跨度溫室溫度場和氣流場的影響。劉煥等[17]利用CFD找到了人工光植物工廠內最佳氣流循環模式。Zhang Y等[18]設計了帶有三排氣孔的空氣管，空氣垂直向下流到作物冠層表面，并借助CFD軟件模擬植物冠層的氣流場，最終確定2根空氣管的方案氣流最為均勻。

本研究擬將植物工廠中的單層栽培架為研究對象，設計管道通風方案，并設置不同進風口風速進行正交試驗，利用CFD仿真軟件對不同方案下的作物冠層表面氣流場進行模擬仿真，分析不同管孔直徑、管孔數量及進風口風速對作物冠層氣流分布的影響，確定最佳管道通風方案，以期為改善植物工廠微環境和緩解作物燒尖現象提供一定的理論參考。

1 材料與方法

1.1 試驗裝置

試驗單層栽培架位于中國農業科學院農業環境與可持續發展研究所植物工廠內，單層栽培架規格(長×寬×高)為1.4 m×0.7 m×0.4 m，采用穿孔聚乙烯管道進行通風，管道直徑為0.05 m，試驗作物為生菜(Lactucasativacv. Tiberius)。本試驗只考慮作物冠層和冠層上部空氣層，燈板到栽培板的距離一般為0.3 m，成熟期的生菜平均株高為0.1 m[19]，作物冠層和冠層上部空氣層的高度分別設為0.1 m和0.2 m，管道架在栽培架外側，管道氣孔水平方向的中心線與空氣層和作物冠層相交切線重合。在Workbench中的Geometry模塊中創建單層植物栽培架模型(圖1)，以模型東北角為原點O(0,0,0)正南方向為X軸正向，垂直向上為Y軸正向，正西方向為Z軸正向。

1.空氣層；2.作物冠層；3.管道1.Air interior; 2.Crop interior; 3.Tube圖1 試驗單層栽培架結構圖Fig.1 Schematic diagram of single shelf cultivation frame

1.2 試驗設置

在溫室內使用穿孔聚乙烯管道來分配空氣已變得越來越普遍，但在穿孔管道設計中，氣孔的直徑、數量、形狀、間距等因素都會影響空氣流通率和氣流速度。在Well C的試驗中得出，氣孔截面積之和與管道口面積之比小于1.5時風速排放更為均勻[20]。為滿足合適的管孔徑比，本研究通過計算設置了3組方案，所有方案的管道長度都為1.4 m，管道其余參數見表1，并在Workbench中的Geometry模塊中創建了3組方案的管道模型(圖2)。通風管道的氣流由進風口進入管道后經氣孔流入作物冠層內部，進風口風速分別設為2、3、4和5 m/s進行模擬。Kitaya等[21]研究表明，作物冠層表面的適宜風速為0.3～1.0 m/s，所以在本研究中，按風速將植物冠層表面劃分為3個區域：停滯區(≤0.3 m/s)、適宜區(0.3～1.0 m/s)和高速區(≥1.0 m/s)。試驗所用泡沫栽培板尺寸(長×寬×高)為0.7 m×0.65 m×0.014 m，栽培密度17株/m2，每個栽培板16個定植孔(4行4列)，定植孔直徑為0.03 m，選用8個定植孔定植生菜，在栽培板上呈三角形分布。

表1 不同管道通風方案的管道基本參數Table 1 Basic parameters of different tube ventilation cases

圖2 不同管道通風方案的管道模型Fig.2 Tube models of different tube ventilation cases

1.3 測點布置

試驗中選取X=0.35 m截面，在Y=0.5、1.0、1.5、2.0、2.5 m高度處各均勻布置4個測點，共20個測點(P1，P2，…，P20)；選取Z=0.4 m截面，在Y=0.5、1.0、1.5、2.0、2.5 m高度處各均勻布置2個測點，共10個測點(P21，P22，…，P30)。2個截面共30個測點測量風速值進行模擬驗證(圖3)。氣流速度測量使用熱線風速儀(Climomaster 6501-BG，日本加野麥克斯公司)，風速測量范圍為0.01～5.00 m/s，風速測量公差為標準值的±2%，工作溫度范圍為-20～70 ℃，相對濕度范圍為2%～98%，氣壓范圍為-5～5 kPa。

P1，P2，…，P30為風速測點P1，P2，…，P30 are the measuring points of air velocity圖3 X=0.35 m截面(a)和Z=0.4 m截面(b)風速測點布置Fig.3 Air velocity measuring points of sections of X=0.35 m (a) and Z=0.4 m (b)

2 數學模型

2.1 CFD建模與網格劃分

利用Workbench中的Geometry模塊和ANSYS Mesh模塊對栽培裝置進行幾何建模和網格劃分，采用FLUENT 19.2對栽培區域進行模擬計算。在Geometry模塊中創建與單層栽培架原型相同尺寸的栽培裝置模型，包括作物冠層、空氣層和管道3個部分。將構建的模型輸入到ANSYS Mesh模塊中，利用Proximity and Curvature方法進行網格劃分(圖4)，管道部分網格進行局部加密，3組方案模型的網格平均偏斜度和最大偏斜度分別為0.732 5和0.74、0.747 5和0.76、0.722 5和0.74，網格數為465 453～502 345(表2)，網格質量均為優[22]。

圖4 栽培裝置網格劃分圖Fig.4 Mesh in single shelf cultivation frame

表2 不同管道通風方案的網格參數Table 2 Mesh quality of different tube ventilation cases

2.2 控制方程

將栽培層內的通風過程視為三維、穩態、不可壓縮的湍流流動過程；湍流模型選用k-ε模型(realizable)；使用有限體積法對控制方程進行離散；用SIMPLE算法(壓力耦合方程組的半隱式方法)計算求解壓力耦合方程；氣流數值計算按照流體力學守恒數學公式計算，即連續性方程、動量方程、能量方程，符合質量、動量、能量守恒定律[23]，公式如下：

連續性方程

(1)

動量方程

(2)

(3)

(4)

能量方程

(5)

式中：ρ為密度，kg/m3；t為時間，s；v為速度矢量；u，υ，ω分別為流體質點速度在x，y，z方向上的分量；μ為黏性系數；S為源項；T為溫度，K；k為導熱系數，W/(m·K)；Cp為比熱容，J/(kg·K)；p為壓力，Pa。

2.3 參數設置及計算方法

將植物栽培裝置模型進行簡化，作物冠層上部空氣層設為空氣介質，作物冠層設為多孔介質，黏性阻力系數和慣性阻力系數分別為25和1.3[19]；氣流由進風口進入管道內部通過氣孔吹入作物冠層。管道進風口設置為速度入口(Velocity-inlet)，進風速度分別為2、3、4、5 m/s；栽培裝置四周氣流出口設置為壓力出口(Pressure-outlet)，壓力為0 Pa，所有材料的熱物性參數見表3。選用SIMPLE算法進行壓力速度耦合運算?？臻g離散設置中，梯度項選用最小二乘法；動量、能量和粘性項都選用一階迎風格式，使其更快收斂，將能量項的松弛因子設置為10-6，動量項和粘性項皆設置為10-3，迭代步數設為300步。

表3 材料的熱物性參數Table 3 Thermal physical parameters of materials

3 CFD仿真結果

3.1 CFD模型驗證

將30個測點的風速模擬值與實測值進行對比(圖5)，平均絕對誤差為0.02 m/s，均方根誤差為0.024，其中誤差最大的測點是第9個測點，實測風速值為0.15 m/s，CFD模擬風速值為0.085 m/s，絕對誤差為0.065 m/s?？傮w上看，CFD風速模擬值與風速實測值吻合度良好。因此本試驗所構建的單層栽培架實際內部氣流場分布情況與模擬情況吻合良好，能夠進行數據分析與和后續的模擬試驗。

圖5 單層栽培架內風速實測值與模擬值對比圖Fig.5 Comparison between measured and simulated air velocity in the single shelf cultivation frame

3.2 6孔通風管道作物冠層氣流分布

對6孔通風管道方案進行模擬，作物冠層氣流分布云圖見圖6。氣孔附近風速最大，作物冠層平面氣孔位置呈6個氣流束，風速沿氣流束逐漸減弱。除6個氣流束所在位置外，作物冠層平面其他區域風速較小，氣流在氣孔與氣孔之間存在較大面積的停滯區，且在冠層平面邊緣管道兩端位置存在部分零風速區域，整體風速并不均勻。隨著進風口風速增大，作物冠層高速區面積也逐漸增大，停滯區面積逐漸減小。在進氣速度為2、3、4和5 m/s時，作物冠層平面風速平均值分別為0.28、0.43、0.56和0.66 m/s，適宜區的面積比分別為32.7%、40.1%、37.9%和34.2%，其中當進風口風速為3 m/s時，適宜區面積比最大為40.1%(表4)。

圖6 6孔通風管道不同進風口風速(v)下作物冠層表面風速分布云圖Fig.6 The contours of canopy surface air velocity under different inlet air velocity of 6 holes tube

表4 6孔通風管道不同進風口風速下作物冠層表面各區域的面積比Table 4 The area percentage of each area at the canopy surface under different inlet air velocity of 6 holes tube

3.3 12孔通風管道作物冠層氣流分布

對12孔通風管道方案進行模擬，結果見圖7：最大風速仍分布在氣孔附近，當進風口風速為 2 m/s 時，作物冠層平面氣孔位置呈12個氣流束且風速較小，與之相比，當進風口風速為3、4和5 m/s時，氣流分布較均勻，隨著進風口風速的逐漸增大，氣孔間的停滯區面積減小，作物冠層適宜區的面積增大，但冠層平面邊緣管道兩端位置仍存在零風速區域。作物冠層平面風速平均值分別為0.34、0.52、0.71和0.89 m/s，適宜區的面積比分別為47.0%、68.4%、62.2%和52.9%，其中當進風口風速為3 m/s時，面積比最大為68.4%(表5)。

圖7 12孔通風管道不同進風口風速(v)下作物冠層表面風速分布云圖Fig.7 Contours of canopy surface air velocity under different inlet air velocity of 12 holes tube

表5 12孔通風管道不同進風口風速下作物冠層表面各區域的面積比Table 5 The area percentage of each area at the canopy surface under different inlet air velocity of 12 holes tube

3.4 25孔通風管道作物冠層氣流分布

25孔通風管道方案模擬結果(圖8)表明，當進風口風速為2 m/s時，作物冠層平面氣孔位置呈較明顯的氣流束，與之相比，當進風口風速為3、4、5 m/s 時，風速分布沒有氣流束，而是呈片狀梯度分布，氣流分布較均勻，隨著進風口風速的逐漸增大，氣孔間停滯區面積降低，冠層平面邊緣管道兩端零風速面積也相應減小。作物冠層平面風速平均值分別為0.36、0.55、0.74和0.94 m/s，適宜區的面積比分別為54.4%、86.1%、85.8%和75.0%，其中當進風口風速為3 m/s時，面積比最大為86.1%(表6)。

圖8 25孔通風管道不同進風口風速(v)下作物冠層表面風速分布云圖Fig.8 Contours of canopy surface air velocity under different inlet air velocity of 25 holes tube

表6 25孔通風管道不同進風口風速下作物冠層表面各區域的面積比Table 6 The area percentage of each area at the canopy surface under different inlet air velocity of 25 holes tube

3.5 氣流場對比分析

由冠層平面風速分布云圖(圖6～8)對比可知，氣孔孔徑大數量少則風速范圍差異大，氣孔孔徑小數量多則風速范圍差異小。當風速為5 m/s時，25孔通風管道作物冠層平面風速平均值最高，為0.94 m/s，適宜區的面積比為75.0%，而6孔和12孔通風管道分別為34.2%和52.9%；當風速為4 m/s 時，25孔通風管道作物冠層平面風速平均值最高，為0.74 m/s，適宜區的面積比為85.8%，而6孔和12孔通風管道分別為37.9%和62.2%；當風速為3 m/s時，25孔通風管道作物冠層平面風速平均值最高，為0.55 m/s，適宜區的面積比為86.1%，而6孔和12孔通風管道分別為40.1%和68.4%；當風速為2 m/s時，25孔通風管道作物冠層平面風速平均值最高，為0.36 m/s，適宜區的面積比為54.4%，而6孔和12孔通風管道分別為32.7%和47.0%。3種管道在進風口風速逐漸降低的情況下，作物冠層平面平均風速也隨之逐漸降低，但適宜區的面積比逐漸上升，在進風口風速為 3 m/s 時均達到最高值，在進風口風速為2 m/s時，適宜區的面積比均下降。綜合考慮，進風口風速為3 m/s 的25孔管道通風方案最佳(圖9)。

圖9 不同管道通風方案在不同進風口風速下作物冠層表面適宜區的面積比Fig.9 Percentage of suitable area at canopy surface of different tube ventilation cases under different air velocity

4 結論與討論

本研究利用CFD模擬了植物工廠內單層栽培裝置作物冠層內部的氣流分布情況，并將模擬值與實測值進行對比，風速平均絕對誤差為0.02 m/s，均方根誤差為0.024，CFD風速模擬值與實測值吻合，模擬效果良好。

基于植物工廠內作物冠層管道通風氣孔的孔徑與數量設計了3種(6孔、12孔和25孔)管道通風模式，4種不同的進風口風速(2、3、4和5 m/s)進行正交試驗。在Hesham的試驗中，傳統通風模式冠層表面的平均風速為0.23 m/s，停滯區面積比為70%[11]。本試驗所有管道通風方案的平均風速均高于傳統通風模式，停滯區面積占比均低于傳統通風模式，即所有管道通風方案的平均風速及氣流均勻性均優于傳統通風模式。當進風口風速為3 m/s時，25孔通風管道冠層平面適宜區的面積比為86.1%，為最佳管道通風方案。

但當進風口風速為4 m/s時，25孔通風管道冠層平面適宜區的面積比為85.8%，停滯區為8.7%，而高速區面積比為5.5%，與進風口為3 m/s時的面積比進行對比，差值分別為0.3%、2.5%和-2.8%；由此可看出，當進風口風速為4 m/s時，冠層平面適宜區和停滯區的百分比都有所減小，但差值較小，雖然高速區的面積比有所增加，但都聚集在氣孔附近，對作物冠層產生影響較小，所以下一步試驗將進風口風速分別設為3和4 m/s進行模擬。

本試驗只研究了在固定管道高度和角度條件下，不同進風口風速、孔徑和孔數對作物冠層平面氣流分布的影響，但在實際生產中LED燈的散熱會對栽培裝置內部的氣流產生一定影響。因此，下一步試驗將嵌入LED燈熱源模型，研究管道與冠層栽培板的高度以及氣孔的傾斜度對作物冠層平面和冠層內部氣流分布的影響，從而更加精確地模擬栽培裝置內部環境，最終確定最優的管道通風方案。