陳柏平崔凡劉波杜云飛王子昌
1.中國礦業(yè)大學(北京) 地球科學與測繪工程學院,北京 100083;2.河北工程大學 地球科學與工程學院,河北邯鄲 056038
煤炭是我國目前乃至未來很長時期內(nèi)能源供給的重要組成。 近年來,國家對煤炭安全愈發(fā)重視,為保障煤炭開采的安全性、高效性,煤炭開采智能化、無人化等已成為未來發(fā)展趨勢[1-2]。 袁亮[3-4]于2017年提出煤炭精準開采體系,并指出精準開采和物聯(lián)網(wǎng)架構(gòu)將是實現(xiàn)未來無人礦山的重要組成部分。 毛善君等[5]在2018年提出透明化礦山的概念,認為透明化礦山是利用鉆探、物探、煤巖層識別等技術(shù)手段,基于統(tǒng)一的數(shù)據(jù)管控平臺和生產(chǎn)實時揭露信息,實現(xiàn)礦山井上下工程、設(shè)備、開采環(huán)境等實時信息的高精度可視化,為礦山的安全生產(chǎn)、培訓和智能開采提供高技術(shù)手段和統(tǒng)一的可視化管控平臺。 因此,建立精準開采和物聯(lián)網(wǎng)架構(gòu)體系,以礦井地質(zhì)、鉆井和測井、綜合地球物理勘探數(shù)據(jù)為基礎(chǔ)的地質(zhì)保障技術(shù)至關(guān)重要,這些數(shù)據(jù)又是建立在地下空間地質(zhì)體架構(gòu)上的[6-7]。 精確的空間地質(zhì)模型,不僅可以為透明化礦山的建設(shè)提供前期基礎(chǔ)地質(zhì)框架設(shè)計依據(jù),還能降低后續(xù)通過綜合物探技術(shù)獲取的動態(tài)、實時監(jiān)測等數(shù)據(jù)進行地質(zhì)模型修正的工作量。
“十三五”期間,隨著煤炭工業(yè)智能、綠色發(fā)展,構(gòu)建精確的空間地質(zhì)模型是煤炭透明、智能、精準開采方案中實現(xiàn)可視化、指導智能開采的重要基礎(chǔ)[8-10]。 近年來,國內(nèi)眾多學者對透明化礦山建設(shè)和煤炭開采做了眾多有意義的研究,并多次提到三維地質(zhì)模型的重要性。 毛善君等[5]提出的透明化礦山構(gòu)建原則中指出,透明化三維模型構(gòu)建是關(guān)鍵技術(shù)之一。 李梅等[11]提到三維地質(zhì)建模是礦山三維可視化系統(tǒng)的重要組成部分。 目前,獲取井下三維地質(zhì)體信息的常用方法是基于鉆井、測井、三維地震等數(shù)據(jù)進行三維地質(zhì)建模。 陳博等[12]基于鉆井和測井資料建立煤體構(gòu)造識別方法,并利用隨機反演方法實現(xiàn)煤體結(jié)構(gòu)的三維可視化。 呂堂杰等[13]結(jié)合鉆井、測井地質(zhì)分析,并基于層序旋回識別提高了單煤層儲層厚度的解釋精確度。 基于鉆井、測井、地質(zhì)資料的建模方法具有較高的縱向分辨率,對井點和周圍區(qū)域的預(yù)測準確性較佳,但存在施工人力、物力等成本高的問題。 同時,面對較大區(qū)域時,由于井數(shù)量有限、地形和構(gòu)造的影響,井間區(qū)域的預(yù)測效果已不能滿足透明化礦山建設(shè)的高分辨率地質(zhì)模型需求。 為解決高精度地質(zhì)模型構(gòu)建問題,近年來眾多學者將鉆井、測井、地震等數(shù)據(jù)進行結(jié)合,利用三維地震橫向連續(xù)性強的特點,以鉆井、測井信息為約束,形成了基于三維地震反演的儲層建模方法。 胡勇等[14]基于地震反演成果進行建模,提高了建模精度;韓冰等[15]利用測井資料建立的高頻地質(zhì)模型和疊后地震阻抗反演建立的低頻模型,實現(xiàn)了地層模型的三維可視化,提高了地質(zhì)模型的空間精度;孫月成等[16]基于地質(zhì)統(tǒng)計學反演方法綜合鉆井、測井、地震等多尺度信息建立高精度三維地質(zhì)模型,進一步提高了氣田三維地質(zhì)模型精度。
近年來,隨著儲層預(yù)測精度要求的提高和三維地震反演技術(shù)的發(fā)展,常規(guī)的地震反演技術(shù)的精度已不能滿足透明化礦山建設(shè)高精度數(shù)據(jù)的需求[17]。 地質(zhì)統(tǒng)計學反演方法綜合了鉆井、測井、地震和地質(zhì)等多方面資料,橫、縱方向上都具有高分辨率的優(yōu)勢,已被廣泛應(yīng)用到油、氣田、海洋等領(lǐng)域[18-20]。 但在煤田領(lǐng)域上,目前應(yīng)用較少,特別是缺乏針對煤系地層環(huán)境特點的反演參數(shù)和井約束條件的研究。 為此,本文通過建立三維地震正演模型,設(shè)置多組地質(zhì)統(tǒng)計學反演參數(shù)和井約束條件進行實驗,分析地質(zhì)統(tǒng)計學反演在煤系地層中的煤層及其厚度反演效果,并探討反演參數(shù)和約束井選取方案,最后將其應(yīng)用到實際煤田的三維巖性地質(zhì)模型構(gòu)建的研究中,以期獲得更加精準的效果,更好地服務(wù)于透明化礦山基礎(chǔ)模型的構(gòu)建。
地質(zhì)統(tǒng)計學反演是基于貝葉斯框架,通過地震、測井等驗證性地質(zhì)信息擬合得到概率密度函數(shù)(Probability Density Function,PDF)來表示儲層的巖性或物性特征的概率分布。 反演得到的計算公式是一個多維函數(shù),較為復雜,需要采用馬爾科夫鏈蒙特卡洛算法(Markov Chain Monte Carlo,MCMC)構(gòu)建馬爾科夫鏈進行隨機采樣和蒙特卡洛模擬,從而獲得需要的高分辨率數(shù)據(jù)體[21-22]。 貝葉斯判別是一種利用概率統(tǒng)計原理對多種信息進行整合并判別的方法,能夠公平地將多種不確定信息源進行統(tǒng)一分析判斷,常被用來在已知先驗信息概率(E)和假設(shè)條件概率(H)的前提下,計算某一事件(X)發(fā)生的后驗概率[23]。 數(shù)學表達式如下:
式中,P(X|H,E)為待求解的后驗函數(shù)(后驗概率分布);P(X|H)為在已知假設(shè)條件(H)下事件X的先驗函數(shù)(先驗概率分布);P(E|X)為似然函數(shù);P(H|E)為近似全概率模型。
本文涉及的先驗信息概率是通過測井的巖性和縱波阻抗信息擬合得到,主要包括概率密度函數(shù)和變程函數(shù)兩個主要的反演參數(shù)。
概率密度函數(shù)表示儲層參數(shù)或?qū)傩栽诳臻g上的分布規(guī)律[24],可表示為
式中,m為巖性;f為縱波阻抗;i為采樣點;P(f)為巖性的概率;P(m1|f1)為頂部樣本點巖性對應(yīng)下阻抗的概率的分布;P(mi|m(i-1),f(i-1),fi)為第i處樣本的巖性對應(yīng)下阻抗的概率的分布。
可以看出,用來擬合縱波阻抗的高斯分布函數(shù)的均值和方差直接影響到先驗函數(shù)的準確性[25]。
變程函數(shù)表示儲層屬性在空間位置上的延展性,可以通過式(3)進行變差分析后擬合得到[18]。
式中,σf為樣本點對應(yīng)的縱波阻抗擬合的方差;Δt為采樣間隔時間;R為通過變差分析得到的變程。
變程具有一定的主觀性,變程越大相關(guān)性越強。 對于縱波阻抗,變程一般采用高斯函數(shù)進行擬合,變程的大小會影響反演儲層參數(shù)的橫、縱向的延伸程度;對于巖性,橫向變程越大,同一巖性的在橫向方向的延展程度越長,因此橫向變程的取值至關(guān)重要[26]。
基于MCMC 算法的地質(zhì)統(tǒng)計學反演過程中,先驗信息的概率密度函數(shù)和變程參數(shù)以及約束井條件對最終的反演結(jié)果至關(guān)重要[27]。 目前地質(zhì)統(tǒng)計學反演方法已被應(yīng)用到油田、氣田領(lǐng)域中,對薄儲層的預(yù)測效果良好。 但在煤田方面,尤其是透明化礦山建設(shè)中煤系地層地質(zhì)建模的需求,反演參數(shù)的設(shè)置和約束井的選擇,目前尚缺少具有針對性的研究。
本文主要基于GeoSoftware 地學軟件平臺中的Jsaon 軟件進行相關(guān)試驗和研究。 煤系地層地質(zhì)統(tǒng)計學反演參數(shù)和約束井條件選擇方法如下:
(1) 建立楔狀煤層縱波阻抗數(shù)值模型,通過波阻抗模型正演技術(shù)[28]獲得地震數(shù)據(jù)。
(2) 從阻抗模型中提取偽測井數(shù)據(jù)來模擬實際測井數(shù)據(jù),包括縱波阻抗和巖性數(shù)據(jù)。
(3) 根據(jù)偽井的縱波阻抗和巖性分布信息,設(shè)置不同的概率密度函數(shù)均值和橫向變程參數(shù)進行地質(zhì)統(tǒng)計學反演,并分析最佳的反演參數(shù)組合。
(4) 在最佳反演參數(shù)條件下,進一步分析不同井約束條件下反演的煤層厚度和分布情況。
(5) 獲得反演參數(shù)和井約束,選取最佳方案。
楔狀煤層縱波阻抗數(shù)值模型如圖1 所示,模型尺寸設(shè)置為101 條線和101 個數(shù)據(jù)道,長度為1 000 m×1 000 m(線道距為10 m×10 m),總時深設(shè)置為1 s,采樣時間間隔為1 ms,由中部楔狀煤層和上部砂巖、下部泥巖構(gòu)成,煤層的厚度從1 ~25 m均勻變化。 砂巖、泥巖和煤層的密度、縱波傳播速度設(shè)置見表1。 圖1 中不同的顏色代表對應(yīng)巖性的縱波阻抗值。 利用縱波阻抗模型正演方法進行正演,正演過程中通過縱波阻抗與子波進行褶積,采用的子波類型為雷克子波,頻率為50 Hz,波長為80 ms,得到的地震數(shù)據(jù)結(jié)果如圖2 所示。
圖2 正演數(shù)據(jù)剖面Fig.2 Forward data section
表1 模型各巖性的密度、縱波速度參數(shù)Table 1 Parameter table of density and p-wave velocity of each rock in the model
圖1 楔狀煤層縱波阻抗模型示意圖Fig.1 Model diagram of p-wave impedance of wedge-shaped coal seam
地質(zhì)統(tǒng)計學反演方法流程如圖3 所示,主要包括數(shù)據(jù)準備、反演參數(shù)設(shè)置、反演3 個階段。 數(shù)據(jù)準備階段中,需要準備優(yōu)化好的三維地震、測井、巖性、子波等數(shù)據(jù)。 反演參數(shù)設(shè)置是影響最終反演結(jié)果的重要環(huán)節(jié)。 本文重點對巖性和縱波阻抗的概率密度、橫向變程以及井約束條件進行了分析。
圖3 地質(zhì)統(tǒng)計學反演流程Fig.3 Geostatistical inversion process
2.2.1 概率密度函數(shù)及變程參數(shù)設(shè)置
基于正演數(shù)據(jù)的地質(zhì)統(tǒng)計學反演中,反演參數(shù)的大小設(shè)定主要針對煤層,其他巖性參數(shù)設(shè)置為表1 縱波阻抗模型參數(shù)的對應(yīng)值,且反演過程中加入井約束。 由于實際地層的密度和聲波的測井數(shù)據(jù)不是定值,一般可以用高斯分布進行擬合,因此在偽井信息統(tǒng)計的基礎(chǔ)上加入高斯白噪聲,使得偽井的縱波阻抗值接近實際分布。 偽井的提取遵循均勻分布原則,從模型平面范圍上均勻選取9 個點進行提取,如圖4 所示。 在實際反演參數(shù)設(shè)置時,一般通過測井、地震數(shù)據(jù)統(tǒng)計擬合得到。 在擬合過程中,由于數(shù)據(jù)自身和人為經(jīng)驗影響,參數(shù)的選取會存在一定不確定性[29]。 為了進一步針對煤層探討反演參數(shù)對反演結(jié)果的影響,設(shè)置了不同概率密度函數(shù)的均值ε 和方差σ2以及橫向變程H進行實驗分析,詳細參數(shù)設(shè)置見表2。 針對概率密度函數(shù),參照表1 中煤的縱波阻抗值設(shè)置3種不同的值。 針對橫向變程,參照阻抗數(shù)值模型尺寸大小設(shè)置4 種不同的值,見表2 中3、7、8、9號實驗。 其中,第3 組實驗參數(shù)與實際加噪后的偽井數(shù)據(jù)一致,理論上第3 組實驗參數(shù)下效果最好。
表2 地質(zhì)統(tǒng)計學反演參數(shù)Table 2 Geostatistical inversion parameter
2.2.2 井約束條件參數(shù)設(shè)置
此外,反演過程采用井約束的數(shù)量和位置條件對反演結(jié)果也有重要影響[30]。 為分析不同數(shù)量和位置條件下井約束對地質(zhì)統(tǒng)計學反演的效果,基于上面第3 組實驗的地質(zhì)統(tǒng)計學反演參數(shù),并加入不同數(shù)量和不同位置的井約束條件進行了分析研究。9 組不同井約束條件的實驗參數(shù)見表3,表中代表約束井位置的編號與圖4(a)對應(yīng)。
表3 不同井約束條件反演參數(shù)Table 3 Parameter of different well constraint conditions inversion
圖4 偽井提取位置和縱波阻抗測井曲線(加噪)Fig.4 Pseudo-well extraction location and p-wave impedance log (imnoised)
表2 不同反演參數(shù)條件下煤層厚度反演及其相對誤差如圖5 和圖6 所示。 其中實驗組1、3、5和實驗組2、4、6 的厚度和相對誤差結(jié)果表明,在其他條件一樣時,均值的增大和減小都會造成最終的反演結(jié)果與實際誤差增大;對比實驗組1、2,實驗組3、4 和實驗組5、6 的厚度和相對誤差可以看出,在其他條件一樣時,方差的選取不當也會造成最終的反演結(jié)果與實際誤差增大;對比實驗組3、7、8、9的厚度和相對誤差可以看出,橫向變程為200 m 的誤差最大,1 200 m 的誤差最小,并隨著變程的增大,反演結(jié)果的相對誤差逐漸降低。
圖5 不同反演參數(shù)條件下煤層厚度地質(zhì)統(tǒng)計學反演平面圖Fig.5 Geostatistical inversion plane of coal seam thickness under different inversion parameters
圖6 不同反演參數(shù)條件下煤層厚度地質(zhì)統(tǒng)計學反演與實際模型厚度的相對誤差平面圖Fig.6 The relative error of the thickness of coal seam thickness and actual model thickness under different inversion parameters
概率密度函數(shù)的均值和方差影響較小,這反映出基于MCMC 法的地質(zhì)統(tǒng)計學反演結(jié)果傾向于地震信息。 在本次實驗條件下,地震數(shù)據(jù)在反演過程中,對概率密度函數(shù)均值和方差的取值具有一定的約束能力,反演結(jié)果趨近地震信息。
橫向變程取值對反演結(jié)果影響較大。 這是由于橫向變程是通過地震信息變差分析得到的,進一步反映地震數(shù)據(jù)在地質(zhì)統(tǒng)計學反演中的重要性。 本次實驗條件下,煤層橫向上是連續(xù)分布的,橫向變差取值越接近實際煤層橫向延展大小,反演結(jié)果越好,但考慮實際煤層受沉積環(huán)境、構(gòu)造等影響造成局部變薄、尖滅現(xiàn)象,實際反演過程中變程取值需要綜合考慮。
表3 不同井約束條件下煤層厚度反演及其相對誤差如圖7 和圖8 所示。 對比實驗1 ~9 結(jié)果可以看出,約束井數(shù)量越多,反演的厚度整體誤差越小,越接近實際值。 約束井的分布對反演結(jié)果存在影響,在約束井附近,反演的厚度相對誤差較小,越近越準確;遠離約束井的位置,準確性降低。 此外,均勻分散選取約束井,可以降低約束井的數(shù)量要求,如實驗2 相比較實驗1,選取的井減少1/3,但由于選取的井分散在整個區(qū)域,實驗2 的反演結(jié)果基本接近實驗1。
圖7 不同井約束條件下煤層厚度地質(zhì)統(tǒng)計學反演平面圖Fig.7 Coal seam thickness geostatistical inversion floor plan under different well constraints
圖8 不同井約束條件下煤層厚度地質(zhì)統(tǒng)計學反演與實際模型厚度相對誤差平面圖Fig.8 The thickness of coal seam thickness geostatistical inversion and actual model thickness relative error plan for different well constraints
綜上分析,在針對煤系地層的地質(zhì)統(tǒng)計學反演中,煤層的概率密度函數(shù)均值應(yīng)接近實際值,方差不易過大,橫向變程應(yīng)考慮實際煤層分布情況進行選取,由于地震波的振幅信息可以一定程度上反映地層的延展情況,可以考慮通過煤層界面的地震數(shù)據(jù)均方根振幅值的平面分布范圍選取合適的橫向變程。 此外,約束井應(yīng)均勻分散在整個研究區(qū)。
以山西省孝義市曙光煤礦井田某工作面區(qū)域為例,基于正演實驗結(jié)果對該研究區(qū)的目標煤系地層,進行地質(zhì)統(tǒng)計學反演并獲取三維巖性模型。 研究區(qū)地層主要有長城系、寒武系、奧陶系、石炭系、二疊系、三疊系及新生界的第三系、第四系,主要目標煤層為二疊系下統(tǒng)山西組的2 號煤層,煤層的厚度為0 ~3.25 m,平均1.30 m,煤層縱波平均速度為2 450 m/s。 研究區(qū)面積約為3.85 km2,疊后地震的采樣間隔為1 ms,線道間距為5 m×5 m,頻帶范圍在10 ~120 Hz,主頻約50 Hz,區(qū)內(nèi)共有5 個分布較均勻的鉆井,其中05 井不參與最終的反演,作為反演結(jié)果的盲井驗證,如圖9 所示。
圖9 研究區(qū)范圍及三維地震數(shù)據(jù)Fig.9 Research area range and three-dimensional seismic data
為保障反演的準確性,前期對地震數(shù)據(jù)、測井數(shù)據(jù)進行了優(yōu)化和校正處理,以確保先驗信息的準確性。 5 個井的測井曲線如圖10 所示,煤層厚度分別為1.51 m、2.08 m、2.2 m、0.98 m、1.24 m,煤層的平均縱波阻抗最低,砂巖最高。 由于煤層縱波阻抗相對周圍巖性較低,在地震數(shù)據(jù)中煤層頂?shù)捉缑姹憩F(xiàn)為強反射軸,通過同相軸追蹤獲取煤層的地震界面,并將其作為約束條件加入到地質(zhì)統(tǒng)計學反演中,進一步提高反演的正確性。
圖10 測井巖性和縱波阻抗曲線Fig.10 Logging lithology and p-wave impedance curves
為消除子波的旁瓣影響,選取的反演范圍為研究區(qū)2 號煤層頂?shù)捉缑嫔舷赂鲾U大20 ms,反演參數(shù)取值見表4。 根據(jù)2 號煤層地震的均方根振幅信息,煤層的橫向變程設(shè)置為600 m,概率密度擬合結(jié)果如圖11 所示。
表4 反演參數(shù)取值Table 4 Value of inversion parameters
圖11 地質(zhì)統(tǒng)計學反演概率密度函數(shù)擬合結(jié)果Fig.11 Fitting results of geostatistical inversion probability density function
2 號煤地質(zhì)統(tǒng)計學反演的厚度如圖12 所示,厚度范圍為0 ~2.6 m。 地質(zhì)統(tǒng)計學反演厚度在測井處的煤層厚度對比如圖13 所示,結(jié)合表5 的井上厚度誤差分析可以看出,地質(zhì)統(tǒng)計學反演的厚度與測井的實際厚度在井點處吻合度較好,反演的相對誤差為1.82% ~12.24% ,05 井(盲井)的誤差最大。 此外,煤層厚度越大,反演的誤差越小。 基于反演結(jié)果對2 號煤及其頂、底板部分巖性進行了三維巖性地質(zhì)體建模(圖14)。 可以看出,地質(zhì)統(tǒng)計學反演在橫向和縱向上都具有很好的分辨能力,能夠精確地刻畫出各巖性厚度及其空間展布,尤其是針對薄煤層也有比較準確的刻畫能力。
圖12 2 號煤地質(zhì)統(tǒng)計學反演的厚度平面圖Fig.12 Thickness plane of no.2 coal from geostatistical inversion
圖13 地質(zhì)統(tǒng)計學反演與測井處的煤層厚度對比剖面Fig.13 Coal seam thickness contrast profile at geostatistical inversion and logging
表5 地質(zhì)統(tǒng)計學反演煤層厚度井點處相對誤差分析Table 5 Analysis of relative error at well point in geostatistical inversion of coal seam thickness
圖14 基于地質(zhì)統(tǒng)計學反演巖性結(jié)果建立的三維地質(zhì)巖性模型Fig.14 A 3d geological lithology model based on geostatistical inversion
(1) 基于楔狀煤層縱波阻抗數(shù)值模型正演數(shù)據(jù)的多組反演實驗結(jié)果對比分析,概率密度函數(shù)、橫向變程以及約束井條件對地質(zhì)統(tǒng)計學反演結(jié)果均有影響,其中,概率密度函數(shù)影響較小,橫向變程和約束井條件影響較大。 井的數(shù)量越多,對全局的約束力和范圍越大,反演的準確性越高。 此外,約束井位置分布均勻,可以適當降低井數(shù)量的要求。
(2) 實際應(yīng)用表明,地質(zhì)統(tǒng)計學反演的厚度與測井數(shù)據(jù)基本一致,相對誤差為1.82% ~12.24%。能夠較好地刻畫出煤系地層的空間展布信息,可以為透明化礦山前期的地質(zhì)模型構(gòu)建提供可靠的巖性體數(shù)據(jù)。