黑水虻生物轉化餐廚垃圾有機肥離散元模型參數標定

宋世圣，孫松林，方芹，彭才望，周婷，朱海英

（湖南農業(yè)大學機電工程學院，長沙 410128）

隨著人口數量的增加和人們飲食習慣的變化，餐廚垃圾的量日益增多。餐廚垃圾中含有大量的有機質，如處理不當易造成環(huán)境破壞，從而影響人們的健康，餐廚垃圾如何有效且潔凈地處理就顯得尤為重要[1-4]。

黑水虻近年來被廣泛關注，其通過蟲體的取食與消化，可將餐廚垃圾轉化為含有高質量蛋白質的生物質[5-7]。現有關于黑水虻生物轉化餐廚垃圾有機肥的研究主要集中在黑水虻幼蟲的生長發(fā)育以及生物轉化餐廚垃圾的轉化率等方面[8-9]，在黑水虻生物轉化餐廚垃圾有機肥的離散元仿真模型方面研究極少，影響了后續(xù)收集、輸送、篩分機械的研究，降低了生產效益。

離散元素法(discrete element method，DEM)是新型數值計算方法，其認為整體是由大量離散個體組成，同時各個個體之間存在靠攏與離散2種形式，可用于分析與求解復雜系統(tǒng)的動力學問題[10-11]。物料離散元參數標定是研究物料與機械相互作用的基礎，韓樹杰等[12]利用漏斗注入法對新疆果園深施散體廄肥進行了離散元參數標定，為機械化深施廄肥離散元仿真提供了理論依據；王韋韋等[13]對玉米秸稈粉料致密成型的離散元模型進行了參數標定，通過模孔壓縮對比試驗驗證了其標定參數的準確性，可為秸稈粉料致密成型的仿真過程提供基礎參考；石辰風等[14]基于“Hertz-Mindlin with JKR Cohesion”模型對4種中藥浸膏粉末進行了參數標定，同時對比了不同硬件配置對仿真計算結果的影響，可以對中藥制藥中的藥粉混合和輸送過程進行精確離散元模擬；彭才望等[15]利用圓筒提升法對黑水虻處理的豬糞有機肥進行了參數標定，為該物料收集機械的設計提供一定理論參考；邢潔潔等[16]對海南熱區(qū)磚紅壤顆粒進行了離散元模型參數標定，通過破土阻力試驗，驗證了其標定參數的可靠性，在離散元軟件中，可以較好地模擬觸土部件與海南熱區(qū)磚紅壤的接觸作用。

目前，關于黑水虻生物轉化餐廚垃圾的研究大多圍繞在轉化效率[17]及生態(tài)效益[18]等方面，關于該物料的轉運、輸送研究較少。而通過對餐廚垃圾有機肥離散元模型參數的標定，可以對機械化收集、轉運物料的過程進行離散元模擬仿真，極大地減少收集機械的設計研究時間，降低勞動力成本，提高生產效益。因此，本研究利用仿真與物理相結合的方法，構造餐廚垃圾有機肥模型，采用圖像處理方法精確測量堆積角，依次通過Plackett-Burman試驗、最陡爬坡試驗、Box-Behnken試驗確定物料的離散元參數數值選擇，可為黑水虻生物轉化餐廚垃圾有機肥在收集、轉運過程中的離散元仿真提供一定理論依據。

1 材料與方法

1.1 試驗材料

黑水虻生物轉化餐廚垃圾有機肥（以下簡稱餐廚垃圾有機肥）取自湖南農業(yè)大學黑水虻科研基地（28.18°N，113.06°E），由含水率70%左右的餐廚垃圾經黑水虻幼蟲處理10～12 d后形成，該餐廚垃圾主要成分包括肉類、青菜、辣椒、土豆、蓮藕、主食等。餐廚垃圾中因含有骨頭等硬物，會影響后續(xù)處理設備的正常運行[17-18]，因此在餐廚垃圾處理前，使用小型粉碎機進行粉碎，黑水虻生物轉化處理后的餐廚垃圾有機肥呈顆粒狀，經烘干法測得其含水率為37.5%。由于在餐廚垃圾中存在辣椒籽等黑水虻無法轉化的雜質，通過預測量發(fā)現，有機肥顆粒粒徑一般在2 mm以下，雜質顆粒粒徑一般在2.6 mm以上，因此使用10目標準篩（篩孔2 mm）對有機肥進行篩分，保留篩下物作為堆積角試驗的材料。將顆粒分作10組，每組200個顆粒，使用數顯游標卡尺（精度0.02 mm）對粒徑進行測量。粒徑分布如圖1所示，走勢同正態(tài)分布，大多數顆粒粒徑在1.1～1.3 mm之間，粒徑均值1.2 mm，標準差為0.25 mm。

圖1 粒徑分布Fig.1 Particle size distribution

1.2 接觸參數測定

1.2.1 靜摩擦系數的測定 靜摩擦系數是物料受到最大靜摩擦力與法向正壓力的比值，可以較好地反映物料與接觸表面的摩擦特性[19]。首先將傾角儀調0，為防止有機肥顆粒較小而產生的滑動摩擦及有機肥顆粒間的黏結性，因此，使用上下無底的方形模具插入餐廚垃圾有機肥堆體，使其形成長4 mm、寬4 mm、高5 mm的有機肥立方體[27]，放置在光滑的鋼板（圖2）上，鋼板下放置1個滑塊，通過推動滑塊調節(jié)鋼板傾斜角度，在有機肥立方體與鋼板產生相對滑動的瞬間，停止推動，立即記錄傾角儀讀數，試驗重復10次，通過式（4）得到餐廚垃圾有機肥與鋼板的靜摩擦系數。同時在測定餐廚垃圾有機肥之間的靜摩擦系數時，只需在鋼板上利用雙面膠粘上0.5 mm厚的有機肥顆粒排即可，粘結時盡可能讓顆粒排布緊密，以免黏結劑對摩擦系數產生影響[25,28]。

圖2 靜摩擦系數測定原理Fig.2 Principle of static friction coefficient measurement

式中，G為有機肥立方體重力（N）；N為斜面對有機肥立方體的支撐力（N）；f為斜面對有機肥立方體的最大靜摩擦力（N）；μ1為靜摩擦系數；θ1為有機肥立方體滑動瞬間的傾角（°）。

1.2.2 滾動摩擦系數的測定 餐廚垃圾有機肥相較于鋼板彈性模量較小，因此在鋼板滾動時，可忽略微觀變形和彈性滯后，只有少許滑動影響試驗結果，且傾斜角度越大，滑動摩擦影響越大[19]。因此在進行滾動摩擦系數測定試驗前，使用有機肥圓柱體進行預滾動，測得有機肥滾動的最小角度為5.7°。

利用直徑為10 mm的圓形環(huán)刀將餐廚垃圾有機肥制成半徑5 mm、高20 mm的圓柱體，由鋼板頂端滾動至鋼板底端，利用相機拍攝有機肥圓柱體在鋼板滾動的全過程，將視頻導入Premiere 2018軟件中，選取不同幀數記錄1、2、3、4點的位置（圖3），1、2距離為L1，3、4距離為L2，由于L1、L2距離較小，可以求得中點a、b的瞬時速度，通過式（5）得到餐廚垃圾有機肥與鋼板之間的滾動摩擦系數[20]，試驗重復10次。有機肥與有機肥之間的滾動摩擦系數測定同1.2.1。

圖3 滾動摩擦系數測定原理Fig.3 Measuring principle of rolling friction coefficient

式中，μ2為滾動摩擦系數；m為有機肥的質量（kg）；h為有機肥下落高度（m）；vb為b點瞬時速度（m·s-1）；va為a點瞬時速度（m·s-1）；L3為a點到b點位移（m）；g為重力加速度9.81（m·s-2）。

1.3 堆積角物理測定

有機肥堆積角試驗中漏斗下端管口直徑18 mm，漏斗下端管口與鋼板距離100 mm。試驗時，利用擋板將漏斗口堵住，將有機肥顆粒緩慢倒入漏斗，待穩(wěn)定后移開擋板，顆粒經漏斗口落在鋼板上形成顆粒堆，待該顆粒堆穩(wěn)定后，利用相機拍攝顆粒堆視圖，將照片導入Photoshop 2018中，在+x及+y方向分別標注角度取平均值即為該顆粒堆的堆積角，重復10次取平均值，測得堆積角為41.34°。

1.4 EDEM仿真模型構建及仿真參數設置

餐廚垃圾有機肥顆粒外觀近似球形，因此將有機肥顆粒仿真模型設置為球形，有機肥顆粒平均粒徑為1.2 mm，設置動態(tài)顆粒工廠，顆粒尺寸以正態(tài)分布生成，粒徑標準差0.25 mm，選擇Cap，Lower設置為0.5，Upper設置為1.5。在SolidWorks中根據物理試驗模型1∶1進行構建，以IGS文件格式導入EDEM 2020，模型材料為鋼，泊松比為0.30，密度為7 860 kg·m-3，剪切模量為7.9×1010Pa[15]，如圖4所示。

圖4 堆積角仿真模型Fig.4 Simulation model of accumulation angle

設置factory在漏斗上方以virtual狀態(tài)生成顆粒，物理有機肥顆粒百粒重在0.511～0.585 g之間，均值為0.534 g。在物理試驗時取100 g有機肥顆粒進行堆積角試驗，為與物理堆積角試驗條件一致，在構建仿真模型時，有機肥顆粒共18 000個，生成速率為9 000個·s-1，生成顆粒時，擋板堵住漏斗口，2.5 s時待顆粒在漏斗中穩(wěn)定后，移開擋板，使顆粒由漏斗口落在鋼板上。固定時間步長為Rayleigh時間步長的22%，Cell Size為3 Rmin，數據保存的時間間隔為0.01 s。由于餐廚垃圾有機肥較易粘結且具有一定濕度，故選取“Hertz-Mindlin with JKR”模型作為此次離散元仿真試驗的接觸模型[21]。

1.5 參數標定試驗設計

根據斜面試驗及已有研究進展[12-16,22-25]設置參數水平（表1），利用Design Expert 10.0軟件設計Plackett-Burman試驗進行優(yōu)化，以相對誤差為評價指標，當相對誤差最小時說明該組參數即為各組參數中的最優(yōu)選擇，從而得到參數的最優(yōu)區(qū)間。根據上述試驗確定參數范圍，以顯著性參數為試驗因素，以堆積角為響應值，設計Box-Behnken試驗。在Design Expert軟件的響應優(yōu)化模塊中，以餐廚垃圾有機肥物理堆積角41.34°為目標值對回歸模型進行尋優(yōu)求解，得到參數的最優(yōu)組合，在該參數組合下進行仿真試驗，重復試驗3次，與物理試驗結果進行對比。

表1 Plackett-Burman試驗參數Table 1 Factors of Plackett-Burman Design

1.6 數據處理

采用Design Expert 10.0軟件進行試驗設計及數據分析，采用Origin 2017進行繪圖，采用Matlab 2018b進行圖像處理。

2 結果與分析

2.1 Plackett-Burman篩選試驗

Plackett-Burman試驗篩選參數結果如表2所示，可以看出，仿真試驗堆積角數值跨度較大，在30.2°～62.85°之間，與物理試驗堆積角的相對誤差較大，因此需進行方差分析篩選顯著性參數。

表2 Plackett-Burman試驗設計及結果Table 2 Design and results of Plackett-Burman Design

方差分析結果如表3所示，可以看出，模型P＜0.05，表明模型顯著。對響應值有正向效應的因素有A、E、G、H，說明在其他因素不變時，隨著該因素的增大，響應值也隨之增大。對響應值有負向效應的因素有B、C、D、F、J、K，說明在其他因素不變時，隨著上述因素的增大，響應值會逐漸減小。D的P值小于0.01，說明D因素對于堆積角的影響極顯著，G、J的P值均小于0.05，說明其對堆積角的影響顯著，其余因素對堆積角的變化均無顯著影響，比較各參數的P值大小，各因素顯著性排序為D、J、G、B、A、H、K、C、E、F。故選取對堆積角影響較大的D、G、J進行后續(xù)試驗。

表3 Plackett-Burman試驗結果顯著性分析Table 3 Significance analysis of Plackett-Burman test results

2.2 最陡爬坡試驗結果分析

由Plackett-Burman試驗篩選得到3個顯著性參數，根據標準化效應在一定步長范圍內增大或減小，可以使仿真的結果快速逼近物理試驗值[25]。在仿真試驗中，改變顯著參數D、G、J的取值，而其他非顯著參數取中間值，即A為0.2、B為7 MPa、C為2 000 kg·m-3、E為0.8、F為0.15、H為0.65、K為0.1 J·m-2。通過表4可以看出，隨著參數G的增大和參數D、J的減小，仿真堆積角與物理堆積角相對誤差的走勢呈現先減小后增大的趨勢，其中最小誤差為4.74%，最大誤差則達到48.16%。而第2組參數組合下得到的堆積角相對誤差最小，因此分別選擇D為0.58、G為0.42、J為0.09作為后續(xù)Box-Behnken試驗的中心參數組合水平，第1與第3組的參數組合分別為低水平和高水平。

表4 最陡爬坡試驗設計與結果Table 4 Design and results of steepest ascent test

2.3 Box-Behnken試驗結果分析

以3個顯著性參數D、G、J為試驗因素，堆積角作為響應值，設計Box-Behnken試驗，其中，5組為中心水平重復試驗，以探究D、G、J對堆積角的影響效果，試驗設計方案及結果如表5所示。對試驗結果進行多元回歸分析，得到堆積角回歸方程（式3）。

表5 Box-Behnken試驗設計及結果Table 5 Design and results of Box-Behnken Design

方差分析結果如表6所示，可以看出，該模型P＜0.01，可以利用該回歸模型對有機肥顆粒堆的堆積角進行準確的預測。失擬項P＞0.05，不顯著，說明此回歸方程中不存在失擬項。各因素顯著性由高到低排序為J＞D＞G，其中D、J對堆積角影響極顯著，G對堆積角影響顯著，各因素交互作用對響應值均無顯著性影響效果，D2、G2對堆積角有極顯著影響，其余因素二次方項對響應值并無顯著影響效果。變異系數為2.92%，說明此試驗具有一定的可靠性。試驗精確度Ap為14.391，說明試驗具有較高的精確度。決定系數R2為0.967 8，校正決定系數為0.926 4，均接近于1，表明該回歸方程擬合程度較好，可以用于堆積角預測分析。

2.4 參數最優(yōu)組合確定及驗證

以餐廚垃圾有機肥物理堆積角41.34°為目標值對回歸模型進行尋優(yōu)求解，得到最優(yōu)組合為：D為0.46、G為0.51、J為0.09；其他參數數值選取依次為：A為0.2、B為7 MPa、C為2 000 kg·m-3、E為0.8、F為0.15、H為0.65、K為0.1 J·m-2，進行仿真試驗，利用EDEM后處理模塊的測量工具測得3次仿真試驗的堆積角分別為41.85°、43.4°、43.28°。3次仿真試驗的平均值為42.84°，與物理試驗堆積角（41.34°）的相對誤差為3.63%，利用Matlab對圖像進行處理，提取輪廓曲線（圖5），可以看出，通過標定后得到的參數進行的仿真試驗堆積邊界與物理試驗較為接近。t檢驗P值為0.094（＞0.05），證明仿真試驗結果與物理試驗結果無顯著性差異。

圖5 物理試驗及仿真模擬對比Fig.5 Comparison of physical test and simulation test

3 討論

目前圍繞如何提高黑水虻轉化餐廚垃圾效率的研究較多，關于機械收集、轉運餐廚垃圾有機肥的研究較少。張杰等[29]對黑水虻資源化處理有機廢棄物進行了闡述，認為黑水虻處理城市餐廚垃圾領域具有較大的應用前景，但在黑水虻對餐廚垃圾生物轉化完成后如何回收地槽或料盤中的有機肥顆粒是目前的難題。本課題組前期設計了斗式取料機[30]及雙向螺旋裝置[31]，一定程度上可以收集、轉運餐廚垃圾有機肥，但在設計過程中由于缺乏離散元模型參數，產生了設計進度慢、試驗勞動強度大等一系列問題。因此本文以堆積角為評價指標，通過對照比較仿真試驗與物理試驗的堆積角數值，得到餐廚垃圾有機肥顆粒的離散元模型參數，仿真試驗與物理試驗的相對誤差僅為3.63%，可以為有機肥收集機械的離散元仿真提供一定數據參考。

離散元法可以精確地模擬散體物料的運動形態(tài)，但農業(yè)物料的含水量往往對其物理特性具有較大影響，而本文僅對含水率為37.5%的餐廚垃圾有機肥進行離散元模型參數標定，后續(xù)將進一步探究不同含水率餐廚垃圾有機肥的離散元參數變化規(guī)律，完善本領域的離散元仿真數據庫。另外需注意，在收集、轉運有機肥顆粒時，物料內部存在一定數量的五齡期黑水虻幼蟲，尚為活體，生物特性較為復雜，蔡靜[32]對活體蠅蛆采用噴灑化學試劑及降溫法使其靜止，從而建立了兩種姿態(tài)的蠅蛆離散元仿真模型，但由于蠅蛆與黑水虻之間的體態(tài)特征和生物特性存在較大差異，因此黑水虻幼蟲與有機肥顆粒之間的接觸參數有待進一步研究。