TY - JOUR
AU - Deng,, Li
AB - Abstract Ground penetrating radar (GPR) is a very important tool for urban investigations. We introduce a finite difference time domain (FDTD) method to solve Maxwell's equations when the dispersive earth medium can be simulated. We give an iterative method for electrical displacement and magnetic field in a Debye medium, and a method to compute electric field from electrical displacement. The physical parameters are used only in calculating the electric field. The examples demonstrate that the electromagnetic wave in the dispersive medium attenuates very fast and the measured signals are very weak. The cavity detection and the tunnel wall inspection are numerically proven to be possible and effective. ground penetrating radar, Maxwell's equations, finite difference time domain, dispersive medium, cavity detection 1 Introduction Ground penetrating radar (GPR) is a geophysical tool using high frequency radio waves to probe inside a medium (Li 1984). It uses wide-band electromagnetic pulses to delineate the structure below a surface or through an opaque object. With the advancement of recent years, GPR has become an important tool in engineering and environmental geophysics; it is widely applied in highway quality inspection and bedrock structure detection, high building foundation appraisal, shallow soil or rock layers structure, underground water level detection and water pollution mapping, and many other aspects. The geological material is diverse in physical parameters, and the electromagnetic attenuation in a subsurface material is stronger than that in the air; therefore the wave propagation in the subsurface material is more complex than that in the air. It is effective to study the radar wave propagation in such a material using finite difference time domain (FDTD) for improving understanding and interpretation of GPR detection (Liu 2004, Liu and Sato 2003, 2005, Debroux 1996; Holliger and Bergmann 1998). As a powerful computational electromagnetic tool, FDTD has been applied in simulating wave propagation, scattering and antenna radiation related to GPR. In previous FDTD development and application, the dielectric constant of the medium was considered as frequency independent (Liu 2004, Liu and Sato 2003, 2005 Debroux 1996; Holliger and Bergmann 1998). As the physical parameters of the medium are frequency dependent, such as the usual rock and soil probed by GPR, the dispersive characteristic of the medium will affect the electromagnetic wave propagation. Luebbers et al introduced a recursive convolution (RC) scheme in FDTD which takes into account the dispersive medium. Later, this approach became the piecewise linear recursive convolution method (PLRC) (Luebbers and Hunsberger 1992, David et al1992, Uno 1998). Teixeira et al expanded perfected matched layer (PML) for lossless media to lossy and dispersive situations and used the PLRC technique, meanwhile a parallel computation code was developed (Teixeira et al1998). To improve the computational efficiency, Fan and Liu (Fan and Liu 1998, Liu and Fan 1999a, 1996) developed a pseudospectral time domain (PSFD) method to compute the electromagnetic wave in the dispersive medium. We developed a novel FDTD method to compute the electromagnetic wave propagation in the dispersive medium by following Sullivan's theory (Sullivan 2000). Several numerical examples for urban applications of GPR are presented in this paper and the simulated results are analysed. 2 FDTD formulation in the dispersive medium Physical parameters related to electromagnetic waves include electric conductivity σ, dielectric constant ε and magnetic permeability μ. Magnetic permeability reflects the relationship between magnetic field strength H and magnetic flux density B, and is always a constant for most earth materials. Electric conductivity shows hardness when current fluxes the medium, for example, a metal has very a large conductivity of one thousand to ten thousand or more. Another important parameter is the dielectric constant, which reflects the relationship between electric flux density D and electric field strength E. We have 1 In alternate electromagnetic fields, electric parameters of rocks are changed with frequency. If the frequency is lower, the polar water molecule can keep up with the change in electric field. As the frequency increases, the polar molecule cannot keep up with the external field because of the inertia of the polarization process. This results in a new current component which differs from the displacement current by 90° and coincides with the convection current. The polar lag produces appendix conductance. This process occurs with the decrease of ε because the polarization cannot reach its maximum. If we think the polarization reaches its maximum exponentially with time constant τ, the complex dielectric constant can be expressed as 2 where ε∞r and ε0r are dielectric constants when the frequency is infinite and 0, τ is the relaxation time. This relationship is called as the Debye formula (Daev 1981). The frequency-dependent electric flux density D is introduced in order to consider the dispersive medium. Maxwell's equations in the dispersive medium can be expressed as 3E and D are normalized as 4 so that material parameters are omitted in computation, and there are some other advantages for constructing the absorbing boundary condition. Then Maxwell's equations become 5 where and H can be expressed as the six scalar equations: 6a 6b 6c 6d 6e 6f After difference approximation, they can be expressed as 7a 7b Here only corresponding equations for (6c) and (6f) are introduced as examples, other equations can be stated similarly. Δt and Δx are time increment and grid size, respectively. Electric and magnetic field space configuration can be found in Uno (1998). The computation for is complicated and requires deduction from (5). When conductivity is included, the frequency-dependent dielectric constant can be expressed as 8 Electric flux density and electric field strength are related by 9 Taking the first term into the time domain is not a problem because it is a simple multiplication. In the second term, Fourier theory tells us that 1/jω in the frequency domain is the integration in the time domain. The inverse transform of the Debye term is the convolution of ((ε0r - ε∞r)/τ) e-(t/τ) and (Sullivan 2000). Equation (9) becomes 10 Then (10) is approximated as a summation in the sampled time domain: 11 where 12 13 superscript n refers to the value at the nth time step. When ε∞r equals ε0r, equation (11) is fitted for the non-dispersive medium. Substituting (12) and (13) into (11), 14 Therefore, can be calculated by iteration. The inhomogeneous medium needs to be considered only as computation. Equation (14) can be expressed as scalar equations. The initial values are all zero in computation. We follow Sullivan's method to expand the PML technique to the lossy dispersive medium. This algorithm is fitted for both dispersive and non-dispersive media. Sullivan's PML adds fictitious dielectric constants to the normalized Maxwell's equations; these fictitious constants have nothing to do with the real physical values. This method has two distinct advantages over the original Bererger's PML. First, since the split of the fields is avoided, much computer resource is saved. Second, the dielectric constant and the conductivity appear only in the electric flux calculation; therefore it is easily fitted for lossy and dispersive media. 3 Numerical examples We developed a 3D FDTD code for electromagnetic simulation fitted for both dispersive and non-dispersive media. Since we focus on propagation and scattering characteristics of radar waves, we use an infinitesimal dipole instead of a large antenna with size and geometry. The source is excited by a Gaussian pulse. First, we compare the wave propagations in dispersive and non-dispersive media. The physical parameters for the dispersive medium are ε∞r = 2, ε0r = 3.5, σ = 0.001 S m-1 for the dispersive medium and ε0r = ε∞r = 3.5, σ = 0.001 S m-1 for the corresponding non-dispersive medium. The grid size is 0.05 m, the time increment is 8.34 × 10-2 ns, and the antenna's central frequency is 200 MHz. There are 15–20 grids per wavelength. The transmitter and receiver are polarized vertically, and are located in a homogeneous medium. A common source mode is adapted here; the offsets are 1, 2, …, 6 m, respectively. The simulated signals for the two cases are illustrated in figure 1(a). The signals are normalized by the maximum amplitude of each trace, which are illustrated in figure 1(b). The main differences for two situations are wavelet waveforms and amplitudes. Signal wavelets become wider as the distance increases for the dispersive medium than for the non-dispersive medium; this mean dispersive medium lowers the resolution. The radar wave attenuates faster in the dispersive than in the non-dispersive medium. Figure 1 Open in new tabDownload slide Simulated signals (a) in a dispersive medium (solid line) and a non-dispersive medium (dashed line) by FDTD, and the maximum amplitude (b) for each trace. Figure 1 Open in new tabDownload slide Simulated signals (a) in a dispersive medium (solid line) and a non-dispersive medium (dashed line) by FDTD, and the maximum amplitude (b) for each trace. A horizontal layer model is considered in the second example as shown in figure 2(a). The physical parameters are ε∞r = 4, ε0r = 7, σ = 0.005 S m-1, τ = 6 ns for medium 2, while the physical parameters for medium 1 change from ε0r = ε∞r = 2, σ = 0.001 S m-1, for the non-dispersive medium to ε∞r = 2, ε0r = 4, σ = 0.001 S m-1, τ = 5 ns for the dispersive medium. The grid size is 0.05 m, the time increment is 8.34 × 10-2 ns, and the antenna's central frequency is 200 MHz. There are 15–20 grids per wavelength. The layer interface is 2.5 m below the air–ground interface. The common source mode is used here. The transmitter marked by ‘T’ and the receivers marked by ‘R’ are illustrated in figure 2(a), the offset between the transmitter and the first receiver is 1 m, so is the distance between receivers. The simulated signals are shown in figure 2(b) for the non-dispersive medium (solid line) and the dispersive medium (dashed line). Each signal is normalized by the maximum amplitude of the first trace. It is noticed that the event of the reflected signals shows a hyperbola. The amplitude of the signal becomes weaker as the distance increases. The reflected signals in the dispersive medium are much more attenuated than those in the non-dispersive medium. Figure 2 Open in new tabDownload slide Horizontal layered model (a) and simulated signals (b) as medium 1 is dispersive (dashed line) and non-dispersive (solid line). Figure 2 Open in new tabDownload slide Horizontal layered model (a) and simulated signals (b) as medium 1 is dispersive (dashed line) and non-dispersive (solid line). We find that the reflection over a shorter distance is smaller than a longer distance and for far offsets the amplitude of the direct wave is smaller than the reflection for the non-dispersive medium. The wave phenomenon is relatively complicated for GPR. Some energy travels above the air/earth interface and attenuates very fast, this energy produces the air wave signals at the receiving antennas. Some energy transmitted by antenna is coupled into the ground, and this part of the energy moves down as a conical beam. As this beam touches the underground boundary, it will be reflected. As the reflected signal reaches the receiving antenna, it induces the reflected signals. Therefore the received reflected signals have different travel paths including the path before and after the reflection. Each travel path has a different attenuation constant. Generally, a signal with a longer travel path has larger attenuation in the same medium. But it is possible that the later reflected signals are stronger than the early ones in certain special situations because of the directivity of the antenna, though the later signal travels a longer distance than the early one. This situation can also be noticed in measured data as shown in figure 5.4 of Conyers' book (Conyers 2004). This problem will be investigated intensively in the future. Figure 3 shows a simulation for cavity detection. A cavity below the surface of the pavement is a potential hazard to human beings, therefore this kind of detection is very important. The section of a typical cavity below the surface is illustrated in figure 3(a). The cavity is filled with air and resulted from leakage of underground water. The top layer of the pavement is asphalt with ε∞r = ε0r = 4, σ = 0.0001 S m-1, the layer below is made of pebbles with ε∞r = ε0r = 7, σ = 0.001 S m-1, below the pebble layer is soil with ε0r = 10, ε∞r = 7, τ = 1.0 ns, σ = 0.001 S m-1. The grid size is 0.05 m, the time increment is 1.18 × 10-1 ns, and the antenna's central frequency is 100 MHz. There are more than 20 grids per wavelength. The simulated radargram is shown in figure 3(b), and a common offset measurement is chosen with 1 m offset. Both the top and the bottom of the cavity result in clear hyperbolic diffractions between which we can find the reflection wave corresponding to the shape of the cavity. This simulation provides us insights into the cavity detection. Figure 3 Open in new tabDownload slide A cavity model (a) showing the cross-section of a cavity below the pavement, and simulated GPR signals (b). Figure 3 Open in new tabDownload slide A cavity model (a) showing the cross-section of a cavity below the pavement, and simulated GPR signals (b). Figure 4 shows a simulation for a tunnel coating quality inspection. Tunnels are often made into rock. The tunnel wall is generally several metres high and very long in the horizontal direction, and with an arced top. The antenna is touching the wall during measurement. The original rock surface after digging is rough and irregular; therefore a concrete coating with a wire network is made on the rough rock surface to reinforce the tunnel. Ground penetrating radar is often applied to detect the gap between the concrete (ε∞r = ε0r = 6, σ = 0.005 S m-1) and the original rock (ε∞r = ε0r = 5, σ = 0.0001 S m-1). The existence of the gap means there is a defect in the tunnel. The grid size is 0.0025 m, the time increment is 5.897 × 10-3 ns, and the antenna's central frequency is 900 MHz. There are more than 70 grids per wavelength. A section of the tunnel wall is shown in figure 4(a). Metal wires are embedded in the concrete. There are two gap zones between the rock and the concrete: one is air filled (left), and the other is filled with water (right) with parameters ε0r = 81, ε0r = 5.5, τ = 0.0109 ns, σ = 0.01 S m-1. The simulated result is shown in figure 4(b). The offset is 10 cm, and the moving step is 2 cm. Metal wires result in strong reflections which interfere with other reflection signals from the two gaps. The reflection from the two gaps can be seen around 4 ns with different colours. Below these reflections, there are many multiples below 4 ns. In order to evaluate the effects from metal bars, we also present the simulated result without metal wires as shown in figure 4(c). The signals become clearer than those in figure 4(b). We can find a multiple for the left gap at around 6 ns, two multiples for the right gap at around 6 and 9 ns. This simulation shows that the detection of the gap is possible even when the metal wire is present, but a complex analysis is necessary. Figure 4 Open in new tabDownload slide A tunnel wall coating model (a) shows the gaps between the original rock and the concrete, and the simulated GPR signals when the metal wires are present (b) or not (c). Figure 4 Open in new tabDownload slide A tunnel wall coating model (a) shows the gaps between the original rock and the concrete, and the simulated GPR signals when the metal wires are present (b) or not (c). The GPR penetration depth is determined by the antenna frequency and the electrical conductivity of the earth materials being profiled. It is inversely proportional to the antenna frequency as the other factors are fixed. It varies from one or two wavelengths in highly lossy earth to tens of wavelengths in dry materials. 4 Conclusions The frequency-dependent dielectric constant appears only in the process of electric fields computation by normalizing electric flux density D and electric field strength E. The relationship between electric flux density D and electric field strength E is determined by the Debye formula, the iteration algorithm for E in the dispersive medium is presented. Several examples related to urban applications are presented. It is concluded that FDTD is a very effective tool for GPR simulation and GPR can be applied to cavity detection and tunnel wall inspection. Acknowledgment This research is supported by NSFC with no 40474042. References Conyers L B . , 2004 , Ground-Penetrating Radar for Archaeology New York Rowman & Littlefield Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Daev D C . , 1981 , High-Frequency Electromagnetic Logging Method (translated by Gen X W, Zhao Y W, Lu D) Beijing Petroleum Industry Press (in Chinese) Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC David F , Kelley R , Luebbers J . , 1992 Piecewise linear recursive convolution for dispersive media using FDTD , IEEE Trans. Antennas Propag. , vol. 40 (pg. 792 - 7 ) OpenURL Placeholder Text WorldCat Debroux R S . , 1996 3D Modeling of the electromagnetic response of geophysical targets using the FDTD method , Geophys. Prospect. , vol. 44 (pg. 457 - 68 )http://dx.doi.org/10.1111/j.1365-2478.1996.tb00157.x0016802513652478 Google Scholar Crossref Search ADS WorldCat Fan G X , Liu Q H . , 1998 A 3D PML-FDTD algorithm for simulating ground-penetrating radar on dispersive earth media Proc. 7th Int. conf. Ground Penetrating Radar (pg. 579 - 84 ) Holliger K , Bergmann T . , 1998 Accurate and efficient FDTD modeling of ground-penetrating radar antenna radiation , Geophys. Res. Lett. , vol. 25 (pg. 3883 - 996 )http://dx.doi.org/10.1029/1998GL90004900948276 Google Scholar Crossref Search ADS WorldCat Li D X . , 1994 , Method and Application of GPR Beijing Geological Press (in Chinese) Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Liu S X . , 2004 Multi-frequency electromagnetic well logging simulation by 2-D cylindrical FDTD , J. Jilin Univ.: Earth Sci. Edition , vol. 34 (pg. 283 - 6 ) (in Chinese) OpenURL Placeholder Text WorldCat Liu S X , Sato M . , 2003 Simulation of borehole radar , J. Jilin Univ.: Earth Sci. Edition , vol. 33 (pg. 545 - 50 ) (in Chinese) OpenURL Placeholder Text WorldCat Liu S X , Sato M . , 2005 Transient radiation from an unloaded, finite dipole antenna in a borehole: experimental and numerical results , Geophysics , vol. 70 (pg. k43 - k51 )http://dx.doi.org/10.1190/1.21060481070485X Google Scholar Crossref Search ADS WorldCat Liu Q H , Fan G X . , 1999a A frequency-dependent PSTD algorithm for general dispersive media , IEEE Microw. Guid. Wave Lett. , vol. 9 (pg. 51 - 3 )http://dx.doi.org/10.1109/75.75504310518207 Google Scholar Crossref Search ADS WorldCat Liu Q H , Fan G X . , 1999b Simulation of GPR in dispersive media using a frequency-dependent PSTD algorithm , IEEE Trans. Geosci. Remote Sens. , vol. 37 (pg. 2317 - 24 )http://dx.doi.org/10.1109/36.78962801962892 Google Scholar Crossref Search ADS WorldCat Luebbers R J , Hunsberger F . , 1992 FDTD for Nth-order dispersive media , IEEE Trans. Antennas Propag. , vol. 40 (pg. 1297 - 301 )http://dx.doi.org/10.1109/8.2027070018926X Google Scholar Crossref Search ADS WorldCat Sullivan D M . , 2000 , Electromagnetic Simulation Using The FDTD Method New York IEEE Press Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Teixeira F L , et al. , 1998 Finite-difference time-domain simulation of ground penetrating radar on dispersive, inhomogeneous, and conductive soils , IEEE Trans. Geosci. Remote Sens. , vol. 36 (pg. 1928 - 37 )http://dx.doi.org/10.1109/36.72936401962892 Google Scholar Crossref Search ADS WorldCat Uno T . , 1998 , Finite Difference Time Domain Method for Electromagnetic Field and Antenna Analysis Japan Corona Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC © 2007 Nanjing Institute of Geophysical Prospecting
TI - FDTD simulations for ground penetrating radar in urban applications
JF - Journal of Geophysics and Engineering
DO - 10.1088/1742-2132/4/3/S04
DA - 2007-09-01
UR - https://www.deepdyve.com/lp/oxford-university-press/fdtd-simulations-for-ground-penetrating-radar-in-urban-applications-fmGp00itB0
SP - 262
VL - 4
IS - 3
DP - DeepDyve
ER -