徐鳳姣,謝興兵,郭全仕,周 磊,王新宇

(1.中國石油化工股份有限公司石油物探技術(shù)研究院,江蘇南京211103;2.中國地質(zhì)大學(xué)地球物理與空間信息學(xué)院,湖北武漢430074;3.長江大學(xué)油氣資源與勘探技術(shù)教育部重點實驗室,湖北武漢430100)

地?zé)豳Y源作為新型接替能源,因其具有可再生、節(jié)能及環(huán)保等特點,越來越受到廣泛重視[1]。據(jù)中國大陸地區(qū)大地?zé)崃鲾?shù)據(jù)匯編顯示,我國地?zé)豳Y源十分豐富,具有東高、中低、南高以及西北低等特征[2-3]。其中,深層(3～10km埋深)干熱巖地?zé)豳Y源總量為20.9×106EJ,相當(dāng)于714.9×1012t標(biāo)準(zhǔn)煤,是極具開發(fā)潛力的戰(zhàn)略接替能源[4]。

應(yīng)用于地?zé)豳Y源勘探的地球物理方法以間接勘探方法為主,主要依據(jù)巖石物性參數(shù)與溫度的相關(guān)性開展研究。目前,常用的方法主要有重力法、磁法、電磁法、以及地震方法等[5-9]。大地電磁法(MT)因具有勘探深度大、工作效率高和成本低廉等優(yōu)點,已成為深層地?zé)峥碧阶畛Ｓ玫姆椒ㄖ籟10-13]。常規(guī)MT資料處理解釋主要基于電磁感應(yīng)理論,然而,地下巖石在低頻電磁場作用下,會同時發(fā)生電磁感應(yīng)和激發(fā)極化現(xiàn)象,此時的巖石電阻率是一個隨頻率變化的復(fù)電阻率[14-18]。多年來,國內(nèi)外學(xué)者在考慮激電效應(yīng)的MT正演方面做了很多研究。1978年,PELTON等[19]基于實驗結(jié)果,提出由激電效應(yīng)引起的復(fù)電阻率隨頻率的變化規(guī)律能用Cole-Cole模型表示。隨后,羅延鐘等[20]基于電子導(dǎo)體激發(fā)極化的電化學(xué)機理,通過研究過電位充、放電特性和頻譜特征,建立了面極化和體極化的等效電路和數(shù)學(xué)模型,從而在理論上證明了采用Cole-Cole模型描述激電復(fù)電阻率頻譜的合理性。2011年,朱占升等[21]采用有限單元法,通過MT二維正演計算非極化均勻空間中存在的不同極化體,發(fā)現(xiàn)高阻極化體對視電阻率值影響強于低阻極化體。2016年,符超等[22]利用視電阻率比值與相位比值,詳細討論了大地電磁場激發(fā)極化效應(yīng)與極化體埋深和厚度的關(guān)系。2018年,張志勇等[23]采用超松弛迭代雙共軛梯度法,開展了同時考慮激電效應(yīng)和電磁效應(yīng)的復(fù)電阻率法二維數(shù)值模擬研究。2020年,熊治濤等[24]采用有限元數(shù)值計算方法,利用Galerkin加權(quán)余量法,實現(xiàn)了二維MT各向異性極化介質(zhì)正演,該研究成果對頻率域電磁法野外勘探具有指導(dǎo)意義。

隨著基于激電效應(yīng)的MT正演理論的不斷完善,野外MT實測資料激電參數(shù)提取取得了突破性進展。2006年,陳清禮等[25]基于巖礦石激電響應(yīng)的Dias模型,開展MT資料真譜參數(shù)反演方法研究,成功從MT資料中提取激電參數(shù),為地質(zhì)解釋提供更多依據(jù)。曹中林等[26]在油氣檢測中進行MT激電效應(yīng)的模擬研究,該方法利用巖石電阻率、極化率、時間常數(shù)以及頻率相關(guān)系數(shù)等參數(shù),同時進行MT資料處理解釋,提高了MT方法的分辨能力和可靠性。HE等[27]和羅衛(wèi)峰等[28]在柴達木盆地開展MT實測資料激發(fā)極化信息提取研究,試驗結(jié)果表明,氣田呈高電阻率、高極化率異常模式,該成果對研究區(qū)含油氣性檢測評價起到了指示性作用。董莉等[29]通過改進自適應(yīng)策略,提出一種基于非均勻統(tǒng)計分布的自適應(yīng)差分進化兩階段最小構(gòu)造反演方法,提高了對微弱激電信息提取的精度,降低了MT資料解釋的多解性。截止目前,MT激電勘探主要應(yīng)用于礦產(chǎn)和油氣勘探領(lǐng)域,針對深層地?zé)豳Y源勘探的MT方法僅限于常規(guī)方法,為提高MT方法在深層地?zé)峥碧街械膽?yīng)用效果,有必要對激電條件下MT勘探在深層熱儲中的響應(yīng)特征進行深入分析。

本文采用有限差分?jǐn)?shù)值計算方法,通過引入Cole-Cole復(fù)電阻率模型,完成了不同激電模型的二維MT正演計算。通過分析不同激電參數(shù)對深層熱儲視電阻率與相位響應(yīng)特征,結(jié)合熱儲巖石溫度與激電參數(shù)關(guān)系,為深層地?zé)犭姶趴碧教峁┬路椒ā?/p>

1 正演理論

1.1 電磁場基本理論

當(dāng)平面入射的電磁波在地球介質(zhì)中傳播時,假定近地表大地電磁場為穩(wěn)態(tài)場,取時諧因子為e-iwt,基于本構(gòu)方程,頻率域Maxwell方程組的旋度方程可以表示為:

(1)

式中:E,H分別為電場強度和磁場強度;ω為角頻率;μ為磁導(dǎo)率;ε為介電常數(shù);σ為電導(dǎo)率。

設(shè)定二維模型構(gòu)造走向沿x方向,傾向沿y方向,即?/?x=0,在笛卡爾坐標(biāo)系中,旋度方程能解耦得到如公式(2)所示的TE、TM極化模式電磁場方程組。

(2)

式中:Ex,Hx分別為電場強度和磁場強度沿x方向的分量。

1.2 有限差分?jǐn)?shù)值模擬

本文采用有限差分?jǐn)?shù)值計算方法進行大地電磁二維正演計算,利用差分方法可以將公式(2)中的亥姆赫茲方程組分別轉(zhuǎn)化為TE和TM極化模式的二階偏微分方程,因此,不需要采用交錯網(wǎng)格,而是分別求取Ex,Hx對應(yīng)的二階偏微分方程。構(gòu)建如圖1所示的矩形剖分網(wǎng)格,其中,y方向的節(jié)點編號用j表示,z方向的節(jié)點編號用k表示。

圖1 網(wǎng)格剖分節(jié)點編號

在忽略位移電流項的條件下,采用三點差分公式,將公式(2)中的微分方程組改寫為差分方程組,則獲得節(jié)點P(j,k)處的差分公式為:

(3)

(4)

dzk+ρ(j,k-1)dyjdzk-1+ρ(j-1,k)·

dyj-1dzk+ρ(j-1,k-1)dyj-1dzk-1]

(5a)

(dyjdzk+dyjdzk-1)

(5b)

(dyjdzk+dyj-1dzk)

(5c)

本文進行正演計算時統(tǒng)一采用第一類邊界條件,4個邊界都距離探測目標(biāo)區(qū)足夠遠,不受地下不均勻體的影響,上邊界在距離地表足夠遠(一般大于100km)的高空中,取場值為常數(shù)(通常取1),下邊界在地下足夠深的地方(探測深度的5倍以上),取值為0,左、右邊界條件分別為一維正演結(jié)果,則本文差分方程求解的邊界條件為:

(6)

式中:ExL(k),HxL(k),ExR(k),HxR(k)分別為一維正演得到的左、右邊界電、磁場值;Ny為沿y方向的第N個節(jié)點編號。

聯(lián)立公式(3)、公式(4)和公式(6),則可得到大型線性方程組:

Kx=b

(7)

式中:K為大型稀疏矩陣;x為待求解的電磁場值;b為矢量邊界條件。

求解公式(7)可得到電磁場分量,結(jié)合公式(8)所示的卡尼亞視電阻率公式,即可計算出地表觀測點對應(yīng)的視電阻率與阻抗相位。

(8)

1.3 基于Cole-Cole模型的二維MT正演

本文采用COLE等[30]提出的Cole-Cole復(fù)電阻率模型,開展基于激電效應(yīng)的二維MT正演研究,其中Cole-Cole復(fù)電阻率數(shù)學(xué)表達式如下:

(9)

式中:ρ(iω)為不同頻率條件下巖石的復(fù)電阻率;ρ為不存在極化時的直流電阻率;m為極化率;τ為時間常數(shù);c為頻率相關(guān)系數(shù)。

基于激電模型的二維MT正演計算時,由于激電效應(yīng)的存在,網(wǎng)格單元的電阻率應(yīng)為復(fù)電阻率。因此,將公式(3)、公式(4)中的平均電導(dǎo)率和平均電阻率,替換成公式(9)中的平均復(fù)電導(dǎo)率和平均復(fù)電阻率。通過改變激電參數(shù)即可實現(xiàn)不同激電模型的二維MT正演。

2 深層熱儲激電響應(yīng)特征分析

2.1 基于激電參數(shù)的深層熱儲MT勘探可行性分析

在針對巖石電阻率與溫度關(guān)系的研究中,國內(nèi)外很多學(xué)者認(rèn)為:以離子導(dǎo)電方式為主的地下巖石,隨著溫度的升高,離子活性增強,巖石內(nèi)部電離程度增加,使得導(dǎo)電離子總數(shù)增多,從而電阻率降低。大量實驗研究發(fā)現(xiàn),巖石電阻率與溫度的經(jīng)驗公式[31-33]為:

(10)

式中:ρ為溫度T時的電阻率;ρi為溫度Ti時的電阻率;α為溫度系數(shù),其值與其巖性和地下水溶液礦化度相關(guān),一般取0.02。

將公式(10)變型,即可得到電阻率變化率與溫度變化的關(guān)系式:

(11)

式中:ερ為電阻率相對變化率;ΔT=T-Ti為溫度差。

計算不同溫度差條件下巖石電阻率相對變化率,結(jié)果如表1所示。從表1中可以看出,當(dāng)溫度差達到20℃時,能夠引起28.57%的電阻率相對變化率,說明電磁方法進行地?zé)峥碧骄哂锌尚行浴?/p>

表1 不同溫度差對應(yīng)電阻率相對變化率

由于地下巖石具有頻散效應(yīng),其電阻率是一個隨頻率變化的值,即復(fù)電阻率。本文引入Cole-Cole復(fù)電阻率模型來表征巖石復(fù)電阻率,聯(lián)立公式(9)和公式(10)得到溫度與巖石激電參數(shù)的關(guān)系式:

(12)

公式(12)說明,電阻率、極化率、時間常數(shù)及頻率相關(guān)系數(shù)能夠反映溫度變化,為基于激電效應(yīng)的MT方法進行深層地?zé)峥碧教峁┝死碚撘罁?jù)。

2.2 算法驗證

本文采用3層層狀介質(zhì)模型(K型)進行二維MT有限差分正演算法驗證,層狀介質(zhì)模型參數(shù)如表2 所示。

表2 層狀介質(zhì)模型參數(shù)

對比二維MT有限差分正演結(jié)果與一維MT解析解結(jié)果發(fā)現(xiàn),兩種方法的計算結(jié)果基本一致(圖2),從而驗證了二維MT有限差分方法的正確性。

圖2 二維MT有限差分正演結(jié)果與一維MT解析解結(jié)果對比a 視電阻率結(jié)果; b 相位結(jié)果

2.3 深層熱儲極化異常響應(yīng)特征研究

為研究深層地?zé)釤醿T激電響應(yīng)特征,本文設(shè)置如圖3、表3所示的激電模型。如圖3所示,在電阻率為100Ω·m的均勻半空間中,存在一個電阻率為10Ω·m,埋深為3km,大小為1km×2km的熱儲極化異常體,其激電參數(shù)如表3所示。本文正演計算頻率范圍為10-3～103Hz,按對數(shù)等間距取61個頻點,地表測點范圍為-5～5km,每0.1km一個測點,共101個測點。

圖3 深層地?zé)釤醿る娔Ｐ?/p>

本文設(shè)置3種不同激電參數(shù)模型進行深層熱儲極化異常體響應(yīng)特征研究。

2.3.1 不同極化率的MT響應(yīng)

如表3中異常體模型一所示,在100Ω·m的均勻半空間中,存在一個電阻率為10Ω·m、時間常數(shù)為10-1、頻率相關(guān)系數(shù)為0.4、極化率分別為0,0.2,0.4,0.6,0.8等5種情況的極化異常體?；贑ole-Cole復(fù)電阻率模型,開展不同極化率條件下的二維MT正演,得到如圖4所示的視電阻率和相位擬斷面圖。

表3 熱儲模型激電參數(shù)

圖4a至圖4d自上而下分別顯示了極化率為0,0.2,0.4,0.6,0.8時,TE、TM極化模型下二維MT正演視電阻率和相位響應(yīng)。由圖4可知,地下介質(zhì)的視電阻率和相位隨極化率參數(shù)的增加而變大,其中TM極化模式響應(yīng)明顯強于TE極化模式,且TM極化模式的視電阻率響應(yīng)大于相位響應(yīng)。基于此,開展TM極化模式下視電阻率相對變化率計算,結(jié)果如圖5 所示,研究發(fā)現(xiàn),當(dāng)極化率差大于0.2時,視電阻率相對變化率最大值高于20%,表明因極化率參數(shù)改變而產(chǎn)生的激電響應(yīng)較強,說明該參數(shù)可作為深層地?zé)峥碧降挠行?shù)。

圖4 不同極化率的二維MT正演結(jié)果擬斷面(自上而下極化率分別為0,0.2,0.4,0.6,0.8)a 視電阻率-TE; b 視電阻率-TM; c 相位-TE; d 相位-TM

圖5 TM極化模式下視電阻率相對變化率(Ⅰ)a 極化率差為0.2; b 極化率差為0.4; c 極化率差為0.6; d 極化率差為0.8

2.3.2 不同時間常數(shù)的MT響應(yīng)

設(shè)置如圖3、表3中異常體模型二所示的二維MT正演模型,在100Ω·m的均勻半空間中,存在一個電阻率為10Ω·m、極化率為0.4、頻率相關(guān)系數(shù)為0.4、時間常數(shù)分別為10-2,10-1,100,101,102等5種情況的極化異常體?；贑ole-Cole復(fù)電阻率模型,采用有限差分算法,開展不同時間常數(shù)條件下二維MT正演,得到如圖6所示的視電阻率和相位擬斷面圖。

圖6清晰反映了時間常數(shù)對大地電磁視電阻率和相位有較大影響,隨著時間常數(shù)的增大,視電阻率和相位有明顯的減小趨勢,其中視電阻率的響應(yīng)大于相位的響應(yīng),且TM極化模式對時間常數(shù)的變化相對敏感。計算TM極化模式下視電阻率相對變化率,結(jié)果如圖7所示,可以發(fā)現(xiàn),當(dāng)時間常數(shù)差大于(102-101)時,視電阻率相對變化率最大值大于20%,表明時間常數(shù)改變引起的激電響應(yīng)較弱。

圖6 不同時間常數(shù)的二維MT正演結(jié)果擬斷面(自上而下時間常數(shù)分別為10-2,10-1,100,101,102)a 視電阻率-TE; b 視電阻率-TM; c 相位-TE; d 相位-TM

圖7 TM極化模式下視電阻率相對變化率(Ⅱ)a 時間常數(shù)差為(10-1-10-2); b 時間常數(shù)差為(100-10-1); c 時間常數(shù)差為(101-100); d 時間常數(shù)差為(102-101)

2.3.3 不同頻率相關(guān)系數(shù)的MT響應(yīng)

在笛卡爾直角坐標(biāo)系下,設(shè)計100Ω·m的均勻半空間中,存在一個電阻率為10Ω·m、極化率為0.4、時間常數(shù)為10-1、頻率相關(guān)系數(shù)分別為0.2,0.4,0.6,0.8等4種情況的極化異常體,如圖3、表3中異常體模型三所示?；贛T正演理論,采用有限差分算法,通過引入Cole-Cole復(fù)電阻率模型,得到不同頻率相關(guān)系數(shù)條件下,二維MT正演視電阻率和相位擬斷面圖,結(jié)果如圖8所示。

由圖8可知,大地電磁視電阻率和相位隨著頻率相關(guān)系數(shù)的增大而增大,其中TM極化模式響應(yīng)特征比TE極化模式的響應(yīng)特征明顯,且TM極化模式條件下,視電阻率的變化更敏感。計算TM極化模式下視電阻率相對變化率,結(jié)果如圖9所示,可以發(fā)現(xiàn),頻率相關(guān)系數(shù)差大于0.6時,視電阻率相對變化率最大值為20%,表明頻率相關(guān)系數(shù)改變引起的激電響應(yīng)也相對較弱。

圖8 不同頻率相關(guān)系數(shù)的二維MT正演結(jié)果擬斷面(從上到下頻率相關(guān)系數(shù)分別為0.2,0.4,0.6,0.8)a 視電阻率-TE; b 視電阻率-TM; c 相位-TE; d 相位-TM

圖9 TM極化模式下視電阻率相對變化率(Ⅲ)a 頻率相關(guān)系數(shù)差為0.2; b 頻率相關(guān)系數(shù)差為0.4; c 頻率相關(guān)系數(shù)差為0.6

3 結(jié)論

基于有限差分方法,同時引入Cole-Cole模型,實現(xiàn)了二維激電介質(zhì)的MT正演。通過設(shè)置不同激電參數(shù)(極化率、時間常數(shù)、頻率相關(guān)系數(shù))模型,模擬計算并分析深層熱儲模型的MT視電阻率和相位響應(yīng)特征,得到如下結(jié)論與建議:

1) 研究不同激電模型二維MT響應(yīng)發(fā)現(xiàn),視電阻率和相位隨極化率、頻率相關(guān)系數(shù)的增大而增大,隨時間常數(shù)的增大而減小;

2) 二維極化介質(zhì)條件下,TM極化模式的響應(yīng)特征明顯強于TE極化模式的響應(yīng)特征,且視電阻率與相位相比,對激電參數(shù)的變化更敏感;

3) 計算TM極化模式下視電阻率相對變化率,結(jié)合溫度與電阻率變化率理論關(guān)系,發(fā)現(xiàn)由極化率參數(shù)引起的激電響應(yīng)較強,而由時間常數(shù)和頻率相關(guān)系數(shù)引起的激電響應(yīng)相對較弱。

綜上所述,在實際深層地?zé)岫SMT勘探時,為提高資料處理解釋的精度與速度,建議采用TM極化模式進行資料反演,并采用電阻率和極化率多參數(shù)聯(lián)合解釋方式,以提高MT方法在深層地?zé)峥碧街械膽?yīng)用能力。