葛蒙蒙，陳桂香，劉文磊，劉超賽

河南工業(yè)大學 土木工程學院，河南 鄭州 450001

倉儲環(huán)境的溫濕度是影響糧食品質(zhì)的主要因素，糧食本身的生命活動與糧倉內(nèi)外環(huán)境間的熱量傳導影響糧堆的溫濕度[1-2]。因此，掌握并準確預測倉儲糧堆的溫度變化是確保糧食安全儲藏的基礎。由于糧食在儲藏過程中受到生物和非生物因素的多重影響，倉儲糧堆溫度變化難以得到準確、有效地預測。

研究糧堆熱濕傳遞規(guī)律方法有試驗和數(shù)值模擬兩種。試驗結果往往是可靠的，但試驗過程會比較煩瑣，甚至需要投入大量的人力、物力資源，而且每次試驗得出的結果不盡相同，具有不可重復性。數(shù)值模擬可以很好地解決試驗方面的不足，它可以準確有效地分析和預測糧倉生態(tài)系統(tǒng)的變化。Thorpe等[3]在分開考慮倉儲糧堆間氣體溫度和糧食本身溫度的情況下，建立了倉儲糧堆溫度變化的數(shù)學模型。Jia等[4]構建了笛卡爾坐標系下的圓筒倉糧堆溫度場模型。Gaston 等[5]利用有限元法，模擬了在外界環(huán)境季節(jié)性變化以及小麥的呼吸作用下，小麥糧堆溫度和水分變化規(guī)律。張前等[6]根據(jù)試驗回歸模型預測倉儲糧堆的溫度變化。李軍軍等[7]依據(jù)有限元法建立了倉儲糧堆溫度場數(shù)學模型。蔣華偉等[8]構建了在單一局部發(fā)熱下的倉儲糧堆溫度場的數(shù)學模型。梁醒培等[9]利用數(shù)值模擬的方法建立了小麥糧堆的溫度場數(shù)學模型，并且研究了小麥在不同季節(jié)的溫度變化。門艷忠[10]利用有限元分析法和傳熱學相關理論模擬了非穩(wěn)態(tài)下稻谷的溫度變化，結果表明糧堆內(nèi)會出現(xiàn)一個相對高溫區(qū)，但由于其數(shù)學模型不夠精細，所以與實際情況有一定的差距，不能很好地應用于實際。劉建英[11]根據(jù)實倉內(nèi)3個高度層的測溫數(shù)據(jù)，推導出糧溫經(jīng)驗公式，但其僅能反映倉內(nèi)局部的溫度，并不能精確地反映倉內(nèi)其他位置的溫度變化。周全申等[12]通過現(xiàn)場實測，統(tǒng)計分析了不同氣候條件對糧堆溫差的影響。閆艷霞等[13]依據(jù)熱力學相關理論與平衡理論等構建了局部發(fā)熱條件下的糧堆溫度場數(shù)學模型，發(fā)現(xiàn)倉儲糧堆溫度總是處于一種動態(tài)變化當中。白忠權[14]在考慮谷物吸附和解吸的情況下，利用數(shù)值模擬與試驗驗證，得出谷物糧堆傳熱傳質(zhì)的變化規(guī)律。

作者利用COMSOL軟件建立靜態(tài)儲糧的平房倉三維物理模型，對糧倉以及糧堆的傳熱過程進行數(shù)值模擬,分析了倉儲稻谷糧堆整體平均糧溫變化，并與實測數(shù)據(jù)進行對比。

1 模型模擬

1.1 平房倉三維物理模型

模擬倉為浙江某高大平房倉，平房倉長29.37 m，寬20.44 m，高8.2 m，糧堆高6 m。墻體采用空心磚墻結構，倉頂采用預應力混凝土板，倉底為混凝土地面。利用COMSOL軟件建立高大平房倉三維物理模型，模擬計算域由糧堆與倉內(nèi)上部空氣層兩個體區(qū)域組成。對物理模型采用非結構網(wǎng)格極細化處理，共劃分為170多萬個網(wǎng)格。糧倉模型及網(wǎng)格劃分如圖1所示。

圖1 高大平房倉三維物理模型及網(wǎng)格劃分Fig.1 Three-dimensional physical model and mesh division diagram of large horizontal warehouse

1.2 模擬參數(shù)

1.2.1 糧堆外界環(huán)境溫度

對于靜態(tài)儲糧，引起糧堆溫度變化的主要因素是外界環(huán)境溫度的變化。將模擬倉2017年外界月平均氣溫變化繪成曲線，結果如圖2所示。

圖2 浙江某高大平房倉2017年外界月平均氣溫變化Fig.2 Average temperature change outside a large horizontal warehouse in Zhejiang province in 2017

對圖2曲線進行近似擬合得出外界溫度變化函數(shù)：

T=-0.606 4(t/2 592 000-7.07)2+197.42，

式中：T為外界環(huán)境平均溫度,K；t為時間,s。

倉壁和糧堆的初始溫度取自2016年12月底的平均溫度，全年環(huán)境溫度變化采用上式所示溫度，每月按30 d進行計算。

1.2.2 糧倉壁面邊界

靜態(tài)儲糧過程中，糧倉壁面的溫度變化會直接影響糧堆溫度。外界環(huán)境通過改變倉壁溫度，進而影響糧堆溫度。在長時間自然儲糧過程中，外部的傳熱主要是熱輻射與對流換熱，在考慮太陽輻射條件時，糧倉倉壁與其倉頂所受太陽輻射量不盡相同，此時模擬計算過程不僅復雜且結果會產(chǎn)生很大誤差。為簡化計算，將倉壁與倉頂所受太陽輻射轉(zhuǎn)化為倉壁溫度，倉底部為絕熱層。

1.2.3 模型參數(shù)

COMSOL模擬靜態(tài)儲糧傳熱過程，研究發(fā)現(xiàn)影響糧堆傳熱計算的主要參數(shù)有稻谷堆的導熱系數(shù)、比熱容、密度與孔隙率等。具體參數(shù)值如表1所示[15]。

表1 COMSOL數(shù)值模擬的相關條件和材料參數(shù)設置Table 1 Relevant conditions and parameter settings of COMSOL numerical simulation

1.3 控制方程

稻谷儲藏期間涉及的方程有質(zhì)量守恒、動量守恒以及能量守恒方程。

質(zhì)量守恒方程：

式中：ρa為空氣密度，kg/m3；t為時間，s；u為風速，m/s。

動量守恒方程：

式中：μ是空氣動力黏度，Pa·s；K為滲透率,m2；g為重力加速度，m/s2；β是空氣膨脹系數(shù),1/K；δij為kronecker函數(shù);T為稻谷糧堆的溫度，K；T0為初始溫度，K；p為壓強，Pa。

能量守恒方程：

式中：ρg為稻谷糧堆的干密度，kg/m3;Cg為稻谷糧堆的比熱容，kJ/(kg·K);Ca為空氣的比熱容，kJ/(kg·K)；Kg為稻谷糧堆的導熱系數(shù)，W/(m·K)；Wg為干基水分；hfg是糧堆的蒸發(fā)潛熱，kJ/kg；qh是呼吸過程釋放的熱量，J/mg；YCO2是糧堆呼吸24 h的CO2的釋放率。

在自然儲藏過程中不存在強制的通風過程，所以不用考慮糧堆與孔隙間的氣流溫差，進而可忽略由溫度變化引起的自然對流傳熱。糧倉中以倉壁與糧堆間的熱傳導為主，其熱量傳導的數(shù)學模型：

Q=-kT，

式中：Q為熱通量密度，kJ/(m2·s)；k為糧堆的熱傳導系數(shù)，W/(m·K)；T為糧堆的溫度，K。

2 結果與分析

2.1 糧倉寬度方向上糧堆的溫度分布

圖3為模擬糧倉在寬度方向上中間垂直截面的溫度分布變化情況。

圖3 糧倉寬度方向中垂面溫度分布變化情況Fig.3 Variation of vertical temperature distribution in the width direction of the warehouse

由圖3可知，稻谷糧堆溫度沿糧倉寬度方向呈現(xiàn)出很明顯的不均勻分布。儲藏前3個月，近壁面糧堆溫度隨外界環(huán)境溫度上升逐漸上升，并逐漸向糧堆內(nèi)部傳遞；儲藏6個月時，隨著外界環(huán)境溫度值的進一步升高，距壁面稍遠處糧堆溫度也隨之上升；儲藏9個月時，此時外界環(huán)境溫度雖然已經(jīng)開始下降，但除近壁面糧堆溫度略有下降外，糧堆內(nèi)部溫度仍有上升趨勢；儲藏12個月時，隨著外界環(huán)境溫度的持續(xù)下降，近壁面糧堆溫度也隨之下降，但由于內(nèi)部熱量不能及時傳遞出去，近壁面會出現(xiàn)相對高溫區(qū)。倉頂空氣層是糧堆上層沒有出現(xiàn)高溫區(qū)的主要原因，是由于空氣熱傳導性弱于倉頂混凝土，所以能夠在一定程度上穩(wěn)定糧堆的溫度。

2.2 糧倉長度方向上糧堆的溫度分布

圖4為模擬糧倉在長度方向上中間垂直截面的溫度分布變化情況。由圖4可知，糧倉長度方向上稻谷糧堆的溫度變化情況與寬度方向上相似，但不同方向上糧堆溫度變化受外界環(huán)境影響程度不同，糧倉長度方向上糧堆高溫帶的范圍和規(guī)模均小于寬度方向。因此，在靜態(tài)儲藏過程中，外界環(huán)境對糧堆溫度在糧倉寬度方向上的影響要大于長度方向。

圖4 糧倉長度方向中垂面溫度分布變化情況Fig.4 Variation of vertical temperature distribution in the length direction of the warehouse

為了更直觀地分析稻谷糧堆內(nèi)部溫度分布及變化情況，選取堆高為3 m平面上沿長度方向的中軸線研究糧堆溫度的變化情況，結果如圖5所示。

圖5 3 m處糧堆中軸線沿糧倉長度方向溫度分布Fig.5 Temperature distribution of the central axis of the grain pile at 3 m along the length of the warehouse

由圖5可知，2—8月，近壁面處稻谷糧堆溫度最高，越往中間，糧堆溫度越低，尤其是6—8月最為明顯，此時外界環(huán)境逐漸達到1年中的最高溫，近壁面處糧溫很快達到全年最高值24.1 ℃。9—12月，此時外界環(huán)境越過1年中的最高溫開始下降，近壁面處高溫糧堆溫度也隨之開始下降，但其積聚的熱量則繼續(xù)向內(nèi)傳遞，這樣就出現(xiàn)糧堆溫度由壁面往內(nèi)迅速上升，達到最大值后逐漸下降的現(xiàn)象。在距糧倉壁面1～3 m處，糧堆產(chǎn)生了明顯的高溫區(qū)，且在12月最為顯著，距糧倉壁面約2 m處高溫帶溫度達到峰值。在距糧倉壁面較遠處的糧堆隨儲藏時間延長，其糧溫有所提高，儲藏12個月的最高糧溫為14.7 ℃。在距壁面4 m左右，無論是長度方向還是寬度方向上，糧溫始終處于一個穩(wěn)定的低溫區(qū)。

2.3 試驗驗證

為驗證構建的研究模型，將倉儲稻谷糧堆整體平均糧溫實測數(shù)據(jù)與模擬結果進行比較分析，試驗時用糧情檢測系統(tǒng)測定高大平房倉稻谷糧堆溫度的變化，將高大平房倉1年內(nèi)稻谷糧堆整體平均溫度模擬結果和實測結果對比，如圖6所示。

圖6 高大平房倉一年內(nèi)糧堆平均溫度變化模擬與實測對比結果Fig.6 Comparison of simulated and measured results of average temperature change of grain pile in a large warehouse within one year

由圖6可知，模擬與實測結果變化趨勢基本一致，糧堆溫度實測值略高于模擬值，最大誤差為0.91 ℃，在誤差允許范圍內(nèi)。

3 結論

利用COMSOL軟件對靜態(tài)儲藏稻谷糧堆1年內(nèi)的溫度變化進行了研究，得出如下結論：在靜態(tài)儲糧過程中糧倉外界環(huán)境尤為重要，糧堆溫度受外界溫度變化而改變。近壁面處糧溫受外界環(huán)境溫度變化影響最為顯著，在距糧倉壁面較遠處的糧堆隨儲藏時間延長，其糧溫有所提高。當外界環(huán)境溫度經(jīng)過由低到高再到低的變化后，在距糧倉壁面1～3 m處，稻谷糧堆產(chǎn)生了明顯的高溫區(qū)。糧倉長度方向上稻谷糧堆的溫度變化情況與寬度方向上相似，但在靜態(tài)儲藏過程中，外界環(huán)境對糧堆溫度在糧倉寬度方向上的影響要大于長度方向。將模擬與實測結果對比可知，兩者變化趨勢基本一致，稻谷糧堆溫度實測值略高于模擬值，最大誤差為0.91 ℃。利用COMSOL對糧堆溫度傳遞進行模擬是一種簡便且有效的途徑，可以為今后的優(yōu)化工作提供更加直觀可靠的研究工具。