鐘汝能, 鄭勤紅, 姚斌, 向泰

(1.云南師范大學(xué)能源與環(huán)境科學(xué)學(xué)院， 昆明 650500； 2.云南師范大學(xué)物理與電子信息學(xué)院， 昆明 650500)

谷物(如稻米、小麥、玉米等)籽粒多為顆粒物質(zhì)，這些物質(zhì)在儲(chǔ)藏、輸送、加工過(guò)程中都涉及到顆粒堆積，在農(nóng)業(yè)機(jī)械化、電氣化的發(fā)展中，由顆粒物質(zhì)和空氣所組成的混合堆積物料宏觀上所表現(xiàn)出來(lái)的等效介電特性是研究者關(guān)注的主要參數(shù)之一。介電特性是指物質(zhì)分子中的束縛電荷對(duì)外加電場(chǎng)的響應(yīng)特性(ε=ε′-jε″)[1]，其中，介電常數(shù)(ε′)表征電磁波能量的儲(chǔ)存，介電損耗因子(ε″)表征電磁波能量的轉(zhuǎn)化。針對(duì)介電特性的基礎(chǔ)研究可為農(nóng)產(chǎn)品的微波加工、無(wú)損檢測(cè)等提供依據(jù)。

前人的研究成果較好地推進(jìn)了農(nóng)業(yè)現(xiàn)代化的發(fā)展，然而，當(dāng)前針對(duì)農(nóng)業(yè)物料電磁效應(yīng)的應(yīng)用研究與市場(chǎng)需求規(guī)模仍然存在差距，相關(guān)基礎(chǔ)研究需要進(jìn)一步深入。將傳統(tǒng)的工業(yè)顆粒填充混合物(復(fù)合材料)仿真模型拓展應(yīng)用到農(nóng)業(yè)工程領(lǐng)域面臨兩個(gè)方面的現(xiàn)實(shí)問(wèn)題：一是與工業(yè)復(fù)合材料中顆粒的隨機(jī)(或均勻)分布形式不同，谷物籽粒在儲(chǔ)藏、加工過(guò)程中以自然落料堆積的形式存在[2]，在相應(yīng)的仿真研究中需要借助離散元技術(shù)、有限元流體技術(shù)等獲得顆粒的自然堆積狀態(tài)，使模擬顆粒的空間分布與谷物真實(shí)存在狀態(tài)相吻合。二是谷物籽料多為非球體顆粒(如稻米、小麥籽粒為橢球狀，玉米籽粒為多結(jié)構(gòu)等)，如果采用傳統(tǒng)的球體模型進(jìn)行模擬仿真，則離散單元與實(shí)際物料差別較大，但是，非球體顆粒模型的本構(gòu)關(guān)系很難建立、數(shù)據(jù)結(jié)構(gòu)煩瑣、接觸判斷復(fù)雜，相關(guān)的技術(shù)方法需要研究者開發(fā)。截至目前，專門針對(duì)農(nóng)業(yè)顆粒物料等效介電特性的模擬仿真研究成果相對(duì)較少，部分學(xué)者所開展的相關(guān)研究也主要集中于堆積狀態(tài)仿真或者局域場(chǎng)模擬計(jì)算中的某一方面，覆蓋農(nóng)業(yè)散粒體的生成、堆放和介電特性計(jì)算等全周期的模擬仿真成果鮮見報(bào)道[20]。

為了節(jié)約研究成本，探究適宜分析堆積型谷物物料等效介電特性的仿真模型(簡(jiǎn)稱為堆積模型)，本研究采用離散元法模擬谷物籽料的自然堆積狀態(tài)，獲得顆?？臻g分布的位置坐標(biāo)和方向分布矩陣數(shù)據(jù)，并基于COMSOL軟件構(gòu)建了谷物籽料/空氣混合物等效介電特性的有限元法模擬模型，采用平均能量法(average energy method，AEM)對(duì)堆積模型的有效性進(jìn)行了驗(yàn)證。

1 原理與方法

1.1 混合物料等效介電特性的仿真原理

當(dāng)混合物中顆粒的尺寸遠(yuǎn)遠(yuǎn)小于電磁波波長(zhǎng)時(shí)，可用準(zhǔn)靜態(tài)條件下的平行板電容器模型分析混合物的等效介電特性[21]。如圖1所示，平行板電容器上、下極板的電勢(shì)φ1和φ2分別為1 和0 V，側(cè)面邊界條件為?φ/?n=0，極板間為顆粒填充混合物。

圖1 顆粒填充混合物模擬模型Fig.1 Simulation model of mixtures filled with particles

在有限元方法中，每個(gè)網(wǎng)格單元所擁有的靜電能如式(1)所示。

δWe(k)=

(1)

式中，k為網(wǎng)格單元的編號(hào),其取值由網(wǎng)格剖分的粗細(xì)程度決定,εk和vk分別表示第k個(gè)單元的介電特性和體積。

混合體所儲(chǔ)存的靜電能計(jì)算如下。

(2)

由電磁理論可知，圖1所示電容器中所存儲(chǔ)的靜電能可等效地表示為式(3)。

(3)

式中，(φ1-φ2)為兩個(gè)極板之間的電勢(shì)差，d為兩個(gè)極板之間的距離，S為極板的面積，ε0=8.85×10-12F·m-1，εeff為混合物等效介電特性。令We=Wc，由式(2)(3)即可獲得電容器中混合物(圖1)的等效介電特性(εeff)。

為準(zhǔn)確模擬計(jì)算各向同性混合物的介電特性，可生成同一條件下的多個(gè)模擬模型，取其所有模型計(jì)算結(jié)果的平均值作為數(shù)值結(jié)果。

(4)

1.2 谷物顆粒落料堆積的離散元模擬

1.2.1模擬顆粒建模 以稻谷籽粒為對(duì)象，采用基于離散元法的EDEM軟件進(jìn)行非球體顆粒建模：將稻谷籽粒的外形簡(jiǎn)化為橢球體,在EDEM中導(dǎo)入擬構(gòu)建顆粒的橢球體結(jié)構(gòu)輪廓，采用不同半徑的球形進(jìn)行“無(wú)限制”間接填充，直到相應(yīng)的輪廓全部填滿；導(dǎo)出全部填充球的位置坐標(biāo)和半徑坐標(biāo)，進(jìn)行格式轉(zhuǎn)換并替換為全局坐標(biāo)數(shù)據(jù)；在EDEM中加載坐標(biāo)數(shù)據(jù)，完成非球形結(jié)構(gòu)(橢球體)顆粒的構(gòu)建。

1.2.2材料參數(shù)設(shè)置 采用立方狀不銹鋼箱體作為物料裝載器件,依據(jù)稻谷籽粒和不繡鋼材料的一般物理和力學(xué)特性[2]，設(shè)置稻谷籽粒和不銹鋼容器的泊松比為0.25和0.29，彈性模量為375 和75 000 MPa，密度為1 350和8 000 kg·m-3，設(shè)定顆粒與顆粒之間的恢復(fù)系數(shù)、靜摩擦系數(shù)和滾動(dòng)摩擦系數(shù)分別為0.6、0.3和0.01；設(shè)定顆粒與容器之間的恢復(fù)系數(shù)、靜摩擦系數(shù)和滾動(dòng)摩擦系數(shù)分別為0.5、0.56和0.02。

1.2.3顆粒落料堆積 設(shè)定籽粒與容器的接觸模型為Hertz Mindlin(收獲后的稻谷籽粒含水率不高,籽粒間的黏附力可以忽略，可將其近似為理想顆粒體)。設(shè)定顆粒下落的自由落體加速度為-9.81 m·s-2，顆粒生成方式為dynamic。指定擬生成顆粒數(shù)量,采用“Simulator”功能完成顆粒的自由落料堆積仿真。仿真結(jié)束后，從“Result”中導(dǎo)出每一個(gè)橢球顆粒的位置坐標(biāo)(Xn,Yn,Zn;n為顆粒編號(hào))和方向向量矩陣數(shù)據(jù)(XXn,XYn,XZn;YYn,YXn,YZn;ZXn,ZYn,ZZn)。

1.3 谷物顆粒物料的有限元模擬模型

依據(jù)方向余弦矩陣?yán)碚?，編制Matlab語(yǔ)言程序，將從EDEM中導(dǎo)出的單個(gè)顆粒所對(duì)應(yīng)的方向向量矩陣換算為笛卡爾全局坐標(biāo)系中X、Y、Z軸對(duì)應(yīng)的單位矢量(in,jn,kn)。同時(shí)，在COMSOL軟件的“App開發(fā)器”中創(chuàng)建“Model method”程序，在“聲明”中定義不同的數(shù)組、變量分別儲(chǔ)存各顆粒的位置坐標(biāo)值、體積比和結(jié)構(gòu)參數(shù)，生成一個(gè)大的立方體代表基體，使用“錄制”方法生成橢球體顆粒的建模代碼并嵌入到主程序中，依次讀取每一個(gè)顆粒的位置坐標(biāo)(Xn,Yn,Zn)和方向向量(in,jn,kn)。執(zhí)行method程序，在指定區(qū)域生成橢球體顆粒及立方體，執(zhí)行布爾運(yùn)算后形成“聯(lián)合體”，完成模擬模型構(gòu)建。

綜合考慮網(wǎng)格剖分的收斂性及計(jì)算效率，設(shè)置網(wǎng)格剖分級(jí)別為“細(xì)化”，采用“中等”優(yōu)化級(jí)別進(jìn)行單元質(zhì)量?jī)?yōu)化。在材料屬性中依次選取“基體區(qū)域”“累積”分別對(duì)基質(zhì)、基體進(jìn)行物理屬性賦值。求解式(2)～(4)，即可得到稻谷顆粒/空氣混合物的等效介電特性數(shù)據(jù)。

1.4 谷物堆積角測(cè)量

1.4.1真實(shí)谷物堆積角測(cè)量 稻谷籽粒購(gòu)自云南呈貢農(nóng)業(yè)種子有限公司，濕基含水率為14.6%，10次測(cè)量后取其長(zhǎng)半軸、短半軸的平均值分別為3.46 和1.45 mm，平均粒徑比為2.39∶1。經(jīng)除雜、去芒和篩分挑選后得到相對(duì)均勻的1 500粒長(zhǎng)粒稻谷籽粒，將其以落料堆積在不銹鋼金屬板上[2]，得到真實(shí)稻谷堆的堆積形態(tài)后，從不同方位分多次測(cè)量稻谷堆底的直徑(d)和高度h，采用Matlab函數(shù)(式5)計(jì)算得到稻谷堆的平均堆積角。

tanθ′=2h/d

(5)

1.4.2模擬谷物堆積角測(cè)量 采用1.2節(jié)所述方法模擬生成橢球體稻谷籽粒在不銹鋼金屬板上自然落料堆積，設(shè)定顆粒尺寸、數(shù)量與實(shí)驗(yàn)樣品的平均值相一致。得到模擬稻谷堆的堆積形態(tài)后，采用Origin軟件的圖像識(shí)別得到單側(cè)稻谷堆的邊緣輪廓曲線，提取輪廓曲線數(shù)據(jù)并得到擬合直線的斜率，計(jì)算模擬谷物堆的堆積角(θ)[2]。

θ=arctan|k|/π)

(6)

1.5 谷物堆積模型及數(shù)值驗(yàn)證

1.6 實(shí)驗(yàn)測(cè)量技術(shù)及測(cè)量參數(shù)

采用矩形波導(dǎo)傳輸反射方法[22]分別測(cè)量并得到不同體積分?jǐn)?shù)條件下稻米物料、小麥物料的等效介電特性數(shù)據(jù)(24 ℃)，其中，稻米顆粒(“Lebonnet”) 的密度、含水率和測(cè)量頻率分別為1.476 g·cm-3、12.2%、11.0 GHz；小麥顆粒(“Scoutland”) 的密度、含水率和測(cè)量頻率分別為1.406 g·cm-3、10.9%、11.7 GHz。

采用無(wú)校準(zhǔn)同軸傳輸反射方法[14]分別測(cè)量并得到不同含水率(濕基)條件下小米顆粒物料的等效介電特性數(shù)據(jù)(24 ℃)，其中小米顆粒含水率分別為2.0%、7.5%、10.0%、15.0%、19.7%，體積密度分別為1.47、1.38、1.35、1.27、1.24 g·cm-3。測(cè)量頻率為2～6 GHz。

1.7 儀器設(shè)備及數(shù)據(jù)處理

采用離散元顆粒分析軟件EDEM 2018(EDEM Solutionsinc.英國(guó)愛丁堡)完成顆粒模型構(gòu)建及堆積仿真，采用多物理場(chǎng)耦合有限元軟件COMSOL Multiphysics 5.3a(COMSOL INC. 瑞典斯德哥爾摩)、數(shù)據(jù)處理軟件Matlab R2012a(Math Works，美國(guó)馬薩諸塞州)和Origin 8.5(Atos Origin, 荷蘭阿姆斯特丹)完成堆積重現(xiàn)、模型構(gòu)建、數(shù)值計(jì)算和數(shù)據(jù)處理分析。使用矢量網(wǎng)絡(luò)分析儀(ZNB20, Rohde & Schwarz Ltd，德國(guó)慕尼黑)、85051B 7 mm/APC-7同軸空氣線(Agilent Technology，馬來(lái)西亞檳城)和鹵素水分測(cè)定儀(DHS-16，常州衡正電子儀器有限公司，精度為5 mg，水分范圍：0.00%～100.00%，溫度范圍：室溫～160 ℃)測(cè)量農(nóng)業(yè)物料的等效介電特性。使用游標(biāo)卡尺(MITUTOYO 500-173，日本香川，精度為0.02 mm)測(cè)量顆粒的規(guī)格參數(shù)。

2 結(jié)果與分析

2.1 谷物顆粒堆積模型的有效性

2.1.1真實(shí)與模擬谷料堆形態(tài)分析 采用真實(shí)稻谷籽粒堆和模擬稻谷籽粒堆的堆積角偏差來(lái)評(píng)價(jià)堆積模型的準(zhǔn)確性。真實(shí)谷粒堆、模擬谷粒堆的側(cè)向堆積狀態(tài)如圖2所示，從堆體的邊緣擴(kuò)散、結(jié)構(gòu)形態(tài)上看，二者的堆面形態(tài)基本吻合。實(shí)驗(yàn)測(cè)量得到真實(shí)谷粒的堆積角為20.20°。

圖2 真實(shí)谷料堆和模擬谷料堆側(cè)向形態(tài)比較Fig.2 Comparison of side view of accumulation form between real rice grains and simulation rice grains

2.1.2真實(shí)與模擬谷料堆誤差分析 采用圖像識(shí)別技術(shù)得到模擬谷料堆的單側(cè)邊緣輪廓曲線、擬合直線，如圖3所示，計(jì)算得到模擬谷粒堆的堆積角為20.11°。對(duì)比分析表明，二者的堆積角誤差小于0.5%，說(shuō)明基于離散元法的谷物顆粒堆積模型是有效性的。

圖3 模擬谷料的單側(cè)輪廓線及線性擬合直線Fig.3 Unilateral contour line and linear fitting of simulation rice grains

2.1.3模擬谷料堆準(zhǔn)確性分析 為了驗(yàn)證有限元法重現(xiàn)顆粒堆積狀態(tài)的準(zhǔn)確性，分別建立相應(yīng)的模型?；陔x散元法的橢球體谷物顆粒在不銹鋼材質(zhì)立方體裝載器件中的落料堆積情況和基于有限元法的橢球形顆粒堆積狀態(tài)重現(xiàn)情況如圖4所示?？梢钥闯觯w粒的大小、位置和方向完全一致，說(shuō)明基于有限元法的顆粒堆積狀態(tài)重現(xiàn)的計(jì)算程序是正確的。

圖4 不同技術(shù)方法形成的顆粒堆積狀態(tài)比較Fig.4 Comparison of particle accumulation by different technical methods

2.2 顆?？臻g分布對(duì)仿真模型數(shù)值的影響

圖5 不同的顆粒空間分布形式對(duì)計(jì)算結(jié)果的影響Fig.5 Effects of different particle spatial distribution forms on calculated results

2.3 堆積模型可行性分析

2.4 堆積模型準(zhǔn)確性驗(yàn)證結(jié)果分析

2.4.1不同體積分?jǐn)?shù)谷物物料驗(yàn)證結(jié)果 分別測(cè)量并獲得不同體積分?jǐn)?shù)條件下稻米物料、小麥物料的等效介電特性數(shù)據(jù)。同時(shí)，依據(jù)ε～ρ的關(guān)系式計(jì)算得到對(duì)應(yīng)籽粒的介電特性數(shù)據(jù)，在此基礎(chǔ)上，采用堆積模型對(duì)稻米物料、小麥物料的等效介電特性進(jìn)行反演計(jì)算，數(shù)值計(jì)算結(jié)果與實(shí)驗(yàn)測(cè)量值[22]、不同介電混合方程計(jì)算結(jié)果的對(duì)比情況如圖7所示?？梢钥闯?，在體積分?jǐn)?shù)42%～83%條件下，無(wú)論是等效介電常數(shù)還是等效介電損耗因子，堆積模型、LLL方程的數(shù)值結(jié)果與實(shí)驗(yàn)測(cè)量值最吻合，最大偏離幅度為-3.12%；采用CRI方程得到的數(shù)值結(jié)果大于測(cè)量值，最大偏離幅度達(dá)到8.36%；采用Lichtenecker(Li)方程得到的數(shù)值結(jié)果小于實(shí)驗(yàn)測(cè)量值，最大偏離幅度達(dá)到-8.97%。說(shuō)明堆積模型的數(shù)值結(jié)果具有一定的準(zhǔn)確性，該模型可用于上述微波頻段、體積分?jǐn)?shù)條件下谷物類農(nóng)產(chǎn)品的介電特性預(yù)測(cè)和分析。

圖7 堆積模型結(jié)果與實(shí)驗(yàn)測(cè)量值和其他公式計(jì)算結(jié)果比較Fig.7 Results comparison of between the model and the measured data and other formulas

圖8 小米堆積模型結(jié)果與實(shí)驗(yàn)測(cè)量值比較Fig.8 Results comparison between the model and the measured data for mill

3 討論