不同二氧化碳含量對馬鈴薯貯藏室內溫度分布的影響

甄琦， 閆廣澤， 塔娜， 趙志勇， 于慧敏

（內蒙古農業大學能源與交通工程學院，呼和浩特 010018）

馬鈴薯作為常見的糧菜飼兼用作物之一，目前已成為中國的第四大主糧，對保障國家糧食安全意義重大，如何延長其保鮮周期是農產品貯藏領域的熱點研究方向[1-4]。半地下式貯藏室作為中國北方地區農戶普遍采用的貯藏設施，在馬鈴薯貯藏過程中起著至關重要的載體作用，其優點是建設成本低、調控溫濕度簡單、操作方便等[5-6]；但該類貯藏室無排氣和通風裝置，導致貯藏環境與馬鈴薯的生理特性出現偏差從而影響馬鈴薯的貯藏質量[7-10]。二氧化碳作為馬鈴薯貯藏呼吸熱的副產物，其溫室效應對密閉貯藏室內溫度的分布有較大影響。已有學者基于二氧化碳的含量問題[11-13]展開了研究，為各類果蔬在二氧化碳影響下的貯藏條件和方法提供了研究基礎。若能平衡二氧化碳對貯藏室的溫室效應與抑制馬鈴薯呼吸作用的關系，可從溫度控制及生理抑制2 個角度提升馬鈴薯貯藏期間的保鮮效果，所以研究不同含量的二氧化碳對貯藏室內溫度分布的影響有重要的意義。

計算流體力學（computational fluid dynamic，CFD）是一種通過計算機解析數學模型流體流動和傳熱問題的可靠工具[14]。目前有較多學者利用CFD 研究冷庫、貯藏室等農業設施室內的氣體流速分布與溫度分布[15-17]。He等[18]建立了11跨塑料溫室大棚的三維CFD 數值模型，研究通風口配置和開口尺寸對溫室小氣候模式及內部氣候時空變化的影響。Zhang 等[19]通過模擬單個貨架生產系統中的生長環境，設計并提出了1 種改進的空氣循環系統，有助于防止生菜生產中發生變質。Echaroj等[20]采用CFD 數值仿真對3種幾何形狀的儲藏室進行了數值模擬，得到了不同幾何形狀儲藏室的氣流速度分布，發現貯藏室幾何形狀對氣流速度分布的影響。目前，氣體含量對果蔬貯藏質量影響的研究主要集中在氣體對果蔬相關生理指標的作用上[21-22]，而二氧化碳氣體造成的溫室效應對密閉環境的貯藏室內溫度變化影響較為顯著，卻鮮有相關的研究。

本研究綜合考慮了室內馬鈴薯堆垛呼吸熱、二氧化碳含量等多方面因素，以我國西部地區典型的半地下式馬鈴薯貯藏室為載體，構建三維CFD 數值計算模型，將仿真計算的結果與試驗測得的數據進行比較，驗證相關數學模型的有效性和準確性，基于該模型討論不同含量二氧化碳貯藏室內溫度分布情況。

1 材料與方法

1.1 試驗材料

本研究以呼和浩特市內蒙古農業大學農學院試驗基地半地下式貯藏室為試驗載體，該貯藏室主要貯藏食用型馬鈴薯，外觀如圖1A 所示。貯藏室為東西走向，整體尺寸為長40.0 m，寬11.5 m，高2.6 m，其內部共分為10 個南北走向的單間貯藏室，每個貯藏室的尺寸為8.0 m×4.0 m×2.6 m，本研究選取1 個貯藏室作為試驗場所。貯藏室的墻體為磚墻，厚度為0.37 m，貯藏室有1.87 m 處于地平面之下，北墻外部堆積有厚度為1.0 m 土層，南墻開2.0 m×1.5 m 的大門，東西居中且距門1.60 m，貯藏室底部距外地表深度為1.56 m。

圖1 半地下式貯藏室外觀、尺寸及設備布置Fig. 1 Appearance， size and equipment layout of the semi-underground storage room

圖2 馬鈴薯貯藏室幾何模型Fig. 2 Diagram of numerical model of potato storage room

以‘冀張薯12 號’作為供試材料，試驗周期為2022 年9 月25 日（入庫）至2023 年3 月2 日（出庫），貯藏室均按照正常的農業生產進行貯藏保溫，半地下式貯藏室的尺寸及設備布置示意圖如圖1B~D所示。

1.2 試驗方法

將澤豐ZFC-0型二氧化碳發生器（體積1.30 m×0.71 m×0.62 m，出氣口直徑1.2 cm，出氣速度2 m·s-1，中國）布置于貯藏室北墻的中心位置，用于保持貯藏室有穩定的二氧化碳氣體含量。馬鈴薯以1.35 m×1.10 m×0.91 m 的近立方體結構堆積在貯藏室中心位置，分別將5 個DigiTH 型溫度傳感器（量程-40~80 ℃，精度≤±2 ℃，分辨率0.01 ℃，中國）布置在馬鈴薯堆垛的5 個外表面（除了底面）的中心位置，分別命名為北側T1，南側T2，西側T3，東側T4，頂部T5，長、寬、高分別距馬鈴薯堆垛邊界約0.675、0.550、0.450 m（尺寸及位置參考圖1~2），傳感器以2 次·h-1的頻率記錄溫度數據，試驗用堆體及設備布置如圖3所示。

圖3 貯藏室內馬鈴薯堆體及設備Fig. 3 Diagram of potato stack and equipment in storage room

考慮到馬鈴薯入庫周期集中在我國北方冬季，通過對外界小型氣象站數據的監控，發現2022 年1 月31 日（晴，氣溫-6~-18 ℃，北風3~4 級）為典型的北方冬季天氣條件，以該天14∶00時的貯藏室內外環境因子數據作為試驗條件，對貯藏室內溫度進行數據分析，并通過數值建模對比不同二氧化碳含量下的貯藏室內溫度的變化。

1.3 半地下式貯藏室CFD建模

此模型建立是以上述貯藏室為參考模型，詳細尺寸為：長8.0 m，寬4.0 m，高2.6 m；南墻開設貯藏室門，其尺寸為2.0 m×1.5 m。在貯藏室北墻放置ZFC-0 型二氧化碳發生器，其尺寸為：長1.3 m，寬0.71 m，高0.62 m，在二氧化碳發生器的中心處設有1 個出風口，直徑為1.2 cm。馬鈴薯堆垛近似繪制為長1.35 m、寬1.10 m、高0.91 m 的幾何體，放置在貯藏室中心位置。考慮到貯藏室常閉狀態，但門與墻體間存在間隙，故出風口位于門與墻體的縫隙處。馬鈴薯貯藏室模圖如圖2所示。

1.4 數學模型及模型假設

因馬鈴薯在貯藏過程中釋放出大量的呼吸熱，將馬鈴薯堆體假設為中心具有內熱源的多孔介質[23]，考慮馬鈴薯堆體與室內空氣的傳熱作用，在貯藏過程中氣溫較低，熱輻射作用不明顯故可以忽略。本研究應用到的數學模型具體如下。

1.4.1內熱源數學模型 內熱源在流場中的傳熱主要是因為內熱源自身散熱在流場中產生熱量的傳遞。在導熱的過程中其熱量傳遞符合導熱微分方程[24]。

其中，ρ表示密度，kg·m-3；c為定壓比熱容，J·kg-1·K-1；ρc為單位體積的物體溫度升高1 K 所需的熱量，J·m-3·K-1；λ為熱導率，W·m-1·K-1；w是內熱源的發熱量，J；T表示時間為t時的瞬時溫度，K；t為時間，s；x、y、z是直角坐標下的不同的方向上的距離，m。

內熱源加熱的時間越長，其熱量擴散情況也會隨著時間變化，經過一段時間的傳遞后會使整個流場內的熱量傳遞達到物理平衡的狀態。

1.4.2多孔介質模型 由于在本研究中，馬鈴薯被視作多孔介質與外界空氣進行對流換熱，需明確多孔介質內的流體流動狀態所產生的影響及其能量、質量、熱量傳遞等狀態，涉及到的低速的流體流動的連續方程（達西定律）如公式（2），對流換熱方程如式（3）[25]。

式中，v為流體流速，m·s-1；k為滲透率系數，%；μ為流體的黏度系數，kg·m-1·s；ρa為馬鈴薯的密度，kg·m-3；ρb為馬鈴薯堆垛內部的空氣密度，kg·m-3；φa為馬鈴薯的孔隙率，%；φb為空氣的孔隙率，%；H為馬鈴薯的比焓，J·kg-1；h為馬鈴薯之間空氣的比焓，J·kg-1；keff為馬鈴薯堆垛傳熱系數，W·m-1·K-1；qr為內熱源釋放的熱量，J。

1.4.3氣體區數學模型 貯藏室內熱質傳遞過程較為復雜，為簡化計算過程，根據文獻[26]對模型做出如下假設：①室內氣體為牛頓流體且不可壓縮；②氣體物性參數為常數；③室外溫度穩定不變，忽略墻壁引起的熱質損失；④庫內氣體在壁面上無滑移。基于以上假設，并忽略由于溫度、含量差所引起的升浮力影響，室內的流場可以簡化為三維穩態、不可壓縮、粘性的湍流流場，選擇k-ε湍流模型，在直角坐標系下，聯立連續性方程、動量方程及能量方程，上述控制方程可以表示為公式（4）。

式中，φ表示通量，當φ為1時，代表該方程表示質量守恒方程，當φ為速度矢量v=[u v w]時，代表該方程表示動量守恒方程，u、v、w是速度矢量v在3 個方向的速度標量，當φ為T時，代表該方程表示能量守恒方程，T為溫度，K；Γ 為與φ相對應的廣義擴散系數;S為與φ相對應的廣義源項。

1.4.4其他條件參數設置 根據文獻[25]設置各項參數。

①邊界條件。入口參數：速度入口，風速設置為2 m·s-1；湍流強度5.88%；水力直徑0.012 m；溫度4 ℃。出口參數：壓力出口，壓力為0 Pa；湍流強度7.04%；水力直徑0.099 m；溫度0 ℃；壁面邊界條件具體設置如表1。

表1 壁面邊界條件參數設置Table 1 Setting parameters of wall boundary conditions

②多孔介質。將馬鈴薯堆體假設為多孔介質，其參數設置如下：孔隙率0.436，滲透率0.734，慣性阻力系數136 239，黏性阻力系數0.119。

③內熱源參數。將幾何模型的馬鈴薯堆垛中心設置為熱源，溫度設置為2 ℃，其材料參數見表2。

表2 材料參數設置Table 2 Material parameter settings

④入口氣體含量條件。入口氣體含量為進氣口氣體含量，參考GB/T 51124—2015《馬鈴薯貯藏設施設計規范》[27]，分別設置0.00%、0.15%、0.30%（體積分數）的二氧化碳氣體作為對比。

1.5 室內數值方法可靠性驗證

因目前尚未有關于半地下式馬鈴薯貯藏室三維數學模型的準確性研究。白通通[28]在蘋果冷藏庫中豎直壁面貼附送風模式下送風溫度與送風速度的研究比較深入和成熟，其團隊在該方面的研究中應用的模型較為準確，因馬鈴薯貯藏過程與蘋果在冷藏過程具有相似性，故采取了此方法開展室內數值方法可靠性的驗證工作。為了驗證數值模擬方法的正確性，采用尹海國等[29]的豎壁貼附送風方式的軸線風速的計算公式，

式中，u(y*)為距送風口距離為y*時的軸線速度，y*為房間高度和對應高度的差值，u0為送風速度，b為條縫風口寬度。

1.6 馬鈴薯堆體等效熱模型有效性的驗證

考慮到馬鈴薯貯藏過程中釋放出大量的呼吸熱，為進一步準確描述貯藏室內馬鈴薯堆體與環境換熱的溫度分布，將堆體虛擬為多孔介質，并在其內部添加內熱源，將馬鈴薯堆體整體等效為1 個發熱體。為研究該數值模型的精確程度，將通過傳感器所采集的貯藏室內馬鈴薯堆表面不同測點的溫度實測值與數值計算獲得的模擬值進行對比，進而驗證馬鈴薯堆體等效熱模型應用的準確性。根據1.2方法中布置的溫度測點位置（北側T1、南側T2、西側T3、東側T4、頂部T5，如圖2 所示），提取相應測點的溫度數據。并通過后處理軟件提取與試驗測點位置對應的數值仿真結果進行對比。

1.7 網格劃分及無關性驗證

在數值模擬的過程中，網格質量和網格數量都會對仿真計算產生極大的影響，網格數量增加會增加仿真計算的精度，但會造成計算資源的極大消耗。在數值仿真計算前進行網格無關性驗證。條件及參數設置不變的情況下，改變計算域內的網格密度，看相對應的網格數變化，將3 種網格數的仿真結果與試驗結果進行對比。

2 結果與分析

2.1 室內數值方法可靠性驗證結果分析

室內數值方法可靠性驗證結果如圖4 所示，通過模擬仿真得到的3 種風速的數據與文獻[28]結果進行比較，得到軸線的速度的均方根誤差分別為3.15%、2.74%、3.35%，證明了此數值方法是正確的，可以較為精準的模擬冷藏庫、貯藏室內的溫度場和流場狀態。

圖4 數值方法驗證結果Fig. 4 Numerical method verification results

2.2 馬鈴薯堆體等效熱模型有效性驗證結果分析

由圖5 可知，測點位置的模擬結果與試驗結果數據基本吻合且變化規律一致，最大相對誤差和平均相對誤差分別為9.77%和8.26%，模型的計算誤差在允許范圍之內，故將馬鈴薯堆體虛擬為具有內熱源模型的多孔介質的等效熱模型，能夠對馬鈴薯堆體發熱情況的溫度分布進行較為精確的描述，確認了馬鈴薯堆體等效熱模型有效性。

圖5 模擬值與實測值對比Fig. 5 Comparison of simulated value and measured value

2.3 網格劃分及無關性驗證分析

網格無關性驗證結果表明，在條件及參數設置不變的情況下，改變計算域內的網格密度，相對應的網格數分別為153 萬、230 萬、420 萬。將3 種網格數的仿真結果與試驗結果相對比，其橫方向截面的溫度和結果如圖6 所示。根據仿真結果，153 萬、230 萬這2 種網格數的結果相近，230 萬網格數的溫度與420 萬網格數的溫度最大誤差大于17%，而153 萬、230 萬的溫度結果的最大誤差為11%，因此可以得出153 萬網格結果較為合適，且與430 萬網格數相比網格數少，計算速度快，因此本研究采用153 萬左右的網格數。劃分網格的尺寸設置為300 mm，增長率為1.2，目標偏度0.8。并對模型的速度入口、壓力出口以及馬鈴薯堆垛位置進行了加密。

圖6 網格無關性驗證結果Fig. 6 Grid independence verification results

2.4 二氧化碳環境因子對貯藏室溫度的影響

通過不同入口氣體含量模型的模擬仿真溫度云圖，分析其對半地下式馬鈴薯貯藏室傳熱的影響。當充入的二氧化碳含量為0.00%、0.15%、0.30%時，貯藏室平均溫度與馬鈴薯堆平均溫度分別為1.34 和1.93 ℃；1.36 和1.93 ℃；1.36 和1.94 ℃。為更清晰展示所獲得的結果，在貯藏室幾何模型中分別進行3 個方向上的取面，X、Y、Z方向截面的立體圖如圖7所示。

圖7 貯藏室不同位置截面Fig. 7 Different positions in the storage room

二氧化碳氣體含量為0.15%的情況下，半地下式貯藏室內溫度分布如圖8 所示。通過添加內熱源模型的方式模擬馬鈴薯堆體在貯藏室中的發熱情況，可以看到在X=4 m 位置馬鈴薯堆體的中心位置溫度最高，熱量在堆體內部進行傳遞，繼而加熱整個馬鈴薯貯藏室。X=7.5 m 處接近出口位置，故在其切面頂部和底部溫度較低（圖8A）。Y=0.3 m 處，由于左側遠離貯藏室出口溫度較高，向右側依次遞減，貯藏室地面（Y=2.3 m處）由于接近貯藏室大門，受到外界較低溫度空氣的影響較為明顯，可以看到靠近大門位置的地面溫度較低（圖8B）。

圖8 0.15% CO2含量下貯藏室內溫度分布Fig. 8 Temperature distribution in storage room under 0.15% CO2 gas content

由圖8C 可知，不同Z切面反映出貯藏室左上方溫度較高，這是由于二氧化碳發生器氣體出口位置靠前，噴出后氣體由于溫度較高向上運動，由于左上方頂墻、左墻及二氧化碳發生器壁面的共同作用，使得較熱的氣體堆積在左上方。貯藏室大部分區域受到馬鈴薯堆體與出口氣體的共同加熱，溫度較均衡，在1.21~1.29 ℃之間，能夠較好地貯藏馬鈴薯。而受到出口溫度的影響，右側地面及右墻的整體溫度較低，該位置影響馬鈴薯貯藏，應合理改進。該結果與實際馬鈴薯貯藏室內溫度分布較為符合，能夠在一定程度上真實反映靜態條件下的馬鈴薯貯藏熱量傳遞情況。

2.5 不同二氧化碳含量對馬鈴薯貯藏過程中傳熱的影響

為研究不同組分的二氧化碳氣體對馬鈴薯貯藏過程中傳熱的影響情況，分別截取二氧化碳含量為0.00%、0.15%、0.30%時，貯藏室在X=4 m、Z=2 m處的溫度分布云圖截面進行對比。從圖9的溫度云圖可以得到貯藏室內的整體溫度分布狀態。云圖中紅色中心區域是馬鈴薯堆垛的存放區，可以看到這個位置的溫度是由馬鈴薯堆垛的中心區域向外部逐漸發散，表明熱源的加入會以輻射的方式對馬鈴薯周邊溫度產生一定的影響。當二氧化碳含量由0.00%提高到0.30%時，馬鈴薯堆體表面的空氣溫度由1.23 ℃提升至1.37 ℃，因受多孔介質模型的影響，堆體內部溫度呈近球形分布，由1.85 ℃逐層遞減至室溫。貯藏室內低于1.05 ℃的溫度區域，隨著二氧化碳含量的提升越來越小，而室內頂部接近1.35 ℃的區域卻越來越大。

圖9 貯藏室在X=4 m處不同CO2氣體含量的溫度分布Fig. 9 Temperature distribution of different CO2 gas content in the storage room at X=4 m

從圖10 可以看出，貯藏室內溫度整體分布由屋頂到地面遞減，北側墻溫度高于南側墻，發熱的馬鈴薯堆對底層空氣有加熱作用，二氧化碳含量的提升擴大了高溫區域的范圍。由二氧化碳發生器輸入的較高溫度（3.85 ℃）的二氧化碳氣體在北墻附近形成了1 個溫度較高的回流區，隨著二氧化碳含量的提升，回流區的中心溫度由0.75 ℃提升至1.28 ℃。由于發生器噴出氣體速度較快，在噴口及斜上方形成了高溫區域，該區域的整體溫度也隨著二氧化碳含量的提升而變大。二氧化碳的噴出對發熱的馬鈴薯堆體溫度分布也造成了影響，接近發生器一側的馬鈴薯堆體附近的整體溫度明顯高于另外一側，這會造成堆體溫度分布不均，影響馬鈴薯的品質。

圖10 貯藏室在Z=2 m處不同CO2氣體含量的溫度分布Fig. 10 Temperature distribution of different CO2 gas content in the storage room at Z=2 m

2.6 貯藏室中心溫度點分析

當二氧化碳氣體的含量分別為0.00%、0.15%、0.30%時，半地下式馬鈴薯貯藏室中心位置的橫向切面溫度如圖11 所示，含有0.30%二氧化碳氣體的半地下式馬鈴薯貯藏室整體各點溫度均高于其他2 個處理的半地下式馬鈴薯貯藏室整體各點溫度。其中，在貯藏室橫切面取點中心位置，3 個二氧化碳組分的貯藏室點溫度均達到最大值，是因為在馬鈴薯堆垛中心加入了熱源，表明熱源對于馬鈴薯貯藏室溫度有一定的影響，可以輻射的方式使周圍溫度升高。

圖11 貯藏室中心橫切面溫度分布Fig. 11 Temperature distribution in the cross-section of the center of the storage room

不同二氧化碳含量下，半地下式馬鈴薯貯藏室中心位置的豎向切面溫度如圖12 所示，在豎切面取點坐標的初始位置，0.15%和0.30%二氧化碳氣體含量所對應的溫度均處于較高狀態，因為此處距離二氧化碳發生器的出氣口最近，而二氧化碳是溫室氣體，且二氧化碳發生器在生成二氧化碳時會產生一定的熱量，從而導致此處溫度略高。在豎切面取點坐標的中心位置，3 組溫度點線圖均達到了最大值，同樣是因為熱源的影響。而在取點坐標末尾處0.15%和0.30%二氧化碳含量所對應的貯藏室溫度取點數據均達到了最低值，一方面是因為此處距離二氧化碳發生器位置較遠，另一方面因為處于中間位置的馬鈴薯堆垛阻擋住了二氧化碳氣體的擴散，致使此處溫度相較于貯藏室整體溫度略偏低。

圖12 貯藏室中心豎切面溫度分布Fig. 12 Temperature distribution of the vertical section in the center of the storage room

3 討論

為了明確半地下式貯藏室貯藏馬鈴薯的過程中二氧化碳含量對室內熱量傳遞的影響，本研究基于實際場地獲得的環境因子數據，利用CFD 數值技術對密閉狀態下貯藏室內二氧化碳含量對室內溫度分布的影響展開研究，通過對數值方法、網格無關性及與實測值的對比，多角度驗證了數值方法在半地下式貯藏室內開展相關研究的有效性，發現在密閉環境下二氧化碳的含量提升，明顯地提高了貯藏室內溫度，長期積累對馬鈴薯的貯藏產生不利影響。周博等[11,30]利用近似的數值模型研究了在氣調庫的快速降氧環節中，庫內氧氣和二氧化碳含量隨時間變化的規律以及氮氣純度對降氧時間的影響，與試驗數據有較好的一致性，確保了模擬方法的精確性，證明了CFD 數值方法在半地下式貯藏室內的密閉環境中研究氣體含量與環境因子的相互作用具有可行性。受馬鈴薯呼吸作用的影響，二氧化碳含量的逐漸升高使貯藏室內溫度逐步提升，馬鈴薯在高溫及缺氧狀況下造成果體二氧化碳積累，無氧呼吸作用增強，導致馬鈴薯塊莖活性降低，內部出現黑心病變，影響馬鈴薯的貯藏品質。

因受條件限制，對照試驗中使用的馬鈴薯堆體體積較小且測點布置數量有限，本研究采用穩態模擬馬鈴薯堆體的發熱以及貯藏室內的環境，若考慮真實情況，馬鈴薯的散熱和室內各項環境因子的變化更趨向于非穩態，如能將時間變化考慮到室內環境的模擬中，能夠獲得更符合實際的室內溫度分布規律。