# Improved Magnetotelluric Zohdy-Oldenburg Direct Inversion

Authors: Hui Cao; Yun Liu; Kun Peng Wang
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

Based on researches of Baiyao (1995) and Zhang and Xu (2001), this paper proposes an improved 2D MT Zohdy-Oldenburg direct inversion method, in the least-square sense, embodying the features of Zohdy’s ratio method and Oldenburg’s difference method, in the condition of rugged topography, with phase information. It bypasses large calculations of the Jacobian matrix and large sparse linear systems of equations and enables direct modifications and comparisons of the model parameters. According to the calculation and analysis of examples, it shows faster convergence and higher precision. In contrast with the conventional linear inversion, the calculation speed of this new method can be increased by more than 10 times.


---

## Body

## 1. Introduction

Calculating Jacobian matrix and solving large systems of linear equations are required in conventional magnetotelluric (referred to later as MT) linear or nonlinear optimized inversion method [1–4]. It is imaginable that the calculation is very complicated and the systems of equations are often ill conditioned. Therefore in this paper we think it is necessary to find a quick and effective inversion method to solve practical problems in engineering projects.Zohdy’s [5] method for the inversion calculation of direct-current resistivity sounding method (1989) converts polar distance (AB/2) to equivalent depth and modifies the model resistivity with the ratio between forward curve and measured curve to eliminate errors caused by zero-order asymptotic inversion (ratio method). This is called the asymptotic iterative inversion. Oldenburg and Ellis [6] in 1991 also proposed an iterative inversion method for MT studies, which modifies the model resistivity and model depth with the difference, calculated from Bostick, of forward result data and measured data of the prediction model inversion (difference method). Zohdy’s method is to modify the model resistivity with ratio method, so the rate of iterative convergence is apparently slower than that of Oldenburg’s method. In Zohdy’s iterative process, sometimes the model modification may be too much that the fitting errors will diverge. To solve this problem, two papers authored, respectively, by Shizhe and Bin (1995) [7, 8] suggest that one correction factor should be added to Zohdy’s method to control the step size which modifies the model resistivity. This inversion method is named as the curve comparison. Baiyao (1995) [9] repeated his study on the methods used by Zohdy and Oldenburg and proposed a least-square-sense-based direct inversion method, which combines Zohdy’s ratio method with Oldenburg’s difference method by filling longitudinal conductivity in the model grid and applies the combination in the inversion of direct-current resistivity sounding method [10] and MT sounding method [11]. On the basis of researches of Shizhe and Bin (1995) [7, 8] and Weidelt’s MT impedance phase and apparent resistivity theory [12], Zhang and Xu [13, 14] added impedance phase information to Zohdy’s ratio method and thus improved the rate of convergence in its iterative process.This shows that the MT direct inversion, bypassing calculating Jacobian matrix and solving large systems of linear equations, directly modifying and providing comparisons for the model parameters, can fully rival Bostick inversion. Therefore, to solve practical problems like slow speed of inversion and poor fitting precision found in engineering projects, this paper, based on previous studies, aims to improve the direct inversion method, namely, to propose an improved 2D MT Zohdy-Oldenburg direct inversion method with phase information.


## 2. Improved Zohdy-Oldenburg Direct Inversion in the Least-Square Sense

Suppose that there areN parameters of a geoelectric prediction model, and they are m=(m1,m2,…,mN)T. Their relationship with surface observed values, Ma, is as follows:
(1)Ma=Mc(m)+O(m).
In this equation, Mc(m) means theoretical observed values obtained from the forward computation of N model parameters, and O(m) means errors between observed values and theoretical values.To obtain the resultO(m)≈0, that is to say, the errors between observed values Ma and theoretical values Mc(m) can reach a certain degree of precision, multiple iterations and modifications of m, parameters of the geoelectric model, are often needed before a group of optimized m are produced. This is the general process of geophysical inversion.During thekth iterative calculation, for the prediction model m(k), we expand (1) to Taylor series and neglect terms which contain higher than second power in Δmj(k)=mj-mj(k), and then we get
(2)Ma=Mc(m(k))+∑j=1N∂Mc(m(k))∂mj(k)Δmj(k).
When there are M observed data, in matrix form, the above equation can be put into
(3)Ma=Mc(k)+A(k)Δm(k),
where A(k) is the M×N matrix, the Jacobian matrix (partial derivative matrix).Establishing the objective function is often needed in conventional linear inversions:(4)φ(m)=∥Ma-(Mc(k)+A(k)Δm(k))∥α,
where α is the matrix norm, and when α is 2, it becomes the least-square method. To make the objective function φ(m) tend to minimum and thus determine the amount of the model parameters’ modification Δm, calculating the Jacobian matrix and solving systems of linear equations are needed. This kind of inversion is burdened with large calculation and often contains ill-conditioned systems of equations.To avoid this problem, we introduce the Zohdy-Oldenburg direct inversion method, which can also be called the asymptotic iterative method. If we replace theA(k) in the Jacobian matrix with B(k), the objective function in (4) becomes
(5)φ(m)=∥Ma-(Mc(k)+B(k)Δm(k))∥α.
And when we use the second-order norm, this turns into a question about the least-square method, and its solution is
(6)Δm(k)=W(k)(Ma-Mc(k)),
where W(k)=(B(k))-1 and Δm(k)=m(k+1)-m(k) and therefore the following equation can be derived:
(7)m(k+1)=m(k)+W(k)(Ma-Mc(k)).To keep the algorithm from converging, a correction factorc(c≤1) will be added to (7) to limit the step size of modification. Thus we can obtain the iterative formula of Zohdy-Oldenburg direct inversion method in the least-square sense.In this case, when it comes to thejth model parameter, the inversion formula for (k+1)th iterative calculation is
(8)mj(k+1)=mj(k)+cWj(k)(Maj-Mcj(k)).
When c=Wj=1, (8) becomes Oldenburg’s asymptotic iterative method, which modifies the model with the difference method.WhenWj(k)=mj(k)/Mcj(k), (8) turns into
(9)mj(k+1)=mj(k)+c(MajMcj(k)-1)mj(k).
If c=1, (9) turns into Zohdy’s asymptotic iterative method, which modifies the model with the ratio method.Hence we can conclude that the improved Zohdy-Oldenburg direct inversion has features of asymptotic iterative method of the two scholars, ensuring that the algorithm in the iterative process can achieve faster convergence and higher stability of data.


## 3. Two-Dimensional Improved Zohdy-Oldenburg Direct Inversion

### 3.1. Two-Dimensional Direct Inversion

In the 2D inverse modeling of MT, we adopt an improved Zohdy-Oldenburg direct inversion method. The geoelectric prediction model’s parameter,m=(m1,m2,…,mN)T, is divided into two parts, the resistivity ρ=(ρ1,ρ2,…,ρM)T and the depth h=(h1,h2,…,hM)T. So (9) can be changed to the MT 2D improved Zohdy-Oldenburg direct inversion iterative formula:
(10)ρi,j(k+1)=ρi,j(k)+c1(ρai,jρci,j(k)-1)ρi,j(k),(11)hi,j(k+1)=hi,j(k)+c2(hai,jhci,j(k)-1)hi,j(k).
Equation (10) is the iterative formula for the improved Zohdy-Oldenburg direct inversion and gives resistivity inversion of the model. In this formula, ρi,j(k) stands for the true resistivity (or called intrinsic resistivity) of kth model at ith sounding point, of jth layer; ρai,j the measured apparent resistivity at ith sounding point, jth frequency point; ρci,j(k) the theoretical apparent resistivity of kth model at ith sounding point, jth frequency point; in this paper the correction factor c1=0.5.Equation (11) is the iterative formula for the improved Oldenburg’s iterative inversion and gives model depth inversion. In this formula, hi,j(k) represents the depth of kth model at ith sounding point, of jth layer; hai,j the apparent depth at ith sounding point, of jth layer (calculated from Bostick’s inversion formula); hci,j(k) the theoretical apparent depth of kth model at ith sounding point, of jth layer (calculated from Bostick’s inversion formula); in this paper the correction factor c2=0.5c1.


### 3.2. Two-Dimensional Direct Inversion with Phase Information

Zhang and Xu [14] and Weidelt [12], 1972, derived the approximate relationship between MT apparent resistivity and impedance phase:
(12)φ(f)≈π4(1+dlgρa(f)dlgf).
In this relationship, f represents the frequency, φ(f) the impedance phase, and ρa(f) the apparent resistivity.IfD(f)=dlgρa(f)/dlgf, then
(13)D(f)≈4πφ(f)-1.
Bostick inversion formula is already known, so we get
(14)ρ(f)=ρa(f)1-D(f)1+D(f),
where ρ(f) stands for the true resistivity in Bostick inversion formula.When we plug (13) into (14), we obtain the inversion formula of impedance phase:
(15)ρ(f)=ρa(f)(π2φ(f)-1).
When small changes of impedance phase bring changes to resistivity, namely, both sides of (15) take a derivative with respect to φ(f), the result is
(16)dρ(f)dφ(f)=ρa(f)π2-1φ(f)2.
Ignoring the local minimum of differential, we find that (16) turns into
(17)Δρ(f)=ρa(f)π2-Δφ(f)φ(f)2.
If we plug (17) into (10), the iterative formula of Zohdy’s direct inversion, we get
(18)ρi,j(k+1)=ρi,j(k)+c1(ρai,jρci,j(k)-1)ρi,j(k)+c1ρai,jπ2φci,j(k)-φai,jφai,j2.
Equation (18) is the improved iterative formula with phase information for Zohdy’s direct inversion. In this equation, φai,j stands for the measured impedance phase at ith sounding point, jth frequency point; φci,j(k) the theoretical impedance phase of kth model at ith sounding point, jth frequency point.The model depth inversion still uses (11), the improved iterative formula for Oldenburg direct inversion. If we put (18) and (11) in the simultaneous form, then we will obtain the 2D MT improved Zohdy-Oldenburg direct inversion with phase information.


### 3.3. Calculation Method of Direct Inversion

The improved Zohdy-Oldenburg direct inversion is based on the direct inversion of Zohdy and Oldenburg, and it constructs the initial model of inversion with Bostick inversion formula and combines the formulas of Zohdy and Oldenburg simultaneously while it inverses the model resistivity and model depth. This method, featuring faster speed, higher precision, and faster convergence, can tackle a series of problems caused by solving large sparse Jacobian matrix.Process of the improved Zohdy-Oldenburg direct inversion is as follows.Step 1.
Use Bostick inversion formula to calculate Bostick true resistivity and layer depth of each electric horizon, and fill 2D initial prediction model in the mesh model and perform a forward computation is as follows.Step 2.
Compare the measured apparent resistivity with the theoretical apparent resistivity obtained from the forward computation to check whether the precision can meet requirements. If the answer is yes, then it is certain that the prediction model is the result of the inversion.Step 3.
If not, modify the model parameters with Zohdy-Oldenburg iterative formula and obtain a new 2D prediction model. Conduct Bostick inversion, and get true resistivity and layer depth of each electric horizon. Refill the new 2D prediction model in the mesh model and again perform a 2D forward computation.Step 4.
Return to Step2.The fitting error between the measured apparent resistivity with the theoretical apparent resistivity will adopt mean square error; that is to say,(19)Error(k)=1NP×N∑i=1NP∑j=1N(ρai,j-ρci,j(k)ρai,j)2,
where Error(k) stands for the kth prediction model’s mean square error between the measured apparent resistivity and the theoretical apparent resistivity, NP the number of measuring points, and N the number of frequencies.


## 3.1. Two-Dimensional Direct Inversion

In the 2D inverse modeling of MT, we adopt an improved Zohdy-Oldenburg direct inversion method. The geoelectric prediction model’s parameter,m=(m1,m2,…,mN)T, is divided into two parts, the resistivity ρ=(ρ1,ρ2,…,ρM)T and the depth h=(h1,h2,…,hM)T. So (9) can be changed to the MT 2D improved Zohdy-Oldenburg direct inversion iterative formula:
(10)ρi,j(k+1)=ρi,j(k)+c1(ρai,jρci,j(k)-1)ρi,j(k),(11)hi,j(k+1)=hi,j(k)+c2(hai,jhci,j(k)-1)hi,j(k).
Equation (10) is the iterative formula for the improved Zohdy-Oldenburg direct inversion and gives resistivity inversion of the model. In this formula, ρi,j(k) stands for the true resistivity (or called intrinsic resistivity) of kth model at ith sounding point, of jth layer; ρai,j the measured apparent resistivity at ith sounding point, jth frequency point; ρci,j(k) the theoretical apparent resistivity of kth model at ith sounding point, jth frequency point; in this paper the correction factor c1=0.5.Equation (11) is the iterative formula for the improved Oldenburg’s iterative inversion and gives model depth inversion. In this formula, hi,j(k) represents the depth of kth model at ith sounding point, of jth layer; hai,j the apparent depth at ith sounding point, of jth layer (calculated from Bostick’s inversion formula); hci,j(k) the theoretical apparent depth of kth model at ith sounding point, of jth layer (calculated from Bostick’s inversion formula); in this paper the correction factor c2=0.5c1.


## 3.2. Two-Dimensional Direct Inversion with Phase Information

Zhang and Xu [14] and Weidelt [12], 1972, derived the approximate relationship between MT apparent resistivity and impedance phase:
(12)φ(f)≈π4(1+dlgρa(f)dlgf).
In this relationship, f represents the frequency, φ(f) the impedance phase, and ρa(f) the apparent resistivity.IfD(f)=dlgρa(f)/dlgf, then
(13)D(f)≈4πφ(f)-1.
Bostick inversion formula is already known, so we get
(14)ρ(f)=ρa(f)1-D(f)1+D(f),
where ρ(f) stands for the true resistivity in Bostick inversion formula.When we plug (13) into (14), we obtain the inversion formula of impedance phase:
(15)ρ(f)=ρa(f)(π2φ(f)-1).
When small changes of impedance phase bring changes to resistivity, namely, both sides of (15) take a derivative with respect to φ(f), the result is
(16)dρ(f)dφ(f)=ρa(f)π2-1φ(f)2.
Ignoring the local minimum of differential, we find that (16) turns into
(17)Δρ(f)=ρa(f)π2-Δφ(f)φ(f)2.
If we plug (17) into (10), the iterative formula of Zohdy’s direct inversion, we get
(18)ρi,j(k+1)=ρi,j(k)+c1(ρai,jρci,j(k)-1)ρi,j(k)+c1ρai,jπ2φci,j(k)-φai,jφai,j2.
Equation (18) is the improved iterative formula with phase information for Zohdy’s direct inversion. In this equation, φai,j stands for the measured impedance phase at ith sounding point, jth frequency point; φci,j(k) the theoretical impedance phase of kth model at ith sounding point, jth frequency point.The model depth inversion still uses (11), the improved iterative formula for Oldenburg direct inversion. If we put (18) and (11) in the simultaneous form, then we will obtain the 2D MT improved Zohdy-Oldenburg direct inversion with phase information.


## 3.3. Calculation Method of Direct Inversion

The improved Zohdy-Oldenburg direct inversion is based on the direct inversion of Zohdy and Oldenburg, and it constructs the initial model of inversion with Bostick inversion formula and combines the formulas of Zohdy and Oldenburg simultaneously while it inverses the model resistivity and model depth. This method, featuring faster speed, higher precision, and faster convergence, can tackle a series of problems caused by solving large sparse Jacobian matrix.Process of the improved Zohdy-Oldenburg direct inversion is as follows.Step 1.
Use Bostick inversion formula to calculate Bostick true resistivity and layer depth of each electric horizon, and fill 2D initial prediction model in the mesh model and perform a forward computation is as follows.Step 2.
Compare the measured apparent resistivity with the theoretical apparent resistivity obtained from the forward computation to check whether the precision can meet requirements. If the answer is yes, then it is certain that the prediction model is the result of the inversion.Step 3.
If not, modify the model parameters with Zohdy-Oldenburg iterative formula and obtain a new 2D prediction model. Conduct Bostick inversion, and get true resistivity and layer depth of each electric horizon. Refill the new 2D prediction model in the mesh model and again perform a 2D forward computation.Step 4.
Return to Step2.The fitting error between the measured apparent resistivity with the theoretical apparent resistivity will adopt mean square error; that is to say,(19)Error(k)=1NP×N∑i=1NP∑j=1N(ρai,j-ρci,j(k)ρai,j)2,
where Error(k) stands for the kth prediction model’s mean square error between the measured apparent resistivity and the theoretical apparent resistivity, NP the number of measuring points, and N the number of frequencies.


## 4. Analysis of Examples

### 4.1. One-Dimensional Layered MT Model

Figure1 is the H-type MT cross-sectional model of horizontally layered earth (3 layers), which is nonproportional and for illustration only. The resistivity and the layer thickness of the first layer are, respectively, 100 Ω·m and 500 m; for the second layer, they are 10 Ω·m and 500 m; the third layer’s resistivity is 1000 Ω·m. The computation domain is 0~4 km, the dot pitch is 100 m, and it totaled 41 measuring points. The range of frequency is 103~10−3 Hz, it totaled 61 frequency points, and the sampling interval is base-10 logarithms.Figure 1
H-type MT cross-sectional model.Figure2(a) is the inversion result of apparent resistivity, in which the horizontal axis represents depth (depth, unit: m) while the vertical axis represents resistivity (ρ, unit: Ω·m); Figure 2(b) is the curve graph of inversion fitting error, in which the horizontal axis represents iteration times and the vertical axis represents fitting error.Inversion result.
(a)
Inversion result of apparent resistivity
(b)
Inversion fitting errorAs can be seen from the inversion result of TE mode and TM mode, the anomalies of H-type stratum are revealed clearly. In the TE mode inversion result, the iteration times are 9, fitting errors percentage is 3.9%, and it takes 623 s; as for the TM mode inversion result, the iteration times are 9, fitting errors percentage is 6.2%, and it takes 410 s. The maximum relative error between these 2 modes is 5.3% and the mean square error is less than 2%.


### 4.2. Rugged Topography 2D Structural Model

Figure3 is the rugged topography 2D structural model. The largest drop in the rugged topography is 500 m, which contains faults and high-resistivity and low-resistivity layers; the measuring points cover 0~20 km, the dot pitch is 200 m, and it totaled 101 measuring points; the frequency is within the common frequency spectrum of a V8 electrical prospecting apparatus, 102.5~10−2.96 Hz, the sampling interval is base-10 logarithms, and it totaled 38 frequency points.Figure 3
Rugged topography 2D structural model.Figure4 shows the cross-sectional view of inversion results, in which the horizontal axis represents measuring point location (X, unit: km) and the vertical axis depth (Depth, unit: km); Figures 4(a) and 4(b) are apparent resistivity inversion results of TE mode and TM mode, respectively, while Figures 4(c) and 4(d) are apparent resistivity inversion results of TE mode and TM mode with phase. Figure 5 is the fitting error curve of inversion iteration, in which the horizontal axis represents number of iterations or iteration times, while the vertical axis represents fitting errors. Figure 5(a) is the comparison curve of fitting error of TE mode inversion iteration with and without phase while Figure 5(b) is the comparison curve of fitting error of TM mode inversion iteration with and without phase.Cross-sectional view of inversion results.
(a)
Cross-sectional view of TE mode inversion result
(b)
Cross-sectional view of TM mode inversion result
(c)
Cross-sectional view of TE mode inversion result (with phase)
(d)
Cross-sectional view of TM mode inversion result (with phase)Fitting error curve of inversion iteration.
(a)
Comparison curves of fitting error of TE mode inversion iteration with and without phase
(b)
Comparison curves of fitting error of TM mode inversion iteration with and without phaseFrom the above cross-sectional view of inversion results, it can be seen that the impact of rugged topography is eliminated, and inversions with topography of TE mode and TM mode are able to reflect locations and anomalies of faults and high-resistivity and low-resistivity layers. But it fails to reflect the second low-resistivity layer in the middle area of the figure, due to the equivalence phenomenon in geophysics.In the TE mode inversion without phase, the iteration times are 15 and the fitting error is 6.2%. In the TE mode inversion with phase, the iteration times are 17 and the fitting error is 5.7%. In the TM mode inversion without phase, the iteration times are 20 and the fitting error is 5.3%. In the TM mode inversion with phase, the iteration times are 25 and the fitting error is 2.8%.


### 4.3. Model with Abnormal Body under Rugged Topography

Figure6 shows the model with abnormal body under rugged topography, which is nonproportional and just for the convenience of easier understanding. In this figure, the drop is 60 m, the background resistivity is 100 Ω·m, and it contains a high-resistivity abnormal body of 500 Ω·m, which is 1600 m long and 150 m high, and it is 500 m away from its upper top surface to the ground surface. The measuring points cover 0~4 km, the dot pitch is 100 m, and it totaled 41 measuring points; the frequency is 104~10−2 Hz, the sampling interval is base-10 logarithms, and it totaled 61 frequency points.Figure 6
Geoelectric cross section model with abnormal body under rugged topography.Figure7 is the cross section of the inversion result. It is a resistivity (or intrinsic resistivity) contour map (ρ, unit: Ω·m), in which the horizontal axis stands for measuring points location (X, unit: km) and the vertical axis depth (Depth, unit: km). Figures 7(a) and 7(b) are apparent resistivity inversion results of TE mode and TM mode, respectively, while Figures 7(c) and 7(d) are apparent resistivity inversion results of TE mode and TM mode with phase. Figure 8 is the fitting error curve of inversion iteration, in which the horizontal axis represents number of iteration times, while the vertical axis represents fitting errors. Figure 8(a) is the comparison curve of fitting error of TE mode inversion iteration with and without phase while Figure 8(b) is the comparison curve of fitting error of TM mode inversion iteration with and without phase.Cross-sectional view of inversion results.
(a)
Cross-sectional view of TE mode inversion result
(b)
Cross-sectional view of TM mode inversion result
(c)
Cross-sectional view of TE mode inversion result (with phase)
(d)
Cross-sectional view of TM mode inversion result (with phase)Fitting error curve of inversion iteration.
(a)
Comparison curves of fitting error of TE mode inversion iteration with and without phase
(b)
Comparison curves of fitting error of TM mode inversion iteration with and without phaseThe key of inversions with topography lies in the reliability of numerical modeling, which is influenced by topography, in forward calculation. As can be observed from the cross-sectional view of the inversion results, the impact of rugged topography is significantly reduced in the TE mode and the TM mode. And since the impact of topography on TE mode is relatively small, the inversion result of abnormal body of TE mode is better than that of TM mode. In the TE mode inversion without phase, the iteration times are 2 and the fitting error is 6.6%; in the TE mode inversion with phase, the iteration times are 22 and the fitting error is 0.25%. In the TM mode inversion without phase, the iteration times are 20 and the fitting error is 0.58%; in the TM mode inversion with phase, the iteration times are 25 and the fitting error is 4.3%.


## 4.1. One-Dimensional Layered MT Model

Figure1 is the H-type MT cross-sectional model of horizontally layered earth (3 layers), which is nonproportional and for illustration only. The resistivity and the layer thickness of the first layer are, respectively, 100 Ω·m and 500 m; for the second layer, they are 10 Ω·m and 500 m; the third layer’s resistivity is 1000 Ω·m. The computation domain is 0~4 km, the dot pitch is 100 m, and it totaled 41 measuring points. The range of frequency is 103~10−3 Hz, it totaled 61 frequency points, and the sampling interval is base-10 logarithms.Figure 1
H-type MT cross-sectional model.Figure2(a) is the inversion result of apparent resistivity, in which the horizontal axis represents depth (depth, unit: m) while the vertical axis represents resistivity (ρ, unit: Ω·m); Figure 2(b) is the curve graph of inversion fitting error, in which the horizontal axis represents iteration times and the vertical axis represents fitting error.Inversion result.
(a)
Inversion result of apparent resistivity
(b)
Inversion fitting errorAs can be seen from the inversion result of TE mode and TM mode, the anomalies of H-type stratum are revealed clearly. In the TE mode inversion result, the iteration times are 9, fitting errors percentage is 3.9%, and it takes 623 s; as for the TM mode inversion result, the iteration times are 9, fitting errors percentage is 6.2%, and it takes 410 s. The maximum relative error between these 2 modes is 5.3% and the mean square error is less than 2%.


## 4.2. Rugged Topography 2D Structural Model

Figure3 is the rugged topography 2D structural model. The largest drop in the rugged topography is 500 m, which contains faults and high-resistivity and low-resistivity layers; the measuring points cover 0~20 km, the dot pitch is 200 m, and it totaled 101 measuring points; the frequency is within the common frequency spectrum of a V8 electrical prospecting apparatus, 102.5~10−2.96 Hz, the sampling interval is base-10 logarithms, and it totaled 38 frequency points.Figure 3
Rugged topography 2D structural model.Figure4 shows the cross-sectional view of inversion results, in which the horizontal axis represents measuring point location (X, unit: km) and the vertical axis depth (Depth, unit: km); Figures 4(a) and 4(b) are apparent resistivity inversion results of TE mode and TM mode, respectively, while Figures 4(c) and 4(d) are apparent resistivity inversion results of TE mode and TM mode with phase. Figure 5 is the fitting error curve of inversion iteration, in which the horizontal axis represents number of iterations or iteration times, while the vertical axis represents fitting errors. Figure 5(a) is the comparison curve of fitting error of TE mode inversion iteration with and without phase while Figure 5(b) is the comparison curve of fitting error of TM mode inversion iteration with and without phase.Cross-sectional view of inversion results.
(a)
Cross-sectional view of TE mode inversion result
(b)
Cross-sectional view of TM mode inversion result
(c)
Cross-sectional view of TE mode inversion result (with phase)
(d)
Cross-sectional view of TM mode inversion result (with phase)Fitting error curve of inversion iteration.
(a)
Comparison curves of fitting error of TE mode inversion iteration with and without phase
(b)
Comparison curves of fitting error of TM mode inversion iteration with and without phaseFrom the above cross-sectional view of inversion results, it can be seen that the impact of rugged topography is eliminated, and inversions with topography of TE mode and TM mode are able to reflect locations and anomalies of faults and high-resistivity and low-resistivity layers. But it fails to reflect the second low-resistivity layer in the middle area of the figure, due to the equivalence phenomenon in geophysics.In the TE mode inversion without phase, the iteration times are 15 and the fitting error is 6.2%. In the TE mode inversion with phase, the iteration times are 17 and the fitting error is 5.7%. In the TM mode inversion without phase, the iteration times are 20 and the fitting error is 5.3%. In the TM mode inversion with phase, the iteration times are 25 and the fitting error is 2.8%.


## 4.3. Model with Abnormal Body under Rugged Topography

Figure6 shows the model with abnormal body under rugged topography, which is nonproportional and just for the convenience of easier understanding. In this figure, the drop is 60 m, the background resistivity is 100 Ω·m, and it contains a high-resistivity abnormal body of 500 Ω·m, which is 1600 m long and 150 m high, and it is 500 m away from its upper top surface to the ground surface. The measuring points cover 0~4 km, the dot pitch is 100 m, and it totaled 41 measuring points; the frequency is 104~10−2 Hz, the sampling interval is base-10 logarithms, and it totaled 61 frequency points.Figure 6
Geoelectric cross section model with abnormal body under rugged topography.Figure7 is the cross section of the inversion result. It is a resistivity (or intrinsic resistivity) contour map (ρ, unit: Ω·m), in which the horizontal axis stands for measuring points location (X, unit: km) and the vertical axis depth (Depth, unit: km). Figures 7(a) and 7(b) are apparent resistivity inversion results of TE mode and TM mode, respectively, while Figures 7(c) and 7(d) are apparent resistivity inversion results of TE mode and TM mode with phase. Figure 8 is the fitting error curve of inversion iteration, in which the horizontal axis represents number of iteration times, while the vertical axis represents fitting errors. Figure 8(a) is the comparison curve of fitting error of TE mode inversion iteration with and without phase while Figure 8(b) is the comparison curve of fitting error of TM mode inversion iteration with and without phase.Cross-sectional view of inversion results.
(a)
Cross-sectional view of TE mode inversion result
(b)
Cross-sectional view of TM mode inversion result
(c)
Cross-sectional view of TE mode inversion result (with phase)
(d)
Cross-sectional view of TM mode inversion result (with phase)Fitting error curve of inversion iteration.
(a)
Comparison curves of fitting error of TE mode inversion iteration with and without phase
(b)
Comparison curves of fitting error of TM mode inversion iteration with and without phaseThe key of inversions with topography lies in the reliability of numerical modeling, which is influenced by topography, in forward calculation. As can be observed from the cross-sectional view of the inversion results, the impact of rugged topography is significantly reduced in the TE mode and the TM mode. And since the impact of topography on TE mode is relatively small, the inversion result of abnormal body of TE mode is better than that of TM mode. In the TE mode inversion without phase, the iteration times are 2 and the fitting error is 6.6%; in the TE mode inversion with phase, the iteration times are 22 and the fitting error is 0.25%. In the TM mode inversion without phase, the iteration times are 20 and the fitting error is 0.58%; in the TM mode inversion with phase, the iteration times are 25 and the fitting error is 4.3%.


## 5. Conclusions

In this paper, the proposed rugged topography 2D MT improved Zohdy-Oldenburg direct inversion is a asymptotic iteration method based on the least square, so the modification direction of this model always points to the descent direction of the objective function. This method embodies the features of both Zohdy’s asymptotic iteration method (ratio method) and Oldenburg’s asymptotic iteration method (difference method). The advantages of faster convergence and higher precision bypass large calculations of the Jacobian matrix and large sparse linear systems of equations and enable direct modifications and comparisons of the model parameters. In contrast with the conventional linear inversion, the calculation speed of the new method in this paper can be increased by more than 10 times.



---

*Source: 174586-2013-12-30.xml*