數(shù)值模擬技術(shù)預(yù)測(cè)風(fēng)機(jī)兩種擺放方式對(duì)冷庫(kù)堆垛貨物的影響

杜子崢 謝 晶 朱進(jìn)林

（1.上海海洋大學(xué)食品學(xué)院，上海 201306；2.上海水產(chǎn)品加工與貯藏工程技術(shù)研究中心，上海 201306）

果蔬等農(nóng)產(chǎn)品在進(jìn)入超市，供給消費(fèi)者前，需在冷庫(kù)中進(jìn)行較長(zhǎng)時(shí)間貯藏，冷庫(kù)內(nèi)部合理的溫度分布可顯著延長(zhǎng)商品貨架期，減少商品損耗，提高經(jīng)濟(jì)效益。為研究冷庫(kù)內(nèi)氣流的分布機(jī)理以改善貨物貯藏條件，國(guó)內(nèi)外研究人員［1－6］對(duì)冷庫(kù)內(nèi)氣流組織分布做了大量的研究，主要通過(guò)計(jì)算流體力學(xué)模擬和實(shí)驗(yàn)驗(yàn)證兩種方式相互驗(yàn)證。初期的研究主要通過(guò)對(duì)冷庫(kù)空庫(kù)模型進(jìn)行數(shù)值模擬，使學(xué)者對(duì)冷庫(kù)內(nèi)部氣流組織的溫度、風(fēng)速、渦旋分布有了初步認(rèn)識(shí)。然而冷庫(kù)主要功能是貯藏貨物，貨物不同的擺放方式、堆垛尺寸、潛在的蒸發(fā)潛熱、呼吸熱均會(huì)顯著影響其內(nèi)部氣流組織，因此僅對(duì)空庫(kù)進(jìn)行模擬研究遠(yuǎn)遠(yuǎn)不夠，必須考慮貨物對(duì)冷庫(kù)氣流組織的影響。劉妍玲等［7］在果蔬擺放方式研究中發(fā)現(xiàn)中間有通道式果蔬擺放方式，氣流組織分布較合理。Chourasia等［8］利用數(shù)值模擬技術(shù)，圍繞冷庫(kù)土豆堆垛寬高比、堆垛間距做了大量的研究，發(fā)現(xiàn)加大貨物垂直間距可顯著改善貨物溫度分布。Delele等［9］在對(duì)側(cè)吹風(fēng)式水果冷庫(kù)加濕系統(tǒng)CFD模擬研究中，建立了單個(gè)堆垛貨架模型，并對(duì)貨架采用多孔介質(zhì)模型模擬，成功預(yù)測(cè)了冷庫(kù)內(nèi)部空氣相對(duì)濕度、速度以及貨架溫度分布，發(fā)現(xiàn)貨物頂部溫度最高，貨物最熱區(qū)域與冷庫(kù)空氣側(cè)存在0.92℃溫差。Tanaka等［10］對(duì)多種貨物裝載方式進(jìn)行數(shù)值模擬研究，結(jié)合試驗(yàn)測(cè)量值，確定了貨物最優(yōu)化堆放方式。但是以上研究并沒(méi)有在冷庫(kù)滿(mǎn)載情況下進(jìn)行。

本研究擬以筆者實(shí)際參觀(guān)走訪(fǎng)的一家中型冷庫(kù)為藍(lán)本，建立模型，利用數(shù)值模擬技術(shù)，在貨物滿(mǎn)載下研究對(duì)吹風(fēng)式和下吹風(fēng)式兩種風(fēng)機(jī)擺放方式對(duì)庫(kù)內(nèi)氣流分布以及堆垛貨物溫度分布的影響，旨在尋找適合中型冷庫(kù)的最佳風(fēng)機(jī)設(shè)置方案，為未來(lái)冷庫(kù)建設(shè)提供一定的參考依據(jù)。

1 模型和方法

1.1 物理模型建立

研究對(duì)象為一間24.3m×21.6m×7.2m的中型冷庫(kù)，庫(kù)內(nèi)貯藏果蔬等農(nóng)產(chǎn)品，配備有貨架以方便貨物堆垛，冷庫(kù)實(shí)際安裝4臺(tái)吊頂式的Kuba GEA冷風(fēng)機(jī)，采用側(cè)面對(duì)吹方式，風(fēng)機(jī)總功率為8.96kW，每臺(tái)風(fēng)機(jī)有3個(gè)直徑為0.45m的圓形風(fēng)口，出風(fēng)溫度為0℃，風(fēng)速12m/s。庫(kù)內(nèi)共有42個(gè)貨架，貨架間距0.5m，離壁面距離1m，貨架底部距離地面0.3m，貨物堆垛尺寸為4.5m×1.2m×5.5m。冷庫(kù)中間設(shè)有寬度為2m走道，以方便人員進(jìn)行堆垛和捯垛操作。為研究風(fēng)機(jī)不同擺放方式對(duì)庫(kù)內(nèi)貨物影響，分別建立兩種模型：① 吊頂風(fēng)機(jī)對(duì)送風(fēng)模型，該模型以實(shí)際冷庫(kù)作為模板，風(fēng)機(jī)設(shè)置在兩側(cè)墻面，結(jié)構(gòu)見(jiàn)圖1（a）；② 風(fēng)機(jī)上回下送風(fēng)模型，風(fēng)機(jī)設(shè)置在冷庫(kù)頂部，且出風(fēng)口下側(cè)處于貨架間隙，模型結(jié)構(gòu)見(jiàn)圖1（b）。由于兩種模型均具有對(duì)稱(chēng)結(jié)構(gòu)，不影響計(jì)算結(jié)果的前提下，為縮短計(jì)算時(shí)間，只對(duì)模型取對(duì)稱(chēng)面一側(cè)進(jìn)行模擬。

圖1 兩種不同風(fēng)機(jī)擺放冷庫(kù)模型Figure 1 Two different cold store models with various fan arrangement

為簡(jiǎn)化計(jì)算，對(duì)模型進(jìn)行如下假設(shè)：

（1）冷庫(kù)壁面假設(shè)為無(wú)滑移恒溫壁面，壁面溫度恒定為3℃，冷庫(kù)內(nèi)部密封良好，不考慮外界滲透作用的影響；

（2）冷庫(kù)氣體為不可壓縮的理想氣體，符合Boussinesq假設(shè)；

（3）貨物側(cè)存在呼吸熱，且呼吸熱為恒定值；

（4）貨物側(cè)為多孔結(jié)構(gòu)，且內(nèi)部氣體是層流流動(dòng)；

（5）不考慮內(nèi)部因潛熱導(dǎo)致的溫度變化。

1.2 堆垛貨物模型方程

冷庫(kù)內(nèi)部堆放貨物被視為多孔介質(zhì)，多孔介質(zhì)模型的動(dòng)量方程是在標(biāo)準(zhǔn)動(dòng)量方程的后面加上動(dòng)量源項(xiàng)。貨物側(cè)流動(dòng)阻力可以參照Darcy-Forchheimer公式：

其中右側(cè)第一項(xiàng)為考慮粘性損失的達(dá)西公式項(xiàng)，右側(cè)第二項(xiàng)為慣性損失項(xiàng)。Darcy-Forchheimer公式通常也可以表示為：

式中：

p——壓強(qiáng)，Pa；

u——速度矢量，m/s；

K——通過(guò)多孔區(qū)域的達(dá)西滲透率，D；

C——慣性阻力系數(shù)。

當(dāng)局部雷諾數(shù)大于1時(shí)，必須考慮慣性損失項(xiàng)。局部雷諾數(shù)公式：

式中：

u——速度矢量，m/s；

ρ——流體密度，kg/m3；

dpe——微粒的有效當(dāng)量直徑，m。

當(dāng)量直徑可由式（4）進(jìn)行計(jì)算：

式中：

V——顆粒的體積，m3。

為簡(jiǎn)化模型，將貨物看做均質(zhì)的體積相同的個(gè)體。在多孔介質(zhì)模型的使用上，模擬多孔板或者管束系統(tǒng)，可以忽略滲透項(xiàng)，只使用慣性損失項(xiàng)。但是對(duì)堆垛貨物忽略滲透項(xiàng)會(huì)導(dǎo)致模擬精度下降，無(wú)法實(shí)現(xiàn)貨物內(nèi)部氣流的準(zhǔn)確預(yù)測(cè)，需考慮堆垛貨物內(nèi)部滲透項(xiàng)系數(shù)，使用ergun方程［11］的填充床模型對(duì)貨物側(cè)進(jìn)行計(jì)算，ergun方程表示為：

式中：

L——填充床厚度，m。

滲透阻力系數(shù)和慣性阻力系數(shù)可表示為：

式中：

ε——多孔介質(zhì)區(qū)域的孔隙率，%；

λ——形狀系數(shù)。

孔隙率及形狀系數(shù)公式：

式中：

A——顆粒表面積，m2；

ρb——貨物側(cè)區(qū)域密度，kg/m3；

ρp——貨物密度，kg/m3。

1.3 湍流模型選擇

工程應(yīng)用中湍流的數(shù)值模擬主要分為三類(lèi)：直接數(shù)值模擬（DNS）、大渦模擬（LES）和基于雷諾平均N—S方程組（RANS）模型。由于計(jì)算機(jī)條件約束，前兩種模擬方式仍在探索階段，RANS模型在工程應(yīng)用中使用最廣泛。冷庫(kù)湍流模型主要使用k—ε湍流模型、SSTk—ω湍流模型以及RSM湍流模型。早期研究一般采用工程通用的k—ε模型，但其預(yù)測(cè)精度仍有待檢驗(yàn)，Delete等［12］對(duì)多種兩方程模型（standardk—ε、RNGk—ε、realizablek—ε、standardk—ε和SSTk—ω湍流模型）預(yù)測(cè)精度進(jìn)行檢測(cè)，上述模型對(duì)平均風(fēng)速預(yù)測(cè)的相對(duì)誤差分別為24.3%，22.4%，23.5%，18.2%，SSTk—ω湍流模型預(yù)測(cè)較準(zhǔn)確，因此，本次模擬采用SSTk—ω湍流模型。

1.4 多孔介質(zhì)區(qū)域參數(shù)設(shè)定

將貨物側(cè)視為多孔介質(zhì)區(qū)域，考慮貨物內(nèi)部滲透性較差以及其幾何尺度比湍流渦的尺度要大，因此必須抑制多孔介質(zhì)區(qū)域湍流影響，使用laminar zone模型對(duì)湍流生成進(jìn)行抑制，多孔介質(zhì)區(qū)域視為各項(xiàng)同性。貨物側(cè)堆垛為有效直徑80 mm的蘋(píng)果，其呼吸熱可依據(jù)公式Sef＝1.91Tf＋4.64（Anno）［13］計(jì)算。蘋(píng)果物性參數(shù)［14］：密度為837kg/m3，熱導(dǎo)率為0.558W/（m·℃），比熱為3 665.6J/（kg3·℃）。蘋(píng)果表面積A可根據(jù)公式A＝8.849－0.835 9dpe＋0.032 5［15］進(jìn)行計(jì)算，貨物側(cè)滲透系數(shù)、慣性系數(shù)以及其他參數(shù)見(jiàn)表1。

表1 多孔介質(zhì)區(qū)域具體參數(shù)設(shè)定Table 1 Specific parameters of the porous medium region

1.5 邊界條件設(shè)定及計(jì)算方法

使用非結(jié)構(gòu)網(wǎng)格對(duì)模型進(jìn)行網(wǎng)格劃分，計(jì)算域網(wǎng)格尺寸為1.4cm，下吹風(fēng)式冷庫(kù)生成3.3×106個(gè)網(wǎng)格，兩種模型分別生成3.25×106和3.26×106個(gè)網(wǎng)格。風(fēng)機(jī)出口采用Velocity inlet速度出口，回風(fēng)邊界設(shè)為outflow，紊流強(qiáng)度和水力直徑定義湍流，湍流強(qiáng)度根據(jù)公式I＝0.16Re－1/8求得為5%，水利直徑為出風(fēng)口特征尺寸0.45m。所有固體表面設(shè)為無(wú)滑移壁面，冷庫(kù)外墻和地面為現(xiàn)場(chǎng)聚氨酯發(fā)泡，其熱導(dǎo)率為0.022W/（m·K），墻體厚度0.15m，外墻和地面采用定熱流密度邊界條件，熱流密度可根據(jù)式（10）進(jìn)行計(jì)算。

式中：

λw——墻體熱導(dǎo)率，W/（m·℃）；

tw——墻體厚度，m；

Tw、T——分別為室外溫度（10℃）和冷庫(kù)內(nèi)部溫度（0℃），℃。

經(jīng)計(jì)算其熱流密度為1.47W/m2。

考慮因浮力驅(qū)動(dòng)的自然對(duì)流的影響，引入Boussinesq假設(shè)。利用fluent對(duì)兩種模型進(jìn)行穩(wěn)態(tài)模擬以獲得冷庫(kù)內(nèi)部氣流相對(duì)穩(wěn)定時(shí)的氣流分布，連續(xù)性殘差設(shè)定為10－3，其余各值設(shè)為10－4。

2 穩(wěn)態(tài)模擬結(jié)果分析

通過(guò)fluent計(jì)算之后，得到兩種風(fēng)機(jī)擺放方式下冷庫(kù)內(nèi)部速度及溫度數(shù)據(jù)，分別對(duì)它們主流截面、貨架底部截面及貨物側(cè)進(jìn)行速度及溫度的對(duì)比分析。

2.1 速度分布及對(duì)比

由于風(fēng)機(jī)設(shè)置形式不同，出風(fēng)主流截面風(fēng)速分布差異顯著（見(jiàn)圖2）。下吹式風(fēng)機(jī)在出風(fēng)口左右兩側(cè)形成對(duì)稱(chēng)的回風(fēng)主流，風(fēng)速為0.6m/s；對(duì)吹式風(fēng)機(jī)回風(fēng)主流出現(xiàn)在風(fēng)機(jī)下側(cè)，風(fēng)速為0.4m/s。風(fēng)機(jī)下吹冷庫(kù)在貨物區(qū)域平均流速高于風(fēng)機(jī)對(duì)吹型冷庫(kù)，風(fēng)機(jī)對(duì)吹冷庫(kù)貨物側(cè)流速主要處于0.1m/s的較低值，原因有兩方面：① 由于側(cè)吹風(fēng)機(jī)在重力和貨物呼吸熱產(chǎn)生的冷熱氣流密度差影響下，流速衰減嚴(yán)重，導(dǎo)致風(fēng)機(jī)遠(yuǎn)端存在大面積送風(fēng)死角；② 由于貨物對(duì)側(cè)吹氣流產(chǎn)生阻礙作用。這一結(jié)果可能導(dǎo)致出風(fēng)主流遠(yuǎn)端貨物區(qū)由于冷量不夠影響換熱，以致部分貨物區(qū)域溫度較高。

圖2 兩種不同風(fēng)機(jī)擺放下吹風(fēng)主流速度分布Figure 2 Comparisons of velocity distribution at the mainstream with two different fan arrangement

圖3為出風(fēng)主流處渦旋分布截面圖。風(fēng)機(jī)下吹庫(kù)在垂直方向近地面處發(fā)現(xiàn)對(duì)稱(chēng)的回流渦旋（圖3（a））。風(fēng)機(jī)對(duì)吹庫(kù)在水平方向靠近風(fēng)機(jī)左右兩側(cè)存在對(duì)稱(chēng)的回流渦旋（圖3（b）），風(fēng)機(jī)對(duì)吹庫(kù)形成的回風(fēng)循環(huán)較短，大部分冷空氣只通過(guò)貨物中部便加入回風(fēng)循環(huán)，并未實(shí)現(xiàn)均勻送風(fēng)。堆垛貯藏中，一般設(shè)立墊倉(cāng)板使貨物與地面間隔一段距離，以實(shí)現(xiàn)貨物底部充分回流換熱并方便堆垛貨物。圖4為兩種風(fēng)機(jī)擺放在冷庫(kù)底部的風(fēng)速分布，風(fēng)機(jī)下吹方式在垂直方向未受到貨物阻礙，且貨架的間隔形成天然風(fēng)道，對(duì)冷空氣進(jìn)行分流，平均流速高于風(fēng)機(jī)對(duì)吹庫(kù)，風(fēng)速呈現(xiàn)放射狀衰減，不存在回流死角，風(fēng)機(jī)對(duì)吹庫(kù)風(fēng)速呈現(xiàn)階梯狀衰減，在貨物底部形成大面積回流死角（見(jiàn)圖4（b）），風(fēng)速僅為0.01～0.02m/s。在同樣風(fēng)量條件下，風(fēng)機(jī)下吹庫(kù)回流更加均勻，有利于冷庫(kù)內(nèi)氣流通暢流動(dòng)。

2.2 溫度分布及對(duì)比

圖3 兩種不同風(fēng)機(jī)擺放下吹風(fēng)主流處渦旋分布Figure 3 Comparisons of vorticity distribution at the mainstream with two different fan arrangement

圖4 兩種不同風(fēng)機(jī)擺放下冷庫(kù)底部流速分布Figure 4 Comparisons of velocity distribution at the bottom of cold store

圖5 兩種不同風(fēng)機(jī)擺放下貨物側(cè)溫度分布三維圖Figure 5 Temperature distributions at the cargo side

圖5為冷庫(kù)貨物溫度分布三維視圖。由圖5可知，下吹風(fēng)方式在貨物區(qū)溫度均勻程度優(yōu)于對(duì)吹風(fēng)方式，風(fēng)機(jī)下吹風(fēng)擺放方式在貨物側(cè)氣流分布不但均勻，且平均流速高于側(cè)吹風(fēng)冷庫(kù)，使得冷空氣能夠較好地與貨物進(jìn)行換熱，下送風(fēng)冷庫(kù)貨物最高溫度在近風(fēng)機(jī)兩側(cè)，溫度為274.8K，對(duì)吹式冷庫(kù)風(fēng)機(jī)遠(yuǎn)端貨物溫度較高且高溫區(qū)面積較大，最高溫度可達(dá)276K，而該區(qū)域恰恰存在上述提到的低流速回流死角，可見(jiàn)在相同風(fēng)量條件下，冷庫(kù)內(nèi)部充分均勻回流對(duì)溫度場(chǎng)合理分布起關(guān)鍵作用，這同時(shí)也印證了Delete等［12］在研究中的發(fā)現(xiàn)。

圖6為兩種風(fēng)機(jī)擺放方式下冷庫(kù)底部溫度分布對(duì)比，下吹風(fēng)式風(fēng)機(jī)在底部溫度較低并且分布均勻（圖6（a）），冷庫(kù)氣流在貨物間隔底部形成兩個(gè)低溫帶，多方向主流回風(fēng)循環(huán)。對(duì)吹式冷庫(kù)在底部存在溫度梯度，風(fēng)機(jī)射流遠(yuǎn)端溫度較高，中間區(qū)域溫度較低（圖6（b））；在射流末端貨物區(qū)域形成了多個(gè)渦旋（圖7），渦旋對(duì)冷庫(kù)氣流產(chǎn)生動(dòng)能耗散作用，使射流末端貨物換熱進(jìn)一步惡化，導(dǎo)致末端貨物溫度較高。

圖6 兩種不同風(fēng)機(jī)擺放下貨物底部溫度分布截面圖Figure 6 Comparisons of temperature distribution at the bottom of cold store

圖7 對(duì)吹風(fēng)冷庫(kù)渦旋速度分布截面圖Figure 7 Vorticity distributions of the airflow in cold store with back blowing fans

3 討論

在本研究中，與風(fēng)機(jī)對(duì)吹方式比較，下吹方式可以顯著改善冷庫(kù)及貨物內(nèi)部的溫度分布。對(duì)吹風(fēng)冷庫(kù)在風(fēng)機(jī)吹風(fēng)主流下部及射流末端貨物區(qū)溫度較高，這是由于對(duì)吹風(fēng)主流在受到熱壓的影響下流速衰減，在射流末端產(chǎn)生多個(gè)渦旋，加上氣流受貨物阻礙形成較短回流循環(huán)，導(dǎo)致流經(jīng)風(fēng)機(jī)遠(yuǎn)端形成大面積送風(fēng)死角，影響貨物換熱。單臺(tái)下吹式風(fēng)機(jī)在庫(kù)內(nèi)可形成多方向主流回風(fēng)循環(huán)，冷庫(kù)內(nèi)回流方向增加無(wú)疑也對(duì)氣流分布起到改善作用。此外，下吹式風(fēng)機(jī)直接對(duì)貨物進(jìn)行送風(fēng)降溫，冷量的利用率高于側(cè)吹式風(fēng)機(jī)，貨物區(qū)域溫度更低也更均勻。在本研究中，冷庫(kù)采用風(fēng)機(jī)下吹方式冷卻效果優(yōu)于對(duì)吹方式。目前，中國(guó)大部分冷庫(kù)采用側(cè)吹風(fēng)型風(fēng)機(jī)，下吹風(fēng)型風(fēng)機(jī)雖在送風(fēng)均勻性上有優(yōu)勢(shì)，但在實(shí)際使用中卻有眾多限制，比如出風(fēng)主流易受貨物阻礙，貨物區(qū)域流速過(guò)高干耗嚴(yán)重等等。

4 結(jié)論

選擇將貨物視為多孔介質(zhì)，即將貨物視為有大量密積微小空隙構(gòu)成的物質(zhì)，通過(guò)計(jì)算貨物內(nèi)部粘性阻力、慣性阻力、孔隙率及其他參數(shù)以實(shí)現(xiàn)貨物內(nèi)部氣流組織準(zhǔn)確預(yù)測(cè)，這種方式的缺點(diǎn)是使建立的冷庫(kù)堆垛貨物模型復(fù)雜化，模擬計(jì)算時(shí)間明顯延長(zhǎng)，收斂速度顯著變慢，但更接近實(shí)際情況。通過(guò)對(duì)模型模擬計(jì)算結(jié)果對(duì)比，得到以下主要結(jié)論：

（1）建立了兩種不同風(fēng)機(jī)擺放方式下冷庫(kù)堆垛貨物模型，并利用多孔介質(zhì)模型對(duì)貨物側(cè)溫度和速度分布進(jìn)行模擬，獲得冷庫(kù)內(nèi)部氣流組織速度和溫度分布，可為優(yōu)化冷庫(kù)貨物擺放和風(fēng)機(jī)設(shè)置提供理論參考。

（2）貨物貯藏過(guò)程中，氣流組織分布均勻性是影響貨物溫度分布的關(guān)鍵因素，對(duì)吹風(fēng)式冷庫(kù)由于受熱壓和重力影響，風(fēng)機(jī)下側(cè)貨物和風(fēng)機(jī)遠(yuǎn)端易形成送風(fēng)死角，溫度較高，下吹式冷庫(kù)在貨物區(qū)域的流速和溫度分布較對(duì)吹式冷庫(kù)均勻。

（3）本研究對(duì)比分析了兩種風(fēng)機(jī)擺放方式下冷庫(kù)內(nèi)部溫度場(chǎng)和速度場(chǎng)分布情況，未來(lái)研究中，可對(duì)比研究?jī)煞N方式下冷庫(kù)相對(duì)濕度的變化。

1 荊有印，王長(zhǎng)海，丁桂艷，等.旋轉(zhuǎn)氣流及貨物對(duì)冷庫(kù)流場(chǎng)影響的研究［J］.冷藏技術(shù)，2008（4）：34～37，41.

2 繆晨，謝晶.冷庫(kù)空氣幕流場(chǎng)的非穩(wěn)態(tài)數(shù)值模擬及驗(yàn)證［J］.農(nóng)業(yè)工程學(xué)報(bào)，2013（7）：246～253.

3 Xie Jing，Qu Xiao-hua，Shi Jun-ye，et al.Effects of design parameters on flow and temperature fields of a cold store by CFD simulation［J］.Journal of Food Engineering，2006，77（2）：355～363.

4 湯毅，謝晶，王金鋒，等.三維計(jì)算流體力學(xué)技術(shù)用于預(yù)測(cè)小型冷庫(kù)內(nèi)氣流分布［C］//中國(guó)制冷學(xué)會(huì).第八屆全國(guó)食品冷藏鏈大會(huì)論文集.北京：中國(guó)制冷學(xué)會(huì)，2012.

5 繆晨，謝晶.空氣幕的研究進(jìn)展［J］.食品與機(jī)械，2012，28（4）：237～262.

6 湯毅，謝晶，王金鋒，等.計(jì)算流體力學(xué)在冷庫(kù)優(yōu)化中的應(yīng)用研究進(jìn)展［J］.食品與機(jī)械，2011，27（5）：186～189.

7 劉妍玲，張巖，王世清，等.果蔬擺放形式對(duì)冷庫(kù)內(nèi)氣流場(chǎng)分布影響的研究［J］.青島農(nóng)業(yè)大學(xué)學(xué)報(bào)（自然科學(xué)版），2008，25（1）：24～27.

8 Chourasia M K，Goswami T K.Simulation of effect of stack dimensions and stacking arrangement on cool-down characteristics of potato in a cold store by computational fluid dynamics［J］.Biosystems Engineering，2007，96（4）：503～515.

9 Delele M A，Schenk A，Tijskens E，et al.Optimization of the humidification of cold stores by pressurized water atomizers based on a multiscale CFD model［J］.Journal of Food Engineering，2009，91（2）：228～239.

10 Tanaka F，Konishi Y，Kuroki Y，et al.The use of CFD to improve the performance of a partially loaded cold store［J］.Journal of Food Process Engineering，2012，35（6）：874～880.

11 Ergun S，Orning A A.Fluid flow through randomly packed columns and fluidized beds［J］.Industrial ＆Engineering Chemistry，1949，41（6）：1 179～1 184.

12 Delele M A，Schenk A，Tijskens E，et al.Optimization of the humidification of cold stores by pressurized water atomizers based on a multiscale CFD model［J］.Journal of Food Engineering，2009，91（2）：228～239.

13 Delele M A，Vorstermans B，Creemers P，et al.CFD model development and validation of a thermonebulisation fungicide fogging system for postharvest storage of fruit［J］.Journal of Food Engineering，2012，108（1）：59～68.

14 Lisowa H，Wujec M，Lis T.Influence of temperature and variety on the thermal properties of apples［J］.International Agrophysics，2002，16（1）：43～52.

15 Schotsmans W，Verlinden B E，Lammertyn J，et al.Factors affecting skin resistance measurements in pipfruit［J］.Postharvest Biology and Technology，2002，25（2）：169～179.