基于CFD的蘋果隔板包裝預冷溫度場研究

朱文穎 史 策 韓 帥 劉 歡 楊信廷

(1.國家農業(yè)信息化工程技術研究中心， 北京 100097； 2.農產品質量安全追溯技術及應用國家工程實驗室， 北京 100097)

0 引言

目前，果品產業(yè)中最大的關注點是采收后水果在供應鏈中的質量損失，據(jù)統(tǒng)計從水果收獲后直至到達消費者手中，損失率高達13%～38%[1]，減少水果的采后損失是確保未來水果產業(yè)食品安全的重要組成部分[2]。溫度調控失誤和物理損傷是導致水果質量損失的重要因素，隔板包裝(包裝箱內部含有隔板)能夠防止水果的擠壓和碰撞，研究隔板包裝的內部溫度場分布具有重要意義[3-4]。

圖1 兩種包裝箱箱體的結構示意圖Fig.1 Structure diagram of two kinds of packaging1.通風箱體 2.開孔 3.隔板 4.蘋果

計算流體力學(Computer fluent dynamic，CFD)數(shù)值模擬方法在研究預冷過程中的溫度場分布發(fā)揮著重要作用。趙春江等[5-6]總結歸納了CFD技術在處理預冷問題時前處理、求解過程和后處理3個階段的具體方法，闡述了近10年中CFD技術在冷鏈物流發(fā)展中的應用，展望了冷鏈物流未來的發(fā)展趨勢。HAN等[7]針對內含兩層蘋果的包裝箱開孔方式進行了模擬分析，比較了包裝箱內部氣流溫度場分布云圖，發(fā)現(xiàn)增加包裝箱開氣孔數(shù)量能夠增強預冷速率，減少預冷時間。文獻[8-10]模擬不同大小橘子在單個包裝盒內部的氣流形式，為水果不規(guī)則排布的溫度分布提供了參考。FERRUA等[11-12]模擬裝有不規(guī)則形狀、不同尺寸草莓的塑料盒在不同開孔方式下盒中溫度分布情況，為不規(guī)則形狀溫度場分布提供了參考。HAN等[13]收集整理了10種常見不同尺寸和開孔方式的包裝箱，得到最優(yōu)預冷包裝箱。DELELE等[14]模擬9種相同尺寸、不同開孔方式的包裝箱在裝載和空載情況下的預冷狀況，對比分析得到開孔數(shù)目、開孔形狀和開孔位置對包裝箱內部氣流溫度的影響。

然而，包裝箱內部結構對溫度場和氣流場分布也具有重要意義[15]。相比普通包裝箱(內部不含隔板的包裝箱)，隔板包裝箱具有較好抗震功能，能夠防止水果的擠壓和碰撞[16-17]，減少水果質量損失。而正確評價隔板包裝箱的預冷效果在預冷研究中同樣至關重要。

為此，本文以蘋果為研究材料，分別建立隔板包裝模型和普通包裝模型，對滿載相同數(shù)量蘋果的兩種包裝箱分別進行CFD模擬計算，得到兩種包裝箱內部溫度氣流分布，對比分析兩種包裝箱內部蘋果的冷卻時間、冷卻速率以及冷卻均勻性，客觀評價隔板包裝箱和普通包裝箱的預冷效果。

1 材料與方法

1.1 物理模型

蘋果在打包過程中，包裝箱設計應兼顧防止蘋果機械損傷和保證冷空氣流通均勻性的作用。文獻[18-19]表明，改變包裝箱通風孔尺寸、形狀位置等因素，水果的冷卻時間、冷卻速率以及冷卻均勻性都會發(fā)生改變。本文以裝有蘋果的包裝箱為研究對象，分為單層隔板包裝和單層普通包裝，采用Icem軟件進行幾何建模，兩種包裝箱結構如圖1所示。包裝箱材料采用市場上常用B楞的瓦楞紙，尺寸(長×寬×高)為55 cm×40 cm×14 cm，箱體厚度為0.5 cm，隔板包裝箱的隔板厚度為0.3 cm，包裝箱內分為12個隔間，每個隔間內裝有一個蘋果，蘋果直徑為12 cm，兩種包裝箱內蘋果擺放位置相同。兩種包裝開孔方式相同，進風口和出風口均為圓形，直徑4 cm，單面開孔率都為3.3%，具體位置如圖1所示。

1.2 網格劃分

本文采用Icem進行網格劃分，為提高計算精度，保證計算的收斂性，整個箱體及內部結構采用六面體結構化網格進行劃分，蘋果與壁面、蘋果與隔板之間留有一定空隙，蘋果、進風口以及出風口等球體或圓形區(qū)域采用“O”型網格技術生成六面體結構化網絡，并進行局部加密。隔板包裝箱最大網格尺寸不超過3 mm，網格數(shù)量為646 443個單元，普通包裝箱最大網格尺寸不超過5 mm，網格數(shù)量為617 422個單元，劃分結果如圖2所示。由于本文采用六面體結構化網格對模型進行劃分，利用正交性和寬高比對網格進行檢測[20]，正交性指標分布從0到1，趨于0則網格質量較差。結果表明2種包裝網格的正交性和寬高比良好，正交性指標均大于0.6，寬高比均控制在1∶1～5∶1之間。

圖2 網格劃分Fig.2 Mesh generation

1.3 數(shù)學模型

建立數(shù)學模型求解兩種包裝箱內部預冷過程的溫度分布情況時，為了對實際模型進行簡化，減少計算時間，對模型做如下假設：兩種包裝箱密封性良好，不存在漏氣現(xiàn)象；空氣為不可壓縮氣體；蘋果、包裝箱和隔板的物理參數(shù)不隨溫度的變化而改變；蘋果是各項同性連續(xù)介質等直徑球體，忽略蘋果產生的呼吸熱和蒸騰熱，忽略蘋果之間的輻射作用；蘋果與紙箱之間采用縮小法建模[21]，縮小系數(shù)為0.93。

模型的材料主要包括空氣、外包裝箱體、隔板以及蘋果，依據(jù)文獻[22-23]，空氣、蘋果、外包裝箱體及隔板的具體熱物理性能參數(shù)如表1所示。

表1 熱物性能參數(shù)Tab.1 Thermo-physical properties

1.4 初始條件與邊界條件

(1)初始條件：當時間t=0時，外包裝箱體、隔板、包裝箱內部空氣以及蘋果的初始溫度T=25℃。

(2)入口邊界：將包裝箱迎風口面的兩個開孔設置為速度入口邊界條件，文獻[24]表明，當送風速度超過2.5 m/s時，制冷時間、冷卻速度以及冷卻均勻性無顯著改變，因此本文入口邊界設置為：風速2.5 m/s、空氣溫度0℃、湍流動能0.04 m2/s2、耗散率12.5 s-1。

(3)出口邊界：將包裝箱背風口面的兩個開孔設置為出口邊界條件，邊界上所有物理量梯度為零。

(4)壁面邊界：蘋果表面、紙箱表面和隔板表面設置為無滑移壁面條件，壁面上速度為零，且垂直于壁面的速度也為零。

1.5 數(shù)值模擬方法

采用基于有限體積法的CFD商用軟件Fluent 15.0進行求解計算和后處理分析，考慮重力的影響，利用非穩(wěn)態(tài)剪切壓力傳輸(Shear stress transport, SST)k-ω湍流模型、基于壓力的分離式求解器進行計算，在空氣流體區(qū)采用的控制方程為質量守恒方程、動量守恒方程以及能量守恒方程，將熱源加載到蘋果區(qū)域控制方程中，采用SIMPLE(Semi-implicit method for pressure-linked equations)算法求解，時間步長設置為6 min，步數(shù)設為100步，共計10 h。

1.6 評估參數(shù)

采用預冷模擬中常用的3個參數(shù)進行模擬評估，分別是冷卻時間、瞬時冷卻速率和瞬時冷卻均勻性。

無量綱數(shù)θ用來評估達到冷卻溫度的農產品占整體農產品的比值，計算公式為

θ=(Tin-T0)/(Tave-T0)

(1)

式中Tin——包裝箱內蘋果初始平均溫度，℃

Tave——所有蘋果平均溫度，℃

T0——進風口預冷溫度，即預期預冷溫度，℃

當θ=7/8時，果品預冷基本達到預期溫度，所用時間稱為7/8冷卻時間(Seven-eight cooling time, SECT)；當θ=1/2時，冷卻速率達到平均水平，所用時間稱為1/2冷卻時間(Half cooling time, HCT)。

瞬時冷卻速率用來評估蘋果在冷卻過程中的降溫快慢，計算公式為

R=(Tt-i-Tt)/i

(2)

式中Tt——t時刻蘋果溫度,℃

Tt-i——t-i時刻蘋果溫度,℃

i——時間步長,min

通常采用溫度變異系數(shù)來評估冷卻均勻性，溫度變異系數(shù)越大，箱體內部蘋果溫度數(shù)據(jù)離散程度越大，溫度分布越不均勻；反之，溫度變異系數(shù)越小，箱體內部蘋果溫度數(shù)據(jù)離散程度越小，溫度分布越均勻。溫度變異系數(shù)計算公式[25]為

(3)

式中n——包裝箱內蘋果總數(shù)量

Tk——蘋果k的溫度

溫度變異系數(shù)通常采用熱力學溫度進行計算。

圖4 相同邊界條件兩種包裝不同位置蘋果內部平均溫度曲線Fig.4 Changing curves of average temperature of apples at different positions with the same boundary condition

2 結果與分析

2.1 冷卻時間

圖3為兩種包裝在相同初始條件和相同邊界條件下，通過CFD模擬得到的蘋果內部平均溫度變化曲線和蘋果表面平均溫度變化曲線。虛線與曲線的交點橫坐標代表7/8冷卻時間。

圖3 相同邊界條件蘋果平均溫度變化曲線Fig.3 Changing curves of average temperature of all apples with the same boundary condition

由圖3可知，隔板包裝和普通包裝的溫度曲線趨勢基本相同，隨著冷卻時間增長，冷卻速率逐漸減慢。兩種包裝蘋果表面平均溫度差值不超過0.67℃，兩種包裝蘋果內部平均溫度差值不超過0.26℃。蘋果表面平均溫度低于蘋果內部平均溫度，隨著冷卻時間增長，蘋果內部受熱傳導作用的影響，表面平均溫度曲線與內部平均溫度曲線逐漸重合。

但是，隔板包裝和普通包裝在箱體內部相同位置的蘋果溫度差別較大，相同包裝箱內部不同位置蘋果的溫度也有所差別，按照位置的不同，將包裝箱內部的蘋果分為3部分：進風口位置蘋果、中間位置蘋果和出風口位置蘋果，內部平均溫度曲線如圖4所示。

如圖4c所示，隔板包裝和普通包裝中間位置蘋果溫度基本相同；而隔板包裝進風口位置蘋果溫度低于普通包裝進風口位置蘋果溫度，如圖4b所示，存在該現(xiàn)象的原因歸結于內部隔板對冷空氣的傳播起到阻礙作用，冷空氣從進風口進入包裝箱后，受到隔板阻礙，在進風口隔間內停留，形成局部低溫區(qū)域，無法傳遞到出風口處，致使隔板包裝在出風口隔間位置相對高溫，如圖4d所示，中間區(qū)域次之。

2.2 冷卻速率

冷卻過程中，兩種包裝在相同初始條件和邊界條件下，蘋果瞬時冷卻速率隨時間的變化曲線如圖5所示。

圖5 相同邊界條件蘋果瞬時冷卻速率曲線Fig.5 Curves of instantaneous cooling rates of all apples with the same boundary condition

由圖5可知，在整個預冷過程中，隔板包裝和普通包裝內部蘋果冷卻速率曲線整體趨勢基本相同，隨著冷卻時間增加冷卻速率逐漸減慢，兩種包裝蘋果表面冷卻曲線和蘋果內部冷卻曲線基本重合。冷卻時間30 min內，蘋果表面冷卻速率遠大于蘋果內部冷卻速率；冷卻時間超過30 min后，冷卻速率下降到0.15℃/min,蘋果表面冷卻速率逐漸慢于蘋果內部冷卻速率；冷卻時間超過300 min后，蘋果表面冷卻速率和蘋果內部冷卻速率基本重合。整個過程中，兩種包裝瞬時冷卻速率差值不超過0.012℃/min。

同理比較包裝箱內部3個位置的蘋果溫度，隔板包裝進風口位置蘋果冷卻速率快于普通包裝進風口位置冷卻速率，隔板包裝出風口位置蘋果冷卻速率慢于普通包裝出風口位置冷卻速率，兩種包裝中間位置蘋果冷卻速率基本相同。

圖7 冷卻過程中不同時間下隔板包裝內的溫度分布云圖Fig.7 Instantaneous static temperature contours of partition packaging with the same boundary condition

比較隔板包裝和普通包裝的冷卻時間和冷卻速率，兩種包裝的整體預冷效果基本一致，均能滿足冷鏈運輸過程中的預冷需求。

2.3 冷卻均勻性

隔板包裝和普通包裝溫度變異系數(shù)隨時間變化的曲線如圖6所示。

圖6 相同邊界條件下溫度變異系數(shù)曲線Fig.6 Curves of instantaneous temperature variation coefficients with the same boundary condition

由圖6可知，隔板包裝和普通包裝的溫度變異系數(shù)均為先增大后減小，隨著冷卻時間增長，包裝箱內蘋果之間溫度差值逐漸減小，溫度變異系數(shù)逐漸減小。對比兩種包裝的冷卻過程，普通包裝的溫度變異系數(shù)小于隔板包裝的溫度變異系數(shù)，普通包裝蘋果冷卻均勻性優(yōu)于隔板包裝蘋果冷卻均勻性。存在該現(xiàn)象的原因歸結于隔板包裝內部冷空氣受到隔板阻礙，隔斷之間冷空氣無法直接熱傳遞，形成溫度差異，隨著冷卻時間增長，兩種包裝內部溫度逐漸趨于零，包裝內部溫度場逐漸均勻，溫度變異系數(shù)逐漸趨于零。

2.4 云圖分析

圖7和圖8給出不同冷卻時間隔板包裝和普通包裝的內部溫度場分布云圖，由圖7和圖8可以明顯看出，在相同冷卻時間，隔板包裝與普通包裝溫度場分布基本相同，冷空氣進入包裝箱后，進風口位置蘋果首先預冷降溫，進風口位置蘋果溫度小于出風口位置蘋果溫度。冷卻3 h后，兩種包裝箱內部所有蘋果的溫度均降至10℃以下；冷卻5 h后，兩種包裝箱內部所有蘋果的溫度降至4℃以下，滿足7/8冷卻時間的冷卻條件；冷卻10 h后，兩種包裝箱內空氣溫度和蘋果溫度基本保持在0℃。

圖8 冷卻過程中不同時間下普通包裝內的溫度分布云圖Fig.8 Instantaneous static temperature contours of ordinary packaging with the same boundary condition

普通包裝內溫度分布云圖如圖8所示，冷空氣從進風口均勻流入，出風口均勻流出，內部溫度分布平緩，沒有隔板阻礙現(xiàn)象和溫度突變現(xiàn)象。相比普通包裝，隔板包裝箱的內部隔板對冷空氣具有明顯阻礙現(xiàn)象，隔板之間存在微小溫度差，如圖7所示。

圖9和圖10給出冷卻時間為120 min隔板包裝和普通包裝內部的風速分布圖。由圖9和圖10能夠明顯看出，兩種包裝內部風速范圍基本相同，均保持在1.0 m/s 左右。普通包裝內部風速分布較為平穩(wěn)均勻，如圖10所示。相比普通包裝，隔板包裝內部隔板對風速傳播具有阻礙作用，隔板處風速約為2.5 m/s，高于箱體內部其他位置風速，如圖9所示。

圖9 冷卻120 min隔板包裝的風速分布圖Fig.9 Velocity distribution of partition packaging at 120 min

圖10 冷卻120 min普通包裝的風速分布圖Fig.10 Velocity distribution of ordinary packaging at 120 min

3 實驗驗證

3.1 驗證方法

為驗證模擬結果的正確性，以裝滿蘋果的隔板包裝箱為研究對象，采用天津商業(yè)大學的壓差預冷實驗臺進行實驗驗證，實驗裝置布置圖如圖11所示。圖中隔板包裝紙箱中裝有12個直徑為8 cm的蘋果，分別編號為1～12，分為進風口位置蘋果(編號為1～4)、中間位置蘋果(編號為5～8)和出風口位置蘋果(編號為9～12)3組，每個蘋果如圖11所示插入2個溫度傳感器，實時監(jiān)測蘋果內部溫度，圖11中“·”為溫度測試點。實驗前，首先確保紙箱內部蘋果平均溫度為25℃，通入風速為2.5 m/s的冷空氣，隔板紙箱內蘋果開始降溫，采用日本Yokogawa溫度記錄儀進行測試記錄，溫度傳感器靈敏度為0.5℃，每隔2 min記錄一次數(shù)據(jù)，記錄時間為3 h，共90個數(shù)據(jù)。

圖11 實驗平臺結構與實驗裝置布署圖Fig.11 Structure diagram of experiment platform and equipment deployment1.壓差風機 2.蒸發(fā)器 3.靜壓室 4.蜂窩網 5.蓋子 6.隔板紙箱 7.插有溫度傳感器的蘋果 8.支架 9.加濕器 10.溫度記錄儀

3.2 結果與分析

讀取實驗中溫度記錄儀的數(shù)據(jù)并進行數(shù)據(jù)處理，計算平均值，風速為2.5 m/s時，3 h內蘋果的溫度模擬值和溫度實驗值曲線如圖12所示。

圖12 溫度模擬值和實驗值比較Fig.12 Comparison of experiment and simulation results of temperature

由圖12可知，3 h內隔板包裝箱內全部蘋果平均溫度模擬值與平均溫度實驗值變化趨勢相同，對比所有溫度模擬值與實驗值得出，溫度最大差值小于1.5℃。模擬值與實驗值存在偏差的可能原因主要歸結于蘋果形狀的簡化、測量儀器自身誤差、外界環(huán)境的簡化以及壓差預冷實驗臺風速溫度的誤差。

風速為2.5 m/s時，3 h內隔板包裝箱內不同位置蘋果的溫度實驗值如圖13所示。

圖13 風速2.5 m/s不同位置蘋果溫度實驗曲線Fig.13 Experimental temperature at different positions with wind speed of 2.5 m/s

由圖13可知，隔板包裝內不同位置蘋果溫度實驗結果和溫度模擬結果趨勢基本相同，區(qū)別在于出風口位置溫度實驗值低于出風口位置溫度模擬值。中間位置和出風口位置蘋果溫度實驗值曲線重合。

采用均方根誤差(Root mean square error,RMSE)和平均相對誤差(Average relative deviation, ARD)判定實驗值與模擬值之間的關系，計算結果見表2。

表2 均方根誤差與平均相對誤差Tab.2 Values of RMSE and ARD

由表2可知，進風口位置模擬值和實驗值偏差小，出風口位置模擬值和實驗值偏差大，存在偏差的原因可歸結于模擬值忽略包裝箱體外界環(huán)境，實驗過程中隔板包裝箱外部溫度逐漸降低，出風口位置蘋果在開孔處與外界冷空氣相接觸，溫度模擬值比溫度實驗值大。

4 結論

(1)隔板包裝和普通包裝內蘋果的冷卻時間和冷卻速率曲線基本一致，冷卻溫度最大差值不超過0.67℃，冷卻速率最大差值不超過0.012℃/min。

(2)冷卻過程中，隔板包裝溫度變異系數(shù)略大于普通包裝溫度變異系數(shù)，隔板包裝進風口位置蘋果溫度低于普通包裝進風口位置蘋果溫度，隔板包裝出風口位置蘋果溫度高于普通包裝出風口位置蘋果溫度，但不影響冷卻效果。

(3)冷卻過程中隔板包裝和普通包裝在相同時間內部的溫度分布范圍基本一致，120 min時隔板包裝和普通包裝內的風速分布趨勢基本相同，基本保持在內部風速約為1.0 m/s。

(4)經過實驗驗證，實驗結果與模擬結果基本一致，溫度最大差值小于1.5℃，全部蘋果的均方根誤差為1.13℃，平均相對誤差為8.2%，實驗驗證了本研究利用CFD模擬的合理性與準確性。