玉米秸稈飼料除塵篩出物離散元數值模擬參數標定

肖子卿 田海清* 張 濤 王 迪 盛 越 李大鵬 劉 飛 周建成

(1.內蒙古農業大學 機電工程學院，呼和浩特 010018；2.內蒙古領異農科技發展有限公司，呼和浩特 010200)

我國秸稈資源豐富，將秸稈資源轉化為飼草料，不僅能解決畜牧養殖飼草料缺乏的問題，而且還能改善生態環境。現有秸稈飼料收獲機可對田間秸稈收獲并破碎為直接可飼喂的段狀或絲狀，但收獲機撿拾的秸稈中摻雜大量塵土，從而導致收獲的秸稈飼料清潔程度較低，適口性較差。因此對收獲后的秸稈飼料進行除塵提高清潔度，對提高飼料利用率意義重大。

目前帶有除塵裝置的飼料收獲機市場所占比例逐年增加。如張開飛等設計的小麥收獲機復合除塵裝置，實現了割臺靜電除塵和噴淋除塵相結合的復合除塵；杜洋研制的秸稈除塵設備，實現了對秸稈的凈化，進而保障秸稈二次加工產物的質量；李嘉誠等設計的新型無塵脫粒機，采用空氣幕防塵技術從源頭對秸稈除塵；類似的還有秸稈收獲機螺旋輸送除塵裝置、具有除塵功能的滾筒篩和帶有除塵結構的飼料粉碎機等。近幾年結合收獲、除塵、打捆為一體的飼料收獲機也開始興起，市場上出現的有順邦4JL-2.2型秸稈飼料收獲機和綠萌機械9YFG-2.2P型秸稈揉碎除塵方草捆打捆機等。帶有除塵裝置的飼料收獲機雖逐漸增加，但對除塵機理的研究還較為少見。研究除塵過程中秸稈飼料與塵土的運動分離機理，對優化設計除塵裝置，提升除塵效率有積極意義。

本研究前期對收獲機收獲后的飼料采用標準篩篩分時發現，所篩出的雜質為塵土與秸稈碎末的混合物。通過試驗方法分析秸稈飼料與塵土運動分離機理較為困難，而離散元法是模擬混合物內部運動狀態的有效方法，國外相關研究在運用此方法對固體顆粒、粉末等混合物的內部運動及參數標定中已取得較多成果。國內也有研究將此方法應用在混合物中存在秸稈時的運動狀態分析及參數標定中，如張濤等為研究玉米秸稈與揉碎機的相互作用及運動規律，標定了玉米秸稈與揉碎機錘片、秸稈間的接觸參數。王韋韋等標定了玉米秸稈粉料與不銹鋼板、秸稈粉料間的接觸參數；馬彥華等標定了苜蓿秸稈與45鋼、苜蓿秸稈間的接觸參數。綜上，秸稈已有較為完整的物理參數與數值模擬參數，但對于塵土與秸稈碎末之間的物理參數與數值模擬參數的研究仍較為少見。

本研究擬采用物理試驗與數值模擬相結合的方法對玉米秸稈飼料的篩出物(塵土與秸稈碎末)進行離散元數值模擬參數標定。利用EDEM軟件對塵土與秸稈碎末運動情況進行數值模擬，以休止角為響應值，分別通過Plackett-Burman試驗、最陡爬坡試驗和Box-Behnken試驗確定較優的塵土與秸稈碎末的數值模擬參數值，以期為塵土與秸稈運動分離機理分析及除塵裝置的優化提供數據支持。

1 材料及方法

1.1 試驗材料

2020年11月進行秸稈收獲作業試驗，試驗地點位于內蒙古自治區呼和浩特市土默特左旗大圪賁村(121°58′ E，36°15′ N)試驗田，試驗所用秸稈飼料和土壤樣品均取自收獲作業后的試驗田。根據GB/T 5917.1—2008 《飼料粉碎粒度測定兩層篩篩分法》選用篩孔直徑為3和7 mm的標準篩對秸稈飼料進行篩分，篩出物包括塵土和秸稈碎末。利用游標卡尺測得秸稈碎末平均直徑為2 mm，平均長度5 mm。利用標準土壤篩對塵土進行篩分，塵土平均粒徑為0.5 mm。

1.2 試驗方法

運用休止角物理試驗與數值模擬相結合的方法對玉米秸稈飼料的篩出物(塵土與秸稈碎末)進行離散元數值模擬參數標定，利用休止角測定儀進行塵土-秸稈碎末休止角測定試驗，測量塵土與秸稈碎末的錐體高度，計算得出休止角。利用EDEM軟件進行塵土-秸稈碎末休止角數值模擬試驗，分別通過Plackett-Burman試驗、最陡爬坡試驗和Box-Behnken試驗確定最優參數組合。應用最優參數組合進行休止角和標準篩除塵率數值模擬，與休止角和標準篩除塵率物理試驗結果相對比，驗證標定參數組合的準確性。

1.3 休止角測定試驗

休止角測定采用FT-104B型休止角測定儀，漏斗上方進料口直徑為135 mm，下方出料口直徑20 mm，高度140 mm。透明環形容器位于漏斗下方50 mm處，直徑為100 mm，高度25 mm。將儀器平置于操作平臺后，將基板調整水平，將透明環形容器置于基板中心位置，使容器的中心與漏斗軸線一致。參考文獻中對休止角測定方法，將篩選出的塵土和秸稈碎末勻速倒入儀器中的漏斗，待其流出停止后120 s，測量錐體高度。試驗重復5次取平均值，錐體平均高度為69.5 mm，計算得到秸稈碎末和塵土顆粒堆物理試驗休止角

α

為41.669°。圖1為休止角測定裝置及試驗材料。

1.漏斗；2.秸稈碎末與土；3.透明環形容器；4.基板1.Funnel; 2.Straw debris and soil;3.Transparent annular container; 4.Substrate圖1 休止角測定裝置及試驗材料Fig.1 Angle of repose measuring device and test material

物理試驗休止角(

α

)計算公式為：

(1)

式中：

h

為錐尖到基板的高度，mm；

L

為錐體底端到基板的高度，mm；

R

為透明環形容器半徑，mm。

1.4 接觸參數測定

秸稈飼料的篩出物中有塵土與秸稈碎末，秸稈碎末中包含秸稈外皮和秸稈內瓤，所以需要測定的接觸參數有塵土-外皮滾動摩擦因數、塵土-外皮靜摩擦因數、塵土-內瓤滾動摩擦因數和塵土-內瓤靜摩擦因數，測定結果作為數值模擬輸入參數選擇依據。

物理試驗選用的秸稈均為飼料收獲機收獲后殘留于田間的漏收秸稈，本研究前期已測得其外皮密度為508 kg/m、內瓤密度40 kg/m、總密度110.9 kg/m、泊松比0.3、剪切模量64 MPa。與秸稈接觸的材料為鋼，密度為7 865 kg/m、泊松比0.3、剪切模量79 000 MPa。

1

.

4

.

1

靜摩擦因數測定試驗使用CNY-1型斜面儀對塵土-外皮和塵土-內瓤靜摩擦因數進行測定。斜面上的物料受3個力作用，分別為重力

G

、支持力

F

和摩擦力

F

。根據力的正交分解得：

(2)

得出：

μ

=tan

θ

(3)

式中：

G

為物料所受重力，N；

m

為物料質量，kg；

g

為重力加速度，m/s；

F

為物料所受支持力，N；

F

為物料所受摩擦力，N；

θ

為斜面儀的傾斜角度,(°)；

μ

為摩擦因數。將取自試驗田的土壤樣品裝于直徑90 mm、高10 mm的透明圓柱容器內，平整其表面并固定于斜面儀上。將秸稈切割為30 mm的規則圓柱體，為防止在測定過程中秸稈發生滾動，將2個規則秸稈段粘在一起并置于土樣之上。調整斜面儀使其指針位于零刻度處，將斜面儀緩慢逆時針轉動直至秸稈在土樣上勻速滑動，記錄斜面儀傾斜角度

θ

，根據式(3)求得靜摩擦因數。重復10次試驗，得出秸稈在土樣上發生滑動的角度為21.801°～34.992°，塵土-外皮的靜摩擦因數為0.4～0.7。去除秸稈外皮，將內瓤切割為長、寬、高均為15 mm的規則正方體，將內瓤置于土樣之上，緩慢逆時針轉動斜面儀直至內瓤方塊在土樣上緩慢滑動，記錄斜面儀傾斜角度

θ

。重復10次試驗，得出秸稈內瓤在土樣上發生滑動的角度為30.964°～41.987°，塵土-內瓤的靜摩擦因數為0.6～0.9。

1

.

4

.

2

滾動摩擦因數測定試驗測定滾動摩擦因數的儀器也為CNY-1型斜面儀，測定方法與靜摩擦因數測定方法類似。將切割好的單個秸稈段置于土樣之上，調整斜面儀使其指針位于零刻度處，緩慢逆時針轉動斜面儀，當秸稈開始滾動時記錄斜面儀傾斜角度

θ

。重復10次試驗，得出秸稈在土樣上發生滾動的角度范圍為11.310°～26.565°，塵土-外皮的滾動摩擦因數為0.2～0.5。去除秸稈外皮，將圓柱形的秸稈內瓤置于土樣之上，將指針轉動到零刻度處時開始緩慢逆時針轉動斜面儀，當秸稈內瓤開始滾動時停止轉動，記錄斜面儀傾斜角度

θ

。重復10次試驗，得出秸稈內瓤在土樣上發生滾動的角度為16.699°～30.964°，塵土-內瓤的滾動摩擦因數為0.3～0.6。

為了保證數值模擬的準確性，根據物理試驗測得的數據，得出測定參數的取值范圍，為后續的數值模擬輸入參數取值提供依據，最后通過數值模擬確定參數的最優取值。

2 秸稈與塵土離散元模型的建立及參數標定

2.1 離散元模型的建立

秸稈飼料經篩分后留下的雜質包括秸稈碎末和塵土，秸稈碎末包含秸稈外皮和秸稈內瓤，外皮和內瓤的長度與直徑差異性不大。因此在EDEM軟件中建立數值模擬模型時將秸稈外皮和秸稈內瓤的數值模擬模型都設置為長5 mm、直徑2 mm的圓柱體，塵土的數值模擬模型設置為直徑0.5 mm的球體(圖2)。

圖2 秸稈外皮、內瓤(a)和塵土(b)數值模擬標定模型Fig.2 Calibration model for numerical simulation of straw skin and pulp (a) and soil particle (b)

導入休止角測定儀模型，模型漏斗進料口直徑為135 mm，出料口直徑20 mm，下方圓盤直徑50 mm、高度25 mm(圖3)。設置顆粒模型參數，結合GEMM數據庫與相關文獻，得出部分參數取值見表1。

1.漏斗；2.秸稈碎末與土；3.圓盤1.Funnel; 2.Straw debris and soil; 3.Disk圖3 休止角數值模擬模型Fig.3 Numerical simulation model of angle of repose

表1 塵土與秸稈碎末休止角數值模擬試驗部分參數取值
Table 1 Some parameter values of dust and straw debris angle of repose numerical simulation test

參數Parameter數值Numericalvalue參數Parameter數值Numericalvalue秸稈外皮-外皮靜摩擦因數Straw skin-straw skin static friction factor0.142秸稈外皮-外皮碰撞恢復系數Straw skin-straw skin restitution coefficient0.485秸稈外皮-外皮滾動摩擦因數Straw skin-straw skin rolling friction factor0.078秸稈外皮-內瓤碰撞恢復系數Straw skin-straw pulp restitution coefficient0.263秸稈外皮-內瓤靜摩擦因數Straw skin-straw pulp static friction factor0.495秸稈內瓤-內瓤碰撞恢復系數Straw pulp-straw pulp restitution coefficient0.348秸稈外皮-內瓤滾動摩擦因數Straw skin-straw pulp rolling friction factor0.166秸稈外皮-鋼碰撞恢復系數Straw skin-steel restitution coefficient0.663秸稈內瓤-內瓤靜摩擦因數Straw pulp-straw pulp static friction factor0.427秸稈內瓤-鋼碰撞恢復系數Straw pulp-steel restitution coefficient0.325秸稈內瓤-內瓤滾動摩擦因數Straw pulp-straw pulp rolling friction factor0.144秸稈外皮-塵土碰撞恢復系數Straw skin-soil restitution coefficient0.213秸稈外皮-鋼靜摩擦因數Straw skin-steel static friction factor0.226秸稈內瓤-塵土碰撞恢復系數Straw pulp-soil restitution coefficient0.135秸稈外皮-鋼滾動摩擦因數Straw skin-steel rolling friction factor0.119塵土泊松比Soil poisson’s ratio0.4 秸稈內瓤-鋼靜摩擦因數Straw pulp-steel static friction factor0.375塵土密度ρ/(kg/m3)Soil density2 600秸稈內瓤-鋼滾動摩擦因數Straw pulp-steel rolling friction factor0.150

本研究需標定的試驗參數取值分別為塵土-鋼滾動摩擦因數、塵土-鋼靜摩擦因數、塵土-外皮滾動摩擦因數、塵土-外皮靜摩擦因數、塵土-內瓤滾動摩擦因數、塵土-內瓤靜摩擦因數、塵土-塵土滾動摩擦因數、塵土-塵土靜摩擦因數、塵土接觸模型JKR表面能系數、塵土-塵土碰撞恢復系數及塵土-鋼碰撞恢復系數。

將篩分后的塵土和秸稈碎末稱重，將秸稈碎末進行人工挑選，剩余的塵土進行稱重，得出塵土所占比率。重復3次試驗，求其平均值，得出塵土占比率為9.34%。塵土具有一定的粘性，在數值模擬過程中選用Hertz-Mindlin with JKR接觸模型，設置生成秸稈外皮和秸稈內瓤顆粒的質量均為15.69 g，生成塵土質量為3.27 g。時間步長5.026×10s，總時間5 s。在漏斗上方建立顆粒工廠，每1 s生成顆粒質量為100 g，5 s后所有顆粒從漏斗下端流出，在下方圓盤上形成穩定的顆粒堆，通過后處理測定休止角。

2.2 數值模擬參數的標定

2

.

2

.

1

顯著性因素的篩選

結合物理試驗所獲得的接觸特性參數以及查閱文獻，確定數值模擬中各參數的取值范圍見表2，利用Design-Expert軟件中的Plackett-Burman進行試驗設計，以數值模擬顆粒堆休止角為響應值，篩選出對數值模擬休止角影響最顯著的因素。

表2 休止角數值模擬影響因素Plackett-Burman試驗參數取值
Table 2 Influencing factors of numerical simulation of angle of repose and values of Plackett Burman test parameters

試驗參數Test parameter低水平Low-level高水平High-level塵土-鋼滾動摩擦因數X1Soil-steel rolling friction factor X10.050.20塵土-鋼靜摩擦因數X2Soil-steel static friction factor X20.51.2塵土-外皮滾動摩擦因數X3Soil-straw skin rolling friction factor X30.20.5塵土-外皮靜摩擦因數X4Soil-straw skin static friction factor X40.40.7塵土-內瓤滾動摩擦因數X5Soil-straw pulp rolling friction factor X50.30.6塵土-內瓤靜摩擦因數X6Soil-straw pulp static friction factor X60.60.9塵土-塵土滾動摩擦因數X7Soil-soil rolling friction factor X70.10.2塵土-塵土靜摩擦因數X8Soil-soil static friction factor X80.321.16

表2(續)

試驗參數Test parameter低水平Low-level高水平High-level塵土接觸模型JKR表面能系數X9JKR surface energy coefficient of soil contact model X9011.25塵土-塵土碰撞恢復系數X10Soil-soil restitution coefficient X100.150.75塵土-鋼碰撞恢復系數X11Soil-steel restitution coefficient X110.20.5

在Plackett-Burman試驗設計中，將表2內11個試驗參數取值的最高值和最低值分別標為+1、-1 水平，設置一個中心點即選取每個試驗參數最高值和最低值的平均值作為0水平，共進行13次數值模擬。數值模擬得出的塵土與秸稈碎末數值模擬模型堆見圖4，利用EDEM軟件中后處理的量尺工具測量顆粒堆的高度，計算得出數值模擬休止角

α

。

H為數值模擬顆粒堆高度；R為圓盤半徑。H is numerical simulation of particle height; R is the radius of disk.圖4 塵土與秸稈碎末數值模擬模型堆Fig.4 Numerical simulation model pile of soil and straw debris

漏斗下方的圓盤半徑

R

為25 mm，已知顆粒堆的高度，利用式(4)求得顆粒堆的休止角(

α

)：

(4)

式中：

H

為數值模擬顆粒堆高度，mm；

R

為圓盤半徑，mm。計算所得的數值模擬休止角(

α

)見表3，利用Design-Expert軟件對結果進行方差分析。影響休止角的各個因素顯著性排序見表4：影響休止角的最顯著性因素有3個，其中塵土接觸模型JKR表面能系數(

X

)對休止角的影響最大，其次為塵土-塵土碰撞恢復系數(

X

)和塵土-塵土滾動摩擦因數(

X

)。在進行后續數值模擬時，其余不具顯著性影響的因素，參數取值均為平均值，分別為塵土-鋼滾動摩擦因數0.125、塵土-鋼靜摩擦因數0.85、塵土-外皮滾動摩擦因數0.35、塵土-外皮靜摩擦因數0.55、塵土-內瓤滾動摩擦因數0.45、塵土-內瓤靜摩擦因數0.75、塵土-塵土靜摩擦因數0.74、塵土-鋼碰撞恢復系數0.35。

表3 休止角數值模擬影響因素Plackett-Burman試驗設計及結果
Table 3 Influencing factors of numerical simulation of angle of repose and Plackett-Burman experimental design and results

序號Serialnumber各因素水平編碼 Coding of each factor levelx1x2x3x4x5x6x7x8x9x10x11休止角α2/(°)Repose angle11-111-1111-1-1-141.3482-1-11-111-1111-142.6143-1111-1-1-11-11130.9644-11-111-1111-1-147.2035-1-1-11-111-111141.348

表3(續)

序號Serialnumber各因素水平編碼 Coding of each factor levelx1x2x3x4x5x6x7x8x9x10x11休止角α2/(°)Repose angle6-1-1-1-1-1-1-1-1-1-1-138.6607-111-1111-1-1-1138.66080000000000040.03091-1-1-11-111-11141.3481011-1-1-11-11-11145.00011111-1-1-11-111-144.421121-1111-1-1-11-1143.8311311-1111-1-1-11-137.954

注：,,…,分別為,,…,的水平編碼。,,…,各因素的物理意義見表2，下同。
Note： , , …, are the horizontal codes of , , …, respectively. The physical meanings of the factors , , …, are shown in Table 2. The same below.

表4 休止角數值模擬影響因素顯著性排序
Table 4 Significance ranking of influencing factors of numerical simulation of angle of repose

因素Factor效應Effect均方和Mean square sum影響率/%Influence ratio排序SortX12.4117.438.764X2-0.832.051.039X3-1.617.783.917X4-1.345.392.718X51.658.134.086X60.080.020.5510X72.5519.499.793X80.601.080.5411X95.92104.9852.731X10-2.6821.5010.802X11-1.8410.155.105

2

.

2

.

2

最陡爬坡試驗通過Plackett-Burman試驗確定了3個對休止角影響最顯著的因素，以數值模擬所測休止角(

α

)和物理試驗休止角(

α

)的相對誤差作為評價指標，利用最陡爬坡試驗可迅速確定各因素的最佳參數取值范圍。最陡爬坡試驗設計及結果見表5：隨著塵土接觸模型JKR表面能系數(

X

)、塵土-塵土碰撞恢復系數(

X

)和塵土-塵土滾動摩擦因數(

X

)的不斷增大，試驗所測休止角先增大后減小，整體呈現一個拋物線形狀。觀察各試驗相對誤差發現，5號試驗相對誤差最小，所以最佳參數取值范圍在試驗5附近。繼續取試驗4、試驗5和試驗6參數分別為低、中和高水平進行后續數值模擬。

表5 休止角顯著因素最陡爬坡試驗設計及結果
Table 5 Design and results of the steepest climbing test for the apparent factors of angle of repose

序號Serial number因素 FactorX7X9X10休止角α2/(°)Repose angle相對誤差/%Relative error10.1000.1537.9548.9120.12 2.250.2740.0303.9330.14 4.500.3943.8315.2040.16 6.750.5145.0008.0050.18 9.000.6341.3480.7860.2011.250.7540.0303.93

2

.

2

.

3

Box

-

Behnken試驗設計及結果分析對最陡爬坡試驗結果分析，以試驗相對誤差為評價指標，確定了試驗最優參數取值范圍。以試驗4、試驗5和試驗6參數分別為低(-1)、中(0)和高水平(1)進行Box-Behnken試驗設計(表6)。試驗中取3個中心點對試驗誤差進行評估，得出Box-Behnken試驗設計表(表7)。Box-Behnken試驗共進行15次數值模擬，當數值模擬完成后，同樣利用EDEM中后處理的量尺工具測量顆粒堆高度，求得數值模擬休止角(

α

)結果見表7。

表6 休止角顯著因素最優參數取值
Table 6 Optimal parameter value of apparent factor of angle of repose

因素Factors低水平Low-level中水平Medium level高水平High-levelX70.160.180.20X96.759.0011.25 X100.510.630.75

表7 休止角顯著因素Box-Behnken試驗設計及結果
Table 7 Box-Behnken test design and results of significant factors of repose angle

序號Serial number因素水平編碼 Coding of each factor levelx7x9x10休止角α2/(°)Repose angle100043.232-1-1038.88310-139.35410142.61501-142.33601142.617-10138.668-11041.32900042.96100-1-139.0011-10-140.031200043.56131-1039.351411045.32150-1140.03

利用Design-Expert軟件對其結果進行多元回歸分析，得出顯著性參數與數值模擬休止角(

α

)的二次多項式方程為：

α

=43

.

28+0

.

97

X

+1

.

79

X

+0

.

40

X

+
0

.

88

X

X

+1

.

16

X

X

-0

.

19

X

X

-
1

.

45

X

-0

.

62

X

-1

.

67

X

(5)

對Box-Behnken試驗進行方差分析，統計結果見表8。塵土-塵土滾動摩擦因數(

X

)、塵土接觸模型JKR表面能系數(

X

)、塵土-塵土滾動摩擦因數與塵土接觸模型JKR表面能系數(

X

X

)、塵土-塵土滾動摩擦因數與塵土-塵土碰撞恢復系數(

X

X

)、塵土-塵土滾動摩擦因數的二次項(

X

)和塵土-塵土碰撞恢復系數的二次項(

X

)對休止角的影響顯著。模型的

P

<0.05，說明休止角與回歸方程關系顯著；失擬項的

P

=0.33>0.05，表明方程擬合較好；變異系數CV=1.18%，表明試驗可靠性較高；決定系數

R

=0.980 7，校正決定系數兩者的數值均接近1，模型可以較準確的反映實際狀況。

表8 休止角顯著因素Box-Behnken試驗設計回歸模型方差分析
Table 8 Box-Behnken experimental design regression model analysis of variance for significant factors of angle of repose

方差來源Source of variation均方和Mean square sum自由度Degree of freedom均方Mean squareP值P value模型 Model60.4096.710.000 9 X77.4917.490.002 5 X925.63125.630.000 1 X101.2811.280.067 8 X7X93.1213.120.015 2 X7X105.3615.360.005 1 X9X100.1410.140.476 2 X727.7417.740.002 3 X921.4111.410.058 8 X10210.33110.330.001 2殘差 Residual1.1950.24失擬項 Lack of fit1.0030.330.222 5純誤差 Pure error0.1820.09總和 Sum61.5814

2

.

2

.

4

回歸模型交互效應分析根據二次多項式模型方差分析結果，塵土-塵土滾動摩擦因數與塵土接觸模型JKR表面能系數(

X

X

)、塵土-塵土滾動摩擦因數與塵土-塵土碰撞恢復系數(

X

X

)的

P

<0.05，這2個交互項對顆粒堆休止角影響顯著。利用Design-Expert軟件對塵土-塵土滾動摩擦因數與塵土接觸模型JKR表面能系數(

X

X

)以及塵土-塵土滾動摩擦因數與塵土-塵土碰撞恢復系數(

X

X

)交互作用的響應曲面進行繪制，交互項對休止角的影響見圖5。通過對比響應面坡度趨勢得出，塵土接觸模型JKR表面能系數(

X

)所對應的響應面坡度趨勢最陡，塵土-塵土滾動摩擦因數(

X

)和塵土-塵土碰撞恢復系數(

X

)所對應的響應面坡度較為平緩。這說明塵土接觸模型JKR表面能系數(

X

)對顆粒堆休止角的影響顯著。

圖5 部分因素交互作用對休止角(α2)的影響Fig.5 Effect of interaction of some factors on angle of repose

3 參數優化及試驗驗證

3.1 休止角驗證試驗

以秸稈碎末和塵土顆粒堆物理試驗休止角

α

=41.668°為目標，利用Design-Expert軟件對顯著性參數與數值模擬休止角(

α

)的二次多項式方程(式(5))進行尋優求解，在軟件提供的眾多參數組合中尋找若干組組合，進行數值模擬。最后選擇與物理試驗休止角相對誤差最小的一組組合為數值模擬最優取值，即塵土-塵土滾動摩擦因數0.182、塵土接觸模型JKR表面能系數9.052和塵土-塵土碰撞恢復系數0.52。為驗證上述參數的準確性，以上述參數再次進行3次重復數值模擬，所得數值模擬休止角

α

分別為40.356°、40.741°和42.332°，平均值為41.143°，與物理試驗顆粒堆休止角

α

的相對誤差為1.26%。

3.2 標準篩除塵率驗證試驗

為了進一步驗證數值模擬參數取值的準確性，進行標準篩除塵率對比驗證試驗。標準篩直徑為450 mm，高度150 mm，篩孔直徑分別為3 mm和7 mm。利用標準篩進行除塵物理試驗，將篩出物稱重得出除塵率，重復6次試驗并取平均值，得出物理試驗除塵率為13.24%。利用三維建模軟件Solidworks繪制標準篩三維模型，模型尺寸與實物尺寸相同，具體結構見圖6。

圖6 標準篩三維模型Fig.6 3D model of standard sieve

將三維模型導入EDEM軟件中，在標準篩上添加正弦平移方向的力，模擬手動搖晃標準篩的運動狀態，頻率為10 Hz，振幅8 mm。在標準篩上方設置顆粒工廠，顆粒分別為段狀的秸稈外皮、段狀的秸稈內瓤、塵土、碎末狀的秸稈外皮和碎末狀的秸稈內瓤。生成段狀的秸稈外皮和內瓤模型，模型直徑為3 mm，長度20 mm，質量88 g，塵土顆粒和秸稈碎末顆粒模型不變，質量分別為2 和21 g。設置時間步長為5.026×10s，總時間8 s。數值模擬結束后，在后處理面板的Setup selections選項中選擇Grid bin group，網格位于標準篩下方，設置網格數為1，尺寸為長220 mm、寬170 mm、高43 mm，在數據顯示窗口中選擇顯示統計顆粒質量為塵土顆粒和秸稈碎末顆粒，最后選擇顯示顆粒屬性即顯示從標準篩中篩出的塵土顆粒和秸稈碎末顆粒質量。數值模擬所篩出的顆粒質量為28.44 g，除塵率為14.29%。標準篩數值模擬除塵率與物理試驗除塵率相對誤差為7.93%，再次驗證了數值模擬參數取值的準確性。

4 結 論

本研究采用物理試驗與數值模擬相結合的方法對玉米秸稈飼料的篩出物(塵土與秸稈碎末)進行離散元數值模擬參數標定，分別采用Plackett-Burman試驗、最陡爬坡試驗和Box-Behnken試驗確定對數值模擬休止角影響最顯著的因素以及參數取值，主要結論如下：

1)利用斜面儀測定塵土-外皮滾動摩擦因數為0.2～0.5，塵土-外皮靜摩擦因數0.4～0.7，塵土-內瓤滾動摩擦因數0.3～0.6，塵土-內瓤靜摩擦因數0.6～0.9。

2)通過Design-Expert軟件進行Plackett-Burman試驗和最陡爬坡試驗得到了3個對顆粒堆休止角影響最顯著因素，分別為塵土接觸模型JKR表面能、塵土-塵土滾動摩擦因數和塵土-塵土恢復系數。

3)以物理試驗顆粒堆休止角為目標值，通過Box-Behnken試驗以及方差分析得到了顯著性參數最優解，分別為塵土接觸模型JKR表面能系數9.052 J/m、塵土-塵土滾動摩擦因數0.182、塵土-塵土恢復系數0.52。對優化后的參數進行數值模擬得到休止角為41.143°。

4)對所得出的試驗數據分別進行休止角驗證試驗和標準篩除塵率驗證試驗，數值模擬休止角與物理試驗休止角相對誤差為1.26%，標準篩數值模擬除塵率與物理試驗除塵率相對誤差為7.93%，進一步驗證了數值模擬參數的可靠性。