耿艷峰，孫 建，王偉亮，程民彪

（中國(guó)石油大學(xué)(華東)控制科學(xué)與工程學(xué)院，青島 266580）

旋轉(zhuǎn)導(dǎo)向鉆井技術(shù)是20 世紀(jì)90 年代發(fā)展起來(lái)的一項(xiàng)自動(dòng)化鉆井技術(shù)[1]。與傳統(tǒng)的滑動(dòng)導(dǎo)向鉆井相比，旋轉(zhuǎn)導(dǎo)向鉆井工具的井眼凈化效果好，井身軌跡控制精度高，位移延伸能力強(qiáng)[1,2]。工具面角參數(shù)是導(dǎo)向鉆井工具中的一項(xiàng)重要的姿態(tài)參數(shù)，表征了鉆頭的鉆進(jìn)方向，實(shí)現(xiàn)工具面角的準(zhǔn)確可靠測(cè)量是提高旋轉(zhuǎn)導(dǎo)向鉆井工具性能的前提。常用的工具面角測(cè)量傳感器包括：加速度計(jì)、磁通門、陀螺儀。大量研究表明單一傳感器往往難以實(shí)現(xiàn)準(zhǔn)確測(cè)量[3-5]，通過(guò)融合不同傳感器數(shù)據(jù)可以顯著提高測(cè)量精度，互補(bǔ)濾波和自適應(yīng)卡爾曼濾波[6,7]兩種數(shù)據(jù)融合方法的應(yīng)用較為常見。鉆井過(guò)程噪聲復(fù)雜，既有鉆具軸向、周向以及扭轉(zhuǎn)振動(dòng)產(chǎn)生的周期性噪聲、又有隨著鉆井深度、壓力增加的趨勢(shì)性噪聲，噪聲分布往往是未知的?？紤]到鉆井過(guò)程參數(shù)范圍基本可以確定，將非周期性噪聲、模型線性化、陀螺儀漂移等影響因素視為幅值有界但分布未知誤差(unknown but bounded error, UBB)進(jìn)行處理具有理論上的可行性。針對(duì)UBB 誤差的集員估計(jì)理論近年來(lái)得到了快速發(fā)展，文獻(xiàn)[8]利用基于橢球集合的擴(kuò)展集員卡爾曼濾波算法實(shí)現(xiàn)車輛位置和速度估計(jì)；文獻(xiàn)[9]提出一種基于中心對(duì)稱多胞體的集員狀態(tài)估計(jì)算法實(shí)現(xiàn)一類包含未知輸入和系統(tǒng)干擾系統(tǒng)的狀態(tài)估計(jì)；文獻(xiàn)[10]提出一種基于觀測(cè)器誤差界分區(qū)估計(jì)的集員估計(jì)方法，實(shí)現(xiàn)一類包含未知但有界干擾和噪聲的廣義系統(tǒng)的狀態(tài)估計(jì)。

高溫、高壓、強(qiáng)振動(dòng)的鉆井環(huán)境極易導(dǎo)致工具面角傳感器出現(xiàn)故障，進(jìn)一步可能引發(fā)其它事故。如果故障無(wú)法在井下進(jìn)行處理，就要進(jìn)行起鉆維修，嚴(yán)重影響了導(dǎo)向鉆井技術(shù)的推廣應(yīng)用。文獻(xiàn)[11]中詳細(xì)分析了旋轉(zhuǎn)導(dǎo)向穩(wěn)定平臺(tái)可能發(fā)生的故障類型，并針對(duì)不同傳感器提出了相應(yīng)的故障檢測(cè)與處理方法，基本解決了加速度計(jì)的故障處理問(wèn)題，但對(duì)于陀螺儀故障缺少相應(yīng)的故障診斷和處理方法。實(shí)際系統(tǒng)中，故障的幅值信息是后續(xù)處理的重要基礎(chǔ)，能夠提供幅值信息的故障估計(jì)問(wèn)題得到越來(lái)越多的學(xué)者關(guān)注。文獻(xiàn)[12]使用無(wú)跡卡爾曼濾波器對(duì)故障參數(shù)進(jìn)行估計(jì)，實(shí)現(xiàn)故障發(fā)生后對(duì)控制器進(jìn)行重構(gòu)；文獻(xiàn)[13]在考慮模型不確定性以及環(huán)境噪聲的條件下，基于H∞理論利用線性矩陣不等式設(shè)計(jì)非線性觀測(cè)器研究了一類非線性系統(tǒng)的故障檢測(cè)、隔離和估計(jì)方法；文獻(xiàn)[14]提出一種新的滑模觀測(cè)器實(shí)現(xiàn)不滿足常見故障估計(jì)條件的系統(tǒng)的故障估計(jì)。

基于上述分析，本文嘗試?yán)米赃m應(yīng)擴(kuò)展集員卡爾曼濾波(adaptive set-membership extended Kalman filtering, ASM-EKF)算法解決導(dǎo)向鉆井工具面角傳感器的故障估計(jì)與處理問(wèn)題。將陀螺儀故障項(xiàng)增廣為狀態(tài)變量，利用加速度計(jì)建立測(cè)量方程來(lái)估計(jì)陀螺儀故障，并對(duì)陀螺儀輸出進(jìn)行校正，實(shí)現(xiàn)陀螺儀故障條件下工具面角的可靠測(cè)量。最后通過(guò)模擬樣機(jī)測(cè)試與仿真驗(yàn)證方案的有效性，為故障條件下工具面角的精確測(cè)量與控制系統(tǒng)的正常運(yùn)行奠定基礎(chǔ)。

1 工具面角測(cè)量系統(tǒng)濾波動(dòng)態(tài)模型

1.1 工具面角測(cè)量方案

以動(dòng)態(tài)指向式旋轉(zhuǎn)導(dǎo)向鉆井工具為例，其工具面角的傳感器布局方案參見文獻(xiàn)[7]，如圖1 所示。采用雙加速度計(jì)和陀螺儀實(shí)現(xiàn)組合測(cè)量，雙加速度計(jì)對(duì)稱布局消除運(yùn)動(dòng)加速度的影響，消除運(yùn)動(dòng)加速度后的各軸加速度測(cè)量值用αx，αy，αz表示，加速度計(jì)的x軸平行于穩(wěn)定平臺(tái)軸線，y軸指向穩(wěn)定平臺(tái)切線，z軸指向穩(wěn)定平臺(tái)軸心。

圖1 穩(wěn)定平臺(tái)傳感器示意圖Fig.1 The diagram of stabilized platform sensors

則工具面角可由式(1)計(jì)算：

陀螺儀敏感軸與穩(wěn)定平臺(tái)軸向平行，用來(lái)測(cè)量穩(wěn)定平臺(tái)對(duì)地轉(zhuǎn)速，陀螺儀測(cè)量工具面角的公式為：

其中：φk是k時(shí)刻的工具面角，ωk為k時(shí)刻的陀螺儀輸出轉(zhuǎn)速，Ts為采樣時(shí)間。

1.2 工具面角組合測(cè)量模型

考慮陀螺儀輸出發(fā)生加性故障，根據(jù)式(2)獲得校正后的工具面角的迭代形式如下：

其中：fk表示陀螺儀的加性故障。

將工具面角和陀螺儀故障作為系統(tǒng)狀態(tài)，陀螺儀輸出作為系統(tǒng)輸入，建立如下的狀態(tài)方程：

根據(jù)式(1)建立雙加速度計(jì)量測(cè)方程，結(jié)合式(4)的狀態(tài)方程，建立工具面角測(cè)量系統(tǒng)的狀態(tài)空間模型：

2 自適應(yīng)擴(kuò)展集員卡爾曼濾波器

2.1 擴(kuò)展集員卡爾曼濾波算法

擴(kuò)展集員卡爾曼濾波算法是在卡爾曼濾波框架下，引入未知但有界誤差項(xiàng)來(lái)提升算法估計(jì)性能。擴(kuò)展卡爾曼濾波算法需要對(duì)非線性模型線性化，將線性化誤差看作未知但有界誤差項(xiàng)描述為橢球集合，并與模型中涉及的未知但有界誤差求取閔可夫斯基和(Minkowski sum)，獲得新的未知但有界誤差項(xiàng)[8,15]。模型(5)的量測(cè)方程為非線性方程，線性化造成的未知但有界誤差項(xiàng)描述為。

對(duì)量測(cè)方程線性化，將式(5)變換為如下形式：

具有未知但有界誤差項(xiàng)和隨機(jī)誤差項(xiàng)的非線性系統(tǒng)的最優(yōu)卡爾曼濾波增益、一步預(yù)測(cè)協(xié)方差以及新息協(xié)方差計(jì)算公式如下[8]：

2.2 UBB 誤差的數(shù)學(xué)描述

由式(5)可知，狀態(tài)方程是線性方程，不存在線性化誤差。因此，狀態(tài)方程的未知但有界誤差考慮陀螺儀的漂移以及模型誤差，將描述為如下橢球集合：

量測(cè)方程存在線性化誤差，將線性化誤差作為未知但有界誤差項(xiàng)考慮到濾波算法中，用來(lái)描述。量測(cè)方程的線性化過(guò)程如下：

陀螺儀具有多種故障形式，當(dāng)發(fā)生突變故障時(shí)，新息協(xié)方差矩陣Pk會(huì)突然增加，導(dǎo)致濾波模型中的測(cè)量UBB 誤差增大，進(jìn)而影響濾波算法收斂速度，甚至導(dǎo)致發(fā)散。在濾波過(guò)程趨于穩(wěn)態(tài)值后，新息矩陣趨于極小值，可以根據(jù)新息矩陣的大小來(lái)判斷是否需要調(diào)節(jié)預(yù)測(cè)值，此處引入如下評(píng)價(jià)指標(biāo)：

式中：ke為測(cè)量誤差矩陣，表示新息協(xié)方差矩陣的F范數(shù)，為測(cè)量誤差矩陣的2 范數(shù)，γ根據(jù)實(shí)際工況進(jìn)行選取，γ值越大，故障估計(jì)的快速性越好，濾波算法的穩(wěn)定性越差。

取nr=1，采用外邊界橢球來(lái)描述拉格朗日余項(xiàng)，同時(shí)保證正定，將進(jìn)行如下描述：

其中：Hesk表示量測(cè)方程的Hessian矩陣，c1∈(0,1]為調(diào)節(jié)因子，若式(14)滿足，則認(rèn)為濾波狀態(tài)穩(wěn)定，c1=1 ；否則根據(jù)實(shí)際工況在取值范圍內(nèi)選取。

2.3 集員Sage-Husa 自適應(yīng)濾波算法

井下工況復(fù)雜多變，導(dǎo)向工具鉆進(jìn)的過(guò)程中，地質(zhì)的變化、穩(wěn)定平臺(tái)運(yùn)動(dòng)狀態(tài)的變化、振動(dòng)強(qiáng)度的變化等情況，都會(huì)導(dǎo)致系統(tǒng)噪聲和測(cè)量噪聲的變化，Sage-Husa 自適應(yīng)濾波算法可以不斷地調(diào)整系統(tǒng)模型參數(shù)和噪聲模型統(tǒng)計(jì)參數(shù)，使濾波性能達(dá)到最優(yōu)，常規(guī)的 Sage-Husa 算法都是針對(duì)傳統(tǒng)卡爾曼濾波模型[16,17]，下面基于集員卡爾曼框架對(duì)該算法進(jìn)行改進(jìn)。

結(jié)合集員卡爾曼濾波框架以及式(7)的數(shù)學(xué)模型可得：

由式(7)和式(16)可得：

對(duì)式(17)左右同乘對(duì)應(yīng)轉(zhuǎn)置，求取數(shù)學(xué)期望可得：

對(duì)式(18)等式右邊項(xiàng)進(jìn)行處理：

由式(9)(18)結(jié)合式(19)～(22)變換，Qk-1替換為Qk得：

根據(jù)式(16)可得：

故測(cè)量噪聲Rk的自適應(yīng)算法為：

為提升算法快速性，提出如下自適應(yīng)方案：

綜上所述，自適應(yīng)擴(kuò)展集員卡爾曼濾波算法的具體步驟如下：

步驟2：根據(jù)式(14)、式(15)更新測(cè)量UBB 誤差；

步驟4：根據(jù)式(8)更新最優(yōu)卡爾曼增益Kk；

步驟5：根據(jù)式(9)更新一步預(yù)測(cè)協(xié)方差Pk/k-1；

步驟6：根據(jù)式(10)更新協(xié)方差陣Pk；

步驟7：根據(jù)式(26)更新過(guò)程噪聲Qk和測(cè)量噪聲Rk

3 仿真與實(shí)驗(yàn)

粘滑運(yùn)動(dòng)是鉆進(jìn)過(guò)程中常見的工況，頻率在0.05～0.5 Hz 之間[18]，根據(jù)穩(wěn)定平臺(tái)測(cè)控系統(tǒng)性能要求，在0.5 Hz，200%粘滑度的工況下，工具面角波動(dòng)范圍應(yīng)小于15 °[11]，這對(duì)算法的跟蹤性能提出了較高的要求。為了測(cè)試更嚴(yán)苛的條件，考慮井斜角為90 °的情況；使用0.5 Hz，波動(dòng)幅值20 °的正弦模擬粘滑工況進(jìn)行測(cè)試。

3.1 典型鉆井工況仿真

假設(shè)導(dǎo)向鉆井工具工作在造斜工況，工具面角固定在183 °；30 s 后發(fā)生粘滑，工具面角在183 °基準(zhǔn)上以頻率為0.5 Hz，幅值為20 °的正弦波動(dòng)；60 s 后改成穩(wěn)斜鉆進(jìn)，工具面角以36 ° /s的速度在0～360 °范圍內(nèi)線性變化，90 s 后繼續(xù)造斜，工具面角維持在183 °持續(xù)到120 s。Baker Hughes 是世界三大油服公司之一，其井下振動(dòng)手冊(cè)[19]將鉆具振動(dòng)分為7 個(gè)等級(jí)，其中處于良好安全運(yùn)行環(huán)境下的鉆具振動(dòng)級(jí)別為0 級(jí)[0, 0.5 g)以及1 級(jí)[0.5 g,1 g)，因此本文選取加速度計(jì)的振動(dòng)噪聲方差為0.5g2，模擬1 級(jí)鉆具振動(dòng)條件下的加速度計(jì)信號(hào)輸出情況；參考陀螺儀手冊(cè)[20]，陀螺儀的速率噪聲密度為陀螺儀噪聲帶寬設(shè)置為1600 Hz，故將陀螺儀測(cè)量噪聲方差設(shè)置為100( ° /s)2，模擬陀螺儀固有噪聲。

為驗(yàn)證ASM-EKF 性能以及陀螺儀故障估計(jì)性能，分別在如表1 所示三種陀螺儀故障狀態(tài)下進(jìn)行仿真。采用雙加速度計(jì)測(cè)量方案、陀螺儀測(cè)量方案以及自適應(yīng)擴(kuò)展集員卡爾曼濾波算法來(lái)對(duì)上述工況以及故障類型下的工具面角的進(jìn)行解算。

表1 仿真故障類型Tab.1 Simulation fault type

圖2 是無(wú)故障情況下三種方案的處理效果，加速度計(jì)解算受噪聲影響非常劇烈，解算工具面角的RMSE 為50.52 °，無(wú)法反映實(shí)際工具面角；陀螺儀解算相對(duì)精確，但陀螺儀解算工具面角值存在0.1 /s° 的漂移，120 s 后存在13.74 °的解算誤差，難以滿足性能要求；ASM-EKF 算法則利用數(shù)據(jù)融合的優(yōu)勢(shì)，精準(zhǔn)地實(shí)現(xiàn)了對(duì)工具面角的估計(jì)，RMSE 為2.55 °。圖3展示了加性故障I 下的解算效果，陀螺儀由于故障原因失去解算能力，ASM-EKF 算法仍然能夠精確的進(jìn)行工具面角估計(jì)，RMSE 為2.49 °，并且準(zhǔn)確地估計(jì)出陀螺儀的故障。

圖2 無(wú)故障仿真效果Fig.2 Fault-free simulation effect

圖3 加性故障I 仿真效果Fig.3 Additive fault I simulation effect

圖4 展示了加性故障II 下的解算效果，突變故障下，ASM-EKF 算法的工具面角估計(jì)值有短暫的偏離，但是很快調(diào)整到實(shí)際值附近波動(dòng)，RMSE 為3.57 °，故障估計(jì)值很快跟蹤上突變故障。

圖4 加性故障II 仿真效果Fig.4 Additive fault II simulation effect

圖5 展示了加性故障III 下的解算效果，驗(yàn)證了ASM-EKF 算法對(duì)于10 °/s 范圍內(nèi)故障下的工具面角仍具有較好的估計(jì)效果，RMSE 為2.55 °，與無(wú)故障時(shí)的精度相當(dāng)。陀螺儀的故障估計(jì)值可以反映故障變化趨勢(shì)，但是相對(duì)于故障幅值具有較大的波動(dòng)，這是由于陀螺儀固有噪聲與故障幅值處于同一量級(jí)所致。此時(shí)故障已在陀螺噪聲范圍內(nèi)，嚴(yán)格意義地講，其已不能按故障進(jìn)行處理。綜上說(shuō)明ASM-EKF 算法具有大幅值陀螺儀加性故障的快速跟蹤估計(jì)能力，可以實(shí)現(xiàn)多種復(fù)雜加性故障下工具面角的可靠測(cè)量。

圖5 加性故障III 仿真效果 Fig.5 Additive fault III simulation effect

3.2 模擬樣機(jī)實(shí)測(cè)結(jié)果

課題組自行開發(fā)的動(dòng)態(tài)指向式旋轉(zhuǎn)導(dǎo)向鉆井工具原理樣機(jī)參見圖6，樣機(jī)采用兩組三軸重力加速度計(jì)MMA8451 以及一個(gè)單軸陀螺儀ADXRS645 測(cè)量工具面角。

圖6 原理樣機(jī)結(jié)構(gòu)圖Fig.6 Principle prototype structure drawing

利用原理樣機(jī)測(cè)試導(dǎo)向鉆井工具四種工作狀態(tài)來(lái)驗(yàn)證該算法的實(shí)用性，同時(shí)分別與文獻(xiàn)[7]中采用的自適應(yīng)擴(kuò)展卡爾曼 (adaptive extended Kalman filtering,AEKF) 算法以及文獻(xiàn)[8]中提出的集員卡爾曼(set-membership extended Kalman filtering, SM-EKF)算法進(jìn)行對(duì)比，驗(yàn)證本文算法的性能。

實(shí)驗(yàn)?zāi)M的四種典型工況：恒工具面角、粘滑工況、穩(wěn)斜工況、調(diào)工具面角。實(shí)驗(yàn)中將原理樣機(jī)固定在振動(dòng)臺(tái)上模擬井下振動(dòng)。振動(dòng)臺(tái)開啟后，振動(dòng)頻率由0Hz 過(guò)渡到50Hz，之后持續(xù)在50Hz 振動(dòng)，振動(dòng)頻率的變化能夠更好地模擬井下實(shí)際振動(dòng)，同時(shí)增加了處理的復(fù)雜性。

陀螺儀故障主要針對(duì)持續(xù)故障以及突變故障兩種情況進(jìn)行模擬，具體故障模擬形式如表2 所示。

表2 實(shí)驗(yàn)故障類型Tab.2 Experimental fault type

利用上述三種算法獲得的工具面角估計(jì)值與陀螺儀故障估計(jì)值如圖7、8 所示。不同故障條件下，不同算法獲得的工具面角估計(jì)值與參考值的均方根誤差(RMSE)對(duì)比如表3 所示。噪聲方差以及未知但有界誤差方差的初始值設(shè)置為[8]：

表3 工具面角估計(jì)均方根誤差Tab.3 Root mean square error of toolface estimation

圖7 加性故障I 實(shí)測(cè)效果Fig.7 Additive fault I measurement effect

從圖7、圖8 中可以看出，對(duì)于加性故障I，三種算法都可以實(shí)現(xiàn)故障的跟蹤估計(jì)，文獻(xiàn)[7]所采用的AEKF 算法始終存在一定的偏離，導(dǎo)致工具面角的估計(jì)也存在一定程度的持續(xù)偏離。對(duì)于加性故障II，AEKF 算法完全失去故障估計(jì)能力，導(dǎo)致工具面角估計(jì)完全偏離，無(wú)法達(dá)到可靠測(cè)量要求；文獻(xiàn)[8]采用的SM-EKF 算法雖然可以快速跟蹤加性故障II，但是波動(dòng)太大，從表3 可以看出，該方法在一定程度上會(huì)造成較大的工具面角估計(jì)誤差，無(wú)法滿足實(shí)際鉆井過(guò)程對(duì)工具面角測(cè)量的要求；而本文提出的ASM-EKF 算法可以快速準(zhǔn)確地實(shí)現(xiàn)陀螺儀故障估計(jì)，在發(fā)生故障情況下，工具面角估計(jì)值的最大均方根誤差為3.841°，滿足實(shí)際鉆井技術(shù)需求。

圖8 加性故障II 實(shí)測(cè)效果Fig.8 Additive fault II measurement effect

4 結(jié) 論

為解決動(dòng)態(tài)指向式旋轉(zhuǎn)導(dǎo)向工具面角測(cè)量陀螺儀故障估計(jì)與處理問(wèn)題，本文提出了一種自適應(yīng)擴(kuò)展集員卡爾曼濾波算法，實(shí)現(xiàn)了在陀螺儀發(fā)生加性故障情況下仍能夠準(zhǔn)確地測(cè)量工具面角。仿真和實(shí)驗(yàn)測(cè)試結(jié)果驗(yàn)證了該方法的有效性。主要結(jié)論如下：

1）獲得了一種可靠測(cè)量工具面角的算法，在陀螺儀無(wú)故障時(shí)可以實(shí)現(xiàn)工具面角的準(zhǔn)確測(cè)量，在陀螺儀發(fā)生加性故障情況時(shí)可以實(shí)時(shí)估計(jì)陀螺儀故障幅值，并對(duì)測(cè)量值進(jìn)行補(bǔ)償修正。持續(xù)故障與突變故障下工具面角測(cè)量的RMSE 分別為2.075 °、3.841 °，滿足井下控制系統(tǒng)對(duì)工具面角測(cè)量精度的需求。

2）可以將線性化、陀螺儀漂移、井下振動(dòng)等未知因素對(duì)井下傳感器測(cè)量數(shù)據(jù)的影響等效成未知但有界誤差項(xiàng)，SM-EKF 算法既有高效跟蹤突變故障的能力、又有較高的濾波精度，可以用于同類井下測(cè)量數(shù)據(jù)的處理。

本文是在實(shí)驗(yàn)室進(jìn)行的模擬測(cè)試，在振動(dòng)強(qiáng)度以及振動(dòng)變化規(guī)律等方面與實(shí)際鉆井環(huán)境仍有差別，需要通過(guò)實(shí)際鉆井環(huán)境進(jìn)行測(cè)試，進(jìn)一步提升該算法的實(shí)用性。