大定源瞬變電磁一維自適應(yīng)正則化反演

徐玉聰，趙 寧，秦 策，王 銳，張召彬

(1. 成都理工大學(xué)地球物理學(xué)院，四川成都 610059；2. 河南理工大學(xué)計(jì)算機(jī)學(xué)院，河南焦作 454000)

?

大定源瞬變電磁一維自適應(yīng)正則化反演

徐玉聰1，趙 寧2，秦 策1，王 銳1，張召彬1

(1. 成都理工大學(xué)地球物理學(xué)院，四川成都 610059；2. 河南理工大學(xué)計(jì)算機(jī)學(xué)院，河南焦作 454000)

自適應(yīng)正則化反演對(duì)初始模型要求較低，直接對(duì)瞬變電磁響應(yīng)數(shù)據(jù)進(jìn)行反演。正則化因子通過計(jì)算每次迭代的數(shù)據(jù)目標(biāo)函數(shù)和模型目標(biāo)函數(shù)自適應(yīng)的得到,從而快速獲得地下的地電結(jié)構(gòu)。本文采用自適應(yīng)正則化反演算法對(duì)大定源回線瞬變電磁一維層狀模型進(jìn)行反演，使用均勻半空間作為初始模型，首次采用模型對(duì)深度的二階導(dǎo)數(shù)極小的模型約束，通過典型理論模型的反演計(jì)算，證明了TEM自適應(yīng)正則化一維反演算法擬合效果好、精度高，由于不用反復(fù)搜索正則化因子，并且收斂速度快，體現(xiàn)了良好的穩(wěn)定性和可靠性。

瞬變電磁 大定源 自適應(yīng)正則化反演 最光滑模型約束

Xu Yu-cong, Zhao Ning, Qin Ce1, Wang Rui, Zhang Zhao-bin. One-dimensional adaptive regularization inversion of transient electromagnetic sounding with a large fixed source[J]. Geology and Exploration, 2015, 51(2)：0360-0365.

1 引言

瞬變電磁法是通過對(duì)觀測(cè)到的二次場(chǎng)響應(yīng)信息進(jìn)行提取分析，從而獲得地下介質(zhì)的地電結(jié)構(gòu)。大定源回線裝置發(fā)射線圈一般采用邊長(zhǎng)數(shù)百米的矩形回線,由于布設(shè)一次發(fā)射回線能夠在多個(gè)點(diǎn)同時(shí)接收,因此該方法工作效率相對(duì)其他裝置要高的多，所以廣泛應(yīng)用于許多工程領(lǐng)域(樸化榮，1990；李貅等，2002)。因此需要尋找一種穩(wěn)定、可靠的反演解釋方法。

一維瞬變電磁反演解釋常用有兩種思路,一種是最優(yōu)化反演,一種是煙圈反演。煙圈理論最早是由M.N.Nabighian等提出的，這是一種近似反演方法,反演速度非常快,而且不依賴初始模型,不過由于反演結(jié)果不能提供與層厚相關(guān)的信息,而且精度達(dá)不到要求,因此,煙圈反演往往只用來做定性分析。目前常用最優(yōu)化反演方法是最小二乘，由于最小二乘只是盡量擬合觀測(cè)數(shù)據(jù)，對(duì)反演模型約束不夠，所以難以獲得好的反演結(jié)果；翁愛華(翁愛華等，2007)采用OCCAM反演法直接對(duì)感應(yīng)電動(dòng)勢(shì)進(jìn)行反演，OCCAM方法可以獲得很好的反演效果，由于OCCAM反演(Constable S C,1987)方法每次迭代要多次正演計(jì)算，所以耗時(shí)較長(zhǎng)。

本文結(jié)合大定源回線瞬變電磁的具體特點(diǎn)，采用陳小斌(陳小斌等，2005)自適應(yīng)正則化反演法，實(shí)現(xiàn)了約束模型下的正則化反演，其中正則化因子通過自適應(yīng)的方式給出，避免了OCCAM反演中為搜索正則因子進(jìn)行的大量正演計(jì)算過程，因此通過極少的迭代次數(shù)和計(jì)算時(shí)間即可獲得較好的的反演結(jié)果。

2 自適應(yīng)正則化反演算法原理

正則化反演是基于盡可能擬合觀測(cè)數(shù)據(jù)，同時(shí)獲得最理想的約束模型的一種算法，其基本思路是把地下介質(zhì)分成多層層狀結(jié)構(gòu)，把每層厚度加入模型約束函數(shù)，即粗糙度矩陣，通過計(jì)算得到每層電阻率，從而構(gòu)建地下地電結(jié)構(gòu)。其總目標(biāo)函數(shù)可歸結(jié)為：

(1)

式中φ(m)為總目標(biāo)函數(shù)；λ為正則化因子；φm(m)為模型約束目標(biāo)函數(shù)；φd(m)為觀測(cè)數(shù)據(jù)目標(biāo)函數(shù)；m為模型向量。

其中(1)式中觀測(cè)數(shù)據(jù)目標(biāo)函數(shù)φd(m)可表示為

(2)

式(2)中Δd是觀測(cè)數(shù)據(jù)與理論響應(yīng)差向量；Wd為數(shù)據(jù)加權(quán)矩陣；F為正演算子；J為經(jīng)過泰勒級(jí)數(shù)展開后正演響應(yīng)對(duì)電阻率的偏導(dǎo)數(shù)矩陣，即雅克比矩陣。

模型約束目標(biāo)函數(shù)φm(m)由(3)式給出

(3)

φm(m)中Rm粗糙度矩陣，m為模型向量；目前φm(m)構(gòu)建方式有三種，即最小約束模型、最平緩約束模型、最光滑約束模型；模型約束目標(biāo)函數(shù)的構(gòu)建可完全轉(zhuǎn)化為粗糙度核矩陣Rm的計(jì)算，本文采用最光滑模型約束(模型對(duì)深度的二階導(dǎo)數(shù))。

其構(gòu)建方式如下：Rm=Rn×n，h0、h1…h(huán)n為模型層厚距離地面垂直深度。

(4)

正則化因子采用陳小斌提出的CMD方案：

(5)

式(5)中k表示第k次迭代反演，通過自適應(yīng)的調(diào)節(jié)正則化因子λ，可以平衡觀測(cè)數(shù)據(jù)目標(biāo)函數(shù)和模型約束目標(biāo)函數(shù)的權(quán)重，從而獲得最佳擬合結(jié)果。

將(1)式對(duì)模型向量求偏導(dǎo)數(shù)，并做線性化處理，根據(jù)極小化原則，可得到反演方程：

(6)

式(6)中Δm為待求的模型修正向量，Jk為當(dāng)前模型的雅克比矩陣；解此線性方程，獲得模型修正量，進(jìn)而得到下一次迭代反演的初始模型mk+1，其中mk+1=mk+Δm。

反演的終止條件定義為：

當(dāng)反演擬合差小于給定RMSinit，或者反演迭代達(dá)到最大迭代次數(shù)，即終止反演，獲得局部最優(yōu)的反演模型。

3 反演算例

為檢驗(yàn)上述反演算法的可靠性，利用層狀介質(zhì)模型的理論正演數(shù)據(jù)加4%的隨機(jī)噪聲作為反演數(shù)據(jù)，對(duì)大定源發(fā)射線框內(nèi)不同測(cè)點(diǎn)數(shù)據(jù)進(jìn)行反演。大回線瞬變電磁裝置參數(shù)設(shè)置如下：發(fā)射線圈為正方形線圈，線圈匝數(shù)設(shè)置為1，邊長(zhǎng)200 m, 發(fā)射電流10 A,接收線圈接收面積為1 m2，匝數(shù)為1。裝置布置如圖1，a1、a2處為測(cè)點(diǎn)，三層和四層模型的反演初始模型，均采用100 Ω·m的均勻半空間。

圖1 大定源裝置布置圖Fig.1 Measurement map of large fixed source transient electromagnetic device

3.1 三層模型反演試算

三層模型中最典型的是H型和K型地電模型，圖2所示是H型模型(各層電阻率分別為300 Ω·m、50 Ω·m、300 Ω·m，層厚為100 m、100 m)測(cè)點(diǎn)a1的反演結(jié)果，反演第3次得到的修正模型已經(jīng)基本反映出了200 m～300 m的相對(duì)高阻層，隨著反演迭代次數(shù)的增加逐漸逼近真實(shí)模型的電阻率和深度；從圖2所示理論數(shù)據(jù)和反演迭代響應(yīng)數(shù)據(jù)曲線對(duì)比，可以看出反演迭代響應(yīng)數(shù)據(jù)快速擬合理論數(shù)據(jù)，在最佳修正模型結(jié)果下的反演迭代數(shù)據(jù)和原始數(shù)據(jù)相對(duì)擬合差下降到1.5e-3,幾乎完全擬合了原始理論數(shù)據(jù)。反演過程中正則化因子和相對(duì)擬合差快速下降，可以看出自適應(yīng)正則化反演算法能夠快速穩(wěn)定的收斂。

圖2 H型模型測(cè)點(diǎn)a1反演結(jié)果以及擬合差、正則化因子變化曲線Fig.2 Inversion results of H model in survey station of a1, with the change curve of misfit and regularization factor

采用相同H型模型測(cè)點(diǎn)a2反演結(jié)果如圖3所示，通過6次反演迭代，獲得了較為理想的反演模型，模型響應(yīng)也擬合得較好。

H型最終反演結(jié)果非常接近真實(shí)模型，低阻層的深度基本相當(dāng)，中間層低電阻率接近真實(shí)值，測(cè)點(diǎn)a1反演結(jié)果第一層和第三層高阻也與真實(shí)值一致；測(cè)點(diǎn)a2接近大定源裝置觀測(cè)的極限位置，由于采用中心回線反演方式，對(duì)第一層高阻層的反應(yīng)不夠準(zhǔn)確。

為了驗(yàn)證該反演算法對(duì)高阻層的反應(yīng)能力，采用三層K型地電結(jié)構(gòu)(各層參數(shù)如下:ρ1=50 Ω·m、ρ2=300 Ω·m、ρ3=50 Ω·m，層厚分別為100 m、100 m)進(jìn)行驗(yàn)證。

圖4所示是K型地電模型測(cè)點(diǎn)a1的反演結(jié)果，由于第一層低阻層的屏蔽作用，反演需要11次迭代才得到了較好的反演結(jié)果，基本反映出了真實(shí)的電阻率和深度，高阻層電阻率基本接近真實(shí)模型。采用相同的地電模型，測(cè)點(diǎn)a2的反演結(jié)果如圖5所示。羅延鐘(羅延鐘等，2003)的正演研究結(jié)果表明，由于低阻層的屏蔽作用，瞬變電磁法對(duì)高阻層響應(yīng)非常微弱，所以也導(dǎo)致對(duì)高阻層的反演需要較多的迭代次數(shù)；強(qiáng)建科(強(qiáng)建科等，2013)采用OCCAM方法對(duì)航空瞬變電磁反演的K型模型結(jié)果也證明，瞬變電磁法對(duì)高阻層的反演不夠靈敏，基本無法反演出真實(shí)的高阻電阻率。

從圖4、圖5所示的響應(yīng)擬合曲線、相對(duì)擬合差和正則化因子變化曲線，可以看出隨著迭代次數(shù)增加，反演過程快速而穩(wěn)定。本文的自適應(yīng)正則化算法，由于采用了最光滑模型約束，對(duì)K型地電模型反演結(jié)果有很好的改善。

圖3 H型模型測(cè)點(diǎn)a2反演結(jié)果以及擬合差、正則化因子變化曲線Fig .3 Inversion results of H model in survey station of a2, with the change curve of misfit and regularization factor

3.2 四層模型反演試算

由于大地構(gòu)造的復(fù)雜性(毛立峰等，2011)，本文設(shè)計(jì)了四層HK型地電模型來驗(yàn)證本文算法的有效性。圖6所示HK型地電模型(其中ρ1=200 Ω·m、ρ2=50 Ω·m、ρ3=250 Ω·m、ρ4=100 Ω·m，層厚分別為100 m、100m、100m)測(cè)點(diǎn)a2處的反演結(jié)果，無論是深度位置還是電阻率值與理論模型都擬合得較好；還是由于低阻層的屏蔽作用，相近的高阻層無法反映出真實(shí)的電阻率值，但是與周圍地層電阻率的差異基本反映了出來。圖6所示響應(yīng)擬合曲線、數(shù)據(jù)擬合情況、正則化因子變化，說明了該反演算法的有效性和對(duì)復(fù)雜地電結(jié)構(gòu)的適應(yīng)性。

圖4 K型模型測(cè)點(diǎn)a1反演結(jié)果以及擬合差、正則化因子變化曲線Fig.4 Inversion results of K model in survey station of a1, with the change curve of misfit and regularization factor

圖5 K型模型測(cè)點(diǎn)a2反演結(jié)果以及擬合差、正則化因子變化曲線Fig .5 Inversion results of K model in survey station of a2, with the change curve of misfit and regularization factor

圖6 四型模型測(cè)點(diǎn)a2反演結(jié)果以及擬合差、正則化因子變化曲線Fig.6 Inversion results of four layers model in survey station of a2, with the change curve of misfit and regularization factor

4 結(jié)論

本文實(shí)現(xiàn)了大定源回線瞬變電磁資料的一維自適應(yīng)正則化反演，正則化因子的自適應(yīng)調(diào)整，是該方法快速而穩(wěn)定的收斂；由于采用最光滑模型約束模型目標(biāo)函數(shù)，對(duì)于瞬變電磁法對(duì)高阻層反應(yīng)不夠靈敏的情況有所改善；對(duì)于復(fù)雜地電結(jié)構(gòu)，該方法也有較好的有效性和適應(yīng)性；并且反演迭代10次左右，就能夠快速穩(wěn)定的收斂，可以得到較為滿意的結(jié)果，驗(yàn)證了該方法對(duì)于一維大定源回線瞬變電磁資料解釋是可行的。

盡管利用自適應(yīng)正則化反演算法進(jìn)行大定源瞬變電磁一維反演具有諸多優(yōu)點(diǎn),但其亦存在不足。就算法本身而言，對(duì)于靠近發(fā)射線框區(qū)域測(cè)點(diǎn)，也采用中心回線反演方法近似處理，對(duì)于淺部地層的反演結(jié)果不夠準(zhǔn)確，同時(shí)瞬變電磁方法本身對(duì)相對(duì)高阻層的反應(yīng)也不夠靈敏，但是本文采用的反演方法基本可以達(dá)到預(yù)期的效果。

Anders V C, Niels B C. 2003.A quantitative appraisal of airborne and ground-based transient electromagnetic(TEM) measure events in Denmark[J].Geophysics,68(2):523-534

Alan Y L, James M, Andrea V. 2010.Breaks in lithology :Interpretation problems when handing 2D structures with a 1D approximation[J].Geophysics, 70(4):179-188

Brodie R, Sambridge M. 2006.A holistic approach to inversion of frequency-domain airborne EM data[J]. Geophysics,71(6): 301-312

Christensen N B, Reid J E, Halkjmr M.2009.Fast laterally smooth inversion of airborne time-domain electromagnetic data[J] .Near Surface Geophysics, 7(5):599-612

Chen Xiao-bin, Zhao Guo-ze, Tang Ji, Zhan Yan, Wang Ji-jun. 2005. An adaptive regularized inversion algorithm for regularized data[J]. Chinese Journal of Geophysics, 48(4): 937-946(in Chinese with English abstract)

Constable S C, Parker R L, and Constable. 1987. Occam’s Inversion: A practical algorithm for generating smooth models from electro sounding data[J]. Geophysics, 5(2):289-300

Deszcz Pan M，F(xiàn)itterman D V, Labson V F. 1998.Reduction of inversion errors in helicopter EM data using auxiliary information[J].Exploration Geophysics, 29(2):142-146

Huang H, Palacky G J. 1991.Damped least-squares inversion of time-domain airborne EM data based on singular value decomposition[J]. Geophysical Prospecting, 39(6):827-844

Li Xiu. 2002. The theory and application of transient electromagnetic sounding. Xian: Shanxi Science and Technology Press:35-42(in Chinese)

Luo Yan-zhong, Zhang Sheng-ye, Wang Wei-ping. 2003. A research on one dimension forward for aerial electromagnetic method in time domain[J] .Chinese Journal of Geophysics,46(5): 719-724(in Chinese with English abstract)

Misac N N. 1992.Ectromagnetic Methods in Applied Geophysics [M].Zhao Jing-xiang tran.Beijing:Geological Publishing House: 195-207(in Chinese)

Misae N.Nabihgina. 1979. Quasi-static transient Response of a conducting half-space approximate presentation [J].Geophysics,Vol.44: 10-14

Mao Li-Feng, Wang Xu-ben, Chen Bin. 2011.Study on an adaptive regularized 1D inversion method of helicopter TEM data [J].Progress in Geophysics, 20 (1):300-305(in Chinese with English abstract)

Piao Hua-rong . 1990. The theory of electromagnetic sounding [M]. Beijing: Geological Publishing House :15-21(in Chinese)

Qiang Jian-ke, Li Yong-xing, Long Jian-bo. 2013.1-d Occam inversion method for airborne transient electromagnetic data[J]. Computing Techniques for Geophysical and Geochemical Exploration, 53(3):01-05(in Chinese with English abstract)

Sattel D. 2005.Inverting airborne (AEM) data with Zohdy ’s method[J].Geophysics, 70(4):77-85

Vallee M A, Smith R S. 2009. Application of Occam’s inversion to airborne time-domain electromagnetic [J].The leading Edge，28(3):284-287

Weng Ai-hua.2007. Occam’s inversion and its application to transient electromagnetic method [J]. Geology and Exploration,43(5):74-76(in Chinese with English abstract)

Yin C, Hodges G. 2007.Simulated annealing for airborne EM inversion [J].Geophysics, 72(4):189-195

[附中文參考文獻(xiàn)]

樸化榮.1990.電磁測(cè)深法原理[M].北京：地質(zhì)出版社:15-21

李 貅．2002.瞬變電磁測(cè)深的理論與應(yīng)用[M]．西安：陜西科學(xué)技術(shù)出版社:35-42

翁愛華. 2007.Occam反演及其在瞬變電磁測(cè)深中的應(yīng)用[J].地質(zhì)與勘探, 43(5)：74-76

米薩克N.納比吉安.1992.勘查地球物理-電磁法[M].趙經(jīng)祥譯.北京：地質(zhì)出版社:195-207

陳小斌，趙國(guó)澤,湯 吉,詹 艷，王繼軍．2005.大地電磁自適應(yīng)正則化反演算法[J]．地球物理學(xué)報(bào),48(4)：937-946

羅延鐘,張勝業(yè),王衛(wèi)平．2003.時(shí)間域航空電磁法一維正演研究[J]．地球物理學(xué)報(bào),46(5)：719-724

強(qiáng)建科，李永興，龍劍波.2013.航空瞬變電磁數(shù)據(jù)一維Occam反演[J]．物探化探計(jì)算技術(shù),53(3)：01-05

毛立峰，王緒本，陳 斌．2011.直升機(jī)航空瞬變自適應(yīng)正則化一維反演方法研究[J]．地球物理學(xué)進(jìn)展，20(1)：300-305

One-Dimensional Adaptive Regularization Inversion of Transient Electromagnetic Sounding with a Large Fixed Source

XU Yu-cong1, ZHAO Ning2, QIN Ce1, WANG Rui1, ZHANG Zhao-bing1

(1.CollegeofGeophysics,ChengduUniversityofTechnology,Chengdu,Sichuan610059； 2.CollegeofComputerScienceandTechnology,HenanPolytechnicUniversity,Jiaozuo,Henan454000)

Adaptive regularization inversion has less requirement of the initial model, and directly inverts transient electromagnetic response data. In the 1D inversion algorithm, regularization factors are adaptively obtained by calculating the data object function and model object function in each iteration, so that the geoelectric structure in the subsurface can be quickly determined. This paper adopts the adaptive regularization inversion to compute an one-dimensional layered model of transient electromagnetic sounding with a large fixed source loop, using a homogeneous half-space model as the initial model. The model object function is employed to the depth of the second derivative of tiny model constraints for the first time. By calculating the typical theory model, the result of inversion fits very well and has high precision. Furthermore, avoiding repeated search for regularization factors, this method has fast convergence and shows good stability and reliability.

TEM, large fixed source, adaptive regularized inversion, smoothest constraint model

2014-09-19；

2014-12-16；[責(zé)任編輯]郝情情。

國(guó)家自然科學(xué)基金重點(diǎn)基金(U1262206)資助。

徐玉聰(1988年-)，男，在讀碩士研究生，主要從事電磁法數(shù)值模擬。E-mail: xuyc1988@foxmail.com。

P618

A

0495-5331(2015)02-0360-06