機(jī)收麻山藥離散元模型構(gòu)建及其仿真參數(shù)標(biāo)定

郝建軍，龍思放，李 浩，賈雅麗，馬志凱，趙建國

機(jī)收麻山藥離散元模型構(gòu)建及其仿真參數(shù)標(biāo)定

郝建軍1，龍思放1，李 浩1，賈雅麗2，馬志凱1，趙建國1

（1. 河北農(nóng)業(yè)大學(xué)機(jī)電工程學(xué)院，保定 071001；2. 河北省農(nóng)業(yè)機(jī)械鑒定監(jiān)理總站，石家莊 050031）

由于麻山藥收獲過程缺乏有效數(shù)值模擬，在很大程度上阻礙了麻山藥收獲機(jī)的設(shè)計(jì)與優(yōu)化。該文測定了麻山藥的密度、長度、徑向尺寸、抗壓、抗彎及抗剪強(qiáng)度，基于離散元法建立了麻山藥雙峰分布模型，并對黏結(jié)參數(shù)進(jìn)行校核；以土壤堆積角為響應(yīng)值，對沙壤土基質(zhì)間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)、滾動摩擦因數(shù)和表面能4個(gè)參數(shù)進(jìn)行標(biāo)定，建立了土壤堆積角與4個(gè)參數(shù)之間的回歸模型并進(jìn)行驗(yàn)證，標(biāo)定了麻山藥與鋼板、沙壤土間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)和滾動摩擦因數(shù)。試驗(yàn)結(jié)果表明，麻山藥雙峰分布模型能夠表征麻山藥的力學(xué)特性，參數(shù)校核得到法向剛度、切向剛度、臨界法向應(yīng)力、臨界切向應(yīng)力及黏結(jié)半徑分別為9.3×105N/m、3.0×106N/m、0.58 MPa、0.14 MPa、3.5 mm；沙壤土基質(zhì)間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)、滾動摩擦因數(shù)及表面能分別為0.42、0.20、0.30、0.40 J/m3，離散元仿真試驗(yàn)后得到的土壤堆積角與試驗(yàn)結(jié)果平均誤差為1.48%；麻山藥與鋼板之間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)及滾動摩擦因數(shù)分別為0.34、0.26、0.049，與沙壤土之間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)及滾動摩擦因數(shù)分別為0.21、0.38、0.075。研究結(jié)果可為麻山藥機(jī)械化收獲及產(chǎn)后加工等仿真試驗(yàn)提供一定的理論參考。

農(nóng)作物；模型；離散元法；麻山藥；參數(shù)標(biāo)定

0 引 言

麻山藥是藥食同源的名貴食材，在中國具有悠久的栽培歷史。其主產(chǎn)于河南、河北、山東、江蘇、湖北、安徽、山西等地。隨著生活水平的提高，富含多種維生素、氨基酸，且兼具食藥兩用的麻山藥及其加工制品備受青睞。麻山藥扎根較深，在地表以下0.6～1.2 m不等，收獲過程中極易劃傷或折斷。在中國，小面積麻山藥種植區(qū)多以人工挖掘收獲為主；較大面積種植區(qū)多采用挖溝機(jī)先在麻山藥兩行間開出一條深溝，然后人工再用鐵鍬挖掘。收獲過程機(jī)械化程度低、費(fèi)時(shí)費(fèi)力，而且極易對麻山藥造成損傷，影響品質(zhì)。

目前針對麻山藥收獲的研究主要集中在物理試驗(yàn)方面，較少涉及數(shù)值模擬研究工作，所以有必要利用便捷的研究方法實(shí)現(xiàn)麻山藥收獲過程可視化、參數(shù)化。EDEM是國際通用的基于離散單元法模擬和分析顆粒系統(tǒng)的CAE軟件，目前已在農(nóng)業(yè)領(lǐng)域得到了廣泛應(yīng)用[1-6]。其中Hertz-Mindlin with bonding黏結(jié)模型能夠模擬農(nóng)作物在收獲、加工過程中的破碎現(xiàn)象，廣泛應(yīng)用于農(nóng)業(yè)物料建模領(lǐng)域[7-8]。國內(nèi)外學(xué)者進(jìn)行參數(shù)標(biāo)定及校核[9-11]，大多是根據(jù)理論公式計(jì)算出BPM模型所需黏結(jié)參數(shù)值，再通過物理試驗(yàn)去驗(yàn)證選定的參數(shù)值，使虛擬試驗(yàn)表觀結(jié)果與物理試驗(yàn)結(jié)果相近似，選擇表觀特征變化最接近的一組參數(shù)作為農(nóng)作物破碎仿真中使用的黏結(jié)參數(shù)。這類方法能夠定性分析物料破碎狀態(tài)，對于物料破碎時(shí)所受外力的定量分析相對不夠精確。本文基于Hertz-Mindlin with bonding黏結(jié)模型建立小白嘴麻山藥雙峰分布模型，并對適宜其種植的沙壤土的物理特性及兩者間的接觸參數(shù)進(jìn)行標(biāo)定，利用物理試驗(yàn)與虛擬試驗(yàn)相結(jié)合的方法校核麻山藥離散元模型以及所標(biāo)定參數(shù)的準(zhǔn)確性。研究結(jié)果可為分析麻山藥在收獲過程中的運(yùn)動、受力狀態(tài)及破損方式提供理論基礎(chǔ)。

1 麻山藥離散元模型建立及參數(shù)校核

1.1 麻山藥力學(xué)特性參數(shù)測定

為保證麻山藥在收獲過程中不受損傷，需分析其在收獲過程中的受力狀態(tài)并建立其離散元模型。為此，要得到麻山藥的三軸尺寸和力學(xué)參數(shù)。本文以河北省安平縣普遍種植的小白嘴麻山藥為研究對象，其物理特性參數(shù)為：長度400～900 mm，徑向尺寸為20～45 mm，平均密度為1 010 kg/m3，平均含水率為80%。采用HY-0580微型機(jī)電子萬能力學(xué)試驗(yàn)機(jī)，通過壓縮、彎曲、剪切試驗(yàn)測量其力學(xué)特性參數(shù)。

麻山藥的彈性模量是求得剪切模量的先決條件。材料在彈性變形階段時(shí)，其應(yīng)力和應(yīng)變成正比關(guān)系，比例系數(shù)則為彈性模量。彈性形變符合胡克定律，壓縮與彎曲彈性模量計(jì)算公式分別為

式中1為壓縮彈性模量，MPa；為載荷，N；為截面面積，mm2；為應(yīng)變；2為彎曲彈性模量，MPa；為載荷，N；為跨距，mm；為彎曲撓度，mm；為截面慣性矩，mm4。經(jīng)試驗(yàn)與計(jì)算得到麻山藥力學(xué)參數(shù)，如表1所示。由表中數(shù)據(jù)可知，麻山藥軸向與徑向彈性模量較為接近，結(jié)構(gòu)較均質(zhì)。

表1 麻山藥力學(xué)試驗(yàn)結(jié)果

1.2 麻山藥根莖結(jié)構(gòu)模型構(gòu)建

組成顆粒群的球形顆粒粒徑分布方式包括單一分布、高斯分布和雙峰分布，不同的顆粒粒徑分布方式會導(dǎo)致顆粒間黏結(jié)鍵分布密度、長短、數(shù)量的不同，如圖1所示。其中單一分布采用球形顆粒等粒徑排布而成，顆粒周圍黏結(jié)鍵分布均勻，對破碎對象的力學(xué)特性不能很好的描述；高斯分布中具有非等粒徑球形顆粒，顆粒間黏結(jié)鍵疏密分布，較單一分布能夠更好的反映破碎對象的力學(xué)特性。而組成雙峰分布的顆粒粒徑服從正態(tài)分布，大顆粒占據(jù)主要空間位置，小顆粒則緊密排布在大顆粒周圍，顆粒具有較高的配位數(shù)，從而提高了顆粒群的填充密度，降低了孔隙率，保證了模型具有很好的黏結(jié)強(qiáng)度[12]。

a. 單一分布a. Mono distributionb. 高斯分布b. Gaussian distributionc. 雙峰分布c. Bimodal distribution

本文使用雙峰分布堆積方式對麻山藥建模（大尺度顆粒：平均半徑3 mm，標(biāo)準(zhǔn)差0.2 mm；小尺度顆粒：平均半徑1 mm，標(biāo)準(zhǔn)差0.2 mm），共9 887個(gè)非等粒徑的顆粒、111 723個(gè)有效黏結(jié)鍵，平均單個(gè)球體配位數(shù)約為11，孔隙率較低，能更好的模擬麻山藥力學(xué)特征，如圖2a所示。參考相關(guān)學(xué)者對研究對象的處理方法[6,13-15]及力學(xué)試驗(yàn)結(jié)果，將麻山藥等效為各向同性結(jié)構(gòu)，泊松比取值0.4[2]。另外，麻山藥根莖表面有直徑小于1 mm、長度0～50 mm不等且與麻山藥主根莖連接、容易斷裂的毛細(xì)根，收獲過程中所表現(xiàn)出與土壤之間的黏結(jié)力可以忽略不計(jì)，在模型中可將其略去。圖2b為添加BPM黏結(jié)鍵后的麻山藥離散元模型正、俯視截面圖。

圖2 麻山藥離散元模型

1.3 黏結(jié)參數(shù)確定及校驗(yàn)

本文根據(jù)麻山藥力學(xué)試驗(yàn)結(jié)果，參考國內(nèi)外相關(guān)學(xué)者研究理論基礎(chǔ)對麻山藥黏結(jié)參數(shù)進(jìn)行標(biāo)定，在離散元軟件中進(jìn)行壓縮和剪切仿真試驗(yàn)，選擇曲線特征變化與實(shí)際力學(xué)試驗(yàn)結(jié)果最為接近的一組參數(shù)，使標(biāo)定結(jié)果更加精確。

1.3.1 BPM黏結(jié)參數(shù)確定

結(jié)合麻山藥力學(xué)特性曲線圖，法向剛度K用式（3）計(jì)算[7]

K=DF/D（3）

式中DF為垂直于麻山藥的載荷增量，N；D為麻山藥產(chǎn)生的形變增量，mm。

前文對麻山藥進(jìn)行力學(xué)試驗(yàn)后得到麻山藥形變與載荷關(guān)系曲線，如圖3a所示。

圖3 麻山藥形變與載荷關(guān)系圖

段為彈性形變階段，隨著負(fù)載的增加，麻山藥形變呈線性增加。截取段進(jìn)行分析，線性關(guān)系如圖3b所示。對形變和負(fù)載進(jìn)行一元線性回歸分析，回歸方程中=99.642，決定系數(shù)2為0.995 6。根據(jù)公式（3）可以得到法向剛度為斜率的103。選取5個(gè)樣品進(jìn)行測量，取的平均值，即得法向剛度K=100 201 N/m。

切向剛度K為法向剛度K的2/3～1，研究中取切向剛度為法向剛度的1倍[16]，故得到法向剛度系數(shù)與切向剛度系數(shù)均為100 201 N/m。由前文軸向壓縮試驗(yàn)可得，平均壓縮臨界載荷為400 N，由抗壓強(qiáng)度公式：/，其中為壓力，為受力面積，為臨界法向應(yīng)力。得：=0.58 MPa。根據(jù)彈性結(jié)構(gòu)體切應(yīng)力計(jì)算公式，得到麻山藥破碎中臨界切向應(yīng)力為：=0.14 MPa。顆粒黏結(jié)半徑一般為顆粒半徑的1.1～1.2倍，本文中取黏結(jié)半徑為3.5 mm。

式中是臨界切向應(yīng)力，Pa；是剪切載荷，N；是受剪截面面積，mm2。

1.3.2 參數(shù)修正試驗(yàn)與結(jié)果

在進(jìn)行虛擬預(yù)試驗(yàn)過程中，虛擬試樣所承受最大壓縮和剪切載荷與實(shí)際值比較存在誤差。因此參考上文理論計(jì)算值數(shù)量級和國內(nèi)學(xué)者研究方法[17]，在EDEM軟件中進(jìn)行虛擬壓縮與剪切修正試驗(yàn)，以最大壓縮載荷及最大剪切力作為響應(yīng)值（分別用1、2表示），對這2個(gè)黏結(jié)參數(shù)進(jìn)行進(jìn)一步修正，結(jié)果表明隨著2個(gè)參數(shù)數(shù)值的增加，麻山藥模型抗壓與抗剪載荷均呈現(xiàn)上升趨勢，且切向剛度值對于麻山藥模型的抗破壞能力影響大于法向剛度值。試驗(yàn)因素編碼值與仿真試驗(yàn)結(jié)果如表2所示。使用State-East公司的Design-Expert V8.0.6.1軟件對試驗(yàn)數(shù)據(jù)進(jìn)行響應(yīng)面分析，結(jié)果如表3、4所示，得到優(yōu)化后的回歸模型分別為

應(yīng)用Design Expert 8.0.6.1軟件對2個(gè)回歸模型以最初設(shè)計(jì)點(diǎn)為起始點(diǎn)求解，使模擬結(jié)果盡可能滿足麻山藥壓縮和剪切力學(xué)試驗(yàn)所得到的力學(xué)曲線特征，最終得到一組最優(yōu)標(biāo)定結(jié)果如表5所示。

表2 仿真試驗(yàn)方案與結(jié)果

表3 Y1模型方差分析

注：**極顯著（<0.01）*顯著（<0.05），下同。

Note: ** is highly significant (<0.01), * is significant (<0.05), the same below.

表4 Y2模型方差分析

表5 BPM黏結(jié)參數(shù)

1.3.3 BPM黏結(jié)參數(shù)校驗(yàn)

利用虛擬標(biāo)定所得的黏結(jié)參數(shù)在EDEM軟件中對長40 mm的仿真試樣進(jìn)行壓縮和剪切試驗(yàn)，將壓縮載荷-形變曲線特征及剪切力-位移曲線特征與實(shí)際力學(xué)試驗(yàn)結(jié)果進(jìn)行比對，校核黏結(jié)參數(shù)。

圖4 麻山藥模型仿真試驗(yàn)

從圖4a看出，麻山藥經(jīng)過軸向壓縮后，最大應(yīng)力集中在表皮中部，斷裂趨勢自左向右傾斜延伸。在虛擬試驗(yàn)中由應(yīng)力云圖也可以看到最大應(yīng)力集中于虛擬試樣表皮中部。兩者破裂形式基本一致。從圖4b看出，麻山藥徑向剪切破壞后，由外表皮至內(nèi)部纖維被切斷。虛擬試驗(yàn)中對試樣進(jìn)行徑向剪切時(shí)，顆粒間黏結(jié)鍵斷裂，最大應(yīng)力集中于刀刃處。兩者切口寬度與刀片厚度相一致，且試樣在虛擬試驗(yàn)中的剪切破壞行為與物理試驗(yàn)基本吻合。由圖5可以看出，仿真試樣的壓縮載荷-形變曲線與相鄰曲線特征基本一致，較平均峰值水平偏低，剪切力-位移曲線峰值稍高于平均水平，原因是法向剛度K與切向剛度K不斷變化的過程中，剪切力隨之變化的靈敏程度較壓縮載荷高，但兩曲線整體的特征趨勢與力學(xué)試驗(yàn)結(jié)果基本相符。綜上所述，經(jīng)過以上步驟所得的麻山藥最優(yōu)黏結(jié)參數(shù)所表現(xiàn)出來的力學(xué)特征值與麻山藥實(shí)際情況較為接近，能夠表征麻山藥的基本力學(xué)特性。

圖5 麻山藥壓縮、剪切試驗(yàn)載荷-位移曲線

2 沙壤土參數(shù)標(biāo)定

2.1 試驗(yàn)材料與方法

采用土壤堆積角試驗(yàn)法對沙壤土基質(zhì)間的接觸參數(shù)進(jìn)行標(biāo)定，仿真模型如圖6所示。由于麻山藥種植于沙壤土中，顆粒之間存在一定的黏附關(guān)系，選用Hertz-Mindlin with JKR接觸模型[18]作為土壤基質(zhì)間的接觸模型。試驗(yàn)中所用土壤取自河北省安平縣麻山藥種植田，土質(zhì)松散、無結(jié)塊現(xiàn)象、粒徑大部分在4 mm以下，屬于沙質(zhì)土壤。稱量計(jì)算可知土樣平均容重1 362 kg/m3。在烘干箱中烘干土樣直至質(zhì)量恒定，稱量計(jì)算得到平均含水率為14.15%，仿真參數(shù)見表6[5-6,19-20]。土壤堆積試驗(yàn)中所用漏斗進(jìn)料口徑為100 mm，出料口徑為16 mm，出料口與堆積底板間的高度差為110 mm。試驗(yàn)時(shí)，漏斗中的土壤顆粒經(jīng)漏斗出料口落于堆積底座上，在側(cè)面對堆積體進(jìn)行拍照，使用MATLAB將獲取的圖像進(jìn)行灰度化、二值化處理，并使用Canny算子進(jìn)行邊緣檢測，再將圖像導(dǎo)入Origin軟件獲取輪廓點(diǎn)坐標(biāo)，得到散點(diǎn)圖之后，用高斯分布對其進(jìn)行擬合，找到圖像拐點(diǎn)。然后取拐點(diǎn)切線與軸所夾銳角為土壤的堆積角[21]，試驗(yàn)過程如圖7a、b、c所示。

1.顆粒工廠 2.漏斗 3.計(jì)算域 4.堆積底板

表6 仿真參數(shù)

圖7 土壤顆粒休止角計(jì)算過程

2.2 參數(shù)標(biāo)定試驗(yàn)與結(jié)果

接觸模型參數(shù)中沒有顯著性影響的參數(shù)并不能基于堆積角來標(biāo)定。本文參照相關(guān)文獻(xiàn)[22-26]的標(biāo)定結(jié)果可知，顆粒間碰撞恢復(fù)系數(shù)、顆粒間靜摩擦因數(shù)、顆粒間滾動摩擦因數(shù)以及顆粒間表面能4個(gè)因素對于堆積角的影響是顯著的，所以標(biāo)定結(jié)果可靠。進(jìn)行堆積角Box-Behnken試驗(yàn)，試驗(yàn)因素為顆粒間碰撞恢復(fù)系數(shù)、顆粒間靜摩擦因數(shù)、顆粒間滾動摩擦因數(shù)和顆粒表面能，分別對應(yīng)、、、。設(shè)3組0水平試驗(yàn)，共進(jìn)行27組試驗(yàn)，試驗(yàn)因素編碼值與仿真試驗(yàn)結(jié)果如表7所示。使用State-East公司的Design-Expert V8.0.6.1軟件對試驗(yàn)數(shù)據(jù)進(jìn)行響應(yīng)面分析，結(jié)果如表8所示，得到回歸模型為

表7 堆積角仿真試驗(yàn)方案與結(jié)果

模型顯著性檢驗(yàn)<0.000 1，決定系數(shù)2=0.990 9，失擬項(xiàng)=0.207 5，說明回歸模型極顯著，失擬項(xiàng)不顯著，擬合程度很高。由表8可知，對土壤堆積角的影響，、、、、、2極顯著，顯著、影響顯著順序?yàn)?>2>>>>。試驗(yàn)因素與堆積角之間存在二次非線性關(guān)系和交互作用。

對回歸模型進(jìn)行優(yōu)化，以實(shí)際堆積角平均值為優(yōu)化目標(biāo)，顆粒間碰撞恢復(fù)系數(shù)、顆粒間靜摩擦因數(shù)、顆粒間滾動摩擦因數(shù)、顆粒間表面能作為變量，建立參數(shù)優(yōu)化數(shù)學(xué)模型，使用MATLAB軟件的優(yōu)化工具箱進(jìn)行求解。結(jié)果表明當(dāng)顆粒間碰撞恢復(fù)系數(shù)為0.42，顆粒間靜摩擦因數(shù)為0.20，顆粒間動摩擦因數(shù)為0.30，顆粒間表面能為0.40J/m3時(shí)，仿真得到土壤堆積角數(shù)值為41.25°與實(shí)際堆積角平均值41.87°基本一致，誤差為1.48%，如圖8所示。

圖8 優(yōu)化解下的堆積角仿真試驗(yàn)與物理試驗(yàn)結(jié)果

表8 通用旋轉(zhuǎn)組合設(shè)計(jì)回歸方程方差分析

3 麻山藥與鋼板及沙壤土間接觸參數(shù)的標(biāo)定

麻山藥根莖較長，在試驗(yàn)過程中截取形態(tài)較為均勻的圓柱形根莖段進(jìn)行試驗(yàn)。在EDEM軟件中使用顆粒填充法建立麻山藥顆粒仿真模型[27], 仿真模型與實(shí)物對比如圖9a所示。將安平縣沙壤土均勻填入長500 mm、寬400 mm、高15 mm的無上蓋鐵盒中，并輕微壓實(shí)，得到沙壤土土床。在EDEM軟件中，使用上一節(jié)所測得沙壤土相關(guān)參數(shù)建立顆粒床，待顆粒完全靜止后用平面輕微壓實(shí)顆粒床，保證壓實(shí)后顆粒床高度與實(shí)際土床高度基本一致，如圖9b所示。

3.1 麻山藥碰撞恢復(fù)系數(shù)測定

碰撞恢復(fù)系數(shù)的定義為碰撞后法向相對分離速度與碰撞前相對接近速度的比值[26,28]，碰撞恢復(fù)系數(shù)的計(jì)算如式

式中H為麻山藥下落高度，m；h為麻山藥彈起高度，m。

為保證麻山藥在下落后不會破裂同時(shí)彈起高度易于分辨，進(jìn)行預(yù)試驗(yàn)確定麻山藥下落最佳高度為700 mm。使用高速攝像機(jī)對試驗(yàn)過程錄像，應(yīng)用卷尺及鋼板尺測量麻山藥與鐵板碰撞后彈起的高度，測量精度為5 mm。試驗(yàn)后，將錄像視頻進(jìn)行0.1倍慢放，確定麻山藥塊彈起高度。分別試驗(yàn)10次后取平均值，得到麻山藥-鋼板之間的平均碰撞恢復(fù)系數(shù)為0.34，麻山藥-沙壤土之間的平均碰撞恢復(fù)系數(shù)為0.21。在EDEM軟件中進(jìn)行仿真驗(yàn)證試驗(yàn)，麻山藥仿真模型下落碰撞鋼板、顆粒床后彈起高度與實(shí)際試驗(yàn)彈起高度平均值誤差僅分別為0.64%、3.7%。

3.2 靜摩擦因數(shù)仿真標(biāo)定

利用斜面滑動法進(jìn)行靜摩擦因數(shù)的測定。使用數(shù)顯角度尺測量傾斜角度，精度為0.05°，試驗(yàn)測得麻山藥在鋼板上產(chǎn)生滑動的平均角度為13.45°，在沙壤土土床上產(chǎn)生滑動的平均角度為20.30°。仿真過程中將滾動摩擦因數(shù)設(shè)為0，待麻山藥顆粒完全靜止后，鋼板與沙壤土土床分別以0.5 °/s的速度定軸轉(zhuǎn)動，直至麻山藥顆粒開始滑動。按一定梯度調(diào)整靜摩擦因數(shù)，得到仿真試驗(yàn)結(jié)果如表9所示。建立擬合方程分別為

式中1、2分別表示麻山藥顆粒在鋼板、沙壤土表面開始產(chǎn)生滑動角度，（°）；1、2分別表示麻山藥-鋼板、麻山藥-沙壤土之間的靜摩擦因數(shù)。

兩擬合方程決定系數(shù)2分別為0.999 9、0.999 3，擬合程度較高。將1=13.45，2=20.30代入式中，得到麻山藥-鋼板間、麻山藥-沙壤土間的靜摩擦因數(shù)分別為1=0.26和2=0.38。進(jìn)行仿真試驗(yàn)，得到仿真試樣在鋼板、顆粒床上開始產(chǎn)生滑動的平均角度分別為13.60°和20.85°，與實(shí)際試驗(yàn)誤差分別為1.1%、2.6%。

表9 靜摩擦仿真試驗(yàn)方案與結(jié)果

3.3 滾動摩擦因數(shù)仿真標(biāo)定

麻山藥在鋼板和沙壤土土床表面的滾動將產(chǎn)生滾動摩擦，試驗(yàn)中麻山藥試樣從一定傾斜角度的鋼制平板上某一位置以初速度為0沿著斜面向下滾動（為減小土壤表面不夠平滑所帶來的試驗(yàn)誤差，2次試驗(yàn)均用鋼板作為傾斜板），最終分別滾落至水平放置的鋼板和沙壤土土床上，滾動一段距離之后靜止，采用卷尺測定麻山藥塊在水平鋼板和沙壤土床表面的滾動距離，測量精度為5 mm，以滾動距離為評價(jià)指標(biāo)，對麻山藥-鋼板和麻山藥-沙壤土間的滾動摩擦因數(shù)進(jìn)行測定。為保證麻山藥在滾動過程中不會停置、產(chǎn)生彈跳以及滾出測試區(qū)域等情況，通過預(yù)滾動試驗(yàn)確定了斜置鋼板傾斜角度為17°，在斜置鋼板上滾動距離分別設(shè)為1=300 mm，2=200 mm。假設(shè)麻山藥試樣為一顆理想的圓柱體[3,23]，將試驗(yàn)過程認(rèn)定為純滾動，不考慮靜摩擦力的影響，由能量守恒定律得

式中為鋼板傾斜角度，(°)；為麻山藥塊在傾斜剛板上的滾動距離，m；1為麻山藥與鋼板間的滾動摩擦因數(shù)，μ為麻山藥與沙壤土間的滾動摩擦因數(shù)；L為麻山藥在鋼板上的滾動距離，m；2為麻山藥在沙壤土表面的滾動距離，m。滾動摩擦試驗(yàn)結(jié)果如表10所示。將表中數(shù)據(jù)代入公式得到麻山藥與鋼板、麻山藥與沙壤土之間的平均滾動摩擦因數(shù)分別為0.104、0.125。

表10 滾動摩擦試驗(yàn)結(jié)果

沿傾斜鋼板面向上距離底端1=300 mm、2=200 mm處進(jìn)行麻山藥顆粒仿真驗(yàn)證試驗(yàn)。仿真過程中發(fā)現(xiàn)，麻山藥在鋼板、土顆粒床上的平均滾動距離分別為549、305 mm，相比實(shí)際試驗(yàn)滾動距離平均值偏小。隨后進(jìn)行滾動摩擦仿真標(biāo)定試驗(yàn)，按一定梯度調(diào)整滾動摩擦因數(shù)，試驗(yàn)結(jié)果如表11所示，建立擬合方程分別為

式中3、4分別表示麻山藥顆粒在鋼板、沙壤土表面的滾動距離，m；3、4分別表示麻山藥顆粒-鋼板、麻山藥顆粒-沙壤土間的滾動摩擦因數(shù)。

表11 麻山藥顆粒滾動摩擦仿真試驗(yàn)方案與結(jié)果

方程決定系數(shù)2分別為0.986 3、0.998 8，擬合方程較為準(zhǔn)確可靠。將實(shí)際試驗(yàn)滾動距離平均值分別代入兩方程中，得到滾動摩擦因數(shù)3=0.049，4=0.075。對得到的參數(shù)進(jìn)行仿真驗(yàn)證，重復(fù)5次試驗(yàn)得到平均滾動距離分別為558、320 mm，與實(shí)際試驗(yàn)滾動距離誤差分別為1.60%、4.91%。

通過仿真標(biāo)定試驗(yàn)得到麻山藥與鋼板、沙壤土間仿真滾動摩擦因數(shù)1=0.049、2=0.075，與實(shí)際試驗(yàn)滾動距離誤差僅為1.60%和4.91%；而根據(jù)能量守恒試驗(yàn)法測得1=0.104、2=0.125，與實(shí)際試驗(yàn)相對誤差較大，分別為112.25%、67%。原因可能是：真實(shí)麻山藥并不是理想的圓柱體，同時(shí)所選取的麻山藥顆粒表皮以及鋼板表面特性與仿真模型有差異；在實(shí)際試驗(yàn)條件下麻山藥重力勢能的損耗并不完全來自滾動摩擦，還包括顆粒-鋼板、沙壤土所產(chǎn)生的靜摩擦以及顆粒在滾動過程中與鋼板、沙壤土碰撞所損耗的能量，而假設(shè)理想試驗(yàn)將這些能量消耗所忽略，故出現(xiàn)以能量守恒法所測值與實(shí)際試驗(yàn)值相對誤差較高的情況。因此，在進(jìn)行滾動摩擦因數(shù)測定時(shí)，要充分考慮試驗(yàn)對象的特殊性，選用適宜的方法進(jìn)行測定。

4 結(jié) 論

1）通過試驗(yàn)測得麻山藥根莖的含水率、密度、彈性模量等本征參數(shù)，發(fā)現(xiàn)麻山藥徑向與軸向彈性模量值相差較小，表明其結(jié)構(gòu)較為均質(zhì)；根據(jù)所測參數(shù)建立了麻山藥離散元雙峰分布模型，并從表觀和力學(xué)性能2個(gè)方面進(jìn)行了驗(yàn)證，表明雙峰分布模型能夠較好的表征麻山藥物理特性。

2）確定了安平縣沙壤土堆積角度，結(jié)合Box-Behnken響應(yīng)面分析試驗(yàn)建立并優(yōu)化了麻山藥種植區(qū)的沙壤土堆積角與顆粒間碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)、滾動摩擦因數(shù)及表面能4個(gè)顯著性參數(shù)之間的回歸模型，得到最優(yōu)參數(shù)組合依次為0.42、0.20、0.30、0.40J/m3，為土壤模型的建立及麻山藥與土壤接觸參數(shù)的標(biāo)定提供參考。

3）在EDEM軟件中建立了麻山藥仿真試樣、鋼板及顆粒床模型，通過將實(shí)際試驗(yàn)與數(shù)值模擬試驗(yàn)相結(jié)合，得到了麻山藥與鋼板的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)、滾動摩擦因數(shù)分別為0.34、0.26、0.049，與沙壤土之間的碰撞恢復(fù)系數(shù)、靜摩擦因數(shù)、滾動摩擦因數(shù)分別為0.21、0.38、0.075。為后續(xù)建立麻山藥-沙壤土復(fù)合體模型提供了準(zhǔn)確的數(shù)值參考，進(jìn)而使麻山藥收獲過程實(shí)現(xiàn)可視化、參數(shù)化；同時(shí)標(biāo)定結(jié)果可為麻山藥產(chǎn)后加工的相關(guān)仿真試驗(yàn)及麻山藥產(chǎn)業(yè)全程機(jī)械化進(jìn)程的推進(jìn)提供理論基礎(chǔ)。

[1] 王國強(qiáng)，郝萬軍，王繼新. 離散單元法及其在EDEM上的實(shí)踐[M].西安：西北工業(yè)大學(xué)出版社，2010.

[2] 高國華，謝海峰，王天寶. 設(shè)施蔬菜收獲機(jī)拉拔力學(xué)性能EDEM仿真與試驗(yàn)[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2017，33(23)：24－31.

Gao Guohua, Xie Haifeng, Wang Tianbao. EDEM simulation and experiment of pullout force of protected vegetable harvester[J]. Transactions of the Chinese Society of Agricultural Engineering, 2017, 33(23): 24－31. (in Chinese with English abstract)

[2] 劉文政，何進(jìn)，李洪文，等. 基于離散元的微型馬鈴薯仿真參數(shù)標(biāo)定[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2018，49(5)：132－142，149.

Liu Wenzheng, He Jin, Li Hongwen, et al. Calibration of simulation parameters for potato minituber based on EDEM[J]. Transactions of the Chinese Society for Agricultural Machinery, 2018, 49(5): 132－142, 149. (in Chinese with English abstract)

[3] 原建博，李驊，吳崇友，等. 基于離散單元法的水稻籽?？焖兕w粒建模研究[J]. 南京農(nóng)業(yè)大學(xué)學(xué)報(bào)，2018，41(6)：187－194.

Yuan Jianbo, Li Hua, Wu Chongyou, et al. Study on apace particle modeling of rice grain basis on the discrete element method[J]. Journal of Nanjing Agricultural University, 2018, 41(6): 187－194. (in Chinese with English abstract)

[4] 張濤，劉飛，趙滿全，等. 大豆種子與排種器接觸物理參數(shù)的測定與離散元仿真標(biāo)定[J]. 中國農(nóng)業(yè)大學(xué)學(xué)報(bào)，2017，22(9)，86－92.

Zhang Tao, Liu Fei, Zhao Manquan, et al. Measurement of physical parameters of contact between soybean seed and seed metering device and discrete element simulation calibration[J]. Journal of China Agricultural University, 2017, 22(9): 86－92. (in Chinese with English abstract)

[5] 程一啟，張兆國，張丹，等. 三七收獲機(jī)組合式挖掘鏟減阻效果研究：基于EDEM[J]. 農(nóng)機(jī)化研究，2019，41(2)：55－60.

Cheng Yiqi, Zhang Zhaoguo, Zhang Dan, et al. Research of resistance reduction effect of panax notoginseng harvest combined digging shovel: Based on EDEM[J]. Journal of Agricultural Mechanization Research, 2019, 41(2): 55－60. (in Chinese with English abstract)

[6] 王笑丹，王洪美，韓云秀，等. 基于離散元法的牛肉咀嚼破碎模型構(gòu)建[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2016，32(4)：228－234.

Wang Xiaodan, Wang Hongmei, Han Yunxiu, et al. Structure of beef chewing model based on discrete element method[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2016, 32(4): 228－234. (in Chinese with English abstract)

[7] 江濤，吳崇友，湯慶，等. 基于ANSYS和EDEM的小麥莖稈切割仿真研究[J]. 江蘇農(nóng)業(yè)科學(xué)，2018，46(17)：231－234.

[8] Quist J, Evertsson C M. Cone crusher modelling andsimulation using DEM[J]. Minerals Engineering, 2016, 85: 92－105.

[9] 徐金明，謝芝蕾，賈海濤. 石灰?guī)r細(xì)觀力學(xué)特性的顆粒流模擬[J]. 巖土力學(xué)，2010，31(增刊2)：390－395.

Xu Jinming, Xie Zhilei, Jia Haitao. Simulation of mesomechanical properties of limestone using particle flow code[J]. Rock and Soil Mechanics, 2010, 31(Supp.2): 390－395. (in Chinese with English abstract)

[10] Quist J. Cone Crusher Modelling and Simulation[D]. Gothenburg: Chalmers University of Technology, 2012.

[11] Mishra B K. A review of computer simulation of tumbling mill by the discrete element method: Part I-contact mechanics[J]. International Journal of Mineral Processing, 2003, 71(1/2/3/4): 73－93.

[12] 張濤，劉飛，趙滿全，等. 基于離散元的排種器排種室內(nèi)玉米種群運(yùn)動規(guī)律[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2016，32(22)：27－35.

Zhang Tao, Liu Fei, Zhao Manquan, et al. Movement law of maize population in seed room of seed metering device based on discrete element method[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2016, 32(22): 27－35. (in Chinese with English abstract)

[13] 賴慶輝，袁海闊，胡子武，等. 滾筒板齒式三七種苗分離裝置結(jié)構(gòu)設(shè)計(jì)與試驗(yàn)[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2018，49(4)：121－129.

Lai Qinghui, Yuan Haikuo, Hu Ziwu, et al. Design and experiment on seedling separation device of panax notoginseng seedlings based on roller zigzag mechanism[J]. Transactions of the Chinese Society for Agricultural Machinery, 2018, 49(4): 121－129. (in Chinese with English abstract)

[14] Yu Y, Fu H, Yu J. DEM-based simulation of the corn threshing process[J]. Advanced Powder Technology, 2015, 26(5): 1400－1409.

[15] 張鋒偉，宋學(xué)鋒，張雪坤，等. 玉米秸稈揉絲破碎過程力學(xué)特性仿真與試驗(yàn)[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2019，35(9)：58－65.

Zhang Fengwei, Song Xuefeng, Zhang Xuekun, et al. Simulation and experiment on mechanical characteristics of kneading and crushing process of corn straw[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2019, 35(9): 58－65. (in Chinese with English abstract)

[16] 霍麗麗，孟海波，田宜水，等. 粉碎秸稈類生物質(zhì)原料物理特性試驗(yàn)[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2012，28(11)：189－194.

Huo Lili, Meng Haibo, Tian Yishui, et al. Experimentalstudy on physical property of smashed crop straw[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2012, 28(11): 189－194. (in Chinese with English abstract)

[17] 劉凡一. 清選裝置中小麥顆粒和短莖稈離散元建模研究[D]. 楊凌：西北農(nóng)林科技大學(xué)，2018.

Liu Fanyi. Discrete Element Modelling of the Wheat Particles and Short Straw in Cleaning Devices[D]. Yangling: Northwest A&F University, 2018. (in Chinese with English abstract)

[18] 姚宗路，歐陽雙平，孟海波，等. 顆粒狀秸稈物料流動特性試驗(yàn)[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2012，43(7)：112－116.

Yao Zonglu, Ouyang Shuangping, Meng Haibo, et al. Flow characterization of biomass particle straw[J]. Transactions of the Chinese Society for Agricultural Machinery, 2012, 43(7): 112－116. (in Chinese with English abstract)

[19] 畢秋實(shí)，王國強(qiáng)，陳立軍，等. 基于離散元-多體動力學(xué)聯(lián)合仿真的機(jī)械式挖掘機(jī)挖掘阻力仿真與試驗(yàn)[J]. 吉林大學(xué)學(xué)報(bào)：工學(xué)版，2019，49(1)：106－116.

Bi Qiushi, Wang Guoqiang, Chen Lijun, et al. Numerical simulation and experiment on excavation resistance of mechanical excavator based on DEM-MBD co-simulation[J]. Journal of Jilin University (Engineering and Technology Edition), 2019, 49(1): 106－116. (in Chinese with English abstract)

[20] 張銳，韓佃雷，吉巧麗，等. 離散元模擬中沙土參數(shù)標(biāo)定方法研究[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2017，48(3)：49－56，343－350.

Zhang Rui, Han Dianlei, Ji Qiaoli, et al. Calibration methods of sandy soil parameters in simulation of discrete element method[J]. Transactions of the Chinese Society for Agricultural Machinery, 2017, 48(3): 49－56, 343－350. (in Chinese with English abstract)

[21] 羅帥，袁巧霞，GOUDA Shaban，等. 基于JKR粘結(jié)模型的蚯蚓糞基質(zhì)離散元法參數(shù)標(biāo)定[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2018，49(4)：343－350.

Luo Shai, Yuan Qiaoxia, GOUDA Shaban, et al. Parameters calibration of vermicomposting nursery substrate with discrete element method based on JKR contact model[J]. Transactions of the Chinese Society for Agricultural Machinery, 2018, 49(4): 343－350. (in Chinese with English abstract)

[22] 夏鵬，李郁，楊公波. 散粒物料堆積角離散元仿真研究[J].起重運(yùn)輸機(jī)械，2015(2)：107－110.

Xia Peng, Li Yu, Yang Gongbo. Study on the discrete element simulation of scattered materials[J]. Hoisting and Conveying Machinery, 2015(2): 107－110. (in Chinese with English abstract)

[23] 韓燕龍，賈富國，唐玉榮，等. 顆粒滾動摩擦系數(shù)對堆積特性的影響[J]. 物理學(xué)報(bào)，2014，63(17)：4501－4507.

Han Yanlong, Jia Fuguo, Tang Yurong, et al. Influence of granular coefficient of rolling friction on accumulation characteristics[J]. Acta Physica Sinica, 2014, 63(17): 4501－4507. (in Chinese with English abstract)

[24] 賈旭光，陳曦，李鑫. 不同粒度和堆載形態(tài)下散體瞬時(shí)自然安息角的實(shí)驗(yàn)研究[J]. 現(xiàn)代礦業(yè)，2015，549(1)：25－27.

Jia Xuguang, Chen Xi, Li Xin. Experiment research on instantaneous natural repose angle of granular slope under different granularity and stack form[J]. Modern Mining, 2015, 549(1): 25－27. (in Chinese with English abstract)

[25] 德頤姆科技有限公司. 含水率：物料參數(shù)中JKR表面能的設(shè)置[Z]. 深圳：EDEMBulkSim官方微博，2018[2018-04-09]. https://weibo.com/ttarticle/p/show?id= 2309404226950502699919.

[26] González-Montellano C, Fuentes J M, Ayuga-Téllez E, et al. Determination of the mechanical propertiesof maize grains and olives required for use in DEM simulations[J]. Journal of Food Engineering, 2012, 111(4): 553－562.

[27] 宋學(xué)鋒. 基于離散元法的玉米秸稈揉絲機(jī)工作過程研究[D]. 蘭州：甘肅農(nóng)業(yè)大學(xué)，2018.

Song Xunfeng. Study on the Operation Process of the Corn Straw Kneading Machine Based on Discrete Element Method[D]. Lanzhou: Gansu Agricultural University, 2018. (in Chinese with English abstract)

[28] Zhou Y C, Xu B H, Yu A B, et al. Numerical investigation of the angle of repose of monosized spheres[J]. Physical Review E, 2001, 64(2Pt1): 021301.

Development of discrete element model and calibration of simulation parameters for mechanically-harvested yam

Hao Jianjun1, Long Sifang1, Li Hao1, Jia Yali2, Ma Zhikai1, Zhao Jianguo1

(1.,,071001,; 2.,050031,)

As the parameters of discrete element model for yam harvesting are difficult to obtain directly from literatures or experiments, this study proposes a parameter calibration method that combines experiments with numerical simulations. Firstly, a physical test was carried out. Compression, shear and bending tests were carried out on yam with the moisture content of 80%. The axial and radial elastic modulus of yam were obtained, which were 4.48 and 4.03 MPa, respectively. The results showed that the axial and radial elastic modulus were close, and the structure of yam was homogeneous. At the same time, the density, length and radial dimension of yam were measured. Based on the discrete element method, the bimodal distribution model of yam was established by the particle filling method. According to the research theories of domestic and foreign scholars, the bonding parameters of yam were calibrated. Through the compression and shear simulation tests in the discrete element software, the calibration parameters were checked and verified to make the calibration results more accurate. The experimental results showed that there were 9 887 non-equal size particles and 111 723 effective bonds in the bimodal distribution model of yam, the coordination number was about 11, the porosity was low, which could better characterize the mechanical properties of yam. The normal stiffness, shear stiffness, critical normal stress, critical shear stress and bonding radius were 9.3×105N/m, 3.0×106N/m, 0.58 MPa, 0.14 MPa, 3.5 mm. Secondly, the contact parameters of sand loam particles were calibrated by the soil accumulation angle test. The particle size of sandy loam was less than 4 mm, and the moisture content was 14.15%. After screening, the particles of sandy loam met the requirement of uniformity. The method of funnel test was used to form a certain stacking angle of sandy loam and the stacking body was photographed on the side. The acquired image was grayed and binarized by MATLAB and the edge was detected by Canny operator. Then the image was imported into Origin software to obtain the contour point coordinates. The scatter map was fitted with Gaussian distribution to find the inflexion point of the image and the sharp angle between the tangent of inflection point andaxis as the accumulation angle of soil was took. Then, several tests were carried out to determine the accumulation angle of sandy loam. Combined with the Box-Behnken response surface analysis test, the regression model between the sand soil accumulation angle and coefficient of restitution, coefficient of static friction, coefficient of rolling friction and surface energy was established and optimized. The optimal parameter combination was 0.42, 0.20, 0.30 and 0.40J/m3. The results showed that the average value of simulation soil accumulation angle was 41.25°, the average value of actual soil accumulation angle was 41.87°, and the error between them was only 1.48%. Finally, the simulation sample of yam, steel plate and particle bed model were established in EDEM software. Through the test of coefficient of restitution, static friction coefficient and coefficient of rolling friction, the coefficient of restitution, static friction coefficient and coefficient of rolling friction between yam-steel plate and yam-sandy loam were 0.34, 0.26, 0.049, 0.21, 0.38 and 0.075, the calibration results was reliable and accurate compared with the actual experimental errors of less than 5%. It provides accurate contact parameters for the subsequent establishment of composite model of yam and sandy loam, thereby we can visualize and parameterize the harvest process of yam; meanwhile, the calibration results can provide theoretical basis for the relevant simulation test of yam post-production processing and the promotion of the whole mechanization process of yam industry.

crops; models; discrete element method; yam; calibration of parameters

郝建軍，龍思放，李 浩，賈雅麗，馬志凱，趙建國. 機(jī)收麻山藥離散元模型構(gòu)建及其仿真參數(shù)標(biāo)定[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2019，35(20)：34－42.doi：10.11975/j.issn.1002-6819.2019.20.005 http://www.tcsae.org

Hao Jianjun, Long Sifang, Li Hao, Jia Yali, Ma Zhikai, Zhao Jianguo. Development of discrete element model and calibration of simulation parameters for mechanically-harvested yam[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2019, 35(20): 34－42. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2019.20.005 http://www.tcsae.org

2019-08-15

2019-09-30

河北省現(xiàn)代農(nóng)業(yè)創(chuàng)新團(tuán)隊(duì)建設(shè)項(xiàng)目（HBCT2018090206）

郝建軍，教授，博士生導(dǎo)師，主要從事農(nóng)機(jī)裝備設(shè)計(jì)與制造研究。Email：hjjpaper@163.com

10.11975/j.issn.1002-6819.2019.20.005

S539; S313

A

1002-6819(2019)-20-0034-09