榮文杰， 瞿林濤， 李嘉鈺， 李寶寬

（1.東北大學(xué) 冶金學(xué)院， 沈陽 110819； 2.礦冶過程自動控制技術(shù)國家重點實驗室， 北京 102628）

滾筒是在冶金工業(yè)中被廣泛應(yīng)用的顆粒物料處理設(shè)備，如回轉(zhuǎn)窯、干燥機(jī)等［1］.以冶煉鎳鐵合金工藝的回轉(zhuǎn)窯為例，在生產(chǎn)過程中，紅土鎳礦、石灰石和煙煤等顆粒狀材料進(jìn)入回轉(zhuǎn)窯后，會隨著窯體轉(zhuǎn)動而不斷翻滾.但由于這些顆粒物料存在形狀、尺寸、密度、表面粗糙度等多方面的差異，顆粒在運動過程中呈現(xiàn)出混合與偏析的現(xiàn)象.不同顆粒的混合程度直接影響物料間的熱質(zhì)傳遞過程，進(jìn)而影響產(chǎn)品質(zhì)量，因此滾筒內(nèi)的顆?；旌吓c偏析過程一直備受研究者們的關(guān)注.

近年來，很多學(xué)者針對滾筒內(nèi)顆?；旌吓c偏析過程開展了物理實驗和數(shù)值模擬研究［2］.與流體混合類似，顆粒混合的3 個主要機(jī)理為對流混合、擴(kuò)散混合和剪切混合［3］.不同物性的顆粒在混合過程中會出現(xiàn)偏析現(xiàn)象，其中顆粒直徑和密度的差異是引起顆粒偏析的主要原因.由顆粒直徑不同而引起的顆粒偏析現(xiàn)象的機(jī)理包括滲透、流動和振動，而浮力是不同密度顆粒偏析的機(jī)理.離散元模型（discrete element model， DEM）是常用的數(shù)值模擬方法，該模型在計算過程中可以追蹤到每個顆粒的位置和速度，從而得到詳細(xì)的顆粒信息.由于它具有能夠預(yù)測實驗難以測量的顆粒信息，以及比連續(xù)性模型更直觀顯示顆粒的優(yōu)勢，故得到了廣泛的應(yīng)用［4］.Zhu 等［5-6］對DEM 的相關(guān)理論發(fā)展及其在顆粒系統(tǒng)的應(yīng)用進(jìn)行了總結(jié)，他們認(rèn)為DEM 是研究顆粒運動行為的有效方法.Yang 等［7］應(yīng)用DEM 對不同轉(zhuǎn)速滾筒內(nèi)同一直徑顆粒的宏觀和微觀行為進(jìn)行研究，展示并對比了多種顆粒流動狀態(tài).胡陳樞等［8］采用DEM 對滾筒內(nèi)不同直徑二元顆粒的運動過程進(jìn)行模擬，提出了表觀混合指數(shù)和動態(tài)混合指數(shù)的概念，并將它們用于計算顆粒的徑向混合程度.Xie 等［9］通過研究滾筒內(nèi)顆粒床的活動區(qū)發(fā)現(xiàn)，隨著弗勞德數(shù)的增大，活動區(qū)的顆粒速度逐漸增加.葛良等［10］采用DEM 模擬波形滾筒內(nèi)顆?；旌线^程發(fā)現(xiàn)，在波形滾筒內(nèi)使用較低的轉(zhuǎn)速將達(dá)到與圓形滾筒較高轉(zhuǎn)速下相同或相近的顆?；旌铣潭?Fernandes等［11］對顆粒的摩擦系數(shù)、恢復(fù)系數(shù)、楊氏模量和泊松比等輸入?yún)?shù)進(jìn)行不確定性的量化研究，他們發(fā)現(xiàn)摩擦系數(shù)是滾筒內(nèi)顆粒流動最主要的影響參數(shù).Zhang 等［12］研究了具有提升板的滾筒內(nèi)顆粒運動情況并與PIV 測量結(jié)果進(jìn)行對比，結(jié)果表明顆粒速度隨滾筒轉(zhuǎn)速的增加及填充比的降低而增加.Dury 等［13］研究了滾筒內(nèi)不同直徑顆粒的徑向偏析與滾筒轉(zhuǎn)速的關(guān)系，提出了描述最終偏析數(shù)量和偏析速度的參數(shù).董彥龍［14］研究了滾筒內(nèi)不同密度顆粒的偏析過程，發(fā)現(xiàn)密度比越大，顆粒偏析越明顯.Xu 等［15］研究了顆粒尺寸、密度和摩擦對滾筒內(nèi)顆?；旌系挠绊懀Y(jié)果表明顆粒尺寸和密度對顆?；旌系挠绊戄^大，但對顆粒摩擦的影響較小.Li 等［16］對橢圓柱形滾筒內(nèi)不同顏色的顆?；旌线^程進(jìn)行研究，他們發(fā)現(xiàn)在橢圓柱形滾筒內(nèi)，轉(zhuǎn)速低時顆粒的混合得到改善，轉(zhuǎn)速較高時顆粒在水平方向上的偏析得到加強.趙永志等［17］通過控制滾筒的轉(zhuǎn)速和填充比分別得到了月亮模式和太陽模式的顆粒徑向偏析.Maione 等［18］還對滾筒內(nèi)圓柱形顆粒的軸向偏析現(xiàn)象進(jìn)行了研究.上述研究展示了均勻或者二元顆粒在滾筒內(nèi)的運動過程和混合及偏析現(xiàn)象，并揭示了不同性質(zhì)顆粒的混合及偏析機(jī)理.

當(dāng)前的研究大多針對均勻顆?；蛘叨w粒，然而工業(yè)中應(yīng)用的顆粒物料往往更加復(fù)雜，一般是多種不同性質(zhì)的顆?；旌显谝黄?因此，本文中針對不同直徑的三元顆粒建立離散元模型，對3 種不同直徑顆粒在同體積和同數(shù)量工況下的顆粒離散分布、顆粒床動態(tài)休止角及顆粒速度分布進(jìn)行探究，以此開展多元顆粒混合過程的模擬計算，為研究更接近實際的滾筒內(nèi)顆?；旌线^程提供參考.

1 數(shù)學(xué)模型及計算條件

1.1 離散元模型

DEM 最早由Cundall 和Strack 在20 世紀(jì)70年代提出［19］，近年來逐漸得到廣泛應(yīng)用.DEM 利用牛頓經(jīng)典力學(xué)方程，對相鄰元之間因接觸產(chǎn)生的作用力進(jìn)行計算.每個單元受力后在合力和合力矩的作用下發(fā)生形變和位移，通過牛頓力學(xué)方程和運動定律的循環(huán)計算，就可以追蹤到每個顆粒的位置和速度.離散單元接觸模型主要有硬球模型和軟球模型.硬球模型假設(shè)顆粒與顆粒之間的碰撞是瞬間發(fā)生的.軟球模型則假設(shè)顆粒與顆粒之間發(fā)生的碰撞是在某個連續(xù)的時間范圍內(nèi)，通過重疊量的大小以及牛頓動力學(xué)定律來計算出顆粒間的接觸力等狀態(tài)參數(shù).在軟球模型模擬顆粒碰撞的過程中，將兩個離散單元顆粒之間的碰撞作用在法向和切向方向都抽象為彈簧、阻尼器及滑動摩擦器的共同作用，通過彈性、阻尼及滑移的力學(xué)機(jī)理考慮顆粒間的相互作用.圖1 為顆粒接觸受力示意圖.

圖1 顆粒接觸示意圖Fig.1 Schematic of two particles in contact

在時間t內(nèi)，顆粒i的位置、線速度和角速度可依據(jù)牛頓定律進(jìn)行計算.

式中：m，I，v，ω和g分別表示質(zhì)量（kg）、轉(zhuǎn)動慣量（kg·m2）、線速度（m／s）、角速度（rad／s）和重力加速度（m／s2）；xi表示顆粒i的位置；Fcnij和Fctij分別表示顆粒法向接觸力和切向接觸力（N）；Fdnij和Fdtij分別表示顆粒法向阻尼力和切向阻尼力（N）；Ttij和Trij分別表示顆粒受到切向力產(chǎn)生的扭矩和滾動摩擦力產(chǎn)生的扭矩（N·m）.DEM 廣泛用于顆粒運動過程計算，詳細(xì)的力和力矩計算公式可以參考文獻(xiàn)［5-9］.

1.2 邊界條件及模擬參數(shù)

本文中根據(jù)Parker 等［20］開展的滾筒內(nèi)顆粒運動實驗建立了三維滾筒離散單元模型.實驗中滾筒長度為650 mm，但為了減少計算量，本文中模擬了一個長度為16 mm 的短滾筒，并為滾筒兩端設(shè)置周期性邊界條件來保證計算結(jié)果的可靠性，該方法已在文獻(xiàn)［7］中得到證實.滾筒內(nèi)放置的顆粒直徑為3 mm，密度為2 500 kg／m3，顆粒之間的接觸參數(shù)及模型其他參數(shù)見表1.滾筒內(nèi)顆粒的運動模式隨著滾筒轉(zhuǎn)速不同一般分為6 種，即隨著轉(zhuǎn)速增加，依次呈現(xiàn)滑動、坍塌、滾動、傾瀉、崩塌和離心.在實際生產(chǎn)中顆粒運動往往處于滾動模式，本文中研究的滾筒內(nèi)顆粒運動在轉(zhuǎn)速42 r／min時處于滾動模式.在計算過程中，顆粒在滾筒內(nèi)的指定位置生成后受重力作用下落到底部，考慮到顆粒初始位置及分布并不影響最終的混合結(jié)果，顆粒生成位置設(shè)定為按空間坐標(biāo)均勻分布.但是，為了保證顆粒位置和速度在初始狀態(tài)是穩(wěn)定的，當(dāng)所有顆粒下落停穩(wěn)后，設(shè)置滾筒開始旋轉(zhuǎn).滾筒轉(zhuǎn)動6 圈后，顆粒運動呈現(xiàn)動態(tài)穩(wěn)定狀態(tài).為了保證計算結(jié)果的可靠性，下文討論的是滾筒旋轉(zhuǎn)到8 圈以后的計算結(jié)果.

表1 模型參數(shù)Table 1 Parameters in model

為了完成模型驗證，將單一直徑顆粒的計算結(jié)果與實驗結(jié)果進(jìn)行對比.在此基礎(chǔ)上，本文中進(jìn)一步探究了直徑為3，4，5 mm 的顆粒在同數(shù)量和同體積兩種工況下的混合過程.主要原因如下：一方面是考慮到該顆粒體系仍然能夠使用原計算網(wǎng)格；另一方面是考慮到滲透原理會導(dǎo)致不同直徑顆粒在混合過程中出現(xiàn)偏析，顆粒直徑差距過大時偏析規(guī)律顯而易見，而三元顆粒在直徑相近情況下的混合與偏析情況更值得研究.表2 列出了不同工況的顆粒數(shù)量.

表2 兩種工況下的顆粒數(shù)量Table 2 Particle numbers under two working conditions

1.3 離散中心及離散半徑

同一直徑顆粒在滾筒內(nèi)會呈離散分布，本文中定義“離散中心”來描述某一時刻同一直徑顆粒在整個離散空間分布的中心位置.以建立的滾筒模型中心位置為三維坐標(biāo)系原點，離散中心坐標(biāo)公式如下所示：

式中：（xi，yi，zi）為某一直徑的顆粒i在三維空間的坐標(biāo)；ni為該直徑顆粒個數(shù)；（）為離散中心在三維空間下的坐標(biāo).

以離散半徑表示某直徑顆粒距離其離散中心的平均值，可用此參數(shù)反映離散空間上顆粒群分布范圍的大小.以R代表離散半徑，其計算公式如下所示：

式中：n為顆?？倐€數(shù).

2 結(jié)果與討論

2.1 模型驗證

Parker 等［20］在滾筒內(nèi)顆粒運動的實驗中采用正電子發(fā)射粒子跟蹤方法（positron emission particle tracking， PEPT），測量了不同轉(zhuǎn)速下顆粒床動態(tài)休止角和顆粒沿顆粒床中心線的切向速度.本文中將轉(zhuǎn)速為42 r／min 的模擬結(jié)果與實驗結(jié)果進(jìn)行對比，結(jié)果如圖2 和圖3 所示.由圖可知，顆粒床動態(tài)休止角和顆粒切向速度結(jié)果分別與實驗值吻合，該結(jié)果驗證了模型的準(zhǔn)確性.

圖2 顆粒床動態(tài)休止角對比Fig.2 Comparison of dynamic angle of repose of particle bed

圖3 顆粒切向速度對比Fig.3 Comparison of particle tangential speed

2.2 顆粒離散分布

根據(jù)前述定義，3 種不同直徑顆粒在同數(shù)量和同體積兩種工況下的離散半徑如圖4 所示.從圖中可以看出，在滾筒旋轉(zhuǎn)至6 圈后，顆粒運動狀態(tài)基本穩(wěn)定.在兩種工況下，3 mm 顆粒離散半徑最小，4 mm 顆粒離散半徑居中，5 mm 顆粒離散半徑最大.在同數(shù)量工況下，3 種顆粒的離散半徑差距相當(dāng)，約為0.005 m；在同體積工況下，5 mm 顆粒與4 mm 顆粒的離散半徑差別明顯縮小，4 mm顆粒與3 mm顆粒的離散半徑差別略微增大.這是由于在同體積工況下，5 mm 顆粒的數(shù)量較少，而3 mm 顆粒的數(shù)量較多.

圖4 兩種工況下離散半徑的變化Fig.4 Revolution of discrete radius under two working conditions

圖5 更直觀地展示了3 種顆粒在同數(shù)量工況下的整體分布及各直徑顆粒的單獨分布.從圖中可以看出：3 mm 顆粒集中分布在顆粒床的內(nèi)部區(qū)域；5 mm 顆粒大部分分布在顆粒床的外部區(qū)域，并幾乎填滿了整個外圍部分；而4 mm 顆粒大部分填充在3 mm 顆粒與5 mm 顆粒之間.3 種顆粒形態(tài)分布大體上可看作為：5 mm 顆粒區(qū)域包裹4 mm顆粒，4 mm 顆粒區(qū)域包裹3 mm 顆粒，形成典型的“大包小”現(xiàn)象.

圖5 同數(shù)量工況下的顆粒分布Fig.5 Spatial distribution of particles under the working condition of the same number

圖6 展示了3 種顆粒在同體積工況下的整體分布及各直徑顆粒的單獨分布.從圖中可以看出：3 mm 顆粒主要集中在顆粒床的內(nèi)部區(qū)域，但由于其數(shù)量更多，有一部分顆粒向外圍進(jìn)行擴(kuò)散；4 mm顆粒大部分分布在3 mm 顆粒群的外部，并在一定程度上擠占了本來5 mm 顆粒所處的外圍區(qū)域；而5 mm 顆粒仍然大部分分布在顆粒床的外部區(qū)域.總的來說，3 種顆粒在空間上所占據(jù)的總體積是相等的，但依然出現(xiàn)了“大包小”的現(xiàn)象.

圖6 同體積工況下的顆粒分布Fig.6 Spatial distribution of particles under the working condition of the same volume

2.3 顆粒床動態(tài)休止角

圖7 示出了顆粒床動態(tài)休止角的計算結(jié)果.由圖可知，在轉(zhuǎn)速、填充比及顆粒接觸參數(shù)一致的情況下，同數(shù)量工況下顆粒床動態(tài)休止角約為27°，同體積工況下顆粒床動態(tài)休止角約為28°，兩種工況下的顆粒床動態(tài)休止角結(jié)果與Parker等［20］實驗測量得到的顆粒床的動態(tài)休止角幾乎是一致的.換言之，當(dāng)顆粒在滾筒內(nèi)處于滾動流態(tài)時，顆粒之間緊密接觸，并且以顆粒床整體的形式一起隨著滾筒的轉(zhuǎn)動而運動.這也進(jìn)一步驗證了Mellmann 等［21］對滾筒內(nèi)混合顆粒運動形態(tài)的研究結(jié)果：顆粒的大小及數(shù)量在固定的操作條件和接觸參數(shù)下并不會影響顆粒混合后的床層形態(tài).

圖7 顆粒床動態(tài)休止角對比Fig.7 Comparison of dynamic angle of repose of particle bed

2.4 顆粒速度分布

圖8 展示了3 種顆粒在同數(shù)量工況下的整體速度分布及各直徑顆粒的單獨速度分布.從圖中可以看出：顆粒床的傾斜表面較為平整，并有一層清晰可見的顆粒群從上往下快速滾落；與這些顆粒相比，其他顆粒則以較慢的速度在床層內(nèi)部和底部循環(huán)進(jìn)行爬升滾落；單一直徑顆粒靠近顆粒床傾斜表面的數(shù)量越多時，它所占據(jù)的大速度顆粒就越多.因此，大部分位于外圍區(qū)域的5 mm 顆粒群含有的大速度顆粒最多，4 mm 顆粒群次之，3 mm顆粒群中幾乎沒有速度較大的顆粒.同時，計算出各顆粒群速度的平均值如下：5 mm 顆粒群0.024 3 m／s；4 mm 顆粒群0.016 2 m／s；3 mm 顆粒群0.013 3 m／s.對于各顆粒群的平均速度而言，5 mm 最大，4 mm次之，3 mm 最小，這說明顆粒直徑與該直徑顆粒的平均速度呈正相關(guān).

圖8 同數(shù)量工況下的顆粒速度分布及其矢量形式Fig.8 Particle velocity distribution with the corresponding vector form under the working conditions of the same number

圖9 展示了3 種顆粒在同體積工況下的整體速度分布及各直徑顆粒的單獨速度分布.從圖中可以看出， 與同數(shù)量工況類似，顆粒床表面有一層清晰可見的顆粒群從上往下快速滾落.3 mm顆粒群整體速度偏低，而5 mm 顆粒群所含有的大速度顆粒變少，甚至沒有4 mm 顆粒群所含有的大速度顆粒多.這是因為5 mm 顆粒群的外圍空間被4 mm顆粒群擠占了一部分.各顆粒群平均速度結(jié)果如下：5 mm 顆粒群為0.019 2 m／s；4 mm顆粒群為0.016 9 m／s；3 mm 顆粒群為0.014 1 m／s.與同數(shù)量工況的結(jié)果類似，顆粒直徑與該直徑顆粒的平均速度呈正相關(guān).

圖9 同體積工況下的整體顆粒速度分布及其矢量形式Fig.9 Particle velocity distribution with the corresponding vector form under the working conditions of same volume

3 種直徑顆粒在同數(shù)量和同體積兩種工況下，全部顆粒的平均速度均為0.017 9 m／s.因此，當(dāng)滾筒的操作條件和顆粒接觸參數(shù)并未改變時，滾筒壁面?zhèn)鬟f給其內(nèi)部顆粒的動量不受顆粒數(shù)量和直徑的影響.

3 結(jié) 論

（1）3 種直徑顆粒在同數(shù)量和同體積兩種工況下，顆粒分布均出現(xiàn)了“大包小”的現(xiàn)象，顆粒的直徑越大，對應(yīng)的離散半徑越大.在同數(shù)量工況下，3 種直徑顆粒的離散半徑差距相當(dāng)；在同體積工況下，5 mm 顆粒的離散半徑與4 mm 顆粒的離散半徑差別明顯縮小，4 mm 顆粒的離散半徑與3 mm顆粒的離散半徑差別略微增大.

（2）顆粒的直徑大小及數(shù)量在操作條件和接觸參數(shù)不變的情況下并不會影響到滾筒內(nèi)顆粒床的動態(tài)休止角大小.

（3）3 種直徑顆粒在同數(shù)量和同體積兩種工況下，顆粒直徑均與該直徑顆粒的平均速度呈正相關(guān)關(guān)系；在同體積和同數(shù)量工況下顆粒的最大速度分別為0.143 m／s 和0.109 m／s；在兩種工況下顆粒整體的平均速度均為0.017 9 m／s.