A study on Parallel Computation Based on Finite Element Forward Modeling of 2D Magnetotelluric ()
1. Introduction
Magnetotelluric sounding method is a natural electromagnetic field exploration method [1] . The natural field contains different frequencies from high to low and electromagnetic wave has different penetration depth with different frequencies [2] , so the method can explore the depth of the target [3] . Because the method has many advantages such as exploring deeply, avoiding the impact of high resistivity body and being sensitive to the low resistivity layer, more and more people do research on it [1] .
With the development of computer technology, 2D magnetotelluric forward modeling is becoming mature; 2D finite element [4] [5] magnetotelluric forward can get more accurate results and solve the rough terrains problem. Zeng has developed the algorithm of magnetotelluric forward modeling [3] in 2008. Xie has developed the algorithm of magnetotelluric forward modeling [6] with terrain in 2007. Liu has developed the algorithm of magnetotelluric forward modeling [7] based on second field. However, it costs more time when the 2D magnetotelluric program computes relatively big grid and the forward algorithm is called for dozens of times; thus, we need to develop an effective algorithm.
We have developed the parallel algorithm in order to solve the problem of forward modeling in the paper. I will introduce the algorithm.
2. Methodology
2D finite element magn (1)
etotelluric forward modeling.
Maxwell equations
(2)
X axis is parallel with the target, Y axis is vertical to the target, Z axis is vertical to the ground. When the angular frequency is and time factor is [3] , the field equation is shown as follows [8] .
TE mode:
(3)
(4)
(5)
TM mode:
(6)
(7)
(8)
is the electric field toward x direction, and are the magnetic field toward y direction and z direction.
I will introduce the theory in TE mode.
In theory of electromagnetic fields the electromagnetic fields is divided into the primary electromagnetic field and second electromagnetic field. The primary electromagnetic field is the field when there is no target in the ground. The second electromagnetic field is the field caused by the target. is the primary electric field of x direction, and is the primary magnetic field of y and z direction. When there is no target in the ground, the equation is shown as follows.
(9)
(10)
(11)
When there is a target in the ground, total field equation minus primary field equation equals second field [9] equation. The resistivity of target minus the background resistivity equals. is the second electric field of, is the second magnetic field of, is the second magnetic field of
(12)
(13)
(14)
Substituting Equations (12) and (13) into Equation (14), we get the helmohoz equation.
(15)
Multiply both side of the equation by v, integrate the equation.
(16)
When there is no target in the ground , we solve the 1D magnetotelluric equation [10] and get.
Divide the ground into M × N grid and every rectangle is divided into 4 triangles. Through isoparametric element transformation of the equation, we put the element into the matrix A [11] . We get the equation A
and solve the equation, we get.. Substituting into equation 3, 4,
5, we get, and the apparent resistivity of every observed point.
3. Experiment
3.1. The Overall Design of Parallel Computation of 2D Magnetotelluric
Through analyzing the algorithm of forward, for example 9 frequencies, because the frequency is different, the matrix A is different in equation Ax = b. We need to solve equation for 9 times. In parallel algorithm (such as 4 processes), every process needs to solve equation for 3 times at most. However, in serial algorithm (only 1 process) the process needs to solve equation for 9 times. The less the times (the process needs to solve equation) are, the program running time is less, the parallel efficiency is higher.
The 9 frequencies are from freq (1) to freq (9). The value is 100, 31.6, 10, 3.16, 1, 0.316, 0.1, 0.0316, 0.01 hz. The higher the frequency is, the less the time of solving the equation is. There are two kinds of processes in the algorithm. One type is the main process, the other type is subprocess. 0 process is the main process. It aims to distribute the tasks, broadcast the global data, gather the result and output the files. The subprocesses aims to receive data from the main process, perform tasks and send the result to the main process. In order to perform more tasks, the main process performs a small task. We will introduce the performing process.
1) MPI_INIT(), Init the parallel environment, the 0 process reads the frequency, model file and the observed data. MPI_Bcast() [12] , it broadcasts the data to the other processes.
2) Table 1, the relevant frequencies of a process, is shown as follows. All processes are assigned to calculate the data of the relevant frequencies separately. We get the primary field and Substitute into Equation (9) and (10) to get,. Through isoparametric element transformation of the equation, we get the matrix A of the relevant frequencies.
Table 1. The relevant frequencies of a process.
3) we solve the equation A, by gauss method and we get the second electric field of every observed point. Substituting into Equation (12) and (13), we get,., ,. The apparent resistivity is, All processes are assigned to get apparent resistivity of the relevant frequencies separately. We get the value of m_mpi in each process. MPI_gatherv() [12] , 0 process gathers m_mpi of other processes in m_all, m_all has
9 m
_mpi and the sequence is 1, 7, 2, 3, 6, 4, 8, 5, 9. We sort the m_all and get the m.
4) 0 process writes the resistivity data to the file, MPI_FINALIZE() finalize the MPI environment.
3.2. The Program Is Developed in the Environment
OS: windows xp 64 bit;
CPU: intel core i7 3.2 GHz support 8 processes;
Memory: 8 GB;
Develop language: fortran;
Compiler: Compaq fortran;
Parallel environment: mpich2;
Command: Mpiexec-np N./mt2d N is the number of processes.
4. Result and Discussions
4.1. Low Resistivity Model
The depth of the target is 1000 m from the ground, the size is 500 m × 250 m. Background resistivity is 100 W×m. The resistivity of the target is 10 W×m. The size of the model’s grid is 100 × 100; the observed position is 32, 34, 36, 38, 40, 42, 44, 46, 48, 50, 52, 54, 56, 58, 60, 62, 64, 66, 68, 70 point on the ground. 9 frequencies are from 0.01 Hz to 100 Hz. We get the forward results from the parallel program and compare the forward result of two programs when the frequency is 10 hz in Figure 1.
In Figure 1, the horizontal axis is the serial number of observed point and the vertical axis is the apparent resistivity. express the result of parallel computation as o,express the result of normal computation as +, the result is the same.
4.2. Validity of the Result
The parallel computation is based on the serial program. According to the characteristics, that program solves equation for different frequency separately; it distributes the tasks to all processes and never changes other algorithm, so the parallel computation result is the same with the serial program’s result. The study compares the forward result of two programs in Figure 1; the result is the same, so it proves the validity of the program.
4.3. Discussions of Parallel Efficiency
In order to evaluate the efficiency of the parallel program for different account of processes, we calculate parallel speedup and parallel efficiency. The running time of serial program divided by the running time of parallel for N processes is Parallel speedup, parallel speedup devided by the number of processes N is parallel efficiency. Table 2 is the statistics of the running time for the algorithm of 2D magnetotelluric forward modeling in TE
Table 2. The statistics of the running time for the algorithm.
mode.
Analysis of Table 2 shows that the efficiency is the highest when the number of processes is 2. The speedup changes slightly and the efficiency declines when the number of processes changes from 6 to 8. Because the communication of the processes occupy more time and both of the second process need to solve equation for 2 times when the number of processes increases from 6 to 8. We can see that the effect of parallel algorithm is very obvious. We will do further study for the efficiency of the algorithm.
5. Conclusion
The computation of 2D finite element magnetotelluric forward for relatively big grid spends much time and the forward algorithm is called for dozens of times, so the key to the problems is parallel computation. The study realizes the parallel algorithm for 2D finite element magnetotelluric forward in the MPI environment. The algorithm is proved correct and efficient. The study lays the foundation for the parallel computation for 2D finite element magnetotelluric forward and inversion.