仵 杰,雪宇超,路 濤,李高仁,陳廣識
(1.西安石油大學(xué) 電子工程學(xué)院,陜西 西安 710065; 2.光電油氣測井與信息檢測教育部重點實驗室,陜西 西安710065;3.中國石油集團(tuán)測井有限公司,陜西 西安 710077; 4.中國石油天然氣集團(tuán) 長慶油田勘探開發(fā)研究院,陜西 西安 710018)
在測井勘探中采集的井孔聲波測量值會受到圍巖影響,也會受到砂巖層、頁巖層引起的各向異性、不平衡地層應(yīng)力以及裂縫的影響,致使聲波測井的解釋變得復(fù)雜,并使巖層力學(xué)特性參數(shù)的估算產(chǎn)生偏差[1-2]。其中地層不平衡應(yīng)力是造成儲層油氣運移的主要原因,確定地層應(yīng)力參數(shù)及分布能有效預(yù)測儲層方位[3]。地層油氣運移會改變應(yīng)力方向上地層的聲學(xué)特性從而引起地層縱波、橫波等聲學(xué)各向異性[4-5]。關(guān)于地層縱波各向異性的研究,國內(nèi)外已做了大量的工作[6-7]。目前,定量分析地層橫波各向異性對聲波測井陣列波形的影響方面的研究較少。偶極聲波測井激發(fā)的彎曲波是描述不同方向、不同強(qiáng)度的地層異常應(yīng)力的關(guān)鍵手段[8]。與均勻地層相比,橫波各向異性地層情況下偶極子聲源激發(fā)的聲場更為復(fù)雜,無法采用解析的方法進(jìn)行模擬,一般采用三維有限元方法對不同地層各向異性下的偶極聲源激發(fā)的反射聲場影響因素進(jìn)行數(shù)值模擬[9]。
本文基于COMSOL仿真軟件實現(xiàn)了各向異性地層下偶極聲波測井陣列波形的數(shù)值模擬,研究不同應(yīng)力方向角及地層橫波各向異性強(qiáng)度對陣列波形的影響規(guī)律,并通過合適算法反演各向異性地層參數(shù), 以期為橫波各向異性地層下偶極聲波測井工程應(yīng)用解釋提供理論依據(jù)。
某油田主要產(chǎn)油層段儲層非均質(zhì)性較強(qiáng),單井產(chǎn)量普遍較低。為改善油氣開發(fā)效果,需要采用人工壓裂等措施。人工壓裂會使儲層及地層應(yīng)力方位發(fā)生改變,同時可能會影響非均質(zhì)儲層產(chǎn)生地層各向異性,從而影響偶極聲波測井?dāng)?shù)據(jù),無法呈現(xiàn)真實儲層信息。若能通過測井響應(yīng)趨勢特征判斷偶極聲波測井探測儲層是否具有各向異性,可以避免不必要的工程損失,提高測井效率。
國內(nèi)外巖石物理實驗證實,當(dāng)?shù)貙哟嬖趹?yīng)力時,沿著應(yīng)力方向的地層變得更致密,從而聲速增大,而垂直方向地層聲速則相應(yīng)減小。為此本文利用COMSOL有限元軟件建立圖1(a)所示的各向異性地層模型,并進(jìn)行網(wǎng)格剖分。將實際地層各向異性簡化為在相互正交的2個方向上存在聲速不同的地層,能夠以較為簡單的方式反映地層應(yīng)力產(chǎn)生的各向異性特性。在聲場有限元數(shù)值模擬中,網(wǎng)格對聲場的解析受聲場中聲波波長控制,可知聲波波長由波速與頻率共同決定,本文確定聲源頻率為4 kHz,即模型網(wǎng)格疏密程度受其聲速影響。圖1(a)中可觀察到模型網(wǎng)格由密到疏依次為井孔、低聲速地層以及高聲速地層,符合剖分規(guī)律。
圖1 正演模型剖分圖及xoz截面示意圖
圖1(b)給出了xoz截面充液井孔計算模型。模型井眼與儀器同軸,沿垂直井軸方向從里到外的介質(zhì)依次為儀器、泥漿濾液、地層,儀器源距(偶極聲源到第一個接收換能器之間的距離)為2.5 m;儀器間距(z方向上相鄰2個接收換能器陣列之間的距離)為0.15 m;井眼直徑為標(biāo)準(zhǔn)8 in(1 in=25.4 mm);數(shù)值模擬中不考慮隔聲體,用正負(fù)偶極子作為橫波聲源。
為了使激勵不易被噪聲干擾,選取時域平滑、頻譜信號分辨率高的高斯函數(shù)為偶極聲源函數(shù),在儀器外徑處加載2個震動相位相反的點聲源模擬偶極子聲源,高斯聲源函數(shù)時域表達(dá)式如下:
G(t)=-2π2f02(t-tp)e-π2f02(t-tp)2。
(1)
式中:t為時間;tp為初振相,tp=1/f0;f0為中心頻率,4 kHz。
偶極子聲波測井儀器收集4分量的偶極子聲波陣列信號,32個接收換能器對稱分布在三維坐標(biāo)系的X和Y方向。4分量波形數(shù)據(jù)中,存在2組同向分量(XX和YY)和2組交叉分量(XY和YX),本文中分量的第一個字母代表發(fā)射偶極子的方向,第二個字母代表接收器的方向。
各向異性地層使沿不同方向傳播的橫波方位及速度產(chǎn)生差異。偶極聲源激發(fā)的滑行橫波穿過井壁進(jìn)入各向異性地層,由于橫波分裂理論,將分為分別沿快慢波面以不同速度傳播的正交極化的快主橫波(FP)和慢主橫波(SP),在實際地層勘探中,偶極聲源偏振方向一般與快慢橫波面有一定夾角θ,如圖2所示。
圖2 橫波分裂:發(fā)射、接收橫波與快慢橫波分量關(guān)系
2.2.1 波形匹配反演算法
本文的各向異性反演方法是以快慢主波波形相似的數(shù)值模擬理論為依據(jù),考慮到4分量32個接收換能器接收到的陣列中快慢主波的各種可能組合建立反演目標(biāo)函數(shù),通過求其全局最優(yōu)解得到各向異性參數(shù)。反演算法采用了快慢波匹配、方位角、各向異性值聯(lián)合反演的方法,目的是能夠較快地得到地層的各向異性參數(shù)及方位角。此算法具有高精度、高效率的優(yōu)勢,通過多次迭代抑制和消除波形數(shù)據(jù)的噪聲干擾,最終求取全局最優(yōu)解。本方法成為現(xiàn)今通過偶極聲波測井響應(yīng)波形信號反演各向異性參數(shù)應(yīng)用最多的方法,其宗旨是通過優(yōu)化算法最小化目標(biāo)函數(shù),即
(2)
式中:T積分時窗取整個波形時間跨度;s1和s2分別為聲源到接收器陣列跨度內(nèi)快、慢橫波的慢度,慢度差δs=s2-s1;θ為方位角;m和n為決定進(jìn)行波形相似匹配的任意2個接收器位置;zm是聲源到第m個接收器之間的距離;δz是相鄰接收換能器之間的距離;N為接收換能器的總數(shù)(N=8);t為積分時間。
針對本文所研究的地層各向異性情況,依托目標(biāo)函數(shù)式(2)所表示的快、慢波形間的數(shù)學(xué)相似性原理對響應(yīng)波形進(jìn)行解析。如圖3所示,首先假設(shè)快、慢波相似,若慢波在時間上向前移動到快波位置與之匹配,兩者應(yīng)該基本重合。
圖3 波形匹配示意圖
圖3所示的是偶極聲波測井儀器4個收發(fā)聲系陣列的其中一個子陣列的排列情況,圖中參數(shù)分別為:T表示發(fā)射聲系;Rn、Rm分別為此子陣列接收聲系中的任意2個,n、m分別表示第幾個接收器;δtmm表示在Rm接收器處快慢橫波的到時時差,δtmn表示Rm接收的快橫波和Rn接收的慢橫波之間的到時時差。
第m個接收器處,波形匹配的表達(dá)式為
(3)
將第m個接收器與第n個接收器處波形進(jìn)行匹配
SPn(t+δtmn)=SPn(t+δszm+s2(n-m)δz)≈
FPm(t)
(4)
通過分析4分量32個接收換能器采集的波形信號,利用算法能夠提取出快橫波傳播方向以及基于偶極聲源偏振方向的偏移角θ和快慢橫波慢度差,其結(jié)果可以表征應(yīng)力方向以及各向異性強(qiáng)度。
2.2.2 能量各向異性算法
地層能量各向異性是表征地層各向異性強(qiáng)度的參數(shù)之一。此算法的內(nèi)核是直接計算交叉分量相對于4分量波形總能量的能量占比,由于能量各向異性結(jié)果受計算誤差的影響很小,其結(jié)果可靠性高,可以高精度地表征地層各向異性強(qiáng)度,同時可以校正經(jīng)波形匹配反演算法計算的陣列各向異性結(jié)果(即經(jīng)陣列波形結(jié)果進(jìn)行波形匹配計算的各向異性強(qiáng)度)。
若某方位角θ上接收到4分量偶極波形分別為XX,XY,YX,YY,定義一個波形的能量為波形各個采樣點值的絕對值之和
(5)
式中:f為波形函數(shù),N為波形取樣點數(shù),E為波形能量。
定義交叉分量與波形總能量的百分比,即交叉分量相對能量為
(6)
定義交叉分量相對能量最大值和最小值,分別為
(7)
定義能量各向異性A為
(8)
本文數(shù)值模擬主要研究基于地應(yīng)力的地層各向異性對偶極聲波測井全波波形的影響。快慢地層分別在應(yīng)力與垂直于應(yīng)力的方向。為清晰分析各向異性參數(shù)對模型的影響,定義引起地層各向異性方向與偶極聲源X偏振方向之間的夾角為θ,定義地層各向異性強(qiáng)度系數(shù)
K=vl/vh。
(9)
式中:vh和vl分別為造成地層橫波各向異性的快慢橫波波速。K值大小與地層各向異性強(qiáng)度呈反比。
數(shù)值模擬采用的井孔流體和各向異性地層快慢橫波聲速參數(shù)的確定既要符合本文正演研究假設(shè),同時也要結(jié)合實際的合理的測井參數(shù)??紤]本文研究的井旁各向異性情況,給定井孔內(nèi)流體和地層快慢聲學(xué)參數(shù),見表1。
表1 井孔流體及各向異性地層聲學(xué)參數(shù)
擬定采用的地層各向異性強(qiáng)度為K=0.73,各向異性方位角分別為θ=30°和90°。數(shù)值模擬結(jié)果如圖4、圖5所示。
圖4 K=0.73,θ=30°各向異性對偶極測井陣列波形的影響
圖5 K=0.73,θ=90°各向異性對偶極測井陣列波形的影響
波形幅度分析:
(1)交叉分量:對于各向同性地層,其交叉分量波形信號強(qiáng)度很低;當(dāng)?shù)貙哟嬖诟飨虍愋郧腋飨虍愋詮?qiáng)度K=0.73,方位角θ分別為30°時,各向異性地層在介質(zhì)分界面發(fā)生耦合作用,彎曲波發(fā)生分裂現(xiàn)象,致使交叉分量能量變高,陣列波形幅度變大。當(dāng)各向異性方位角θ=90°時,主應(yīng)力方向與Y發(fā)射方向一致,彎曲波分裂程度很低,其交叉分量波形能量與各向同性地層波形振幅在同一數(shù)量級。
(2)同向分量:各向異性強(qiáng)度相同條件下,2種不同各向異性方位角對同向分量波形的影響為首波波形振幅均衰減明顯,且全波波形相位均向時間零刻度方向發(fā)生偏移,但各向異性方位角變化對同相分量陣列波形影響不大。
FFT是離散傅里葉變換的簡稱。對聲波數(shù)據(jù)流式信號的分析通常分為2個方面,一方面是其在時間域上連續(xù)的波形特征分析,例如振幅以及振型到時和結(jié)束等;另一方面是頻譜分析,從頻域角度解析波形得到更多有用的信息,例如可以得到波形在一個時間區(qū)段上的不同頻率下的能量分布和幅值分布,進(jìn)行波形解析后可以得到所研究波形的主要能量及其幅度所在頻率范圍。對上節(jié)不同各向異性方位角下偶極聲波測井響應(yīng)波形進(jìn)行頻譜變換得到頻譜圖,如圖6—圖7所示。
圖6 θ=30°時頻譜圖
圖7 θ=90°時頻譜圖
波形能量頻譜分析:
(1)各向異性方位角θ=30°時:由圖6(a)(b)可知,同相分量波形能量在聲源中心頻率處波形幅度衰減最大,同向分量波形能量衰減接近70%,結(jié)果與同向分量波形衰減趨勢吻合,XX分量波形平均能量衰減接近43%,YY分量波形能量平均衰減近55%,表明各向異性方位角θ=30°造成橫波分裂后,XX分量的能量大于YY分量,結(jié)論與各向異性地層橫波分裂理論相符; 交叉分量波形能量受各向異性影響漲幅較大,波形平均能量增長近5倍。
(2)各向異性方位角θ=90°時:由圖7(a)(b)可知,當(dāng)主應(yīng)力方向與Y發(fā)射方向同向,橫向各向異性產(chǎn)生的快慢波面分別坐落于儀器的Y和X接收方向,此時橫向各向異性造成橫波分裂程度最低,同向分量波形能量衰減減少,XX分量平均能量衰減近41%,YY分量平均能量衰減近33%;交叉分量平均能量與各向同性地層中測得的波形的能量相近。
各向異性地層所表示出聲學(xué)特性的差異決定了各向異性強(qiáng)度。本節(jié)模擬分析3種各向異性方向情況下不同各向異性強(qiáng)度對偶極聲波測井全波波形的影響。計算得到的波形如圖8—圖9所示。
圖8 θ=30°時3種各向異性強(qiáng)度對波形影響趨勢
圖9 θ=90°時3種各向異性強(qiáng)度對波形影響趨勢
波形分析如下:
(1)同向分量:由圖8(a)、圖9(a)可知,各向異性地層對波形幅度影響顯著,同向分量波形幅度明顯減小,波形幅度減少程度隨著各向異性強(qiáng)度系數(shù)K的變小而變大。
(2)交叉分量:交叉分量是各向異性地層不同橫波波速地層分界面耦合的產(chǎn)物,通常分界面兩邊地層聲學(xué)性質(zhì)差異越大,其波形能量越強(qiáng),如圖8(b)所示。機(jī)理在于本文正演模型地層的特殊性,地層各向異性分界面始終處于交叉分量聲波從發(fā)射到接收的途徑中,聲波通過2種不同聲學(xué)特性介質(zhì)分界面時發(fā)生反射與折射,導(dǎo)致交叉分量波形能量發(fā)生衰減,且聲波衰減程度隨分界面兩邊地層橫波波速差異變大而變強(qiáng)。圖9(b)中,3種橫波各向異性強(qiáng)度造成交叉分量波形幅度很小且不隨系數(shù)K變化,結(jié)果產(chǎn)生了逆常理性現(xiàn)象,原因在于此時快慢波面分別坐落于偶極聲源偏振方向,造成彎曲波分裂程度很低,致使交叉分量波形振幅極低。而無各向異性的均勻地層則幾乎不存在交叉分量信息,如圖8(b),圖9(b)交叉分量第一道波形所示。
本文將有限元數(shù)值模擬正演得到的4分量陣列波形響應(yīng)數(shù)據(jù)作為式(2)的輸入,以各向異性地層中橫波分裂為具有相似特征的快、慢橫波為理論基礎(chǔ),使用模擬退火算法得到目標(biāo)函數(shù)在求解域中的全局最小值,反演出地層各向異性參數(shù)包括慢橫波的慢度s2、快慢橫波之間的差異δs以及最大應(yīng)力方位角θ,此時快橫波慢度s1=s2-δs。根據(jù)本文各向異性模型的正演計算陣列波形進(jìn)行參數(shù)反演,反演目的是為了驗證正演結(jié)果的正確性。 為了便于比較分析反演效果,選取預(yù)設(shè)模型參數(shù),見表2。
表2 參與反演過程的正演橫波各向異性參數(shù)
將各向異性方位角θ=30°、90°對應(yīng)的正演陣列波形結(jié)果作為反演對象,設(shè)置內(nèi)部迭代次數(shù)N=100進(jìn)行迭代計算。由于模擬退火算法每次迭代在全局中搜尋最小值,其目標(biāo)函數(shù)一般呈波動式逐步減小,缺點為收斂速度較慢;若給定合理反演參數(shù)范圍,可以加速反演收斂速度,同時得到正確反演結(jié)果。本次給定反演角度范圍為0~90°,得到反演結(jié)果見表3、表4。
表3 模型主應(yīng)力方向角θ=30°時反演結(jié)果對比
表4 模型主應(yīng)力方向角θ=90°時反演結(jié)果對比
反演結(jié)果分析如下:
(1)模擬退火優(yōu)化算法實現(xiàn)了偶極聲波測井模型各向異性參數(shù)的高精度反演,得到了與預(yù)設(shè)模型各向異性參數(shù)吻合的結(jié)果。
(2)各向異性參數(shù)變化的同時,對同一種各向異性方位角進(jìn)行多次反演,其結(jié)果與正演參數(shù)基本一致,經(jīng)計算反演相對誤差小于5%,在數(shù)值計算誤差允許范圍內(nèi)。
本文各向異性地層模型描述的地層狀態(tài)較為復(fù)雜,造成地層不同深度處的各向異性描述困難,通過波形匹配反演算法輸出快慢橫波慢度s1和s2,已知接收器陣列各向異性強(qiáng)度計算公式為2(s2-s1)/(s1+s2),計算主應(yīng)力方向角θ=30°、90°下2種工況下的陣列各向異性與能量各向異性,對比結(jié)果見表5。
表5 陣列各向異性與能量各向異性對比
由波形匹配算法思想可知,陣列各向異性的結(jié)果由算法多次迭代產(chǎn)生,算法迭代存在誤差的累積,其結(jié)果可能會發(fā)生偏差;而能量各向異性表征地層各向異性強(qiáng)度,因能量各向異性算法簡潔高效,故結(jié)果可靠性更高。將3種各向異性結(jié)果進(jìn)行對比,陣列各向異性結(jié)果略大于能量各向異性結(jié)果,數(shù)值上的差異可能是由于誤差迭代產(chǎn)生的,其結(jié)果可以正確表征地層的各向異性強(qiáng)度。
(1)地層各向異性方向角對偶極測井波形具有影響,地層各向異性使同相分量波形幅度及能量降低,各向異性方向角的改變對同相分量波形影響不顯著;地層出現(xiàn)一定角度異常應(yīng)力會使交叉分量波形振幅和能量產(chǎn)生明顯增長,特殊情況下譬如各向異性方向與聲源偏振方向一致,交叉分量波形振幅及能量達(dá)到最小。
(2)地層各向異性強(qiáng)度對偶極測井波形影響較大,同向分量波形振幅的變化程度隨地層各向異性強(qiáng)度系數(shù)變大而顯著減?。桓鶕?jù)聲波在介質(zhì)分界面的反射定律和折射定律解釋了交叉分量波形產(chǎn)生的逆常理性趨勢。
(3)精確反演地層各向異性參數(shù)驗證了正演模擬的正確性;相較于傳統(tǒng)的陣列各向異性反演方法,能量各向異性反演方法能夠更準(zhǔn)確地反映儲層各向異性強(qiáng)度。