垂直彎管段高壓密相煤粉氣力輸送特性的數值模擬

裴 宇，熊源泉，周海軍

（東南大學能源熱轉換及其過程測控教育部重點實驗室，南京210096）

高壓超濃相煤粉氣力輸送是大型氣流床高壓氣化的關鍵技術之一，由于受到目前試驗方法和條件的限制，還不能實時準確地測量管道中氣固兩相的流動參數，因此對高壓密相煤粉氣力輸送過程進行數值模擬研究有利于輸送系統的設計與優化.

類比稠密氣體分子運動論，Savage等、Jenkins等、Shahinpour 等、Lun 等、Johnson 等和Syamlal發展了顆粒動力學理論，引入顆粒碰撞恢復系數的概念和稠密氣體的Boltzmann輸運方程，修正了經典的分子運動論對兩顆粒碰撞的處理，推導出顆粒相應力的表達式和適合于描述顆粒流動的封閉方程組.顆粒動力學理論的引入是近幾十年來雙流體模型最重要的發展.在高壓密相模擬中，水平管和垂直管的研究取得了一定進展［1-2］，但彎管方面的模擬研究鮮有報道.Mcglinchey等［3］在顆粒動力學基礎上提出了常壓下90°彎管的計算模型（以下簡稱Mcglinchey彎管模型），該模型中忽略了摩擦應力，而且固相壁面采用的是無滑移條件，模擬的壓降結果誤差在10%～90%內.Ma等［4］運用Mcglinchey彎管模型研究了粒徑對水平轉垂直向上彎管壓降的影響，模擬得到的整體管道壓降與試驗較吻合，但就彎管段而言，管道壓降誤差高達50%.Vashisth等［5］采用鏡面系數和顆粒與壁面碰撞恢復系數來表示壁面剪切和擬熱流條件，對Mcglinchey彎管模型進行了修正，模擬獲得的彎管固相體積分數和速度分布云圖與實際較吻合，但沒有對壓降進行分析.

在本文高壓密相煤粉氣力輸送管道內，離心加速度垂直指向管道壁面外側，而且遠大于重力加速度，是垂直彎管氣力輸送過程中影響氣固流動的最主要因素.在輸送過程中，離心力的作用導致彎管壁面外側煤粉堆積形成堆積層，并以顆粒群的形式運動，而彎管壁面內側顆粒呈現懸浮流動.懸浮流動時，固相應力主要來源于顆粒之間的瞬時碰撞和平移產生的動量傳遞，這在顆粒動力學理論中已闡述清楚；在顆粒堆積層內，顆粒間除了瞬時碰撞和平移外，由于離心力作用，沉積層內的顆粒相互擠壓且長時間接觸產生了摩擦應力，并通過摩擦應力進行顆粒間動量和能量的傳遞與耗散.由于離心力較大，堆積層頂部體積分數接近自然堆積體積分數，顆粒與顆粒以及顆粒與壁面都保持長時間接觸，摩擦應力成為堆積層內固相應力和壁面剪切應力的重要部分，必須加以充分考慮.

以往文獻在運用運動顆粒動力學解決彎管問題時，在顆粒內和顆粒與壁面處忽略了摩擦應力的影響，造成數值模擬產生較大的誤差.為此，筆者在顆粒動力學的基礎上考慮了摩擦應力的影響，通過一體化管道對垂直向上轉水平彎管（其曲率半徑為0.2m，管道內徑為10mm，固相速度為5～8m／s）的高壓密相煤粉氣力輸送進行模擬，并結合多組試驗工況，驗證了模型的準確性，同時分析了垂直彎管內氣固流動情況以及補充風體積流量對輸送特性的影響.

1 數學模型

1.1 控制方程

連續方程

動量方程

式中：φ為各相體積分數；ρ為密度；v為速度矢量；τ為應力張量；pg為氣相壓力；Fsg表示氣固相間作用力；g為重力加速度；下標g、s分別表示氣相、固相.

在密相煤粉氣力輸送過程中氣固相間作用構成復雜，由于氣固速度差而產生的曳力是相間作用的主導因素，所以本文模型中的氣固相間作用只考慮曳力.

式中：β為氣固相間的曳力系數，大量參考文獻采用Gidaspow 提出的曳力模型［6］（以下簡稱Gidaspow模型）.

當φg≤0.8時，β采用Ergun提出的適用于密相的曳力模型：

式中：μg為氣體黏性系數；ds為顆粒直徑.

當φg＞0.8時，β采用Wen等提出的適用于稀相的模型：

式中：CD為流體阻力系數.

但φg＞0.8（即φs≤0.2）時并不完全是稀相流動，而且Gidaspow 模型在φg＝0.8 時不具有連續性，因此Lu等［7］對Gidaspow 模型進行了修正，提出Huilin-Gidaspow 模型：

1.2 顆粒動力學理論

Savage等認為顆粒相內的應力可以看成是動力學應力與摩擦應力之和，即，其中為顆粒動力學應力張量，為摩擦應力.

式中：為顆粒動力學應力張量中的固相壓力；為顆粒相黏性系數；λks為顆粒相體積黏性系數；Ss為顆粒相速度變形率張量；I為單位張量.

固相壓力由顆粒運動引起的動力效應和顆粒碰撞引起的碰撞效應組成：

顆粒相黏性系數包括顆粒相碰撞黏度和動力黏度2部分：

式中：Θs為顆粒擬溫度；ess為顆粒間碰撞恢復系數；g0，ss為顆粒的徑向分布函數.

顆粒相體積黏性系數表示顆粒相對壓縮的懸浮阻力，其表達式為：

徑向分布函數為：

顆粒擬溫度的輸運方程如下：

式中：qs為擬熱流密度，在顆粒相內，表示沿顆粒擬溫度Θs梯度方向的擬熱流密度，其中kΘs為擴散系數；表示顆粒相內由剪切應力產生的顆粒脈動；γΘs表示顆粒與顆粒之間碰撞導致的顆粒脈動的耗散；φsg表示氣體與顆粒間脈動的交換［1-2，7］.

Johnson等的摩擦應力模型認為固相體積分數達到0.5后摩擦應力才起作用，但高壓密相煤粉氣力輸送中，固相體積分數低于0.5時也存在摩擦應力.因此，筆者認為固相體積分數達到0.1后摩擦應力開始起作用，并根據實際情況修正經驗參數，具體表達式如下：

式中：pf為摩擦應力；φ為顆粒間的摩擦角；經驗參數F＝0.1，r＝2，s＝5，φs，min＝0.1，φs，max＝0.55.

在氣固模擬中氣固湍流輸運方程較統一，具體可參考文獻［1］、文獻［7］和文獻［8］.

1.3 邊界條件

1.3.1 入口和出口條件

對于氣相，入口處徑向速度分布為光滑管充分發展湍流流動，出口處設置為壓力出口.

式中：Ug為表觀氣速；φs，in為入口處固相體積分數；D為管道半徑.

對于固相，入口處徑向速度為均勻分布，出口處設置為壓力出口.

其中，qm，s為固相質量流量.

1.3.2 壁面條件的修正

對于氣相，壁面采用無滑移邊界條件.對于固相，壁面剪切應力僅考慮顆粒與壁面通過碰撞產生的顆粒動力學應力.

式中：usw為顆粒與壁面的滑移速度，下標w 表示壁面；φ′為鏡面系數，與管道材料、顆粒類型和壁面加工條件等因素有關，本文彎管段的鏡面系數設定為0.05，水平段和垂直段的鏡面系數設定為1×10-5.

式中：esw為顆粒與壁面碰撞恢復系數，取值為0.6.

Jenkins等［9］認為顆粒流動分為2種區域：快速剪切區域和緩慢剪切區域.快速剪切區域中固相體積分數不高，固相應力主要來源于顆粒碰撞與平移導致的動量傳遞；緩慢剪切區域中固相體積分數高，單個顆粒與周圍接觸面緊密接觸產生相互作用，固相應力主要來源于顆粒與周圍接觸面長時間接觸而產生的靜摩擦應力.離心加速度是影響流動的最主要因素，本文中垂直彎管的離心加速度遠大于重力加速度，這使得顆粒與壁面的接觸壓力較大，此時摩擦應力對壁面剪切應力的影響不容忽視.Johnson等［10］在流化床模型中將摩擦應力表示為·tanδ，其中pf為摩擦正應力即摩擦壓力，δ為顆粒與壁面的摩擦角.引用該表達式，可得到考慮壁面摩擦應力的壁面剪切應力表達式：

此時，壁面擬熱流密度可采用式（24）進行計算.

2 模擬條件

如圖1所示，模擬中采用一體化管道，管道內徑為10mm，其中垂直管長2m、水平管長4m、垂直彎管曲率半徑為0.2m，管道截面網格中軸向采用六面體cooper，軸間距為2mm，共計56.8萬個網格.

圖1 管道示意圖和截面網格Fig.1 Piping schematic and the cross-section grid

3 試驗裝置

高壓密相煤粉氣力輸送試驗裝置見圖2.高壓氮氣經緩沖罐分為流化風、充壓風和補充風3路，發料罐內的煤粉經流化風、充壓風驅動進入輸送管道，在發料罐出口引入補充風調節輸送過程中的固氣比（即固相與氣相質量比），充壓風用來維持發料罐內的壓力.進入受料罐的含粉氣流經布袋除塵器分離煤粉后放空，受料罐壓力由電動調節閥調節.在輸送管的水平段、垂直段和垂直彎管段均設有差壓變送器來測量壓差.在高壓發料罐上安裝電子秤，測量煤粉質量流量.試驗裝置上壓力、流量、壓差和煤粉質量流率數據均采用計算機數據采集系統采集.輸送管道內徑為10mm，輸送距離為53m.

圖2 試驗系統和發料罐示意圖Fig.2 Schematic diagram of the experimental system and supply tank

4 結果與討論

選取內蒙褐煤的輸送試驗進行模擬，試驗中保持背壓2.5 MPa、總輸送壓差0.5 MPa和流化風體積流量0.4m3／h不變，通過改變補充風體積流量來改變輸送固氣比；模擬中驗證了不同固氣比下模型的準確度，其中褐煤全水分質量分數為5.45%，密度為1 496kg／m3，平均粒徑為204.3μm，選取管程中間段的6m 為幾何模型，具體參數設置見表1.

表1 模型參數設置Tab.1 Model parameter settings

4.1 考慮壁面摩擦應力對模擬結果的影響

分別采用未考慮摩擦應力和考慮摩擦應力的壁面條件對4組工況進行模擬對比，壓降模擬結果見表2.

表2 考慮摩擦應力對垂直彎管壓降的影響Tab.2 Effect of friction stress on pressure drop of the vertical elbow

注：1）表示未考慮摩擦應力時模擬得到的垂直彎管壓降與試驗值的誤差；2）表示考慮摩擦應力時模擬得到的垂直彎管壓降與試驗值的誤差.

由表2可知，4組工況下采用未考慮摩擦應力的壁面條件模擬得到的壓降與試驗值間的誤差均較大，遠超出合理范圍，采用考慮摩擦應力的壁面條件模擬得到的壓降均大于采用未考慮摩擦應力的壁面條件模擬得到的壓降，計算準確性得到較大改善，基本達到了模擬精度的要求.

針對是否考慮摩擦應力對壁面剪切應力的影響較大，以補充風體積流量為0.6m3／h的工況為例，圖3給出了45°截面壁面剪切應力.從圖3 可以看出，未修正前的壁面剪切應力只考慮顆粒與壁面的碰撞，最大值為180Pa；而修正后的壁面剪切應力引入了摩擦應力的影響，最大值增大到240Pa，壁面的流動阻力增大，相應的垂直彎管壓降從9.03 kPa增大到12.25kPa.

圖3 考慮摩擦應力對45°截面壁面剪切應力的影響Fig.3 Effect of friction stress on the wall shear stress of 45°crosssection

圖4給出了45°截面沿徑向分布的固相速度和固相湍動能.從圖4（a）可以看出，由于修正后的壁面條件考慮了顆粒與壁面長時間接觸而產生的摩擦應力，從而增大了壁面剪切應力，增大了壁面流動阻力，造成顆粒密集區域速度減小.對于氣相而言，顆粒密集區域流動阻力增大，速度減小，而在顆粒稀疏區域加速帶動固相速度增大.從圖4（b）可以看出，由于壁面是固相湍動能產生和耗散的重要部分，考慮摩擦應力后，顆粒在壁面處的剪切應力增大，壁面擬熱流密度中的產生項隨之增大，壁面擬熱流密度也隨之增大，壁面產生的湍動能經顆粒間傳遞造成截面顆粒密集區域的湍動能增大，而在顆粒稀疏區域，由于氣固速度增大，顆粒間曳力增大，氣固湍動能交換更頻繁，固相湍動能隨之增大，因此考慮壁面摩擦應力后彎管內的固相湍動能均有所增大.

圖4 考慮摩擦應力對45°截面沿徑向分布的固相速度和固相湍動能的影響Fig.4 Effect of friction stress on the radial solid velocity distribution and turbulent kinetic energy of 45°cross-section

4.2 彎管整體流動過程的分析

以補充風體積流量為0.6m3／h工況為例，對彎管整體流動過程進行分析.圖5給出了沿彎管軸向的固相體積分數分布云圖.圖6給出了從彎管入口到中點每隔5°各截面固相體積分數云圖.從圖5和圖6可以看出，顆粒相在入口段中心體積分數達到最大值，進入彎管后由于離心力的作用，顆粒相逐漸向彎管外壁移動形成高體積分數區，內壁處低體積分數區逐漸擴大.

在到達15°截面之前，靠近外壁面的部分顆粒開始與壁面發生碰撞，反彈后以一定的動量沖向未發生碰撞的顆粒群，此時顆粒群依舊保持垂直向上的速度且煤粉質量較大，反彈的顆粒無法突破顆粒群，會隨著顆粒群再次向上運動，所以在15°截面之前外壁面的高體積分數區和內壁面的低體積分數區一直在擴大.在15°截面之后，顆粒群與外壁面發生碰撞并開始向內壁面反彈，由于低體積分數區中心氣流速度較大，攜帶能力較強，反彈到中心的顆粒被氣流迅速帶走，而反彈到壁面的顆粒開始聚集并向低體積分數區中心移動；在15°截面到30°截面之間，從圖6（d）～圖6（g）可以清晰地看出，低體積分數區兩側的固相體積分數開始升高，顆粒向低體積分數區中心發展.當發展到45°截面時，高體積分數區顆粒與壁面的碰撞、顆粒與載氣的曳力和離心力已經達到平衡，流形基本穩定，低體積分數區形成了倒“V”形的穩定區域.

圖5 沿彎管軸向的固相體積分數分布云圖Fig.5 Solid concentration contours along axial direction of elbow

圖6 垂直彎管各截面固相體積分數云圖Fig.6 Solid concentration cloud on various sections of vertical elbow

4.3 補充風體積流量對垂直彎管壓降的影響

補充風體積流量增大導致輸送管道補充風入口處的壓力升高，減小了發料罐與補充風入口處的壓差，可視為輸送物料運動的能量減少，所以煤粉質量流量與體積通量均逐漸減小.另一方面，發料罐出煤量減少，用于填充煤粉流出后罐內形成空隙的沖壓風體積流量減小，因此隨補充風體積流量的增大，表觀氣速并非線性增大.補充風體積流量較小時，補充風體積流量增大對出煤量的影響不大，充壓風體積流量變化不大，此時隨補充風體積流量的增大，表觀氣速增大較快；補充風體積流量較大時，補充風體積流量增大對出煤量的影響較大，出煤量和補充風體積流量明顯減小，此時隨著補充風體積流量的增大，表觀氣速增大趨勢減緩（見表3）.采用考慮摩擦應力的壁面條件模擬這4組工況，模擬中其他條件保持不變，改變入口處氣相平均速度、固相平均速度和煤粉體積分數.

圖7給出了垂直彎管壓降隨補充風體積流量的變化.從圖7可以看出，修正后的垂直彎管段模擬所得壓降與試驗值的誤差減小至20%以內，基本符合模擬要求.輸送中垂直彎管段壓降由氣相壓降和固相壓降2部分組成，其中固相壓降占主導地位，固相壓降又分為動能壓降、重力壓降和摩擦壓降，其表達式如下：

表3 不同補充風體積流量下對應的輸送參數Tab.3 Conveying parameters at different supplementary air flows

式中：Gs為煤粉質量通量；H為高度變化量；f為固相摩擦壓損系數；Us為固相表觀速度.

圖7 垂直彎管壓降隨補充風體積流量的變化Fig.7 Changes of vertical elbow pressure drop with supplementary air flow

4.4 補充風體積流量對垂直彎管流動特性的影響

選取垂直彎管段中點45°截面進行分析，各工況下的固相體積分數分布云圖見圖8.從圖8可以看出，隨著補充風體積流量的增大，截面煤粉質量通量減小，輸送固氣比降低，在云圖中表現為體積分數分布變窄；補充風體積流量增大造成表觀氣速隨之增大，從而增強了載氣與顆粒間的曳力作用，煤粉更容易被攜帶，高體積分數區逐漸減小.

圖8 不同補充風體積流量下垂直彎管段45°截面的固相體積分數云圖Fig.8 Solid concentration cloud on 45°cross-section at different supplementary air flows

圖9給出了各工況下垂直彎管段45°截面的固相湍動能和固相擬溫度徑向分布.從9（a）可以看出，垂直彎管段截面存在2個高固相湍動能區域：貼近壁面外側的高體積分數區和貼近壁面內側的低體積分數區.

隨著補充風體積流量的增大，表觀氣速增大，截面煤粉質量通量減小，在高體積分數區，補充風體積流量不大時，貼近壁面處的固相體積分數最大值不變（見圖8），此時固相速度增大導致顆粒與壁面的滑移速度增大，顆粒與壁面的剪切應力、壁面擬熱流密度qs｜w和固相湍動能增大，當補充風體積流量達到1.0m3／h時，貼近壁面處的固相體積分數最大值開始減小，此時顆粒與壁面的碰撞概率和摩擦應力均減小，壁面擬熱流密度qs｜w和固相湍動能也開始減小.

圖9 各工況下垂直彎管段45°截面固相湍動能和固相擬溫度徑向分布Fig.9 Solid-phase turbulent kinetic energy and granular temperature distribution on 45°cross-section under different conditions

在低體積分數區，氣固相間脈動傳遞強烈，表觀氣速增大對固相湍動能的增強效果更加明顯，在補充風體積流量達到1.0m3／h時，截面煤粉質量通量急劇減小，造成內壁面低體積分數區煤粉體積分數進一步減小，此時該區域湍流強度有所減弱.固相擬溫度曲線也存在著與固相湍動能相似的變化規律.

5 結 論

（1）考慮顆粒與壁面長時間接觸而產生的摩擦應力后，模擬所得的垂直彎管壓降與試驗值的誤差減小至20%以內，相比采用未考慮壁面摩擦應力時模擬所得的結果，計算精度得到較大改善.

（2）隨著補充風體積流量的增大，垂直彎管壓降先增大后減小，表觀氣速增大，垂直彎管壁面外側煤粉堆積減少，低體積分數區范圍增大，固相湍動能和固相擬溫度均有所增大，當補充風體積流量達到1.0m3／h時，由于固相體積分數急劇下降，固相湍動能和固相擬溫度均減小.

［1］蒲文灝，趙長遂.水平管加壓密相煤粉氣力輸送沙丘流的數值模擬［J］.動力工程學報，2011，31（2）：120-126.

PU Wenhao，ZHAO Changsui.Numerical simulation on dune flow in pressurized dense phase pneumatic conveying of pulverized coal in horizontal pipe［J］.Journal of Chinese Society of Power Engineering，2011，31（2）：120-126.

［2］蒲文灝，趙長遂，熊源泉，等.垂直管密相輸送的數值模擬［J］.動力工程，2008，28（1）：95-99.

PU Wenhao，ZHAO Changsui，XIONG Yuanquan，etal.Numerical simulation of the dense phase conveying in vertical pipe［J］.Journal of Power Engineering，2008，28（1）：95-99.

［3］MCGLINCHEY D，COWELL A，KNIGHT E A，et al.Bend pressure drop predictions using the Euler-Euler model in dense phase pneumatic conveying［J］.Particulate Science and Technology，2007，25（6）：495-506.

［4］MA A C，WILLIAMS K C，ZHOU J M，etal.Numerical study on pressure prediction and its main influence factors in pneumatic conveyors［J］.Chemical Engineering Science，2010，65（23）：6247-6258.

［5］VASHISTH S，GRACE J R.Simulation of granular transport of Geldart type-A，-B，and -D particles through a 90°elbow［J］.Industrial ＆ Engineering Chemistry Research，2012，51（4）：2030-2047.

［6］GIDASPOW D，JUNG J，SINGH R K.Hydrodynamics of fluidization using kinetic：an emerging paradigm 2002flour-daniel lecture［J］.Powder Technology，2004，148（2）：123-141.

［7］LU H L，GIDASPOW D，BOUILLARD J，etal.Hydrodynamic simulation of gas-solid flow in a riser using kinetic theory of granular flow［J］.Chemical Engineering Journal，2003，95（1／2／3）：1-13.

［8］PU W H，ZHAO C S，XIONG Y Q，etal.Threedimensional numerical simulation of dense pneumatic conveying of pulverized coal in a vertical pipe at high pressure［J］.Chemical Engineering ＆ Technology，2008，31（2）：215-223.

［9］JENKINS J T，SAVAGE S B.A theory for the rapid flow of identical，smooth，nearly，elastic spherical particles［J］.J Fluid Mech，1983，130（5）：187-202.

［10］JOHNSON P C，JACKSON R.Frictional-collisional constitutive relations for granular materials，with application to plane shearing［J］.J Fluid Mech，1987，176（9）：67-93.