Three-Dimensional Numerical Visualization And Simulation of Multiphysics in Point Source Dc Magnetometric Resistivity Method

Magnetometric Resistivity(MMR) Method, which can use the power supply method of the traditional apparent resistivity method to measure the magnetic field. At present, the application research abroad is relatively extensive, but the domestic(China) research on the application of the MMR method is very few, and it is not even well-known. This paper is based on the MMR theoretical method under the point source DC condition, combined with the abnormal potential method and the modified Biot-Savart law, and using the three-dimensional numerical calculation method of the finite difference method to calculate the abnormal potential field, electric field, and magnetic field on the matlab2018a platform. Calculate, realize the multi-physics simulation in the electromagnetic field through Matlab platform programming, and verify the correctness of the algorithm by the spherical anomaly model with an analytical solution. Through the visual simulation of the threedimensional data volume of the multi-physics field in the electromagnetic field, we can better understand the response mechanism of the electromagnetic field under DC conditions and grasp their three-dimensional spatial distribution rules. It is hoped that the research in this article can help the research of MMR Personnel better use this method for exploration.


Introduction
The magnetometric resistivity(MMR) method is a type of electromagnetic exploration system that measures the magnetic field excited by artificial source non-inductive (direct current or pseudo direct current) current, rather than potential difference [1,2]. When a point source direct current is introduced into the earth, the traditional apparent resistivity method measures the potential difference. However, on the conductivity interface containing three-dimensional anomalous bodies, due to the change of the conductivity on the discontinuous interface, an abnormal magnetic field will be caused. For abnormal magnetic fields, the MMR method needs to be studied. The MMR method under DC conditions has a lot of applied research abroad [3][4][5][6], but there are few kinds of research on the application of this method in China. As early as 1958, Stefanescu proposed a modified Biot-Savart law based on horizontal terrain [7], which can calculate the magnetic field by potential or current density. In 1978, Professor R.N. Edwards also proved it [8].
Compared with the abstract field concept and the corresponding mathematical equations, the visual representation of the three-dimensional field can give researchers an intuitive experience. However, due to the constraints of the computer and other technical conditions at the time, the corresponding three-dimensional visual simulation of the electromagnetic field rarely appears. But now, with the rapid development of computer technology, the Matlab platform is a powerful computing software integrating calculation and visualization, so this article uses Matlab's powerful computing functions and data visualization functions, as well as simplified algorithms, which is very simple to realize the three-dimensional visualization simulation of the relevant electromagnetic field excited by the point source current in the MMR method (the visualization also combines the VTK format data). In China, in 1991, Chinese researcher Yang Jianwen proved that the revised Biot-Saffar law is also applicable to undulating terrain [9]. In 1993, Professor Ruan Baiyao gave an arbitrary twodimensional geoelectric cross-section point current source underground conduction current magnetic field fraction expression [10]. In 2017, Professor Weng Aihua analyzed the basic principles of the magnetoelectric method, the development status at home and abroad, and its application prospects [11], making the magnetoelectric method more familiar to the Chinese people.
Based on previous theoretical studies, this paper calculates the abnormal potential by the numerical calculation method of the abnormal potential method and the finite difference method [12][13][14] and calculates the magnetic field by the integral equation: modified Biot-Savart law. At the same time, the three-dimensional visual simulation processing of the relevant electromagnetic field in the MMR method under the point source condition was programmed on the Matlab platform. That is, the abnormal potential Ua is calculated by the abnormal potential method, and then the three components of the abnormal electric field E in the X, Y, and Z-axis directions are calculated. The electric field data is used to obtain the abnormal magnetic field H through the integrated Biot-Savart law. Finally, use the visualization toolkit in Matlab(Also combined with "Paraview") to visualize the potential field, electric field, and magnetic field. At the same time, a spherical anomaly model with an analytical solution was used to verify the accuracy of the algorithm. From the simulation results, we can see that the distribution of the electromagnetic field is mainly concentrated near the anomalous body, and their response in areas far from the boundary is relatively small, and the areas where the anomalous electromagnetic field and the abnormal potential field respond are very beautifully symmetrical.

Abnormal Potential Method
In the three-dimensional Cartesian coordinate system, the field source in the homogeneous semiinfinite space under the ground is a point current source whose coordinates are on the ground( 0 , 0 , 0 ) and the current intensity is I. The solution space model is shown in Figure 1, (1) and (2) regions Respectively represent the anomalous body space and background space of the conductivity distribution in the borehole and surrounding rock. Using the abnormal potential method, the potential is divided into the background potential 0 and the abnormal potential , namely = 0 + , represents the actual conductivity σ minus the background conductivity Rates 0 , namely = σ − 0 , T0, T1...T5 respectively represent the upper boundary, left boundary, front boundary, right boundary, lower boundary, and back boundary of the space model. According to electromagnetic theory, through derivation, the differential equation and its definite solution conditions for the abnormal potential can be obtained: Where = ( , , ), = ( , , ), 0 , and are known. When the point source is on the ground, wa=2π, and when the point source is underground, wa=4π, so only is an unknown number. Solving the equations can solve the abnormal potential 0 .
The abnormal electric field E is:

Modified Biot-Savart law
The revised Biot-Savart law, that is, the magnetic field at any point 0 ( 0 , 0 , 0 ) in space is: Where V represents the entire solution space, that is, the sum of space area (1) and space area (2) in the model diagram, ( , , ), ( , , ) is the value of any point in space Abnormal potential and conductivity. By expanding the above formula, the integral expressions of , , can be obtained, especially when performing full-space calculations, you need to pay attention to the "singularity". For the detailed process, please refer to the corresponding literature [9,10].
From the equation, we can see that the magnetic field generated by the current is closely connected to the potential generated by the current on the interface of the electrical inhomogeneity. In this way, the calculation problem of the conduction current magnetic field is transformed into the calculation of the potential. Therefore, some basic theories and numerical calculation methods in traditional electrical exploration can be used to study the calculation of the conduction current magnetic field.

Finite Difference Method for Solving Abnormal Potential
Before performing numerical calculations, we need to discretize the solution space (using grid nodes). At the same time, we can use a finer grid to discretize the abnormal body and the space around it and use sparseness for the solution area outside the boundary of the abnormal body. Use (i, j, k) instead of (x, y, z) to number the nodes, and the discrete space shape is shown in Figure 2 ( Figure 2 is a schematic diagram of the three-dimensional space 20m*20m*20m grid discrete ).

Difference scheme of abnormal potential method
Through the finite difference discretization of the differential equation (1) that the abnormal potential satisfies, the linear equation is obtained: In the numerical calculation, we need to pay attention to the calculation of the node potential within the boundary and the node potential on the boundary. This article only introduces the main equations used. When the node (i, j, k) is within the boundary, the matrices [ ], [ ] and{ } in the linear equation (4)  When the node (i, j, k) is located on the ground T0 boundary, because the ground node is a half unit, there is no (i, j, k-1) node in the space, and the field source is located on the ground, on the ground boundary, According to the principle of mirror image method, the current intensity of the corresponding node current source will be doubled.

Accuracy verification
We use the classical spherical anomaly with potential analytical solution to analyze the relative error of the numerical solution and the analytical solution, regardless of the magnetic field or the electric field, we are based on the abnormal potential obtained by the abnormal potential method in combination with Equation (2) and Equation (3). Therefore, we only need to analyze the relative error between the total potential numerical solution obtained by the abnormal potential and the potential obtained by the analytical solution(and the magnetic field is only generated on the discontinuous electrical interface, so the most important thing is that we need to analyze the accuracy of the potential near the abnormal body. If the potential obtained by the abnormal potential method and the analytical solution is used to obtain the relative error of the potential to meet our needs, It can be said that the algorithm is correct), to verify the accuracy of the algorithm. It can be seen from Figure 3 that the relative error is 0-0.45%, and the maximum error is less than 0.45%, which shows that the algorithm is correct and has high accuracy.

Simulation of abnormal potential
We use Section 3 "Finite Difference Method to Solve Abnormal Potential" to calculate the threedimensional abnormal potential (unit: V), and then combine the visualization toolkit in matlab2018a and the corresponding visualization software (such as "Paraview", etc.) to perform the three-dimensional data volume visual processing. The obtained static abnormal potential field is shown in Figure 4. The red cube is the outline of the abnormal body. The volume diagram of the abnormal potential presents its three-dimensional spatial distribution and its changing law. The point source is located directly above the cubic anomaly. The response area near the body is mainly concentrated on the top and ground of the abnormal body, and the response of the abnormal potential on the top surface of the abnormal body is negative, and the response on the bottom surface is positive. Figure 4. The range of the space is 20m*20m*20m, the background resistivity is 300ohm·m, the size of the abnormal body is 2m*2m*2m, the resistivity of the abnormal body is ρ=30ohm·m, the point source is located directly above the abnormal body, the current intensity I=3A, the conductance The rate σ=1/ρ; the red cube is the contour of the abnormal body, and the Volume diagram is the threedimensional abnormal potential distribution.

Simulation effect of electromagnetic field induced by abnormal body in MMR method
From the abnormal potential obtained in section 4.2, combine equation (2) and equation (3) to obtain the electric field E and the magnetic field H. Through the simulation image, as shown in Figure  5-6 (6a: remove the abnormal body, observe the electromagnetic field in the abnormal body The internal spatial distribution and its changing law; 6b: Amplify the electromagnetic field near the anomalous body and observe the distribution of electromagnetic and magnetic fields around and within the anomaly). The direction of the magnetic field is perpendicular to the electric field, and the magnetic field is mainly distributed around the boundary of the anomaly, and the electric field passes through the anomaly Inside, in the entire spatial distribution, the length of the arrow reflects the magnitude of the magnetic field or electric field value. Through the simulation diagram, we can see the three-dimensional spatial distribution and direction of the abnormal electric or magnetic field (or according to the current density j=σE, you can see The distribution of the "induced" current density). Turning abstract physical concepts into visual images can give researchers a more intuitive representation of the physical field. As shown in Figure 7, when only describing the vector diagram of the abnormal electromagnetic field in the XY direction, we can intuitively understand the electromagnetic field changes on the ground or near the surface, and it can also help MMR prospectors better arrange the MMR acquisition device. green arrow: abnormal magnetic field vector direction; volume diagram: three-dimensional electric potential distribution diagram, color scale reflects the magnitude of electric potential value, you can see its three views in Appendix B Figure A1, Figure A2 and Figure A3.  Usually, we measure the electric field and magnetic field components, so it is meaningful to fully understand the spatial distribution of each component of the electromagnetic field and its changing law. In section 4.3.1 and section 4.3.2 of this article, we visualize the three-dimensional specific shape of each component of the electric field and magnetic field.

Electric field induced by anomalous body
According to equation (2), the gradient of the three-dimensional abnormal potential is obtained, and the three-dimensional data volume of each component of the abnormal electric field E (unit V/m) is obtained. At the same time, the visualization toolkit of matlab2018a (combined with "Paraview", etc.) is used to perform three-dimensional visual simulation of the electric field. As shown in Figure  8 (after extracting the abnormal potential field and magnetic field, analyze the spatial distribution of the electric field and its changing law separately), the arrow size of the electric field vector represents the relative magnitude of the electric field value in the space. In Figure 9 (9a: represents the threedimensional contour map of the electric field components Ex, Ey, Ez; 9b: represents the threedimensional contour map of Ex, Ey, Ez; 9c: represents the center of the abnormal body and the data volume Ex, Ey, Ez. For smooth slices, the slice position can be extracted from the coordinates of the slice.). We use three-dimensional contours to draw the specific shape of the three-dimensional data volume of each component of the electric field. The simulation image looks like a "ball of wool", it is not smooth. However, it vividly presents the continuous change of the electric field component in space. At the same time, comprehensively in Figure 9, we can intuitively see that the abnormal distribution of Ex and Ey is at an angle of 90 , and we can roughly judge by Ez Out of abnormal body shape.

Magnetic field induced by anomalous body
From equation (3), the components of the magnetic field H (unit: 10 -8 A/m) in the X-axis, Y-axis, and Z-axis directions are obtained through the integration of the modified Biot-Savart law, through the visualization toolkit in matlab2018a (Combined with "Paraview", etc.), the visualization of the three-dimensional vector direction of the magnetic field H is realized. As shown in Figure 10, combined with Figure 5, we can see that the magnetic field is mainly distributed near the boundary of the anomaly, that is, near the discontinuous conductivity interface. The size of the arrow represents the relative size of the magnetic field. As shown in Figure 11 (11a: represents the three-dimensional contour map of the electric field components Hx, Hy, Hz; 11b: represents the top view of the threedimensional contour of Hx, Hy, Hz; 11c: represents the data volume Hx, Hy, Hz through the center of the abnormal body For smooth slices, the position of the slice can be seen from the coordinates of the slice), we also use the form of "wool ball" to represent the continuous distribution of magnetic field components in space. Hx and Hy are similar in shape to electric fields Ex and Ey, but the magnetic field can better reflect the boundary of the discontinuous decomposition surface of anomalous conductivity, and Hz can better explain this problem. The shape of the anomalous body and its boundary is clearer than that of the electric field. Therefore, by measuring the magnetic field, we can better describe the shape and location of the underground anomaly. This also shows the superiority of the MMR method. The disadvantage is that the magnetic field value is too small, which is several orders of magnitude different from the electric field value. High-precision measuring instrument It is susceptible to interference during measurement, but this problem is being solved [5].

Conclusions
In this paper, through matlab2018a platform programming, the simulation and visualization of the related potential field, magnetic field, and electric field induced by the abnormal body in the MMR method under the point source DC condition are realized. Through simulation, we visualized the physical relationship between electromagnetic fields, and intuitively expressed the physical relationship between the potential field and electric field, and between electric field and magnetic field.
Through Figure 4-11, we can see that the magnetic resistivity method is more advantageous in identifying the size and shape of anomalous bodies. At the same time, we can also conclude that when the anomalous resistivity is low relative to the background, the response of the abnormal potential near the point source is negative, and the response on the side away from the point source is positive; the abnormal electric field Ex, Ey is at There will be four abnormal response areas around the abnormal body, and positive and negative abnormalities appear in pairs. Ez will have three abnormal response areas, the negative abnormal response is in the middle, and the positive abnormal response and the negative abnormal area are continuously distributed in the upper, middle, and lower parts; abnormal magnetic field Hx, Hy has two abnormal response areas around the abnormal body, Hz has 8 abnormal response areas, and the positive and negative anomalies are paired and distributed continuously.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix B
The three views of Figure 5