2D DC Resistivity modeling for anisotropic subsurface with variable topography based on mimetic finite-difference method

Abstract

The Direct Current (DC) resistivity method is a classical geophysical method to obtain subsurface geoelectrical images. This technique is utilized in groundwater, mineral mapping, subsurface pollution monitoring, saltwater intrusion and other civil engineering applications, where a reliable data analysis requires a versatile and robust forward modeling algorithm. The present study develops a 2-dimensional (2D) DC resistivity forward modeling algorithm employing mimetic finite difference methods (MFDM). The MFDM preserves the valuable properties of the continuum governing partial difference equation in discrete space, leading to a better representation of actual electrical potential by simulated potential. This study presents the first application of MFDM for DC resistivity modeling. The accuracy of the developed scheme is benchmarked utilizing analytical responses of a dyke model and two-layer anisotropic models. Since there are no analytical solutions for the variable topography cases, the accuracy of the scheme is demonstrated by comparing the solution with the published responses. A three-layer model is used to examine the stability of the devised algorithm by incorporating non-orthogonal grids. Non-orthogonal grids are produced by randomly varying the nodal coordinate of orthogonal grids. The observed error trends show that the algorithm is highly stable with regard to grid distortion and can accurately simulate complicated models involving topography and anisotropic subsurface. Furthermore, the numerical computation time analysis reveals that the developed algorithm is computationally stable to grid distortion. To efficiently accommodate the 3D character of the source in a 2D DC resistivity modeling, a new space domain approach is devised in this study. The developed algorithm is valuable in the case of long current-potential electrode spacing, including the case of high resistivity contrast and anisotropic subsurface. A half-space model was employed to examine the limitations of various wavenumber schemes. In the case of wavenumber-based modeling techniques, these wavenumber schemes are used in inverse cosine transform needed for space domain computation. It was observed that after a particular offset, all the wavenumber schemes deviated from the analytical result except the Gauss quadrature method with 120 wavenumbers, suggesting the long offset data requires simulation of large wavenumbers. Consequently, the wavenumber scheme becomes computationally expensive in the case of long offset simulation. Motivated by this analysis, a space domain modeling algorithm that utilizes a new boundary condition applicable to the plane that passes through the source position is developed. The proposed approach is computationally competitive with the wavenumber domain approach. It is likely to be even more efficient in case of large offsets as a small number of grids are sufficient to discretize the space in the strike direction. Extensive numerical simulations are carried out to demonstrate that the developed method is reliable and versatile for deep imaging surveys having variable topography and anisotropic subsurface, including tilted transversely isotropic cases. The construction of an algorithm based on a modified boundary condition, which aids in overcoming the wavenumber problem and gives accurate solutions for huge offsets, can be used as an essential tool to study responses for geophysical models till large offsets. Hence, we use the developed algorithm to obtain azimuthal apparent resistivity curves for various two-layered models, including isotropic, tri-axial anisotropic, and tilted transversely isotropic (TTI) models. The simulated azimuthal apparent resistivity plots provide insight into the scenarios when these plots behave like an isotropic case, even in the case of anisotropic subsurface. Further, the sensitivity curves are generated by taking the derivatives of apparent resistivity values with respect to the parameters that govern anisotropy. It is found that for a 2D case, the DC data generally shows sensitivity to all four parameters governing the anisotropy, which include the three principal resistivity values and one angle defining the angle of the tilted symmetry of anisotropy. However, when the azimuthal apparent resistivity plot evolves circularly, the principal resistivity along the profile direction becomes insensitive to the observed data.