基于溫濕度場耦合的糧堆離散測點溫度場重現分析

尹 君 吳子丹 吳曉明 張忠杰 趙會義 李福君 吳文福

（吉林大學生物與農業工程學院1，長春 130025）

（國家糧食局科學研究院2，北京 100037 ）

儲糧生態系統主要包括非生物子生態系統和生物子生態系統兩部分，其中糧堆包含于生物子生態系統中。糧堆溫度對糧堆內的水分遷移、微生物滋長和蟲害活動等都有十分重要的影響，是決定儲糧生態系統安全的關鍵因子之一[1]，因此準確掌握糧堆溫度場分布、通過溫濕度場的耦合分析預測糧堆溫度場變化趨勢，是預判儲糧安全狀態的重要方法之一。目前，糧堆溫度場的研究方法主要有試驗研究、數值模擬研究及半試驗半仿真研究等。由于試驗研究存在成本投入高、樣本量大、試驗結果重復性差等弊端，許多學者開始采用數值模擬方法或半試驗半仿真方法，研究糧堆溫度場分布及變化。Thorpe等[2]將糧食籽粒溫度與糧堆孔隙間氣體溫度分開考慮，構建了儲糧過程中糧堆溫度場變化的數學模型。Gastón等[3]和 Iguaz等[4]依據糧食水分等溫吸附方程和熱質局部平衡原理，視糧堆為多孔介質，建立了糧堆熱濕耦合模型。Xu等[5]研究了糧堆孔隙率對儲糧通風系統均勻性的影響，并建立了相關數學模型。Jian等[6]利用長期實測數據探明了加拿大北部地區金屬鋼板倉小麥溫度隨外界環境變化規律和建立了糧堆水分遷移、再分配模型。Jia等[7]在笛卡爾坐標系下建立了二維圓筒倉糧堆溫度場模型。Rehman等[8]、吳子丹等[9]、李鐵盤等[10]利用短期實測數據運用多元回歸方法構建了回歸模型，并實現了對糧堆各測點溫度短期（3 d）預測。張前等[11]利用試驗回歸模型預測糧層內的溫度變化。李軍軍等[12]采用有限元分析法，構建糧堆溫度場變化數學模型，在此基礎上開發了空倉預警系統。前人運用CFD模擬技術，模擬分析不同倉型的準靜態儲藏和機械通風2種狀態下糧堆溫度場的分布及變化趨勢[13-18]。但是多數的模擬研究尚未實現對復雜的糧堆精確仿真，難以準確預測預警糧情的演變。吳子丹等[17]成功開發了智能糧情檢測、低劑量環流熏蒸、智能通風和高效谷物冷卻4項技術為一體的“四合一”儲糧新技術，并廣泛應用于“十五”以來使用國債資金建設的550×108kg的淺圓倉、高大平房倉；其中智能糧情檢測技術中采用了將糧堆歷史溫度場變化數據進行回歸建模，并用實時檢測數據加以自動修正的方法，實現了糧溫異常判斷和預警[18]，后來又提出糧堆多場耦合理論，其中基于溫濕度場耦合理論的糧情預測研究是重點之一[1]，這也為儲糧生態系統理論體系提出了新的研究思路。

經實倉實踐表明：糧堆高度較高、柱狀的糧倉如淺圓倉、立筒倉，容易在冬春季節發生倉壁和糧面結露而導致霉變。本試驗通過對我國北方一個典型的淺圓倉進行模擬研究，利用溫度傳感器陣列采集不同季節、不同點位的糧堆溫度，應用糧溫擬合算法和Matlab模擬軟件構建糧堆溫度場模型，將實測數值和計算數值與時空域相關聯，重現糧堆溫度場變化，并采用溫濕度場耦合理論進行云圖分析，試圖找出導致淺圓倉發生結露的機理和征兆，為預判儲糧安全提供理論依據。

1 糧堆溫度場重現與溫濕度場耦合

1.1 淺圓倉物理模型

研究對象為中央儲備糧天津東麗直屬庫21號淺圓倉（圖1），倉內直徑為22.8 m，倉檐高為16.3 m，倉頂高為22.5 m，糧堆高度為15.0 m。

圖1 淺圓倉物理模型

1.2 儲糧糧堆及測試儀器

本倉所儲為河北產的混合小麥，水分為11.8%，儲量為5 502.8 t。選用由天津明倫電子科技有限公司生產的糧情測控系統對儲糧糧堆離散點進行糧溫監測[20]。選用由上海精宏實驗設備有限公司生產的DHG型電熱恒溫干燥箱對采取105℃恒重法對抽取的樣品進行水分檢測[21]。

1.3 分布式溫度測試系統

在淺圓倉的xoy平行剖面上，共分3圈測溫點，由內而外3圈測溫點半徑和測點數分別為3.7、7.4、11.1 m；3點、8點和12點（圖2）。即最外圈測溫點距倉壁0.3 m，每層共有23個測溫點。其每圈測溫點布置如圖2所示。全倉共有8層，上、下層測溫點分別距糧面、糧底0.5 m，其余6層測溫點間距為2.0 m，全倉測溫點共計184（23×8）個。

圖2 xoy平行剖面上測溫點布置圖

1.4 糧堆溫度場重現的擬合算法

由數值計算方法可知，已知n個離散數據（xi，f（xi））i＝1，2，3…n，在[a，b]上滿足

式中Pn（x）稱為牛頓插值函數，本試驗研究選用的是二次牛頓插值法，余項R2（x）已滿足試驗精度要求。再運用Matlab模擬軟件重現糧堆各個剖面的等溫線云圖。

1.5 糧堆溫濕度場耦合算法

基于糧食平衡濕度原理，將糧堆任一點的溫度值、糧食水分和小麥糧種擬合系數值代入WU（CAE）模型[17]之一，耦合得到糧堆各點的糧食平衡絕對濕度值：

式中：EAH為糧堆內任一點的水汽分壓值／Pa，M為糧食濕基含水率／%，t為該點的糧食溫度／℃，A1、A2、B1、B2、D為不同糧種解吸和吸附不同狀態下的擬合系數，本研究選用小麥的擬合系數值。

再代入WU（CAE）模型的又一表達式，計算糧堆任一點水汽分壓達到飽和狀態時的溫度，即露點方程：

式中：DPT為糧堆內任一點的露點溫度／℃。

2 結果與討論

2.1 不同季節糧堆溫度場云圖分析

2.1.1 夏季糧堆不同方向溫度場云圖分析

由圖3可知，糧堆受外界環境的影響很大。同一時刻不同高度糧堆溫度場分布略有差異。夏季，糧堆表層平均溫度最高，糧堆底層平均溫度次之，糧堆中心區域平均溫度最低。倉頂是受太陽輻射、空氣對流、熱傳導影響最大的位置，故靠近表層的糧溫隨外界變化最為劇烈，糧堆中心區域受外界影響速度最緩慢、幅度最小。

分析圖3、圖4可知，由于糧食是熱的不良導體，小麥導熱系數僅為0.15W／（m·℃），盡管受夏季高溫影響的糧堆表層最高溫度已達到近30℃。而高熱區域只出現在近糧表面不足1 m處，在糧堆內部仍然保持一個巨大的“冷芯”。按照糧堆溫度濕度場耦合原理，由 WU（CAE）模型[17]可測算得到熱區域與冷區域的平衡水汽分壓分別為：當含水量為11.8%的小麥在熱區域30℃時，水汽分壓EAH＝184.0 Pa（解吸），當含水量為11.8%的小麥在冷區域邊緣3℃時，水汽分壓EAH＝33.4 Pa（吸附）。由于熱區域水汽分壓遠遠大于冷區域，據局部熱質平衡原理和熱濕轉移規律可知，熱量、水分將沿著溫度梯度方向進行遷移，即由糧堆表面向糧堆內部逐漸擴散。一方面，糧堆表層雖然溫度較高，但隨水分散失其含水量會降低1%～2%，又根據糧食儲藏穩定性的溫度、水分雙因子負相關原理可知，此時糧堆表層的糧食是安全的。另一方面，由于“冷芯”的體積遠大于糧堆表面的熱區體積，即V冷＞＞V熱，由糧堆表層傳入的有限水分被巨大“冷芯”吸收，冷區域的平均水分上升幅度極小，在低溫環境下糧堆內部也處于安全狀態。從溫度場云圖分析可以看出，如果夏季糧堆中存有一個巨大的“冷芯”，即使隨外界溫度上升糧堆表層溫度略有升高，但糧堆整體仍可處于安全狀態。據此分析，夏季低溫儲糧的關鍵是盡可能保持糧堆的“冷芯”，而那種要求全倉、包括易受外界影響的糧堆表層也降低溫度的做法，既耗費制冷能源、又無任何實際意義。

圖3 2012年8月30日21號淺圓倉xoy平行剖面（Z＝0.5、8.5、14.5 m）溫度場云圖

圖4 2012年8月30日21號淺圓倉東西向y＝0垂面（XZ）溫度場云圖

2.1.2 秋、冬季節糧堆不同方向溫度場云圖分析

隨著一個夏季的熱量傳入和積累，糧堆上部出現溫度較高區域，“冷芯”有所減小，但是對于5 502.8 t的大糧堆而言，總體仍保持了較大的冷區域體積。圖5和圖6展示了秋、冬季節糧堆中垂面溫度場分布變化。在圖6中糧堆出現了3個不同溫度區域，靠近倉底受地溫影響出現厚度不足0.5m、跨度小于8m、均溫約在13℃的偏暖區域，此區域的熱量、水分均隨著溫度梯度向上進行遷移，但是由于此區域上部的冷區域體積較大，故其傳出的熱量、水分對冷區域的影響可忽略。另外，在糧堆上層出現一個由積熱產生的厚度約2 m、跨度約18 m、均溫約在13℃、最高溫度達到17℃的暖區域。此區域的熱量在向外散發的同時，其熱量、水分也向糧堆內部冷區域進行擴散、遷移，但是下方的冷區域體積遠遠大于暖區域體積，加之糧食的比熱容較大[小麥為1.86 kJ／（kg·℃）]，所吸收的能量分布于冷區域各點后，并未使冷區域發生較大的變化，故此時糧堆仍處于安全狀態。

圖5 2011年9月20日21號淺圓倉東西向y＝0垂面（XZ）溫度場云圖

2.1.3 冬末春初季節糧堆不同方向溫度場云圖分析

隨著秋冬季外界氣溫不斷降低，糧堆表層溫度受到頂層空氣的影響會逐漸下降，致使糧堆表層出現冷區域，來年春季糧堆出現新的溫度分層變化。圖7～圖9展示了冬末春初糧堆中垂面和不同高度水平剖面溫度場分布圖。糧堆出現了4個不同溫度段的區域，靠近倉底是厚度不足1m、跨度小于15m、均溫約在12℃的相對偏暖區域，如上所述，此區域傳出的熱量、水分對其上部厚度約8 m、跨度約22 m的大“冷芯”而言影響可忽略。糧堆表層1 m之下存在一個厚度約3 m、跨度約18 m、均溫約在12℃、最高溫度達到14℃的相對熱區域，此“熱芯”是夏季積熱形成的。對此熱區域下方的巨大冷區域而言，其溫度、水分的變化量較小也可忽略。值得關注的是，在此熱區域上方還有一個小的冷區域。即糧堆表層冷區域體積（厚度不足1 m、跨度約20 m）比此熱區域體積小很多，熱區域的糧食水分隨著溫度梯度方向向表層傳遞，然而此冷區域的吸熱能力很強，在2月外溫較低時節，熱量不斷通過倉壁散失到外界，而水分無法穿透倉壁故滯留于糧堆表層和倉壁內，水分逐步積累最終將形成結露。

由 WU（CAE）模型[17]可以測算出，該熱區域的糧食露點溫度為5.6℃（吸附），由云圖顯示距離熱區域僅1 m的糧堆表層糧溫僅-0.9～3.6℃，已經低于熱區域的糧食露點，尤其是靠近糧堆表層中部、表層邊部方向由于溫度梯度最大，為最早可能發生結露位置，屬于危險區域。

圖6 2011年12月22日21號淺圓倉東西向y＝0垂面（XZ）溫度場云圖

圖7 2013年2月8日21號淺圓倉東西向y＝0垂面（XZ）溫度場云圖

圖8 2013年2月8日21號淺圓倉南北向x＝0垂面（YZ）溫度場云圖

圖9 2013年2月8日21號淺圓倉xoy平行剖面（Z＝0.5、8.5、12.5、14.5 m）溫度場云圖

由圖9亦知，糧堆中有4個不同溫度區域時，“冷芯”區域較大、糧溫較均衡、分布較勻，故穩定性較好。但是糧堆上部的熱區域與糧堆表層冷區域互相影響較大，此時外界環境溫度較低，其向糧堆表層傳遞的熱量和水分對糧堆表層影響更為突出，水分隨熱量不斷遷移到表層和倉壁，逐漸累積而升高，故在冬春季之際，近糧堆表層下方厚度小于1 m的冷區域易出現結露。據此分析，在秋冬季節，如果糧堆中存在“熱芯”，是糧堆可能結露的明顯征兆。

2.2 機械通風作業后糧堆溫度場分布云圖分析

為了消除糧堆安全隱患，可借助冬季外界低溫環境的有利時機，通常于當年12月至次年1月進行間歇式機械通風作業，消除可能導致糧堆結露的“熱芯”和溫差。圖10所示為同一糧倉在機械通風作業結束后的春季糧堆不同高度水平剖面溫度場分布云圖。各個高度糧堆溫度場分布基本相同，整個糧堆形成了糧溫均衡的全冷區域。隨著氣溫逐步上升，這個均衡的冷區域將再次被周邊的熱糧包圍，形成新的“冷芯”。由此分析，冬季適時通風降溫，不僅可以消除當期的儲糧隱患，而且為全年安全儲糧制造了“冷芯”條件。

圖10 春季（3月8日）經冬季通風的21號淺圓倉

3 結論

利用溫度傳感器陣列測得不同季節糧堆溫度真實值，運用擬合算法和Matlab模擬軟件重現糧堆溫度場在不同季節、不同方位的變化，通過云圖分析，得出以下結論：

3.1 在夏季高溫的外界環境影響下，糧堆頂部、底部、倉壁之處糧堆溫度上升，但是如果糧堆中心存在大體積的“冷芯”作用，仍可確保糧堆安全度夏。

3.2 在秋冬之際，較高的糧堆易出現“熱-冷”區域分層，糧堆較大熱區域在散熱過程中易形成外冷內熱的“熱芯”，由“熱芯”驅動的濕熱轉移是導致結露發生的內因。

3.3 在冬春交替之際，糧堆中如果出現明顯的體積較大的“熱芯”，在溫度梯度最大的方位，即靠近倉壁、糧堆表層處易發生結露。

3.4 選擇冬季有利時機進行機械通風調控，驅散糧堆“熱芯”，使糧堆總體處于低溫均衡狀態，消除了濕熱轉移的驅動源，可以有效防控糧堆結露現象的發生；同時也為安全度夏創造了“冷芯”條件。

3.5 借助構建糧堆溫度場模型，重現糧堆溫度場分布和變化，運用糧堆溫度濕度場耦合理論進行云圖分析，可提前預測結露，為預測較長時期糧堆的儲藏安全穩定性預警提供了新方法。但是，本研究提出的特定溫度濕度場耦合條件導致結露的概念，仍需要進行定量研究才能準確預報發生結露的時間和相應的位置；對于不同的糧堆大小和形狀、不同的糧種和雜質分布、不同的氣候類型，本方法的運用方式也需要進一步研究完善。

[1]吳子丹，趙會義，曹陽，等.糧食儲藏生態系統的仿真技術應用研究進展 [J].糧油食品科技，2014，22（1）：1-6

[2]Thorpe G R.Modelling ecosystems in ventilated conical bottomed farm grain silos[J].Ecological Modelling，1997，94：255-286

[3]Gastón A，Abalone R，Bartosik R E，et al.Mathematicalmodelling of heat and moisture transfer of wheat stored in plastic bags（silobags）[J].Biosystems Engineering，2009，104：72-85

[4]Iguaz A，Arroqui C，Esnoz A，et al.Modeling and Validation of Heat Transfer in Stored Rough Rice without aeration[J].Biosystems Engineering，2004，88（4）：429-439

[5]Xu S，Digvir S Jayas，White N D，et al.Momentum-diffusivemodel for gas transfer in granularmedia[J].Journal of Stored Products Research，2002，38（5）：455-462

[6]Jian F J，Digvir S J，White N D G.Temperature fluctuations and moisture migration in wheat stored for 15 months in a metal silo in Canada[J].Stored Products Research，2009，45：82-90

[7]Jia C C，Sun DW，Cao CW.Computer simulation of temperature changes in a wheat storage bin[J].Journal of Stored Products Research，2001，37：165-177

[8]Rehman Z U，Habib F，Zafar S I.Nutritional changes in maize（Zeamays）during storage at three temperatures[J].Food Chemistry，2002，77：197-201

[9]吳子丹，楊萬生.我國糧油科學技術學科的發展與展望[J].糧油食品科技，2011，19（3）：1-4

[10]李鐵盤，宋友林.糧庫溫度分布函數及應用 [J].鄭州工程學院學報，2004，25（3）：77-79

[11]張前，周永杰，辛立勇，等.高大平房倉儲糧溫度變化規律及數學模型研究[J].糧食儲藏，2003，32（6）：25-30

[12]李軍軍，傅宏.基于糧堆溫度場數學模型的空倉預警模塊研究 [J].軟件導刊，2012，11（4）：93-94

[13]王遠成，張忠杰，吳子丹，等.計算流體力學技術在糧食儲藏中的應用 [J].中國糧油學報，2012，27（5）：86-91

[14]彭威，張忠杰，任廣躍，等.倉儲糧堆溫度場CFD模擬應用研究 [J].糧油食品科技，2011（6）：5-8

[15]張忠杰，李瓊，楊德勇，等.準靜態倉儲糧堆溫度場的CFD模擬 [J].中國糧油學報，2010，25（4）：45-50

[16]陳桂香，岳龍飛，林淼達，等.糧堆內熱濕耦合傳遞數值模擬與試驗驗證 [J].中國糧油學報，2014，29（1）：77-82

[17]吳子丹.綠色生態低碳儲糧新技術 [M].北京：中國科學技術出版社，2011

[18]吳子丹，曹陽，趙會義.糧情預測系統、糧情控制系統及糧堆安全狀態判別方法：中國，102156007A[P].2013-2-20

[19]歐陽曙光，呂青青，付樂樂，等.應用Matlab繪制活性炭吸附等溫線和孔徑分布圖 [J].實驗室研究與探索，2012，31（9）：26-28

[20]GB／T 29890—2013，糧油儲藏技術規范[S]

[21]GB／T 5497—1985，糧食、油料檢驗 水分測定法 [S].