楊金顯，趙 淳

(河南理工大學(xué) 電氣工程與自動化學(xué)院，河南 焦作 454000)

隨鉆測量MWD(measurement while drilling)是指鉆機(jī)在鉆進(jìn)時連續(xù)不斷地檢測有關(guān)鉆孔或鉆頭的信息[1]。三軸加速度計是MWD核心測量器件，實際工作中不斷向控制系統(tǒng)提供井斜信息，并對工具面角的解算提供重要信息。由于加速度計距離鉆頭很近，其輸出信號含大量振動干擾，導(dǎo)致隨鉆測量結(jié)果精度差，因此需要建立隨鉆非線性振動誤差模型并予以補(bǔ)償。

目前在隨鉆加速度計振動誤差信號處理上，常見方法有以下幾種。孫云濤等[2]在井斜動態(tài)測量時對重力加速度固有頻帶外的振動干擾采用低通濾波器處理，實際上振動加速度的頻帶范圍較寬，在固有頻帶內(nèi)仍有振動干擾，且對姿態(tài)估計實時性有很大影響[3]。徐寶昌等[4]提出在旋轉(zhuǎn)導(dǎo)向鉆進(jìn)中建立鉆具振動加速度噪聲模型，采用改進(jìn)UKF濾波去除有色噪聲，但該方法主要考慮了鉆具軸向振動影響，對于復(fù)雜化振動干擾將出現(xiàn)較大的解算誤差。王濤等[5]提出了適用于隨鉆測量振動條件的加速度計36階Kalman濾波器，測量精度提高，但需要對各項器件誤差進(jìn)行分離、辨識與補(bǔ)償。Xue等[6]建立一種連續(xù)測量數(shù)學(xué)模型來減少加速度計動態(tài)測量誤差，在鉆頭不同運(yùn)動模式下使用不同信號處理方法。Yang等[7]提出一種借助磁和慣性數(shù)據(jù)提取鉆具重力信息的方法，振動干擾顯著時精度不高。

考慮到加速度計振動誤差模型的辨識難點在于復(fù)雜振動干擾，不確定激勵與響應(yīng)的關(guān)系，人為構(gòu)造模型難以描述過程特征。因此，本文旨在通過黑箱模型估計振動誤差，基于改進(jìn)Hammerstein(以下簡稱H)模型建立非線性振動函數(shù)，利用較少的輸入輸出辨識數(shù)據(jù)，實現(xiàn)對加速度計振動誤差建模，具有簡單、易辨識、計算量少、較好反應(yīng)過程特征等優(yōu)點。

1 振動誤差建模

1.1 建模方法

鉆柱振動受鉆井參數(shù)、鉆具組合、井深結(jié)構(gòu)、鉆柱結(jié)構(gòu)等因素影響[8]。鉆具振動主要分為縱向振動、橫向振動和扭轉(zhuǎn)振動[9]。實際隨鉆工作中，振動的表現(xiàn)形式更加多樣[10]。Lobo等[11]證實由于不確定性和環(huán)境因素，用隨機(jī)模型建立鉆柱振動模型。Tian等[12]對鉆柱粘滑振動特性建立扭轉(zhuǎn)動力學(xué)模型。Mohammadzadeh等[13]研究鉆柱的全耦合非線性振動動態(tài)模型。鉆柱振動的非線性特性復(fù)雜，振動沖擊對鉆柱的全部干擾無從得知，人為構(gòu)造模型難以描述過程特征，因此選擇H黑箱模型分析振動特性。目前，H模型因其特殊結(jié)構(gòu)將非線性控制問題簡化為線性模型預(yù)測控制問題，廣泛應(yīng)用于非線性系統(tǒng)識別、信號處理和控制工程實踐當(dāng)中[14]，如多傳感器系統(tǒng)、電力系統(tǒng)和神經(jīng)系統(tǒng)等。

井斜測量時，隨鉆測量系統(tǒng)緊鄰鉆頭安裝[15]。MWD加速度計安裝于鉆桿內(nèi)，各加速度計的安裝軸向如圖1所示。

圖1 加速度計安裝示意圖Fig.1 Installation diagram of the accelerometer

隨鉆系統(tǒng)進(jìn)行作業(yè)時，實際鉆進(jìn)軌跡與設(shè)計軌跡會出現(xiàn)偏差，因此傳感器數(shù)據(jù)準(zhǔn)確對軌跡修正至關(guān)重要。鉆進(jìn)軌跡的傾斜角(傾角、井斜角)指的是井眼軸線的切線與鉛垂直線的夾角，完全由加速度解算得到?；趫D1三軸加速度計安裝坐標(biāo)系，建立圖2隨鉆鉆進(jìn)圖。

圖2 隨鉆鉆進(jìn)圖Fig.2 Diagram of drilling while drilling

鉆進(jìn)過程中，鉆具的強(qiáng)烈振動使三軸加速度計的輸出包含大量振動干擾，導(dǎo)致姿態(tài)解算精度不高。為了從含有大量振動信號的原始信號中得到所需的有用數(shù)據(jù)，需要對三軸加速度信號進(jìn)行處理。

鉆柱振動信號是典型的非平穩(wěn)隨機(jī)信號[16]。加速度傳感器受到的振動干擾也是隨機(jī)信號，顯著特點是非周期性和瞬時值不易預(yù)測，但其統(tǒng)計特性是有規(guī)律的[17]。為能在井底復(fù)雜環(huán)境下，加速度計仍能得到精確的井斜角，需要對隨機(jī)振動特性進(jìn)行誤差分析。本文基于H模型，從中間變量不可測的角度，視振動輸入為加速度計的一種有色噪聲，表示為線性動態(tài)模塊與非線性靜態(tài)模塊串聯(lián)的形式，然后分別進(jìn)行參數(shù)辨識[18]。由于H模型辨識參數(shù)少，且較好的解決了時變性、非線性問題，提高了辨識精度，該模型如圖3所示。

圖3 隨鉆振動的H模型Fig.3 H model of vibration while drilling

其中，加速度計隨鉆振動H模型由非線性環(huán)節(jié)、線性環(huán)節(jié)和隨機(jī)振動干擾信號組成；非線性部分N(.)以uN真實加速度輸入作為輸入信號，非線性環(huán)節(jié)輸出N(uN)，表示含非線性誤差的有效加速度；線性部分L(.)以N(uN)為輸入信號，xN為加速度計受到離心力等影響的動態(tài)轉(zhuǎn)動加速度，尚未受到振動干擾，是未知中間變量；ξN為隨鉆隨機(jī)振動干擾下輸入的有色噪聲；yN為加速度計隨機(jī)振動實際輸出。

圖3(a)中，需要將待辨識模型中的有色噪聲ξN白化為白噪聲εN。非線性環(huán)節(jié)輸出N(uN)作為簡化模型圖3(b)的輸入信號，xN為過程模型輸出，εN為零均值白噪聲，ξN為有色噪聲，yN為含噪輸出加速度，A(z-1)，B(z-1)，C(z-1)單位后移算子z-1多項式。

從圖3可知，H模型輸出可表示為

yN=xN+ξN

(1)

(2)

(3)

式中：βi為非線性環(huán)節(jié)基函數(shù)g(.)的系數(shù)；r為加速度計噪聲模型中非線性環(huán)節(jié)的階數(shù)。

過程模型

(4)

噪聲模型

(5)

式中：A(z-1)=1+a1z-1+a2z-2+…+anaz-na；B(z-1)=b1z-1+b2z-2+…+bnbz-nb；C(z-1)=1+c1z-1+c2z-2+…+cncz-nc。

將式(3)代入式(2)得到

(6)

式中：φN=(-yN-1,…,-yN-Na,g1(uN-1),…,g1(uN-Nb)；g2(uN-1),…,g2(uN-nb),…,gm(uN-nb)；εN-1,…,εN-nc)T；θ=(a1,…,ana,β1b1,…,β1bnb,β2b1,…,β2bnb,…,βrb1,…,βrbnb,c1,…,cnc)T。

圖4 加速度計信號處理Fig.4 Accelerometer signal processing

具體算法步驟：

步驟1利用測量數(shù)據(jù)集進(jìn)行模型參數(shù)辨識，得到非線性結(jié)構(gòu)參數(shù)初值c1,…,cnc、線性結(jié)構(gòu)參數(shù)初值a1,…,ana以及b1,…,bnb；

步驟2采樣當(dāng)前系統(tǒng)輸出，得到最新的誤差Δδ；

步驟3通過與條件閾值Δλ比較更新調(diào)節(jié)系數(shù)γ；

步驟4得到最新的H誤差模型；

步驟5計算加速度值；

步驟6下一采樣點，返回步驟2繼續(xù)循環(huán)。

1.2 參數(shù)辨識

定義參數(shù)向量

θ=[θaθβbθc]T

(7)

式中：θa=[a1,a2,…,ana]T；θβb=[β1b1,…,β1bnb,β2b1,…,βrbnb]T；θc=[c1,c2,…,cnc]T；nβb=nβ*nb，na=nb，n=na+nβb+nc。

定義信息向量

(8)

式中，φNs=[-yN-1,…,-yN-na,g1(uN-1),…,g1(uN-nb),…,gr(uN-nb)]T；φNn=[ξN-1,…,ξN-nc]T。

(9)

根據(jù)貝葉斯公式，參數(shù)的后驗概率密度函數(shù)為

(10)

(11)

(12)

綜上，給定隨鉆隨機(jī)振動時的實測輸入、輸出數(shù)據(jù)集DN，參數(shù)的后驗概率密度

(13)

(14)

(15)

(16)

從上式可以得出，當(dāng)前時刻的噪聲方差可近似等于上一時刻的噪聲方差

(17)

估計參數(shù)向量

(18)

增益矩陣

(19)

協(xié)方差矩陣

(20)

(21)

經(jīng)過適當(dāng)?shù)木仃囎儞Q，得到βi的表達(dá)式

(22)

2 模型修正

(23)

考慮到實際隨鉆測量振動環(huán)境的強(qiáng)時變性、強(qiáng)非線性，還有隨著上述隨鉆振動模型參數(shù)遞推辨識算法對時變參數(shù)追蹤能力的增強(qiáng)，也會造成對噪聲愈加敏感，因此提出引入條件閾值Δλ、調(diào)節(jié)系數(shù)γ，將隨鉆系統(tǒng)的鉆頭運(yùn)動狀態(tài)與調(diào)節(jié)系數(shù)關(guān)聯(lián)起來，對隨機(jī)振動模型進(jìn)行修正，提高辨識靈敏度

(24)

(25)

式(25)表征鉆頭的運(yùn)動狀態(tài)與調(diào)節(jié)系數(shù)之間的反比關(guān)系。將初始調(diào)節(jié)參數(shù)與鉆頭運(yùn)動狀態(tài)關(guān)聯(lián)起來，可以拓寬調(diào)節(jié)系數(shù)的變化范圍，使參數(shù)遞推辨識過程更細(xì)致地跟隨鉆頭的隨機(jī)振動過程，從而建立更加精確的隨機(jī)振動誤差模型。

(1)Δδ<Δλ，認(rèn)為當(dāng)前隨鉆振動模型已經(jīng)比較準(zhǔn)確，故保持模型不變。

(2)Δλ<Δδ，認(rèn)為當(dāng)前隨鉆振動模型無法準(zhǔn)確描述樣本信息，需要自適應(yīng)更新參數(shù)。

(26)

式中：初始調(diào)節(jié)參數(shù)γ0與鉆頭運(yùn)動狀態(tài)的關(guān)系可以用分段函數(shù)表示，δ1，δ2，Δγ1，Δγ2是與鉆頭運(yùn)動狀態(tài)相關(guān)的常值。Δδ/Δλ≤1表征鉆頭位于靜止?fàn)顟B(tài)，1<Δδ/Δλ≤δ1，表征鉆頭工作位于低動態(tài)，δ1<Δδ/Δλ≤δ2，表征鉆頭處于高動態(tài)，實際隨鉆工作中鉆頭的轉(zhuǎn)速通常處于低動態(tài)。Δδ為當(dāng)前時刻系統(tǒng)實際輸出值與估計的殘差，Δλ為條件閾值，γ0為默認(rèn)調(diào)節(jié)系數(shù)，算法推導(dǎo)流程如圖5。

圖5 算法流程Fig.5 Algorithm flow

參數(shù)向量

(27)

增益矩陣

(28)

協(xié)方差矩陣

(29)

3 實驗驗證及分析

3.1 鉆進(jìn)模擬實驗

為驗證算法性能，在實驗室條件下進(jìn)行垂直鉆進(jìn)模擬實驗。實驗場地條件為北緯35°，東經(jīng)113°，地球自轉(zhuǎn)角速度15°/h，地球重力加速度為9.8 m/s2。實驗所需設(shè)備有某型工程用小型鉆機(jī)，MEMS隨鉆測量模塊和模擬地質(zhì)層等。將隨鉆測量核心模塊(已標(biāo)定)固定在小型鉆機(jī)鉆桿上，隨鉆桿一起轉(zhuǎn)動，如圖6所示。根據(jù)隨鉆工作時加速度計的隨機(jī)振動誤差特性，選取該模型線性環(huán)節(jié)和非線性環(huán)節(jié)階數(shù)分別為r=2，na=nb=nc=3，垂直鉆進(jìn)實驗之前，調(diào)整鉆桿垂直并固定，測量靜態(tài)初始加速度；鉆進(jìn)實驗中，以轉(zhuǎn)速360°/s及恒定鉆壓進(jìn)行垂直鉆進(jìn)，以200 Hz采樣頻率采集加速度數(shù)據(jù)作為初始數(shù)據(jù)集，將測量數(shù)據(jù)集DN按照1.2節(jié)分析流程估計模型參數(shù)，最后不斷采樣當(dāng)前系統(tǒng)輸出，更新模型參數(shù)；實驗結(jié)束后，選取其中1 000組數(shù)據(jù)進(jìn)行算法處理。

圖6 鉆進(jìn)模擬Fig.6 Drilling simulation

由前文分析可知，真實隨鉆加速度的獲取，依賴振動加速度ξ的確定。實際應(yīng)用中，已知振動誤差模型階數(shù)na=3，nc=3，由初始數(shù)據(jù)集確定初始參數(shù)后，不斷采樣更新模型參數(shù)向量θN，分解得到參數(shù){c1,c2,…,cn}，根據(jù)噪聲模型計算得到當(dāng)前時刻的振動加速度ξN=C(z-1)/A(z-1)εN，算法處理前后加速度曲線的對比如圖7所示。由圖可知，隨機(jī)振動下加速度計的輸出誤差具有顯著波動性和隨機(jī)性，變化范圍始終保持在一定區(qū)間內(nèi)，但趨勢項并不明顯，所建立的隨鉆隨機(jī)振動模型呈現(xiàn)出和原始序列相同趨勢，擬合效果較為理想。

圖7 模擬鉆進(jìn)加速度曲線Fig.7 Vibration acceleration curves

為進(jìn)一步精確驗證所得模型的有效性，將隨機(jī)振動狀態(tài)下的實際輸出與基于非線性模型的預(yù)測輸出對比，利用標(biāo)準(zhǔn)均方根誤差(NRSME)適應(yīng)度值評價模型預(yù)測輸出與實測輸出的擬合程度，如表1所示。

表1 各軸向數(shù)據(jù)NRSME適應(yīng)度值Tab.1 NRSME fitness value of each axial data

(30)

通過NRSME適應(yīng)度值表可看出基于H模型的自適應(yīng)算法預(yù)測誤差較小，擬合程度較高，該模型能較好描述隨鉆隨機(jī)振動情況下加速度計的振動動態(tài)特性，驗證了建模方法的有效性。

3.2 鉆進(jìn)實驗

為進(jìn)一步證明基于H模型的隨鉆振動誤差模型的可行性，選擇焦作某煤礦的垂直鉆進(jìn)實驗測試數(shù)據(jù)進(jìn)行實驗驗證，部分實鉆數(shù)據(jù)如表2所示，高精度測井儀器測量的實際軌跡曲線和測量軌跡曲線對比如圖8所示?？紤]隨鉆加速度計振動誤差特性，選取該模型線性環(huán)節(jié)和非線性環(huán)節(jié)階數(shù)取值同前，垂直鉆進(jìn)開始時采集少量輸入輸出數(shù)據(jù)作為初始數(shù)據(jù)集，用測量數(shù)據(jù)集辨識模型參數(shù)，并采樣系統(tǒng)當(dāng)前輸出，在線更新加速度計的隨鉆振動誤差模型。

表2 部分實鉆加速度數(shù)據(jù)Tab.2 Partially drilled data

圖8 鉆進(jìn)軌跡Fig.8 Drilling trajectory

通過前述流程計算相應(yīng)的振動加速度，對比基于改進(jìn)H模型的隨鉆振動誤差模型和原始垂直鉆進(jìn)的三軸加速度計輸出，產(chǎn)生了如圖9所示的加速度擬合曲線。從圖中可看出，基于改進(jìn)H模型得出的預(yù)測序列，與原始序列趨勢相同，擬合效果較為理想，實鉆加速度數(shù)據(jù)算法處理后的NRSME適應(yīng)度值表如表3所示。

圖9 實鉆加速度曲線Fig.9 True drilling acceleration curves

表3 各軸向數(shù)據(jù)NRSME適應(yīng)度值Tab.3 NRSME fitness value of each axial data

從表3可直觀看出基于改進(jìn)H模型的自適應(yīng)算法預(yù)測誤差較小，擬合程度較高，說明自適應(yīng)非線性H模型能夠較好地描述隨鉆隨機(jī)振動情況下加速度計的振動動態(tài)特性，驗證了建模方法的有效性。

最后，隨鉆加速度計振動誤差建模成功后，對測量數(shù)據(jù)的處理是否符合實際的軌跡曲線角度，通過計算算法補(bǔ)償后的井斜角與實際軌跡角度的對比圖來說明補(bǔ)償效果，如圖10所示。

圖10 實鉆井斜角曲線圖Fig.10 Actual drilling angle curve

從圖8～10可知，基于改進(jìn)H模型的自適應(yīng)算法預(yù)測誤差準(zhǔn)確，數(shù)據(jù)處理后的測量值與實際軌跡曲線較為一致，實驗證明自適應(yīng)非線性H模型能夠較好地描述隨鉆隨機(jī)振動情況下加速度計的振動動態(tài)特性，驗證了建模方法的有效性。

4 結(jié) 論

目前對隨鉆測量過程中鉆具的實際振動進(jìn)行測試還是有一定困難的，針對隨鉆測量時鉆具振動使加速度測量信號含大量振動干擾，測量精度不高的問題，提出一種建立非線性隨鉆振動誤差模型的方法來獲取真實隨鉆加速度。本文設(shè)計一種改進(jìn)的基于自適應(yīng)H模型算法，建立加速度隨機(jī)振動誤差模型，利用遞推貝葉斯算法辨識參數(shù)，引入條件閾值和調(diào)節(jié)系數(shù)與鉆頭運(yùn)動狀態(tài)關(guān)聯(lián)起來，提高辨識靈敏度，最后設(shè)計實驗研究MWD真實振動環(huán)境下輸出加速度的隨機(jī)振動誤差規(guī)律，實驗結(jié)果表明很好地實現(xiàn)對加速度計隨機(jī)振動誤差序列的擬合估計，說明建模方法是有效的，對工程應(yīng)用中隨鉆隨機(jī)振動誤差的建模與實時補(bǔ)償具有重要的參考價值。