陳清華蘇國用孫美華姜闊勝劉 萍

(1. 安徽理工大學(xué)機械工程學(xué)院，安徽 淮南 232001；2. 智能礦山技術(shù)與裝備安徽省重點實驗室，安徽 淮南 232001)

中國是農(nóng)業(yè)和糧食消費大國，糧食儲存對保障國家安全和社會穩(wěn)定意義重大，如何實現(xiàn)糧食的安全儲藏一直是研究熱點之一[1-2]。研究[3]表明冷卻干燥通風(fēng)可以控制糧堆溫濕度，進而實現(xiàn)糧食的長期安全存放，而糧食作為典型的多孔介質(zhì)，干燥通風(fēng)過程中其內(nèi)部的傳熱傳質(zhì)規(guī)律較為復(fù)雜。導(dǎo)熱系數(shù)、熱擴散率等作為糧食的重要熱物性參數(shù)，由于可直接表征糧食傳遞熱量的能力[4-5]，其準(zhǔn)確獲取對于指導(dǎo)和優(yōu)化糧食干燥和安全儲藏工藝至關(guān)重要。由于多孔介質(zhì)傳熱機理較為復(fù)雜，試驗測試仍是獲取其熱物性參數(shù)的主要手段[6]。例如谷和平等[7]和龔紅菊等[8]分別運用正規(guī)狀況法及Dickerson圓桶瞬態(tài)熱流法，對顆粒狀堆積物、稻谷的導(dǎo)熱系數(shù)和熱擴散率進行了試驗測試。但現(xiàn)有測試方法中，較少的對試驗過程中糧堆內(nèi)部熱量的側(cè)向流動進行考慮，從而使測試精度受到一定的影響。文獻[9]提出基于非穩(wěn)態(tài)導(dǎo)熱乘積法測顆粒狀材料熱擴散系數(shù)，優(yōu)點是消除了側(cè)向熱流誤差帶來的影響，因而更適用于松散物料熱物性測試。在此基礎(chǔ)上，岳高偉等[10]利用高低溫交變濕熱試驗箱設(shè)定溫度邊界(-50～100 ℃)，將應(yīng)用范圍擴展至低溫環(huán)境。由于非穩(wěn)態(tài)導(dǎo)熱乘積法模型只能直接測算熱擴散系數(shù)a，如需獲取比熱容Cp和導(dǎo)熱系數(shù)λ，仍需借助磁力攪拌水卡計[11]等測試手段，大大限制了其推廣應(yīng)用性。通過易于獲取的數(shù)據(jù)信息結(jié)合明確的數(shù)學(xué)模型，反演或估計被研究對象的多個參數(shù)被稱為多宗量傳熱反問題[12]，目前已在多個領(lǐng)域得到應(yīng)用并被認(rèn)可。本試驗擬基于多維非穩(wěn)態(tài)導(dǎo)熱乘積法數(shù)學(xué)模型，并結(jié)合參數(shù)估計法[13-14]對糧食顆粒的導(dǎo)熱系數(shù)λ、比熱容Cp等進行反演計算，本測試方法試驗操作難度較低，便于獲取高效準(zhǔn)確的糧粒熱物性參數(shù)，進而為制定糧食儲藏技術(shù)方案提供科學(xué)依據(jù)。

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

在某一刻τ將一直徑為D，高為δ，上端面絕熱，初始溫度為t0的短圓柱體，放入溫度為tw的恒溫環(huán)境中。設(shè)其溫度分布與角度無關(guān)，則其導(dǎo)熱微分方程見式(1)[15]。

(1)

式中：

T——溫度，℃；

r——徑向變量，m；

x——軸向變量，m；

α——熱擴散系數(shù)，m2/s。

顯然該問題為二維非穩(wěn)態(tài)熱傳導(dǎo)問題，一般采用數(shù)值解法，但其步驟復(fù)雜且不易編程實現(xiàn)。介于此，試驗中引用多維非穩(wěn)態(tài)導(dǎo)熱問題的乘積解法[15]。如圖1所示，一短圓柱體由直徑為r的無限長圓柱和一塊厚度為x的無限大平板垂直相貫而成，則有：

T(R,x,τ)=[T(x,τ)]×[T(r,τ)]。

(2)

式(1)的解析解：

(3)

式(3)中包含無窮級數(shù)，顯然難以直接求解，根據(jù)試算取級數(shù)前6項即可滿足精度要求，式(3)可改寫為：

(4)

式中：

θ(r,x,τ)——短圓柱過余溫度，℃；

θ(x,τ)——無限大平板過余溫度，℃；

θ(r,τ)——無限長圓柱過余溫度，℃；

θ0——初始過余溫度，θ=t-tw，其中t為任意時刻試樣溫度，℃；

tw——恒溫箱內(nèi)部溫度，℃；

L——無限大平板厚度，m；

R——無限長圓柱半徑，m；

μn——特征值；

F0——傅里葉數(shù)，F(xiàn)0=ατ/y2，其中α為熱擴散系數(shù)(m2/s)；y為特征長度；對于無限大平板y=L，對于無限長圓柱體y=R。

根據(jù)熱擴散率公式a=λ/(ρCp)，顯然如果已知試樣密度ρ、導(dǎo)熱系數(shù)λ、比熱容Cp及厚度δ等物理參數(shù)，可計算得到短圓柱體內(nèi)任意一點x處，τ時刻的無量綱過余溫度θ(r,x,τ)，此為傳熱正問題。

圖1 非穩(wěn)態(tài)導(dǎo)熱乘積法原理圖

2 熱物性參數(shù)估計靈敏度分析

為研究試樣內(nèi)任意一點處的溫度Y(τ,η)變化規(guī)律，某一時刻τ將試樣放入溫度為tw的環(huán)境中，溫度擾動記為u(τ)?？紤]試樣的導(dǎo)熱系數(shù)λ、比熱容Cp等熱物性參數(shù)，以及環(huán)境溫度tw、時間τ、空間位置等因素都會影響u(τ)的值，故將這些參數(shù)構(gòu)成一個向量：

η=(η1,η2,η3,……,ηm)T，

(5)

式中：

m——參數(shù)個數(shù)。

顯然，在試驗過程中對于一個確定的溫度測點，式(5)中未知參數(shù)η1，η2為試樣的λ和Cp。取試樣中同一高度處，距中心軸均為r的k個離散測點處的溫度測量值記為Yi(τ,η)(i=1,……,k)，然后給出一組λ和Cp的估計值，根據(jù)式(4)計算出的過余溫度記為θ(τ,η)，各離散測點i的過余溫度記為θi(τ,η)(i=1，……，k)，并進行對比：

(6)

由于測量誤差以及參數(shù)估計具有隨機性，Yi(τ,η)與θi(τ,η)之間必然存在偏差，即ε(η)>0。為此，參數(shù)估計的目標(biāo)是基于式(4)迭代，使偏差ε(η)→min。參數(shù)向量η中的λ和Cp值要能同時以足夠的精度估計出來，必須滿足在最小二乘估計意義下，在參數(shù)估計的時間區(qū)間內(nèi)，參數(shù)的靈敏度線性無關(guān)。本試驗中靈敏度是指θ(τ,η)對參數(shù)η的一階偏導(dǎo)。同時，考慮到方程式的復(fù)雜性，直接求偏導(dǎo)數(shù)較為困難，采用二階中心差商進行計算：

Xij=[θ(τi,η1,……,ηj+△ηj,……,ηm)-θ(τi,η1,……,ηj-△ηj,……,ηm)]/(2△ηj)，

(7)

式中：

Xij——τi時刻θ(τi,η)對參數(shù)ηj的靈敏度系數(shù)。

通常取△ηj=0.000 1ηj。

只有在λ和Cp線性無關(guān)時，2個參數(shù)才能同時估計，因此在估計之前進行靈敏度相關(guān)性分析。分別選取玉米和稻谷進行研究。根據(jù)參考文獻[16]和[17]分別選定玉米的相關(guān)參數(shù)為：密度750 kg/m3，導(dǎo)熱系0.176 5 W/(m·℃)，比熱容1 500 J/(kg·℃)，含水率17%；稻谷參數(shù)為：密度1 480 kg/m3，導(dǎo)熱系數(shù)0.099 5 W/(m·℃)，比熱容911.9 J/(kg·℃)。

由圖2、3可以看出，玉米和稻谷的導(dǎo)熱系數(shù)和比熱容參數(shù)靈敏度線性均不相關(guān)，但比熱容Cp靈敏度均非常小(僅為10-3數(shù)量級)，若直接進行參數(shù)估計，反演得到的比熱容與實際相差較大[14]。

圖2 玉米試樣導(dǎo)熱系數(shù)及比熱容靈敏度值

圖3 稻谷試樣導(dǎo)熱系數(shù)及比熱容靈敏度值

3 關(guān)鍵技術(shù)

基于多維非穩(wěn)態(tài)導(dǎo)熱乘積法原理，設(shè)計測糧食熱物性參數(shù)的測試系統(tǒng)，需要解決的關(guān)鍵技術(shù)有測試系統(tǒng)構(gòu)建、多維非穩(wěn)態(tài)傳熱條件、參數(shù)估計算法等。

3.1 測試系統(tǒng)設(shè)計

如圖4所示，系統(tǒng)主要由恒溫箱、信號采集與傳輸模塊、信號發(fā)射端和接收端、試樣筒、溫度傳感器和微型計算機等構(gòu)成。為了避免在恒溫箱上開孔影響恒溫箱效果，以及增加試驗操作難度，系統(tǒng)采用深圳蜂匯公司型號為Z-0004的ZigBee無線數(shù)據(jù)采集模塊實現(xiàn)信號采集與傳輸。恒溫箱溫度范圍為0～100 ℃，控溫精度±0.1 ℃；測溫?zé)犭娕季取?.2 ℃；ZigBee數(shù)據(jù)采集模塊的射頻芯片CC2530的接收靈敏度為-97 dBm，有效保證了溫度的實時傳輸。系統(tǒng)軟件通過LabView編程建立基本松散物體參數(shù)輸入?yún)^(qū)、采集溫度變化顯示區(qū)、控制按鈕區(qū)和計算結(jié)果區(qū)等，便于對測試結(jié)果的觀察以及對不同糧食顆粒測算參數(shù)設(shè)置。

1. 采集信號接收端 2. 微型計算機 3. 信號采集與傳輸單元 4. 絕熱端蓋 5. 恒溫箱 6. 試樣 7. 試樣筒體 8. 溫度測點 9. 采集信號發(fā)射端

圖4 測試裝置簡圖

Figure 4 Diagram of testing system

3.2 非穩(wěn)態(tài)導(dǎo)熱乘積物理邊界的實現(xiàn)

根據(jù)非穩(wěn)態(tài)導(dǎo)熱乘積法原理及物理模型要求，試樣為圓柱形，同時上表面絕熱，下表面和圓柱面為恒溫邊界。此處利用恒溫箱實現(xiàn)恒溫環(huán)境，并設(shè)計圓柱形試樣筒盛放試樣，為保證試樣放入恒溫箱后較短時間內(nèi)與環(huán)境溫度保持一致，筒體采用黃銅制作，在內(nèi)壁設(shè)置溫度測點，待溫度與環(huán)境溫度一樣時開始試驗測試，以此實現(xiàn)圓柱面和下表面的恒溫邊界條件。同時，試樣筒上端蓋為絕熱棉制作，形成試樣上絕熱表面。

3.3 比熱容Cp靈敏度較低問題的處理

針對試驗選用的材料Cp靈敏度較低，比熱容不可直接估算，采用隨機共軛梯度法進行求解[18]，由于式(4)中的F0=aτ/y2，只有一個未知數(shù)a，因此，首先反演估算a，然后利用a對估計的Cp值進行修正，具體流程(見圖5)：

① 反演估計a；

② 估計材料λ及Cp的原始估計值；

③ 根據(jù)式(4)、(5)、(6)求得滿足精度要求的a及λ；

④ 根據(jù)a估計值對原始估計值Cp進行修正，得到最終值Cp*。

為避免在測試過程中，因外界干擾因素造成的個別數(shù)據(jù)波動對估算結(jié)果造成影響，本試驗首先計算出每個特定時刻τ對應(yīng)的一組數(shù)值λ和Cp*，然后對計算出的一系列物性參數(shù)求加權(quán)平均作為最終參數(shù)的估計值。

圖5 參數(shù)估計流程圖

3.4 溫度測試準(zhǔn)確性及可靠性

當(dāng)材料質(zhì)地均勻，滿足各項同性的性質(zhì)時，溫度測點的布置對物性參數(shù)測算結(jié)果的影響不大。因此，試驗前對稻谷顆粒進行篩選，裝入箱體時進行壓實處理。同時，在離圓柱側(cè)面和底面相同距離的不同位置處，安裝多個溫度傳感器，并取加權(quán)平均作為最終測試結(jié)果。此外，系統(tǒng)采取無線信號發(fā)射與接收的方式，有效保證了試驗箱體內(nèi)部為恒溫邊界條件，減少外界環(huán)境溫度對測試結(jié)果的影響。

4 熱物性測試試驗分析

4.1 糧食顆粒熱物性參數(shù)測試

選取淮南地區(qū)含水率為10.5%的皖稻121。首先利用孔徑為2 mm的篩網(wǎng)對稻谷進行篩選，然后裝入直徑為400 mm，高為80 mm的短圓柱筒試樣盒中，填充密度580 kg/m3，環(huán)境溫度20 ℃，恒溫箱溫度保持90 ℃。在離試樣筒壁面和底面均為40 mm的位置處布置2個溫度測點，取平均值作為最終測試值，依次編號為測點1和測點2。某一時刻τ突然將圓柱體試樣盒放入恒溫箱內(nèi)，同時觀察銅板內(nèi)壁面溫度變化，當(dāng)銅板內(nèi)壁溫度穩(wěn)定且趨近于恒溫箱體內(nèi)部溫度時，開始采集稻谷內(nèi)部測點溫度。為消除試驗初期，測試溫度的不穩(wěn)定對測算結(jié)果的影響，選取100 s以后的采集數(shù)據(jù)為分析對象，試驗共進行3 600 s，溫度采樣間隔10 s，各測點溫度曲線見圖6。

熱擴散率α初始猜測值為2.02×10-7m2/s，導(dǎo)熱系數(shù)λ和Cp初始猜測值分別為0.099 0 W/(m·℃)和823 J/(kg·℃)，取參數(shù)估計結(jié)果見表1。

圖6 測點1、測點2的溫升曲線

可以看出熱擴散率和導(dǎo)熱系數(shù)估計結(jié)果與文獻[17]中的測算結(jié)果較為接近，但比熱容誤差較大(>10%)，而修正后的比熱容則較為吻合。為進一步驗證測算結(jié)果準(zhǔn)確性，對同一種試樣在相同條件下進行多次試驗，通過相對偏差△η衡量測算值彼此接近的程度。

(8)

式中：

由表2可以看出，導(dǎo)熱系數(shù)與熱擴散率參數(shù)反演結(jié)果以及修正后的比熱容Cp*與對應(yīng)平均值的最大相對偏差均小于7%，即試驗滿足可重復(fù)的要求。

表1 參數(shù)估計結(jié)果

表2 參數(shù)估計試驗可重復(fù)精度分析

4.2 理論計算與數(shù)值仿真結(jié)果對比

將表2中獲取的熱物性參數(shù)的平均值代入式(4)中，計算τ在100～3 600 s時的溫度理論值，同時在Fluent中仿真模擬稻谷內(nèi)部的溫度場變化。數(shù)據(jù)曲線對比見圖7。

由圖7可以看出，稻谷的仿真曲線、理論計算曲線與實測溫升曲線三者較為吻合，進一步驗證了參數(shù)估計的準(zhǔn)確性。值得注意的是，試驗至約2 600 s后，實測溫升與仿真和理論計算溫升開始出現(xiàn)偏差，且隨時間的延長，偏差具有變大趨勢。原因可能是熱量傳遞至試驗箱上端絕熱邊界后因無法及時轉(zhuǎn)移，造成熱量積聚效應(yīng)，進而對測點溫度產(chǎn)生影響，使實測溫升偏大，且隨著熱量積聚越來越多，影響也就越來越明顯。從而為盡量消除此因素的影響，選取100～2 500 s 時的溫升為最終有效測算數(shù)據(jù)。

4.3 3種糧食顆粒熱物性參數(shù)測試分析

為進一步驗證測試系統(tǒng)的穩(wěn)定性與準(zhǔn)確性，分別選取材料Y兩優(yōu)900水稻、聯(lián)創(chuàng)11號和源育15玉米作為試樣，含水率分別為13.1%，5.0%，12.5%，填充密度分別為600，720， 800 kg/m3。

圖7 實測、理論與數(shù)值仿真溫升曲線對比

由表3數(shù)據(jù)可以看出，測得的熱導(dǎo)率λ與文獻[16]和[17]數(shù)據(jù)較為接近，其中“源育15”號玉米的導(dǎo)熱系數(shù)與文獻數(shù)據(jù)相對誤差最大，為4.64%，分析原因可能是不同品種的玉米本身的熱物性參數(shù)存在差異，且玉米粒與稻谷相比顆粒直徑較大，存入試驗箱時，玉米顆粒之間的孔隙較大，測溫過程溫度波動較大，最終導(dǎo)致測算結(jié)果與參考值有偏差，但仍小于5%。符合測試精度要求。

表3 3種糧食試樣測試結(jié)果綜合分析

5 結(jié)論

(1) 采用恒溫箱結(jié)合黃銅短圓柱試樣盒構(gòu)建的恒溫邊界模型便于控制實現(xiàn)，且運用非穩(wěn)態(tài)導(dǎo)熱乘積法結(jié)合參數(shù)估計法對糧食顆粒的熱物性參數(shù)進行測試在原理與技術(shù)上均是可行的。

(2) 采用ZigBee無線數(shù)據(jù)采集方式，數(shù)據(jù)采集更方便，且不受試驗箱體結(jié)構(gòu)的限制，降低了試驗操作的難度。

(3) 糧粒的比熱容Cp的靈敏度系數(shù)較小，需先估計再修正，才能得到較為準(zhǔn)確的估計值。利用估算的熱物性參數(shù)值通過理論計算與數(shù)值仿真溫升變化曲線與實測溫度變化曲

線一致，進一步說明了參數(shù)估計的準(zhǔn)確性；對包括皖稻121在內(nèi)的4種糧食進行熱物性測算，結(jié)果與相關(guān)文獻吻合，能夠滿足實際應(yīng)用要求。