基于CFD-DEM耦合模型的巖屑運(yùn)移數(shù)值模擬分析

胡金帥，張光偉，李峻嶺，陳雨，閆豐平

（西安石油大學(xué)機(jī)械工程學(xué)院，陜西 西安 710065）

0 引言

在大位移水平井鉆進(jìn)過程中，由于重力作用，巖屑顆粒極易在水平段形成巖屑床。過高的巖屑床會導(dǎo)致鉆桿高摩阻、高阻力矩，嚴(yán)重時會引發(fā)鉆井安全事故，因此巖屑顆粒運(yùn)移問題已經(jīng)成為制約井下巖屑運(yùn)移效率與井眼清潔程度的關(guān)鍵［1］。巖屑在環(huán)空流場中的運(yùn)移是個復(fù)雜的物理過程，影響巖屑運(yùn)移的因素很多。為了獲取影響巖屑顆粒運(yùn)移的相關(guān)規(guī)律，國內(nèi)外學(xué)者進(jìn)行了大量的研究。Nguyen等［2-3］將環(huán)空巖屑顆粒運(yùn)動形式采用3層巖屑運(yùn)移模型進(jìn)行研究，即純流體層、混合懸浮層、巖屑床層；汪海閣等［4-7］針對不同工況對巖屑運(yùn)移的多種影響因素進(jìn)行了研究；宋洵成等［8-9］通過分析巖屑受力，并考慮巖屑顆粒之間產(chǎn)生的相互作用力，分析了巖屑顆粒的運(yùn)移影響規(guī)律。目前，不同學(xué)者廣泛采用的Mixture模型忽略了顆粒與顆粒之間的碰撞和顆粒與流體之間的相對運(yùn)動，模擬結(jié)果誤差較大，不適合具有大顆粒的工況；而DPM模型能夠較為準(zhǔn)確地模擬顆粒在較低體積分?jǐn)?shù)時的運(yùn)動，但隨著顆粒體積分?jǐn)?shù)增大，顆粒運(yùn)動和顆粒間的作用變得越來越重要，模擬結(jié)果誤差也會增大，導(dǎo)致結(jié)果最終失真。巖屑的所有運(yùn)動都是巖屑顆粒受力后所表現(xiàn)的現(xiàn)象，因此利用CFD-DEM耦合方法［10-11］，流體運(yùn)動由Euler框架下的連續(xù)方法描述，在考慮巖屑顆粒之間碰撞及顆粒-流體作用的基礎(chǔ)上運(yùn)用DEM方法，可以建立一種更加準(zhǔn)確的數(shù)值模型來研究顆粒群在環(huán)空的流動。

鑒于此，本文基于CFD-DEM耦合模型，采用離散元求解器計算巖屑顆粒在環(huán)空流場中的受力情況和運(yùn)動軌跡，充分考慮各種接觸與碰撞，以準(zhǔn)確模擬和分析巖屑顆粒在井下環(huán)空固液兩相復(fù)雜流動狀態(tài)中的流動規(guī)律。通過對巖屑顆粒的數(shù)值模擬，研究不同工況對巖屑運(yùn)移規(guī)律的影響，得出提高巖屑運(yùn)移效率與井眼清潔程度的方法，為水平井巖屑顆粒運(yùn)移研究提供了理論設(shè)計指導(dǎo)。

1 巖屑運(yùn)移數(shù)學(xué)模型

1.1 流體運(yùn)動的控制方程

1.1.1 連續(xù)性方程

根據(jù)質(zhì)量守恒原理，可以得到井底流場控制方程的連續(xù)性方程［1］。在井眼環(huán)空流場中，單位時間內(nèi)通過井底流場單元體積的質(zhì)量,等于單元體積內(nèi)質(zhì)量的變化量，即：

式中：φ為固體顆粒集合體的孔隙度；ρl為鉆井液密度，kg/m3；t為時間，s；▽為拉普拉斯算子；ul為鉆井液運(yùn)動速度，m/s。

1.1.2 動量方程

根據(jù)牛頓第二定律，可以得到流場控制方程的動量方程。單位時間內(nèi)，單元體積內(nèi)動量的增量來自巖屑顆粒各個受力所引起的動量變化量，即：

式中：up為巖屑運(yùn)動速度，m/s；μ為流體動力黏性系數(shù)；p為流體壓力，Pa；Fp為周向壓力差，N；g為重力加速度，m/s2。

1.2 巖屑顆粒運(yùn)動模型

1.2.1 巖屑顆粒受力分析

鉆進(jìn)時，巖屑運(yùn)動軌跡由其作用力平衡方程求得。單個巖屑顆粒遵循牛頓第二定律，且只考慮巖屑有效重力、拖拽力、壓力梯度力、旋轉(zhuǎn)科氏力、附加質(zhì)量力以及旋轉(zhuǎn)離心力［7］。

1.2.2 巖屑顆粒碰撞模型

在井眼環(huán)空中，巖屑隨著旋轉(zhuǎn)流體做圓周運(yùn)動，由于其大小不一，在運(yùn)動過程中必然發(fā)生接觸碰撞，故引入DEM控制方程。巖屑無滑移碰撞模型采用Hertz-Mindlin非線性接觸模型［12］。

巖屑顆粒受到的抗轉(zhuǎn)矩Nij（下標(biāo)i，j代表不同的顆粒）可分為彈簧轉(zhuǎn)矩和黏性阻尼轉(zhuǎn)矩，表示為

式中：kr為滾動剛度，N/m；θp為巖屑雙顆粒的相對滾動角速度，rad/s；Δt為時間增量，s。

其中：

式中：Rr為巖屑雙顆粒滾動有效半徑，m；Fn為巖屑顆粒法向應(yīng)力，Pa；μr為巖屑滾動阻力系數(shù)；ri，rj分別為2個顆粒的半徑，m；為彈簧轉(zhuǎn)矩極限，N·m。

1.3 CFD-DEM耦合模型

CFD-DEM耦合模型計算的關(guān)鍵在于鉆井液與巖屑顆粒固液多物理場的數(shù)據(jù)交換，巖屑顆粒在空間中的動態(tài)位置變化引起巖屑顆粒的體積分?jǐn)?shù)與動量力源交換的變化。根據(jù)巖屑顆粒體積分?jǐn)?shù)和流體動量力源的變化，本文進(jìn)行相應(yīng)的CFD-DEM耦合收斂分析。CFD數(shù)值求解采用SIMPLE算法，對鉆井液的流體運(yùn)動控制方程進(jìn)行求解，數(shù)值求解迭代至收斂后，得到巖屑顆粒與鉆井液的相關(guān)速度，以及巖屑顆粒在鉆井液中所受的力。將求解數(shù)據(jù)傳輸?shù)紼DEM中，進(jìn)行離散元求解，并將巖屑顆粒運(yùn)動方程代入，EDEM迭代計算完成后，得到離散相巖屑顆粒的空間位置、速度、受力等數(shù)據(jù)。將數(shù)據(jù)再次傳遞回CFD求解器，并更新計算域信息進(jìn)入下一個循環(huán)時間。

2 物理模型

2.1 井眼環(huán)空幾何模型構(gòu)建

根據(jù)相關(guān)鉆井設(shè)計資料，選取環(huán)空管道的長度為6 m，井眼和鉆桿的直徑分別為130，76 mm。使用workbench完成計算流場域構(gòu)建和網(wǎng)格劃分，為保證仿真模擬的精準(zhǔn)性，網(wǎng)格劃分采用六面體結(jié)構(gòu)化網(wǎng)格。在EDEM的顆粒工廠中建立巖屑顆粒模型，巖屑顆粒呈球狀。井下環(huán)空三維建模及網(wǎng)格劃分示意如圖1所示，網(wǎng)格數(shù)量為976 960個。

圖1 井下環(huán)空三維建模及網(wǎng)格劃分示意

2.2 模型參數(shù)與邊界條件設(shè)定

計算模型為非定常計算，總體時間5 s，流體計算時間步長設(shè)置為0.000 1 s。CFD-DEM耦合模型中顆粒計算時間步長要求不能大于流體計算時間步長，參考瑞利時間步長，顆粒計算時間步長控制在瑞利時間步長的60%以內(nèi)。湍流流場的計算采用二階迎風(fēng)格式，湍流模型為RNG k-ε模型，為了減少計算，不考慮井筒溫度等因素［13］。重力加速度為9.81 m/s2，重力方向?yàn)榇怪毕蛳?。鉆井液的入口定義為速度入口邊界條件，出口定義為自由出口邊界；設(shè)置鉆井液和巖屑顆粒從入口進(jìn)入，鉆井液設(shè)定為水，密度為1 000 kg/m3，黏度為0.9 mPa·s，選擇能夠更好地描述水基鉆井液流變模式的赫謝爾-巴爾克萊模式；巖屑粒徑為1.5 mm，密度為2 600 kg/m3，泊松比為0.26，彈性模量為15 000 MPa，顆粒工廠的顆粒生成速度為5 000個/s。井壁壁面采用無滑移邊界條件，鉆桿設(shè)定為剛體旋轉(zhuǎn)壁面，壁面附近采用標(biāo)準(zhǔn)壁函數(shù)法。

2.3 模型數(shù)值模擬對比驗(yàn)證

為了驗(yàn)證本文計算模型的準(zhǔn)確性，以文獻(xiàn)［14］中的實(shí)驗(yàn)為例，采用本文CFD-DEM耦合模型計算巖屑顆粒運(yùn)移速度的變化，并與實(shí)驗(yàn)結(jié)果進(jìn)行對比。結(jié)果表明，采用CFD-DEM耦合模型計算的巖屑顆粒運(yùn)移速度與實(shí)驗(yàn)結(jié)果基本一致。2種方法的巖屑顆粒運(yùn)移速度誤差不超過5%，從而證明本文所建立的CFD-DEM耦合模型的可行性。

3 模擬分析

3.1 不同鉆桿偏心度下巖屑顆粒運(yùn)移規(guī)律

在井下，鉆桿中心軸線會與井眼的中心軸線發(fā)生位移偏離，產(chǎn)生不同的偏心度。環(huán)空流道分布不對稱，將使得鉆井液環(huán)空流速及巖屑顆粒運(yùn)移速度發(fā)生改變。保持鉆桿轉(zhuǎn)速為8 rad/s，鉆井液入口流速為3.0 m/s，巖屑粒徑為1.5 mm，模擬得出距鉆井液入口的距離（Z）為3 m處不同鉆桿偏心度下鉆井液環(huán)空流速的分布（見圖2）。不同鉆桿偏心度下巖屑顆粒運(yùn)移速度見圖3。

圖2 不同鉆桿偏心度下鉆井液環(huán)空流速的分布

鉆桿偏心度會影響鉆井液在環(huán)空流道中的流速。由圖2可知：鉆井液在環(huán)空尺寸較大處的流速大，在環(huán)空尺寸較小處的流速??；在偏心度為0時，流道中的鉆井液流速在軸向上分布均勻，流速均值為2 m/s；在偏心度為0.75時，環(huán)空尺寸較小處的流速降至0.17 m/s。由圖3可知，巖屑顆粒在通過環(huán)空流道時，隨著偏心度的增大，巖屑顆粒的運(yùn)移能力不斷削弱。當(dāng)偏心度達(dá)到0.75時，環(huán)空尺寸較大處巖屑顆粒運(yùn)移速度分布均勻，環(huán)空尺寸較小處巖屑顆粒運(yùn)移速度較小，最終巖屑顆粒平均運(yùn)移速度降至0.50 m/s。流速的降低導(dǎo)致巖屑顆粒的運(yùn)移受到抑制，會產(chǎn)生巖屑床的堆積，在軸向上的巖屑顆粒運(yùn)移距離逐漸減小。同時，在偏心度的增加過程中，鉆桿會逐漸接近環(huán)空流道中的巖屑床，從而帶動巖屑顆粒產(chǎn)生相應(yīng)的切向速度，使得巖屑床厚度減小，巖屑顆粒進(jìn)入懸浮層，提升了巖屑顆粒運(yùn)移效率。由于偏心度的增加，鉆桿對巖屑床的攪拌作用更加明顯，從而降低了井筒內(nèi)巖屑床的厚度［15］，但最終巖屑床的堆積速度會超過鉆桿帶動巖屑的速度。位于環(huán)空尺寸較小處的巖屑顆粒大量堆積，難以進(jìn)入懸浮層。即使增加鉆井液環(huán)空流速，巖屑床堆積的巖屑顆粒也很難運(yùn)移。

圖3 不同鉆桿偏心度下巖屑顆粒運(yùn)移速度

3.2 不同鉆桿轉(zhuǎn)速下巖屑顆粒運(yùn)移規(guī)律

在鉆井過程中，鉆桿轉(zhuǎn)速也會影響井下環(huán)空流道中的巖屑運(yùn)移。鉆桿的旋轉(zhuǎn)會帶動鉆井液產(chǎn)生漩渦流，在入口處可以觀察到巖屑顆粒產(chǎn)生漩渦狀態(tài)的分布。鉆桿和巖屑之間的機(jī)械刮削和攪動作用，使得鉆桿旋轉(zhuǎn)對巖屑床造成了破壞［16］。為探究巖屑顆粒在不同鉆桿轉(zhuǎn)速下的運(yùn)移規(guī)律，保持鉆桿偏心度為0.75，鉆井液入口流速為3.0 m/s，巖屑粒徑為1.5 mm，模擬得出Z=3 m處不同鉆桿轉(zhuǎn)速下鉆井液環(huán)空流速的分布（見圖4）。不同鉆桿轉(zhuǎn)速下巖屑顆粒運(yùn)移速度見圖5。

鉆桿的旋轉(zhuǎn)帶動鉆井液產(chǎn)生環(huán)空旋轉(zhuǎn)運(yùn)動，環(huán)空旋轉(zhuǎn)運(yùn)動產(chǎn)生的切向加速度隨著環(huán)空徑向距離的增大而減小，在靠近鉆桿處最大，靠近井壁處最小。從圖4、圖5可知：在鉆桿轉(zhuǎn)速為6 rad/s時，位于懸浮層中心部位的鉆井液環(huán)空流速最高可達(dá)2.4 m/s，巖屑顆粒運(yùn)移速度在t=5 s時僅為0.50 m/s；隨著鉆桿轉(zhuǎn)速逐漸增至12 rad/s，巖屑顆粒運(yùn)移速度在t=5 s時可達(dá)到0.62 m/s。從圖4看出：鉆桿的旋轉(zhuǎn)帶動巖屑床的巖屑顆粒獲得相應(yīng)的舉升力，使得巖屑顆粒進(jìn)入懸浮層；隨著鉆桿轉(zhuǎn)速的增大，環(huán)空尺寸較小處巖屑顆粒的沉降會大大減弱。因此，鉆桿旋轉(zhuǎn)有利于巖屑顆粒的運(yùn)移，隨著鉆桿轉(zhuǎn)速的增大，鉆井液與巖屑顆粒的運(yùn)移速度增大，也有利于減弱巖屑床的堆積效應(yīng)。

圖4 不同鉆桿轉(zhuǎn)速下鉆井液環(huán)空流速的分布

圖5 不同鉆桿轉(zhuǎn)速下巖屑顆粒運(yùn)移速度

3.3 鉆井液不同入口流速下巖屑顆粒運(yùn)移規(guī)律

在實(shí)際鉆井過程中，鉆井液環(huán)空流速的大小通過排量進(jìn)行調(diào)節(jié)，根據(jù)鉆井液最小排量來確定攜帶巖屑顆粒運(yùn)移所需環(huán)空流速的最低值。環(huán)空流速較低時，巖屑在重力作用下會在環(huán)空下側(cè)堆積，形成巖屑床［17］，因此，需要進(jìn)一步研究鉆井液入口流速對巖屑運(yùn)移規(guī)律的影響。保持偏心度為0.75，鉆桿轉(zhuǎn)速為8 rad/s，巖屑粒徑為1.5 mm，模擬得出Z=3 m處鉆井液不同入口流速下的環(huán)空流速分布（見圖6）。鉆井液不同入口流速下巖屑顆粒運(yùn)移速度見圖7。

圖6 鉆井液不同入口流速下的環(huán)空流速分布

圖7 鉆井液不同入口流速下巖屑顆粒運(yùn)移速度

從圖6可以看出：鉆井液入口流速的升高將使得環(huán)空流速增大，產(chǎn)生水平向推力帶動巖屑顆粒運(yùn)動，使得處于懸浮層的巖屑向前運(yùn)移的距離增大；落入巖屑床表面的巖屑也在入口流速增大產(chǎn)生的水平向推力與鉆桿旋轉(zhuǎn)帶來的舉升力的共同影響下進(jìn)入懸浮層，懸浮層的顆粒增多，加大了顆粒間的碰撞，為巖屑的運(yùn)移提供了能量，使得巖屑床厚度降低。結(jié)合圖7可知：當(dāng)鉆井液入口流速為1.5 m/s，在t=5 s時巖屑顆粒運(yùn)移速度僅為0.31 m/s；隨著鉆井液入口流速增大至2.5 m/s時，在t=5 s時巖屑顆粒運(yùn)移速度增大至1.18 m/s，相應(yīng)的巖屑水平向運(yùn)移距離也會增大。因此，鉆井液入口流速的增加，有利于環(huán)空巖屑顆粒的運(yùn)移，減小巖屑床厚度，提高巖屑顆粒的運(yùn)移效率。

3.4 不同巖屑粒徑下巖屑顆粒運(yùn)移規(guī)律

對密度相同的巖屑，由于小粒徑顆粒比表面積大，自身重力不超過其受到的黏滯阻力，能夠穩(wěn)定懸浮于鉆井液中；大粒徑顆粒生成后會迅速沉降至井壁，形成穩(wěn)定的巖屑床。在環(huán)空流道中，巖屑粒徑的大小也會影響巖屑運(yùn)移：粒徑小的顆粒更易沖刷到懸浮層，顆粒之間的碰撞不會影響巖屑運(yùn)移；粒徑較大的顆粒由于重力的作用難以升至懸浮層，顆粒之間或顆粒與壁面的碰撞會影響顆粒的高效運(yùn)移。保持偏心度為0.75，鉆桿轉(zhuǎn)速為8 rad/s，鉆井液入口流速為3.0 m/s，模擬得出Z=3 m處不同巖屑粒徑下鉆井液環(huán)空流速的分布（見圖8）。不同巖屑粒徑下巖屑顆粒運(yùn)移速度見圖9。

圖8 不同巖屑粒徑下鉆井液環(huán)空流速分布

由圖8可知，巖屑粒徑對鉆井液環(huán)空流速的改變影響不大，對巖屑的沉降效果影響較為明顯。當(dāng)巖屑粒徑從0.5 mm增大至2.5 mm時，巖屑沉降效果比較顯著，重力方向下最底部速度為0的區(qū)域逐漸變大，表明巖屑粒徑的增大會導(dǎo)致巖屑顆粒的堆積增加。由圖9可知：當(dāng)巖屑粒徑為0.5 mm，在t=5 s時巖屑顆粒運(yùn)移速度為1.02 m/s；而當(dāng)巖屑粒徑增大至2.5 mm，巖屑顆粒運(yùn)移速度在t=5 s時僅為0.38 m/s。這表明在運(yùn)移過程中，巖屑粒徑增大會導(dǎo)致巖屑顆粒運(yùn)移速度的降低，粒徑較大會明顯影響巖屑顆粒的運(yùn)移效果。巖屑粒徑的增大，導(dǎo)致巖屑顆粒自身的重力增大，鉆桿的旋轉(zhuǎn)產(chǎn)生的徑向舉升力不足，不能帶動大顆粒巖屑運(yùn)移至懸浮層，進(jìn)而導(dǎo)致巖屑床的形成。在較小的環(huán)空中，巖屑的堆積使得鉆井液增速，加快了顆粒之間的碰撞，導(dǎo)致更多的巖屑堆積，從而產(chǎn)生環(huán)空堵塞，甚至導(dǎo)致井下復(fù)雜情況的發(fā)生。

圖9 不同巖屑粒徑下巖屑顆粒運(yùn)移速度

4 結(jié)論

1）通過CFD-DEM的仿真與實(shí)驗(yàn)對比，驗(yàn)證了該耦合模型的可行性。相對于普通的CFD方法，該模型更能夠模擬實(shí)際工況，為后續(xù)的理論研究和實(shí)驗(yàn)研究提供參考，但對于非球形狀的巖屑顆粒缺乏準(zhǔn)確的模擬仿真分析，仍需進(jìn)一步的研究。

2）鉆桿偏心度的增大將導(dǎo)致懸浮層巖屑顆粒平均運(yùn)移速度降低，其中偏心度大于0.50時，巖屑顆粒因自身重力的作用更易形成巖屑床。巖屑粒徑的增大將使得沉降速度增大，從而導(dǎo)致巖屑顆粒運(yùn)移速度降低，巖屑粒徑超過1.5 mm時，巖屑顆粒運(yùn)移速度降低較為明顯。

3）提高鉆桿轉(zhuǎn)速，對巖屑顆粒運(yùn)移速度的提升效果較小，但會使得巖屑滯留量減少。鉆桿旋轉(zhuǎn)的反復(fù)碾磨也可以減小巖屑粒徑，有利于破壞巖屑床，從而增大巖屑顆粒升至懸浮層的概率。增大鉆井液入口流速，能夠使得巖屑顆粒運(yùn)移速度同比例增大。為了減少巖屑床的形成，增大井眼清潔程度，應(yīng)增大鉆桿轉(zhuǎn)速與鉆井液入口流速，減小巖屑粒徑，從而提高巖屑顆粒運(yùn)移效率。