基于序貫高斯模擬法的水利工程壩址區(qū)速度建模

李偉, 周懷來*, 劉興業(yè), 王元君, 馬慧蓮

(1.成都理工大學(xué)地球物理學(xué)院, 成都 610059; 2.成都理工大學(xué)地球勘探與信息技術(shù)教育部重點實驗室, 成都 610059)

全球?qū)λこ淘O(shè)施的建設(shè)需求的不斷增加,勘探力度逐步加大。中國擁有眾多河流,發(fā)展前景良好。在追求發(fā)展速度的同時,越來越注重水利等工程項目的質(zhì)量。水庫作為水利工程中的基礎(chǔ)設(shè)施之一,為區(qū)域居民的生活和生產(chǎn)提供了重要的水資源保障[1]。在修建水庫之前,壩址的選擇至關(guān)重要,必須對影響水庫后期投入使用的各種因素進(jìn)行必要的分析,也包括地質(zhì)災(zāi)害分析工作等[2],在全面符合條件的情況下才能進(jìn)行規(guī)劃設(shè)計和施工,以確保水庫的正常投入使用。在壩址的選擇中,分析壩址區(qū)的速度和巖性是重要的兩個參考因素[3]。通常,壩址區(qū)巖層速度的橫向變化范圍大,僅依靠測井資料的垂向速度信息很難作為依據(jù),需要利用其他條件加以約束[4]。通過建立速度模型或巖性模型可用于分析壩址區(qū)地下結(jié)構(gòu)的穩(wěn)定性,然而這對建模成像精度提出了更高的要求[5]。

水庫堤壩進(jìn)行選址勘察時,會對速度及地下巖石的性質(zhì)進(jìn)行重點勘察,速度分析是必不可少的一步。巖石越堅硬致密,測井上的速度響應(yīng)相對越高,而在巖石松散或孔隙較大部位,速度呈現(xiàn)異常低值,這對水利工程的修建以及后期的運行形成垮塌的風(fēng)險隱患。分析不同類型巖石的速度對堤壩類型和選址的影響也是重要的一步[6]。對壩址區(qū)建立地下巖石的速度模型,是分析壩址穩(wěn)定性的一種方法。然而這對建模精度有很高的要求,常規(guī)的建模方法很難實現(xiàn),此問題待解決。為得到更高精度的近地表速度模型,學(xué)者們開展了許多的研究。趙玲芝等[7]運用多信息融合的方法建立近地表速度模型,將近地表調(diào)查資料與大炮初值層析反演結(jié)果進(jìn)行協(xié)克里格插值技術(shù)融合,得到更高精度的速度模型并在靜校正及疊前深度偏移中得到了很好的效果。劉靈等[8]提出一種相控速度建模方法,綜合地震速度、測井速度、地震屬性及沉積相等多種信息,通過多屬性融合及神經(jīng)網(wǎng)絡(luò)聚類分析來建立低頻相控速度模型,并用于約束反演結(jié)果,極大提高了速度建模的精度及分辨率。

水利工程設(shè)施的建設(shè)對地下速度建模精度要求嚴(yán)格[9],而常規(guī)的速度建模思路是利用測井資料得到數(shù)據(jù)進(jìn)行插值,多為線性插值方法[10],對于水利壩址區(qū)此類復(fù)雜地層結(jié)構(gòu),常規(guī)插值方法難以解決問題,此外井?dāng)?shù)據(jù)只包含縱向速度信息,對于井之間的橫向大范圍空白區(qū)域缺乏約束條件。綜上所述常規(guī)插值方法比較簡單[11],建模效果不佳,需尋找高精度的速度建模方法。為此,有學(xué)者提出了利用地震層位約束來進(jìn)行高精度建模的方法。在地震勘探領(lǐng)域,孟會杰等[12]為解決地震偏移成像精度不高的問題,應(yīng)用層位約束的網(wǎng)格層析建模方法,通過不斷迭代、反演最終建立地下深度域高精度速度模型,并消除了實際工區(qū)中的構(gòu)造成像假象。后來,有學(xué)者在對地質(zhì)統(tǒng)計學(xué)建模的研究中發(fā)現(xiàn),相比于常規(guī)線性插值建模,基于地質(zhì)統(tǒng)計學(xué)理論的建模方法在效果上更好,精度更高,能夠充分利用變量的空間變化特征[13]。

為提高建模精度,準(zhǔn)確分析水利壩址區(qū)內(nèi)部構(gòu)造特征,服務(wù)于壩址的選取及水庫的建設(shè),鑒于前人專家學(xué)者的研究,現(xiàn)提出一種基于地質(zhì)統(tǒng)計學(xué)理論的地震與測井相結(jié)合的速度建模方法,對常規(guī)的測井插值建模的方法以及地質(zhì)統(tǒng)計學(xué)建模方法進(jìn)行完善和補充,應(yīng)用地質(zhì)統(tǒng)計學(xué)序貫高斯模擬隨機(jī)建模,以地震層位信息作為橫向約束條件聯(lián)合測井資料建立速度模型,為壩址區(qū)巖石穩(wěn)定性提供數(shù)據(jù)支撐,對選壩提供重要參考。同時,在前人對速度建模方法研究的基礎(chǔ)之上,創(chuàng)新性地將本文提出的基于地質(zhì)統(tǒng)計學(xué)理論的建模方法應(yīng)用于水利工程領(lǐng)域的壩址選取工作上。

1 方法原理

巖石速度是壩址選取的重要影響因素,因此以速度作為關(guān)鍵切入點對水利壩址區(qū)建立速度模型,進(jìn)而評估待選壩址區(qū)的穩(wěn)定性。常規(guī)速度建模方法較為簡單,模型精度無法滿足水利工程建設(shè)的要求,而使用基于地質(zhì)統(tǒng)計學(xué)理論的地震與測井相結(jié)合的方法建立高精度速度模型,旨在解決壩址區(qū)速度建模精度低的問題,在高精度建模結(jié)果上圈定有利于修建壩址的區(qū)域,更加準(zhǔn)確且具有說服力,為選壩修壩提供重要依據(jù)。

地質(zhì)統(tǒng)計學(xué)作為一門發(fā)展于20世紀(jì)60年代的數(shù)學(xué)地質(zhì)邊緣學(xué)科,包含空間相關(guān)性分析、克里金插值以及隨機(jī)模擬三部分[14]。將研究區(qū)域在未知位置的函數(shù)值作為一組相關(guān)隨機(jī)變量,通過隨機(jī)變量理論分析來預(yù)測未知點的信息。相較于線性的插值方法有獨特的優(yōu)勢特征,且超越了普通的插值問題[15],這也是利用地質(zhì)統(tǒng)計學(xué)進(jìn)行速度建模的理論依據(jù)。

地質(zhì)統(tǒng)計學(xué)提出區(qū)域化變量[16]這一概念,能用空間分布來表達(dá)具有結(jié)構(gòu)性的隨機(jī)變量。同時,為度量區(qū)域化變量的空間變異性,提出了變差函數(shù)的概念,它兼顧區(qū)域化變量的隨機(jī)變化和空間結(jié)構(gòu)性規(guī)律[17],表示空間中數(shù)據(jù)方差隨著兩點之間相對距離的變化而在特定方向上的變異性[18]。假設(shè)隨機(jī)變量ξ滿足本征假設(shè),則在任一方向α,相對距離為|h|的兩個隨機(jī)變量ξ(u)和ξ(u+h)的增量的方差稱為變差函數(shù)γ(h),表示為

(1)

式(1)中:u為空間中的任意點;h為兩點之間的距離,在地質(zhì)統(tǒng)計學(xué)中稱為滯后距;E為數(shù)學(xué)期望。在平穩(wěn)假設(shè)條件下變差函數(shù)γ(h)的變化只與滯后距h有關(guān),而γ(h)隨h的變化稱為變差函數(shù)曲線,如圖1所示。圖1中c0表示滯后距h很小時,兩點之間具有一定的變異性。a表示變程,指區(qū)域化變量在空間中的相關(guān)范圍,變程越大說明觀測數(shù)據(jù)在更大范圍內(nèi)相關(guān)。超出變程范圍后,兩點之間無關(guān);在變程內(nèi),變量越近,變異性越小。c為基臺值,代表變量在空間中總變異性的大小。在速度建模過程中,對變差函數(shù)參數(shù)進(jìn)行分析是至關(guān)重要的一步,其結(jié)果直接影響建模效果,所以需要對參數(shù)進(jìn)行細(xì)致的分析。

圖1 變差函數(shù)曲線Fig.1 Variation function curve

通常實際觀測獲取的樣品數(shù)目有限,例如在速度建模過程中能提供的鉆井曲線信息較少,而通過有限觀測值經(jīng)式(1)計算獲得的變差函數(shù)作為實驗變差函數(shù)γ*(h),表示為

(2)

式(2)中:N(h)為滯后距為h的點對個數(shù);ξ為隨機(jī)變量;ui為空間中的任意點。

以h為橫坐標(biāo),γ*(h)為縱坐標(biāo)即可得到實驗變差函數(shù)圖。實際應(yīng)用中,可供利用的點對數(shù)量越多,獲取的實驗變差函數(shù)就具有越強(qiáng)的代表性。

實際的地質(zhì)統(tǒng)計學(xué)研究用到的都是理論變差函數(shù),因此需要在獲得實驗變差函數(shù)之后進(jìn)行變差函數(shù)擬合以得到理論變差函數(shù)的參數(shù)[19]。其中球狀模型、指數(shù)模型、冪模型及高斯模型等是一些基本理論變差函數(shù)。根據(jù)實際原始數(shù)據(jù)特點選取相應(yīng)的變差函數(shù)模型,不同模型的插值結(jié)果也會不同。本文研究使用球型變差函數(shù)模型,其靈活性大、數(shù)據(jù)適應(yīng)性強(qiáng),模型在原點處呈線性變化,在變程處達(dá)到基臺值。圖2為球型、高斯型、指數(shù)型這三類常用理論模型的變差函數(shù)曲線。

圖2 常用變差函數(shù)模型Fig.2 Common variation function model

使用球形模型進(jìn)行變量空間估計時,需定義球體變差函數(shù)的長軸、短軸和中長軸的變程以及塊金值、拱高和基臺值等。同時定義球體的空間走向,作為變差函數(shù)的方位參數(shù)。本文的速度建模方法將測井得到的速度數(shù)據(jù)視作地質(zhì)統(tǒng)計學(xué)中的區(qū)域化變量,對其進(jìn)行變差函數(shù)分析,研究其數(shù)據(jù)特征選擇球型變差函數(shù)模型,求取最佳的變差函數(shù)參數(shù)。

確定變差函數(shù)重要參數(shù)后,即可進(jìn)行地質(zhì)統(tǒng)計學(xué)序貫高斯模擬。它是隨機(jī)模擬的一種方法,以已知信息為約束,使用隨機(jī)函數(shù)建立多個不同地質(zhì)模型的方法。在模擬過程中考慮區(qū)域化變量隨機(jī)性的同時,還考慮了其空間相關(guān)性。且要求模擬結(jié)果忠實于觀測點的值[20]。假設(shè)隨機(jī)變量符合高斯分布,結(jié)合高斯概率理論與序貫?zāi)M算法對連續(xù)分布的區(qū)域化變量進(jìn)行模擬,以像元為單位按照選定的隨機(jī)路徑進(jìn)行,并且將已經(jīng)模擬過的節(jié)點的狀態(tài)值也直接加入條件數(shù)據(jù)中。模擬步驟如圖3所示[21]。序貫高斯模擬方法建立速度模型,相比于常規(guī)的建模方法在原理上有明顯優(yōu)勢,它能充分利用實際觀測數(shù)據(jù)的空間信息,并對其進(jìn)行變差函數(shù)分析以此建立更為精準(zhǔn)的速度模型,為壩區(qū)提供高精度的地下速度信息,作為實際選壩的主要參考。

序貫高斯模擬可將地震層位作為約束條件加入模擬結(jié)果進(jìn)行隨機(jī)建模。而無層位約束的模擬結(jié)果無法反映地下巖層的橫向變化特征,會得到一個相對光滑、連續(xù)的模型結(jié)果。水利壩址區(qū)地下構(gòu)造因其速度在橫向變化范圍大,若無橫向方向的約束將缺失建模的精度,最終模型結(jié)果無說服性。因此文章將地震層位作為約束進(jìn)行速度建模,提高建模的精度與準(zhǔn)確性。

序貫高斯模擬法在多次模擬時會生成多個等概率的模型,其屬于條件模擬方法,因此模擬結(jié)果忠于觀測點處的數(shù)據(jù),而對于非觀測點區(qū)域,模型之間差異較大?？蓪Χ啻文M得到的速度模型之間有明顯差異的區(qū)域做不確定性分析,選壩時則應(yīng)避開此類位置。

2 建模步驟

(1)對測井速度曲線進(jìn)行粗化處理,提高垂向采樣使數(shù)據(jù)更加連續(xù),同時進(jìn)行井曲線的標(biāo)準(zhǔn)化,使其更加符合地質(zhì)規(guī)律,便于后續(xù)的速度插值與建模。

(2)分析變差函數(shù)參數(shù)。根據(jù)已有數(shù)據(jù)特征和精度,選取球形模型變差模型進(jìn)行建模,確定球形模型主、次、垂向變程,以及最優(yōu)的塊金值、拱高和基臺值等重要參數(shù)。

(3)對選取的變差函數(shù)參數(shù)進(jìn)行序貫高斯模擬,首先進(jìn)行無層位約束模擬,再以地震層位作為橫向約束進(jìn)行模擬,完成不同層位的變差函數(shù)分析和分層速度建模。

(4)使用協(xié)克里金法將上述有、無層位約束情況下建立的速度模型進(jìn)行加權(quán)[22],增強(qiáng)模型的縱向連續(xù)性以及準(zhǔn)確性,形成最終的速度模型。并與常規(guī)的建模方法進(jìn)行對比、評價。

(5)將最終得到的模型使用Garden公式計算出密度,得出反射系數(shù),與雷克子波進(jìn)行褶積得到合成地震記錄,與實際地震記錄相對比,評價速度建模的效果。

3 應(yīng)用實例

在對以上方法研究的基礎(chǔ)上,對位于四川省蒙頂山某一處的近地表鉆井資料進(jìn)行實際數(shù)據(jù)測試,該區(qū)域地層速度橫向變化范圍大,難以獲得精確的地下速度分布,因此需在高精度速度建模的結(jié)果上進(jìn)行選址分析。工區(qū)范圍內(nèi)分布11口可用井(W1～W11),井位信息如圖4所示,井位分布比較均勻。首先對工區(qū)內(nèi)井曲線做批量標(biāo)準(zhǔn)化處理,得到的速度直方圖如圖5所示。圖6為井W6～W10 5口井標(biāo)準(zhǔn)化后的結(jié)果。此外對井速度數(shù)據(jù)做“粗化”處理,取速度值的算術(shù)平均對測井曲線進(jìn)行重采樣,增強(qiáng)垂向上的連續(xù)性。

圖4 壩區(qū)井位信息Fig.4 Well location information of dam area

圖5 井曲線標(biāo)準(zhǔn)化前后直方圖Fig.5 Histogram of well curve before and after standardization

圖6 井曲線(W6～W10)標(biāo)準(zhǔn)化前后結(jié)果Fig.6 Well curve (W6～W10 ) results before and after standardization

完成上述兩步驟后進(jìn)行序貫高斯模擬,若以地震層位作為橫向約束,則需要對每一層變差函數(shù)進(jìn)行分析。表1為有層位無層位約束的變差函數(shù)參數(shù)的分析結(jié)果。本文變差函數(shù)模型選擇球形模型,由于水平方向變差函數(shù)參數(shù)對模擬結(jié)果影響很小,因此將其視為各向同性,只對垂向變差函數(shù)參數(shù)進(jìn)行分析。

表1 有層位、無層位約束變差函數(shù)參數(shù)Table 1 Variation function parameters with horizon constraints and without horizon constraint

為分析不同變差函數(shù)對模擬結(jié)果的影響,下文使用無層位約束的模擬方式并通過設(shè)置不同的塊金值n(0、0.15、0.30、0.45)及垂向變程v(5、25、45、65)進(jìn)行測試,如圖7和圖8所示。

圖7 不同塊金值的模擬結(jié)果Fig.7 The simulation results of different nugget values

圖8 不同垂直變程的模擬結(jié)果Fig.8 The simulation section results of different vertical ranges

通過分析,塊金值越大,變量的隨機(jī)性表現(xiàn)得越強(qiáng)烈,當(dāng)塊金值為0時,模擬結(jié)果平滑,能較好體現(xiàn)出研究區(qū)域的速度變化趨勢,隨機(jī)性較小;對于垂向變差函數(shù)分析結(jié)果,不同垂直變程對模擬結(jié)果的差別較大。取較小的垂直變程能夠提高垂向分辨率,但對區(qū)域整體趨勢刻畫不足;若取較大的垂直變程,能體現(xiàn)出區(qū)域整體的變化趨勢,卻缺失垂向分辨率。因此選取最佳垂向變差參數(shù)才能達(dá)到較好的模擬結(jié)果,此處選取25作為最佳垂直變程。

完成變差函數(shù)分析后,使用序貫高斯模擬方法得到速度模型。首先圖9(a)展示無地震層位約束的速度模型過井W6～W10的剖面,圖中5條黑色橫線為地震層位,由上至下依次命名為S1～S5(S1、S5分別對應(yīng)t=0 ms和t=90 ms時刻)。從圖中得出由于沒有層位的約束,在分層處速度橫向變化特征未能突出。圖9(b)為有地震層位約束后的模擬結(jié)果,在對每個層位的井?dāng)?shù)據(jù)都做了獨立的變差函數(shù)分析后,得到了相比無層位約束更高精度的結(jié)果,同時能反映出橫向速度變化特征,但其速度模擬結(jié)果在局部層位分界處之間變化劇烈,與常規(guī)地質(zhì)規(guī)律認(rèn)識不符。

圖9 無層位約束與有層位約束模擬結(jié)果過井W6～W10剖面Fig.9 Simulation section results without and with horizon constraint cross well W6～W10

將無層位約束的速度模擬結(jié)果與有層位約束的模擬結(jié)果進(jìn)行協(xié)同克里金加權(quán)插值,加權(quán)步驟中將有層位約束的結(jié)果設(shè)置為主變量,無層位約束的模擬結(jié)果為次變量,通過分析確定權(quán)值后進(jìn)行插值最終得到圖10所示的速度模型。對模擬結(jié)果進(jìn)行分析,模型既能體現(xiàn)出速度的橫向變化的復(fù)雜性,也符合地質(zhì)認(rèn)識,且相較于上文結(jié)果精度進(jìn)一步提升,為水庫的壩址選取分析以及安全性評價提供重要的數(shù)據(jù)支撐。

圖10 協(xié)同克里金序貫高斯結(jié)果過井W6～W10剖面Fig.10 Co-kriging sequential gaussian section results cross well W6～W10

本文研究中使用的序貫高斯模擬屬于隨機(jī)模擬方法,多次模擬能得到多個等概率的模擬結(jié)果,但在鉆井位置均能忠實于觀測數(shù)據(jù)。多次模擬結(jié)果的差異體現(xiàn)了建模結(jié)果的不確定性。每次模擬后可對速度值差異明顯的區(qū)域進(jìn)行分析和取舍,以此縮小范圍確定最佳的選壩位置。圖11展示使用協(xié)克里金序貫高斯再次模擬兩次的結(jié)果。

圖11 兩次協(xié)克里金序貫高斯模擬結(jié)果過井W6～W10剖面Fig.11 Twice co-kriging sequential gaussian simulation section results cross well W6～W10

結(jié)合圖10,在剖面上紅色區(qū)域圈定出了一處三次模擬結(jié)果未發(fā)生明顯變化且速度響應(yīng)較高的區(qū)域,說明該位置的模擬結(jié)果響應(yīng)更趨近于真實地下構(gòu)造。在水庫修建中,基于安全考慮選壩位置在速度高值響應(yīng)特征的巖層中最為合適,而此分析結(jié)果為壩址選取提供了一項重要的參考依據(jù)。

下文分析常規(guī)速度建模方法效果。圖12(a)為克里金插值法得到的速度模型,其結(jié)果得到的是一種相對較平滑的結(jié)果,局部的速度異常點未得到突出,只能用作反映全局特征。其次克里金插值是一種確定性方法,每次插值得到的模型相同,無法反映模型之間的差異性以及不確定性。

圖12 克里金插值與井速度散點插值結(jié)果過井W6～W10剖面Fig.12 Kriging interpolation and well velocity scatter interpolation section results cross well W6～W10

圖12(b)是井速度散點直接插值為速度體后的結(jié)果,圖中的黑色豎線為鉆井提供的速度觀測值。分析后得出插值效果并不理想,速度呈明顯條狀分布,且插值結(jié)果不完整,模型出現(xiàn)大片空白區(qū)域以及速度值為0的異常位置,不宜作為最終速度模型的參考。

通過以上分析,使用地震層位約束的地質(zhì)統(tǒng)計學(xué)協(xié)克里金序貫高斯模擬得到的速度模型精度最佳,能較好地展示壩址區(qū)地下巖石的速度信息。圖13為t=63.60 ms時刻的速度模型平面圖,通過分析速度在平面上的分布特征,發(fā)現(xiàn)速度高值響應(yīng)大多分布于壩區(qū)靠南區(qū)域,圖13中紅色實線為上文展示的過井剖面圖在此平面的投影,藍(lán)色框線部分圈定了一處有利修建水庫的壩址區(qū),區(qū)域內(nèi)速度值較高,同時在過井剖面上為高值特征響應(yīng),其結(jié)果為水利工程壩址區(qū)的選址提供了重要的參考依據(jù)。

圖13 平面上t=63.60 ms時刻協(xié)克里金序貫高斯法速度模型Fig.13 Co-kriging sequential gaussian velocity model on the plane at time t=63.60 ms

以圖10速度模型作為最終結(jié)果,利用Garden公式計算密度及反射系數(shù),并與主頻為50 Hz的雷克子波進(jìn)行褶積得到合成地震記錄,檢驗建模效果,如圖14(a)所示。

圖14 合成地震記錄與實際地震記錄過井W6～W10剖面Fig.14 Synthetic seismograms and the actual seismic record section cross well W6～W10

將計算得到的合成記錄與實際地震記錄圖14(b)進(jìn)行對比,從中分析得出:合成記錄與實際資料的吻合程度較好,地震層位對應(yīng)地震剖面上的同相軸部分連續(xù)性較好,通過本文的方法建立的速度模型對壩址區(qū)的穩(wěn)定性分析及選址都能夠提供重要的依據(jù),達(dá)到了速度建模的目的及實際效果。

本文提出的基于地質(zhì)統(tǒng)計學(xué)理論的高精度速度建模方法圈定預(yù)測出的實際區(qū)域有利壩址區(qū),已被當(dāng)?shù)叵嚓P(guān)水利施工單位認(rèn)證,在實際現(xiàn)場勘測中,區(qū)域內(nèi)地下巖性普遍致密,速度呈現(xiàn)高值響應(yīng)。壩址位置選取正確,建設(shè)方案已被采取,現(xiàn)正進(jìn)入施工階段。

4 結(jié)論

(1)地層的速度信息是影響壩址區(qū)選擇以及后期水利工程安全投入使用的重要因素。因此,對速度進(jìn)行分析是必要的一步,有必要運用近地表速度模型的精細(xì)構(gòu)建技術(shù)來確立精確近地表速度模型。

(2)針對水利工程壩址的選取,本文提出基于地質(zhì)統(tǒng)計學(xué)隨機(jī)模擬的速度建模方法,并引入地震層位橫向約束建模過程。根據(jù)壩址區(qū)近地表速度分布的特點選取合理的變差函數(shù)模型是確保模型精度的重要前提,若參數(shù)選取不當(dāng)則得到的速度模型無法用來評估壩址區(qū)的穩(wěn)定性以及后期投入使用。利用測井?dāng)?shù)據(jù)作為硬約束,地震層位作為橫向約束的地質(zhì)統(tǒng)計學(xué)序貫高斯模擬方法建立的速度模型符合地質(zhì)規(guī)律。

(3)多次模擬可以得到多個等概率的速度模型,不同模型之間的差異可以反映建模結(jié)果的不確定性,為壩址選擇的風(fēng)險決策提供依據(jù)。將提出的方法應(yīng)用于實際數(shù)據(jù)資料上并取得了較好的成果,依據(jù)建立的速度模型分析并圈定出了有利于修建水利工程的壩址區(qū),圈定方案被采用,實際正在施工當(dāng)中,為水庫的選壩的合理性及后期的投入使用提供了重要的數(shù)據(jù)支撐與保障。