甘德清 閆澤鵬 薛振林 劉志義 孫海寬

（1.華北理工大學(xué)礦業(yè)工程學(xué)院，河北唐山063200；2.河北省礦業(yè)開發(fā)與安全技術(shù)重點(diǎn)實(shí)驗(yàn)室，河北唐山063009）

近年來，淺部礦產(chǎn)資源逐漸減少，導(dǎo)致礦山開采深度不斷增加，伴隨而來的“三高”（高溫、高壓、高井深）問題日趨嚴(yán)重［1］。充填采礦法的有效應(yīng)用有助于解決深部采場(chǎng)的高地應(yīng)力問題，并以其綠色、環(huán)保、安全等特點(diǎn)，逐漸成為地下礦山的主要開采方法［2-3］。與其他充填方法相比，高濃度全尾砂充填工藝具有成本低、易于實(shí)現(xiàn)機(jī)械化的優(yōu)點(diǎn)［4］。管道輸送是高濃度全尾砂充填工藝的關(guān)鍵環(huán)節(jié)，是確保系統(tǒng)穩(wěn)定、高效運(yùn)行的前提［5］。管道輸送過程中常發(fā)生堵管、磨損等事故，限制了充填礦山的生產(chǎn)能力，有研究認(rèn)為管道阻力是影響充填管道穩(wěn)定性及管路堵塞的主要原因［6］。

影響管道輸送阻力的因素眾多，由于受到測(cè)量手段和試驗(yàn)裝置限制，常規(guī)的環(huán)管試驗(yàn)受到多種條件制約，粗略的測(cè)量值很難作為論證依據(jù)［7］。隨著計(jì)算機(jī)技術(shù)的發(fā)展，數(shù)值模擬（數(shù)值計(jì)算）成為研究流體流動(dòng)的重要手段［8］。甘德清等［9］通過Fluent軟件模擬得出料漿管道阻力損失與管徑呈反比。溫震江等［10］通過模擬驗(yàn)證得出，相同濃度時(shí)管道阻力損失與流量呈正相關(guān)。楊波等［11］通過ANSYS/FLOTRAN模擬料漿在充填管道內(nèi)的流動(dòng)特征，得出入口速度增大時(shí)，濃度對(duì)管道阻力的影響極小。謝翠麗等［12］應(yīng)用Fluent軟件對(duì)煤水兩相流動(dòng)過程中的表觀滑移現(xiàn)象進(jìn)行了研究，發(fā)現(xiàn)固相體積分?jǐn)?shù)較低時(shí)，出現(xiàn)超前滑移現(xiàn)象。

由于高濃度料漿屬于非牛頓體范疇，因而其在管道輸送過程中必然存在表觀滑移現(xiàn)象［13］，忽略壁面滑移效應(yīng)對(duì)管道阻力損失的影響得出的模擬結(jié)果，與充填料漿管道輸送的實(shí)際情況會(huì)有一定的偏差。為此，本研究基于Comsol軟件建立考慮壁面滑移效應(yīng)的管道輸送數(shù)值模型，分析管徑、灰砂比及濃度對(duì)管道輸送阻力的影響規(guī)律，并以某礦山為例，給出最佳的充填管路技術(shù)參數(shù)。

1 數(shù)學(xué)模型

1.1 模型假設(shè)

管道內(nèi)部主流區(qū)與滑移層接觸部分的流體速度分布是連續(xù)的，壁面滑移速度是剪切流動(dòng)區(qū)和滑移層連接處的速度與管道內(nèi)壁的速度差［14］。為便于分析，本研究對(duì)漿體在管內(nèi)的流動(dòng)行為進(jìn)行假設(shè)：①管道內(nèi)的流體為均質(zhì)體，滑移層內(nèi)流體在輸送管道內(nèi)各區(qū)域滑移層厚度相同；②管道輸送過程中不存在熱交換；③滑移層的狀態(tài)不會(huì)影響主流區(qū)的漿體。

1.2 理論模型

（1）滑移流動(dòng)的基本方程。根據(jù)管流中剪切速率的定義γ=4vR及管流基本公式，得出管徑為R的管道中料漿產(chǎn)生壁滑移時(shí)的基本流動(dòng)方程［15］為

式中，R為管道半徑，m；v為漿體平均流速，m/s；τ為壁面剪切應(yīng)力，Pa；τw為臨界剪切應(yīng)力，Pa；vslip為壁面滑移速度，m/s。

（2）黏性耗散。管道輸送過程存在的黏性耗散量是一個(gè)與速度梯度相關(guān)的物理量，可表示為

式中，φ為黏性耗散量；η為流體黏度，Pa/s；u、v、w分別為X、Y、Z軸方向的速度分量，m/s；T為料漿溫度，℃為耗散率平均值，%。

1.3 控制方程

對(duì)于不可壓縮的流體，Comsol軟件中控制流體流動(dòng)的基本方程主要有能量守恒方程、質(zhì)量方程和狀態(tài)方程，控制流動(dòng)的變量有流體壓力P、密度ρ、溫度T和流速u。

Navier-Stokes動(dòng)量守恒方程（忽略外部體積力）可表示為

式中，u0為料漿初始速度，m/s；ρ為料漿密度，kg/m3；g為重力加速度，9.8 m2/s；u、v、w分別為X、Y、Z軸的3個(gè)速度分量，m/s。

質(zhì)量守恒方程可表示為

式中，ρ為料漿密度，kg/m3；u為料漿速度，m/s；t為料漿輸送時(shí)間，s。

狀態(tài)方程可表示為

式中，ρ為料漿密度，kg/m3；P為流體壓力，Pa；T為溫度，℃。

2 模擬試驗(yàn)設(shè)計(jì)

2.1 試驗(yàn)材料

本研究全尾砂樣品取自河北某礦選礦廠，從底堆處取料，經(jīng)烘干后備用測(cè)定其密度為2 850 kg/m3。采用激光粒度儀對(duì)尾砂顆粒的粒級(jí)組成進(jìn)行測(cè)定，得到其中位粒度d60=17.53 μm。根據(jù)礦山充填開采經(jīng)驗(yàn)，以325#礦渣硅酸鹽水泥為膠凝劑，經(jīng)測(cè)定其密度為3 080 kg/m3。

2.2 方案設(shè)計(jì)

根據(jù)國內(nèi)礦山自流充填系統(tǒng)的發(fā)展情況［16］，選擇充填管路倍線5，管徑為100 mm、150 mm和200 mm。為了更加科學(xué)合理地得到試驗(yàn)數(shù)據(jù)，同時(shí)盡可能減小試驗(yàn)工作量，以灰砂比、料漿濃度及管徑為主要影響因素設(shè)計(jì)了考慮三因素三水平的正交試驗(yàn)，試驗(yàn)參數(shù)取值如表1所示。

根據(jù)雷諾數(shù)Re計(jì)算公式［17］，分別計(jì)算了不同條件下的雷諾數(shù)，發(fā)現(xiàn)雷諾數(shù)均小于2 300時(shí)，管內(nèi)流動(dòng)狀態(tài)屬于層流。設(shè)定邊界條件“入口—充分發(fā)展的流動(dòng)—平均速度—2 m/s”；Z方向設(shè)置重力為體積力；由于設(shè)置重力為體積力，設(shè)置管道出口邊界條件為0或?qū)恿髁鞒鰰?huì)導(dǎo)致模擬計(jì)算不收斂，需要設(shè)置出口為“開放邊界”且無黏滯應(yīng)力并加入壓力積分為0的約束；設(shè)置壁面條件為滑移速度（介于無滑移和滑移之間可以較好地表征壁面滑移效應(yīng)）。

2.3 幾何模型建立及網(wǎng)格剖分

本研究基于Comsol軟件三維空間維度對(duì)充填管路進(jìn)行建模，彎管曲率半徑為0.5 m，由于充填倍線為5，設(shè)置豎直管長(zhǎng)度為1.5 m，水平管長(zhǎng)度為9.5 m。Comsol軟件自帶的網(wǎng)格剖分工具有很強(qiáng)的編輯能力，用戶可根據(jù)需求調(diào)整網(wǎng)格的形狀、大小等。考慮到滑移層極薄，為使計(jì)算結(jié)果更精確，設(shè)置模型的邊界層數(shù)為5，邊界層拉伸因子為2，最終剖分單元數(shù)目為350 623。建模效果如圖1所示。

3 數(shù)值模型可靠性驗(yàn)證與結(jié)果分析

3.1 壁面滑移數(shù)值模型可靠性驗(yàn)證

吳愛祥等［18］經(jīng)物理試驗(yàn)研究，提出的壁面滑移作用下的管道阻力公式為

式中，ΔP為管道壓力損失，Pa；D和L分別為管道直徑和管道長(zhǎng)度，m；μp為高濃度料漿黏度，Pa/s；τy為料漿屈服應(yīng)力，N；δ為滑移層厚度，取δ=5 μm；μslip為滑移層黏度（黏度極低將其看作常溫狀態(tài)下的水，取μslip=1.05 mPa·s）。

文獻(xiàn)［19］提出的管道阻力損失計(jì)算的常用公式如下

式中，ij為充填料漿水力坡度，mH2O/m；i0為清水水力坡度，mH2O/m；g為重力加速度，9.8 m2/s；D為管道直徑，m；v為料漿流速，m/s；ρj和ρ0分別為料漿和清水的密度，kg/m3。

為驗(yàn)證Comsol數(shù)值模型在充填料漿管道輸送阻力損失研究方面的可靠性，將上述試驗(yàn)條件代入式（7）和式（8）進(jìn)行計(jì)算，對(duì)比分析數(shù)值模擬結(jié)果和公式計(jì)算結(jié)果，如表2所示。

由表2可知：式（8）計(jì)算結(jié)果與數(shù)值模型的計(jì)算結(jié)果相對(duì)誤差為48.1%～386.1%，誤差極大。原因是該式在多次試驗(yàn)的基礎(chǔ)上所得，試驗(yàn)條件對(duì)公式的適用性影響較大，且該式未考慮壁面滑移對(duì)管內(nèi)料漿流動(dòng)的影響；式（7）計(jì)算結(jié)果與數(shù)值模擬計(jì)算結(jié)果的相對(duì)誤差在20%以內(nèi)，說明本研究試驗(yàn)所建立的考慮壁面滑移效應(yīng)的數(shù)值模型相對(duì)可靠。

3.2 數(shù)值模擬結(jié)果分析

將數(shù)值模擬數(shù)據(jù)導(dǎo)入DPS數(shù)據(jù)處理軟件中，采用極差分析法研究料漿濃度、灰砂比、管徑對(duì)阻力損失的影響程度，并建立了阻力損失關(guān)于三者的回歸方程。試驗(yàn)結(jié)果如表3所示。

由表3可知：考慮壁面滑移效應(yīng)時(shí)，各因素對(duì)管道阻力的影響程度依次為管徑＞質(zhì)量濃度＞灰砂比。通過回歸擬合建立的阻力損失與濃度、灰砂比、管徑的函數(shù)關(guān)系式為

注：K1，K2，K3為每個(gè)因素各個(gè)水平下的極差總和為每個(gè)因素各個(gè)水平下的極差均值。

式中，Y為阻力損失，Pa/m；x1為質(zhì)量濃度，%；x2為尾砂質(zhì)量與水泥質(zhì)量的比值；x3為充填管徑，mm。回歸方程的相關(guān)系數(shù)R=0.999 99，調(diào)整后Ra=0.999 97，可見回歸方程較為顯著，方程可靠程度高。

4 充填料漿管道輸送阻力損失變化特征

4.1 管徑對(duì)輸送阻力的影響

管徑與管道阻力損失之間的關(guān)系如圖2所示，可以看出管徑與阻力損失呈負(fù)相關(guān)。圖2（a）表明管徑增大時(shí)，灰砂比對(duì)阻力損失影響的顯著程度減小，3條曲線逐漸發(fā)散；圖2（b）表明，灰砂比不變，濃度增大時(shí)，料漿流動(dòng)的“黏性”增大，輸送阻力越大。

高濃度料漿中尾砂顆粒遷移主要受到壁面剪切作用力影響，管徑增大，壁面剪切作用力減小，顆粒的遷移運(yùn)動(dòng)變緩，滑移效應(yīng)減弱，對(duì)摩擦阻力的抵消作用減小。管徑較小時(shí)，滑移層的減阻效果明顯，可以看出，管徑為100～150 mm時(shí)，管道輸送阻力快速降低；當(dāng)管徑增大到一定值時(shí)，料漿中的尾砂顆粒所受到的剪切作用力較弱，致其無法進(jìn)行遷移運(yùn)動(dòng)，不能形成覆蓋管道內(nèi)壁的滑移層，此時(shí)滑移作用效果較弱，隨著管徑增大，管道輸送阻力緩慢減小。

4.2 濃度對(duì)輸送阻力的影響

料漿濃度與管道阻力損失之間的關(guān)系曲線如圖3所示。由圖3可知：料漿濃度與管道阻力損失基本呈正相關(guān)。圖3（a）表明，管徑不變料漿濃度增加時(shí)，管道輸送阻力逐漸增大，同時(shí)可以看出灰砂比對(duì)輸送阻力影響較大，水泥含量越高，阻力損失值越大。圖3（b）反映灰砂比一定時(shí)濃度對(duì)阻力損失的影響規(guī)律，管徑越大，管道阻力損失越小，同時(shí)隨著濃度增大，圖中3條曲線開始聚攏，說明濃度變化導(dǎo)致輸送阻力的差異越來越小。

濃度較低時(shí)，料漿管道輸送過程中滑移層的形成主因?yàn)殪o態(tài)壁面損耗效應(yīng)［19］，滑移層主要由水構(gòu)成。較小的剪切力作用下就會(huì)引發(fā)滑移流動(dòng)，產(chǎn)生滑移效應(yīng)，對(duì)管道形成“潤(rùn)滑”作用，其與料漿輸送過程所受到的部分摩擦阻力相互抵消，隨著濃度增大，輸送阻力增長(zhǎng)速率加快；濃度較大時(shí)，滑移層由尾砂顆粒的遷移運(yùn)動(dòng)控制，尾砂顆粒需要在較大的剪切應(yīng)力作用下才會(huì)發(fā)生遷移運(yùn)動(dòng)，且濃度增大時(shí)，料漿流動(dòng)性變差，摩擦阻力增大，二者相互作用下，管道阻力損失快速增長(zhǎng)。結(jié)合圖3（b）分析發(fā)現(xiàn)，料漿濃度超過70%時(shí)阻力損失增長(zhǎng)過快，因此本研究推薦最大輸送濃度不應(yīng)超過70%。

4.3 灰砂比對(duì)輸送管道阻力的影響

砂灰比與管道阻力損失之間的關(guān)系如圖4所示，可以看出砂灰比與管道阻力損失呈負(fù)相關(guān)。圖4（a）表明，管徑一定時(shí)，管道阻力損失隨著砂灰比增大而呈線性減??；圖4（b）表明料漿濃度一定時(shí)，料漿中的水泥含量降低，管道阻力損失快速下降，當(dāng)水泥含量低于一定值后，管道阻力損失減小的速率快速減小。

水泥遇水發(fā)生水化反應(yīng)的過程會(huì)大量消耗料漿中的自由水，使得高濃度料漿濃度增大，表觀黏度及屈服應(yīng)力增加，導(dǎo)致料漿的流動(dòng)性變差，料漿的黏性耗散量增大。同時(shí)由于滑移層的組成成分主要為自由水，當(dāng)大量的自由水被消耗時(shí)，不能形成可以填充整個(gè)管道粗糙內(nèi)壁的滑移層，此時(shí)管道輸送阻力較大。當(dāng)水泥含量下降（砂灰比增大）時(shí)，水化反應(yīng)需水量下降，料漿流動(dòng)性較好，可以形成有效的滑移層，此時(shí)管道輸送阻力快速下降；當(dāng)水泥含量較少時(shí)，其水化反應(yīng)消耗的水對(duì)滑移層結(jié)構(gòu)的影響較小，管道輸送阻力降低速率減小。由圖4（b）可知：料漿在管道輸送過程中，灰砂比與阻力損失的關(guān)系曲線存在明顯轉(zhuǎn)折，當(dāng)灰砂比超過1∶8時(shí)，阻力損失的降低速率明顯變小，曲線趨于水平，灰砂比低于1∶10時(shí)會(huì)影響充填體強(qiáng)度的形成，因此本研究推薦的最佳灰砂比范圍為1∶8～1∶10。

5 工程實(shí)例驗(yàn)證

冀東地區(qū)某礦為地下礦山，由于生產(chǎn)需要建有1座地面充填站，標(biāo)高+110 m。該礦山采用全尾砂膠結(jié)充填技術(shù)進(jìn)行充填作業(yè)。充填選用的尾砂來自選廠濃縮池，經(jīng)濃縮后在攪拌桶中與水泥混合制備成濃度為70%～73%、灰砂比為1∶4～1∶6的充填料漿。充填料漿經(jīng)砂漿泵加壓后通過充填鉆孔（鉆孔長(zhǎng)度110 m）進(jìn)入地下采場(chǎng)。其中，充填管路內(nèi)徑180 mm，砂漿流量為100～140 m3/h，其中充填系統(tǒng)水平距離為2 100 m，豎直距離為500 m，坡面夾角為30°。該礦山充填系統(tǒng)如圖5所示。為準(zhǔn)確計(jì)算該礦山充填系統(tǒng)的輸送阻力，以充填鉆孔底端為0 m標(biāo)高，根據(jù)礦山實(shí)際充填管路，建立了模型尺寸和充填管路尺寸比為1∶1的幾何模型。

根據(jù)礦山充填料漿實(shí)際配比參數(shù)選擇灰砂比為1∶6，濃度分別為72%、73%的料漿進(jìn)行模擬試驗(yàn)。同時(shí)根據(jù)前述分析所得的參數(shù)最優(yōu)取值范圍，選擇灰砂比為1∶8，濃度為66%、68%、70%的料漿進(jìn)行對(duì)比驗(yàn)證，不同料漿所對(duì)應(yīng)的流變參數(shù)如表4所示。

選擇漿體初始速度為1.5 m/s進(jìn)行模擬試驗(yàn)，結(jié)果如圖6所示。由圖6可知：當(dāng)采用礦山實(shí)際充填參數(shù)時(shí)，管道輸送阻力較大。礦山充填管道的輸送阻力隨著料漿濃度增加而快速增大，當(dāng)料漿濃度小于70%時(shí)，可大大減小料漿管道輸送過程中的能量損耗，同時(shí)也會(huì)降低堵管、爆管的風(fēng)險(xiǎn)?？紤]到礦山充填結(jié)構(gòu)體的強(qiáng)度需求，為有效提高礦山生產(chǎn)效率和減小滲水量，充填料漿濃度不宜低于66%，因此，得到本試驗(yàn)條件下料漿的最佳輸送濃度為66%和68%。灰砂比由1∶6降低為1∶8時(shí)充填管道的輸送阻力減小，但是水泥含量過低時(shí)會(huì)影響礦山充填體早期強(qiáng)度形成。當(dāng)灰砂比大于1∶10時(shí)，水泥含量變化對(duì)充填體強(qiáng)度的影響較小，同時(shí)基于管道輸送阻力的變化情況得出該試驗(yàn)條件下的最佳灰砂比為1∶8。

6 結(jié) 論

（1）通過Comsol軟件建立了考慮壁面滑移效應(yīng)的管道輸送模型，將數(shù)值模擬結(jié)果與長(zhǎng)沙礦冶研究院的經(jīng)驗(yàn)公式和考慮滑移效應(yīng)的推導(dǎo)公式進(jìn)行比較，相對(duì)誤差分別為48.1%～386.1%和4.8%～16.7%，證明本研究所建立的數(shù)值模型用于計(jì)算考慮壁面滑移效應(yīng)的阻力損失具有可靠性。

（2）基于壁面滑移效應(yīng)分析管道阻力損失的變化特征，認(rèn)為阻力損失與灰砂比、濃度呈正相關(guān)關(guān)系，與管徑呈負(fù)相關(guān)關(guān)系。壁面滑移效應(yīng)是影響管道阻力損失的重要條件，研究成果為提高充填料漿管道輸送的穩(wěn)定性提供了理論依據(jù)。

（3）通過極差分析，發(fā)現(xiàn)各因素對(duì)料漿管道阻力損失的影響程度為管徑＞料漿濃度＞灰砂比。以冀東地區(qū)某充填礦山為研究背景進(jìn)行模擬試驗(yàn)，得到最佳充填參數(shù)為濃度66%、68%，灰砂比1∶8，對(duì)于提高礦山輸送管路的穩(wěn)定性有一定的參考意義。