基于EDEM 的洛稻籽粒堆積仿真研究

孔令今楊開敏王夢月

(山東建筑大學 熱能工程學院，山東 濟南250101)

0 引言

在糧食儲藏過程中，稻谷等糧食籽粒的堆積進料是首要問題。 深入研究并掌握糧食籽粒堆積形成糧堆的機制，是簡化籽粒模型構(gòu)建過程和縮短籽粒堆積仿真模擬時間的重要措施。 在使用軟件模擬糧食籽粒的堆積過程時，首要解決的是如何進行非球形糧食籽粒的填充，但是在構(gòu)建非球形籽粒模型時，大多采用自動填充模式進行非球形籽粒的填充。 王美美等[1]對玉米籽粒進行自動填充，填充的過程需要引入二次開發(fā)系統(tǒng)，填充非球形籽粒所需的球形籽粒數(shù)量為351 個；原建博等[2]通過非接觸式三維激光掃描法重構(gòu)水稻籽粒，自動填充籽粒模型的填充球數(shù)為203 個，填充過程耗時為1 h；劉彩玲等[3]通過三維激光掃描等技術(shù)重構(gòu)水稻籽粒，填充顆粒模型的球數(shù)為87 個。 但是，在實際的計算過程中，上述方法由于EDEM 的仿真模擬計算需要考慮到每一個球體之間的接觸力和每一個球體與壁面的接觸力，導致用數(shù)量多的球形籽粒填充非球形籽粒模型時計算的過程非常復雜，大幅增加了計算時間。

在仿真模擬過程中，多采用對比實驗與堆積角測量以驗證籽粒仿真模型各項參數(shù)的準確性。 目前，實驗法測籽粒物料堆堆積角的主要方法有排除法[4-5]、注入法[6-8]和傾斜法[9]。 在大型糧食儲備倉庫中，進料時一般采用皮帶傳送，將糧食籽粒以拋物線形式送入糧倉，糧食籽粒在糧倉內(nèi)的堆積過程為自由落體式下落過程，而注入法最符合實際的籽粒堆積過程，因此模擬實驗中常采用注入法，使洛稻通過漏斗自上部流入，在底部形成圓錐狀堆積角。

文章選取洛稻998 籽粒測量其長、寬、厚、球度等尺寸參數(shù)，統(tǒng)計出糧食籽粒粒徑的分布規(guī)律等；利用軟件EDEM 基于離散元法實現(xiàn)糧食籽粒堆積的仿真模擬，對非球形籽粒進行手動填充，同時實驗采用注入法測得籽粒物料堆的堆積角，將結(jié)果進行對比分析，驗證軟件EDEM 以手動填充方式進行非球形洛稻籽粒填充的可行性，以期能夠減少籽粒模型構(gòu)建單體的數(shù)量，簡化籽粒構(gòu)建過程和縮短仿真模擬時間，為以后測量籽粒堆積體的孔隙率和探究料堆的傳熱傳質(zhì)等提供尺寸參數(shù)和籽粒模型，可以更好地指導實際的工程。

1 基于離散元法的籽粒模型建立

1.1 籽粒的理論接觸模型

對球形籽粒建立接觸模型，如圖1 所示。

圖1 球形籽粒接觸模型圖

籽粒i在任意時刻t的控制方程[10]由式(1)和(2)表示為

式中:mi為籽粒i的質(zhì)量，kg；Ii為轉(zhuǎn)動慣量，kg·m2；vi為速度，m/s；wi為角速度，rad/s；ki為與籽粒i發(fā)生接觸的顆粒的數(shù)量；Tij為扭矩，N·m；mig為籽粒i的重力，N；fcontact，ij為接觸力，N；fdamp，ij為阻力，N。

1.2 離散元法的接觸模型方程

離散元法DEM(Discrete Element Method) 是Cundull 和Strack 于1971 年提出的一種模擬籽粒材料靜態(tài)與動態(tài)行為的研究方法[10]。 離散元法將求解空間離散為離散單元陣，各個單元之間是不連續(xù)、可分離的，即一個單元可以與其相鄰的單元接觸或分開[11]，可將分析對象看成充分多的離散介質(zhì)，每個籽?；驂K體為1 個單元。 DEM 的基本思路是利用籽粒接觸模型計算相互接觸單元間的受力，并利用牛頓第二運動定律求解籽粒的運動參量[12]，即全過程中的每一時刻各籽粒間的相互作用，接觸力和每個離散單元的運動狀態(tài)。

引入時間步長的概念，遞歸地求解離散單元的運動過程，交替反復運算，實現(xiàn)籽粒運動情況的預測[2]。 對于非球形籽粒的接觸模型可選用軟球模型，其將籽粒間的法向力簡化為彈簧、阻尼器和滑動器，引入彈性系數(shù)和阻尼系數(shù)參量，不考慮籽粒表面變形，依據(jù)籽粒間發(fā)現(xiàn)重疊和切向位移計算接觸力[3]。

接觸模型方程[13]如下:

將籽粒接觸過程中的振動運動進行法向和切向分解，法向振動運動方程式Fn由式(3)表示為

籽粒接觸過程的切向振動運動表現(xiàn)為切向滑動與籽粒的滾動，其切向分量Fs方程式、籽粒所受的外力矩M分別由式(4)和(5)表示為

式中:mi，j為籽粒i、j的等效質(zhì)量，kg；Ii，j為籽粒i、j的等效轉(zhuǎn)動慣量，kg·m2；r為旋轉(zhuǎn)半徑，m；μn、μs分別為籽粒的法向和切向相對位移，m；θ為籽粒自身的旋轉(zhuǎn)角度，rad；Kn、Ks分別為接觸模型中的法向及切向彈性系數(shù)；cn、cs分別為接觸模型中的法向及切向阻尼系數(shù)。

2 洛稻998 籽粒堆積的實驗測定

2.1 實驗測量儀器

實驗儀器包括500 g 托盤電子秤(精度為0.01 g)、0～150 mm 電子游標卡尺(精度為0.02 mm)、量角器(精度為0.01°)、不銹鋼漏斗(底徑為15 mm)、不銹鋼托盤(200 mm×200 mm)。

2.2 洛稻998 籽粒尺寸參數(shù)的測量

實驗隨機選取100 粒洛稻998 籽粒，重約3.06 g，利用游標卡尺對籽粒的長度L、寬度W、厚度T進行測量[14]，計算出球度，如圖2 所示。 為避免隨機誤差，可測量100 次籽粒尺寸，求得平均值和標準差，獲得數(shù)據(jù)見表1。

圖2 籽粒尺寸參數(shù)測量圖

表1 洛稻998 尺寸參數(shù)測定表/mm

對長度、寬度和厚度所測數(shù)據(jù)進行數(shù)據(jù)統(tǒng)計，并依據(jù)測量數(shù)據(jù)繪制洛稻998 籽粒的幾何尺寸分布曲線，如圖3 所示。 洛稻998 籽粒長度、寬度、厚度尺寸分別集中分布在7.2～7.4、3.3～3.4 和2.3～2.5 mm之間，尺寸參數(shù)的標準差較小，數(shù)據(jù)集中性較強，尺寸參數(shù)的比例分布符合正態(tài)分布規(guī)律。

圖3 洛稻998 種子幾何尺寸分布圖

王磊[15]在對類球形蔬菜種子進行了研究，得出了球度SP的計算方程由式(6)表示為

將所求的洛稻998 籽粒尺寸參數(shù)帶入式(6)，計算得到球度為53.24%，因此在EDEM 中建立模型時將籽粒視為橢球型。

張桂花等[7]在對包衣稻種物理特性的研究中，將稻種的3 軸尺寸測得后以3 軸平均值作為單一粒徑。 將洛稻998 籽粒的尺寸參數(shù)進行計算，得到單一粒徑為4.40 mm。

綜上所述，實驗測得的洛稻998 籽粒尺寸可以作為下一步EDEM 仿真的數(shù)據(jù)依據(jù)。

2.3 洛稻998 籽粒料堆的堆積角測定

實驗中參考田曉光等[6]做法，制作固定底面的不銹鋼漏斗裝置，漏斗底直徑為15 mm，約為顆粒單一粒徑的3 倍以上，防止在堆積過程中發(fā)生堵塞等對實驗結(jié)果產(chǎn)生影響。 設置不銹鋼漏斗的底端和不銹鋼托盤的距離為100 mm，盡量真實地模擬糧食籽粒在糧倉堆積過程中的自由下落過程，以保證測得的堆積角數(shù)據(jù)可靠。

李海偉[16]根據(jù)JB/T 9014.7—1999《連續(xù)輸送設備散狀物料堆積角的測定》[17]對散狀物料的堆積角進行了測定，堆積角定義為當散狀物料自由堆積在水平面上，且保持堆形面與水平面的最大夾角。 籽粒在糧倉實際進料過程中，底部保持靜止，糧食籽粒自由下落，靜態(tài)堆積角θi測量方程由式(7)表示為

式中:h為料堆的測量高度，mm；d1為料堆測量上底底徑，mm；d2為料堆的測量底徑，mm。

隨機選取洛稻998籽粒共計2000粒，重約為61.27 g。 將籽粒傾倒入不銹鋼漏斗中，使籽粒在重力的作用下運動，最終形成的料堆靜態(tài)堆積角如圖4 所示。

圖4 料堆靜態(tài)堆積角圖

重復實驗10 次，計算料堆的測量高度和料堆的測量底徑并取平均值，得到的統(tǒng)計數(shù)據(jù)見表2。

表2 料堆測量的高度和底徑表/mm

測量料堆的平均高度h、上底底徑d1、底徑d2分別為4.05、3.98 和13.93 mm，將值代入到式(7)可得實驗測得的堆積角為39.15°。

3 基于EDEM 的籽粒堆積仿真

3.1 非球形籽粒單體建模

基于上述實驗所測得的尺寸參數(shù)可知，洛稻998 的外形可視為橢球型。 在EDEM 中用球形籽粒進行手動填充，分別選取7、11、16 個球形籽粒的3種工況，進行洛稻998 的非球形籽粒單體的填充，便于在測試堆積角過程中進行對比試驗，填充模型如圖5 所示。

由圖5 可知，非球形籽粒填充單體模型的外形接近于洛稻998 籽粒，在模擬過程中的單一顆粒質(zhì)量為0.03 g，實際籽粒的質(zhì)量約為0.030635 g，在密度參數(shù)設定相同的情況下，非球形籽粒模型的體積與實際籽粒體積大致相同。

圖5 非球形籽粒填充模型圖

3.2 非球形籽粒堆積模型建立

在實際糧食籽粒的堆積過程中，糧食籽粒之間會存在靜摩擦力和滑動摩擦力，籽粒與壁面之間也同樣存在這兩種摩擦力；籽粒與籽粒之間或者籽粒與壁面之間發(fā)生碰撞，會發(fā)生形變。 實驗測定中的籽?！蚜：妥蚜！诿嬷g的滾動摩擦力在EDEM 仿真模擬時以滾動摩擦系數(shù)的設置體現(xiàn)；當籽粒從漏斗底部直徑為15 mm 的漏口自由落下時，其與不銹鋼托盤之間的滾動摩擦力和彈性模量會對籽粒堆積角的角度產(chǎn)生影響，在EDEM 仿真過程中同樣要設置彈性模量等參數(shù)，才能夠更加準確地模擬實際糧食籽粒在進料過程中的堆積現(xiàn)象，使仿真模擬的結(jié)果符合工程實際。

在EDEM 軟件仿真構(gòu)建時，設置漏斗材質(zhì)為不銹鋼材質(zhì)，內(nèi)部構(gòu)建虛擬的靜態(tài)籽粒工廠，顆?？倲?shù)設為2000 粒，底部托盤材質(zhì)為不銹鋼材質(zhì)。 洛稻998-洛稻998、洛稻998-不銹鋼的恢復系數(shù)、靜摩擦系數(shù)和滾動摩擦系數(shù)等設置[18-22]見表3。

表3 基于EDEM 的籽粒堆積仿真模型所需參數(shù)表

對于洛稻998 的籽粒物理模型方程，選用軟件中的Hertz-Mindlin(no slip)模型作為顆粒接觸模型，建立籽粒堆積模型如圖6 所示。 籽粒物堆積角仿真模擬實驗采用Euler 模型，時間步長設為Rayleigh 時間步長的30%，計算單元格設為最小單元格的3 倍，總的仿真模擬時間設為6 s，其中重力加速度為9.81 m/s2，數(shù)據(jù)保存時間設為每0.01 s 保存一次，待仿真模擬過程完成后，進行堆積角的測量。

圖6 籽粒堆積模型圖

3.3 仿真模型堆積角的測量與分析

由7、11 和16 個球形籽粒構(gòu)建的單體洛稻998進行堆積模型仿真模擬，如圖7 所示，仿真模擬堆積完成的時間分別為2.68、2.84、2.96 s，實際的堆積時間分別為0.57、1.6 和2.7 h。

根據(jù)實驗數(shù)據(jù)可知，在最大限度不改變非球形籽粒外形尺寸的前提條件下，隨著堆積所需球形籽粒個數(shù)的減少，堆積所需的時間約可縮短為采用自動填充法用多個球形顆粒單體填充籽粒模型進行堆積所需時間的70%。

糧堆處于靜態(tài)，將EDEM 后處理中的剪切(Clipping)和工具(Tools)中的量角器(Protractor)相結(jié)合測量堆積角的角度，以XOZ 面為基準建立切面，利用量角器進行角度的測量，如圖8 所示，不斷改變x軸和z軸的數(shù)值，使測得的角度更加接近實際值。 測量10 次后取平均值，測得的堆積角數(shù)據(jù)統(tǒng)計結(jié)果見表4。

圖7 3 種單體籽粒堆積圖

由表4 可知，由7、11、16 個球形籽粒構(gòu)建的單體洛稻模型在相同的參數(shù)條件下堆積完成后測得的堆積角分別為42.82°、42.08°、41.88°，3 種工況測得的仿真堆積角數(shù)值上差別很小，與實驗注入法測得的堆積角39.15°對比，偏差分別為9.37%、7.48%、6.97%，數(shù)據(jù)之間的偏差均＜10%，從而在EDEM 中驗證了手動填充非球形籽粒具有可行性。

圖8 堆積角測量圖

表4 3 種單體模型堆積建模測得的堆積角統(tǒng)計表/(°)

4 結(jié)論

根據(jù)上述研究結(jié)果，可得出主要結(jié)論如下:

(1) 大量實測數(shù)據(jù)得到了洛稻籽粒平均尺寸參數(shù)，通過繪制長度、寬度、厚度的比例分布曲線，驗證其符合正態(tài)分布，并基于此參數(shù)進行EDEM 非球形籽粒手動填充，填充的非球形籽粒模型體積與實際籽粒體積相近，進而驗證了所測得的籽粒尺寸參數(shù)正確性。

(2) 球形籽粒構(gòu)建的3 種單體洛稻998 模型進行堆積仿真模擬形成的靜態(tài)堆積角與實驗注入法測得的堆積角相比較，計算數(shù)據(jù)偏差＜10%，說明實際糧食籽粒的下落堆積過程可以通過EDEM 仿真模擬實現(xiàn)，驗證了手動填充非球形籽粒具有可行性。

(3) 在保證非球形籽粒尺寸參數(shù)和接觸參數(shù)等參數(shù)與實際測量值相同的情況下，EDEM 用較少的球形籽粒填充出來的非球形籽?？梢阅M實際的糧食進料過程，能夠簡化非球形籽粒的填充過程，堆積所需的時間約可縮短為采用自動填充法所需時間的70%。