基于堆積試驗(yàn)的黑水虻離散元仿真參數(shù)標(biāo)定與分析

彭才望，周 婷，孫松林，謝燁林，魏 源

(湖南農(nóng)業(yè)大學(xué) 機(jī)電工程學(xué)院，湖南 長(zhǎng)沙 410128)

近年來(lái)，隨著規(guī)?；B(yǎng)豬場(chǎng)的迅速發(fā)展，豬場(chǎng)的糞污排放量逐漸增多，對(duì)豬場(chǎng)周?chē)拇髿?、地下水等造成了?yán)重的環(huán)境污染。同時(shí)，給人體的生活健康帶來(lái)一定的危害。黑水虻幼蟲(chóng)(黑水虻)在畜禽糞便“肥料化利用”過(guò)程中作用明顯，能實(shí)現(xiàn)減污、資源循環(huán)利用。近年來(lái)針對(duì)黑水虻的研究主要集中在黑水虻對(duì)畜禽糞便和餐廚垃圾的轉(zhuǎn)化率、黑水虻幼蟲(chóng)營(yíng)養(yǎng)成分、黑水虻自身生長(zhǎng)參數(shù)變化等方面。而黑水虻處理畜禽糞便后續(xù)分離環(huán)節(jié)中篩分機(jī)械應(yīng)用程度低，分離效率低、成本高，降低了規(guī)?；a(chǎn)率。

黑水虻分離是指在黑水虻幼蟲(chóng)養(yǎng)殖成熟后期，通過(guò)一定工程技術(shù)方法將畜禽糞便等養(yǎng)殖基料中的黑水虻排泄物、分泌物以及養(yǎng)殖基料殘?jiān)蛛x，獲得黑水虻產(chǎn)品，是黑水虻處理畜禽糞便實(shí)現(xiàn)資源化利用的關(guān)鍵環(huán)節(jié)之一。黑水虻的物性參數(shù)與接觸參數(shù)是實(shí)現(xiàn)機(jī)械化分離黑水虻的關(guān)鍵之一。堆積角作為散體物料的固有屬性，常被用作散體物料顆粒參數(shù)標(biāo)定，也可用來(lái)研究散體顆粒的復(fù)雜力學(xué)特性。離散元法可以模擬顆粒間微觀(guān)與宏觀(guān)變形，在農(nóng)業(yè)機(jī)械、農(nóng)業(yè)物料領(lǐng)域得到廣泛的應(yīng)用。王國(guó)強(qiáng)等闡述了通過(guò)離散元堆積角試驗(yàn)與實(shí)際試驗(yàn)值對(duì)比驗(yàn)證，獲取離散元參數(shù)。林嘉聰?shù)韧ㄟ^(guò)堆積角試驗(yàn)與離散元仿真堆積角虛擬試驗(yàn)相結(jié)合的方法，量化了蚯蚓糞的滾動(dòng)摩擦因數(shù)。徐莉等對(duì)蠅蛆與固體畜禽糞便混合物分離運(yùn)動(dòng)進(jìn)行了離散元仿真研究，量化了蠅蛆的物性參數(shù)與接觸參數(shù)并揭示了蠅蛆與振動(dòng)篩模型間的運(yùn)動(dòng)學(xué)規(guī)律，為同屬于資源性昆蟲(chóng)的黑水虻仿真研究提供了借鑒依據(jù)。但是，黑水虻與蠅蛆在生物特性與形態(tài)方面存在一定差異，目前應(yīng)用離散元方法對(duì)黑水虻仿真模型的研究鮮有報(bào)道。

本文以黑水虻幼蟲(chóng)為研究對(duì)象，分析黑水虻的生物特性，利用EDEM仿真軟件并選取系統(tǒng)中的“Hertz-Mindlin with JKR”接觸模型，進(jìn)行堆積角仿真試驗(yàn)。首先，應(yīng)用Plackett-Burman設(shè)計(jì)試驗(yàn)篩選對(duì)黑水虻堆積角影響顯著的參數(shù)。然后，進(jìn)一步以堆積角作為響應(yīng)值，通過(guò)三因素二次回歸正交旋轉(zhuǎn)組合響應(yīng)面試驗(yàn)確定出最優(yōu)參數(shù)組合，并將標(biāo)定后的最優(yōu)參數(shù)進(jìn)行仿真，對(duì)比驗(yàn)證堆積角以確定出黑水虻的組合參數(shù)。以期為機(jī)械分離黑水虻的離散元仿真模擬提供基本參數(shù)，并為黑水虻在蟲(chóng)沙分離環(huán)節(jié)中的瞬間運(yùn)動(dòng)狀態(tài)研究奠定基礎(chǔ)。

1 材料與方法

1.1 基本參數(shù)

試驗(yàn)采用湖南大湘農(nóng)環(huán)境科技股份有限公司提供的黑水虻樣品，該黑水虻是8日齡黑水虻幼蟲(chóng)投放到含水率70%～80%的新鮮豬糞中養(yǎng)殖8 d后，且在化蛹前的狀態(tài)，黑水虻取食豬糞中營(yíng)養(yǎng)物質(zhì)獲得生長(zhǎng)并進(jìn)行生物轉(zhuǎn)化。試驗(yàn)時(shí)，因生物應(yīng)激性使黑水虻的外形容易發(fā)生變化，形成伸直(外形自然伸展)與蜷曲(外形受應(yīng)激而萎縮)兩種狀態(tài)的體型，導(dǎo)致與仿真試驗(yàn)結(jié)果存在一定的差異。因此，為得到堆積試驗(yàn)自然狀態(tài)的黑水虻，堆積試驗(yàn)結(jié)束瞬間，噴灑化學(xué)試劑使黑水虻保持瞬間不變狀態(tài)，克服非物性接觸參數(shù)帶來(lái)的不利影響。然后，每次隨機(jī)選取預(yù)蛹前的50條黑水虻幼蟲(chóng)，重復(fù)試驗(yàn)，測(cè)得伸直狀態(tài)黑水虻體長(zhǎng)23.3 mm，標(biāo)準(zhǔn)差0.70 mm，蜷曲狀態(tài)黑水虻體長(zhǎng)19.4 mm，標(biāo)準(zhǔn)差0.58 mm；蟲(chóng)寬4.9 mm，標(biāo)準(zhǔn)差0.34 mm；蟲(chóng)體重為0.177 g。伸展?fàn)顟B(tài)和蜷曲狀態(tài)的黑水虻數(shù)量分布比例為3∶2，黑水虻如圖1所示。

1.2 黑水虻堆積物理試驗(yàn)

黑水虻堆積角通過(guò)箱體抽板法試驗(yàn)測(cè)得。測(cè)量裝置(不銹鋼材料)由箱體(長(zhǎng)100 mm、寬100 mm、高150 mm)、擋板(長(zhǎng)100 mm、高150 mm)和底板(長(zhǎng)200 mm、寬200 mm)組成。試驗(yàn)時(shí)將黑水虻置于箱體中，緩慢向上抽提擋板，設(shè)定擋板向上抽提速度為50 mm·s，失去擋板壁面的支撐后，黑水虻滑出箱體，因黑水虻自身具備生物應(yīng)激性，堆積落地后，黑水虻自身外表產(chǎn)生非自身接觸參數(shù)導(dǎo)致的微觀(guān)變化，影響堆積角真實(shí)測(cè)量。因此，堆積穩(wěn)定瞬間利用高速攝影機(jī)在垂直鋼板底平面的方向連續(xù)拍照，堆積的坡面與底板的水平面形成的夾角即為堆積角。選取黑水虻自然堆積形態(tài)較穩(wěn)定、應(yīng)激性干擾小的圖像導(dǎo)入CAD軟件中進(jìn)行堆積角標(biāo)注。為降低生物應(yīng)激性帶來(lái)的誤差并較真實(shí)反映黑水虻堆積過(guò)程中的自然狀態(tài)，換用新的黑水虻樣品重復(fù)試驗(yàn)5次，取其堆積角平均值34.66°，標(biāo)準(zhǔn)差為0.65°，如圖2所示。

圖1 黑水虻幼蟲(chóng)Fig.1 Black soldier fly larva

2 仿真模型

2.1 黑水虻仿真模型

依據(jù)文獻(xiàn)[21，23-25]，將黑水虻顆粒間的接觸設(shè)定為軟球模型，顆粒接觸形式如圖3所示。顆粒間的接觸方式簡(jiǎn)化為彈簧阻尼系統(tǒng)，接觸模型簡(jiǎn)化如圖4所示。試驗(yàn)過(guò)程中顆粒間所發(fā)生的碰撞力簡(jiǎn)化為軟球模型間的接觸合力，通過(guò)牛頓運(yùn)動(dòng)定律對(duì)、顆粒間的運(yùn)動(dòng)與相互作用進(jìn)行計(jì)算。接觸合力可分解為顆粒與顆粒碰撞時(shí)的接觸法向力與接觸切向力。軟球模型之間的合力計(jì)算公式為：

=∑(+)。

(1)

式中—顆粒和間的接觸合力(N)；—顆粒和間的接觸法向力(N)；—顆粒和間的接觸切向力(N)。

黑水虻顆粒模型在堆積碰撞運(yùn)動(dòng)試驗(yàn)中，可將發(fā)生的碰撞設(shè)定為顆粒和之間的碰撞，顆粒從圖3虛線(xiàn)位置運(yùn)動(dòng)到實(shí)線(xiàn)位置，運(yùn)動(dòng)形式包括法向運(yùn)動(dòng)、切向運(yùn)動(dòng)以及顆粒間的滾動(dòng)運(yùn)動(dòng)，不考慮顆粒和之間的自身變形，、分別為顆粒間碰撞接觸時(shí)重疊處的位移量和切向位移。根據(jù)黑水虻顆粒軟球模型及其受力與運(yùn)動(dòng)分析情況以及黑水虻處理畜禽糞便過(guò)程中養(yǎng)殖基料中存在一定的濕度，黑水虻表面存在一定的黏性，依據(jù)文獻(xiàn)[21]，在EDEM軟件中應(yīng)用Hertz-Mindlin with JKR模型進(jìn)行黑水虻堆積角仿真試驗(yàn)。

1，底板；2，黑水虻幼蟲(chóng)；3，擋板；4，箱體；5，底座板。1， Bottom plate; 2， Black soldier fly larva; 3， Baffle; 4， Box; 5， Base plate.圖2 箱體抽板法試驗(yàn)Fig.2 Test of box side plate lifting method

圖3 顆粒接觸形式Fig.3 Form of particle contact

黑水虻因自身生物特性導(dǎo)致體長(zhǎng)變化，呈現(xiàn)伸展與蜷曲兩種狀態(tài)。為建立合乎實(shí)際的黑水虻模型，在不影響分析結(jié)果的前提下，按照?qǐng)D4測(cè)得黑水虻的基本參數(shù)，采用16個(gè)球顆粒兩兩重疊組合的方式，通過(guò)調(diào)整圓球坐標(biāo)位置，使各個(gè)圓球球心均在同一水平面上，建立黑水虻的仿真模型如圖5所示。

2.2 仿真參數(shù)設(shè)置

EDEM仿真試驗(yàn)需要設(shè)置黑水虻顆粒與邊界相關(guān)參數(shù)。其中邊界為不銹鋼材料，參數(shù)唯一確定。黑水虻顆粒因自身生物特性，實(shí)際物理試驗(yàn)測(cè)量比較困難，離散元參數(shù)值難以確定唯一性，且當(dāng)前缺乏黑水虻離散元仿真與試驗(yàn)研究相關(guān)的文獻(xiàn)。因此，結(jié)合文獻(xiàn)[23-25]中蠅蛆與蛆料振動(dòng)篩分離散元仿真參數(shù)，確定了本研究中各仿真參數(shù)的變化范圍如表1所示。

圖4 顆粒接觸模型簡(jiǎn)化Fig.4 Model of particle contact

圖5 黑水虻仿真模型Fig.5 Simulation model of black soldier fly

2.3 試驗(yàn)?zāi)Ｐ?/h3>

在EDEM軟件中導(dǎo)入與箱體抽板法物理試驗(yàn)尺寸相同的箱體模型(長(zhǎng)100 mm、寬100 mm、高150 mm)，箱體底部添加與箱體垂直的水平面底板(長(zhǎng)200 mm、寬200 mm)，如圖6所示。箱體頂部設(shè)置一個(gè)虛擬平面作為顆粒工廠(chǎng)，顆粒生成的速度為0.20 kg·s，生成時(shí)間為2.50 s，黑水虻顆?？傎|(zhì)量0.5 kg，將伸直與蜷曲兩種狀態(tài)的黑水虻顆粒模型按數(shù)量分布設(shè)置3∶2，以動(dòng)態(tài)隨機(jī)方式生成并落入箱體內(nèi)部。在Simulator中，瑞利時(shí)間步長(zhǎng)設(shè)置為22%，仿真中網(wǎng)格尺寸取2倍最小球形單元尺寸。為了保證黑水虻顆粒快速穩(wěn)定，設(shè)定顆粒初始下落速度為0.5 m·s，待黑水虻顆粒充滿(mǎn)箱體且顆粒群穩(wěn)定后，擋板垂直向上以50 mm·s的速度提升，失去擋板壁面的支撐，顆粒群滑出箱體并從底部逐漸流出，待黑水虻顆粒群穩(wěn)定后，其坡面與底平面所形成的夾角即為黑水虻堆積角，后處理中在垂直顆粒堆坡面方向截取圖像并導(dǎo)入CAD軟件中進(jìn)行堆積角標(biāo)注。

表1 仿真參數(shù)設(shè)置

1，底板；2，黑水虻；3，擋板；4，箱體。1， Bottom plate; 2， Black soldier fly larva; 3， Baffle; 4， Box.圖6 箱體抽板法堆積角仿真Fig.6 Simulation of box side plate lifting method

3 仿真試驗(yàn)

3.1 試驗(yàn)設(shè)計(jì)

應(yīng)用Design-Expert軟件進(jìn)行Plackett-Burman(PB)試驗(yàn)設(shè)計(jì)，篩選對(duì)黑水虻堆積角影響顯著的參數(shù)。根據(jù)EDEM仿真時(shí)用到涉及黑水虻顆粒相關(guān)的參數(shù)，選取表1中影響黑水虻堆積角的10個(gè)參數(shù)變量作為試驗(yàn)因素，并預(yù)留1個(gè)虛擬變量()用作誤差分析。根據(jù)表1給定的參數(shù)建立黑水虻顆粒模型開(kāi)展堆積角仿真試驗(yàn)，其中～分別選擇上限、下限兩個(gè)數(shù)值作為高、低2個(gè)水平(1、-1)，試驗(yàn)中設(shè)3 個(gè)中心點(diǎn)，根據(jù)未知參數(shù)數(shù)量進(jìn)行15次試驗(yàn)。應(yīng)用Design-Expert軟件進(jìn)行PB試驗(yàn)設(shè)計(jì)及結(jié)果如表2所示，其中～為因素水平值，為空白列。

對(duì)表2所示的Plackett-Burman試驗(yàn)結(jié)果進(jìn)行分析，得到10個(gè)因素對(duì)黑水虻堆積角貢獻(xiàn)率如表3所示。由表3可知，10個(gè)因素中對(duì)黑水虻堆積角影響較為顯著(<0.05)的3個(gè)因素依次為：(黑水虻間靜摩擦系數(shù))、(黑水虻間滾動(dòng)摩擦系數(shù))、(黑水虻剪切模量)，且這3個(gè)因素對(duì)堆積角均為正效應(yīng)，即堆積角隨著該3個(gè)因素增大而增大。其余因素影響不顯著(>0.05)。因此，在后續(xù)最陡爬坡試驗(yàn)和二次回歸正交旋轉(zhuǎn)組合響應(yīng)面試驗(yàn)中選取這3個(gè)參數(shù)作為試驗(yàn)因素進(jìn)一步分析，其余非顯著性參數(shù)選取中間水平值。

表2 試驗(yàn)設(shè)計(jì)及結(jié)果

表3 試驗(yàn)結(jié)果方差分析

3.2 最陡爬坡試驗(yàn)

根據(jù)PB試驗(yàn)篩選出的3個(gè)顯著性參數(shù)為、、，為了使這3個(gè)顯著性參數(shù)達(dá)到最優(yōu)區(qū)域，基于PB試驗(yàn)結(jié)果中顯著因素的正效應(yīng)按確定的固定步長(zhǎng)逐漸增加，其余非顯著性參數(shù)選取中間水平值，進(jìn)行最陡爬坡試驗(yàn)，并計(jì)算黑水虻仿真堆積角與實(shí)際堆積角的相對(duì)誤差，試驗(yàn)方案及結(jié)果如表4所示。

由表4可知，隨著3個(gè)試驗(yàn)參數(shù)取值的增加，仿真得到的黑水虻堆積角與物理試驗(yàn)堆積角的相對(duì)誤差先減小后增大，其中在2號(hào)試驗(yàn)對(duì)應(yīng)的參數(shù)組合下，得到二者相對(duì)誤差最小，說(shuō)明最優(yōu)的參數(shù)水平范圍在2號(hào)試驗(yàn)附近。因此，以1、2、3號(hào)試驗(yàn)中各參數(shù)值作為二次回歸正交旋轉(zhuǎn)組合試驗(yàn)的-1、0、1水平的因素進(jìn)行后續(xù)的響應(yīng)面試驗(yàn)。

表4 最陡爬坡試驗(yàn)設(shè)計(jì)及結(jié)果

3.3 響應(yīng)面試驗(yàn)

根據(jù)PB試驗(yàn)篩選出對(duì)堆積角的顯著影響參數(shù)后，以、、為試驗(yàn)參數(shù)，結(jié)合最陡爬坡試驗(yàn)結(jié)果確定各因素試驗(yàn)水平，黑水虻堆積角為試驗(yàn)指標(biāo)，非顯著性試驗(yàn)參數(shù)延用表1各因素試驗(yàn)水平并選取表中間水平值。利用Design-Expert軟件進(jìn)行三因素二次回歸正交旋轉(zhuǎn)組合試驗(yàn)并進(jìn)行尋優(yōu)，仿真試驗(yàn)因素編碼如表5所示，試驗(yàn)設(shè)計(jì)及仿真結(jié)果如表6所示。

表5 仿真因素試驗(yàn)編碼

表6 試驗(yàn)方案與結(jié)果

表7 回歸方程方差分析

(2)

在保證模型顯著、失擬項(xiàng)不顯著的情況下，剔除模型不顯著的回歸項(xiàng)，對(duì)模型進(jìn)行優(yōu)化，回歸方程為

(3)

分析優(yōu)化后的模型，根據(jù)模型的值(<0.000 1)，與模型的失擬項(xiàng)(=0.605 2)可知，模型優(yōu)化可靠。

由表7所示模型方差分析結(jié)果可知，黑水虻剪切模量和黑水虻間滾動(dòng)摩擦系數(shù)的交互項(xiàng)()對(duì)黑水虻堆積角影響顯著(<0.05)。當(dāng)黑水虻間靜摩擦系數(shù)為0.40時(shí)，利用Design-Expert 8.0.6仿真軟件繪制黑水虻剪切模量和黑水虻間滾動(dòng)摩擦系數(shù)交互作用的響應(yīng)面圖，直觀(guān)地看到此兩個(gè)參數(shù)之間的交互效應(yīng)，如圖7所示?？梢钥闯觯逊e角隨黑水虻間滾動(dòng)摩擦系數(shù)變化的曲面變化明顯，而隨黑水虻剪切模量變化的曲面坡度較為平緩。隨著黑水虻剪切模量的增大，堆積角先緩慢減小再緩慢增大。當(dāng)黑水虻剪切模量為定值時(shí)，堆積角隨黑水虻間滾動(dòng)摩擦系數(shù)的增加而顯著增大。

圖7 黑水虻剪切模量與黑水虻間滾動(dòng)摩擦系數(shù)交互作用Fig.7 The interaction between the shear modulus of black soldier fly and the rolling friction coefficient of black soldier fly

4 優(yōu)化驗(yàn)證試驗(yàn)

為驗(yàn)證黑水虻顆粒組合參數(shù)仿真得到的堆積角與實(shí)際堆積角之間的相對(duì)誤差，采用標(biāo)定后的黑水虻顆粒參數(shù)，結(jié)合表1設(shè)置其他因素水平，進(jìn)行堆積角仿真試驗(yàn)與實(shí)際堆積角試驗(yàn)結(jié)果對(duì)比，驗(yàn)證優(yōu)化參數(shù)下的仿真試驗(yàn)結(jié)果與實(shí)際情況的接近程度?；贒esign-Expert8.0.6軟件的優(yōu)化模塊，以黑水虻實(shí)際堆積角為響應(yīng)目標(biāo)值對(duì)回歸方程進(jìn)行求解并尋優(yōu)。目標(biāo)與約束條件方程組為：

(4)

得到多組最佳參數(shù)組合，選擇其中一組參數(shù)得到的仿真結(jié)果與物理試驗(yàn)堆積角相對(duì)誤差最小。同時(shí)，考慮黑水虻模型是由多個(gè)球形重疊組成，模型外表面與活體黑水虻表面相比，表面存在更多的曲面與凹凸，增加了模型顆粒間的摩擦阻力，流動(dòng)性降低。因此，選取最優(yōu)參數(shù)時(shí)，黑水虻間靜摩擦系數(shù)與黑水虻間滾動(dòng)摩擦系數(shù)不宜取最大，以提高仿真結(jié)果與實(shí)際結(jié)果的吻合程度。即，黑水虻剪切模量8.67 MPa、黑水虻間靜摩擦系數(shù)0.43、黑水虻間滾動(dòng)摩擦系數(shù)0.32，其余非顯著性因素選擇中間水平值(黑水虻泊松比0.3、黑水虻密度350 kg·m、黑水虻間碰撞恢復(fù)系數(shù)0.2、黑水虻-不銹鋼碰撞恢復(fù)系數(shù)0.3、黑水虻-不銹鋼靜摩擦系數(shù)0.6、黑水虻-不銹鋼滾動(dòng)摩擦系數(shù)0.45、JKR表面能0.35 J·m)，其余設(shè)置不變，應(yīng)用EDEM軟件進(jìn)行3次重復(fù)試驗(yàn)，得到堆積角均值為35.84°，與箱體抽板法實(shí)際測(cè)得堆積角34.66°的相對(duì)誤差為3.40%，無(wú)顯著差異，如圖8所示。試驗(yàn)結(jié)果表明3個(gè)顯著性參數(shù)最優(yōu)值準(zhǔn)確可靠。

圖8 仿真試驗(yàn)與物理試驗(yàn)對(duì)比Fig.8 Comparison of simulation and physical tests

5 結(jié)論

(1)結(jié)合黑水虻生物特性，以物理堆積試驗(yàn)瞬間獲得的黑水虻作為研究對(duì)象，采用“Hertz-Mindlin with JKR”作為EDEM離散元仿真模擬的接觸模型，應(yīng)用抽板法進(jìn)行堆積角試驗(yàn)，標(biāo)定了黑水虻在接觸過(guò)程中的物理性能參數(shù)。

(2)應(yīng)用Plackett-Burman試驗(yàn)設(shè)計(jì)篩選出對(duì)黑水虻堆積角有顯著影響的參數(shù)，包括：黑水虻剪切模量、黑水虻間的靜摩擦系數(shù)、黑水虻間滾動(dòng)摩擦系數(shù)。結(jié)合抽板法堆積物理試驗(yàn)結(jié)果，通過(guò)最陡爬坡試驗(yàn)確定了3個(gè)顯著性參數(shù)的最優(yōu)值區(qū)域。

(3)基于三因素二次回歸正交旋轉(zhuǎn)組合響應(yīng)面試驗(yàn)，建立了顯著性參數(shù)與堆積角間的二次回歸模型，以堆積角為響應(yīng)值，對(duì)優(yōu)化后的回歸模型進(jìn)行尋優(yōu)，確定了EDEM仿真試驗(yàn)中最佳接觸參數(shù)組合為黑水虻剪切模量為8.67 MPa、黑水虻間靜摩擦系數(shù)為0.43、黑水虻間滾動(dòng)摩擦系數(shù)為0.32，以最優(yōu)參數(shù)組合得到的仿真試驗(yàn)堆積角為35.84°，與實(shí)際測(cè)量堆積角相對(duì)誤差為3.40%，無(wú)顯著性差異，表明標(biāo)定黑水虻離散元仿真模型參數(shù)準(zhǔn)確可靠。