3D Marine Controlled-Source Electromagnetic Numerical Simulation Considering Terrain and Static Effect

Chunying Gu${}^{\mathbf{\text{1}}}$, Suyi Li${}^{\mathbf{\text{1}}}$, Wanyue Zhang${}^{\mathbf{\text{1}}}$, and Silun Peng${}^{\mathbf{\text{2,*}}}$

¹College of Instrumentation and Electrical Engineering
Jilin University, Changchun 130061, China
cygu20@mails.jlu.edu.cn, lsy@jlu.edu.cn, wanyue22@mails.jlu.edu.cn

²College of Automotive Engineering
Jilin University, Changchun 130022, China
${}^{*}$pengsilun@jlu.edu.cn

Submitted On: April 2, 2024; Accepted On: July 13, 2024

ABSTRACT

The marine controlled-source electromagnetic (CSEM) method is an important geophysical method for the exploration of marine hydrocarbon resources. In marine CSEM forward modeling the uplift terrain, such as submarine hills and seamounts, the static effects caused by polymetallic nodules and hydrothermal sulfide are ignored which can lead to the deviation of marine CSEM data. To improve the accuracy of data processing and interpretation, this study realizes efficient 3D numerical simulation considering submarine uplift terrain and static effect based on the finite difference (FD) algorithm in a fictitious wave domain. First, based on the correspondence principle between the fictitious wave domain and a real diffusive domain, we derive the FD electromagnetic field iteration equations of the fictitious wave domain, and apply the fictitious emission source and the boundary condition of complex frequency shifted-perfectly matched layer (CFS-PML). We use the inverse transformation method to convert the electromagnetic response to the diffusive domain. Then, we carry out simulations on typical 1D and 3D reservoir models to verify the correctness and effectiveness of the algorithm. Furthermore, we design an uplift terrain model and a static effect model and study the influence of parameters such as top width, bottom width, height and volume of uplift terrain on the CSEM field response characteristics through the forward modeling of multiple models and discuss the influence of parameters such as electrical conductivity, length, width, thickness, depth and volume of a shallow anomaly on the marine CSEM response. Finally, we analyze the characteristics and rules of electromagnetic field propagation of uplift terrain and static effect, which provides theoretical guidance for the design of marine CSEM exploration systems.

Index Terms: 3D numerical simulation, fictitious wave domain, marine controlled-source electromagnetic (CSEM), static effect, uplift terrain.

I. INTRODUCTION

The marine controlled-source electromagnetic (CSEM) is an important technology for deep-sea oil exploration, marine geological investigation and marine environment survey [1]. The marine CSEM method usually uses an antenna located tens of meters above the seabed to emit low-frequency electromagnetic signals, and a receiver located on the seabed continuously recording electromagnetic signals. By analyzing the collected electromagnetic signals, the resistivity distribution characteristics and electrical structure of the submarine medium can be inferred [1]. Marine CSEM technology can effectively identify high-resistivity oil reservoirs, reduce exploration risk and improve drilling success rate, which has been widely used in offshore oil exploration. After decades of vigorous development, it has achieved considerable success and economic benefits in instrument and equipment research and development, data processing and interpretation, and commercial exploration applications [2, 3]. Nowadays, the inversion interpretation of marine CSEM data is stepping into the 3D practical stage [4]. However, due to its complexity, 3D inversion technology is not mature yet, and the calculation accuracy and efficiency of the inversion algorithm needs to be improved. Moreover, 3D inversion studies considering seabed static effects are lacking [5]. 3D numerical simulation is the premise and theoretical basis of measured data inversion, and the accuracy and efficiency of inversion interpretation mainly depend on a numerical simulation algorithm. Therefore, it is of great theoretical significance and practical value to study the high-efficiency and high-precision 3D marine CSEM numerical simulation method.

With the improvement of computer performance, 3D numerical simulation technology has developed and advanced rapidly [6, 7]. Marine CSEM numerical simulation algorithms mainly include the integral equation (IE) method [8, 9], FD method [10, 11], finite volume (FV) method [12, 13] and finite element (FE) method [14, 15]. The FD method is currently the more used algorithm in marine electromagnetic forward modeling and inversion because of the advantages of simple calculation, strong applicability, high computational efficiency and sufficient computational accuracy [16]. To improve the computational efficiency, in recent years scholars have proposed a numerical simulation method based on wave field transformation technology, which uses correspondence between the fictitious wave domain and the real diffusive domain to convert the governing equation of the electromagnetic field in the diffusive domain to the fictitious wave domain, and then calculate the electromagnetic response in the fictitious wave domain. Finally, the electromagnetic response in the diffusive domain is obtained through inverse transformation [17–19]. This transformation method makes the electromagnetic field fluctuate in the ”fictitious wave domain” and slows down the propagation speed of the electromagnetic wave. Thus, it is easy to meet the stability conditions of the FD time-domain (FDTD) calculation, obtain a larger fictitious time step, and ensure that limited computing resources can be used to achieve efficient numerical calculation. Therefore, this method is widely used in the field of geophysical numerical simulation at present [20–22].

In marine CSEM exploration, the submarine hills and seamounts shaped by plate tectonics and seafloor spreading are usually uplifted, and the small polymetallic nodules and hydrothermal sulfides that may exist in the shallow part of the seafloor will produce static effect. If 3D numerical simulation ignores the effects of uplift terrain and static effect, it will cause large errors in data interpretation, which may further lead to the waste of data collection and an increase in exploration costs. For the forward modeling of submarine terrain, scholars have made some scientific achievements. Yutaka [23] uses the FD method to investigate the effect of 2D seafloor topography on 3D reservoir electromagnetic response and correct the topography by using the comparative method. Yang et al. [24] analyzes the influence of 2D submarine slope topography on the CSEM response of 3D high-resistivity reservoir based on the FV method. Yan and Han [25] realizes frequency-domain marine CSEM 3D forward modeling based on the FD method and discusses the influence of 2D mountain topography on electromagnetic response and the effect of topographic correction method. Yang et al. [26] uses the goal-oriented adaptive FE method to realize 3D marine CSEM modeling with anisotropy and topography. Therefore, it is important to study the impact of the topographic effect of 3D marine CSEM. For the forward modeling of static effect caused by small rocks and minerals on shallow surfaces, researchers have studied much in the numerical simulation of magnetotellurics and controlled-source audio-frequency magnetotellurics [27–31], but the 3D numerical simulation of marine CSEM considering static effect is very rare. In general, the complexity of the geological structure simulated by the FE method will increase, and the volume of calculations will also increase. A single inversion requires a large number of 3D simulations, and the FE method requires a large amount of memory, so calculation efficiency is limited and the calculation cost is greatly increased. The FD method can take into account calculation accuracy and efficiency. It is a practical method to simulate the 3D submarine uplift terrain and static effect.

In this paper, we carry out a numerical simulation study of the uplift terrain and static effect of 3D marine CSEM using the FD method in fictitious wave domain. We first derive the transformation formula between the governing equations in diffusive domain and the governing equations in fictitious wave domain, construct the Maxwell’s equations in the fictitious wave domain, and introduce the fictitious emission source and the boundary condition of complex frequency shifted-perfectly matched layer (CFS-PML). We use the FD method to solve the CSEM response and the inverse transformation method to recover the electromagnetic response in the real diffusive domain, so as to realize the 3D marine CSEM numerical simulation based on the FD method in the fictitious wave domain. To verify the correctness and effectiveness of the algorithm, we then establish a three-dimensional model and a one-dimensional reservoir model, calculate the electromagnetic response of the above models, and compare with the analytic solutions in the frequency domain. Finally, through the electromagnetic response calculation of the 3D reservoir model of uplift terrain and static effect model, we analyze the influence characteristics of uplift terrain and static effect on electromagnetic response, which can provide reference for other scholars’ subsequent research.

II. FICTITIOUS WAVE DOMAIN FD METHOD

To analyze and understand the influence characteristics of submarine terrain and static effect on marine CSEM propagation, and to achieve efficient and precise calculation of electromagnetic response of a complex reservoir model, this study proposes a 3D forward modeling method of submarine uplift terrain and a static effect based on the FD method in the fictitious wave domain. This method is mainly divided into three parts (as shown in Fig. 1). The first part is the construction of the governing equations in the fictitious wave domain, the second part is the solution of the electromagnetic response in the fictitious wave domain, including the discretization of the FD in the fictitious wave domain, the application of the fictitious emission source and the CFS-PML boundary condition, and the third part is the recovery of the electromagnetic response in the diffusive domain. A frame diagram of the method is shown in Fig. 1.

Figure 1: Frame diagram of the FD method in the fictitious wave domain.

A. Construction of the governing equations in the fictitious wave domain

We first design a submarine reservoir model according to the distribution characteristics of submarine topography and mineral resources. Next, we set the exploration simulation parameters, such as emission frequency and offset, according to the actual exploration situation. Finally, we transform the electromagnetic diffusive equations in the diffusive domain into electromagnetic wave equations in the fictitious wave domain through mathematical transformation. When the displacement current is neglected, the frequency domain expressions of quasi-static diffusive Maxwell’s equations are:

	$-\nabla \times \mathbf{\text{H}}(\omega )+\sigma \mathbf{\text{E}}(\omega )=-\mathbf{\text{J}}(\omega ),$		(1)
	$\nabla \times \mathbf{\text{E}}(\omega )-i\omega \mu \mathbf{\text{H}}(\omega )=-\mathbf{\text{K}}(\omega ),$		(2)

where $\mathbf{E}$ and $\mathbf{H}$ are electric and magnetic vector fields, respectively. $\mathbf{J}$ is the electric current density. $\mathbf{K}$ is the magnetic current density. $\omega$ is the angular frequency in the diffusion domain. $\sigma$ is the conductivity. $\mu$ is the magnetic permeability.

To transform Maxwell’s equations in the diffusive domain into the fictitious wave domain, angular frequency ${\omega }^{\prime }$ and dielectric permittivity ${\epsilon }^{\prime }$ in the fictitious wave domain need to be defined as follows:

	${\omega }^{\prime }$	$=(i+1)\sqrt{\omega {\omega }_0},$		(3)
	${\epsilon }^{\prime }$	$={\displaystyle \frac{\sigma }{2{\omega }_0}},$		(4)

where ${\omega }_0=2\pi f_0$ is the scale parameter. The choice of ${\omega }_0$ is in principle arbitrary, we use f${}_{\mathrm{0}}\mathrm{=}$ 1 Hz.

By multiplying both sides of equation (1) by $\sqrt{-i\omega /2{\omega }_0}$ at the same time and substituting fictitious angular frequency and fictitious dielectric permittivity into equations (1) and (2), using the transformation relationship between fictitious wave domain and real diffusive domain, the frequency domain expressions of Maxwell’s equations in fictitious wave domain can be obtained:

	$-\nabla \times {\mathbf{H}}^{\prime }\left({\omega }^{\prime }\right)-i{\omega }^{\prime }{\epsilon }^{\prime }{\mathbf{E}}^{\prime }\left({\mathit{\omega }}^{\prime }\right)$	$=-{\mathbf{J}}^{\prime }\left({\mathit{\omega }}^{\prime }\right),$		(5)
	$\nabla \times {\mathbf{E}}^{\prime }\left({\mathit{\omega }}^{\prime }\right)-i{\omega }^{\prime }\mu {\mathbf{H}}^{\prime }\left({\mathit{\omega }}^{\prime }\right)$	$=-{\mathbf{K}}^{\prime }\left({\mathit{\omega }}^{\prime }\right),$		(6)

where ${\mathbf{E}}^{\prime }$ and ${\mathbf{H}}^{\prime }$ are electric and magnetic vector fields in the fictitious wave domain. ${\mathbf{J}}^{\prime }$ and ${\mathbf{K}}^{\prime }$ are electric current and magnetic current densities in the fictitious wave domain.

Through the inverse Fourier transform, the time-domain expression of Maxwell’s equations in the fictitious wave domain can be obtained:

	$-\nabla \times {\mathbf{H}}^{\prime }\left({\mathbf{t}}^{\prime }\right)+{\epsilon }^{\prime }{\partial }_{t^{\prime }}{\mathbf{E}}^{\prime }\left({\mathbf{t}}^{\prime }\right)$	$=-{\mathbf{J}}^{\prime }\left({\mathbf{t}}^{\prime }\right),$		(7)
	$\nabla \times {\mathbf{E}}^{\prime }\left({\mathbf{t}}^{\prime }\right)+\mu {\partial }_{t^{\prime }}{\mathbf{H}}^{\prime }\left({\mathbf{t}}^{\prime }\right)$	$=-{\mathbf{K}}^{\prime }\left({\mathbf{t}}^{\prime }\right).$		(8)

B. Solution of electromagnetic responses in the fictitious wave domain

The FDTD method is widely used in practical engineering. The form is simple and can simulate large scale models numerically. Compared with the traditional FDTD method in the diffusive domain, the FDTD method in the fictitious wave domain is easier to meet the stability conditions of the iterative calculation, and a larger fictitious time step can be adopted to reduce the number of iterative and improve the computational efficiency [11]. In the specific solution, we introduce the fictitious emission source to avoid calculation of the initial electromagnetic field and the assumption of flat terrain in the diffusive domain. We apply CFS-PML absorption boundary conditions, which can efficiently absorb low frequency and late electromagnetic waves, reduce memory consumption, and further save computing time.

(a) Discretization of FD in the fictitious wave domain

After constructing the governing equations in the fictitious wave domain, it is necessary to solve the governing equations. Considering the accuracy and efficiency of the calculation, this study uses the FDTD method to solve Maxwell’s equations. Yee grid [32] and the explicit Du Fort-Frankel method [33] are used to discretize the fictitious wave domain electric and magnetic fields in space and time, respectively, and the iterative formulas of 3D FD in the fictitious wave domain are obtained. Taking component $E_x^{\prime }$ as an example, the FD expression of the $n+1$ times iteration is given as:

	$E_{x_{i+1/2,j,k}}^{\prime n+1}=E_{x_{i+1/2,j,k}}^{\prime n}+\mathrm{\Delta }{t}^{\prime }{\displaystyle \frac{2{\omega }_0}{\sigma }}[{\partial }_y^{-}{H}_{z_{i+1/2,j+1/2,k}}^{\prime n+1/2}$
	$\mathrm{\hspace{1em}\hspace{1em}}\mathit{\hspace{1em}\hspace{1em}}-{\partial }_z^{-}{H}_{y_{i+1/2,j,k+1/2}}^{\prime n+1/2}]-\mathrm{\Delta }t{}^{\prime }{\displaystyle \frac{2{\omega }_0}{\sigma }}J{}_x{}^{\prime },$		(9)

where $E{}^{\prime }$ and $H{}^{\prime }$ are electric and magnetic fields responses in the fictitious wave domain. The electric current density in the fictitious wave domain is$J{}^{\prime }$. The sampling time in the fictitious wave domain is $\mathrm{\Delta }t{}^{\prime }$. To ensure the stability of the iterative calculation of time-domain FD method, the time step $\mathrm{\Delta }t{}^{\prime }$ is required to satisfy Courant-Friedrichs-Lewy (CFL) condition:

$$\mathrm{\Delta }{t}^{\prime }\le \frac{1}{c_{\max }\sqrt{\frac{1}{\mathrm{\Delta }{x}^2}+\frac{1}{\mathrm{\Delta }{y}^2}+\frac{1}{\mathrm{\Delta }{z}^2}}},$$

(10)

where $\mathrm{\Delta }x$, $\mathrm{\Delta }y$ and $\mathrm{\Delta }z$ are the cell grid dimensions. The maximum propagation speed of an electromagnetic wave in the fictitious wave domain is $c_{\max }=\sqrt{2{\omega }_0/\mu {\sigma }_{\min }}$.

(b) Application of the fictitious emission source

The selection of emission source signal in the virtual wave domain is flexible and not limited to a single function. However, in the selection process, it is necessary to ensure that the transmitted signal contains a large amount of spectrum, and the form should be easy to calculate. Therefore, in this study, the first derivative of Gaussian pulse is used as the fictitious emission source in the process of 3D numerical simulation. Adopting the above fictitious emission source can not only avoid the singularity of the source and the complex calculation of the initial electromagnetic field, but also break the assumption of flat terrain. The time-domain and frequency-domain expressions of the emission source are:

	$J(t)$	$=-2\beta (t-t_0)\sqrt{{\displaystyle \frac{\beta }{\pi }}}{e}^{-\beta {(t-t_0)}^2},$		(11)
	$J(\omega )$	$=i\omega e^{-\frac{{\omega }^2}{4\beta }}{e}^{-i\omega t_0},$		(12)

where $\beta =\pi f_{\max }^2$, $t_0=\pi \mathrm{/}{f}_{\max }$. The maximum frequency of the electromagnetic field transmission in the fictitious wave domain is$f_{\max }$. In the application of the emission source, the current direction in the fictitious wave domain is consistent with the inline direction in the real diffusive domain. The fictitious emission source is only related to the $x$ and $y$ components of the electric field, and the electric field component in the z-direction is zero.

(c) Introduction of the CFS-PML boundary condition

The CFS-PML absorbing boundary condition has a good absorption effect on electromagnetic wave propagation, which can greatly relieve the memory pressure on the computer and improve computational efficiency. Therefore, considering high efficiency computation, this paper introduces the CFS-PML boundary condition in the fictitious wave domain. Taking the electric field component $E_x^{\prime }$ as an example, its FD calculation expression in the fictitious wave domain is:

	$E_{x_{i+1/2,j,k}}^{\prime n+1}$	$=E_{x_{i+1/2,j,k}}^{\prime n}+$
		$\mathrm{\Delta }{t}^{\prime }{\displaystyle \frac{2{\omega }_0}{\sigma }}\left[{\partial }_y^{-1}{H}_{z_{i+1/2,j+1/2,k}}^{\prime n+1/2}-{\partial }_z^{-1}{H}_{y_{i+1/2,j,k+1/2}}^{\prime n+1/2}\right]$
		$+\mathrm{\Delta }{t}^{\prime }{\displaystyle \frac{2{\omega }_0}{\sigma }}\left({\mathrm{\Psi }}_{e_{x{y}_{i+1/2,j,k}}}^{n+1/2}-{\mathrm{\Psi }}_{e_{x{z}_{i+1/2,j,k}}}^{n+1/2}\right),$		(13)
	${\mathrm{\Psi }}_{e_{x{y}_{i+1/2,j,k}}}^{n+1/2}$	$=e^{-\left(\frac{{\sigma }_y}{k_y}+{\alpha }_y\right)\frac{\mathrm{\Delta }{t}^{\prime }}{{\epsilon }^{\prime }}}{\mathrm{\Psi }}_{e_{x{z}_{i+1/2,j,k}}}^{n-1/2}+$
		${\displaystyle \frac{{\sigma }_y}{\left({\sigma }_y{k}_y+k_y^2{\alpha }_y\right)}}\left(e^{-\left(\frac{{\sigma }_y}{k_y}+{\alpha }_y\right)\frac{\mathrm{\Delta }{t}^{\prime }}{{\epsilon }^{\prime }}}-1\right)$
		$\cdot \left(H_{z_{i+1/2,j+1/2,k}}^{\prime n+1/2}-H_{y_{i+1/2,j-1/2,k}}^{\prime n+1/2}\right)/{\mathrm{\Delta }}_y,$		(14)
	${\mathrm{\Psi }}_{e_{x{z}_{i+1/2,j,k}}}^{n+1/2}$	$=e^{-\left(\frac{{\sigma }_z}{k_z}+{\alpha }_z\right)\frac{\mathrm{\Delta }{t}^{\prime }}{{\epsilon }^{\prime }}}{\mathrm{\Psi }}_{e_{x{z}_{i+1/2,j,k}}}^{n-1/2}+$
		${\displaystyle \frac{{\sigma }_z}{\left({\sigma }_z{k}_z+k_z^2{\alpha }_z\right)}}\left(e^{-\left(\frac{{\sigma }_z}{k_z}+{\alpha }_z\right)\frac{\mathrm{\Delta }{t}^{\prime }}{{\epsilon }^{\prime }}}-1\right)$
		$\cdot \left(H_{y_{i+1/2,j+1/2,k}}^{\prime n+1/2}-H_{y_{i+1/2,j-1/2,k}}^{\prime n+1/2}\right)/{\mathrm{\Delta }}_z,$		(15)

where ${\mathrm{\Psi }}_{e_{x{y}_{i+1/2,j,k}}}^{n+1/2}$ and ${\mathrm{\Psi }}_{e_{x{z}_{i+1/2,j,k}}}^{n+1/2}$ are two discrete variables. ${\sigma }_i$, $k_i$ and ${\alpha }_i(i=x,y,z)$ are boundary parameters, which depend on the relative positions of the node and the target area, but are not fixed values:

	${\sigma }_i$	$={\sigma }_{\max }{\displaystyle \frac{{\left|k-k_0\right|}^m}{d^m}},i=x,y,z,$		(16)
	$k_i$	$=1+(k_{\max }-1){\displaystyle \frac{{\left|k-k_0\right|}^m}{d^m}},i=x,y,z,$		(17)
	${\alpha }_i$	$={\alpha }_{\max }{\displaystyle \frac{{\left|k-k_0\right|}^m}{d^m}},i=x,y,z,$		(18)

where $\left|k-k_0\right|$ is the distance from the calculation area to the absorption boundary. The polynomial parameter is $m$. The boundary layer thickness is $d$.

Figure 2: Comparison of the electromagnetic response from this article with semi-analytical solution (semi-AS) for the 1D oil reservoir model: (a) amplitude curve and (b) relative error curve.

C. Recovery of electromagnetic response in the diffusive domain

After solving the time-domain electromagnetic response in the fictitious wave domain, it is transformed into the frequency-domain response result in the diffusive domain. The conversion relationship is asfollows:

	$E_x(\omega )$	$=\sqrt{{\displaystyle \frac{-i\omega }{2{\omega }_0}}}{\displaystyle {\int }_0^T}{E}_x^{\prime }(t^{\prime })e^{-\sqrt{\omega {\omega }_0}{t}^{\prime }}{e}^{i\sqrt{\omega {\omega }_0}{t}^{\prime }}\mathit{d}{t}^{\prime },$		(19)
	$J_x(\omega )$	$={\displaystyle {\int }_0^T}{J}^{\prime }(t^{\prime })e^{-\sqrt{\omega {\omega }_0}{t}^{\prime }}{e}^{i\sqrt{\omega {\omega }_0}{t}^{\prime }}\mathit{d}{t}^{\prime },$		(20)

where $T=n\mathrm{\Delta }{t}^{\prime }$ is the FD computation time in the fictitious wave domain. The sampling time in the fictitious wave domain is$\mathrm{\Delta }{t}^{\prime }$. The above electromagnetic response is excited by the emission source in the fictitious wave domain, not by the real emission source. To obtain the real electromagnetic response in the diffusive domain, it is necessary to construct Green’s function or the impulse response function of the system according to different working devices. The marine CSEM emission source can be equivalent to a point source or a dipole source. The expression of its diffusive frequency domain Green’s function is:

$$G_x(\omega )=\frac{E_x(\omega )}{J_x(\omega )}.$$

(21)

Through the derivation of Green’s function formula in the full space diffusion frequency domain, we can see that Green’s function has the following properties. Green’s function in the diffusive domain becomes independent of the scale parameter${\omega }_0$and the spectral width or maximum frequency that can be recovered in the diffusive domain is independent of the spectral width or maximum frequency in the fictitious wavedomain [11].

We can use Green’s function to recover the electromagnetic response in the real diffusive domain, obtain the marine CSEM response of the submarine reservoir model, and then analyze and understand the electromagnetic field propagation characteristics and rules of different models, so as to provide a theoretical basis for subsequent marine CSEM exploration under complex geological conditions.

III. ACCURACY VERIFICATION

A. 1D oil reservoir model

We chose a 1D oil reservoir model to verify the accuracy of our algorithm against the semi-analytical solutions given by Kong [34]. The first layer of the 1D model is a seawater layer with a thickness of 1 km and a conductivity of 3.2 S/m. The second layer is a seabed sedimentary layer with a thickness of 1 km and a conductivity of 1 S/m. The third layer is an oil layer with a thickness of 0.2 km and a conductivity of 0.01 S/m. The fourth layer is a downward extending sedimentary layer with a conductivity of 1 S/m. The size of the target area is 20$\times$20$\times$11 km. The transmitter and receivers are laid on the surface of the seabed. The 1D oil reservoir model adopts 100$\times$100$\times$100 m discrete grid. The simulation results of the FWD-FD method are compared with the semi-analytical solution, and the calculation results when the transmitting frequency is 0.1 Hz are shown in Fig. 2. It can be seen from Fig. 2 that the numerical solution of the electromagnetic field obtained by the FWD-FD method is basically consistent with the semi-analytical solution. When the offset is between 2 and 9 km, the amplitude relative error is less than 2.8%, which indicates that the proposed algorithm has high accuracy.

Figure 3: (a) Comparison of the electromagnetic response of different algorithms for the 3D oil reservoir model and (b) snapshot of the induced current system in the fictitious wave field.

Figure 4: 3D oil reservoir model under uplift terrain and comparison model: (a) flat terrain without oil reservoir model (M00), (b) flat terrain 3D oil reservoir model (M01), (c) uplift terrain without oil reservoir model (M02) and (d) uplift terrain 3D oil reservoir model (M03).

B. 3D oil reservoir model

A 3D oil reservoir model as shown in the upper left corner of Fig. 1 is constructed. The thickness of seawater layer is 1 km and the conductivity is 3.2 S/m. The conductivity of seabed sedimentary layer is 1 S/m. The oil and gas reservoir is located 1 km below the seabed, with a size of 4$\times$4$\times$0.1 km and a conductivity of 0.01 S/m. The size of the target area is 20$\times$20$\times$11 km. The transmitter is located directly above the oil and gas reservoir. The receiver is laid on the surface of the seabed. The 3D reservoir model adopts 100$\times$100$\times$100 m discrete grid. The numerical simulation results of the fictitious wave domain FD (FWD-FD) method and frequency domain FV method are compared. The amplitude curve of the electromagnetic response when the emission frequency is 0.1 Hz is shown in Fig. 3, and a snapshot of the induced current system in the fictitious wave domain is shown in Fig. 4. It can be seen from Fig. 3 that the calculation results of the above two algorithms are relatively consistent, which further verifies the effectiveness of the calculation program of the FWD-FDmethod.

The number of grids for 3D simulation of the frequency domain FV method is 72$\times$134$\times$64, occupying 3.5 G memory. The maximum memory usage by the MUMPS direct solver for solving the coefficient matrix is 56.2 G, and the computation takes 48 min [35]. In this article, the number of grids for FD numerical simulation in the fictitious wave domain is 220$\times$220$\times$130, occupying 1.2 G memory. The number of iterations is 6000, the maximum memory usage during the computation is 2.1 G. The whole workflow of forward modelling mainly includes FWD-FD, time-frequency domain transformation of electromagnetic response and Green’s function transformation, whose running time is 22 min, 4 min 12 s, 31 s, respectively, and the total execution time of the whole program is 27 min. The machine is configured as follows: CPU is Intel Core i5-7500 and main frequency is 3.4 GHz. Memory is 16 GB. It can be seen that the time domain difference method does not need to solve large linear equations, so it computes a much larger number of grids than the frequency domain FV method on a computer with the same memory. After solving the time domain electromagnetic response in the fictitious wave domain, all frequency domain electromagnetic responses less than$f_{\max }$ can be calculated quickly, and there is no need to simulate again. When the electromagnetic responses of multiple frequencies need to be calculated, the time domain FD method in the fictitious wave domain can save a lot of computation time, which proves the high efficiency of the proposed method.

Figure 5: Electromagnetic response of the uplift terrain model at a frequency of 0.1 Hz, the relative effect of topography and the relative anomaly curves of oil reservoirs: (a) electromagnetic response curve, (b) topography relative influence curve and (c) reservoir relative anomaly curve.

IV. NUMERICAL EXPERIMENTS

In this section, uplift terrains with different parameters are studied first, followed by a static effect model under different conditions. Finally, a synthetic model including terrain and static effect is designed, and its responses are analyzed.

A. Uplift terrain model

To study the influence characteristics and laws of submarine uplift terrain on the marine CSEM response, we designed a 3D oil reservoir model containing uplift terrain and a comparison model (Fig. 4). In Fig. 4, M00 is a flat terrain without oil reservoir model, M01 is a flat terrain 3D oil reservoir model, M02 is an uplift terrain without oil reservoir model and M03 is an uplift terrain 3D oil reservoir model. We simulate the 3D marine CSEM response of the above model using the FD method in the fictitious wave domain, and calculate the relative influence of topography and the relative anomaly of oil reservoir. We first set the top width of the uplift terrain as 0.4 km, bottom width as 2 km, height as 0.2 km, center coordinate as (4000,0,1900), electrical conductivity as 1 S/m. There is a 3D oil reservoir directly below the uplift, the size of which is 4$\times$4$\times$0.2 km, the central coordinate is (4000,0,3100), and the conductivity is 0.01 S/m. The thickness of seawater layer is 2 km and the conductivity is 3.2 S/m. The conductivity of submarine sediment is 1 S/m. The size of the target area is 16$\times$16$\times$8 km, and the receiver is laid on the surface of the seabed along the uplift terrain.

Figure 6: Electromagnetic response curves of uplift terrain with different top widths, bottom widths and heights: (a) response curves with different top widths, (b) response curves with different bottom widths and (c) response curves with different heights.

Figure 7: Amplitude normalization curves of electromagnetic response for different top widths, bottom widths and heights of uplift terrain: (a) amplitude normalization curves for different top widths, (b) amplitude normalization curves for different bottom widths and (c) amplitude normalization curves for different heights.

The flat terrain without oil reservoir model, flat terrain 3D oil reservoir model, uplift terrain without oil reservoir model and uplift terrain 3D oil reservoir model adopt a 100$\times$100$\times$100 m discrete grid. The electromagnetic response calculation results of the four models at a frequency of 0.1 Hz, the relative effect of topography and the relative anomaly curves of oil reservoirs are shown in Fig. 5. It can be seen from the curves of M00 and M01 and their relative anomalies that the electromagnetic anomalies of oil reservoirs are obvious in the flat terrain 3D model. It can be observed from the curves of M00 and M02, M01 and M03 as well as the relative influence curves that the uplift terrain has a great influence on the electric field component. In addition, the M00 and M03 curves and their relative anomaly curves show that the uplift terrain will also weaken the electromagnetic anomalies generated by the reservoir. Therefore, when using the amplitude versus offset curve to qualitatively analyze potential oil reservoirs, it is necessary to consider the influence of uplift terrain on the CSEM response.

Let us now change the top width, bottom width and height of the uplift terrain, where the top width t is set to 600 m, 800 m and 1000 m, the bottom width b is set to 3000 m, 4000 m and 5000 m, and the height h is set to 300 m, 400 m and 500 m. The above models adopt 100$\times$100$\times$100 m discrete grid. The calculated electromagnetic responses at a frequency of 0.1 Hz are shown in Fig. 6, and the amplitude normalization values of the electromagnetic response curves and typical model electromagnetic response curves at different top widths, bottom widths and heights are shown in Fig. 7.

As can be seen from Fig. 6, the electromagnetic response curves of different top widths, bottom widths and heights are not distinct from those of typical terrain, indicating that the changes of top width, bottom width and height of uplift terrain have a weak influence on the electromagnetic response of the model. However, as can be seen from the amplitude normalization curve in Fig. 7, the uplift height has the greatest influence on the calculated results, followed by the bottom width and the top width. With the increase of height, bottom width and top width, the amplitude normalization values decrease gradually.

Figure 8: Electromagnetic response curves and amplitude normalization curves of different volumes in the uplift terrain model: (a) electromagnetic response curves for different volumes of uplift terrain and (b) amplitude normalization curves of electromagnetic response for different volumes of uplift terrain.

Finally, we change the volume of the uplift terrain, and set the top width, bottom width and height as (t, b, h) (400 m, 2000 m, 200 m), (600 m, 3000 m, 300 m), (800 m, 4000 m, 400 m) and (1000 m, 5000 m, 500 m). The grid size and scaling frequency adopted by the model are consistent with those in Fig. 5. The electromagnetic response calculation results when the frequency is 0.1 Hz are shown in Fig. 8 (a), and the amplitude normalization values of the electromagnetic response curves at different volumes of uplift terrain and typical terrain are shown in Fig. 8 (b).

As can be seen from Fig. 8 (a), there are differences between the electromagnetic response curves of uplift terrain with different volumes and those of typical terrain, indicating that the change of uplift terrain volume has a strong influence on the electromagnetic response of the model. However, it can be seen from the amplitude normalization curves in Fig. 8 (b) that with the an increase of uplift terrain volume, the amplitude normalization values gradually decrease and, compared with the amplitude normalization curves of the electromagnetic response with different top width, bottom width and height in Fig. 7, the decrease trend of the amplitude normalization values are more significant, indicating that the volume change is more obvious than the single size change in the topographic effect.

B. Static effect model

To study the influence characteristics of static effects on the marine CSEM response, we designed a static effect model (M04) containing shallow anomalies (such as polymetallic nodules, hydrothermal sulfides) with deep oil reservoir (Fig. 9), where there is a shallow abnormal body B directly above the oil reservoir A. Reservoir A is 4$\times$4$\times$0.2 km in size, central coordinate is (4000,0,3100), and conductivity is 0.01 S/m. Abnormal body B is 0.5$\times$0.5$\times$0.1 km, and central coordinate is (4000,0,2050). The thickness of seawater layer is 2 km and conductivity is 3.2 S/m. The conductivity of sediments is 0.5 S/m. The size of the target area is 16$\times$16$\times$8 km, and the transmitter and receivers are arranged on the seabed surface.

Figure 9: Static effect model (M04).

We first set the conductivity of abnormal body B as $2\mathrm{S}/\mathrm{m},0.2\mathrm{S}/\mathrm{m},0.1\mathrm{S}/\mathrm{m},0.02\mathrm{S}/\mathrm{m},0.01\mathrm{S}/\mathrm{m},0.008$ $\mathrm{S}/\mathrm{m}$, $0.00625\mathrm{S}/\mathrm{m},0.005\mathrm{S}/\mathrm{m}$ and $0.002\mathrm{S}/\mathrm{m}$, respectively. The corresponding resistivity is $0.5\mathrm{\Omega }\cdot \mathrm{m},5\mathrm{\Omega }\cdot \mathrm{m},10\mathrm{\Omega }\cdot \mathrm{m}$, $50\mathrm{\Omega }\cdot \mathrm{m},100\mathrm{\Omega }\cdot \mathrm{m},125\mathrm{\Omega }\cdot \mathrm{m},160\mathrm{\Omega }\cdot \mathrm{m},200\mathrm{\Omega }\cdot \mathrm{m}$ and $500\mathrm{\Omega }\cdot \mathrm{m}$. The no abnormal body and oil reservoir model, reservoir model and abnormal body and reservoir model adopt $100\times 100\times 100\mathrm{m}$ discrete grid. The electromagnetic response calculated by the above model at a frequency of 0.1 Hz is shown in Fig. 10 .

Figure 10: Electromagnetic response curves of static effect model and comparison model at 0.1 Hz.

Figure 11: Electromagnetic response curves and amplitude normalization curves of abnormal body B at different lengths, widths, thicknesses and depths in the static effect model: (a) electromagnetic response curves of abnormal body B at different lengths, widths, thicknesses and depths and (b) amplitude normalization curves of electromagnetic responses.

Figure 12: Electromagnetic response curves and amplitude normalization curves of abnormal body B at different lengths, widths, thicknesses and depths in the static effect model: (a) electromagnetic response curves of abnormal body B at different volumes and (b) amplitude normalization curves of electromagnetic responses of abnormal body B with different volumes.

As can be seen from Fig. 10, the curves of the $0.5\mathrm{\Omega }\cdot \mathrm{m},5\mathrm{\Omega }\cdot \mathrm{m},10\mathrm{\Omega }\cdot \mathrm{m},50\mathrm{\Omega }\cdot \mathrm{m}$ and $100\mathrm{\Omega }\cdot \mathrm{m}$ abnormal body and reservoir models nearly overlap with those of the reservoir models, indicating that the static effect of the above models is not obvious. However, the electromagnetic response calculation results of $125\mathrm{\Omega }\cdot \mathrm{m}$, $160\mathrm{\Omega }\cdot \mathrm{m},200\mathrm{\Omega }\cdot \mathrm{m}$ and $500\mathrm{\Omega }\cdot \mathrm{m}$ models of abnormal body and reservoir show that when the resistivity of the abnormal body is larger than that of the oil reservoir, the amplitude of electromagnetic response curve decreases gradually with an increase of resistivity, and the deviation from the electromagnetic response curve of the reservoir model is larger. The static effect is more obvious.

Next, we set the conductivity of abnormal body B as $0.005\mathrm{S}/\mathrm{m}$, that is, the resistivity is $200\mathrm{\Omega }\cdot \mathrm{m}$. The center coordinates remain unchanged, and the length, width, thickness and depth of the abnormal body are changed, respectively. The length 1 of the abnormal body in the $x$ direction is set as $1000\mathrm{m},1500\mathrm{m},2000\mathrm{m},2500\mathrm{m}$ and 3000 m . The width w of the abnormal body in the $y$ direction is set as $1000\mathrm{m},1500\mathrm{m},2000\mathrm{m},2500\mathrm{m}$ and 3000 m . The thickness k of the abnormal body in the $z$ direction is set to $200\mathrm{m},300\mathrm{m},400\mathrm{m},500\mathrm{m}$ and 600 m . The depth d is set to $100\mathrm{m},200\mathrm{m},300\mathrm{m},400\mathrm{m}$ and 500 m . The above models adopt $100\times 100\times 100\mathrm{m}$ discrete grid. The electromagnetic response calculation results at an emission frequency of 0.1 Hz are shown in Fig. 11 (a). The amplitude normalization curves of the electromagnetic response curves for different lengths, widths, thicknesses and depths of the abnormal body B and the electromagnetic response curves for the original size of the abnormal body B are shown in Fig. 11 (b).

It can be seen from Fig. 11 (a) that the electromagnetic response curves of abnormal body B with different lengths, widths, thicknesses and depths almost overlap with those of abnormal body B with original dimensions, indicating that the changes of abnormal body length, width, thickness and depth have a very weak influence on the static effect of the model. In addition, it can be seen from the amplitude normalization curves in Fig. 11 (b) that the width of the abnormal body has the greatest influence on the calculated results, followed by the thickness of the abnormal body, the length of the abnormal body and the depth of the abnormal body. With the increase of width, thickness, length and depth, the amplitude normalization value increases gradually.

Finally, we assume that the conductivity and central coordinates of abnormal body B are constant, and only the volume of abnormal body is changed. The dimensions of abnormal body B along the x, y and z directions are set as (l, w, k) (500 m, 500 m, 100 m), (1000 m, 1000 m, 200 m), (1500 m, 1500 m, 300 m), (2000 m, 2000 m, 400 m), (2500 m, 2500 m, 500 m), (3000 m, 3000 m, 600 m) and (4000 m, 4000 m, 800 m). The grid size and scaling frequency adopted by the model are consistent with those in Fig. 10. The calculated electromagnetic response at the emission frequency of 0.1 Hz is shown in Fig. 12 (a), and the amplitude normalization curves of the electromagnetic response curves for different volumes of the abnormal body B and the electromagnetic response curves for the original volume of the abnormal body B are shown in Fig. 12 (b).

Figure 12 (a) shows that the electromagnetic response curves of abnormal body B with different volumes almost overlap with those of abnormal body B with original volumes, indicating that the change of abnormal body volume also has a weak influence on the static effect of the model. However, it is not difficult to find from the amplitude normalization curves in Fig. 12 (b) that with the increase of the volume of abnormal body B, the amplitude normalization value gradually increases and, compared with the amplitude normalization curves of electromagnetic responses at different lengths, widths, thicknesses and depths in Fig. 11 (b), the increase trend of the amplitude normalization value is stronger, indicating that the static effect caused by volume change is more obvious than that caused by single-direction extension.

C. Synthetic model of uplift terrain and static effect

To study the influence characteristics of terrain and static effects on the marine CSEM response, we designed a synthetic model (M05) including terrain and static effect (Fig. 13). The central coordinates and conductivity of the reservoir, terrain and abnormal body are the same as that shown in M03 and M04.

We changed the volume of the uplift terrain and set the volume parameters, mesh size and scaling frequency consistent with Fig. 8. The electromagnetic response at the emission frequency of 0.1 Hz and the amplitude normalization of the response curve for different volumes of the uplifted terrain and the electromagnetic response curve for the original volume of the uplifted terrain are shown in Fig. 14.

Figure 13: Synthetic model of uplift terrain and static effect (M05).

Figure 14: Electromagnetic response curves and amplitude normalization curves of different volume uplift terrain in the synthetic model: (a) electromagnetic response curves of uplift terrain at different volumes and (b) amplitude normalization curves of electromagnetic responses of uplift terrain with different volumes.

Figure 15: Electromagnetic response curves and amplitude normalization curves of different volume abnormal body in the synthetic model: (a) electromagnetic response curves of abnormal body at different volumes and (b) amplitude normalization curves of electromagnetic responses of abnormal body with different volumes.

Figure 16: Electromagnetic response curves and amplitude normalization curves of the terrain and abnormal body at different volume in the synthetic model: (a) electromagnetic response curves of the terrain and abnormal body at different volumes and (b) amplitude normalization curves of electromagnetic responses of the terrain and abnormal body with different volumes.

Then, we changed the volume of the abnormal body and set the volume parameters of abnormal body B, mesh size and scaling frequency consistent with Fig. 12. Electromagnetic response at the emission frequency of 0.1 Hz and amplitude normalization of the response curve for different volumes of abnormal body B and the electromagnetic response curve for the original volume of the abnormal body B are shown in Fig. 15. Finally, we changed both the volumes of uplifted terrain and the abnormal body and set the volume parameters of uplifted terrain and abnormal body B, mesh size and scaling frequency consistent with Figs. 8 and 12. The electromagnetic response at the emission frequency of 0.1 Hz and the amplitude normalization of the response curve for different volumes of the terrain and abnormal body B and the electromagnetic response curve for the original volume of the terrain and abnormal body B are shown in Fig. 16.

It can be seen from Figs. 14 and 15 that the change of uplift terrain volume has a strong influence on the electromagnetic response of the synthetic model, and the amplitude normalization value gradually decreases with an increase of uplift terrain volume. The influence of the volume of abnormal body B on the static effect of the synthetic model is weak, but the amplitude normalization curve shows that the amplitude normalization value increases gradually with an increase of the volume of abnormal body B. The above results are similar to our previous conclusions. In addition, it can be seen from Fig. 16 that the influence of the simultaneous changes of terrain and abnormal body volume on the electromagnetic response of the synthetic model is equivalent to the superposition of the effects of the two changes separately.

V. CONCLUSION

In view of the uplift terrain such as submarine hills and seamounts shaped by plate tectonics and seafloor spreading and the static effect caused by polymetallic nodules and hydrothermal sulfide, this study adopts the FD method in fictitious wave domain to realize the 3D numerical simulation study of marine CSEM for the uplift terrain and static effect model. We compare the electromagnetic response of the whole-space model and the one-dimensional reservoir model with the analytical solution to verify the correctness and effectiveness of the algorithm.

We analyze and study the electromagnetic response through the forward simulation of the uplift terrain 3D reservoir model and the static effect model. The results show that the uplift terrain and the static effect produced by shallow abnormal body have different degrees of influence on the marine CSEM field. The uplift terrain distorts the electric field components and weakens the electromagnetic anomalies generated by oil reservoirs. In addition, compared with the top width, bottom width and height, the change of uplift terrain volume has a stronger influence on the electromagnetic response of the model. With an increase of the resistivity of the abnormal body, the influence of the static effect on the marine CSEM response increases gradually. Meanwhile, compared with length, width, thickness and depth, the static effect caused by a change of abnormal body volume is more obvious. Furthermore, results of the synthetic model corroborate the above conclusions, and the effect of the simultaneous change of terrain and abnormal body volume on the electromagnetic response of the synthetic model is equivalent to the superposition of the effect when the volume changes alone. Therefore, it is necessary to consider the effects of seafloor uplift terrain and static effects in marine CSEM exploration.

ACKNOWLEDGMENT

The authors would like to thank the anonymous reviewers and the editor for their valuable comments that greatly improved this article. This research was financially supported by the National Key Research and Development Program of China (Grant Number 2022YFC2807904).

REFERENCES

[1] S. Constable, “Ten years of marine CSEM for hydrocarbon exploration,” Geophysics, vol. 75, no. 5, pp. 75A67-75A81, Oct. 2010.

[2] S. Constable, “Review paper: Instrumentation for marine magnetotelluric and controlled source electromagnetic sounding,” Geophys. Prospect., vol. 61 Suppl, pp. 505-532, Jan. 2013.

[3] J. Liu, T. Guo, B. Wang, and Z. Guo, “Review of marine electromagnetic methods for hydrocarbon exploration” [in Chinese], Geophys. Prospect. Petro., vol. 60, no. 4, pp. 527-538, July 2021.

[4] E. Attias, K. Weitemeyer, S. Holz, S. Naif, T.A. Minshull, A.I. Best, A. Haroon, M. Jegen-Kulcsar, and C. Berndt, “High-resolution resistivity imaging of marine gas hydrate structures by combined inversion of CSEM towed and ocean-bottom receiver data,” Geophys. J. Int., vol. 214, no. 3, pp. 1701-1714, Sep. 2018.

[5] M. P. Miensopust, “Application of 3-D electromagnetic inversion in practice: Challenges, pitfalls and solution approaches,” Surv. Geophys., vol. 38, no. 5, pp. 869-933, Sep. 2017.

[6] Q. Wu, Y. Zhang, S. Guan, D. Li, and Y. Ji, “The 3D modeling of GATEM in fractured random media based on FDTD,” Appl. Comput. Electrom., vol. 36, no. 11, pp. 1401-1406, Nov. 2021.

[7] Y. Inoue and H. Asai, “Efficient electromagnetic simulation including thin structures by using multi-GPU HIE-FDTD method,” Appl. Comput. Electrom., vol. 33, no. 2, pp. 212-215, Feb. 2018.

[8] D. Yoon, M. S. Zhdanov, J. Mattsson, H. Z. Cai, and A. Gribenko, “A hybrid finite-difference and integral-equation method for modeling and inversion of marine controlled-source electromagnetic data,” Geophysics, vol. 81, no. 5, pp. E323-E336, Oct. 2016.

[9] J. Tang, F. Zhou, Z. Ren, X. Xiao, L. Qiu, C. Chen, and H. Chen, “3-D forward modeling of the controlled-source electromagnetic problem based on the integral equation method with an unstructured grid,” Chin. J. Geophys., vol. 61, no. 4, pp. 1549-1552, Apr. 2018.

[10] R. Streich, “3D finite-difference frequency-domain modeling of controlled-source electromagnetic data: Direct solution and optimization for high accuracy,” Geophysics, vol. 74, no. 5, pp. F95-F105, Oct. 2009.

[11] R. Mittet, “High-order finite-difference simulations of marine CSEM surveys using a correspondence principle for wave and diffusion fields,” Geophysics, vol. 75, no. 1, pp. F33-F50, Feb. 2010.

[12] H. Jahandari and C. G. Farquharson, “Finite-volume modelling of geophysical electromagnetic data on unstructured grids using potentials,” Geophys. J. Int., vol. 202, no. 3, pp. 1859-1876, July 2015.

[13] R. Peng, X. Hu, J. Li, and Y. Liu, “3D inversion of frequency-domain marine CSEM data in VTI media,” Chin. J. Geophys., vol. 62, no. 6, pp. 2165-2175, June 2019.

[14] V. Puzyrev, J. Koldan, J. de la Puente, G. Houzeaux, M. Vazquez, and J. M. Cela, “A parallel finite-element method for 3-D controlled-source electromagnetic forward modelling,” Geophys. J. Int., vol. 193, no. 2, pp. 678-693, May 2013.

[15] H. Chen, T. Li, B. Xiong, H. Wang, R. Zhang, and S. P. Li, “3D MCSEM modeling using an edge-based finite element method based on an unstructured grid and incremental model,” Chin. J. Geophys., vol. 61, no. 6, pp. 2560-2577, June 2018.

[16] R. U. Börner, “Numerical modelling in geo-electromagnetics: Advances and challenges,” Surv. Geophys., vol. 31, no. 2, pp. 225-245, Oct. 2009.

[17] K. H. Lee, G. Liu, and H. F. Morrison, “A new approach to modeling the electromagnetic response of conductive media,” Geophysics, vol. 54, no. 9, pp. 1180-1192, Sep. 1989.

[18] A. T. de Hoop, “A general correspondence principle for time domain electromagnetic wave and diffusion fields,” Geophys. J. Int., vol. 127, no. 3, pp. 757-761, Dec. 1996.

[19] F. A. Maaø, “Fast finite-difference time-domain modeling for marine-subsurface electromagnetic problems,” Geophysics, vol. 72, no. 2, pp. A19-A23, Feb. 2007.

[20] R. Mittet, “Seismic wave propagation concepts applied to the interpretation of marine controlled-source electromagnetics,” Geophysics, vol. 80, no. 2, pp. E63-E81, Feb. 2015.

[21] J. Lu and Y. G. Li, “3-D marine CSEM modeling in fictitious wave domain,” Chin. J. Geophys., vol. 62, no. 8, pp. 3189-3198, Aug. 2019.

[22] Y. Ji, X. Meng, and G. Ren, “Shallow water three-dimensional transient electromagnetic modelling by using fictitious wave field methods,” Appl. Comput. Electrom., vol. 35, no. 1, pp. 72-81, Jan. 2020.

[23] S. Yutaka, “Bathymetric effects and corrections in marine CSEM data,” Geophysics, vol. 76, no. 3, pp. F139-F146, Apr. 2011.

[24] B. Yang, Y. Xu, Z. He, and W. B. Sun, “3D frequency-domain modeling of marine controlled source electromagnetic responses with topography using finite volume method,” Chin. J. Geophys., vol. 55, no. 4, pp. 1390-1399, Apr. 2012.

[25] B. Yan and B. Han, “Bathymetric analysis and corrections for 3-D marine controlled-source electromagnetic field” [in Chinese], Comput. Tech. Geophys. Geochem Explor., vol. 39, no. 1, pp. 9-16, Jan. 2017.

[26] X. Yang, M. Yue, D. Hu, Y. Li, and X. Wu, “Goal-oriented 3-D time-domain marine CSEM modeling with anisotropy and topography,” IEEE Trans. Geosci. Remote Sens., vol. 60, no. 1, pp. 1-14, Jan. 2022.

[27] T. V. Carlos and X. B. Francis, “Principles of spatial surface electric field filtering in magnetotellurics: Electromagnetic array profiling (EMAP),” Geophysics, vol. 57, no. 4, pp. 603-622, Apr. 1992.

[28] J. Tang and J. He, Controlled-Source Audio-Frequency Magnetotelluric Method and its Application. China: Central South University Press, 2005.

[29] L. Xia, “Response characteristic analysis and recognition and correction of static effect of MT sounding,” M.S. thesis, School of Geophysics and Information Technology, China University of Geosciences (Beijing), Beijing, China, June 2016.

[30] Y. Li and M. Kris, “Denoising multicomponent CSEM data with equivalent source processing techniques,” Geophysics, vol. 78, no. 3, pp. E125-E135, Apr. 2013.

[31] Y.Q. Yang, “Research on magnetotelluric static effect suppression and its application based on wavelet threshold,” M.S. thesis, Faculty of Land Resources Engineering, Kunming University of Science and Technology, Kunming, Yunnan, China, May 2021.

[32] Y. Kane, “Numerical solution of initial boundary value problems involving maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag., vol. 14, no. 3, pp. 302-307, May 1966.

[33] E. C. Du Fort and S. P. Frankel, “Stability conditions in the numerical treatment of parabolic differential equations,” Math. Tables Aids Comput., vol. 7, no. 43, pp. 135-152, May 1953.

[34] F. Kong, Electromagnetic Field Calculation and MATLAB Implementation of Layered Media. China: Jiangsu Phoenix Science and Technology Press, 2016.

[35] B. Han, X. Hu, Y. Huang, R. Peng, J. Li, and J. Cai, “3-D frequency-domain CSEM modeling using a parallel direct solver,” Chin. J. Geophys., vol. 58, no. 8, pp. 2812-2826, Aug. 2015.

BIOGRAPHIES

Chunying Gu (Graduate Student Member, IEEE) received the M.S. degree from the College of Mechanical and Vehicle Engineering, Hunan University, Changsha, China, in 2011. She is currently pursuing the Ph.D. degree in detection technology and automation from Jilin University, Changchun, China. Her research interests include 3D marine controlled-source electromagnetic (CSEM) modeling and data processing.

Suyi Li (Member, IEEE) received her M.Sc. and Ph.D. degrees both from Jilin University in 2002 and 2009, respectively. From 2008 to 2009, she studied at the University of Illinois at Urbana-Champaign as a joint Ph.D. Now she is a professor in Jilin University. Her main research interests include computer applications and digital signal processing.

Wanyue Zhang (Graduate Student Member, IEEE) received the B.S. degree from the College of Instrumentation and Electrical Engineering, Jilin University, Changchun, China, in 2022, where she is currently pursuing the M.S. degree. She is primarily engaged in marine controlled-source electromagnetic (CSEM) data processing.

Silun Peng received the Ph.D. degree from the College of Automotive Engineering, Jilin University, Changchun, China, in 2014. He is currently a deputy senior engineer with Jilin University. His research interests include, but are not limited to, computer applications, digital signal processing and hardware system design.

ACES JOURNAL, Vol. 39, No. 6, 505–519
doi: 10.13052/2024.ACES.J.390604
© 2024 River Publishers

ABSTRACT

I. INTRODUCTION

II. FICTITIOUS WAVE DOMAIN FD METHOD

A. Construction of the governing equations in the fictitious wave domain

B. Solution of electromagnetic responses in the fictitious wave domain

C. Recovery of electromagnetic response in the diffusive domain

III. ACCURACY VERIFICATION

A. 1D oil reservoir model

B. 3D oil reservoir model

IV. NUMERICAL EXPERIMENTS

A. Uplift terrain model

B. Static effect model

C. Synthetic model of uplift terrain and static effect

V. CONCLUSION

ACKNOWLEDGMENT

REFERENCES

BIOGRAPHIES