基于EDEM-Fluent的氣動式魚塘投飼機性能優(yōu)化

趙思琪 丁為民，2 趙三琴，2 顧家冰，2 張建凱

(1.南京農業(yè)大學工學院， 南京 210031； 2.江蘇省智能化農業(yè)裝備重點實驗室， 南京 210031)

0 引言

隨著大規(guī)模高密度水產養(yǎng)殖模式的發(fā)展，飼料已成為養(yǎng)殖過程主要的成本因素，約占總養(yǎng)殖成本的60%[1-2]，投飼機作為投飼工序的核心裝備，其作業(yè)性能直接影響飼料利用效率和養(yǎng)殖效益。

針對當前大面積高密度養(yǎng)殖水域的作業(yè)需求，氣動式投飼技術得以迅速發(fā)展[3-5]，借助氣力輸送技術實現飼料遠距離輸送，經拋料機構高速轉動將飼料360°拋撒至水面，完成投飼過程。在實際應用中，由于給料箱與拋料結構相距較遠，現有氣動式投飼機以偏置式進料、離心式拋撒為主，但普遍存在拋料不均勻問題[6-8]，研究表明飼料分布不均勻不僅會降低魚類攝食效率，而且會增加魚群體間競爭強度，影響飼料利用率和魚群同等生長[9-10]。為提高氣動式投飼機拋料均勻性，國外專家學者通過建立投飼機拋料機構數學模型，分析風速、顆粒大小、拋料角度對飼料分布的影響，以優(yōu)化拋料性能，但難以從根本上解決不均勻性問題[11-12]；目前國內相關研究主要集中在整機結構設計和控制系統(tǒng)研究[13-15]，對影響拋料性能的主要因素及優(yōu)化改進方案研究較少。近年來，隨著計算機技術的發(fā)展，離散元素法(Discrete element method，DEM)及其數值模擬仿真軟件EDEM在農業(yè)顆粒物料運動分析方面得到廣泛應用[16-17]，為研究離散飼料顆粒與拋料機構相互作用等提供了良好的平臺和手段。同時，Fluent作為一種有效的三維流場分析軟件，包含豐富而先進的湍流模型，可以精確模擬分析無粘流、層流、湍流等流場，實現對流場內不同空間壓力、濃度、速度等物理量的研究，為研究氣力輸送過程提供了良好的平臺基礎[18]。因此，本文以傳統(tǒng)氣動式投飼機為研究對象，應用EDEM-Fluent耦合方法[19-20]建立氣動式投飼機拋料機構的運動學仿真模型，研究飼料顆粒在拋料作業(yè)中運動過程，分析影響投飼機拋料不均勻的主要原因，在此基礎上提出相應優(yōu)化改進方案，探討各因素及其相互作用對拋料均勻性的影響規(guī)律，以得到提高拋料性能的最優(yōu)設計參數組合，并對改進效果進行仿真和場地試驗驗證。

1 氣動式投飼機結構與工作原理

如圖1所示，氣動式投飼機結構主要包括料箱、拋料盤、進料管、進料室、交流電動機、支架等部件。

圖1 氣動式投飼機結構示意圖Fig.1 Structure diagrams of pneumatic feeder1.拋料盤 2.進料管 3.進料室 4.交流電動機 5.支撐底座 6.支架 7.浮球 8.拋料葉片 9.動力軸 10.導流板

關鍵部件為拋料盤和進料室，拋料盤由上下圓板和葉片組成，圓板半徑R1為270 mm，葉片長、寬為195、72 mm，截面為矩形，共12條，徑向繞動力軸均布在拋料盤內，作業(yè)時拋料盤以角速度ω高速旋轉，對飼料起到負壓吸送和離心拋撒作用；進料室由腔體和導流板組成，腔體主要起到對氣固兩相流的混合和過渡輸送作用，進料室內徑R2為50 mm，進料室高度h為110 mm，腔體底部進料口為圓面，進料管內徑R3為35 mm；導流板主要起到引流和均分飼料的作用，將飼料流均勻往左右兩個半腔引流，有利于飼料在腔內均勻混合，由寬度30 mm、長度70 mm的“L”型角鐵制作，安置在料腔底部進料口中間位置。

拋料工作原理：拋料盤高速旋轉使投飼機和輸料管內部產生負壓吸力，將前端料箱供給的飼料在負壓氣力作用下沿輸料管進入進料室內，經導流板均分引流后沿室壁進入腔體內部混合，進而在氣力作用下向上輸送至拋料盤，在水平拋料盤高速旋轉作用下沿徑向葉片導軌做加速運動，至料盤邊緣以一定速度拋出后做類平拋運動，最后落入水面完成投飼過程。

2 顆粒飼料運動特性分析

2.1 飼料顆粒在離心盤上的運動

假定飼料顆粒為剛性質點，質量為m，因拋料盤內氣流對飼料質點的作用力遠小于離心力，在此忽略不計[21]。則飼料顆粒模型在拋料盤內受力分析如圖2所示。

圖2 飼料顆粒模型受力分析Fig.2 Force analysis of feed particle model

顆粒隨拋料盤做繞軸轉動和沿葉片滑移的合成運動，由圖可知，顆粒所受離心力F1為

F1=mω2l

(1)

式中l(wèi)——飼料顆粒與拋料盤中心線垂直距離，mm

顆粒所受科氏力Fk為

(2)

葉片側面對顆粒的摩擦力Ff為

(3)

式中μ——飼料顆粒與拋料葉片之間摩擦因數

顆粒重力引起的摩擦力F2為

F2=μmg

(4)

則沿OA方向，由牛頓第二定律得

(5)

式中r0——拋料葉片末端與拋料盤中心線垂直距離，mm

解式(5)得飼料顆粒在離心盤內的運動方程為

(6)

(7)

(8)

(9)

(10)

2.2 飼料顆粒脫離拋料盤后的運動

飼料顆粒以一定速度脫離拋料盤后，經類平拋運動最終落入水面，在空氣域中受重力、空氣阻力等作用。因飼料顆粒較小，在無邊界的空氣流體中可以將飼料顆粒的圓柱形狀近似于球面進行分析[22]，則在空氣中運動時，空氣阻力Fa為

(11)

式中Cd——空氣阻力系數

ρ1——空氣密度，kg/m3

rp——飼料顆粒當量半徑，mm

v——飼料顆粒相對空氣速度，m/s

以離心盤水平面為xy面，以垂直盤方向與重力同向為z軸，則飼料顆粒在空中3個方向運動方程為

(12)

(13)

(14)

(15)

式中vx——沿x方向速度分量，m/s

vy——沿y方向速度分量，m/s

vz——沿z方向速度分量，m/s

通過對飼料在拋料盤上的運動分析(式(6))可知，飼料顆粒的拋撒速度取決于拋料盤轉速、拋料盤半徑和飼料顆粒與拋料盤初始接觸位置，而飼料顆粒的拋撒速度決定了投飼機的最大投飼距離、投飼面積等基本作業(yè)參數，為保證在進行投飼均勻性能改進的同時投飼機的基本作業(yè)參數不受影響，以滿足實際工況及投飼機標準[24]中有關投飼機基本作業(yè)參數的要求，在EDEM-Fluent仿真時應設置邊界條件拋料盤轉速和半徑保持不變(符合實際作業(yè)工況)。飼料顆粒隨拋料盤做繞軸轉動和沿葉片滑移的合成運動，運動分析知其在脫離拋料盤后應做斜拋曲線運動，對飼料脫離拋料盤的運動軌跡分析可與EDEM-Fluent仿真結果中飼料顆粒的軌跡進行對比驗證。

3 拋料過程仿真分析

3.1 三維仿真模型構建與物料特性參數選擇

為便于模擬和分析，去除圖1中與飼料運動過程接觸無關的部件，應用Pro/E三維軟件進行實體建模并以.igs格式導入EDEM中，如圖3所示。研究表明在氣動式投飼機拋料過程中，飼料顆粒大小對飼料分布規(guī)律性沒有顯著性影響[10-11]，因此，本文以國內池塘養(yǎng)殖常用的3 mm柱形顆粒飼料為研究對象進行氣動式投飼機拋料性能分析研究，參照文獻[16-17,25-27]中相關研究及試驗方法，對飼料顆粒物理特性參數進行給定和測定，確定仿真參數如表1所示。EDEM-Fluent選擇Eulerian-Eulerian雙歐拉模型進行耦合，采用Fluent軟件中k-ε湍流模型，由于飼料顆粒之間無黏連、顆粒表面光滑，在此采用EDEM軟件中Hertz-Mindlin無滑動接觸模型進行仿真。根據前面對飼料顆粒的運動學分析，為確保投飼機的相關基本作業(yè)參數不受影響，根據實測值設定拋料盤轉速為2 800 r/min，入口風速根據實測值設定為24 m/s，設置顆粒工廠以2 000顆/s(即最大投飼能力)的速率生成飼料顆粒，進行飼料顆粒運動仿真分析，參數設置貼近實際魚塘作業(yè)，以提高模擬仿真的準確性。

圖3 拋料機構仿真圖Fig.3 Simulation diagrams of parabolic mechanism

表1 物料特性參數Tab.1 Material properties of particles

3.2 仿真試驗設計與試驗指標

根據離心式拋撒機械性能評價標準[22,28]，使用飼料周向分布變異系數作為拋撒均勻性能評價指標。應用EDEM后處理Selection模塊，以拋料盤中心為圓心，半徑為22 m(場地試驗得氣動投飼機最大拋撒距離)的圓面上均勻分成12個圓心角為30°的扇形計算網格(圖4)，將落入每個計算網格內的飼料總質量進行統(tǒng)計，計算飼料周向分布變異系數

(16)

其中

(17)

(18)

式中Cv——周向分布變異系數，%

S——標準差，g

n——計算網格數量

M0——收集域內網格收集飼料質量平均值，g

Mi——第i個計算網格收集飼料顆?？傎|量，g

圖4 EDEM中計算網格布置圖Fig.4 Layout diagram of computational grid in EDEM

3.3 拋料性能虛擬仿真結果與分析

運用EDEM-Fluent耦合方法對投飼機拋料性能進行虛擬仿真試驗，以確定造成拋料不均勻的主要原因。為保證仿真的連續(xù)性和可靠性，設置固定時間步長為2×10-5s，為Rayleigh步長的18%，總時間為1 s。如圖5所示，以飼料顆粒流為研究對象，設置其以流線型(stream)顯示，圖中彩色流線表示粒子運動軌跡，不同顏色變化代表運動過程速度變化，表示顆粒在拋料過程的運動狀態(tài)，可以發(fā)現其在拋料盤內和脫離拋料盤后的運動軌跡與飼料顆粒理論運動特性分析相一致，在脫離拋料盤后呈現斜拋曲線運動軌跡。同時，從圖5和圖6可知，在進料腔內氣流對飼料顆粒作用明顯，飼料顆粒在氣流作用下沿進料腔左側壁(進料口對側)被帶入拋料盤，呈明顯偏置分布狀態(tài)，這與氣流流場分布一致，進入拋料盤后氣流場作用減弱，主要受拋料盤高速旋轉離心作用，隨拋料盤做繞軸轉動和沿葉片滑移的合運動，這與實際作業(yè)過程一致。圖7為飼料顆粒拋撒分布圖，飼料顆粒在整個投飼區(qū)周向呈現倒“V”形的分布，周向分布變異系數Cv為62.16%，在0°～60°和270°～360°扇形區(qū)域分布最為集中，占總顆粒量的64.72%，這與實際情況相符合，同時與飼料顆粒和氣流場在進料腔內的分布一致，均呈現明顯偏向一側。分析可知保證飼料顆粒能夠均勻輸送至拋料盤內，是解決拋撒不均勻的關鍵，而進料腔和導流板結構是影響腔體內氣流場分布和飼料顆粒分布的關鍵結構。

圖5 拋料過程性能仿真Fig.5 Performance simulation diagram of throwing process

圖6 拋料過程進料室氣流場仿真Fig.6 Performance simulation diagrams of air flow field in feeding chamber during throwing process

圖7 飼料顆粒周向分布Fig.7 Pellet circumferential distribution

4 結構優(yōu)化設計與仿真

4.1 結構優(yōu)化設計及分析

圖8 氣動式投飼機拋料結構改進示意圖Fig.8 Structure diagram for improvement of pneumatic feeder

為從根本上解決拋料不均勻問題，根據3.3節(jié)分析，對進料腔體和導流板結構進行改進設計(圖8)。將進料腔設計為由外腔和內腔兩部分組成，實現內腔進料口與外腔進料口錯位布置，以優(yōu)化腔體內氣流場分布和提高氣固有效混合時間和空間，飼料顆粒在氣流帶動下首先進入外腔體氣固混合，然后經內腔體底部開口從兩側進入內腔(圖8紅色箭頭)并在氣流作用下沿內腔壁送入拋料盤，為不影響結構安裝及機器作業(yè)，設計外腔體半徑Ro為75 mm，外腔體高度ho為80 mm，內腔體尺寸與第1節(jié)所述相同。同時改進導流板為圓弧形，左右對置各一個，使氣固混合流可以沿著內壁切線方向進入外腔體，并只能從內腔體輸送至拋料盤，有利于改善腔內氣固分布狀態(tài)，為不影響入口進料，設計半徑與改進前進料室內徑R2相同，導流板圓心角α為90°。

圖9和圖10為改進后拋料過程仿真結果，可以發(fā)現混合室內飼料顆粒和氣流場在整個內腔體基本均勻分布，不存在偏置現象，飼料顆粒在氣流作用下可以較均勻的從內腔輸送至拋料盤內，完成拋撒作業(yè)。圖11為改進后飼料顆粒周向分布情況，各個區(qū)域內飼料顆粒質量分數在8%上下浮動，各區(qū)域間的分布差異性明顯降低，周向分布變異系數Cv為16.49%，同改進前相比降低48.23個百分點，可以更好地滿足機具作業(yè)需求，其中出現的浮動主要是由料腔內形成的兩個小的湍流團引起。

圖9 改進后拋料過程性能仿真Fig.9 Performance simulation diagram of throwing process after improvement

圖10 改進后拋料過程進料室氣流場仿真Fig.10 Performance simulation diagrams of air flow field in feeding chamber during throwing process after improvement

圖11 改進后飼料顆粒分布Fig.11 Pellet distribution after improvement

4.2 正交回歸旋轉中心仿真優(yōu)化

在對拋料結構優(yōu)化改進的基礎上，研究關鍵設計參數：外腔體半徑、外腔體高度和導流板圓心角對拋料性能的影響，并通過分析優(yōu)化出投飼機拋料機構最優(yōu)設計參數組合。根據實際安裝和作業(yè)需求，設置試驗因素取值范圍分別為：外腔體半徑70～80 mm；外腔體高度70～100 mm；導流板圓心角80°～100°。以外腔體半徑Ro、外腔體高度ho和導流板圓心角α為影響因素，飼料顆粒周向分布變異系數Cv為評價指標，采用Design-Expert 8.06軟件設計三元二次回歸正交旋轉組合試驗，試驗因素編碼見表2，仿真方案和仿真試驗結果如表3所示。

表2 試驗因素編碼Tab.2 Coding of factors

4.2.1試驗結果方差分析

對試驗結果進行多元回歸擬合分析，得到周向分布變異系數方差分析如表4所示，二次回歸模型高度顯著(P<0.000 1)，失擬項(P=0.076 7>0.05)不顯著，回歸方程不失擬，說明可以用回歸方程代替試驗真實點對試驗結果進行分析，方程較好地表征了周向分布變異系數與外腔體半徑、外腔體高度和導流板圓心角的相關關系。經逐步回歸法剔除不顯著項后，得到Cv對外腔體半徑、外腔體高度和導流板圓心角的二次多元回歸方程

表3 試驗方案和結果Tab.3 Test plan and results

圖12 周向分布變異系數的雙因素交互響應曲面Fig.12 Response surfaces of two factors for variation coefficient of circumferential distribution

(19)

4.2.2響應曲面分析

由表4知，Roho和hoα對周向分布變異系數影響顯著，Roα對其影響不顯著，在此對Roho和hoα交互作用對變異系數的影響效應進行分析。圖12a為α取90°時，外腔體半徑Ro與外腔體高度ho交互作用響應曲面圖，由圖知外腔體高度ho在70～85 mm范圍時隨著ho的增大變異系數逐漸降低，在85～100 mm范圍內時對變異系數的影響規(guī)律相反，但影響程度減弱；外腔體半徑Ro在70～75 mm范圍時隨著Ro的增加周向分布變異系數逐漸降低，在75～80 mm范圍內影響規(guī)律相反。當外腔體半徑Ro在74～78 mm，外腔體高度ho在80～90 mm之間時，周向變異系數Cv均較小。

表4 周向分布變異系數方差分析Tab.4 Variance analysis of variation coefficient for circumferential distribution

注：** 表示影響高度顯著(P<0.01)；*表示影響顯著(P<0.05)。

圖12b為Ro取75 mm時，外腔體高度ho與導流板圓心角α交互作用響應曲面圖，分析圖知外腔體高度ho在70～85 mm范圍時隨著ho的增大變異系數逐漸降低，在85～100 mm范圍內時對變異系數的影響規(guī)律相反，但影響程度減弱；導流板圓心角α在80°～90°范圍時隨著α的增加周向變異系數逐漸降低，在90°～100°范圍內影響規(guī)律相反。當外腔體高度ho在80～90 mm，α在85°～95°范圍時周向變異系數Cv較小。

4.2.3最優(yōu)參數組合確定

為確定影響飼料周向分布變異系數的最優(yōu)組合設計參數，利用Design-Expert的中心組合響應曲面設計功能進行參數優(yōu)化，以表4中各因素范圍為約束條件，以周向分布變異系數回歸方程(式(19))為目標函數，求解其最小值。參數優(yōu)化結果如圖13所示，圖中黃色區(qū)域表示優(yōu)化范圍，優(yōu)化后的最優(yōu)參數組合：外腔體半徑Ro為77.01 mm(取77.00 mm)，外腔體高度ho為85.41 mm(取85.40 mm)，導流板圓心角α為93.17°(取93.20°)，此時周向分布變異系數Cv最小，其理論結果為14.13%。

圖13 中心組合參數優(yōu)化結果Fig.13 Parameter optimization result with central composite design

5 試驗驗證

根據優(yōu)化分析所得最優(yōu)設計參數組合重新制作加工拋料機構(圖14)，并同改進前拋料機構在南京農業(yè)大學工學院農機試驗場地進行定點拋料性能試驗，如圖15所示。飼料選擇與仿真條件一致的榮威1038柱形顆粒漁用飼料，含水率為12.3%，按照現有對離心式撒肥機械性能試驗標準[24，28]，將以投飼機拋料盤中心為圓心，半徑為22 m的圓收集域等分成12份圓心角為30°的扇形區(qū)域，將360 mm×260 mm×60 mm 的收集盒在以投飼機拋料盤中心為圓心，半徑分別為3、6、9、12、15、18、21 m的圓上均勻布置，其中半徑為3 m和6 m的圓上每隔10°布置一個收集盒，半徑為12、15、18、21 m的圓上每隔5°布置一個收集盒(圖15)。試驗結束收集每個區(qū)域內所有收集盒飼料總質量，計算得飼料周向分布變異系數Cv，改進前后各進行3次重復試驗。得到改進前變異系數分別為：64.28%、66.57%、65.90%，均值為65.58%，平均相對誤差為5.51%；優(yōu)化改進后變異系數分別為：14.79%、15.32%、15.13%，均值為15.08%，平均相對誤差為6.72%，試驗驗證結果與仿真預測結果基本一致，說明了軟件優(yōu)化參數的準確性和可行性。同改進前相比，周向分布變異系數降低50.50個百分點，表明優(yōu)化改進方案的有效性。

圖14 投飼機進料腔體改進前后實物圖Fig.14 Pictures of feeder before and after improvement

圖15 試驗現場圖及收集盒布置示意圖(0°～90°)Fig.15 Diagram of test site and distribution map of collection box (0°～90°)

6 結論

(1)利用EDEM-Fluent耦合方法對氣動式投飼機作業(yè)過程進行模擬仿真，分析知進料腔和導流板結構參數是影響拋料性能的關鍵因素。

(2)在拋料結構優(yōu)化改進的基礎上，設計三元二次回歸正交旋轉組合試驗和響應面分析，結果表明各影響因素對飼料周向分布變異系數的影響程度從大到小為：導流板弧度角α、外腔體半徑Ro、外腔體高度ho。優(yōu)化分析得，當Ro為77.00 mm、ho為85.40 mm、α為93.20°時，周向分布變異系數Cv最小，其理論結果為14.13%。

(3)場地試驗驗證表明，最優(yōu)參數組合下周向分布變異系數為15.08%，改進前周向分布變異系數為65.58%，與軟件模擬仿真結果基本一致，表明回歸模型能夠較好地模擬拋料顆粒周向分布情況，借助EDEM-Fluent耦合法分析優(yōu)化氣動式投飼機拋料性能是可行的。