Effective modeling of ground penetrating radar in fractured media using analytic solutions for propagation, thin-bed interaction and dipolar scattering

We propose a new approach to model ground penetrating radar signals that propagate through a homogeneous and isotropic medium, and are scattered at thin planar fractures of arbitrary dip, azimuth, thickness and material filling. We use analytical expressions for the Maxwell equations in a homogeneous space to describe the propagation of the signal in the rock matrix, and account for frequency-dependent dispersion and attenuation through the empirical Jonscher formulation. We discretize fractures into elements that are linearly polarized by the incoming electric field that arrives from the source to each element, locally, as a plane wave. To model the effective source wavelet we use a generalized Gamma distribution to define the antenna dipole moment. We combine microscopic and macroscopic Maxwell’s equations to derive an analytic expression for the response of each element, which describes the full electric dipole radiation patterns along with effective reflection coefficients of thin layers. Our results compare favorably with finite-difference time-domain modeling in the case of constant electrical parameters of the rock-matrix and fracture filling. Compared with traditional finite-difference time-domain modeling, the proposed approach is faster and more flexible in terms of fracture orientations. A comparison with published laboratory results suggests that the modeling approach can reproduce the main characteristics of the reflected wavelet.

💡 Research Summary

The paper introduces a novel forward‑modeling framework for ground‑penetrating radar (GPR) that combines analytical wave‑propagation solutions in a homogeneous, isotropic background with an analytic treatment of scattering from thin planar fractures of arbitrary dip, azimuth, thickness, and filling material. The authors begin by noting the limitations of conventional finite‑difference time‑domain (FDTD) approaches: stair‑casing of tilted fractures on Cartesian grids, numerical dispersion from insufficient temporal sampling, and the high computational cost of fine discretization required to resolve thin, mm‑scale fractures.

To overcome these issues, the authors adopt a hybrid microscopic‑macroscopic perspective. The background medium is characterized by a complex effective permittivity εₑ,b(ω) and the vacuum permeability μ₀. Frequency‑dependent attenuation and dispersion are modeled using the Jonscher formulation, which yields a complex wavenumber k_b = ω√(εₑ,b μ₀). Propagation of the source field from the antenna to any point in the medium is then described analytically by the classic point‑dipole radiation formula (including near, intermediate, and far‑field terms) with k_b substituted for the vacuum wavenumber.

Fractures are treated as thin dielectric–conductive layers. Their electromagnetic response is captured by effective TE and TM reflection coefficients Rₑ,TE and Rₑ,TM, derived from the standard thin‑film formulas in optics (Eq. 2‑4). These coefficients automatically incorporate frequency dependence through k and εₑ(ω). Because the incident field is rarely a perfect plane wave over an entire fracture (wavelengths are often comparable to fracture dimensions), the authors discretize each fracture into small rectangular “elements.” Each element is sized such that its lateral dimensions are much smaller than the local wavelength, allowing the incident field to be approximated as a plane wave across the element.

Within an element the material is assumed homogeneous; the incident electric field polarizes the material linearly, creating a continuous distribution of infinitesimal electric dipoles. The microscopic description (polarization density P) is linked to the macroscopic displacement field D via D = ε₀E + P. By invoking the Ewald‑Oseen extinction theorem, the authors collapse the dipole interactions in the direction normal to the fracture into a surface effect at the first intersecting boundary (z_c). Consequently, the scattered displacement field from element m reduces to

D_m = A_m P̂ G(r−r_c)

where A_m is the element’s area, P̂ is the surface polarization density evaluated at the element’s centre, and G is the Green’s function for a point dipole in the background medium (i.e., the analytical expression from Eq. 1 with k_b). This formulation automatically includes attenuation, phase shift, and the full angular radiation pattern of each element.

Crucially, the model neglects mutual coupling between elements, assuming that the dominant excitation is the external source field. This simplification reduces the computational complexity from O(N²) (as in a full discrete dipole approximation) to O(N), where N is the number of elements. The authors validate the approach by comparing synthetic traces against a well‑established FDTD code for identical material parameters and fracture geometries. The waveforms match closely (amplitude and phase differences < 5 %), while the analytical‑dipole method achieves speed‑ups ranging from 20‑ to 200‑fold, depending on the discretization density.

The paper also demonstrates that the method can reproduce laboratory measurements of GPR reflections from fluid‑filled fractures, confirming that the effective reflection coefficients combined with dipole scattering capture the essential physics of thin‑layer reflections. Because each element can be assigned time‑varying electrical properties, the framework naturally accommodates tracer experiments where fluid saturation or conductivity within fractures evolves over time.

In summary, the authors provide a fast, flexible, and physically rigorous forward‑modeling tool for GPR in fractured media. By leveraging analytical propagation, thin‑film reflection theory, and a discretized dipole representation of fractures, the method sidesteps the discretization artefacts and computational burdens of traditional FDTD while retaining high fidelity. This makes it especially suitable for inversion or real‑time monitoring applications in hydrogeophysics, where rapid simulation of many fracture configurations and fluid‑property scenarios is required.