李紫璇 張菲菲 祝鈺明 吳墨染 劉泳敬 于 琛

(1.長江大學(xué)石油工程學(xué)院 2.油氣鉆井技術(shù)國家工程實驗室防漏堵漏研究室 3.油氣鉆采工程湖北省重點實驗室4.青海油田分公司 5.渤海鉆探第一鉆井分公司 6.渤海鉆探第三鉆井分公司 7.渤海鉆探工程技術(shù)研究院)

0 引言

統(tǒng)計數(shù)據(jù)表明[1]，鉆進過程中大約有70%的時間損失與鉆具阻卡有關(guān)。如不能及時發(fā)現(xiàn)，一旦發(fā)生卡鉆事故，往往需要停鉆處理，將會花費大量的時間及處理成本，并有可能引起更為嚴重的次生危害[2]。目前解決此類問題最有效的方法是對鉆井過程中的卡鉆風(fēng)險進行準確預(yù)測分析，并對風(fēng)險的處理提供合理的決策方案以規(guī)避卡鉆，縮短非作業(yè)時間[3]。

目前，卡鉆預(yù)警方法多以基于大數(shù)據(jù)的統(tǒng)計分析預(yù)警方法和基于鉆井模型的分析模型為主。前者主要分為判別分析法和模式識別法。該方法以收集大量有效數(shù)據(jù)為前提，通過對卡鉆發(fā)生原因、卡鉆時的作業(yè)狀況等進行分析，利用神經(jīng)網(wǎng)絡(luò)、多元統(tǒng)計分析、模糊邏輯及層次分析法等方法建立模型或?qū)＜蚁到y(tǒng)判斷是否發(fā)生卡鉆[4－7]。要想實現(xiàn)實時卡鉆預(yù)測，還需建立時間序列ARMA 模型[8]。用于統(tǒng)計分析方法的數(shù)據(jù)多來自于地層信息和鄰井的卡鉆事故，若缺乏鄰井資料，則難以建立有效的卡鉆檢測系統(tǒng)和方法。同時，收集到的數(shù)據(jù)往往受地域限制，使得利用統(tǒng)計分析的卡鉆預(yù)警模型難以進行大范圍推廣[9]?？ㄣ@預(yù)警方法使用的物理模型多為摩阻扭矩模型，使用摩阻扭矩模型對摩阻系數(shù)進行反算，將摩阻系數(shù)隨深度的異常變化作為卡鉆風(fēng)險監(jiān)測的指標[10]；分析對比扭矩預(yù)測值和扭矩實測值，對鉆井過程中的異常進行監(jiān)測[11－12]。卡鉆發(fā)生時會導(dǎo)致鉆井參數(shù)出現(xiàn)異常變化，但傳統(tǒng)物理模型所需計算量大且難以完全反映卡鉆發(fā)生時的變化規(guī)律[13]，在實際應(yīng)用中誤報率高、準確率低、實時性較差，且難以對發(fā)生的卡鉆進行分類[14]。

本文將瞬態(tài)巖屑運移模型、摩阻扭矩模型和機器學(xué)習(xí)相結(jié)合，研究基于實時錄井?dāng)?shù)據(jù)的卡鉆風(fēng)險監(jiān)測及預(yù)警技術(shù)，利用瞬態(tài)巖屑運移模型對巖屑床分布和井眼清潔狀況進行實時監(jiān)測；將瞬態(tài)巖屑運移模型的計算結(jié)果作為輸入，建立考慮井筒巖屑分布對鉆桿受力影響的改進摩阻扭矩模型；利用貝葉斯優(yōu)化訓(xùn)練鉆井模型，使模型能夠不斷進行自我調(diào)整以適應(yīng)新的工況，即使改變井況或井型，也能不斷提高模型精度，且模型訓(xùn)練速度快，僅需當(dāng)前井的錄井?dāng)?shù)據(jù)便能完成模型訓(xùn)練。通過分析對比混合模型實時預(yù)測結(jié)果和實測數(shù)據(jù)，并進行時序數(shù)據(jù)分析，最終得到卡鉆風(fēng)險指數(shù)。

1 鉆井模型

1.1 實時井眼清潔監(jiān)測模型

井眼清潔不充分是造成卡鉆的常見原因之一，約占比例的50%[3]。在鉆進過程中對整個井眼的巖屑分布情況進行準確且實時地監(jiān)測是有效預(yù)測卡鉆事故的關(guān)鍵。

采用瞬態(tài)巖屑運移模型對巖屑沿井眼軌跡分布進行一維數(shù)值模擬[15]，將巖屑運移過程在時間域上進行離散化，通過實時計算每個離散單元的巖屑堆積量和環(huán)流壓力梯度實現(xiàn)井眼清潔實時監(jiān)測。模型公式基于質(zhì)量和動量守恒方程，公式為:

式中:角標α和β分別為流型中的不同部分(也可稱為層)；n為給定流型的總層數(shù)；A為層的橫截面積，mm2；C為巖屑密度，g/cm3；U 為流速，m/s；S為層的周長，mm；p 為流動壓力，Pa；ρ為層的平均密度，g/cm3；τ 為切應(yīng)力，Pa；Ex為不同層之間的巖屑體積交換率；g 為重力加速度，m/s2。

在數(shù)值模擬中使用基于有限體積的半隱式差分方法來求解方程。在鉆井過程中，鉆井參數(shù)不斷變化，巖屑的流動形態(tài)也會隨之改變，為使模型能同時處理時間域和空間域中巖屑流型的變化，使用基于瞬態(tài)力學(xué)模型的廣義流型[16]，在每次計算前對離散單元的巖屑流型進行判斷，根據(jù)巖屑流型建立相應(yīng)的模型方程組。所有流型的每個層都有自己的質(zhì)量守恒和動量守恒方程，給定流型的方程總數(shù)將取決于其總層數(shù)。模型計算過程中，首先根據(jù)初始狀態(tài)通過動量守恒方程求解下一個時間步的速度和壓力，將收斂的速度和壓力帶入質(zhì)量守恒方程求得新的懸浮巖屑密度和沉積巖屑體積，最后把結(jié)果帶入下一個時間步，依次進行計算。模型建立和求解的詳細過程見文獻[17]。

1.2 實時摩阻扭矩計算模型

隨著井筒中巖屑床的堆積，井筒內(nèi)井眼清潔狀況變差，鉆桿扭矩和起下鉆時的大鉤載荷增加明顯，卡鉆風(fēng)險也隨之增加[18－19]。若要模擬巖屑床對鉆桿運動的精確影響，需要對整個井眼中鉆桿在環(huán)空中的位置進行復(fù)雜的有限元模擬。為便于實時計算，在實際應(yīng)用中，通過使用一個附加力來表示巖屑床和鉆桿之間的復(fù)雜相互作用:起下鉆時，附加力施加在鉆桿軸向上，表示巖屑床對鉆桿軸向作用力的影響；鉆桿旋轉(zhuǎn)時，附加力與鉆桿壁曲面相切，方向與旋轉(zhuǎn)方向相反。

為便于描述空間井眼軌跡，將空間直角坐標系與弗萊納坐標系相結(jié)合，如圖1 所示。通過對鉆桿單元體進行受力分析，建立旋轉(zhuǎn)鉆桿的力和力矩平衡方程，方程如下[20]:

圖1 鉆桿單元體在3D 井眼中的受力圖Fig.1 Force diagram of drill pipe unit body in 3D wellbore

式中:Ft為軸向力，N；wp為單位鉆桿在鉆井液中重力，N；wc為接觸力，N；tz、nz和bz分別為z方向上單位切向、法向和副法向向量；θ為沿平面法線與接觸力間的夾角，(°) ； Fc為巖屑床施加的附加力，N；φ為井斜角，(°) ；μ為摩阻系數(shù)，無量綱；rp為鉆桿半徑，m；Mt為鉆桿旋轉(zhuǎn)所需的軸向扭矩，N·m；ds為鉆桿單元的長度，m；κ為曲率，m－1。

通過上述方程，可得Fc等于鉆桿上的附加力。

除了每個單元體的附加摩擦力需要基于局部巖屑床的高度來計算之外，用于解改進摩阻扭矩模型的算法與解標準摩阻扭矩模型的算法相似。

實時井眼清潔監(jiān)測模型與摩阻扭矩模型結(jié)合應(yīng)用時，以實時操作參數(shù)作為實時井眼清潔監(jiān)測模型的輸入，預(yù)測巖屑床在井眼中的實時動態(tài)分布。摩阻扭矩模型將巖屑床的分布預(yù)測結(jié)果作為輸入，計算出每個離散點的附加力，最后得出鉆桿上的摩阻扭矩值。為進一步研究巖屑床對鉆桿的影響，該研究利用試驗設(shè)備[21]，通過一系列試驗數(shù)據(jù)，探究附加力與巖屑床高的關(guān)系，試驗參數(shù)如表1 所示。

表1 試驗測試參數(shù)Table 1 Experimental test matrix

為了便于研究相關(guān)性，分別對巖屑床高和巖屑床高對扭矩的影響進行標準化處理，如下所示:

式中:hc為巖屑床高度，mm；D為井眼直徑，mm；hs為標準化巖屑床高度，無量綱；Ts為標準化巖屑床高hs下的扭矩測量值，kN·m；Tc為相同井眼結(jié)構(gòu)下清潔井眼的扭矩計算值，kN·m；TR為由巖屑引起的扭矩變化率，無量綱。

經(jīng)標準化處理后，測試數(shù)據(jù)無量綱巖屑床高和扭矩變化率的關(guān)系如圖2 所示。

由圖2 可知，標準化巖屑床高和扭矩變化率的關(guān)系大致呈指數(shù)形式。因此，在數(shù)據(jù)集的基礎(chǔ)上發(fā)展指數(shù)函數(shù)相關(guān)性，公式如下:

式中:a和b均為無量綱參數(shù)，圖2 所示曲線中a和b分別為10.2 和2.025。

圖2 標準化巖屑床高和扭矩變化率的關(guān)系圖Fig.2 Relationship between standardized cuttings bed height and torque change rate

由于試驗是在模擬直井的設(shè)施上進行的，扭矩與井眼和鉆桿間的接觸力、井眼和鉆桿間的摩阻系數(shù)成正比。因接觸力在試驗過程中保持不變，故扭矩的變化只與摩阻系數(shù)的變化有關(guān)。又因參數(shù)是無量綱的，故公式(10) 可轉(zhuǎn)換為表示標準巖屑床高與摩阻系數(shù)變化率之間的關(guān)系，如公式(11)所示。

式中:fR為摩阻系數(shù)的變化率，可以用公式(12)定義。

式中:fs為標準化巖屑床高下的摩阻系數(shù)；fc為沒有巖屑時的摩阻系數(shù)。

將公式(11) 應(yīng)用于前文所述的摩阻扭矩計算模型中，用于計算有巖屑井的扭矩值。與試驗測試不同，在實際鉆井應(yīng)用中需要對參數(shù)進行改變。雖然由試驗數(shù)據(jù)得出的相關(guān)性可能不適用于所有情況，但是摩阻系數(shù)和無量綱巖屑床高的變化趨勢應(yīng)是相同的(即a和b的值可能在不同的情況下有所不同，但公式(11) 仍然有效)。

2 機器學(xué)習(xí)和模型訓(xùn)練

在摩阻扭矩模型計算過程中，鉆桿被分成若干離散單元，由于巖屑床高沿鉆桿發(fā)生變化，所以每個單元體都需要一個特定的摩阻系數(shù)。前文所述摩阻系數(shù)與巖屑床高呈指數(shù)關(guān)系，但a和b的值可能在不同的工況下有所不同，為使模型應(yīng)用范圍更廣，更好地適用于現(xiàn)場實際，需要以歷史鉆井?dāng)?shù)據(jù)為輸入，使用參數(shù)估值的方法對模型進行訓(xùn)練。

通過現(xiàn)場數(shù)據(jù)訓(xùn)練改進摩阻扭矩模型的過程，本質(zhì)上是對考慮井眼清潔影響的摩阻扭矩模型進行參數(shù)估值，可以視為尋找最優(yōu)參數(shù)集(即摩阻系數(shù)沿井眼的分布情況)，使得井口扭矩與轉(zhuǎn)盤扭矩差值最小，是一個非線性優(yōu)化問題。貝葉斯優(yōu)化算法可用較少的評估次數(shù)求得復(fù)雜目標函數(shù)的最優(yōu)解，是較為常用的優(yōu)化算法，其關(guān)鍵在于概率代理模型和采集函數(shù)的建立。

訓(xùn)練模型優(yōu)化過程需要最小化扭矩差值，選用納什效率系數(shù)(NSE) 作為模型的評價指標，表示模型擬合方差占總方差的百分比[22]。因極值對NSE影響較大，本文選用lgNSE作為模型評價指標，可有效降低極值對模型總體擬合效果的影響，如公式(13) 所示[23]。式中:yt和分別代表時間步為t時的扭矩實測值和計算值，代表實測值的平均值。

計算結(jié)果越接近于1，參數(shù)優(yōu)化效果越好。貝葉斯優(yōu)化的目的是使lgNSE最大化，以找出最有可能的摩阻系數(shù)分布情況。

使用概率代理模型替代評估代價高昂的復(fù)雜目標函數(shù)，根據(jù)貝葉斯定理不斷地進行迭代以增加信息量，對先驗知識進行修正[24]。高斯過程函數(shù)具有高度的靈活性、可擴展性和可分析性，對線性和非線性關(guān)系均具有良好的替代性[25]，因而成為貝葉斯優(yōu)化中應(yīng)用最廣泛的概率代理模型。高斯過程是多元高斯概率分布的范化，由均值函數(shù)m(x) 和半正定的協(xié)方差函數(shù)k(x，x′) 構(gòu)成，如公式(14)所示。

式中:g(x) 為多元高斯概率分布，GP為高斯過程。

在使用過程中，為方便計算，假定先驗均值函數(shù)為0，該設(shè)定對后驗分布的準確性幾乎沒有影響[26]，協(xié)方差函數(shù)采用的是一類高度靈活的Matérn 協(xié)方差函數(shù)簇，可產(chǎn)生二階可微的樣本函數(shù)，如公式(15) 所示。

式中:v為平滑參數(shù)，l為尺度參數(shù)，Kv為第二類變形貝塞爾函數(shù)。

采集函數(shù)根據(jù)概率代理模型的后驗結(jié)果構(gòu)造主動選擇策略，被用于確定貝葉斯優(yōu)化的方向，選擇采集函數(shù)最大值點作為下一個評估點。采集函數(shù)有多種類型，如提升概率(PI)、期望提升量(EI)、置信上邊界(UCB) 等，其中EI能將提升概率和提升量整合起來，所需參數(shù)少且能兼顧全局搜索和局部變化，如公式(16) 所示。

式中:y為當(dāng)前最優(yōu)函數(shù)值，φ為標準正態(tài)分布累計密度函數(shù)，μt(x) 為均值，σt(x) 為方差，ξ為平衡參數(shù)，無因次。

本文使用貝葉斯優(yōu)化對改進摩阻扭矩模型進行訓(xùn)練的流程圖如圖3 所示。圖3 中錄井?dāng)?shù)據(jù)為實時錄井?dāng)?shù)據(jù)，井眼數(shù)據(jù)為井身結(jié)構(gòu)、井眼信息、鉆具組合和鉆井液參數(shù)等靜態(tài)數(shù)據(jù)。

圖3 使用貝葉斯優(yōu)化訓(xùn)練摩阻扭矩模型流程圖Fig.3 Flow chart of torque and drag model training by Bayesian optimization

3 實時卡鉆風(fēng)險系數(shù)分析

本文所用實時卡鉆風(fēng)險預(yù)測方法通過分析3 個不同來源的扭矩數(shù)據(jù):考慮巖屑影響的扭矩預(yù)測值、不考慮巖屑影響的扭矩預(yù)測值以及離底自由旋轉(zhuǎn)扭矩的測量值；生成兩個獨立的卡鉆風(fēng)險指數(shù)，即與井眼清潔相關(guān)的卡鉆風(fēng)險指數(shù)和與其他因素相關(guān)的卡鉆風(fēng)險指數(shù)。由于摩阻扭矩模型已通過數(shù)據(jù)驅(qū)動模型以歷史數(shù)據(jù)為基礎(chǔ)做好了調(diào)試，所以實測值與預(yù)測值出現(xiàn)偏差表示井下異常情況的出現(xiàn)。本文選用時序數(shù)據(jù)分析的方法對實時數(shù)據(jù)的異常趨勢進行檢測。首先計算相對偏差值，計算公式如下:

式中:Tc為考慮巖屑影響的扭矩預(yù)測值，kN·m；Tn為不考慮巖屑影響的扭矩預(yù)測值，kN·m；Tr為離底自由旋轉(zhuǎn)扭矩的測量值，kN·m；Rc為考慮巖屑影響與不考慮巖屑影響的扭矩預(yù)測值之間的相對偏差，無量綱；Rm為實測值與考慮巖屑影響的扭矩預(yù)測值之間的相對偏差，無量綱。

在趨勢分析中，移動平均值偏差方法被用于追蹤兩個相對偏差值隨時間的變化趨勢，移動平均偏差值為正表示數(shù)據(jù)處于增長趨勢，為負表示數(shù)據(jù)處于降低趨勢[15]，如公式(19)所示。

最后，計算出表示卡鉆事故發(fā)生概率的風(fēng)險指數(shù)，為實時監(jiān)測提供更直觀的監(jiān)測結(jié)果和警報。最終卡鉆風(fēng)險指數(shù)將被轉(zhuǎn)換成介于0 到1 之間的值，0 表示沒有卡鉆風(fēng)險，1 表示鉆桿有極大可能會卡死，公式為:

4 案例分析

選用鉆井實時數(shù)據(jù)對模型進行驗證。圖4 為一段時長為290 min 的卡鉆發(fā)生前的實際井錄井?dāng)?shù)據(jù)，該數(shù)據(jù)每10 s 記錄一次?？ㄣ@發(fā)生在鉆進過程中，井斜角約為56°。如圖4 所示，在早上5:20 左右，當(dāng)其他參數(shù)基本保持不變時，扭矩從15 kN·m 增加到18 kN·m，且增長速度較快。

圖4 卡鉆發(fā)生之前的實時鉆井?dāng)?shù)據(jù)Fig.4 Real-time drilling data before pipe sticking

為更好地解釋扭矩的變化，將扭矩與井深相對應(yīng)，如圖5 所示。從2 255 m 鉆到2 450 m 期間扭矩測量值持續(xù)增加?？紤]巖屑影響的扭矩預(yù)測值與不考慮巖屑影響的扭矩預(yù)測值之間的偏差持續(xù)增加，表明井眼中的巖屑堆積在阻礙鉆進。在鉆到大約2 410 m (早上5:20) 時，扭矩從15 kN·m 迅速增加到18 kN·m，井下出現(xiàn)異常狀況。由于鉆井隊沒有注意到扭矩的迅速增加，繼續(xù)向下鉆進，在鉆桿卡死前又持續(xù)鉆了40 m。在后續(xù)鉆的這40 m 中，扭矩依然在增加。

圖5 扭矩隨深度變化圖Fig.5 Change of torque with depth

將實時卡鉆風(fēng)險分析應(yīng)用于鉆后分析和算法驗證?？ㄣ@事故發(fā)生前5 h 的風(fēng)險指數(shù)如圖6 所示。

圖6 實時卡鉆風(fēng)險分析結(jié)果Fig.6 Analysis results of real-time pipe sticking risks

由圖6 可知，與井眼清潔相關(guān)的卡鉆風(fēng)險指數(shù)在2:05 左右已經(jīng)達到0.6，并在事發(fā)前上升到0.7 以上。5:50 前，與其他因素相關(guān)的卡鉆風(fēng)險指數(shù)處于極低水平(0.2～0.4 之間)，表明目前引起卡鉆的主要因素與井眼清潔相關(guān)。當(dāng)扭矩突然增大，與其他因素相關(guān)的卡鉆風(fēng)險系數(shù)突然增加到0.8 以上，表明發(fā)生卡鉆事故的可能性很大。通過實例分析，該模型能夠捕捉到卡鉆是否發(fā)生，同時能夠?qū)σl(fā)生的卡鉆進行分類，并在事故發(fā)生前給出預(yù)警信號。因此，如果將該工具應(yīng)用于鉆井監(jiān)測系統(tǒng)中，將有助于避免卡鉆事故的發(fā)生。

5 結(jié)論

(1) 研究了一種實時卡鉆預(yù)警技術(shù)，該技術(shù)將鉆井模型和機器學(xué)習(xí)相結(jié)合，通過分析實時鉆井?dāng)?shù)據(jù)對卡鉆進行實時監(jiān)測，并區(qū)分即將發(fā)生的卡鉆是否與井眼清潔有關(guān)。

(2) 將摩阻扭矩模型和瞬態(tài)巖屑運移模型相耦合，充分考慮了巖屑床對鉆桿的影響。

(3) 使用貝葉斯優(yōu)化訓(xùn)練鉆井模型，使鉆井模型能夠適應(yīng)操作參數(shù)的變化不斷進行自我調(diào)整的同時，也能不斷提高模型預(yù)測精度。

(4) 利用時序數(shù)據(jù)分析監(jiān)測實時數(shù)據(jù)的異常趨勢，通過計算卡鉆風(fēng)險指數(shù)，為卡鉆事件提供直接的實時監(jiān)控和警報，可幫助工程技術(shù)人員及時采取有效措施避免卡鉆事故發(fā)生。