黃國嬌　巴　晶　錢　衛(wèi)

(河海大學地球科學與工程學院地球探測研究所，江蘇南京 211100)

1　引言

水力壓裂技術是非常規(guī)油氣資源(如致密油氣、頁巖油氣)開發(fā)過程中重要的增產(chǎn)措施，該技術通過將流體高壓注入地層，使巖石破裂以提高地層滲透率[1]。微地震監(jiān)測技術通過觀測、定位及分析微地震事件監(jiān)測生產(chǎn)活動及地下狀態(tài)，可用于水力壓裂效果分析[2-6]，在儲層開發(fā)監(jiān)測中起到重要作用。目前，微震監(jiān)測技術應用的主要困難在于： ①水力壓裂誘導的微地震事件通常信噪比較低，且觀測數(shù)據(jù)量大，初至旅行時難拾?。?②即使能夠較準確地拾取初至旅行時(旅行時信息包含了地震事件的發(fā)震時刻和旅行時信息)，微地震事件的定位仍然存在誤差，原因是缺乏精確的地下速度模型。微地震反演中采用的速度模型一般根據(jù)地震、聲波測井、VSP等資料建立[6]，并采用壓裂前的射孔資料校正速度模型[7,8]。但是，射孔反演速度模型需已知射孔位置和發(fā)震時刻，而射孔事件的準確發(fā)震時刻往往是未知的。再者，射孔數(shù)據(jù)一般在壓裂前采集，而在水力壓裂過程中，高壓流體會產(chǎn)生裂縫，新增裂縫和孔隙壓力的改變都會引起地層速度的變化[9]。因此如果采用射孔資料校正后的速度模型進行微地震定位，可能存在較大誤差[10]。

針對微地震初至拾取困難的問題，Song等[11]建立了一個微地震強信號的波形數(shù)據(jù)庫，通過分析強弱微地震信號的波形相關性提取弱微震事件的初至； 宋維琪等[12]提出了獨立分量分析與壓縮感知聯(lián)合理論的微地震弱信號的提取方法； 王晨龍等[13]和李振春等[14]采用地震干涉逆時定位算法，不需要拾取初至信息； 王鵬等[15]提出基于熵值的微震事件識別和恢復方法，以解決微地震監(jiān)測數(shù)據(jù)信噪比低的問題； 李穩(wěn)等[16]將稀疏碼收縮法和小波變換相結合處理井下微地震數(shù)據(jù)，以保證有效地震信號的識別和提取； 曹俊海等[17]采用多道匹配跟蹤算法提高微地震信號的信噪比，進而提高微地震信號識別的準確性。

為得到較準確的速度模型，減小速度模型對定位結果的影響，Grechka等[18]同時反演了介質的各向異性參數(shù)和微地震位置，但假設地下為均勻介質； Li等[19]將雙差定位法推廣至同時確定水平層狀各向異性參數(shù)及微地震定位； 譚玉陽等[20]提出基于初至旅行時差的微地震速度模型反演方法，可在射孔起始時間未知的情況下反演地層速度模型。

在井下微地震監(jiān)測中，微地震事件與檢波器之間的距離通常在幾十至數(shù)百米的范圍內，因此將地下介質假設為水平層狀模型[7，21]，提出了地面與井中觀測旅行時聯(lián)合的微地震定位和速度結構同時反演方法： 首先推導了旅行時關于速度、界面深度以及震源參數(shù)的一階偏導數(shù)，給出了Jacobi矩陣計算方法； 然后用最短路徑算法[22]追蹤射線路徑和計算旅行時；再采用共軛梯度算法求解帶約束的阻尼最小二乘問題，更新地下層狀介質的結構(層速度和層界面深度)以及微地震定位(包括震源位置和發(fā)震時刻)； 最后討論了方法的抗噪聲能力。

2　方法原理

2.1　反演問題

一般而言，基于旅行時數(shù)據(jù)的非線性同時反演問題，可歸結為下列帶約束的阻尼最小二乘最優(yōu)化問題，其反演問題的L2范數(shù)目標函數(shù)為

(1)

式中： δt=(δt1，δt2，…，δtM)是微地震事件的實際觀測旅行時與理論數(shù)據(jù)的旅行時殘差，M則為拾取的初至數(shù)量；δm=(λ1δV，λ2δD，λ3δX)，其中δV和δD分別是地下介質各層速度和地層界面深度的相對改變量，δX為微地震事件的震源坐標和發(fā)震時刻的相對改變量，λ1、λ2、λ3為各類參數(shù)的更新權重系數(shù)；A是Jacobi矩陣；μ是阻尼因子；a和b分別是反演解中更新各類參數(shù)更新量的上、下邊界值。

式(1)解的一般形式[23]為

[μCm+ATCdA]Zmδm=ATCdδt

(2)

式中：Cm、Cd為模型和數(shù)據(jù)空間的協(xié)方差矩陣；Zm為分段約束算子，即

(3)

式(2)是非線性問題局部線性化的基本反演公式，其解具有局域解的特征。為了得到具有實際物理意義的解，可采用Cm、Cd及Zm的先驗信息進行約束。因此，選擇不同的組合(Cm、Cd、Zm)可得到不同形式的反演解。

本文采用Tarantola等[24]提出的廣義帶約束的阻尼最小二乘反演解，即Cm、Cd分別取模型和數(shù)據(jù)空間的協(xié)方差矩陣的逆，則解為

(4)

上述矩陣方程可基于共軛梯度法采用迭代算法求解[25]，其關鍵問題是如何求解具有偏導性質的Jacobi矩陣。

2.2　Jacobi矩陣元素計算

當?shù)叵陆橘|結構(層速度和層深度)和微地震事件同時反演時，Jacobi偏導矩陣包含了三部分：一是旅行時關于速度變化的導數(shù)；二是旅行時關于層深度變化的導數(shù)； 三是旅行時關于微地震震源參數(shù)的導數(shù)。

(5)

式中： Δti,j代表第i個震源第j條射線的旅行時擾動；vk和dk分別是第k層內的速度分布和第k層界面的深度； Δvk和Δdk則分別是相應速度和層界面深度的擾動量；N為射線穿過的層數(shù)；K是與射線相交的界面數(shù)； Δxi,k(k=1，2，3，4) 是第i個微地震事件分別在x、y、z坐標及發(fā)震時刻的擾動量。

式(5)中第一項，即旅行時關于速度的一階偏導數(shù)為

(6)

式中L為射線在第k層內的長度。

式(5)中的第二項，即旅行時關于層界面深度變化的偏導數(shù)參照Rawlison[26]的方法根據(jù)圖1分析并推導。當層界面的深度變化Δd時(從圖1中實線到虛線位置)，引起的旅行時變化為

(7)

則旅行時關于層界面深度變化的偏導數(shù)為

(8)

式中：向量wj、wj+1和wz如圖1所示；vj和vj+1則分別是層界面兩側的速度；e和f分別為界面擾動后界面上、下兩側射線長度的變化量。

式(5)中的第三項，即旅行時關于微地震事件參數(shù)(位置和發(fā)震時刻)變化的偏導數(shù)，可由以下步驟計算。首先，旅行時(t=ta-t0)關于發(fā)震時刻的導數(shù)為

(9)

其次，旅行時關于微地震事件位置變化的偏導數(shù)的計算，可參見圖2所示。若微地震震源的深度擾動值為Δzh，則旅行時關于震源位置xh、yh、zh的偏導數(shù)[26]為

(10)

式中：θh如圖2所示；vh為震源處的速度值；φh、ψh分別為射線水平面投影后與x和y軸的夾角。

圖2　旅行時關于震源位置變化偏導數(shù)計算示意圖

3　數(shù)值模擬實例

為驗證該算法的有效性及實用性，本文選用一個三維模型同時進行地下介質的層析成像及微地震事件的定位試驗。如圖3所示，模型尺寸為200m×200m×50m，網(wǎng)格尺寸為4m×4m×5m。模型是四層的水平層狀介質，速度從地表向下依次為2.0、2.4、2.8和3.2km/s，界面的深度分別為16、30和40m。為模擬地面和井中微地震，采用實際微地震中常用的地表“米”字型和井中聯(lián)合觀測方式(如圖3中三角形所示)，地表采集空間間隔為20m，井中采集空間間隔為10m，共46個檢波器。設定了4個虛擬微地震事件，其空間坐標分別為S1(20m，40m，42m)、S2(100m，100m，45m)、S3(150m，180m，48m)、S4(170m，30m，26m)，激發(fā)時刻分別為10、15、20和5ms。

圖3　三維模型及微地震觀測系統(tǒng)三角形代表觀測點，圓點代表微地震事件

在理想條件下(正確的速度模型、無噪聲)采用該反演算法進行微地震事件S1、S2和S3定位。定位結果如圖4所示，即使初始震源參數(shù)與真實值相差較大，反演的震源位置(圖4a～圖4c)和發(fā)震時刻(圖4d)也都與理論參數(shù)完全一致。

然而實際數(shù)據(jù)處理中不可能滿足理想條件，即地下介質的速度結構(層速度和層界面位置)未知、數(shù)據(jù)含噪聲等。下面給出本文反演算法同時進行速度模型反演和微地震定位的結果，并討論其對噪聲的敏感性。首先選用兩種不同的微地震事件組合測試該反演算法同時進行速度模型反演和微地震定位。第一種采用三個微地震事件(S1、S2和S3)的初至P波旅行時，第二種采用四個微地震事件(S1、S2、S3和S4)的初至P波旅行時，同時反演的結果如圖5(速度結構反演結果)和圖6(微地震定位結果)所示。

圖4　理想條件下的微地震定位結果

從圖5 可以看出，當僅采用三個微地震事件的旅行時數(shù)據(jù)時，第三層的速度結構(層界面和層速度)層析的結果(圖5中藍線)與真實模型(圖5中紅線)差別較大；而當采用四個微地震事件的旅行時數(shù)據(jù)時，層析結果明顯改善，無論是層速度值還是層界面位置(圖5中綠線)都接近真實值。出現(xiàn)這種情況是因為增加的微地震事件S4位于第二層介質中，與微地震事件S1、S2和S3位置有較大區(qū)別，射線覆蓋及射線交錯概率有所增加，提高了層析成像的分辨率。所以當速度結構反演結果得到改善時，相應的微地震定位結果精度也有所提高，特別是震源深度(圖6b和6c)的定位結果。因此，當射線的角度覆蓋較好時，同時進行速度層析成像和微地震定位是可行的，反演精度較高。

圖5　不同微地震事件組合同時反演的速度結構

圖6　不同微地震事件組合同時反演的微地震定位結果

考慮拾取初至旅行時可能存在誤差，在理論的初至P波旅行時數(shù)據(jù)中加入幅度為0.4ms的隨機誤差，四個事件同時反演的結果如圖7和圖8所示。結果顯示：加入隨機噪聲后，對微地震定位精度影響不大，且微地震震源水平位置的定位結果(圖8a)明顯好于震源深度的定位(圖8b和圖8c)。表1給出了具體的定位誤差，可見未加隨機噪聲之前水平定位誤差在幾厘米的范圍內，加入隨機噪聲后水平定位誤差在十厘米左右，而深度的定位誤差為1～2m，但仍在可接受的誤差范圍內。速度結構層析的結果變化較大(圖7)，層速度和層界面深度都存在一定程度的偏差。這是由于當速度結構層析和微地震定位同時進行時，兩者的反演結果實際上是一定程度的折中，并且層速度和層界面的深度在反演時也存在平衡的問題。總體而言，微地震定位對旅行時中的隨機噪聲不敏感。

圖7　噪聲敏感性試驗速度層析結果

x誤差/my誤差/mz誤差/mt0誤差/ms無噪加噪無噪加噪無噪加噪無噪加噪S10．002-0．013-0．026-0．1301．3991．975-0．077-0．109S2-0．0010．1680．0010．0500．8240．341-0．053-0．002S30．0140．118-0．0640．0240．8541．014-0．055-0．040S4-0．092-0．0620．1430．145-0．1690．233-0．103-0．227

圖8　噪聲敏感性試驗微地震定位結果

上述的地面和井中聯(lián)合反演結果均迭代了20次左右。

4　結束語

為降低速度模型估計不準對微地震定位精度的影響，本文提出了對水平層狀介質同時進行速度層析成像和微地震定位的反演算法。首先采用最短路徑算法計算初至波的射線路徑和旅行時，導出旅行時關于層速度、層界面深度及震源參數(shù)的偏導數(shù)，然后結合共軛梯度算法求解帶約束的阻尼最小二乘問題，發(fā)展了一種同時反演速度結構(層速度和層界面)和震源參數(shù)(震源位置和發(fā)震時刻)的方法。地面和井中聯(lián)合觀測條件下的數(shù)值模擬結果表明：該反演算法能有效反演地下介質的速度結構并進行微地震定位，而且微地震定位結果對隨機噪聲不敏感。

本文基于模型正演無噪和加入隨機噪聲的數(shù)據(jù)，驗證了所提出方法的有效性。在實際生產(chǎn)中，還存在地層起伏不平、介質非均質性強、數(shù)據(jù)信噪比低等諸多因素，因此還有必要結合實際數(shù)據(jù)與現(xiàn)場地質情況，應用、測試方法的可靠性和穩(wěn)健性，提高實際定位精度。

[1]張光明.水平井水力壓裂數(shù)值模擬研究[學位論文].安徽合肥：中國科學技術大學，2010.

Zhang Guangming.A Numerical Simulation Study on Hydraulic Fracturing of Horizontal Well[D].University of Science and Technology of China,Hefei,Anhui,2010.

[2]Rutledge J T and Phillips W S.Hydraulic stimulation of natural fractures as revealed by induced micro-earthquakes,Carthage Cotton Valley gas field,east Texas.Geophysics,2003,68(2)：441-452.

[3]Maxwell S C,Rutledge J,Jones R et al.Petroleum reservoir characterization using downhole microseismic monitoring.Geophysics,2010,75(5)：A129-A137.

[4]容嬌君,李彥鵬,徐剛等.微地震裂縫檢測技術應用實例.石油地球物理勘探,2015,50(5)：919-924.

Rong Jiaojun,Li Yanpeng,Xu Gang et al.Fracture detection with microseismic.OGP,2015,50(5)：919-924.

[5]Li J L,Kuleli S,Zhang H J et al.Focal mechanism determination of induced microearthquakes in an oil field using full waveforms from shallow and deep seismic networks.Geophysics,2011,76(6)：WC87-WC101.

[6]張山,劉清林,趙群等.微地震監(jiān)測技術在油田開發(fā)中的應用.石油物探,2002,41(2)：226-231.

Zhang Shan,Liu Qinglin,Zhao Qun et al.Application of microseismic monitoring technology in development of oil field.GPP,2002,41(2)：226-231.

[7]Warpinski N R,Sullivan R B,Uhi J et al.Improved microseismic fracture mapping using perforation timing measurements for velocity calibration.Proceedings of the SPE Annual Technical Conference and Exhibition,Denver,Colorado,2003,84488.

[8]Van Dok R,Fuller B,Engelbrechet L et al.Seismic anisotropy in microseismic event location analysis.The Leading Edge,2011,30(7)：766-770.

[9]Block L V,Cheng C H,Fehler M C et al.Seismic imaging using microearthquakes induced by hydraulic fracturing.Geophysics,1994,59(1)：102-112.

[10]張曉林,張峰,李向陽等.水力壓裂對速度場及微地震定位的影響.地球物理學報,2013,56(10)：3552-3560.

Zhang Xiaolin,Zhang Feng,Li Xiangyang et al.The influence of hydraulic fracturing on velocity and microseismic location.Chinese Journal of Geophysics,2013,56(10)：3552-3560.

[11]Song F X,Kuleli H S,Toksoz M N et al.An improved method for hydrofracture-induced microseismic event detection and phase picking.Geophysics,2010,75(6)：A47-A52.

[12]宋維琪,李艷清,劉磊.獨立分量分析與壓縮感知微地震弱信號提取方法.石油地球物理勘探,2017,52(5)： 984-989,1041.

Song Weiqi,Li Yanqing,Liu Lei.Microseismic weak signal extraction based on the independent component analysis and compressive sensing.OGP,2017,52(5)：984-989,1041.

[13]王晨龍,程玖兵,尹陳等.地面與井中觀測條件下的微地震干涉逆時定位算法.地球物理學報,2013,56(9)：3184-3196.

Wang Chenlong,Cheng Jiubing,Yin Chen et al.Microseismic events location of surface and borehole observation with reverse-time focusing using interfero-metry technique.Chinese Journal of Geophysics,2013,56(9)：3184-3196.

[14]李振春,盛冠群,王維波等.井地聯(lián)合觀測多分量微地震逆時干涉定位.石油地球物理勘探,2014,49(4)：661-666,671.

Li Zhenchun,Sheng Guanqun,Wang Weibo et al.Time-reverse microseismic hypocenter location with interferometric imaging condition based on surface and downhole multi-components.OGP,2014,49(4)：661-666,671.

[15]王鵬,常旭,王一博.基于時頻稀疏性分析法的低信噪比微地震事件識別與恢復.地球物理學報,2014,57(8)：2656-2665.

Wang Peng,Chang Xu,Wang Yibo.Automatic event detection and event recovery in low SNR microseismic signals based on time-frequency sparseness.Chinese Journal of Geophysics,2014,57(8)：2656-2665.

[16]李穩(wěn),劉伊克,劉保金.基于稀疏分布特征下的井下微地震信號識別與提取方法.地球物理學報,2016,59(10)：3869-3882.

Li Wen,Liu Yike,Liu Baojin.Downhole microseismic signal recognition and extraction based on sparse distribution features.Chinese Journal of Geophysics,2016,59(10)：3869-3882.

[17]曹俊海,顧漢明,尚新民.基于局部相關譜約束的多道匹配追蹤算法識別微地震信號.石油地球物理勘探,2017,52(4)：704-714.

Cao Junhai,Gu Hanming,Shang Xinmin.Microseismic signal identification with multichannel matching pursuit based on local coherence spectrum constraint.OGP,2017,52(4)：704-714.

[18]Grechka V,Singh P and Das I.Estimation of effective anisotropy simultaneously with locations of microseismic events.Geophysics,2011,76(6)：WC141-WC155.

[19]Li J L,Zhang H J,Rodi W L et al.Joint microseismic location and anisotropic tomography using differential arrival times and differential back azimuths.Geophy-sical Journal International,2013,195(3)：1917-1931.

[20]譚玉陽,何川,張洪亮.基于初至旅行時差的微地震速度模型反演.石油地球物理勘探,2015,50(1)：54-60.

Tan Yuyang,He Chuan,Zhang Hongliang.Microseismic velocity model inversion based on moveouts of first arrivals.OGP,2015,50(1)：54-60.

[21]Rutledge J T,Phillips W S and Schuessler B K.Reservoir characterization using oil-production-induced microseismicity,Clinton County,Kentucky.Tectonophysics,1998,289(1)：129-152.

[22]Bai C Y,Tang X P and Zhao R.2D/3D multiply transmitted,converted and reflected arrivals in complex layered media with the modified shortest path me-thod.Geophysical Journal International,2009,179(1)：201-214.

[23]Menke W.Geophysical Data Analysis：Discrete In-verse Theory.Academic Press Inc.,San Diego,California,1984,40-126.

[24]Tarantola A,Valette B.Generalized non-linear inverse problems solved using the least squares criterion.Reviews of Geophysics Space Physics,1982,20(2)：219-232.

[25]Zhou B,Sinadinovski C,Greenhalgh S A.Iterative algorithm for the damped minimum norm,least-squares and constrained problem in seismic tomography.Exploration Geophysics,1992,23(3)：497-505.

[26]Rawlison N and Sambridge M S.Seismic traveltime tomography of the crust and lithosphere.Advances in Geophysics,2003,46：81-197.