核磁共振與瞬變電磁三維聯(lián)合解釋方法

李貅， 劉文韜， 智慶全， 趙威

1 長安大學(xué)地質(zhì)工程與測繪學(xué)院， 西安 710054 2 中國地質(zhì)科學(xué)院地球物理地球化學(xué)勘查研究所， 河北廊坊 065000

?

核磁共振與瞬變電磁三維聯(lián)合解釋方法

李貅1， 劉文韜1， 智慶全2， 趙威1

1 長安大學(xué)地質(zhì)工程與測繪學(xué)院， 西安 710054 2 中國地質(zhì)科學(xué)院地球物理地球化學(xué)勘查研究所， 河北廊坊 065000

傳統(tǒng)核磁共振地下含水量解釋多采用基于均勻半空間或?qū)訝顚?dǎo)電模型的一維反演，分層給出地下含水信息.然而，這些方法忽略了地下復(fù)雜電阻率分布信息對結(jié)果的影響，也不能很好地反映局部三維含水構(gòu)造.本文從三維電介質(zhì)中核磁共振響應(yīng)的正演理論出發(fā)，提出首先利用瞬變電磁數(shù)據(jù)進(jìn)行基于等效導(dǎo)電平面法的快速電阻率成像，然后將成像結(jié)果作為核磁共振三維反演的電性模型，進(jìn)行聯(lián)合解釋.激發(fā)磁場的分布采用有限元法直接求解，通過引入偽δ源實(shí)現(xiàn)電流源的加載，并強(qiáng)加散度條件排除了三維磁場模擬中“弱解”的影響.針對核磁共振靈敏度矩陣的病態(tài)性和數(shù)據(jù)中存在的干擾信號，提出考慮罰項(xiàng)的非線性擬合目標(biāo)函數(shù)，利用線性化方法進(jìn)行核磁共振反演.模型數(shù)據(jù)表明該方法能較準(zhǔn)確反映地下三維含水構(gòu)造，實(shí)測算例進(jìn)一步證明了方法的有效性.本研究將促使核磁共振方法在巖溶、裂隙水、孤立水體等復(fù)雜水文地質(zhì)條件及隧道、礦井災(zāi)害水源探測等方面得到有效應(yīng)用.

核磁共振； 瞬變電磁； 等效導(dǎo)電平面； 三維聯(lián)合解釋

1 引言

核磁共振(Nuclear Magnetic Resonance，NMR)，作為目前世界上唯一直接探測地下水體的地球物理方法(林君等, 2011;潘玉玲和張昌達(dá), 2000)，因具有直接確定含水體的空間分布和水量的能力，近年來得到了廣泛的應(yīng)用(潘玉玲等, 2000;Chalikakis et al.,2008;Legchenko et al.,2010;Lubczynski and Roy,2005;李振宇等, 2003;林君等, 2011;萬樂和潘玉玲, 1999).目前核磁共振一維正反演理論已比較成熟(戴苗等, 2009;林婷婷等,2013;Legchenkoand Shushakov, 1998; Guillen and Legchenko, 2002; Weichman et al., 2002; 翁愛華等, 2007)，即假設(shè)地下水呈層狀分布，分層給出地下含水信息.對于近似水平層狀分布的地下水，目前解釋方法能取得較好效果，然而對于地下巖溶、裂隙水、孤立水體等局部三維含水構(gòu)造，由于不滿足層狀分布，強(qiáng)行將地下含水構(gòu)造當(dāng)作一維層狀模型來處理往往會產(chǎn)生較大誤差，甚至得不到合理的結(jié)果(陳斌等,2014).為拓展核磁共振找水的應(yīng)用范圍，實(shí)現(xiàn)復(fù)雜水文地質(zhì)環(huán)境下的有效探測，應(yīng)開展NMR技術(shù)的三維正反演研究.

對于這一問題，Warsa和Grandis(2012)對均勻半空間下三維含水構(gòu)造進(jìn)行了NMR正演，Warsa等(2014)又在之前研究基礎(chǔ)上進(jìn)行了二維、三維含水構(gòu)造的反演，模型算例表明三維反演能夠?qū)崿F(xiàn)對三維含水構(gòu)造的較好分辨.此外，國內(nèi)王鵬(2014)也對均勻半空間下三維含水體進(jìn)行了NMR正演，并對幾種常見三維含水構(gòu)造響應(yīng)特征進(jìn)行了一定討論.但上述研究均基于均勻半空間模型，沒有考慮地下電阻率信息的影響，而實(shí)際中地下電阻率信息是影響核磁共振探測效果的一個重要因素(翁愛華等,2007;Braun and Yaramanci, 2011).傳統(tǒng)的NMR數(shù)據(jù)反演方法，通常依據(jù)地質(zhì)情況估計一個等效半空間模型，這導(dǎo)致反演結(jié)果受人為因素影響(萬玲等, 2013)，對于復(fù)雜三維導(dǎo)電模型，這種等效將嚴(yán)重偏離真實(shí)地質(zhì)構(gòu)造.因此，要進(jìn)行核磁共振的精確解釋，必須首先獲得合理的電性模型.

國際上相關(guān)研究起步于2002年，Hertrich和Yaramanci(2002)曾利用模擬退火方法實(shí)現(xiàn)了核磁共振方法和垂直電測深方法的聯(lián)合反演，利用電測深法獲得了電性模型，但電測深法為了獲取與核磁共振相同深度的電阻率一般需要敷設(shè)很長的供電導(dǎo)線，且對地形等工作條件要求較高，不利于野外施工.Braun等(2005)，Braun和Yaramanci(2008)嘗試了利用核磁共振信號的相位信息進(jìn)行反演獲得電性分布的方法，并在理論模型試驗(yàn)上取得了良好的效果.但在實(shí)際工作中，由于儀器和環(huán)境干擾的影響，相位信息往往信噪比更低，利用相位信息的反演結(jié)果可信度不高.為了克服這些困難，有些學(xué)者著手研究了瞬變電磁和核磁共振方法的聯(lián)合解釋技術(shù).

瞬變電磁法(Transient Electromagnetic，TEM)是一種電磁感應(yīng)方法，它利用不接地回線或接地線源發(fā)送一次脈沖電磁場，導(dǎo)電介質(zhì)內(nèi)部受感應(yīng)產(chǎn)生渦旋電流，在一次脈沖電磁場間歇期間，渦旋電流產(chǎn)生的二次磁場不會隨一次場消失而立即消失，利用線圈或接地電極觀測二次磁場，研究其與時間的變化關(guān)系，從而確定介質(zhì)的電性分布結(jié)構(gòu)及空間形態(tài).瞬變電磁法具有場源靈活、方法多樣以及穩(wěn)定高效等優(yōu)點(diǎn)，相比于電測深法而言，瞬變電磁法不需要敷設(shè)較長的導(dǎo)線即可達(dá)到相當(dāng)?shù)奶綔y深度，對工作環(huán)境要求相對較低.Legchenko(2009)研究了核磁共振和瞬變電磁的聯(lián)合探測方法，結(jié)合兩種方法各自的優(yōu)勢，較為準(zhǔn)確地對地下含水層進(jìn)行了圈定.萬玲等(2013)聯(lián)用瞬變電磁數(shù)據(jù)和核磁共振數(shù)據(jù)進(jìn)行了一維聯(lián)合反演，模型測試和實(shí)測數(shù)據(jù)處理結(jié)果表明，聯(lián)合反演方法對于含水層位置和界面的圈定更為準(zhǔn)確.但聯(lián)合反演方法以瞬變電磁法和核磁共振法的正演為基礎(chǔ)，而三維瞬變電磁正演相比核磁共振而言，需要時間遠(yuǎn)遠(yuǎn)為多.因此目前來看，在三維核磁共振解釋中直接進(jìn)行聯(lián)合反演，仍是不現(xiàn)實(shí)的.

為了綜合核磁共振和瞬變電磁法各自的優(yōu)勢，實(shí)現(xiàn)核磁共振和瞬變電磁三維聯(lián)合解釋，并避開瞬變電磁三維正演這一極為耗費(fèi)時間成本的過程，本文提出，在聯(lián)合解釋工作中，可首先依賴基于導(dǎo)電平面法的快速電阻率成像技術(shù)從瞬變電磁數(shù)據(jù)中提取電性信息，以此作為電性模型，利用三維有限元法計算介質(zhì)中各處激發(fā)磁場值，進(jìn)而實(shí)現(xiàn)核磁共振三維正演；然后進(jìn)行核磁共振三維反演，求取含水率分布；考慮靈敏度矩陣的奇異性和數(shù)據(jù)包含誤差，引入帶有罰項(xiàng)的最小二乘目標(biāo)函數(shù).

2 三維電介質(zhì)下NMR正演

2.1 NMR信號初始振幅的計算

(1)

(2)其中，M0為氫原子核磁化強(qiáng)度；γ=2.675×108s-1T-1為氫質(zhì)子磁旋比；n(r)為空間r(x,y,z)位置處的含水率；ρ(r)為空間r(x,y,z)位置處電阻率；q=I0·τ為發(fā)射脈沖矩；K3D(q,ρ(r),r)為三維核函數(shù)，代表接收線圈對地下r(x,y,z)處含水體的靈敏度大小；B⊥為激發(fā)磁場BT垂直于地磁場B0的分量；E0為NMR信號初始振幅，其大小與被激發(fā)的地下氫質(zhì)子數(shù)量成正比.

建立如圖1a直角坐標(biāo)系，X軸指向地理正北方向，Y軸指向地理正東，Z軸垂直向下，D為地磁偏角，I為地磁傾角，則B⊥可利用BT進(jìn)行兩次坐標(biāo)系旋轉(zhuǎn)的方法求得(林君等，2013).記BT在新坐標(biāo)系下為

B″=B″xex+B″yey+B″zez,

(3)

B″與BT關(guān)系可用旋轉(zhuǎn)矩陣RI和RD表示

B″=RIRDBT,

(4)

(5)此時X方向與地磁場方向重合，于是B⊥可表示為

(6)

(7)其中，Δxi,Δyj,Δzk分別為點(diǎn)r(x,y,z)位置剖分單元在X,Y,Z方向上長度；N,M,L分別為X,Y,Z方向剖分單元個數(shù).

選擇一系列遞增的脈沖矩q，觀測初始振幅E0隨脈沖矩q的變化，可反映由淺到深不同深度含水量分布.根據(jù)(1)—(7)式可知，在觀測參數(shù)確定的情況下，計算響應(yīng)信號初始振幅E0的關(guān)鍵是首先獲得地下激發(fā)磁場的分布.

2.2 激發(fā)磁場分布的求取

在僅存在外加電性源、時諧情況下，假定諧變因子為eiω t，電磁場所滿足的Maxwell方程組為

(8)

式中，ε為介質(zhì)的介電常數(shù)，σ是電導(dǎo)率，μ是磁導(dǎo)率，ω為激勵電性源的圓頻率，i為虛數(shù)單位，t為時間，E和H分別為電場和磁場矢量，Js為外加電流源的電流密度.利用節(jié)點(diǎn)有限元法求解(8)式時，為了規(guī)避計算電流密度的旋度或電場法向不連續(xù)的問題，常先計算電場或矢量勢，然后通過求導(dǎo)來獲取磁場值.然而，由于求導(dǎo)運(yùn)算，計算得出的磁場值將比直接用有限元求解的量精度低.本文擬直接利用有限元法計算磁場值，并引入偽δ源來解決電流密度的旋度計算問題.由于不存在磁荷和面電流密度，磁場在內(nèi)部電性界面上，同時具有法向和切向連續(xù)性.為通用起見，在外邊界S2上給定磁場的第三類邊界條件為

(9)

(10)

其中，k2=-iωμσ.

根據(jù)變分原理，通過求解泛函(11)的駐點(diǎn)，能夠獲得滿足方程(10)和邊界條件(9)的磁場解

(11)

(12)

對求解區(qū)域進(jìn)行剖分以離散(12)式，假定剖分單元內(nèi)部物性均勻，并設(shè)單元節(jié)點(diǎn)i上的磁場三個分量分別為Hxi,Hyi,Hzi，則對剖分后的所有單元，求解變分問題(12)的駐點(diǎn)等價于求解方程組

(13)

其中，n為剖分單元個數(shù)，e為單元編號.

(14)令R=rot(J)，于是在每個小單元內(nèi)有線性方程組

圖1 坐標(biāo)系示意圖(a)坐標(biāo)系與地理方位關(guān)系示意圖；(b)網(wǎng)格剖分示意圖.Fig.1 Diagram of coordinates(a) Relationship between coordinates and the geographic orientation; (b) Mesh subdivision.

圖2 均勻半空間模型激發(fā)磁場強(qiáng)度有限元解和數(shù)字濾波解比較(a) Hx實(shí)部; (b) Hz實(shí)部;(c) Hx虛部;(d) Hz虛部.Fig.2 Comparison ofFEM solutions and digital filtering solutions of excitation magnetic field for homogenous half space model(a) Real of Hx; (b) Real of Hz; (c) Imag of Hx; (d) Imag of Hz.

(15)

源的加載采用如(16)式Herrmann(1979)提出的偽δ函數(shù)來等效場源的作用，將源分布于一定范圍內(nèi)，避免在源點(diǎn)處產(chǎn)生奇異性，直接進(jìn)行總場的求解.細(xì)導(dǎo)線源的效應(yīng)，可利用偽δ函數(shù)進(jìn)行等效.如空間(x0,y0,z0)處電流強(qiáng)度為I、沿I方向的導(dǎo)線源，其電流密度的分布可寫為Jy=Iδs(x-x0)δs(z-z0).式(15)中對電流密度求旋度的運(yùn)算，可轉(zhuǎn)換為對偽δ函數(shù)的導(dǎo)數(shù)運(yùn)算.將四條等效導(dǎo)線源沿著一致電流方向圍成一個方形，即完成對NMR工作中常用的方形回線源加載.

(16)

式中，τ是控制源分布寬度和幅值的參數(shù)，x為某一空間維度的坐標(biāo)，x0為源所在位置即奇異點(diǎn)的坐標(biāo)值.

采用上述方法，計算了邊長100m發(fā)射線圈在電阻率為100Ωm的均勻半空間上產(chǎn)生的磁場，并在地表沿x軸布置一條測線，與數(shù)字濾波解進(jìn)行對比，計算結(jié)果如圖2所示，可見除導(dǎo)線附近位置，有限元計算結(jié)果與數(shù)字濾波解基本一致.最后，激發(fā)磁場(磁感應(yīng)強(qiáng)度)可采用(17)式求得

BT=μH.

(17)

3 介質(zhì)導(dǎo)電性對核函數(shù)的影響分析

核磁共振響應(yīng)核函數(shù)代表接收線圈對地下含水體靈敏度的大小，為研究介質(zhì)電阻率對核函數(shù)的影響，本文對核磁共振響應(yīng)三維核函數(shù)進(jìn)行了仿真.采用100m單匝方形同一線圈，源中心位于地面(0,0)處，發(fā)射脈沖矩范圍0.05～15A·s，地磁場強(qiáng)度48000nT，地磁傾角60°N，地磁偏角0°，拉莫爾頻率2044Hz.首先分別對電阻率為5、10、50、100、300Ωm和1000Ωm均勻半空間模型進(jìn)行計算，分別記為a—f模型，圖3所示為其三維核函數(shù)仿真結(jié)果，圖中僅列出一個脈沖矩(q=2.5 A·s)的計算結(jié)果.從圖中可以看出介質(zhì)電阻率大小直接影響核函數(shù)分布，介質(zhì)電阻率減小時，相同脈沖矩所能激發(fā)的影響范圍相應(yīng)減小，這是因?yàn)殡姶挪ㄔ诘碗娮杪式橘|(zhì)中衰減快，電磁波能量難以穿透到深部.在當(dāng)前線圈邊長(l=100 m)下，對于電阻率小于50 Ωm低阻介質(zhì)，其對核函數(shù)分布的影響尤為明顯；而后隨著介質(zhì)電阻率增大，其對核函數(shù)的影響逐漸減弱；當(dāng)電阻率大于100 Ωm時，其對核函數(shù)分布的影響從圖上幾乎難以分辨.

這種特性也直觀體現(xiàn)在核磁共振信號初始振幅E0上.假設(shè)地下有一水平方向50 m×50 m，厚20 m的含水體，含水體中心位于(0，0，30 m)處，含水率30%，圖4a所示即為不同介質(zhì)電阻率下同一線圈采集得到的核磁共振初始振幅E0的仿真結(jié)果.當(dāng)均勻半空間電阻率小于50 Ωm時，不同介質(zhì)電阻率下初始振幅信號差別明顯，隨著介質(zhì)電阻率減小，初始振幅E0最大值減小明顯，最大值出現(xiàn)對應(yīng)脈沖矩q逐漸增大；隨著介質(zhì)電阻率增大，介質(zhì)電阻率對初始振幅信號的影響逐漸減小，當(dāng)介質(zhì)電阻率達(dá)300 Ωm，與1000 Ωm情況進(jìn)行對比，相同脈沖矩對應(yīng)的初始振幅之差最大值不足1.5 nV.

可以看出，低阻背景對核磁共振信號影響明顯，為進(jìn)一步討論局部導(dǎo)電介質(zhì)對核磁共振響應(yīng)核函數(shù)的影響大小，在100 Ωm均勻半空間介質(zhì)中，設(shè)計三種不同埋深和規(guī)模的低阻導(dǎo)電體模型，低阻體電阻率均為5 Ωm，分別記為g—i模型.模型g：低阻導(dǎo)電體尺寸：長50 m，寬50 m，厚20 m；導(dǎo)電體中心位置：(0,0,30 m).模型h：低阻導(dǎo)電體尺寸：長20 m，寬20 m，厚10 m；導(dǎo)電體中心位置：(0,0,30 m).模型i：低阻導(dǎo)電體尺寸：長50 m，寬50 m，厚20 m；導(dǎo)電體中心位置：(0,0,70 m).圖5所示為其三維核函數(shù)仿真結(jié)果，為便于比較也列出脈沖矩q=2.5 A·s的計算結(jié)果，與100 Ωm均勻半空間模型對比可以看出，對于淺層低阻異常體，當(dāng)異常體的規(guī)模較大(橫向尺寸為線圈邊長1/2)時，其對核函數(shù)的影響較大，而當(dāng)異常體的規(guī)模較小(橫向尺寸為線圈邊長1/5)時，其對核函數(shù)的影響較小.對于深層低阻異常體，即使異常體較大規(guī)模(橫向尺寸為線圈邊長1/2)，其對核函數(shù)的影響也較小.采用上述相同的含水體，對這三種局部異常體模型的初始振幅E0進(jìn)行仿真計算，結(jié)果如圖4b，與100 Ωm均勻半空間模型進(jìn)行對比，同樣的特征也體現(xiàn)在初始振幅上，淺層規(guī)模較大的局部低阻異常體對初始振幅信號影響顯著，而規(guī)模較小或位于深部的低阻異常體對初始振幅信號影響較小，這是因?yàn)楹舜殴舱裥盘柺莵碜灶l域激發(fā)磁場對于宏觀磁矩的“扳倒”作用，而當(dāng)局部異常體的規(guī)模較小或埋藏較深時，對于頻域激發(fā)磁場的影響不大.綜上，核磁共振響應(yīng)對低阻背景和規(guī)模較大的淺層低阻異常體表現(xiàn)較為敏感，而對高阻背景、規(guī)模較小局部低阻異常體或深部的局部異常體表現(xiàn)不敏感.

圖3 均勻半空間三維核函數(shù)分布 (X-Z方向切片)介質(zhì)電阻率為 a-5 Ωm； b-10 Ωm；c-50 Ωm； d-100 Ωm； e-300 Ωm； f-1000 Ωm.Fig.3 3D kernel function distribution for homogenous half space model (X-Z slice)Resistivity a-5 Ωm； b-10 Ωm； c-50 Ωm； d-100 Ωm； e-300 Ωm； f-1000 Ωm.

圖4 不同介質(zhì)電阻率下E0-q曲線(含水體參數(shù)：長50 m,寬50 m, 厚20 m, 含水率30%，中心點(diǎn)位置(0,0,30 m))(a) 不同均勻半空間模型計算結(jié)果對比； (b) 100 Ωm背景介質(zhì)中不同尺寸低阻異常體模型計算結(jié)果對比.Fig.4 E0-q curves of different resistivity(water-bearing parameters: length 50 m, width 50 m, thickness 20 m, water content 30%, central position (0,0,30 m))(a)Results of different homogenous half space models；(b) Results of differentsizes of low resistivity anomaly under 100 Ωm resistivity background.

圖5 不同尺寸低阻異常體模型三維核函數(shù)分布(X-Z方向切片，背景電阻率100 Ωm，低阻體電阻率5 Ωm)g低阻體尺寸：長50 m,寬50 m,厚20 m；異常體中心：(0,0,30 m);h低阻體尺寸：長20 m,寬20 m,厚10 m；異常體中心：(0,0,30 m); i低阻體尺寸：長50 m,寬50 m,厚20 m;異常體中心：(0,0,70 m).Fig.5 3-D kernel function distribution for different sizes of low resistivity anomaly (X-Z slice, background resistivity 100 Ωm, anomalous resistivity 5 Ωm)g-Anomalous model size: length 50 m, width 50 m, thickness 20 m; central position (0,0,30 m); h-Anomalous model size: length 20 m, width 20 m, thickness 10 m; central position (0,0,30 m); i-Anomalous model size: length 50 m, width 50 m, thickness 20 m; central position (0,0,70 m).

4 NMR-TEM聯(lián)合解釋方法

前節(jié)指出，低阻背景和淺部規(guī)模較大的低阻異常體對核磁共振響應(yīng)核函數(shù)分布較大影響，核磁共振解釋中應(yīng)考慮介質(zhì)的電性分布的影響.獲取地下介質(zhì)的電性分布可依靠電磁法反演技術(shù)，但目前而言電磁法三維正演均需較長的時間，若在核磁共振解釋中引入三維電磁法反演，將嚴(yán)重影響解釋工作的效率.考慮核磁共振響應(yīng)對規(guī)模不大的局部導(dǎo)電異常體不敏感，并無必要在核磁共振解釋中獲取完全精確的電性分布.本文采用瞬變電磁法對勘探區(qū)域進(jìn)行電性分布探測，利用基于導(dǎo)電平面的快速電阻率成像技術(shù)近似地獲取勘探區(qū)域的電性分布，然后以此作為電性模型進(jìn)行核磁共振三維反演，獲得地下含水量分布信息.

4.1 等效導(dǎo)電平面快速成像技術(shù)基本原理

等效導(dǎo)電平面法是一種瞬變電磁探測的近似解釋方法，它是根據(jù)視縱向電導(dǎo)曲線進(jìn)行快速、直觀電性介質(zhì)成像的(Tartaras et al., 2000).在瞬變電磁法中，信號接收裝置是由多匝導(dǎo)線繞成的線圈或探頭，接收到的是二次磁場隨時間的變化率.我們將某一時刻觀測的二次場響應(yīng)信號等同地認(rèn)為是由某一深度的一個虛擬源(導(dǎo)電薄板)產(chǎn)生的信號，另一時刻觀測到的感應(yīng)信號又可以等同地認(rèn)為是另外一個深度的另外一個虛擬源所產(chǎn)生的信號，這樣，不同時間延遲的觀測數(shù)據(jù)就會對應(yīng)不同深度的導(dǎo)電薄板.依據(jù)等效導(dǎo)電平面理論有

(18)

其中，S為縱向電導(dǎo)，K為與發(fā)射接收有關(guān)的常數(shù)，I為發(fā)射電流大小，a為發(fā)射回線半徑，F(xiàn)為一個與時間和電導(dǎo)都有關(guān)的函數(shù).F的表達(dá)式為

(19)

h為薄板深度，t為接收時刻.縱向電導(dǎo)S可依(20)式求取

(20)

顯見，(20)式是關(guān)于縱向電導(dǎo)S的非線性方程，求解該方程即得縱向電導(dǎo)S.根據(jù)縱向電導(dǎo)的定義，可得電阻率值ρτ(h)=dh/dS，對應(yīng)深度為

4.2 核磁共振三維反演方法

如圖1b，當(dāng)按照設(shè)定的研究尺度，將研究空間離散化為網(wǎng)格單元，各網(wǎng)格內(nèi)部物性均勻分布時，核磁共振初始振幅計算公式可用(7)式表示，為便于書寫,(7)式也可簡寫為：

(21)

其中，N為剖分單元的個數(shù)，Kj為第j個單元的核函數(shù)，nj為第j個單元的含水率值.

在多個脈沖矩激發(fā)下，(21)式可以矩陣形式表示：

K·n=E0,

(22)

φ(n)=‖Kn-E‖L2+α‖n‖L2,

(23)

其中，‖‖L2為二范數(shù)算子，α為非負(fù)正則化參數(shù).使目標(biāo)函數(shù)(23)最小化即可獲得矩陣方程(22)在最小范數(shù)罰項(xiàng)下的最優(yōu)解.目標(biāo)函數(shù)極小值本文采用信賴域方法進(jìn)行求解，反演初值依據(jù)電性模型假設(shè)電阻率與含水率線性反相關(guān)近似給出.

4.3 TEM-NMR聯(lián)合解釋流程

聯(lián)合解釋流程框圖如圖6所示.其具體操作步驟如下：

(1) 導(dǎo)入NMR、TEM實(shí)測數(shù)據(jù)及工作參數(shù)裝置；

(2) 對TEM數(shù)據(jù)采用等效導(dǎo)電平面法獲得電性模型，并在電阻率與含水率的線性反相關(guān)的近似假設(shè)下，給出含水率反演初始模型；

(3) 基于獲得的電性模型，采用三維有限元方法計算B⊥，從而實(shí)現(xiàn)NMR的正演計算；

(4) 對目標(biāo)函數(shù)采用信賴域方法進(jìn)行迭代反演，直到目標(biāo)函數(shù)值小于給定值或達(dá)到迭代次數(shù)，輸出反演結(jié)果.

圖6 TEM-NMR聯(lián)合解釋流程框圖Fig.6 Flow of TEM-NMR joint interpretation

5 仿真模型算例

5.1 均勻半空間電性模型

為驗(yàn)證核磁共振三維反演方法的有效性，首先對均勻半空間電性模型下核磁共振數(shù)據(jù)進(jìn)行核磁共振單獨(dú)反演. NMR數(shù)據(jù)在如下條件下獲得：邊長100 m的單匝方形線圈，采用收發(fā)共線圈方式；均勻半空間電阻率200 Ωm，地磁場強(qiáng)度48000 nT，地磁傾角60°N，地磁偏角0°，拉莫爾頻率2044 Hz，最大激發(fā)脈沖矩15 A·s；含水體模型(圖7a)參數(shù)：長60 m，寬60 m，厚25 m，中心位于(0，0,27.5 m)，含水率為30%；地面布置測點(diǎn)25個，測點(diǎn)沿X，Y方向間距均為50 m，對X為-100～100 m，Y為-100～100 m的范圍進(jìn)行核磁共振探測，仿真初始振幅曲線如圖7e中實(shí)線(以X=0 m測線上5個點(diǎn)為例).

在XYZ方向以10 m×10 m×5 m網(wǎng)格間隔對探測區(qū)域進(jìn)行離散，分別以200 Ωm，50 Ωm，20 Ωm的半空間電性模型，采用NMR數(shù)據(jù)進(jìn)行單獨(dú)反演，結(jié)果分別如圖7b、7c、7d所示.可見在反演所采用的電性模型與真實(shí)電性情況較相符時，反演結(jié)果與實(shí)際含水體邊界和含水量吻合較好(圖7b)，說明在給定合理電性模型前提下，本文核磁共振三維反演方法可以取得較好的反演效果.另一方面，從圖7c、7d反演結(jié)果可以看出，電性模型影響反演結(jié)果的準(zhǔn)確度，反演時電性模型與真實(shí)電性情況相差越大，反演結(jié)果與真實(shí)含水邊界及水量偏離越大.

5.2 單個含水體模型

考慮含水體的電性差異，把5.1節(jié)含水體電阻率改為10 Ωm，保持其他模型參數(shù)和觀測參數(shù)不變，核磁共振初始振幅如圖8a中實(shí)線所示，以X=0 m測線上5個點(diǎn)為例；瞬變電磁采用邊長30 m的方形發(fā)射線圈，中心回線工作方式，地面布置測點(diǎn)49個，測點(diǎn)X，Y方向間距均為30 m，對X為-90～90 m，Y為-90～90 m的范圍進(jìn)行觀測，瞬變電磁仿真數(shù)據(jù)采用三維時域有限差分方法(孫懷鳳等，2013)獲得.對瞬變電磁仿真數(shù)據(jù)進(jìn)行等效導(dǎo)電平面成像，結(jié)果如圖8a所示，與含水體模型(圖7a)對比可以看出成像結(jié)果基本反映出低阻含水體輪廓且視電阻率值與真實(shí)電阻率值較為接近，可為核磁共振三維反演提供合理的電性模型.

在XYZ方向以10 m×10 m×5 m網(wǎng)格間隔對探測區(qū)域進(jìn)行離散.首先，采用均勻半空間電性模型進(jìn)行核磁共振單獨(dú)反演，電性模型給定為200 Ωm的半空間，反演擬合曲線如圖8d虛線所示，可見反演含水率計算的初始振幅與真實(shí)初始振幅擬合得很好，但從含水率反演結(jié)果來看(圖8b)，雖基本反映出真實(shí)含水體輪廓，但含水體位置偏深，且淺部出現(xiàn)一個假含水薄層，局部含水率與真實(shí)含水率偏差較大.

對于相同的離散網(wǎng)格，采用電阻率成像結(jié)果作為電性模型進(jìn)行核磁共振反演，反演擬合曲線(圖8d中點(diǎn)劃線)與真實(shí)曲線同樣吻合很好，含水率反演結(jié)果來看(圖8c)，反演含水體位置與真實(shí)含水體基本一致，與直接采用半空間作為電性模型進(jìn)行反演的結(jié)果對比，可以看出反演精度有較大改善.

5.3 多個含水體模型

為進(jìn)一步驗(yàn)證聯(lián)合解釋方法進(jìn)行含水率三維解

圖7 均勻半空間模型核磁共振單獨(dú)反演結(jié)果(a)含水體模型；(b)200 Ωm半空間反演電性模型下反演結(jié)果；(c)50 Ωm半空間反演電性模型下反演結(jié)果；(d)20 Ωm半空間反演電性模型下反演結(jié)果；(e)初始振幅曲線對比.Fig.7 Single inversion results of NMR data based on homogenous half space model(a) Water-bearing model； (b) Inversion resultsunder the half space electric model of 200 Ωm； (c) Inversion results under the half space electric model of 50 Ωm； (d) Inversion results under the half space electric model of 20 Ωm； (e) Initial amplitude curves comparison.

圖8 單個含水體模型及其反演結(jié)果(a) 電阻率成像結(jié)果; (b) 半空間電性模型含水率反演結(jié)果; (c) 成像電阻率模型含水率反演結(jié)果; (d) 初始振幅曲線對比.Fig.8 Modeling of single water bearing body and inversion results(a) Resistivity image result; (b) Inversion result of water content based on homogenous half space model; (c) Inversion result of water content based on resistivity image model; (d) Initial amplitude curvescomparison.

圖9 多個含水體模型設(shè)定與反演結(jié)果(a) 含水體模型圖; (b) 電阻率成像結(jié)果; (c) 成像電阻率模型含水率反演結(jié)果; (d) 初始振幅曲線對比.Fig.9 Modeling of multi-water bearing bodies and inversion results(a) Water-bearing model; (b) Resistivity image result; (c) Inversion result of water content based on resistivity image model; (d) Initial amplitude curvescomparison.

圖10 (a)電阻率成像結(jié)果;(b)核磁共振反演結(jié)果Fig.10 (a) Result of resistivity imaging; (b) Result of NMR inversion

釋的效果，設(shè)計一個更為復(fù)雜的地下含水體模型(圖9a)，模型參數(shù)如表1，地磁場參數(shù)、探測方式及參數(shù)保持不變，瞬變電磁仿真數(shù)據(jù)等效導(dǎo)電平面成像結(jié)果如圖9b所示，從圖中可以看出視電阻率對含水體Ⅰ刻畫最清楚，含水體Ⅱ也有較明顯的反映，而深部含水體在圖中反映不明顯，這說明瞬變電磁法隨著目標(biāo)深度增加探測能力逐漸減弱.然而，第3節(jié)中指出局部深層低阻體對核磁共振響應(yīng)核函數(shù)的影響很小，對核磁共振初始振幅信號的影響幾可忽略，因此該成像結(jié)果仍可作為核磁共振反演合理的電性模型.

以成像結(jié)果作為電性模型進(jìn)行核磁共振初始振幅反演結(jié)果如圖9c所示，含水體Ⅰ的反演結(jié)果能夠準(zhǔn)確反映原始模型，說明合理的電性模型可明顯提高反演精度.受電性模型影響，含水體Ⅱ和含水體Ⅲ反演效果比含水體Ⅰ稍差，但也較準(zhǔn)確反映出真實(shí)含水體位置和水量，且在垂向上兩個含水體能清晰可辨，再次說明了聯(lián)合解釋方法的有效性.

表1 多個含水體模型Table 1 Modeling of multi-water bearing bodies

6 實(shí)測資料試驗(yàn)

為驗(yàn)證本文提出聯(lián)合解釋技術(shù)對于實(shí)測數(shù)據(jù)的有效性，我們以雅礱江錦屏水電站輔助涵洞隧道的探測結(jié)果為實(shí)例進(jìn)行模型測試.錦屏水電樞紐工程輔助洞位于四川省涼山彝族自治州的木里、鹽源、冕寧三縣交界處的雅礱江干流錦屏大河灣上，地形起伏較大.為探測隧道開挖前方充水裂隙的存在，采用中心回線裝置進(jìn)行了瞬變電磁探測，并利用等效導(dǎo)電平面成像技術(shù)進(jìn)行解釋，如圖10a所示(Z軸繪圖縮放比3∶1).根據(jù)成像結(jié)果可以判斷，在掌子面前方20～50 m處范圍內(nèi)存在低阻異常，可能為充水裂隙發(fā)育帶.以此電性模型為基礎(chǔ)進(jìn)行核磁共振反演，并勾繪含水率為10%的等值面如圖10b.由圖可見，核磁共振反演的充水帶范圍和低電阻率異常具有較好的對應(yīng)關(guān)系.同時表明含水率較高的充水帶位于35～45 m范圍，而20～30 m的低阻異常充水率并不高，可能僅僅是局部電性異常的反映.后期的掘進(jìn)結(jié)果表明，聯(lián)合解釋方法對于含水帶的圈定是準(zhǔn)確的.

7 結(jié)論

(1) 核磁共振響應(yīng)的計算，需要預(yù)先求得頻域激發(fā)磁場的分布.本文從諧變磁場所滿足的變分問題出發(fā)，通過求變分的駐點(diǎn)得到了有限元方程.引入偽δ解決了細(xì)導(dǎo)線中電流密度的旋度計算問題，實(shí)現(xiàn)了電流源的直接加載，避開了二次場算法.對于三維電磁場有限元求解中可能存在的“偽解”問題，采用了在變分中強(qiáng)加散度罰項(xiàng)的方法，以減小“偽解”對計算結(jié)果的影響.與數(shù)字濾波解對比表明，文中三維頻率域正演方法具有較高的精度.

(2) 核磁共振響應(yīng)核函數(shù)代表接收線圈對地下含水體靈敏度大小，通過模型計算研究了介質(zhì)電阻率對核磁共振響應(yīng)核函數(shù)分布的影響規(guī)律.計算結(jié)果表明，核磁共振響應(yīng)核函數(shù)對高阻介質(zhì)不敏感，但低阻背景或淺部較大尺度的低阻體對核磁共振響應(yīng)核函數(shù)的影響是顯著的.因此，進(jìn)行核磁共振解釋工作時，應(yīng)考慮電性分布的影響.

(3) 在分析介質(zhì)電阻率對核磁共振響應(yīng)核函數(shù)影響的基礎(chǔ)上，提出了利用瞬變電磁信號進(jìn)行電阻率成像并將其作為核磁共振三維反演的電性模型的聯(lián)合解釋方法.考慮到實(shí)際核磁共振數(shù)據(jù)包含干擾信號和空間離散程度帶來的誤差，采用包含罰項(xiàng)的最小二乘目標(biāo)函數(shù)，利用線性化方法進(jìn)行核磁共振反演.

(4) 利用仿真模型和實(shí)際觀測數(shù)據(jù)進(jìn)行了瞬變電磁和核磁共振的聯(lián)合解釋，均勻半空間模型驗(yàn)證了核磁共振三維反演方法的有效性，單個含水體模型與多個含水體模型的解釋結(jié)果均表明，聯(lián)合解釋方法能更好地反映出地下含水率分布，相比傳統(tǒng)單獨(dú)核磁共振反演優(yōu)勢明顯.實(shí)測資料測試表明本文所提出聯(lián)合解釋方法綜合利用了瞬變電磁法的快速成像和核磁共振方法的直接找水特點(diǎn)，低阻異常帶和含水帶互相印證,實(shí)現(xiàn)了對地下含水體較為精確的定位.

Braun M, Hertrich M, Yaramanci U. 2005. Study on complex inversion of magnetic resonance sounding signals.NearSurfaceGeophysics, 3(3): 155-163, doi: 10.3997/1873-0604.2005011.Braun M, Yaramanci U. 2008. Inversion of resistivity in magnetic resonance sounding.JournalofAppliedGeophysics, 66(3-4): 151-164, doi: 10.1016/j.jappgeo.2007.12.004.

Braun M, Yaramanci U. 2011. Evaluation of the influence of 2-D electrical resistivity on magnetic resonance sounding.JournalofEnvironmental&EngineeringGeophysics, 16(3): 95-103.

Chalikakis K, Nielsen M R, Legchenko A. 2008. MRS applicability for a study of glacial sedimentary aquifers in Central Jutland, Denmark.JournalofAppliedGeophysics, 66(3-4): 176-187, doi: 10.1016/j.jappgeo.2007.11.005.

Chen B, Hu X Y, Liu D H, et al. 2014. The development history and new progress of magnetic resonance sounding technique.ProgressinGeophysics(in Chinese), 29(2): 650-659, doi: 10.6038/pg20140224.

Dai M, Hu X Y, Wu H B, et al. 2009. Inversion of surface nuclear magnetic resonance.ChineseJ.Geophys. (in Chinese), 52(10): 2676-2682, doi: 10.3969/j.issn.0001-5733.2009.10.028.Guillen A, Legchenko A. 2002. Inversion of surface nuclear magnetic resonance data by an adapted Monte Carlo method applied to water resource characterization.JournalofAppliedGeophysics, 50(1-2): 193-205, doi: 10.1016/S0926-9851(02)00139-8.

Herrmann R B. 1979. SH-wave generation by dislocation sources—a numerical study.BulletinoftheSeismologicalSocietyofAmerica, 69(1): 1-15.

Hertrich M, Yaramanci U. 2002. Joint inversion of surface nuclear magnetic resonance and vertical electrical sounding.JournalofAppliedGeophysics, 50(1-2): 179-191, doi: 10.1016/S0926-9851(02)00138-6.

Legchenko A, Ezersky M, Camerlynck C, et al. 2009. Joint use of TEM and MRS methods in a complex geological setting.ComptesRendusGeoscience, 341(10-11): 908-917, doi: 10.1016/j.crte.2009.07.013.

Legchenko A, Vouillamoz J M, Roy J. 2010. Application of the magnetic resonance sounding method to the investigation of aquifers in the presence of magnetic materials.Geophysics, 75(6):L91-L100, doi: 10.1190/1.3494596.

Legchenko A V, Shushakov O A. 1998. Inversion of surface NMR data.Geophysics, 63(1): 75-84, doi: 10.1190/1.444329.

Li Z Y, Pan Y L, Zhang B, et al. 2003. Using NMR method research the hydrogeology problems and practical examples.Hydrogeology&EngineeringGeology(in Chinese), 30(4): 50-54, doi: 10.3969/j.issn.1000-3665.2003.04.011.

Lin J, Duan Q M, Wang Y J, et al. 2011. Theory and Design of Magnetic Resonance Sounding Instrument for Groundwater Detection and Its Applications (in Chinese). Beijing: Science Press.

Lin J, Jiang C D, Lin T T, et al. 2013. Underground magnetic resonance sounding (UMRS) for detection of disastrous water in mining and tunneling.ChineseJ.Geophys. (in Chinese), 56(11): 3619-3628, doi: 10.6038/cjg20131103.

Lin T T, Hui F, Jiang C D, et al. 2013. Layered multi-exponential inversion method on surface magnetic resonance sounding dataset.ChineseJ.Geophys. (in Chinese), 56(8): 2849-2861, doi: 10.6038/cjg20130833.

Lubczynski M, Roy J. 2005. MRS contribution to hydrogeological system parametrization.NearSurfaceGeophysics, 3(3): 131-139, doi: 10.3997/1873-0604.2005009.

Pan Y L, Zhang C D. 2000. The Theory and Methods of Water Detecting (in Chinese). Wuhan: China University of Geosciences Press.Pan Y L, Li Z Y, Wan L, et al. 2000. Detecting bedrock fissure water with nuclear magnetic resonance (NMR) method.CTTheoryandApplications(in Chinese), 9(1): 22-25.

Sun H F, Li X, Li S C, et al. 2013. Three-dimensional FDTD modeling of TEM excited by a loop source considering ramp time.ChineseJ.Geophys. (in Chinese), 56(3): 1049-1064, doi: 10.6038/cjg20130333.

Tartaras E, Zhdanov M S, Wada K, et al. 2000. Fast imaging of TDEM data based on S-inversion.JournalofAppliedGeophysics, 43(1): 15-32, doi: 10.1016/S0926-9851(99)00030-0. Wan L, Pan Y T. 1999. Detecting of karst water in China with nuclear magnetic resonance (NMR) method.CTTheoryandApplications(in Chinese), 8(3): 15-19.

Wan L, Lin T T, Lin J, et al. 2013. Joint inversion of MRS and TEM data based on adaptive genetic algorithm.ChineseJ.Geophys. (in Chinese), 56(11): 3728-3740, doi: 10.6038/cjg20131114.

Wang P. 2014. The MRS response of 3D model in uniformity medium.CTTheoryandApplications(in Chinese), 23(4): 569-578.

Warsa, Grandis H. 2012. Sensitivity study of 3-D modeling for multi-D inversion of surface NMR. ∥ International Conference on Physics and Its Applications: (ICPAP 2011).Bandung, Indonesia: AIP Publishing, 1454: 130-133, doi: 10.1063/1.4730704.

Warsa, Grandis H, Parnadi W W, et al. 2014. Multi-dimensional inversion modeling of surface nuclear magnetic resonance (SNMR) data for groundwater exploration.JournalofEngineeringandTechnologicalSciences, 46(2): 123-140. Weichman P B, Lun D R, Ritzwoller M H, et al. 2002. Study of surface nuclear magnetic resonance inverse problems.JournalofAppliedGeophysics, 50(1-2): 129-147, doi: 10.1016/S0926-9851(02)00135-0.

Weng A H, Wang X Q, Liu G X, et al. 2007. Nonlinear inversion of surface nuclear magnetic resonance over electrically conductive medium.ChineseJ.Geophys. (in Chinese), 50(3): 890-896, doi: 10.3321/j.issn:0001-5733.2007.03.030.

附中文參考文獻(xiàn)陳斌, 胡祥云, 劉道涵等. 2014. 磁共振測深技術(shù)的發(fā)展歷程與新進(jìn)展. 地球物理學(xué)進(jìn)展, 29(2): 650-659, doi: 10.6038/pg20140224. 戴苗, 胡祥云, 吳海波等. 2009. 地面核磁共振找水反演. 地球物理學(xué)報, 52(10): 2676-2682, doi: 10.3969/j.issn.0001-5733.2009.10.028.

李振宇, 潘玉玲, 張兵等. 2003. 利用核磁共振方法研究水文地質(zhì)問題及應(yīng)用實(shí)例. 水文地質(zhì)工程地質(zhì), 30(4): 50-54, doi: 10.3969/j.issn.1000-3665.2003.04.011.

林君, 段清明, 王應(yīng)吉等. 2011. 核磁共振找水儀原理與應(yīng)用. 北京: 科學(xué)出版社. 林君, 蔣川東, 林婷婷等. 2013. 地下工程災(zāi)害水源的磁共振探測研究.地球物理學(xué)報, 56(11): 3619-3628, doi: 10.6038/cjg20131103. 林婷婷, 慧芳, 蔣川東等. 2013. 分層多指數(shù)磁共振弛豫信號反演方法研究. 地球物理學(xué)報, 56(8): 2849-2861, doi: 10.6038/cjg20130833.

潘玉玲, 張昌達(dá). 2000. 地面核磁共振找水理論和方法. 武漢: 中國地質(zhì)大學(xué)出版社.

潘玉玲, 李振宇, 萬樂等. 2000. 利用核磁共振方法探查基巖裂隙水. CT理論與應(yīng)用研究, 9(1): 22-25.

孫懷鳳, 李貅, 李術(shù)才等. 2013. 考慮關(guān)斷時間的回線源激發(fā)TEM三維時域有限差分正演. 地球物理學(xué)報, 56(3): 1049-1064, doi: 10.6038/cjg20130333.

萬樂, 潘玉玲. 1999. 利用核磁共振方法探查巖溶水. CT理論與應(yīng)用研究, 8(3): 15-19.

萬玲, 林婷婷, 林君等. 2013. 基于自適應(yīng)遺傳算法的MRS-TEM聯(lián)合反演方法研究. 地球物理學(xué)報, 56(11): 3728-3740, doi:10.6038/cjg20131114.

王鵬. 2014. 均勻介質(zhì)中MRS方法三維模型的核磁共振響應(yīng). CT理論與應(yīng)用研究, 23(4): 569-578.

翁愛華, 王雪秋, 劉國興等. 2007. 導(dǎo)電性影響的地面核磁共振反演. 地球物理學(xué)報, 50(3): 890-896, doi: 10.3321/j.issn:0001-5733.2007.03.030.

(本文編輯 何燕)

Three-dimensional joint interpretation of nuclear magnetic resonance and transient electromagnetic data

LI Xiu1, LIU Wen-Tao1, ZHI Qing-Quan2, ZHAO Wei1

1SchoolofGeologyEngineeringandGeomatics,Chang′anUniversity,Xi′an710054,China2InstituteofGeophysicalandGeochemicalExploration,ChineseAcademyofGeosciences,HebeiLangfang065000,China

The traditional interpretation methods of nuclear magnetic resonance (NMR) usually use 1D inversion based on homogeneous half space or layered conductive earth model and the information of water amount in each layer will be derived. However, these methods usually ignore the influence of the inhomogeneous underground conductivity to the interpretation results. Meanwhile, these results are easily restricted to local water-bearing structures. For the purpose of avoiding the problems aforementioned, a novel joint interpretation method combining transient electromagnetic (TEM) and NMR data is proposed in this paper, which provides an effective and practical interpretation method for the detection of complex hydrogeological conditions.This paper proposes a novel joint interpretation method combining transient electromagnetic (TEM) and NMR data based on the forward modeling theory of NMR method in 3D conductive media. Firstly, a fast resistivity image is derived from floating plate interpretation results of TEM data. Afterwards, the resistivity imaging data will be used as the electrical model in 3D NMR inversion. The distribution of transmitting magnetic field is calculated by FE method directly. The transmitting current source is loaded by introducing a pseudo-delta source, and the influence of “weak solution” in three dimensional simulation of magnetic field is eliminated by a force introduction of divergence condition. In order to improve the illness of sensitivity matrix and eliminate the noise, a non-linear fitting objective function considering the penalty terms is proposed and the linearized inversion method is used in the inversion process.The results of homogenous half-space model show that solutions of excitation magnetic field are almost the same as FEM solutions and digital filtering solutions and the influence of “weak solution” is eliminated. Results of the distribution of the NMR kernel function in different conductive models show that: (1) The background conductive of a model affects the distribution of the kernel function, especially when it is a low conductive one, and the lower of the conductive, the smaller range of the influence of excitation at the same pulse moment will be; (2) Local low resistive bodies with reasonably scale also have a significant effect on the distribution of the kernel function, lower resistivity and shallower buried depth causing greater influence, vice versa; Several conclusions can be drawn from the inversion results: (1) smaller variation of the geo-electrical model from the real case will lead to more accurate inversion results; (2) The inversion result can deviate from reality when using a homogenous half-space as the electrical model directly while the low conductive body have a certain scale, and the accuracy of the inversion result will see a significant improvement with the introduction of the floating plate interpretation method based on TEM data; (3) Though the accuracy of the floating plate interpretation decreases with increased depth, the impact of low conductive bodies on the distribution of the kernel function also decreases with that. Therefore, the fast resistivity image result based on TEM data, in spite of being an approximate result, can still offer a reasonable electrical model for NMR inversion.This paper proposes a novel joint interpretation method combining transient electromagnetic (TEM) and NMR data based on the forward modeling theory of NMR method in 3D conductive media. Results of the 3D distribution of kernel function show that low resistive bodies have a significant effect on the distribution of the kernel function, and the media conductivity should be considered during the NMR inversion. The floating plate interpretation based on TEM data can provide a reasonable geo-electrical model, and the accuracy of the NMR inversion will be improved significantly after introducing joint interpretation method combining transient electromagnetic (TEM) and NMR data.

Nuclear magnetic resonance; Transient electromagnetic; Equivalent conductive plane; 3D joint interpretation

國家重點(diǎn)基礎(chǔ)研究發(fā)展計劃項(xiàng)目(973項(xiàng)目)(2013CB036002)、國家重大科研設(shè)備研制項(xiàng)目(ZDYZ2012-1-05-04)和國家自然科學(xué)基金重點(diǎn)項(xiàng)目(51139004)聯(lián)合資助.

李貅,男,1958年生,教授,博士生導(dǎo)師,主要從事瞬變電磁場的理論與應(yīng)用方面的研究.E-mail:Lixiu@chd.edu.cn

10.6038/cjg20150810.

10.6038/cjg20150810

P631

2014-12-14，2015-07-08收修定稿

李貅， 劉文韜， 智慶全等. 2015. 核磁共振與瞬變電磁三維聯(lián)合解釋方法.地球物理學(xué)報，58(8):2730-2744,

Li X, Liu W T, Zhi Q Q, et al. 2015. Three-dimensional joint interpretation of nuclear magnetic resonance and transient electromagnetic data.ChineseJ.Geophys. (in Chinese),58(8):2730-2744,doi:10.6038/cjg20150810.