王 柯 王遠(yuǎn)成 劉濟(jì)洲 尹 君 尉堯方 余 海

(山東建筑大學(xué)熱能工程學(xué)院1，濟(jì)南 250101) (國家糧食和物資儲備局科學(xué)研究院2，北京 100037)

世界各國因儲糧害蟲對糧食造成的損失非常嚴(yán)重，一般損失率為10%～40%[1]，所以研究儲糧通風(fēng)過程中害蟲數(shù)量的變化至為重要。小麥在儲藏期間的重要害蟲主要有谷蠹、米象和玉米象等，溫度、濕度、水分等生態(tài)因素會影響儲糧害蟲種群生長及分布。目前在實(shí)際儲糧過程中控制害蟲增長的方法主要有采用殺蟲劑、熏蒸劑、通風(fēng)降溫等方式，其中，最常見的方法是采用機(jī)械通風(fēng)，即利用通風(fēng)機(jī)將冷空氣送入糧堆內(nèi)部，與糧堆進(jìn)行熱濕傳遞，降低糧堆的溫度，從而抑制害蟲的繁衍，降低糧食在儲藏期間的損耗。

利用計算機(jī)編程來預(yù)測實(shí)倉通風(fēng)過程中的溫度、濕度、水分以及生物場等變化規(guī)律，是近年來興起的一種既科學(xué)又高效的研究方法。國內(nèi)外研究人員借助計算流體力學(xué)模擬軟件對儲糧通風(fēng)過程的研究已經(jīng)取得豐碩的成果。Hunter等[2]設(shè)計了一種通風(fēng)降溫系統(tǒng)，可以形成抑制害蟲增長的儲糧環(huán)境。Thorpe等[3]提出了儲糧通風(fēng)中有關(guān)糧粒內(nèi)部的熱濕耦合與傳遞、殺蟲劑衰減的有限差分?jǐn)?shù)學(xué)模型，得出糧食初始溫度、水分含量以及通風(fēng)時長等因素對殺蟲劑的衰減率的影響。朱振剛等[4]研究了圓形糧倉的通風(fēng)過程，對比得出立式通風(fēng)形式更有利于糧食儲藏。王遠(yuǎn)成等[5-7]對國內(nèi)外儲糧系統(tǒng)的數(shù)學(xué)模型、數(shù)值模擬方法進(jìn)行了綜述，并模擬預(yù)測了高大圓筒倉內(nèi)溫度場、水分場等。

利用商業(yè)CFD(Computational Fluid Dynamics)軟件模擬儲糧通風(fēng)過程的文獻(xiàn)已經(jīng)屢見不鮮。由于各種條件的限制，糧食倉儲企業(yè)使用商用CFD軟件還有諸多不便。本研究基于多孔介質(zhì)熱濕耦合傳遞理論，建立了淺圓倉的糧堆內(nèi)部熱濕傳遞和流動的數(shù)學(xué)模型，采用Fortran語言編程模擬了儲糧通風(fēng)后糧堆內(nèi)部溫度、水分含量、殺蟲劑的濃度衰減規(guī)律以及谷蠹的分布情況，為優(yōu)化儲糧通風(fēng)工藝，降低糧食在儲藏期間的損耗，提高經(jīng)濟(jì)效益提供參考。

1 儲糧通風(fēng)多場耦合模型及其實(shí)驗(yàn)驗(yàn)證

1.1 實(shí)驗(yàn)倉和實(shí)驗(yàn)工況

儲糧通風(fēng)實(shí)驗(yàn)倉結(jié)構(gòu)及尺寸如圖1所示，其中圖1a、圖1c分別為糧倉的正視圖和俯視圖。筒倉高度為4.5，直徑為4.2 m。通風(fēng)道采用十字型通風(fēng)道，其中十字型通風(fēng)口長度為1.2 m。

注：單位：mm。圖1 倉型結(jié)構(gòu)圖及溫度傳感器布置圖

實(shí)驗(yàn)糧種為小麥，初始溫度為40 ℃，初始濕基狀態(tài)下含水量為12%，噸糧通風(fēng)量為49.68 m3/h· t。每天通風(fēng)6 h，共通風(fēng)10 d，利用布置在糧堆內(nèi)部的探頭測出糧堆溫度的變化，水分變化采用埋袋取樣進(jìn)行水分化驗(yàn)。

1.2 數(shù)學(xué)模型

定性和定量分析儲糧生態(tài)系統(tǒng)內(nèi)的各種生物和非生物因素之間的相互作用是研究儲糧生態(tài)系統(tǒng)的關(guān)鍵。數(shù)學(xué)模型是定性和定量分析必不可少的工具?；趧恿?，質(zhì)量和能量方守恒方程以及多孔介質(zhì)熱濕耦合傳遞理論，構(gòu)建了淺圓倉儲糧通風(fēng)過程中物理和生物參數(shù)變化的數(shù)學(xué)模型。通過坐標(biāo)轉(zhuǎn)換將錐形底筒倉映射到柱坐標(biāo)下的空間變量(ξ,φ,η)中，其中ξ代表半徑方向，φ代表角度方向，η代表高度方向。

1.2.1 流動方程

圓柱坐標(biāo)下通風(fēng)時糧倉內(nèi)部流動如式(1)表示：

(1)

式中：Vr、Vφ、Vz分別為半徑方向、圓周方向和高度方向的速度場；p為糧堆的壓力/Pa；ξ為半徑方向；φ為角度方向；η為高度方向；α、β為假設(shè)的瞬態(tài)因子；α1、α2、α3、α4為等比的經(jīng)驗(yàn)常數(shù)；RP為瞬變因子。

1.2.2 能量守恒方程

根據(jù)局部熱平衡理論，糧堆內(nèi)部的能量守恒方程可由式(2)表示：

(2)

式中：T為糧堆中溫度分布/K；W是糧堆中水分的含量/kg；H為空氣中的絕對濕度/kg/kg；Keff為有效熱擴(kuò)散量/J/kg·k；Hw為糧堆吸濕而放出的熱量/J/kg；ρb、ρa(bǔ)分別為糧堆、干空氣的密度/kg/m3；hvap為自由水放出的潛熱/J/kg；hs為糧堆吸收的熱量/J/kg；Qr為糧堆累積耗氧量/kg。ca、(c1)r、(c1)σ和(c2)σ分別為干空氣、水蒸氣、液態(tài)水和糧食的比熱容/J/kg·K；c2v為潛熱方程中的特定常數(shù)/J/kg·K；εr為糧堆的孔隙率；m為干物質(zhì)質(zhì)量/kg。

1.2.3 水分守恒方程

根據(jù)質(zhì)量守恒原理，糧堆內(nèi)部的質(zhì)量守恒方程可由式(3)表示：

(3)

式中：Deff為有效擴(kuò)散率/m2/s；dm/dt為干物質(zhì)損耗率。

1.2.4 殺蟲劑濃度衰減

在儲糧通風(fēng)中加入殺蟲劑，殺蟲劑會附著在糧粒上從而抵抗害蟲造成的損耗。由于糧食顆粒屬于吸濕性多孔介質(zhì)，當(dāng)殺蟲劑作用于糧粒時會在較短時間內(nèi)產(chǎn)生分解，致使自身藥效減小。殺蟲劑濃度衰減過程由式(4)表示[3]。

Cpes=Cie(-1.386 3(rh)time×10B(T-30)/t*)

(4)

式中：Ci為殺蟲劑的初始濃度；B和t*為殺蟲劑的特性參數(shù)，對于谷蠹來說，分別是0.036和8.467E06。rh和T分別為糧堆內(nèi)部的溫濕度，從方程(4)可以看出殺蟲劑濃度的衰減是與糧堆內(nèi)部的溫濕度相關(guān)的，因此，可以看出糧堆內(nèi)部溫度、濕度或水分是殺蟲劑濃度(單向)耦合的。

1.2.5 儲糧害蟲生長模型

儲糧害蟲的生長及分布與谷物的含水量和溫度有關(guān)。Thorpe[3]研究了昆蟲種群每周的增長速率。谷蠹生長率見式(5)。

(5)

式中：Rs為谷蠹的生長率；Cr、tr分別為谷蠹的經(jīng)驗(yàn)常數(shù)，其值分別為0.043 5、13.0，Tweb為濕球溫度/K。

從熱力學(xué)角度，濕球溫度是與干球溫度和空氣濕度相關(guān)的，即糧堆溫濕度是害蟲生長密切相關(guān)的、相互耦合的。

1.3 數(shù)學(xué)模型的驗(yàn)證

1.3.1 儲糧通風(fēng)實(shí)驗(yàn)

實(shí)驗(yàn)糧種為小麥，通風(fēng)時間為每天0∶00—6∶00，共通風(fēng)10 d。實(shí)驗(yàn)送風(fēng)溫度為20 ℃，相對濕度55.5%。在通風(fēng)期間，采集設(shè)置在糧堆內(nèi)部的溫濕度傳感器檢測糧溫的數(shù)據(jù)，通風(fēng)過程水分采用埋袋取樣進(jìn)行水分化驗(yàn)。同時利用程序進(jìn)行通風(fēng)工況的模擬，然后比較兩者結(jié)果，驗(yàn)證模型的準(zhǔn)確性。本次實(shí)驗(yàn)共設(shè)置11個測點(diǎn)，Z1、Z2、Z3布置在筒倉軸心處，在東南、東北、西南、西北不同高度各設(shè)有2個測點(diǎn)，NE1、SE1、NW1、SW1設(shè)置在下部；NE2、SE2、NW2、SW2設(shè)置在上部，具體傳感器布置見圖1b。糧庫內(nèi)配備一套精確的糧情監(jiān)測軟件，可以將某時刻測點(diǎn)的溫濕度輸入電腦端。本次實(shí)驗(yàn)每隔4 h采集一次溫度，但由于糧食水分含量變化緩慢，所以每天取樣一次測定水分含量。

1.3.2 模擬的初始條件及程序運(yùn)行

模擬狀態(tài)下稻谷糧堆溫度為40 ℃，初始濕基狀態(tài)下含水量為12%，送風(fēng)溫度為20 ℃，相對濕度55.5%，噸糧通風(fēng)量為49.68 m3/h·t，均與實(shí)驗(yàn)保持一致。首先對三維筒倉進(jìn)行網(wǎng)格劃分，在ξ方向，筒倉間距是不均勻的，在其他兩個方向上，網(wǎng)格間距等分。最終選用知21×20×21的網(wǎng)格劃分方案，可以在滿足精確度的同時減少不必要的運(yùn)算量。模擬中通風(fēng)方案采用間歇式，每天通風(fēng)6 h共通風(fēng)10 d，通風(fēng)狀態(tài)下循環(huán)時間步長為90 s，非通風(fēng)狀態(tài)下循環(huán)時間步長為900 s。將筒倉尺寸、糧食的各項(xiàng)物性參數(shù)以及害蟲增長的經(jīng)驗(yàn)常數(shù)等定義為全局變量，設(shè)立在子程序中，主程序運(yùn)行時調(diào)用子程序即可運(yùn)行。

1.3.3 數(shù)值模擬與實(shí)驗(yàn)結(jié)果比較

將實(shí)驗(yàn)實(shí)測的溫度數(shù)據(jù)與模擬數(shù)據(jù)導(dǎo)入Origin軟件進(jìn)行處理后，選取其中三個測點(diǎn)呈現(xiàn)結(jié)果見圖2。

圖2 部分測點(diǎn)溫度的實(shí)驗(yàn)數(shù)據(jù)與模擬數(shù)據(jù)對比

模擬水分?jǐn)?shù)據(jù)與實(shí)測水分?jǐn)?shù)據(jù)的對比見表1。

表1 模擬水分?jǐn)?shù)據(jù)與實(shí)測水分?jǐn)?shù)據(jù)的對比

綜合圖2及表1可知，探頭所測出的糧堆溫度與模擬結(jié)果近乎相似，水分的實(shí)測數(shù)據(jù)與模擬數(shù)據(jù)誤差不超過0.5%，相對較小，驗(yàn)證了數(shù)學(xué)模型和這套程序的準(zhǔn)確性。

2 不同通風(fēng)工況的模擬結(jié)果與對比分析

2.1 不同通風(fēng)工況的初始條件

為了進(jìn)一步優(yōu)化儲糧通風(fēng)工藝，將筒倉高度設(shè)為2.7 m，直徑為2.23 m。通風(fēng)道采用十字型加環(huán)型組合式通風(fēng)道，其中十字型通風(fēng)口長度及環(huán)型通風(fēng)口半徑均為0.8 m，優(yōu)化后的儲糧倉俯視圖如圖3所示。模擬初始溫度為30 ℃，初始濕基狀態(tài)下水分含量為12%，噸糧通風(fēng)量為5.65 m3/h·t，其他參數(shù)設(shè)置同實(shí)驗(yàn)工況。分別模擬工況一即每天通風(fēng)10 h、工況二即全天24 h通風(fēng)兩種工況，通風(fēng)總時間均為10 d，對比分析溫度場、水分含量、儲糧害蟲的分布以及殺蟲劑濃度衰減的結(jié)果，為以后糧食的安全儲藏和通風(fēng)設(shè)計提供參考。

注：單位：mm。圖3 優(yōu)化后的筒倉俯視圖

2.2 不同工況下糧堆溫度及水分含量的變化與分析

自然儲藏時，糧倉內(nèi)糧堆的溫度變化較緩慢，但對糧倉進(jìn)行通風(fēng)時，溫度會發(fā)生較為明顯的變化，以筒倉半徑為橫坐標(biāo)，筒倉高度為縱坐標(biāo)繪制不同工況的糧堆溫度分布如圖4a、圖4c所示。糧食初溫為30 ℃，兩種工況下均通風(fēng)10 d后，通風(fēng)道附近最低溫度降至13.4 ℃，相比初溫將低16.6 ℃，下降較多，說明通風(fēng)降溫效果較好。東西兩側(cè)壁面的溫度分布并不對稱，這是因?yàn)樵诓煌较虻耐矀}壁面所接受到的太陽輻射不同，東、北、西和南四個方向的壁面溫度徑向分布如圖4e、圖4f所示?？梢钥闯鲇捎诮邮芰烁邚?qiáng)度的輻射，南向壁面溫度最高，兩種工況下壁面溫度徑向分布規(guī)律相似，但由于工況二是連續(xù)通風(fēng)，距筒倉中心較近處，壁面溫度下降幅度高于工況一。在筒倉中心區(qū)域，溫度從下向上逐漸增大，平均溫度約為17 ℃。對比圖4a與圖4c可知工況二條件下糧堆溫度分布更加均勻，在筒倉中上部區(qū)域溫度略高于工況一，下部區(qū)域低于工況一，但相差不大。且底部通風(fēng)籠附近溫度均升高，可能與糧粒呼吸作用產(chǎn)熱有關(guān)。

由圖4b、圖4d可知，糧堆的水分含量近似對稱分布，但工況二條件下水分含量整體低于工況一，在錐形倉底部，工況二糧堆水分含量降至11.3%，比初始含水量低約0.7%，比工況一低約1%。同時，由于種子的呼吸作用的影響，底部通風(fēng)道周圍糧堆水分含量升高。由于組合式通風(fēng)道通風(fēng)面積較大，在整個通風(fēng)過程中，水分整體變化較小，說明該風(fēng)道形式通風(fēng)保水效果較好。

圖4 不同工況糧堆溫度與水分含量分布等高線圖

2.3 不同工況下儲糧害蟲的分布與對比分析

影響儲糧害蟲生長繁殖的因素很多，包括溫度、濕度、光線等環(huán)境因素和糧粒本身的生物因素，各環(huán)境因素共同作用對儲糧生物場造成綜合性影響，且儲糧害蟲正常的生長、繁殖必須依賴一定范圍的溫度、濕度條件。對儲糧倉進(jìn)行通風(fēng)能有效控制昆蟲種群生長，本次主要模擬了谷蠹在糧倉內(nèi)的增長量及其分布，如圖5所示。對比圖5a、圖5b可知，兩種工況下谷蠹在筒倉內(nèi)的分布規(guī)律相似，筒倉中心區(qū)域均有分層現(xiàn)象，但工況二條件下谷蠹分層更明顯。由于谷蠹具有喜溫性，所以主要分布在壁面處及上層區(qū)域。在筒倉中心區(qū)域，工況二條件下同高度處谷蠹的增長量高于工況一，但相差不是很大。當(dāng)筒倉高度小于1 m時，工況二的谷蠹增長量小于工況一，這是由于連續(xù)通風(fēng)導(dǎo)致該區(qū)域糧堆溫度及水分含量都低于工況一。

圖5 不同工況下谷蠹數(shù)量分布圖

2.4 不同工況下殺蟲劑濃度衰減的變化與分析

對糧倉通風(fēng)降溫能夠有效緩解儲糧害蟲的生長，但簡單的機(jī)械通風(fēng)很難殺死害蟲，防治儲糧害蟲的最直接有效的方法便是使用殺蟲劑。圖6a、圖6b分別是兩種工況下殺蟲劑濃度衰減等高線圖，可以看出在通風(fēng)道附近，靠近風(fēng)口位置，殺蟲劑濃度較高，最高可達(dá)97.6%和97.7%。在壁面附近及筒倉中心區(qū)域，從下向上，殺蟲劑濃度逐漸降低，這是由于筒倉中心區(qū)域的氣流溫度較高，影響殺蟲劑的降解，所以在這個區(qū)域的殺蟲劑濃度較低，且工況二條件下殺蟲劑濃度衰減分層更明顯，整體濃度水平低于工況一。另外，太陽輻射也會影響殺蟲劑濃度的衰減，由圖6可知，壁面處殺蟲劑濃度下降較快，這是由于太陽輻射主要作用于筒倉壁面，致使壁面溫度較高，加速了殺蟲劑的降解。

圖6 不同工況下殺蟲劑濃度衰減等高線圖

3 結(jié)論

本研究基于多孔介質(zhì)熱濕耦合傳遞理論，建立了淺圓倉的糧堆內(nèi)部熱濕傳遞和流動的數(shù)學(xué)模型，采用Fortran語言編程模擬了儲糧通風(fēng)后糧堆內(nèi)部溫度、水分含量、殺蟲劑的濃度衰減規(guī)律以及谷蠹的分布情況。

基于本文建立的儲糧通風(fēng)數(shù)學(xué)模型預(yù)測的溫度和水分值與實(shí)驗(yàn)測定數(shù)據(jù)最大相差分別為1.1 ℃與0.5%，初步證明數(shù)學(xué)模型的可靠性。由于受到太陽輻射的影響，淺圓倉內(nèi)部糧堆的水分含量近似對稱分布，糧倉不同方向壁面的溫度分布并不對稱，朝陽一面貼壁處糧溫高于被陰面。儲糧害蟲在糧倉內(nèi)的數(shù)量分布與溫度、水分等因素有關(guān)。儲糧害蟲在壁面處分布較多，在筒倉中心區(qū)域，害蟲數(shù)量分布出現(xiàn)分層現(xiàn)象。殺蟲劑濃度衰減受溫度的影響，溫度高會影響殺蟲劑的降解，導(dǎo)致殺蟲劑濃度較低。