基于氣固兩相DEM模型的相對旋轉(zhuǎn)式喂料頭的研究

孫曉霞， 孟文俊， 牛雪梅， 袁媛

(1.太原科技大學(xué) 機(jī)械工程學(xué)院，山西,太原 030024; 2.山西能源學(xué)院，山西,太原 030024)

隨著垂直螺旋輸送機(jī)的廣泛應(yīng)用，作為其重要組成部分的取料裝置也有了各種形式，其中最重要的一種形式就是相對旋轉(zhuǎn)式喂料頭(簡稱喂料頭)，其喂料量和喂料效率關(guān)系到輸送機(jī)的輸送效率[1-2]. 在國際上，垂直螺旋輸送機(jī)及其喂料頭的技術(shù)相對先進(jìn)，提出了如多片導(dǎo)料擋板式反向旋轉(zhuǎn)喂料頭、雙口螺線式喂料頭、螺旋擋板式喂料頭、傾斜翼板式喂料頭等多種結(jié)構(gòu)的喂料頭[3].對喂料頭的研制是通過試驗(yàn)，測得不同尺寸的喂料頭模型在不同的喂料頭轉(zhuǎn)速的情況下的充填率、功率及生產(chǎn)率，由試驗(yàn)曲線描述取料段性能指標(biāo)隨喂料頭結(jié)構(gòu)參數(shù)的變化規(guī)律，確定較佳的設(shè)計變量. 相對于國外先進(jìn)的技術(shù)優(yōu)勢，喂料頭在我國目前尚未有統(tǒng)一的標(biāo)準(zhǔn)，每個廠家都有各自的加工標(biāo)準(zhǔn)和輸送量確定的標(biāo)準(zhǔn).

喂料頭的喂料效率和喂料量直接影響到垂直螺旋輸送機(jī)的輸送效率，而對喂料頭的喂料量產(chǎn)生影響的主要因素就是喂料頭的喂料曲面. 因此，喂料曲面的投影、展開平面以及特征曲線的形狀及方程也是各有特點(diǎn). 顆粒進(jìn)入喂料頭后，由于料壁內(nèi)螺旋的高速旋轉(zhuǎn)，會產(chǎn)生氣流影響顆粒的運(yùn)動，所以對于喂料頭中的顆粒，不僅要考慮顆粒與顆粒、喂料頭之間的相互作用，還要考慮顆粒與氣體之間的相互作用. 顆粒在喂料頭中的運(yùn)動比較復(fù)雜，受力也比較復(fù)雜，難以在理論上深入探討顆粒在喂料頭中的運(yùn)動過程[4].

近年來，隨著計算機(jī)技術(shù)的快速發(fā)展，數(shù)值模擬已逐漸成為與試驗(yàn)手段相輔相成的研究工具. 宋斌等[5-6]采用DEM模型及EDEM軟件仿真對喂料頭進(jìn)行了研究，分析了垂直螺旋輸送機(jī)內(nèi)顆粒在不同條件下的運(yùn)動情況，得出影響顆粒運(yùn)動的因素. 但是喂料頭中需要考慮氣體對顆粒的作用，所以需考慮一種數(shù)值模型，不僅需要考慮顆粒與顆粒、喂料頭之間的相互作用，還要考慮顆粒與氣體之間的相互作用[7]. 文中采用E/L法對垂直螺旋輸送機(jī)進(jìn)行模擬，由于顆粒密度較高，顆粒之間的碰撞不能忽略，而離散單元法(DEM)是顆粒碰撞模型中最應(yīng)用廣泛的一種. 對于模擬喂料頭這種復(fù)雜的顆粒流系統(tǒng)，建立氣固兩相DEM模型，采用離散單元法和計算流體力學(xué)(CFD)相結(jié)合的方式，可以跟蹤每個顆粒的運(yùn)動，獲得大量的微觀信息. 在DEM模型及其仿真的方法基礎(chǔ)上，結(jié)合理論分析與氣固兩相DEM模型，及DEM+CFD耦合仿真分析顆粒的喂料過程，考察顆粒在不同曲面喂料頭中的運(yùn)動狀態(tài)及輸送量情況.

1 氣固兩相DEM模型及其仿真方法

在喂料頭中，因?yàn)轭w粒所受外力及摩擦力的復(fù)雜性和不確定性，難以用傳統(tǒng)的牛頓第二定律的方法來確定顆粒的速度，而且，顆粒還受到空氣的影響，因此，下面采用氣固兩相DEM模型及其仿真分析方法來觀察顆粒的輸送量情況.

1.1 氣固兩相DEM模型

1.1.1 顆粒固相控制方程

在氣固兩相DEM模型中，顆粒的運(yùn)動通過求解牛頓第二定律獲得，

ma=Fg+Fp-p+Fp-g+Fge,

(1)

式中:a為固相顆粒對于喂料葉片的加速度;Fg為重力;方向垂直向下;Fp-p為顆粒與顆粒、喂料葉片之間的作用力，方向不確定;Fp-g為氣體與顆粒之間的作用力，方向不確定;Fge為離心力.

1.1.2 氣相控制方程

在氣固兩相DEM模型中，氣體的運(yùn)動規(guī)律由N-S方程描述，設(shè)εg為氣體空隙率，vg為氣體表觀速度，則氣相控制方程為

(2)

(3)

其中，氣相控制方程中動量方程的流固耦合項(xiàng)為

Fp=β(vg-vp),

(4)

式中β為氣固相間曳力系數(shù)，其模型采用Gidaspow曳力模型. 曳力是氣固相間最主要的作用力.可以看出，氣體對顆粒的曳力Fd與顆粒對氣體的力Fp是一對相互作用力.

1.2 接觸曳力模型

文中采用DEM模型及支持此模型的仿真軟件EDEM軟件進(jìn)行仿真，其提供的顆粒之間的接觸模型有Hertz-Mindlin(no slip)模型等[8]. 文中仿真的顆粒都是無黏結(jié)性小顆粒，根據(jù)實(shí)際情況選擇Hertz-Mindlin 無滑動接觸模型，該模型能高效計算兩個顆粒接觸時的接觸作用力，由法向力和切向力以及滑動摩擦力和滾動摩擦力組成. 法向力分量由Hertzian的接觸理論計算，切向力由Middlin-Deresiewicz計算，它們都有阻尼分量，其阻尼系數(shù)和恢復(fù)系數(shù)有關(guān). 庫倫摩擦定律決定滑動摩擦力，接觸獨(dú)立定向恒轉(zhuǎn)矩模型決定滾動摩擦力.

在顆粒的運(yùn)動過程中，顆粒受到來自氣體的曳力，而氣體又因?yàn)轭w粒的存在對其流動產(chǎn)生影響，因此兩者之間存在緊密的耦合作用，曳力影響了顆粒的輸送，因此曳力模型的準(zhǔn)確性是非常關(guān)鍵的.

1.3 氣固兩相耦合

氣固兩相DEM模型模擬程序由空氣流場Fluent計算和顆粒動力學(xué)EDEM計算兩部分構(gòu)成，模擬時，顆粒場與空氣流場的求解交叉進(jìn)行. EDEM+Fluent耦合求解過程如圖1所示. 在一個時間步長內(nèi)，F(xiàn)luent通過UDF和EDEM中的API接口相結(jié)合，先在Fluent中進(jìn)行流場計算，當(dāng)?shù)諗繒r，觸動EDEM開始當(dāng)前時步的顆粒計算，并把Fluent中的流場數(shù)據(jù)傳遞給EDEM，根據(jù)Gidaspow曳力模型計算氣體與顆粒間的相互作用力，再根據(jù)Hertz-Mindlin(no slip)接觸模型用以計算顆粒與顆粒、喂料葉片、螺旋葉片之間的相互作用力，然后，根據(jù)牛頓動力學(xué)方程給顆粒定位. 當(dāng)以上步驟在EDEM中完成后，再把氣體與顆粒間的作用力反饋回Fluent中，開始新一代的迭代.

圖1 EDEM+Fluent軟件耦合求解過程Fig.1 Coupling and solving process of EDEM with Fluent

2 喂料葉片曲面

2.1 曲面曲線方程

當(dāng)顆粒進(jìn)入喂料頭處的窗口時，隨著螺旋徑向尺寸的逐漸減小，迫使顆粒進(jìn)入螺旋輸送管底部中心附近或?yàn)⒙湓诼菪~片上，螺旋軸的高速旋轉(zhuǎn)使顆粒向上輸送，螺旋擋板沿高度的傾斜有效地防止了顆粒離心作用引起的濺出. 圖2是文中研究的喂料頭結(jié)構(gòu)圖.

圖2 喂料頭結(jié)構(gòu)Fig.2 Structure of feeder head

喂料頭的輸送量由進(jìn)入的顆粒和流出的顆粒之差組成，顆粒進(jìn)入的越多，流出的越少，輸送量越多，輸送量的多少是由螺旋的轉(zhuǎn)速、喂料頭的轉(zhuǎn)速和喂料頭的結(jié)構(gòu)決定的，而喂料頭的結(jié)構(gòu)主要考慮喂料導(dǎo)擋葉的曲面尺寸[8].

利用反向工程，得到某種高效的喂料頭的表面數(shù)據(jù)，并利用這些數(shù)據(jù)得到了實(shí)物原型的三維模型以及喂料曲面的展開圖，文中的喂料頭的數(shù)據(jù)都是源于宋斌等[5]、翟曉晨[6]基于反向工程的數(shù)據(jù)獲得的.

喂料頭單個導(dǎo)擋葉曲面俯視簡圖如圖3所示.

圖3 喂料頭單個導(dǎo)擋葉曲面俯視圖Fig.3 Sketch for top-view of single feeder surface of guide blade

要擬合的曲線由俯視狀態(tài)下的點(diǎn)云數(shù)據(jù)擬合獲得，通過測量、分析，在曲線上與中心O的連線達(dá)到最大值的點(diǎn)處將邊緣曲線分成q1、q2兩部分，曲線q1的r值變化較大，曲線q2的r值變化較小，分別對這兩部分曲線進(jìn)行擬合.

以夾角θ為自變量，r為相應(yīng)的因變量，利用Matlab軟件對能夠反映導(dǎo)擋葉喂料曲面的特征曲線進(jìn)行了方程擬合，得到的曲線為一個阿基米德螺旋線，設(shè)該曲線為曲線3，曲線方程為r=a+bθ，其中r為極徑，θ為極角(弧度).

經(jīng)擬合后q1的擬合曲線為

r=0.358+0.106θ， 0≤θ<π/3,

(5)

q2的擬合曲線為

r=0.478+0.025θ， 0≤θ≤π/2.

(6)

以上是在喂料口半徑為0.358 m時的擬合曲線，其喂料面積和該曲線是成正比的，所以喂料口半徑也是和該曲線是成正比的，設(shè)置喂料口半徑為Rf，則

q1的擬合曲線為

r=Rf+0.3Rfθ， 0≤θ<π/3,

(7)

q2的擬合曲線為

r=1.34Rf+0.07Rfθ， 0≤θ≤π/2.

(8)

當(dāng)顆粒在喂料葉片中運(yùn)動時，其有效速度是朝向喂料頭中心的速度vn1，它的朝向中心的受力和運(yùn)動速度由喂料葉片的阿基米德螺旋線決定，而且喂料頭的轉(zhuǎn)速不變，顆粒朝向中心的有效速度就不會改變. 根據(jù)關(guān)系式

t=θ/ωf,

(9)

式中:t為喂料頭旋轉(zhuǎn)的時間;θ為喂料頭轉(zhuǎn)過的角度(弧度);ωf為喂料頭角速度.

顆粒流向喂料頭中心的速度為

vnl=(r-a)/t.

(10)

將式(9)(10)聯(lián)立，并將其簡化得

vnl=bωf,

(11)

式中:a，b為阿基米德螺旋線的系數(shù)，a為θ=0時的極徑;b為每旋轉(zhuǎn)1°時極徑的增加量. 從式(11)可以看出，流向中心的有效速度vn1只與螺旋線的參數(shù)b和喂料頭的轉(zhuǎn)速ωf有關(guān).

下面改變阿基米德螺旋線b的值，觀察不同導(dǎo)擋葉曲線的喂料狀態(tài). 以漸開線在曲線3的基礎(chǔ)上逐漸增大，b分別增大10和20 mm，得到曲線2和曲線1，以漸開線在曲線3的基礎(chǔ)上逐漸減小，b分別減小10和20 mm，得到曲線4和曲線5. 即漸開線從大到小分別為曲線1～5. 圖4分別為各曲線的q1的曲線圖和q2的曲線圖.

圖4 喂料頭導(dǎo)擋葉片各曲線的曲線圖Fig.4 Curve graph of all kinds of guide blade of feeder head

因此，導(dǎo)擋葉俯視圖的各曲線q1、q2的方程分別為

曲線1:

(12)

曲線2:

(13)

曲線3:

(14)

曲線4:

(15)

曲線5:

(16)

喂料頭單個導(dǎo)擋葉左視如圖5所示.

圖5 單個導(dǎo)擋葉的左視圖Fig.5 Left optic of single feeder surface of guide blade

曲線s和曲線m分別為喂料頭導(dǎo)擋葉的上下曲線，根據(jù)反向工程得到曲線s和曲線m的數(shù)據(jù)，并進(jìn)行擬合.

可得s曲線擬合函數(shù)為

s=0.000 47θ+0.159.

(17)

式中:θ為喂料頭轉(zhuǎn)過的角度;s為曲線s離底端的距離.

可得m曲線擬合函數(shù)為

m=-0.000 007 6θ2+0.002 29θ+0.007 35.

(18)

設(shè)導(dǎo)擋葉的曲線m和曲線s不變，分別與前面q1、q2的曲線1～5的5個曲線組成不同的5個曲面，下面分析這5個曲面下的導(dǎo)擋葉的喂料頭的顆粒情況.

2.2 模型建立

由于喂料頭里面的顆粒還受到空氣的影響，采用氣固兩相DEM模型及EDEM+Fluent耦合仿真方法分析喂料頭內(nèi)的顆粒狀態(tài). 建立模型，取螺旋管徑為0.315 m，喂料口直徑為0.34 m，喂料窗口個數(shù)3，將上述在Solidwork中建立的模型導(dǎo)入到EDEM軟件中，然后再結(jié)合Fluent軟件進(jìn)行EDEM+Fluent耦合仿真分析，得到在考慮氣體時喂料頭內(nèi)的顆粒狀態(tài).

不同時刻的喂料情況如圖6所示.

圖6 喂料頭的喂料過程Fig.6 Feeding states of feeder head

從圖6中可以看出，當(dāng)運(yùn)行到第2.0 s時，喂料頭全部埋入顆粒堆中，顆粒開始沿著螺旋向上輸送，此后開始穩(wěn)定向上輸送顆粒. 顆粒在喂料頭中的運(yùn)動軌跡如圖7所示.

由圖7可以看出，顆粒在喂料頭中的運(yùn)動非常復(fù)雜，顆粒的運(yùn)動速度方向有朝向中心，進(jìn)入喂料頭的；也有離開中心，從喂料頭中濺出的. 所以喂料頭的輸送量由兩部分組成，即進(jìn)入喂料頭的顆粒量，濺出喂料頭的顆粒量.

圖7 顆粒運(yùn)動軌跡圖Fig.7 Movement trace of bulk materials

顆粒與喂料頭之間的靜摩擦因數(shù)是指顆粒與喂料頭之間摩擦力與法向正壓力的比值，它不僅與喂料頭的材料有關(guān)，而且與顆粒的表面形狀和粗糙度有關(guān)，該摩擦因數(shù)對喂料頭中顆粒的運(yùn)動狀態(tài)起很大的作用，所以文中采用不同的靜摩擦因數(shù)μ進(jìn)行分析，各種仿真工況如表1所示.

表1 喂料頭數(shù)值模擬參數(shù)表Tab.1 Various simulation working condition of feeding head

如表1所示，分析6種不同靜摩擦因數(shù)的顆粒的喂料狀態(tài)，喂料頭轉(zhuǎn)速ω分別為10，20，30，40，50，60，70，80，90 r/min，顆粒半徑r=2 mm，堆積密度ρ=0.7 t/m3，喂料頭轉(zhuǎn)速ωf為0～90 r/min. 螺旋管壁里面的螺旋轉(zhuǎn)速根據(jù)文獻(xiàn)[12]中的最佳螺旋轉(zhuǎn)速定義得出.

2.3 最佳曲面分析

利用氣固兩相DEM模型及EDEM+Fluent耦合仿真分析，分別考察以上5種曲面，并比較不同摩擦因數(shù)下不同曲面的輸送量的對比情況.

當(dāng)在不同的顆粒與喂料頭的摩擦因數(shù)下，各曲面在不同轉(zhuǎn)速下的輸送量分別如圖8所示.

由圖8中可以看出，隨著顆粒與喂料頭摩擦因數(shù)的增大，推薦使用漸開線更大些的曲面，因?yàn)槟Σ烈驍?shù)增大后，顆粒受到喂料葉片的摩擦力越大，顆粒從喂料口濺出量減小，所以曲面漸開線大些好，進(jìn)料量大，但是曲面漸開線不是越大越好，太大了也容易導(dǎo)致顆粒從喂料口濺出，最大漸開線的曲面1的輸送量并不比曲面2大，在任何摩擦因數(shù)下都不建議使用曲面1. 同理，隨著顆粒與喂料頭摩擦因數(shù)的減小，推薦使用漸開線更小些的曲面，因?yàn)槟Σ烈驍?shù)減少后，顆粒受到喂料葉片的摩擦力越小，顆粒從喂料口的濺出量增大，所以曲面漸開線小些好，但是曲面漸開線不是越小越好，太小了顆粒進(jìn)入喂料口的量減小，所以最小漸開線的曲面5的輸送量并不比曲面4大，在任何摩擦因數(shù)下都不建議使用曲面5. 不同摩擦因數(shù)下推薦使用不同的導(dǎo)擋葉曲面. 當(dāng)μ≤0.2時,推薦使用曲面4, 當(dāng)0.2<μ≤0.5時,推薦使用曲面3; 當(dāng)μ>0.5時,推薦使用曲面2.

圖8 不同摩擦因數(shù)下的5種曲面的輸送量Fig.8 Throughput of five kinds of curve surface in different friction coefficient

3 不同曲線過料面積及輸送量計算

3.1 過料面積

由前面顆粒進(jìn)料的動力學(xué)分析，可知

Qt=Qin-Qout,

(19)

式中:Qt為喂料頭的輸送量;Qin為喂料量，即通過喂料窗戶進(jìn)入喂料頭里面的顆粒量;Qout為濺出量，即隨著喂料頭的轉(zhuǎn)動，從喂料窗口濺出喂料頭的顆粒輸送量，它與導(dǎo)擋葉的漸開度和顆粒與喂料葉片的摩擦因數(shù)有關(guān).

喂料量Qin表示為

Qin=60V,

(20)

式中V為每分鐘喂料的體積.

V=2πRfDfnωf，

(21)

式中:Df為喂料口面積;Rf為喂料口半徑;n為喂料進(jìn)口個數(shù);ωf為喂料頭轉(zhuǎn)速.

由此可得，喂料頭的喂料量為

Qin=120πnDfRfωf，

(22)

單個喂料口的面積Df如圖9中的ABCD所示，垂直輸送管的結(jié)構(gòu)中心線為DE，q1曲線的終點(diǎn)A點(diǎn)與中心E的連線為AE，連線AD與喂料頭的m曲線組成一個截面，該截面被定義為喂料頭上的過料面，可近似地認(rèn)為是一個梯形.

其中，DE為喂料口半徑，AE為角度為65°時q1曲線的值，DE和AE之間的角度相差45°，AD的值可由AE和DE求出.BC線段為曲線m上兩點(diǎn)之間的線段，B點(diǎn)為角度為60°時曲線m的值，即為梯形中AB的值，C點(diǎn)為角度為120°時曲線m的值，即為梯形中CD的值. 前面已知q曲線、m曲線的方程，從而可求出梯形ABCD中各線段的值，并求出梯形ABCD的面積，即單個喂料面的面積.

圖9 過料面積Fig.9 Single material area

由q1、q2曲線方程(12)～(16)、s曲線方程(17)、m曲線方程(18)，可得出不同的5種曲線的過料面積，如表2所示.

表2 各曲線的過料面積Tab.2 Material area of all kinds of curve

3.2 輸送量公式

設(shè)Qt=kQin，則分別計算不同摩擦因數(shù)下的喂料量Qin，并得出不同摩擦因數(shù)下的輸送量仿真值Qt.

由前面分析可得，當(dāng)摩擦因數(shù)為0.1，0.2時，選擇曲面4，并由表3得出喂料面積為0.018 078 m2，當(dāng)摩擦因數(shù)為0.3，0.4，0.5時，選擇曲面3，并得出喂料面積為0.018 257 m2，當(dāng)摩擦因數(shù)為0.6時，選擇曲面2，并得出喂料面積為0.018 648 m2，然后由式(22)，可計算出各喂料頭轉(zhuǎn)速下的喂料量Qin，并與前面所得的仿真輸送量值進(jìn)行比較，得到不同摩擦因數(shù)下的喂料量與輸送量比較圖，如圖10所示. 由圖10可以看出，隨著摩擦因數(shù)的增大，輸送量和喂料量之間的差距越來越小，即顆粒的濺出量隨著顆粒與喂料頭之間的摩擦因數(shù)的增大而減小.

圖10 不同摩擦因數(shù)下的喂料量與輸送量對比圖Fig.10 Comparison chart between feeding quantity and transport quantity under different friction coefficient

從式(22)中可得出，在同一摩擦因數(shù)和同一曲面下，喂料量和喂料頭轉(zhuǎn)速成線性關(guān)系，設(shè)Qin=aωf，則根據(jù)式(22)，可得a=120πRfDfn. 從圖10中可以看出，輸送量和喂料頭轉(zhuǎn)速的關(guān)系也基本成線性關(guān)系，設(shè)Qt=bωf，對圖13中的每個仿真輸送量值進(jìn)行擬合可得出b. 設(shè)k=b/a，則不同摩擦因數(shù)下的k便可以計算出，計算結(jié)果如圖11所示.

擬合函數(shù)為

k=-0.075μ2+0.273μ+0.782.

(23)

所以喂料頭的輸送量為

Qt=(-0.075μ2+0.273μ+0.782)×

(120πnDfRfωf).

(24)

已知垂直螺旋輸送機(jī)的輸送量時，由式(24)，可得喂料頭的轉(zhuǎn)速為

ωf=

(25)

圖11 不同摩擦因數(shù)下的k值Fig.11 Value of k under different friction coefficient

4 結(jié) 論

針對相對旋轉(zhuǎn)式喂料頭進(jìn)行工作機(jī)理研究和動力學(xué)分析，分析了其氣固兩相DEM模型及其仿真方法，以某種高效的喂料頭為研究對象，研究喂料頭導(dǎo)擋葉5種曲面的曲線方程，并利用EDEM+Fluent耦合仿真方法分析喂料頭內(nèi)的顆粒狀態(tài)，比較不同摩擦因數(shù)下不同曲面的輸送量的對比情況，得到不同摩擦因數(shù)下推薦使用的喂料頭導(dǎo)擋葉曲面. 進(jìn)行喂料頭內(nèi)輸送量和喂料量的分析，計算導(dǎo)擋葉5種曲面的喂料窗口的過料面積，進(jìn)而推導(dǎo)出喂料頭的輸送量公式，得到喂料頭的轉(zhuǎn)速公式，對喂料頭的設(shè)計有非常重要的指導(dǎo)意義.