戴峰， 郭亮， 徐奴文,2*， 樊義林, 徐劍， 姜鵬

1 四川大學(xué) 水力學(xué)與山區(qū)河流開發(fā)保護國家重點實驗室， 成都 610065 2 山東大學(xué) 土建與水利學(xué)院， 濟南 250061 3 中國長江三峽集團公司， 北京 100038

?

基于異向波速模型的微震定位改進

戴峰1， 郭亮1， 徐奴文1,2*， 樊義林3, 徐劍3， 姜鵬1

1 四川大學(xué) 水力學(xué)與山區(qū)河流開發(fā)保護國家重點實驗室， 成都 610065 2 山東大學(xué) 土建與水利學(xué)院， 濟南 250061 3 中國長江三峽集團公司， 北京 100038

針對波速分層的區(qū)域巖體，在異向波速模型的基礎(chǔ)上，對垂向上的應(yīng)力波按巖體波速值大小作分段區(qū)別，推導(dǎo)震源應(yīng)力波走時關(guān)系式，建立分層速度定位目標函數(shù)，基于此提出一種由參數(shù)準備、層速度反演、微震定位三個模塊組成的分層速度定位模型SV，并采用遺傳算法進行優(yōu)化求解.然后，對分層速度定位模型在已構(gòu)建微震監(jiān)測系統(tǒng)的白鶴灘水電站左岸巖質(zhì)邊坡進行驗證.微震事件重定位結(jié)果表明，分層速度定位模型定位微震事件的最大、最小和平均偏離層內(nèi)錯動帶程度指標較單一速度模型分別降低了57.17%、36.51%和57.35%，證明了定位模型在波速分層的區(qū)域巖體微震定位應(yīng)用中比單一速度定位模型更加合理可靠.

微震監(jiān)測； 異向波速模型； 分層速度定位模型； 遺傳算法； 巖質(zhì)邊坡

The SV model is proposed based on the anisotropic velocity model, which simplifies the propagation path of stress waves into a straight line, which consists of horizontal and vertical velocities, that is, one for underground sensors and the other for surface sensors. Taking into account that the vertical stress waves pass through some rock formations at different velocities, the model is constructed part by part according to the values of each interval velocity. The first thing is to derive the nonlinear equation of stress wave travel time by spatial analytic geometry′s operation and analysis. The target function of SV is established by using L2 norm of the time residual between the measured value and calculated value of the travel time difference of each two sensors. Secondly, based on the established target function, the SV is proposed, which consists of parameters-preparation module, interval velocities inversion module and seismic location module. In the parameters-preparation module, the regional rock mass is divided into some wave-velocity layers according to engineering geological data. And, using spatial analytic geometry′s operation and analysis, the analytical expressions for layer interfaces are derived. In the interval velocities inversion module, the target function to solve velocities is figured out according to the form of the target function of SV. Then we collect some blast events or rockburst with known coordinate values. Using the genetic algorithm (GA), we can get optimization solutions of velocities. In the seismic location module, using the determined velocities, the target functions to locate MS events are figured out according to the form of the target function of SV as well. It should be emphasized that during the optimal solving, for every wave-velocity layer, the concrete forms of target function and searching regions are different from each other. Using GA, the optimization solutions of a MS event for each wave-velocity can be worked out, and the minimum of optimization solutions is considered to be the location of one MS event.

In order to verify the effectiveness and reliability of SV, we applied it to a slope engineering on the left bank of the Baihetan hydropower station, where a MS monitoring system (ESG, from Canada) has been installed. Firstly, the slope was divided into 11 wave-velocity layers, which included 4 kinds of velocities, according to the sonic wave testing given by engineering geological investigation data. And, the analytical expressions for layer interfaces were worked out. Secondly, we collected 4 blast events that occurred at 610 m-elevation irrigation and drainage galleries. Using GA, we got optimization solutions of velocities, which are 4500 m·s-1, 4900 m·s-1, 3750 m·s-1, 5000 m·s-1, respectively, through minimizing the target function. And then, 11 target functions to locate MS events for each wave-velocity were figured out by using the determined velocities and parameters of layer interfaces. The MS events were collected from 1 to 18, April, 2015. The locations of these MS events were determined by using the simplified single-velocity model, of which the velocity is 4600 m·s-1(approximate mean value of 4 velocities above). Finally, those MS events, which were scattered around two dislocation bands (LS331, LS337), were relocated by using the SV. The relocation results show that the maximum, minimum and average deviations from dislocation bands decrease by 57.17%, 36.51% and 57.35%, respectively. It is clear that SV is more convenient and reliable than the single velocity location model aimed at those regional rock masses with stratified velocities.

The improved location model (SV), which is based on the anisotropic velocity model, takes into account the impact from the different velocities of rock formations with different lithological characters. Compared with the single velocity model, the accuracy of MS events location aimed at the regional rock masses with stratified velocities is improved obviously. However, like all MS events location models, it requires an accurate velocity-distribution of regional rock masses. Further study on improving the accuracy of MS events location can be achieved through considering velocity partitions inside each rock formation and finding a stable location model, which is more suited to the calculation of the actual stress wave traveltimes.

1 引言

巖石在外界應(yīng)力作用下，其內(nèi)部將產(chǎn)生局部彈塑性能集中現(xiàn)象，當(dāng)能量積聚到某一臨界值之后，就會引起巖體微裂隙的產(chǎn)生與擴展，微裂隙的產(chǎn)生與擴展伴隨著彈性波或者應(yīng)力波的釋放并在周圍巖體內(nèi)快速傳播，這種彈性波在地質(zhì)上稱為微震(microseismic/acoustic emission)(Xu et al.，2011).微震現(xiàn)象是由美國礦業(yè)局Obert(1977)發(fā)現(xiàn)，后來Kaiser(1950)對該現(xiàn)象命名為Kaiser效應(yīng).微震監(jiān)測技術(shù)是通過在監(jiān)測區(qū)域巖體布置檢波器或者傳感器來收集微震信號，然后分析該信號以圈定巖體損傷區(qū)域，從而為區(qū)域巖體穩(wěn)定性評價提供有效依據(jù).如今，微震監(jiān)測分析技術(shù)作為工程建設(shè)安全預(yù)測預(yù)警方法被廣泛應(yīng)用到礦山開采(姜福興等，2006；陸菜平等，2010)、石油開采(宋維琪等，2015)、水利水電工程(Xu et al.，2011；張伯虎等，2012)等行業(yè).

為了使微震監(jiān)測技術(shù)更好地發(fā)揮預(yù)測預(yù)警作用，滿足工程需求，定位精度的提高至關(guān)重要.微震定位問題屬于地球物理反演問題，它根據(jù)一次微震事件在各傳感器的震相到時來反演微震的發(fā)生時刻和空間位置(Spencer and Gubbins，1980；Pavlis and Booker，1980；劉福田，1984；劉杰等，2003；高永濤等，2015).但是，在實際工程應(yīng)用中，由于區(qū)域巖體是非均勻各向異性介質(zhì)，同時巖體中存在不同走向和傾角的節(jié)理斷層，微震震源發(fā)出的應(yīng)力波在傳播過程中能量、頻率和波速改變，這都給微震定位精度帶來了較大的影響.為了實現(xiàn)高精度微震定位，需要采用合理的傳感器空間布置方案、有效的微震波形識別技術(shù)以及合適的微震定位方法.

微震定位方法主要有幾何作圖法、相對定位法、空間域定位法、線性定位法和非線性定位法.幾何作圖法是一種直觀的定位方法，依據(jù)走時方程，通過二維或者三維作圖確定震源位置(廉超等，2006).相對定位法又稱為主事件定位法，是根據(jù)一個震源坐標較精確的主事件來計算其附近其他微震事件的坐標(Spence，1980).空間域定位法是由距離殘差代替到時殘差的一種定位方法.微震震相到時和震源參數(shù)的關(guān)系是一種非線性關(guān)系.線性定位法是通過求導(dǎo)、級數(shù)展開方法將該非線性問題轉(zhuǎn)化為線性問題的定位方法，例如Inglada(1928)法和Usbm(1970)法.這類方法在到時拾取精度不高的情況下定位效果較差.非線性定位法是利用迭代或優(yōu)化求解方式的定位方法統(tǒng)稱，例如Thurber(1985)用非線性牛頓迭代法進行定位.隨著計算機計算能力不斷發(fā)展，不同的非線性優(yōu)化算法工具被應(yīng)用到定位中，用于算法優(yōu)化的目標函數(shù)也出現(xiàn)了不同的形式.趙珠和曾融生(1987)在定位中采用了單純形法.但是單純形法采用單向搜索，形式簡單，容易陷入局部最優(yōu)解.針對這一缺點，全局最優(yōu)化算法被應(yīng)用到定位中，王洪體等(2006)介紹了一種基于浮點遺傳算法的定位方法，具有穩(wěn)定收斂、精度高、時效性強的特點；張華等(2014)基于遺傳算法對廠礦區(qū)內(nèi)的人工爆破進行定位，得到了很好的效果；針對速度模型給不準和聯(lián)合定位法震源位置、發(fā)震時刻和微震傳播速度相互關(guān)聯(lián)，解不唯一的問題，陳炳瑞等(2009)提出微震震源定位分層處理方法，并采用粒子群算法進行優(yōu)化求解；董隴軍等(2011)針對單一速度模型提出了無需預(yù)先測速的微震震源定位方法，利用非線性擬合工具進行求解；李健等(2014)針對單一速度模型，提出了無需測速的定位方法，采用單純形法進行最優(yōu)化求解.上述定位算法多是建立在單一速度模型的基礎(chǔ)上.雖然實際巖體中不同區(qū)域和方向上的波速值不盡相同，但是單一速度模型因不考慮太多未知量和關(guān)系，計算穩(wěn)定性強，是最常見的一種速度模型.

近幾年，一些專家和學(xué)者在尋找適應(yīng)各種區(qū)域巖體微震定位的速度模型上做出了不同探索.一般情形下，某個微震定位算法速度模型的優(yōu)劣取決于兩點：(1)是否可以較好地考慮區(qū)域巖體的波速特性，且能夠給出恰當(dāng)?shù)膽?yīng)力波傳播路徑；(2)是否簡潔，從而利于算法穩(wěn)定實現(xiàn)(Feng et al.，2015).射線追蹤理論能夠找到復(fù)雜巖體中應(yīng)力波傳播路徑，從而逐漸被應(yīng)用到波速差異大、帶空洞的區(qū)域巖體微震定位中(Bai and Greenhalgh，2006；Bai et al.，2010；Collins et al.，2014；Hosseini et al.，2015).但是，復(fù)雜巖體的速度模型建立及射線追蹤的實現(xiàn)對計算機要求相當(dāng)高，計算耗時大且不利于算法穩(wěn)定地實現(xiàn).為此，一些專家和學(xué)者針對不同波速分布特征的區(qū)域巖體分別對應(yīng)力波傳播特性作了相應(yīng)簡化處理，在提高定位精度的同時保證了計算的穩(wěn)定性.其中，F(xiàn)eng等(2015)針對隧洞開挖工程中圍巖波速特征，提出按隧洞斷面對傳感器進行波速值分組的定位速度模型，考慮了開挖隧洞中各斷面波速值差異，同時采用粒子群算法對所建立的目標函數(shù)實現(xiàn)穩(wěn)定求解；鞏思園等(2012)針對煤礦上覆巖層多以層狀形式賦存，P波在垂向和水平向的速度及傳播路徑差異大的情況，提出建立“異向波速模型”，模型求解選用遺傳算法與CMAES (Covariance Matrix Adaptation Evolution Strategy，協(xié)方差矩陣自適應(yīng)進化策略)算法結(jié)合的混合算法，經(jīng)驗證定位效果優(yōu)于單一速度模型.

本文綜合分析各種非線性定位方法，在“異向波速模型”基礎(chǔ)之上，將垂向上的應(yīng)力波按各巖層波速值大小進行分段處理，提出一種簡化的分層速度定位模型SV.該定位模型在考慮了不同巖層應(yīng)力波波速值大小不同的基礎(chǔ)上，對復(fù)雜的層間應(yīng)力波傳播路徑進行了簡化處理，采用具有全局尋優(yōu)特性的遺傳算法對建立的分層速度目標函數(shù)穩(wěn)定求解.相比單一速度定位模型，在波速分層的區(qū)域巖體微震定位中有利于定位精度的提高.

2 傳統(tǒng)單一速度定位模型

微震發(fā)生時向周圍巖體釋放縱波(P波)和橫波(S波)，前者傳播速度比后者快.當(dāng)傳播距離足夠遠時，在波形圖中能夠?qū)波與S波起跳時刻區(qū)分開來.然而，對于微震監(jiān)測區(qū)域來說，由于監(jiān)測區(qū)域不大，往往S波的起跳點疊加在P波中，要精確拾取S波初至?xí)r刻很困難.因此，通常僅拾取P波初至?xí)r刻用于定位.

對于單一速度模型，第i個傳感器接收P波初至?xí)r刻ti與震源參數(shù)(x0，y0，z0，t0)之間的關(guān)系由(1)式表述：

(i=1，2，…，m)

(1)

式中：li為震源與第i個傳感器之間的空間距離；(xi，yi，zi)為第i傳感器空間坐標，m是監(jiān)測區(qū)域內(nèi)接收到P波的傳感器數(shù)量；V為監(jiān)測區(qū)域P波波速值.

理論上，微震事件發(fā)生時刻與走時之和t0+Δti應(yīng)同ti相等，即

(2)

但是，在實際工程應(yīng)用中，由于受到監(jiān)測儀器以及人為因素的影響，微震信號起跳時刻的拾取存在誤差，則ξi值不為零.在這樣的情況下，為了實現(xiàn)微震震源定位，需使用非線性優(yōu)化工具進行最優(yōu)化求解來逼近真值.通常單一速度定位模型采用的非線性優(yōu)化目標函數(shù)有如下三種形式：

因變量為到時的目標函數(shù)：

(3a)

因變量為到時差的目標函數(shù)：

(3b)

因變量為到時差商的目標函數(shù)：

或

(3c)

運用非線性優(yōu)化方法，得到使上述目標函數(shù)盡量小的震源參數(shù)逼近解作為微震震源參數(shù)的最終計算結(jié)果.其中，第一種目標函數(shù)(3a)能將震源發(fā)生時刻一起反演；第二、三種目標函數(shù)(3b)、(3c)不能反演震源發(fā)生時刻；第一、二種目標函數(shù)(3a)、(3b)含有巖體波速參數(shù)，在優(yōu)化求解過程中可以與震源參數(shù)一起反演，從而實現(xiàn)均勻波速模型下無需預(yù)先測速的微震定位(董隴軍等，2011)；第三種目標函數(shù)(3c)不含巖體波速參數(shù)，從而能夠?qū)崿F(xiàn)均勻波速模型下的無需測速微震定位(李健等，2014).

3 異向波速模型的改進

在實際工程中，微震監(jiān)測區(qū)域巖體并不是均勻介質(zhì)，不同區(qū)域的應(yīng)力波波速不同.因此，應(yīng)力波射線是按最小走時路徑傳播，而并非幾何最短路徑.在許多工程中，巖體呈層狀分布，例如層狀礦層、層狀邊坡.

在層狀煤礦微震監(jiān)測中，震源應(yīng)力波傳播到礦區(qū)地面的檢波器比到井下巷道經(jīng)過的巖體介質(zhì)更為復(fù)雜，導(dǎo)致P波在垂向上速度變化很大.針對煤礦中這種層狀賦存和離層帶的特點，鞏思園等(2012)提出對安裝在井下巷道和礦區(qū)地面的檢波器進行區(qū)別對待，采用兩種速度建立異向波速模型：針對井下的檢波器，巖體波速較為均勻，P波傳播最小走時路徑近似為最短路徑傳播，速度為Vi=Vpu；針對礦區(qū)地表的檢波器，將實際最小走時路徑(折線)簡化為最短路徑，速度為Vi=Vpg.本文在此異向波速模型基礎(chǔ)之上，將垂向上的應(yīng)力波按各巖層波速值大小進行分段處理，提出一種簡化的分層速度定位模型SV，以期使微震定位精度得到提升.

3.1 應(yīng)力波走時關(guān)系式推導(dǎo)

如圖1所示，假定微震監(jiān)測區(qū)域巖體由四種不同波速的平行巖層組成，自下而上依次編號為①、②、③和④類，對應(yīng)的層速度值分別表示為V1、V2、V3和V4.實際地震波在不同介質(zhì)間傳播時是遵從費曼原理按折線傳播，鞏思園等(2012)將礦區(qū)地表檢波器接收到微震應(yīng)力波(層面垂向上的應(yīng)力波)的傳播路徑近似為直線，并將垂向上應(yīng)力波波速值視作區(qū)別于水平向的均值.本文考慮到垂向應(yīng)力波路徑穿過的各巖層波速值不同這一特征，將對應(yīng)在各層巖體中的走時由相應(yīng)的空間距離以及波速值表示出來，然后求和得到應(yīng)力波總走時.具體地，以3號傳感器接收到的在①類波速巖層發(fā)生的微震事件應(yīng)力波為例，推導(dǎo)其走時關(guān)系式如下：

通過震源與編號3傳感器的空間直線PS可用參數(shù)式方程組表達為：

(4)

式中：k為直線參數(shù)；(x0，y0，z0)為震源(P)的空間坐標；(x3，y3，z3)為編號3傳感器(S)的空間坐標.

由于①類和②類波速巖體交界面、②類和③類波速巖體交界面以及③類和④類波速巖體交界面互相平行，其空間平面一般解析式分別表示為：

(5)

式中：A，B，C是巖體交界面法向向量的三個分量值.D1,2、D2,3、D3,4分別為相應(yīng)交界面一般式方程中常數(shù)項.

經(jīng)過推導(dǎo)，編號3傳感器接收到發(fā)生在①類波速巖層中的微震事件應(yīng)力波走時Δt3：

(6)式中，α3和β3定義為空間距離參數(shù)，由(7)式計算得到：

(7)

在實際工程中，為了達到良好的監(jiān)測效果，監(jiān)測網(wǎng)絡(luò)呈網(wǎng)狀布置，各傳感器位于不同波速的巖層中.這樣，同一個微震事件的應(yīng)力波傳播到各個傳感器所經(jīng)過的巖層個數(shù)和巖層類別不一樣.由前文可以類比得到，某微震事件應(yīng)力波傳播到編號i傳感器的走時Δti可以由(8)式表示為：

(8)

其中，空間距離參數(shù)αi和βi為：

(9)

圖1 速度分層巖體中應(yīng)力波傳播特征Fig.1 Propagation characteristics of stress waves in rock mass with stratified velocities

式(8)和(9)中：從震源發(fā)出的應(yīng)力波傳播到編號i傳感器所經(jīng)過的波速層依次編號為1，2，…，n；層速度依次為Vi1，Vi2，…，Vin；應(yīng)力波所穿過不同波速層交界面的數(shù)學(xué)解析式中常數(shù)項依次為Di(1,2)，Di(2,3)，…，Di(n-1,n).

(10)3.2 目標函數(shù)建立

基于前文速度分層區(qū)域巖體的應(yīng)力波走時關(guān)系，可以構(gòu)建因變量為到時和因變量為到時差兩種形式的分層速度定位目標函數(shù)如下：

因變量為到時的分層速度定位目標函數(shù)：

(11a)

因變量為到時差的分層速度定位目標函數(shù)：

(11b)

上述等式右邊含有參數(shù)：傳感器空間坐標(xi，yi，zi)、分層模型空間幾何參數(shù)(A，B，C，D)、微震信號初至?xí)r刻ti、層速度值Vi、微震震源參數(shù)(x0，y0，z0，t0).本文采用因變量為到時差的目標函數(shù)(11b).

通過上述推導(dǎo)，可以看出，如果微震震源發(fā)生在不同的巖層，那么對應(yīng)目標函數(shù)的具體表達形式也不同.因此，在利用遺傳算法求解微震空間坐標的過程中，需要分波速層進行搜索，具體波速層η對應(yīng)的搜索可行域Ωη按如下線性規(guī)劃不等式加以限定：

(12)

其中：η為波速層編號；Aη，Bη，Cη，Dη為邊界面解析式方程參數(shù)；γ為邊界面?zhèn)€數(shù).

3.3 定位模型構(gòu)建

在分層速度定位模型中，采用遺傳算法進行分層速度定位目標函數(shù)最優(yōu)化求解，屬于目標函數(shù)最小化問題.但是，在遺傳算法中，種群是向適應(yīng)度值高的方向進化.因此，需要建立目標函數(shù)與適應(yīng)度函數(shù)值之間的映射關(guān)系.本文采用的適應(yīng)度函數(shù)H與目標函數(shù)f的映射關(guān)系如下：

(13)

式中：cmax是當(dāng)前所有或最近k代種群中f的最大值，并且cmax隨遺傳代數(shù)不斷變化.圖2是目標函數(shù)的遺傳算法優(yōu)化流程.在優(yōu)化求解過程中，當(dāng)遺傳算法達到以下兩個條件之一則終止算法：達到最大遺傳代數(shù)N(條件A)；遺傳代種群距離小于閾值Δd(條件B).對每個波速層求解時，如果因達到條件A而終止算法，則認為該波速層的震源優(yōu)化解不存在；如果因達到條件B而終止，則認為該波速層的解存在，同時如果終止種群的最優(yōu)個體適應(yīng)度值和平均適應(yīng)度值分別大于其閾值H和h，則認為解是可行的.終止種群中的最優(yōu)個體作為該波速層的震源優(yōu)化解，同時記錄其對應(yīng)的目標函數(shù)值.從解存在的各個波速層中，選取目標函數(shù)值最小者對應(yīng)的震源優(yōu)化解作為最終定位結(jié)果，同時其目標函數(shù)值作為誤差值來評價該定位結(jié)果的可靠性.

圖2 遺傳算法非線性目標函數(shù)優(yōu)化流程圖Fig.2 Flow diagram of the nonlinear objective function optimization using GA

在此基礎(chǔ)上，如圖3所示，本文構(gòu)建分層速度定位模型分為三個模塊：參數(shù)準備、層速度反演、微震定位：

(1) 模塊Ⅰ：參數(shù)準備.確定監(jiān)測區(qū)域巖體范圍，建立坐標系和分層速度幾何模型.然后劃定波速層，并計算不同波速層交界面解析式方程；確定各巖層遺傳算法搜索可行域Ω；根據(jù)監(jiān)測區(qū)域巖體范圍，優(yōu)化布置傳感器網(wǎng)絡(luò)，得到傳感器坐標；綜合計算速度和精度要求，設(shè)定遺傳算法優(yōu)化終止參數(shù).

圖3 分層速度定位模型流程圖Fig.3 Flow chart of stratified velocity location model

(2) 模塊Ⅱ：層速度反演.收集現(xiàn)場巖爆或施工爆破事件數(shù)據(jù)，拾取相應(yīng)各通道傳感器采集到的波形起跳時刻.根據(jù)巖爆或者施工爆破事件空間坐標(xb，yb，zb)判斷其所在波速層編號為k，由巖爆或者施工爆破事件空間坐標、分層模型邊界面解析式參數(shù)和傳感器坐標，按公式(11b)生成巖層k(波速層k)目標函數(shù)并采用遺傳算法計算層速度組合優(yōu)化解，將其用于模塊Ⅲ微震定位.

(3) 模塊Ⅲ：微震定位.微震系統(tǒng)采集現(xiàn)場微震事件波形信號，拾取波形起跳時刻.由分層模型邊界面解析式參數(shù)、傳感器坐標和模塊Ⅱ中求解得到的層速度組合值，按公式(11b)依次生成各巖層(波速層)目標函數(shù)并采用遺傳算法在可行域Ωη中搜索得到相應(yīng)微震震源參數(shù)優(yōu)化解(xη0，yη0，zη0)及相應(yīng)目標函數(shù)值fη.在滿足解存在條件的各波速層中選取目標函數(shù)值最小者對應(yīng)的微震震源優(yōu)化解作為最終微震震源參數(shù)計算值.

實際工程應(yīng)用中，微震監(jiān)測初期階段，根據(jù)地質(zhì)資料對區(qū)域巖體進行巖層劃分，并在建議的巖體波速值附近，利用模塊Ⅱ進行層速度組合最優(yōu)化求解.隨著對區(qū)域巖體的施工開挖和支護開展，巖體質(zhì)量不斷改變、應(yīng)力不斷調(diào)整，實際層速度組合也隨之改變.為了保證微震源定位精度，在監(jiān)測期間，應(yīng)定期對波速層重新劃分并進行現(xiàn)場人工爆破試驗對速度組合值進行更新計算.

4 工程應(yīng)用

白鶴灘水電站位于金沙江下游，左岸邊坡為順向坡，傾向河流上游南東方向，傾角較緩；巖體內(nèi)部發(fā)育主要斷層F14、F17及多個小斷層、層內(nèi)層間錯動帶.左岸拱肩槽建基面開挖高程為600～834 m，規(guī)模巨大.巖體開挖卸荷必然引起斷層及層間、層內(nèi)錯動帶高應(yīng)力區(qū)裂隙萌生、發(fā)育及擴展，內(nèi)部微破裂往往是巖質(zhì)邊坡外觀變形和失穩(wěn)破壞的前兆.

圖4 左岸邊坡傳感器陣列布置圖Fig.4 Layout of sensor array on slope of left bank

如圖4所示，為了實現(xiàn)對壩基邊坡及其坡內(nèi)廊道開挖卸荷過程圍巖損傷的實時監(jiān)測，左岸邊坡安裝了加拿大ESG微震監(jiān)測系統(tǒng)，系統(tǒng)由HNAS主機處理系統(tǒng)、Paladin(臺站)采集系統(tǒng)、單軸加速度傳感器和三維可視化軟件四部分組成.其中，HNAS通過網(wǎng)絡(luò)控制各個臺站(Paladin)，設(shè)置信號增益、采樣頻率和觸發(fā)閾值等參數(shù)，從而實現(xiàn)全天候連續(xù)采集；Paladin為廣泛使用的24位數(shù)據(jù)模數(shù)轉(zhuǎn)換系統(tǒng)；單軸加速度傳感器靈敏度為0.7Volts、頻率響應(yīng)范圍為50 Hz～5 kHz.白鶴灘左岸邊坡安裝了3套Paladin采集系統(tǒng)，共18通道單軸加速度傳感器陣列分別布置在610 m高程的1#灌排洞及施工連接洞、660 m高程的2#灌排洞和750 m高程的4#灌排洞及抗力體排水洞內(nèi).

4.1 左岸邊坡分層速度定位模型建立

4.1.1 波速層劃分

在實際工程中，按照地層巖石物性將區(qū)域巖體分成若干個厚度不等的巖層，各層之間的波速值各不相同，這種波速的分層同地層的地質(zhì)年代、巖性上的分層基本上是一致的.但地質(zhì)分層比速度分層細，在通常情況下地質(zhì)年代不相同但巖性相同的一些地層可以成為一個速度層，這在地震地層學(xué)中被稱為層序“穿時現(xiàn)象”.

如圖5所示，白鶴灘左岸邊坡為順層巖質(zhì)邊坡，各巖層產(chǎn)狀近似為N45°E，SE∠20°.開挖過后的壩基拱肩槽巖體巖質(zhì)堅硬，呈微新狀態(tài)，無卸荷或弱卸荷.根據(jù)壩基勘探資料，其巖體波速分類列于表1中.如圖6所示，根據(jù)傳感器網(wǎng)絡(luò)布置范圍和邊坡坡面形狀，確定模型尺寸和邊界面；結(jié)合巖性和巖石類別，將模型區(qū)域巖體進行波速分層處理：共分為11層，4種類別，各波速交替出現(xiàn)，相應(yīng)厚度由地質(zhì)勘探圖紙資料換算得到.

4.1.2 層速度值確定

現(xiàn)場收集到610 m高程灌排廊道的4個開挖爆破數(shù)據(jù)，采用模塊Ⅱ在波速建議值附近用遺傳算法搜尋層速度組合最優(yōu)解，計算結(jié)果列于表2中.取由4個爆破事件計算得到的層速度優(yōu)化值的近似平均值V1=4500 m·s-1，V2=4900 m·s-1，V3=3750 m·s-1，V4=5000 m·s-1作為分層速度定位模型的層速度組合值并用于模塊Ⅲ進行微震定位.

4.2 單一速度定位模型和分層速度定位模型對比

本文利用ESG微震監(jiān)測系統(tǒng)對左岸邊坡微震事件進行實時采集，拾取各有效通道信號起跳時刻，然后采用Matlab語言編程來分別實現(xiàn)單一速度已知模型和分層速度定位模型定位.

2015年4月1—18日期間，白鶴灘壩址區(qū)邊坡坡外部分未進行開挖，工程施工集中在750 m高程置換洞和610 m高程灌排廊道的開挖.邊坡層內(nèi)錯動帶LS331和LS337受到坡內(nèi)廊道開挖卸荷擾動，內(nèi)部微破裂沿其分布帶附近萌生和發(fā)展，是主要的巖體損傷區(qū).單一速度已知模型采用4個波速區(qū)波速的近似均值4600 m·s-1，分層速度模型采用前文所述模型.首先，利用單一速度已知模型對2015年4月1—18日期間微震監(jiān)測系統(tǒng)采集的微震事件進行定位，圈定層內(nèi)錯動帶LS331和LS337附近區(qū)域的微震事件；然后，利用分層速度模型對上述圈定事件進行重新定位.圖7給出了2015年4月6日18∶59∶55.234″事件的波形圖，其中有12個通道監(jiān)測到了有效的微震波形(DA_001等為通道編號).

表1 巖體波速分類Table 1 Classification of rock mass velocities

圖5 左岸邊坡拱軸線地質(zhì)剖面圖Fig.5 Geologic cross section along the arch axis at slope of left bank

時間爆破孔底坐標(m)層速度優(yōu)化值(m·s-1)NEUV1V2V3V42015-04-1815∶253012794.600589164.126613.3434563.164857.443641.865032.972015-04-1903∶023012750.942589220.404607.0454533.674708.463863.044957.352015-04-1904∶233012793.056589166.741612.9794512.954881.003616.945050.292015-04-1913∶013012848.356589028.321617.9494450.425053.863938.784942.95

圖6 分層速度定位幾何模型Fig.6 Geometric model of SV

單一速度模型和分層速度模型的定位結(jié)果以小球的形式分別畫在圖8a和8b中，其中小球的顏色表示微震事件發(fā)生所在的時間段.從圖8中定位結(jié)果圖對比來看，單一速度已知模型定位的微震事件分散在層內(nèi)錯動帶LS331和LS337附近，無明顯的分布規(guī)律；分層速度模型定位的微震事件更接近層內(nèi)錯動帶LS331和LS337，并呈條帶狀分布.通過計算某個微震事件到兩個層內(nèi)錯動帶的距離，然后選擇較小的距離值作為評價該微震事件偏離層內(nèi)錯動帶程度大小的指標.該指標值越小，表明微震事件的空間位置越靠近層內(nèi)錯動帶；反之，則表明越遠離層內(nèi)錯動帶.計算結(jié)果表明，單一速度模型定位的微震事件中，最大偏離程度指標為32.505 m，最小為1.364 m，平均為9.134 m；分層速度模型定位的微震事件中，最大偏離程度指標為13.923 m，最小為0.866 m，平均為3.895 m.由此可以看出，采用分層速度定位模型，最大、最小和平均偏離程度指標分別降低了57.17%、36.51%和57.35%，這說明白鶴灘分層速度定位模型的微震事件定位結(jié)果較單一速度模型更加合理可靠.需要指出的是，精確的速度層劃分是提高該算法定位精度的關(guān)鍵.白鶴灘分層速度定位模型中，對層速度組合進行更加精細的劃分能夠進一步提高定位精度.

圖7 2015年4月6日18∶59∶55.234″事件波形Fig.7 Waveforms of an event at 18∶59∶55.234″ on April 6, 2015

5 結(jié)論與展望

(1) 本文針對層狀區(qū)域巖體，在異向波速模型的基礎(chǔ)上，對垂向上應(yīng)力波按波速值大小進行分段區(qū)別，推導(dǎo)了分層速度應(yīng)力波走時關(guān)系式，建立了分層定位目標函數(shù)，并基于此提出了簡化的分層速度定位模型SV.該模型包括參數(shù)準備、層速度反演、微震定位三個模塊.

(2) 在進行分層速度目標函數(shù)優(yōu)化求解時，采用了具有全局最優(yōu)、不依賴優(yōu)化初值和強魯棒性特點的遺傳算法作為最優(yōu)化工具，使目標函數(shù)優(yōu)化求解過程穩(wěn)定性和全局最優(yōu)的要求得到了保證.

(3) 將該速度模型應(yīng)用到白鶴灘水電站左岸邊坡，通過定位效果對比，分層速度定位模型定位微震事件的最大、最小和平均偏離層內(nèi)錯動帶程度指標較單一速度模型分別降低了57.17%、36.51%和57.35%.

(4) 本文提出的改進定位模型考慮了因各巖層應(yīng)力波波速值大小不同而對走時計算的影響，同時能采用遺傳算法對建立的目標函數(shù)實現(xiàn)穩(wěn)定求解.相比單一速度模型，定位精度得到了提高.然而，受到巖體應(yīng)力松弛和風(fēng)化作用影響，同層巖體各部位波速存在差異，囿于層速度劃分無法達到絕對精細，勢必在一定程度上影響微震事件的定位精度.為了提高定位精度，后期需要考慮對層內(nèi)速度進行分區(qū)，在尋找能夠穩(wěn)定實現(xiàn)、更加符合應(yīng)力波實際走時計算的定位模型方面進行改進和完善.

Bai C Y, Greenhalgh S. 2006. 3D local earthquake hypocenter determination with an irregular shortest-path method.BulletinoftheSeismologicalSocietyofAmerica, 96(6): 2257-2268.

Bai C Y, Zhao R, Greenhalgh S. 2010. Rapid 3-D earthquake location using a hybrid global-local inversion approach.PureandAppliedGeophysics, 167(11): 1377-1387.

Chen B R, Feng X T, Li S L, et al. 2009. Microseism source location with hierarchical strategy based on particle swarm optimization.ChineseJournalofRockMechanicsandEngineering(in Chinese), 28(4): 740-749.Collins D S, Pinnock I, Toya Y, et al. 2014. Seismic event location and source mechanism accounting for complex block geology and voids.∥ 48th US Rock Mechanics/Geomechanics Symposium. Minneapolis, Minnesota: American Rock Mechanics Association. Dong L J, Li X B, Tang L Z, et al. 2011. Mathematical functions and parameters for microseismic source location without pre-measuring speed.ChineseJournalofRockMechanicsandEngineering(in Chinese), 30(10): 2057-2067.

Feng G L, Feng X T, Chen B R, et al. 2015. Sectional velocity model for microseismic source location in tunnels.TunnellingandUndergroundSpaceTechnology, 45: 73-83.

Gao Y T, Wu Q L, Wu S C, et al. 2015. Source parameters inversion based on minimum error principle.JournalofCentralSouthUniversity(ScienceandTechnology) (in Chinese), 46(8): 3054-3060.Gong S Y, Dou L M, Ma X P, et al. 2012. Study on the construction and solution technique of anisotropic velocity model in the location of coal mine tremor.ChineseJournalofGeophysics(in Chinese), 55(5): 1757-1763, doi: 10.6038/j.issn.0001-5733.2012.05.033.Hosseini Z, Collins D, Shumila V, et al. 2015. Induced Microseismic Monitoring in Salt Caverns.∥ 49th US Rock Mechanics/Geomechanics Symposium. San Francisco, California: American Rock Mechanics Association. Inglada V. 1928. The calculation of the stove coordinates of a Nahbebens.GerlandsBeitrageZurGeophysik, (19): 73-98.

Jiang F X, Yang S H, Cheng Y H, et al. 2006. A study on microseismic monitoring of rock burst in coal mine.ChineseJournalofGeophysics(in Chinese), 49(5): 1511-1516.

Kaiser J. 1950. Studies on the occurrence of noise in the tensile test [Ph. D. thesis]. Munich: Technical University of Munich.

Leighton F, Blake W. 1970. Rock noise source location techniques. Washington: Bureau of Mines.

Li J, Gao Y T, Xie Y L, et al. 2014. Improvement of microseism locating based on simplex method without velocity measuring.ChineseJournalofRockMechanicsandEngineering(in Chinese), 33(7): 1336-1345.

Liu F Y. 1984. Simultaneous inversion of earthquake hypocenters and velocity structure (Ⅰ)——theory and method.ChineseJ.Geophys. (ActaGeophysicaSinica) (in Chinese), 27(2): 167-175.

Liu J, Zheng S H, Wong Y L. 2003. The inversion of non-elasticity coefficient, source parameters, site response using genetic algorithms.ActaSeismologicaSinica(in Chinese), 25(2): 211-218.

Lian C, Li S L, Dong M, et al. 2006. Method of spherical surface shearing in earthquake location.JournalofGeodesyandGeodynamics(in Chinese), 26(2): 99-103.

Lu C P, Dou L M, Wang Y F, et al. 2010. Microseismic effect of coal materials rockburst failure induced by hard roof.ChineseJ.Geophys. (in Chinese), 53(2): 450-456, doi: 10.3969/j.issn.0001-5733.2010.02.024.

Obert L. 1977. The microseismic method: discovery and early history.∥Proceedings of the First Conference on Acoustic Emission/Microseismic Activity in Geologic Structures and Materials. Clausthal-Zellerfeld: Trans. Tech. Publications, 11-12. Pavlis G L, Booker J R. 1980. The mixed discrete-continuous inverse problem: Application to the simultaneous determination of earthquake hypocenters and velocity structure.JournalofGeophysicalResearch:SolidEarth, 85(B9): 4801-4810. Song W Q, Xu B B, Yu Z C, et al. 2015. Reconstruction of the micro-seismic source vector field and fissure interpretation based on the anisotropy analysis.ChineseJ.Geophys. (in Chinese), 58(2): 656-663, doi: 10.6038/cjg20150226.Spencer C, Gubbins D. 1980. Travel-time inversion for simultaneous earthquake location and velocity structure determination in laterally varying media.GeophysicalJournalInternational, 63(1): 95-116. Spence W. 1980. Relative epicenter determination using P-wave arrival-time differences.BulletinoftheSeismologicalSocietyofAmerica, 70(1): 171-183.

Thurber C H. 1985. Nonlinear earthquake location: theory and examples.BulletinoftheSeismologicalSocietyofAmerica, 75(3): 779-790.

Wang H T, Chen Y, Zhuang C T. 2006. A local earthquake locating method based on genetic algorithm.Earthquake(in Chinese), 26(2): 12-18.

Xu N W, Tang C A, Li L C, et al. 2011. Microseismic monitoring

and stability analysis of the left bank slope in Jinping first stage hydropower station in southwestern China.InternationalJournalofRockMechanicsandMiningSciences, 48(6): 950-963.

Zhang B H, Deng J H, Zhou Z H, et al. 2012. Analysis of monitoring microseism in areas controlled by faults near powerhouse in Dagangshan hydropower station.RockandSoilMechanics(in Chinese), 33(Suppl. 2): 213-218.

Zhang H, Zhang Z L, Yao H, et al. 2014. The study on precise positioning of earthquake in Dachang mining area based on genetic algorithm.ChineseJournalofEngineeringGeophysics(in Chinese), 11(2): 260-265.

Zhao Z, Zeng R S. 1987. A method for improving the determination of earthquake hypocenters.ChineseJ.Geophys. (in Chinese), 30(4): 379-388.

附中文參考文獻

陳炳瑞, 馮夏庭, 李庶林等. 2009. 基于粒子群算法的巖體微震源分層定位方法. 巖石力學(xué)與工程學(xué)報, 28(4): 740-749.

董隴軍, 李夕兵, 唐禮忠等. 2011. 無需預(yù)先測速的微震震源定位的數(shù)學(xué)形式及震源參數(shù)確定. 巖石力學(xué)與工程學(xué)報, 30(10): 2057-2067.

高永濤, 吳慶良, 吳順川等. 2015. 基于誤差最小原理的微震震源參數(shù)反演. 中南大學(xué)學(xué)報, 46(8): 3054-3060.

鞏思園, 竇林名, 馬小平等. 2012. 煤礦礦震定位中異向波速模型的構(gòu)建與求解. 地球物理學(xué)報, 55(5): 1757-1763, doi: 10.6038/j.issn.0001-5733.2012.05.033.

姜福興, 楊淑華, 成云海等. 2006. 煤礦沖擊地壓的微地震監(jiān)測研究. 地球物理學(xué)報, 49(5): 1511-1516.

李健, 高永濤, 謝玉玲等. 2014. 基于無需測速的單純形法微地震定位改進研究. 巖石力學(xué)與工程學(xué)報, 33(7): 1336-1345.

劉福田. 1984. 震源位置和速度結(jié)構(gòu)的聯(lián)合反演(Ⅰ)——理論和方法. 地球物理學(xué)報, 27(2): 167-175.

劉杰, 鄭斯華, 黃玉龍. 2003. 利用遺傳算法反演非彈性衰減系數(shù)、震源參數(shù)和場地響應(yīng). 地震學(xué)報, 25(2): 211-218.

廉超, 李勝樂, 董曼等. 2006. 球面交切法地震定位. 大地測量與地球動力學(xué), 26(2): 99-103.

陸菜平, 竇林名, 王耀峰等. 2010. 堅硬頂板誘發(fā)煤體沖擊破壞的微震效應(yīng). 地球物理學(xué)報, 53(2): 450-456, doi: 10.3969/j.issn.0001-5733.2010.02.024.

宋維琪, 徐奔奔, 喻志超等. 2015. 基于各向異性分析的微地震震源矢量場重建和裂縫解釋. 地球物理學(xué)報, 58(2): 656-663, doi: 10.6038/cjg20150226.

王洪體, 陳陽, 莊燦濤. 2006. 基于浮點遺傳算法的近震定位方法. 地震, 26(2): 12-18.

張伯虎, 鄧建輝, 周志輝等. 2012. 大崗山水電站廠房斷層控制區(qū)域微震監(jiān)測分析. 巖土力學(xué), 33(增2): 213-218.

張華, 張忠利, 姚宏等. 2014. 基于遺傳算法的大廠礦區(qū)地震定位研究. 工程地球物理學(xué)報, 11(2): 260-265.

趙珠, 曾融生. 1987. 一種修定震源參數(shù)的方法. 地球物理學(xué)報, 30(4): 379-388.

(本文編輯 何燕)

Improvement of microseismic location based on an anisotropic velocity model

DAI Feng1， GUO Liang1， XU Nu-Wen1,2*， FAN Yi-Lin3， XU Jian3， JIANG Peng1

1StateKeyLaboratoryofHydraulicsandMountainRiverEngineering,SichuanUniversity,Chengdu610065，China2SchoolofCivilandHydraulicsEngineering,ShandongUniversity,Ji′nan250061,China3ChinaThreeGorgesCorporation，Beijing100038，China

As a kind of three-dimensional detection of seismic events, microseismic (MS) monitoring techniques have been widely used in the world for many years to assess the stability of surrounding rock. The location of MS events is foundation of MS monitoring. Resulting from anisotropy and inhomogeneity of regional rock masses, the velocity of stress waves may be different in different directions and areas. So, it is difficult to figure out a reasonably accurate location. Aimed at the regional rock masses with stratified velocities, a simplified and stable stratified velocity location model (SV) is proposed.

Microseismic monitoring; Anisotropic velocity model; Stratified velocity location model (SV); Genetic algorithm (GA); Rock slope

10.6038/cjg20160914.Dai F， Guo L， Xu N W, et al. 2016. Improvement of microseismic location based on an anisotropic velocity model.ChineseJ.Geophys. (in Chinese),59(9):3291-3301,doi:10.6038/cjg20160914.

國家重點基礎(chǔ)研究發(fā)展計劃(973計劃)項目(2015CB057903)；國家自然科學(xué)基金項目(51374149，51679158)資助.

戴峰，男，1978年生，博士，教授，2010年博士畢業(yè)于加拿大多倫多大學(xué)土木工程系，主要從事巖石動力學(xué)與工程研究.

*通訊作者 徐奴文，男，1981年生，博士，副研究員，2005年畢業(yè)于武漢科技大學(xué)土木工程系，主要從事巖土工程災(zāi)害機理與微震監(jiān)測方面的研究. E-mail：xunuwen@scu.edu.cn

10.6038/cjg20160914

P315

2015-10-09，2016-01-28收修定稿

戴峰， 郭亮， 徐奴文等. 2016. 基于異向波速模型的微震定位改進. 地球物理學(xué)報，59(9):3291-3301,