基于離散元的包膜肥料Bonding模型參數(shù)標(biāo)定

都 鑫 劉彩玲 姜 萌 袁 昊 戴 磊 李方林

(中國農(nóng)業(yè)大學(xué)工學(xué)院， 北京 100083)

0 引言

肥料在農(nóng)業(yè)生產(chǎn)中有著不可替代的作用，中國化肥使用量大但肥料利用率低[1]。為保證農(nóng)業(yè)生產(chǎn)的可持續(xù)發(fā)展，養(yǎng)分效率高、環(huán)境友好型的控釋肥料得到廣泛關(guān)注[2-6]。控釋肥料采用聚合物包膜，可定量控制肥料養(yǎng)分釋放量和釋放期，使養(yǎng)分供應(yīng)與作物各生育期需肥規(guī)律吻合，節(jié)肥增效作用顯著[7]。

包膜肥料的養(yǎng)分釋放特性與包膜層材料、結(jié)構(gòu)緊密相關(guān)[8]，排肥器在排施肥料顆粒過程中會對部分肥料顆粒造成不同程度的機(jī)械損傷，導(dǎo)致包膜層破壞進(jìn)而影響肥料養(yǎng)分釋放特性。為設(shè)計(jì)適用于包膜肥料無損排施的排肥器，考慮包膜肥料與排肥器、包膜肥料間的相互擠壓力、摩擦力等因素，改進(jìn)排肥方式、結(jié)構(gòu)參數(shù)等。其中肥料顆粒間的相互作用力難以直接測得，擬采用離散元方法對此過程進(jìn)行數(shù)值模擬[9]。

自CUNDALL等[10]提出適應(yīng)于土力學(xué)的離散元方法后，該方法已廣泛應(yīng)用于散粒體物料的混料、攪拌、輸送、充填、破碎等過程[11-14]。POTYONDY等[15]提出粘合粒子模型(BPM)，將許多小顆粒通過Bond鍵粘結(jié)形成被破碎的物體，在外力作用下Bond鍵發(fā)生形變，當(dāng)力或力矩達(dá)到一定程度時(shí)，粘結(jié)鍵斷裂，小顆粒互相分開發(fā)生破碎。目前物料破碎數(shù)值模擬被廣泛應(yīng)用于巖石鉆探、土壤挖掘、制沙、磨粉等[16-20]研究中，但在農(nóng)業(yè)領(lǐng)域中應(yīng)用較少，尚未見到包膜肥料破碎數(shù)值模擬方面的相關(guān)報(bào)道。

為保證數(shù)值模擬可真實(shí)反映肥料顆粒擠壓破碎情況，需對包膜肥料顆粒破碎模型進(jìn)行參數(shù)標(biāo)定。王云霞等[21]通過堆積角試驗(yàn)標(biāo)定玉米種子的種間靜摩擦因素和滾動摩擦因數(shù)；鹿芳媛等[22]基于摩擦角試驗(yàn)標(biāo)定了不同含水率下的水稻芽種主要接觸參數(shù)；劉彩玲等[23]基于休止角試驗(yàn)系統(tǒng)標(biāo)定了大顆粒尿素的仿真邊界參數(shù)；袁全春等[24]基于JKR模型標(biāo)定了有機(jī)肥散體顆粒離散元參數(shù)；于慶旭等[25]通過碰撞彈跳試驗(yàn)、斜面滑移試驗(yàn)和斜面滾動試驗(yàn)測定三七種子與ABS塑料間接觸參數(shù)，基于響應(yīng)面優(yōu)化方法得到三七種子間接觸參數(shù)。

本文在國內(nèi)外專家學(xué)者的研究基礎(chǔ)上[26-30]，以包膜肥料顆粒為研究對象，采用單軸壓縮試驗(yàn)測得包膜肥料顆粒的極限破碎位移和極限破碎載荷，通過建立與真實(shí)試驗(yàn)相同的仿真模型，利用Placket-Burman試驗(yàn)、Steepest ascent試驗(yàn)和Box-Behnken試驗(yàn)逐步減少因素?cái)?shù)量、縮小因素范圍，最終優(yōu)化得到與真實(shí)試驗(yàn)結(jié)果最接近的Bonding模型參數(shù)組合；為進(jìn)一步驗(yàn)證最優(yōu)參數(shù)組合的正確性，利用外槽輪排肥器排施包膜肥料顆粒，并統(tǒng)計(jì)對比仿真和實(shí)際條件下肥料顆粒破碎率，以期為研究包膜肥料顆粒機(jī)械破碎機(jī)理、優(yōu)化設(shè)計(jì)無損排施的新型排肥器提供參考。

1 包膜肥料顆粒Bonding模型建立

1.1 Bonding模型

POTYONDY等[15]設(shè)計(jì)了Bonding接觸模型以模擬物料破碎。在Bonding接觸模型中，被破碎的物料(顆粒、塊體等)由N個(gè)較小的顆粒構(gòu)建組合而成，小顆粒之間通過Bond鍵固定在一起。Bond鍵具有與有限元方法相近的力學(xué)特性，在外力作用下會發(fā)生形變(拉伸、壓縮、扭轉(zhuǎn))。當(dāng)Bond鍵形變產(chǎn)生的力或力矩達(dá)到一定程度時(shí)，粘結(jié)鍵斷裂，小顆粒互相分開發(fā)生破碎。Bonding接觸模型假設(shè)兩個(gè)顆粒之間的粘結(jié)鍵是一個(gè)虛擬的平板(圓柱)，如圖1所示。

圖1 Bonding模型原理與顆粒粘結(jié)示意圖Fig.1 Bonding model and particle bonding schematic

當(dāng)兩個(gè)顆粒發(fā)生相對運(yùn)動時(shí)，平板(圓柱)承受拉伸、彎曲和剪切，鍵產(chǎn)生力的計(jì)算公式為

(1)

其中

(2)

式中 δFn——法向力，N

δFt——切向力，N

vn——法向速度，m/s

vt——切向速度，m/s

Sn——法向剛度，N/m

St——切向剛度，N/m

δt——仿真時(shí)間步長，s

A——顆粒間的接觸面積，m2

δMn——法向力矩，N·m

δMt——切向力矩，N·m

ωn——法向角速度，rad/s

ωt——切向角速度，rad/s

J——顆粒慣性矩，mm4

RB——粘結(jié)半徑，mm

當(dāng)力和力矩達(dá)到極限值或Fraction顆粒間的距離大于設(shè)定的接觸半徑時(shí)，粘結(jié)鍵發(fā)生斷裂，法向及切向剪切力計(jì)算公式為

(3)

式中σmax——法向剪切力，N

τmax——切向剪切力，N

1.2 肥料顆粒離散元建模

選用山東農(nóng)泱生物科技有限公司生產(chǎn)的貝可施包膜控釋肥料(圖2a)，肥料含水率0.88%，真實(shí)密度1.46 g/cm3，平均三軸尺寸4.08 mm×3.97 mm×3.89 mm，等效直徑為3.98 mm，球形率為0.975，選取長、寬、高與等效直徑相近的肥料顆粒建立輪廓模型。為研究肥料顆粒破碎特性，采用粘接顆粒模型(Bonded particle method, BPM)建立肥料顆粒離散元模型(圖2b)，通過Bond鍵將多個(gè)直徑相等的小球形填充顆粒粘結(jié)成一個(gè)大顆粒。

圖2 肥料顆粒及離散元模型Fig.2 Fertilizer particles and discrete element model

填充顆粒直徑越小，粘結(jié)成一個(gè)大顆粒所需的顆粒數(shù)就越多。如圖3a所示，填充的球形顆粒半徑為0.4 mm，球形顆粒數(shù)量為68，粘結(jié)鍵個(gè)數(shù)為534；如圖3b所示，填充的球形顆粒半徑為0.2 mm，球形顆粒數(shù)量為582，粘結(jié)鍵個(gè)數(shù)為3 611；如圖3c所示，填充的球形顆粒半徑為0.1 mm，球形顆粒數(shù)量為5 581，粘結(jié)鍵個(gè)數(shù)為39 119。

圖3 不同直徑球形顆粒填充圖Fig.3 Diagrams of particles with different diameters

由于小球形顆粒間彼此相互獨(dú)立，在仿真模擬中需要對每個(gè)小球形顆粒的受力、運(yùn)動情況進(jìn)行計(jì)算，因此，顆粒數(shù)量越多，仿真計(jì)算時(shí)間越長，綜合考慮后選取球形顆粒半徑為0.2 mm的粘結(jié)模型進(jìn)行Bonding模型參數(shù)標(biāo)定。

1.3 單軸壓縮試驗(yàn)

1.3.1真實(shí)試驗(yàn)

選取直徑為3.98 mm的肥料顆粒，采用萬能試驗(yàn)儀對肥料顆粒進(jìn)行單軸壓縮試驗(yàn)，以速度0.1 mm/s進(jìn)行加載，直至試樣出現(xiàn)破壞后停止，單軸壓縮試驗(yàn)中試件加載過程如圖4所示。試驗(yàn)重復(fù)20次取平均值，得到肥料顆粒的位移載荷曲線如圖5所示。

圖4 萬能試驗(yàn)儀Fig.4 Universal tester1.壓板 2.底板 3.控制面板

圖5 力-位移曲線Fig.5 Force-displacement curve

由圖5可以看出，肥料顆粒力-位移曲線分為3部分，第1部分為試驗(yàn)開始階段，試樣局部發(fā)生棱角破壞，但由于肥料顆粒球形度較高，力先快速增大后變緩；第2部分為主體發(fā)生破壞的過程，當(dāng)所施加的載荷超過試樣所能承受的極限時(shí)，試樣主體發(fā)生破壞，力達(dá)到最大值，此時(shí)為破壞點(diǎn)；第3部分隨位移增加，力持續(xù)降低，試樣破碎成塊狀。肥料顆粒加載位移為0.45 mm時(shí)加載應(yīng)力最大值為58.668 N，即20粒直徑為3.98 mm的包膜肥料顆粒平均極限破碎載荷為58.668 N，此時(shí)對應(yīng)的極限破碎位移為0.45 mm。

1.3.2仿真試驗(yàn)

采用離散元軟件EDEM 2018建立與實(shí)際試驗(yàn)一致的單軸壓縮仿真模擬試驗(yàn)，如圖6所示。仿真中用到的包膜肥料顆粒本構(gòu)參數(shù)、包膜肥料間和包膜肥料與ABS間的接觸參數(shù)由實(shí)驗(yàn)室標(biāo)定得到，ABS的本構(gòu)參數(shù)如表1所示[23]。EDEM中設(shè)置仿真總時(shí)間為8 s，時(shí)間步長為7.5×10-7s，每隔0.001 s保存1次數(shù)據(jù)，網(wǎng)格尺寸取3.5倍最小球形單元尺寸。

表1 材料物理特性參數(shù)Tab.1 Physical characteristic parameters of material

圖6 單軸壓縮試驗(yàn)仿真圖Fig.6 Simulation diagram of uniaxial compression test

2 Bonding模型參數(shù)標(biāo)定

2.1 Plackett-Burman試驗(yàn)

Plackett-Burman試驗(yàn)通過考察目標(biāo)響應(yīng)與各因子間的關(guān)系，比較各個(gè)因子2水平間的差異性來確定因子顯著性。以極限破碎位移Y1和極限破碎

載荷Y2為目標(biāo)，以單位面積法向剛度X1、單位面積切向剛度X2、臨界法向應(yīng)力X3、臨界切向應(yīng)力X4和粘結(jié)半徑X5為Plackett-Burman試驗(yàn)的5個(gè)因子，每個(gè)因子選取2個(gè)水平。試驗(yàn)選用N=11的Plackett-Burman設(shè)計(jì)，預(yù)留6個(gè)虛擬項(xiàng)作誤差分析，試驗(yàn)因素水平如表2所示。

表2 Plackett-Burman試驗(yàn)因素水平Tab.2 Plackett-Burman test factors and levels

Plackett-Burman試驗(yàn)方案與結(jié)果如表3所示，利用Design-Expert 8.0.6軟件對各因素效應(yīng)分別進(jìn)行方差分析和t檢驗(yàn)，結(jié)果如表4和圖7所示。

表3 Plackett-Burman試驗(yàn)結(jié)果Tab.3 Plackett-Burman test results

由表4和圖7可知，各因素對極限破碎位移Y1影響由大到小為X5、X2、X1、X4、X3，其中X5、X2、X1的P值小于0.05，表明因素對目標(biāo)影響顯著，由帕累托圖的t-value檢驗(yàn)也可得出同樣結(jié)論，且知X5、X1對目標(biāo)的效應(yīng)值為正，X2對目標(biāo)的效應(yīng)值為負(fù)。同理，可知各因素對極限破碎載荷Y2影響由大到小為X5、X1、X4、X2、X3，其中X5、X1、X4的P值小于0.05，表明因素對目標(biāo)影響顯著，由Pareto Chart的t-value檢驗(yàn)也可得出同樣結(jié)論，且知X5、X1、X4對目標(biāo)的效應(yīng)值均為正。

圖7 帕累托圖Fig.7 Pareto Charts

表4 Plackett-Burman試驗(yàn)結(jié)果顯著性分析Tab.4 Analysis of significance of parameters in Plackett-Burman test

各因素與極限破碎位移Y1和極限破碎載荷Y2的模型為

(4)

2.2 Steepest.ascent試驗(yàn)

根據(jù)Plackett-Burman試驗(yàn)結(jié)果，X3對結(jié)果影響不顯著，因此將X3取為5×105Pa。選取因素X1、X2、X4、X5初值分別為6.40×109N/m3、1.92×1010N/m3、8.00×105Pa、0.16 mm，步長分別為1.60×109N/m3、-3.20×109N/m3、2.00×105Pa、0.02 mm，進(jìn)行Steepest ascent試驗(yàn)，以進(jìn)一步找到逼近真實(shí)值的參數(shù)組合，Steepest ascent試驗(yàn)方案和結(jié)果如表5所示。

表5 最速上升試驗(yàn)方案及結(jié)果Tab.5 Steepest ascent test program and results

由1.3.1節(jié)可知，肥料顆粒實(shí)際單軸壓縮極限破碎位移和極限破碎位移載荷分別為0.45 mm、58.668 N，由表5可知，仿真試驗(yàn)中極限破碎位移和載荷與真實(shí)值間的誤差先減小后增大，第3組試驗(yàn)的極限破碎位移Y1和極限破碎載荷Y2與真實(shí)值最接近。

2.3 Box-Behnken試驗(yàn)

為得到最接近真實(shí)值的極限破碎位移Y1和載荷Y2，根據(jù)Steepest ascent試驗(yàn)結(jié)果，選取第2和第4組試驗(yàn)的因素組合為BBD試驗(yàn)的高低水平，試驗(yàn)方案及結(jié)果如表6所示，方差分析結(jié)果如表7所示。

由表7可知，Y1和Y2的二次回歸模型均高度顯著，且失擬項(xiàng)不顯著，表明模型精度高。X1、X2、X4、X5均顯著且對Y1的影響由大到小為X5、X1、X2、X4。X1、X2、X4、X5均顯著且對Y2的影響由大到小為X5、X2、X1、X4。

表7 方差分析Tab.7 Analysis of variance

為準(zhǔn)確求得各因素最優(yōu)參數(shù)組合，需盡量縮小約束條件范圍，已知Y1和Y2目標(biāo)值分別為0.45 mm、58.668 N，由表6可知，試驗(yàn)方案2、13、15和19的Y1和Y2最接近實(shí)際值，且X5均為0.20 mm，根據(jù)方差結(jié)果可知X4對Y1和Y2的影響程度最小，所以設(shè)定約束條件范圍為

(5)

利用Design-Expert 8.0.6軟件進(jìn)行優(yōu)化求解，得到X1、X2、X4、X5分別為9.978 25×109N/m3、1.033 75×1010N/m3、1.20×106Pa、0.2 mm時(shí)，Y1、Y2分別為0.450 mm、58.668 N，對優(yōu)化求解得到的參數(shù)組合進(jìn)行仿真驗(yàn)證，得到Y(jié)1、Y2仿真值分別為0.451 mm、58.993 N，與真實(shí)值的相對誤差分別為0.222%、0.554%。

3 肥料顆粒破碎率驗(yàn)證試驗(yàn)

3.1 試驗(yàn)條件

為進(jìn)一步驗(yàn)證優(yōu)化得到的肥料顆粒Bonding模型參數(shù)組合的可靠性，選取目前播種機(jī)上常用的外槽輪排肥器進(jìn)行排肥性能試驗(yàn)，試驗(yàn)肥料為前文包膜肥料顆粒并過篩(孔徑3 mm)后備用，如圖8a所示。每次試驗(yàn)稱量1 500 g粒徑大于3 mm的完整肥料顆粒，經(jīng)外槽輪排肥器完全排出后再次過篩(孔徑3 mm)，統(tǒng)計(jì)粒徑小于3 mm的肥料質(zhì)量，在不同排肥輪轉(zhuǎn)速下進(jìn)行試驗(yàn)，每組試驗(yàn)重復(fù)3次取平均值。采用三維軟件SolidWorks 2018對排肥器進(jìn)行建模并保存為.IGS文件，將模型導(dǎo)入離散元軟件EDEM 2018中，如圖8b所示。在與真實(shí)試驗(yàn)相同的排肥輪轉(zhuǎn)速、槽輪有效工作長度條件下進(jìn)行虛擬仿真試驗(yàn)，各條件下肥料顆粒破碎率η計(jì)算公式為

圖8 肥料顆粒破碎試驗(yàn)Fig.8 Fertilizer particles crushing test1.肥料 2.接肥盒 3.步進(jìn)電機(jī) 4.仿真肥料顆粒 5.外殼 6.排肥輪 7.排肥舌

(6)

式中m1——粒徑小于3 mm的肥料質(zhì)量，g

m0——每次試驗(yàn)肥料總質(zhì)量，取1 500 g

3.2 試驗(yàn)結(jié)果

在排肥輪轉(zhuǎn)速分別為10、15、20、25、30 r/min時(shí)得到肥料顆粒真實(shí)破碎率和仿真破碎率如圖9所示。由圖9可以看出，肥料顆粒仿真破碎率與真實(shí)情況相比誤差為9.42%～11.40%，說明本文標(biāo)定得到的Bonding模型參數(shù)組合準(zhǔn)確性高。

圖9 不同轉(zhuǎn)速下肥料顆粒破碎率Fig.9 Fertilizer particles crushing rate at different rotational speeds

肥料顆粒在排肥輪和排肥舌的擠壓、摩擦作用下，顆粒間粘結(jié)鍵斷裂，整個(gè)顆粒從中間破碎為兩大塊，顆粒破碎仿真過程如圖10所示。由圖10d～10f可以看出，在排肥輪和肥料剛接觸時(shí)，肥料除受擠壓力外，還伴隨滾動，說明實(shí)際情況下肥料顆粒包膜材料還會發(fā)生研磨現(xiàn)象，此過程不會破碎肥料顆粒，但在一定程度上破壞包膜結(jié)構(gòu)，進(jìn)而影響肥料養(yǎng)分釋放過程。

圖10 排肥過程中肥料顆粒破碎仿真過程Fig.10 Fertilizer granules crushing simulation

4 結(jié)論

(1)通過單軸壓縮試驗(yàn)得到包膜肥料顆粒的極限破碎位移和極限破碎載荷分別為0.45 mm、58.668 N。

(2)Placket-Burman試驗(yàn)結(jié)果表明，粘結(jié)半徑、單位面積切向剛度和單位面積法向剛度對極限破碎位移影響顯著，粘結(jié)半徑、單位面積法向剛度和臨界切向應(yīng)力對極限破碎載荷影響顯著。

(3)Box-Behnken試驗(yàn)結(jié)果表明，粘結(jié)半徑、單位面積法向剛度、單位面積切向剛度、臨界法向應(yīng)力、臨界切向應(yīng)力和粘結(jié)半徑均對極限破碎位移和極限破碎載荷影響顯著。優(yōu)化求解得到當(dāng)單位法向剛度、單位面積切向剛度、臨界切向應(yīng)力和粘結(jié)半徑分別為9.978 25×109N/m3、1.033 75×1010N/m3、1.20×106Pa、0.2 mm時(shí)，極限破碎位移和極限破碎載荷分別為0.450 mm、58.668 N。最優(yōu)參數(shù)條件下得到極限破碎位移和極限破碎載荷與真實(shí)值的相對誤差分別為0.222%、0.554%。

(4)外槽輪排肥器排肥結(jié)果顯示，在排肥輪轉(zhuǎn)速分別為10、15、20、25、30 r/min時(shí)肥料顆粒真實(shí)破碎率和仿真破碎率相對誤差不大于11.40%，說明本文標(biāo)定得到的Bonding模型參數(shù)組合準(zhǔn)確性高。