徐慧星

(江西環(huán)境工程職業(yè)學(xué)院，江西 贛州 341000)

0 引言

自人類(lèi)社會(huì)開(kāi)始生產(chǎn)糧食以來(lái)，就有對(duì)谷物進(jìn)行儲(chǔ)藏的習(xí)慣。據(jù)國(guó)家糧食局統(tǒng)計(jì)測(cè)算，我國(guó)糧食產(chǎn)后僅儲(chǔ)藏、運(yùn)輸、加工等環(huán)節(jié)損失浪費(fèi)量達(dá)350億kg以上。其中，農(nóng)戶儲(chǔ)糧損失比例高達(dá)8%左右，每年因蟲(chóng)霉鼠雀造成損失200億kg以上。由于儲(chǔ)存條件差、損失大，損失損耗逾75億kg。我國(guó)糧食供求緊張，近年來(lái)糧食進(jìn)口量逐年增加，每年進(jìn)口的谷物和大豆約500億kg，而每年又白白地?fù)p失浪費(fèi)上千億斤糧食。

目前，糧倉(cāng)大多采用溫度傳感器網(wǎng)絡(luò)，對(duì)其內(nèi)部溫度場(chǎng)進(jìn)行實(shí)時(shí)監(jiān)測(cè)，大大增加了儲(chǔ)存費(fèi)用，而且隨著年限的增加，安全隱患比較大，容易導(dǎo)致巨大損失。因此，根據(jù)傳熱學(xué)理論，采用有限元法對(duì)糧倉(cāng)瞬時(shí)溫度場(chǎng)進(jìn)行數(shù)學(xué)模擬分析，不僅可減少溫度傳感器網(wǎng)絡(luò)監(jiān)測(cè)和維護(hù)費(fèi)用，對(duì)防止糧食霉變、節(jié)能減排具有十分重要的現(xiàn)實(shí)意義。

1 傳熱學(xué)有限元方程概述

由傳熱學(xué)原理可知，在一定溫差下，物體內(nèi)熱量會(huì)從高溫地向低溫地傳遞。在正常情況下，熱傳遞中，物體內(nèi)溫度T和位置、時(shí)間成函數(shù)關(guān)系，即

T=T(x,y,z,t)

(1)

物體內(nèi)某時(shí)刻各點(diǎn)溫度值組成的集合稱(chēng)為溫度場(chǎng)。由于溫度同性的連續(xù)介質(zhì)中，溫度場(chǎng)與邊界、初始二者條件方不具方向性，因此溫度場(chǎng)為標(biāo)量場(chǎng)。

各向方程為

(2)

在指定物體溫度邊界ST上的控制方程為

(3)

在指定物體熱流邊界Sq上的條件方程為

(4)

在指定物體對(duì)流邊界Sq上的條件方程為

(5)

初始條件為

T=f(x,y,z,t0)

(6)

(7)

取式(7)一階偏微分為0，并用矩陣表示，則有

(8)

采用有限元計(jì)算方式，可以得到式(8)的溫度場(chǎng)有限元基本方程為

([Kk]+[Kc]){T}=
{Rb}-[C]{T}+{Rq}+[Kc]{Tc}

(9)

其中，[Kk]為物體全部熱傳導(dǎo)的矩陣；{T}為未知節(jié)點(diǎn)溫度的矩陣；[Rb]為與物體內(nèi)部本身熱源等效的節(jié)點(diǎn)熱流矩陣；[C]為物體全部熱容矩陣；{T}為總體節(jié)點(diǎn)溫度的熱流矩陣；{Rq}為指定物體熱流邊界產(chǎn)生的等效節(jié)點(diǎn)熱流矩陣；[Kc]為總體對(duì)流矩陣；{Te}為指物體熱流邊界的節(jié)點(diǎn)環(huán)境溫度矩陣。

2 確定糧食熱物理性參數(shù)

模擬糧食儲(chǔ)藏倉(cāng)的溫度場(chǎng)數(shù)值變化，首先需要了解糧食導(dǎo)熱系數(shù)、比熱、重力密度等熱物性參數(shù)，以上參數(shù)采用特定試驗(yàn)很容易獲取。但是，由于糧食的熱物性參數(shù)變化與環(huán)境溫濕度有關(guān)，可供參考比熱參數(shù)很不容找到，且不同地方和不同時(shí)候因季節(jié)和雨水多少的變化，糧食的比熱參數(shù)更是相差甚遠(yuǎn)。

目前，物體比熱容主要有非穩(wěn)態(tài)熱流法、微量熱法、比較量熱法及混合法等測(cè)量法，結(jié)合糧食熱物性參數(shù)特性，本文采用混合法進(jìn)行比熱參數(shù)的測(cè)量?；旌戏ㄒ话悴捎昧繜崞鳙@取物質(zhì)比熱參數(shù)，測(cè)量過(guò)程中，只需事先在量熱器中加入些比例適合的熱水，然后加入比例適合的試樣和熱水混合，混合法所用的主要儀器是量熱器。測(cè)量時(shí)，在量熱器中預(yù)先加入一定量的熱水，再將一定溫度和質(zhì)量的試樣倒入量熱器中并和熱水充分混合?；旌戏崞胶夥匠虨?/p>

C1M1ΔT1=C2M2ΔT2

(10)

其中，C1和C2分別為水和試樣的比熱容；M1和M2分別為水和試樣的質(zhì)量；T1和T2分別為水和試樣的質(zhì)量的溫度變化量。

根據(jù)式(8)可以得出試樣的比熱容為

(11)

由于在進(jìn)行混合法熱平衡過(guò)程中會(huì)有散熱現(xiàn)象出現(xiàn)，因此需要對(duì)式(9)進(jìn)行一定的改進(jìn)才可以正常使用，改進(jìn)主要包括兩部分內(nèi)容：

1)系統(tǒng)散熱情況下熱容量計(jì)算方式。系統(tǒng)散熱對(duì)熱平衡影響較大，假設(shè)測(cè)試環(huán)境溫度不變，系統(tǒng)散熱量可表示為

ΔQ=CZMZKΔt

(12)

其中，△t為整個(gè)熱平衡所需要的時(shí)間(s)；△Q為整個(gè)熱平衡熱量外散總和(J)；K為散熱系數(shù)(K/s)；CZMZ為系統(tǒng)熱容總和(J/K)。

量熱器熱容量需提前測(cè)定，具體辦法為：采用冷熱水混合法。加入質(zhì)量、溫度為M21、T21的冷水和M22、T22的熱水，記下熱平衡所需時(shí)間Tp。假設(shè)水的比熱容為Cs，則系統(tǒng)散熱的熱容量可以表示為

2)比熱容計(jì)算方式。結(jié)合文章前文內(nèi)容可以求出比熱容方程式為

(C0M0+C1M1)ΔT1=CMΔT2+(C0M0+C1M1)Kt3

(14)

因此，倉(cāng)內(nèi)糧食的比熱容可以表示為

(15)

其中，t3為加入試樣后達(dá)到熱平衡花費(fèi)的時(shí)間。

3 糧倉(cāng)溫度試驗(yàn)測(cè)試方法

在試驗(yàn)中，將溫度傳感器固定在糧倉(cāng)溫度待測(cè)點(diǎn)位置，糧倉(cāng)溫度傳感器位置圖如圖1所示。

圖1 糧倉(cāng)溫度傳感器位置圖Fig.1 The position diagram of granary temperature sensor

圖1中數(shù)值1～39為39個(gè)溫度測(cè)試點(diǎn)。試驗(yàn)測(cè)試中選用某大豆種植區(qū)大豆，大豆和糧倉(cāng)的熱物性參數(shù)如表1所示。為了盡量降低大豆的空氣相對(duì)濕度，減少試驗(yàn)誤差，在試驗(yàn)前大豆經(jīng)過(guò)曬干和烘干處理，并多次進(jìn)行翻動(dòng)，確保整個(gè)大豆含水量和溫度一致。在大豆入倉(cāng)時(shí)，保持大豆含水量和溫度的均勻一致。入倉(cāng)時(shí)，水稻的平均含水量為12.5 %，溫度為28℃。

試驗(yàn)開(kāi)始時(shí)，開(kāi)啟糧食烘干機(jī)，由烘干機(jī)風(fēng)向正對(duì)面對(duì)準(zhǔn)布置測(cè)溫點(diǎn)的半截面，然后倒入大豆，啟動(dòng)烘干機(jī)給倉(cāng)內(nèi)供給熱風(fēng)，測(cè)出糧倉(cāng)內(nèi)熱風(fēng)的風(fēng)速約為0.6m/s。在測(cè)試中，前12h烘干機(jī)一直工作，之后關(guān)閉烘干機(jī)，記錄整個(gè)過(guò)程中溫度傳感器網(wǎng)絡(luò)所測(cè)溫度值。

表1 大豆和糧倉(cāng)的熱物性參數(shù)Table 1 The thermophysical parameters of soybean and granary

4 糧食儲(chǔ)藏溫度場(chǎng)有限元模擬分析

4.1 建立有限元數(shù)學(xué)模型

為了充分利用糧倉(cāng)分布溫度傳感器網(wǎng)絡(luò)測(cè)試的實(shí)時(shí)數(shù)據(jù)，在建立數(shù)學(xué)模型時(shí)，選取溫度傳感器所在地為模型節(jié)點(diǎn)。糧倉(cāng)有限元三維模型如圖 2 所示。

圖2 糧倉(cāng)有限元三維模型Fig.2 The finite element 3D model of granary

1)確定大豆熱物性參數(shù)。在確定糧倉(cāng)溫度數(shù)學(xué)模型中，采用表1中大豆的熱物性參數(shù)為測(cè)試值，即：重力密度為1 200kg/m3，導(dǎo)熱系數(shù)為0.173W/(m·K)，比熱為2 538J/(kg·K)。

2)設(shè)定求解條件。確定糧倉(cāng)溫度數(shù)學(xué)模型，進(jìn)行瞬態(tài)溫度場(chǎng)有限元模擬，確定初始和邊界條件也是重要環(huán)節(jié)。

(1)初始條件：進(jìn)行糧倉(cāng)內(nèi)部瞬態(tài)溫度場(chǎng)模擬，首先需要定義各節(jié)點(diǎn)溫度初始值，取2015年12月16日00點(diǎn)時(shí)刻各溫度傳感器網(wǎng)絡(luò)節(jié)點(diǎn)溫度值為初始條件。

(2)邊界條件：根據(jù)有限元模擬需求，設(shè)定第一類(lèi)邊界條件來(lái)邊界條件，即按照指定邊界設(shè)定溫度值確定邊界條件，也就是說(shuō)，邊界上的節(jié)點(diǎn)溫度隨著傳感器實(shí)際測(cè)量時(shí)變化。

(3)設(shè)定迭代時(shí)間和步長(zhǎng)：總時(shí)間為1年，總共包括 8 760h，步長(zhǎng)為12h。

糧倉(cāng)溫度數(shù)學(xué)模型建模流程如圖3所示。

圖3 糧倉(cāng)溫度數(shù)學(xué)模型建模流程圖Fig.3 The flow chart of mathematical model of granary temperature

4.2 有限元數(shù)學(xué)模型結(jié)果分析

為了更好地分析糧食倉(cāng)內(nèi)溫度變化，分別給出了6 000h(2016年8月8日)、8 000h(2016年10月31日)的糧食溫度云圖，如圖4和圖5所示。

從有限元數(shù)學(xué)模型整體中取出一部分的局部區(qū)域可以觀察其內(nèi)部溫度變化情況。取模型中部區(qū)域 (稱(chēng)為局部A)，能夠比較清晰地看出糧倉(cāng)內(nèi)部溫度變化情況，局部A 6 000h(2016年8月8日) 計(jì)算結(jié)果云圖如圖6所示。

圖4 6 000h糧倉(cāng)溫度云圖Fig.4 The 6 000h granary temperature nephogram

圖5 8 000h糧倉(cāng)溫度云圖Fig.5 The 8 000h granary temperature nephogram

圖6 局部A 6 000h計(jì)算結(jié)果云圖Fig.6 The local A 6 000h calculation results nephogram

取糧堆的下部分(稱(chēng)為局部B)，可以看出糧倉(cāng)內(nèi)部某一截面的溫度變化情況，局部B 8 000h(2016年10月31日) 計(jì)算結(jié)果云圖如圖7所示。

圖7 局部B 8 000h計(jì)算結(jié)果云圖Fig.7 The local B 8 000h calculation results nephogram

從上面4幅有限元模擬模型計(jì)算結(jié)果云圖可以看出：夏天天氣溫度高，糧倉(cāng)溫度呈現(xiàn)“冷心熱皮”狀態(tài)，即糧倉(cāng)四周溫度高，中心地區(qū)溫度低；在糧倉(cāng)四周溫度變化比較明顯，而在糧倉(cāng)中心地區(qū)溫度變化不明顯，這符合大豆是不良導(dǎo)體的物理性，說(shuō)明有限元模擬模型計(jì)算是真實(shí)有效的，采用有限元法求解糧倉(cāng)溫度是可行的。

5 結(jié)論

采用有限元法對(duì)大豆糧倉(cāng)建立了有限元溫度模型，計(jì)算模擬出了糧倉(cāng)1年內(nèi)的溫度情況。結(jié)果表明：采用有限元法可以直觀地觀察糧倉(cāng)具體時(shí)間的溫度分布情況，用有限元法求解糧倉(cāng)溫度是可行的。模型計(jì)算結(jié)果表明：外界溫度變化時(shí)，糧倉(cāng)四周溫度變化比較明顯，而在糧倉(cāng)中心地區(qū)溫度變化不太明顯，這符合大豆是不良導(dǎo)體的物理性，也證明有限元模擬模型計(jì)算是真實(shí)有效的。