CN115902875B - Ground penetrating radar forward modeling method and system based on LOD-FDTD - Google Patents

Ground penetrating radar forward modeling method and system based on LOD-FDTD Download PDF

Info

Publication number
CN115902875B
CN115902875B CN202211164835.1A CN202211164835A CN115902875B CN 115902875 B CN115902875 B CN 115902875B CN 202211164835 A CN202211164835 A CN 202211164835A CN 115902875 B CN115902875 B CN 115902875B
Authority
CN
China
Prior art keywords
fdtd
detection target
lod
target information
fine
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202211164835.1A
Other languages
Chinese (zh)
Other versions
CN115902875A (en
Inventor
谢国大
宋子衡
侯桂林
潘攀
黄志祥
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Anhui University
Original Assignee
Anhui University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Anhui University filed Critical Anhui University
Priority to CN202211164835.1A priority Critical patent/CN115902875B/en
Publication of CN115902875A publication Critical patent/CN115902875A/en
Application granted granted Critical
Publication of CN115902875B publication Critical patent/CN115902875B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources

Abstract

The invention relates to a ground penetrating radar forward modeling method and system based on LOD-FDTD. Processing the multistage Debye dispersion medium by adopting a Z transformation technology to obtain an iterative formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium; performing numerical calculation on the fine grid area in a fine grid form by adopting an LOD-FDTD method to obtain detection target information; performing numerical calculation on the coarse grid area in the form of coarse grids by adopting an FDTD method to obtain detection target information; determining the whole detection target information according to the two detection target information; and carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation to obtain the detection target information after the calculation domain is reduced. The invention can solve the problems of low calculation efficiency, high calculation complexity and low calculation accuracy.

Description

Ground penetrating radar forward modeling method and system based on LOD-FDTD
Technical Field
The invention relates to the field of forward modeling of ground penetrating radar, in particular to a forward modeling method and system of ground penetrating radar based on LOD-FDTD.
Background
Ground penetrating radar (Ground Penetrating Radar, GPR) is a non-destructive measurement technique in geophysical prospecting methods that uses electromagnetic waves to locate objects or interfaces buried in visually opaque objects or subsurface media. The GPR consists of an antenna, a control unit, a display unit, a storage unit, other accessory equipment and the like. The geological radar transmits high-frequency electromagnetic waves through a transmitting antenna, and receives the reflected electromagnetic waves through a receiving antenna. By means of the reflected electromagnetic wave data transmitted back by the receiving antenna, specific parameter characteristics are extracted, and the spatial positions and the distribution of different media can be deduced by combining certain data processing and analysis means. When the GPR method is applied to detection, the attribute of a detected target can be deduced and analyzed mainly by comparing the transmitted data, and the reference data for comparison is mainly obtained by two methods of in-situ sampling and computer simulation. The former is to rely on experiments to complete data acquisition, and the implementation is complex. The computer simulation can be used for conveniently obtaining the comparison data through numerical calculation by establishing related environment models including soil, rock, underground water distribution, abnormal targets with high dielectric properties and the like, and the method can save the consumption of manpower, material resources and financial resources in the test and system improvement process and accelerate the progress of theoretical and experimental research. This computer simulation method, which may be referred to as GPR forward modeling, is one of the main research contents of geophysics. In terms of GPR forward modeling, FDTD is a numerical calculation method that should be relatively widespread at present. This is because the FDTD method has the characteristics of simple concept, flexibility, strong versatility, and the like. Meanwhile, the detected targets and the background soil form a detection scene, and the detection scene has the characteristics of non-uniformity, electric characteristic dispersion characteristic, anisotropy and the like. Fortunately, the electromagnetic properties of these structures can be easily simulated precisely using the FDTD method based on Yee cell format, so FDTD has long been used in the study of GPR systems. Giannopoulos et al developed the FDTD method-based ground penetrating radar forward software "GPRMax"; in the literature of forward modeling and ground penetrating radar imaging in tunnel lining leakage detection and grouting evaluation, shaari carries out related research on the target polarity of a ground penetrating radar antenna; although the FDTD method has many advantages in GPR forward modeling, for an electromagnetic model containing a complex medium with a high dielectric constant and a dispersion characteristic and an electrically small structure, the FDTD method needs to use finer fine grids to simulate the propagation characteristic of electromagnetic waves in a region with severe information changes such as amplitude, phase and the like, and meanwhile, the numerical method is guaranteed to have lower numerical dispersion error and anisotropic error. The adoption of the fine mesh subdivision mode can greatly increase the memory occupation of the computer. In addition, the FDTD method adopts an explicit differential strategy, and the time step is limited by the grid size, namely, the CFL stability condition needs to be met, so that a GPR numerical simulation system adopting fine grid subdivision is poor in simulation efficiency, and the calculation advantage of the FDTD method is greatly influenced.
In order to solve the problem, the sub-grid technology adopting multi-scale grid subdivision is an effective means for solving the problem, and the adopted strategies are as follows: the area with severe electromagnetic field amplitude variation caused by medium attribute and non-interface non-uniformity adopts fine mesh subdivision, and other areas adopt coarse mesh subdivision, so that over-sampling of the whole calculation space is avoided, the calculation efficiency can be greatly improved, and the calculation memory can be reduced. The essence of the subgrid technique is the exchange of information in time-space of electromagnetic field components between coarse and fine grids. And because the cell sizes of the coarse and fine grids are different, each calculation region needs to meet the corresponding time stability condition. In order to obtain stable FDTD solutions, the calculation schemes of early subgrid techniques mainly include the following two types: 1) The whole calculation area comprises a coarse grid area and a fine grid area, and the numerical calculation is carried out by adopting a uniform time step which meets the stability condition of the fine grid. 2) The coarse grid and the fine grid areas are respectively subjected to numerical calculation by adopting time steps meeting respective CFL conditions. The first type of method can avoid time interpolation processing between thick and thin grids, only needs to perform spatial interpolation processing, and is easy to realize numerical values, but the whole calculation area adopts smaller time step, so that the calculation efficiency is low. Compared with the first method, the second method has higher calculation efficiency, but needs to perform double interpolation operation in time and space, so the numerical implementation is complex. In order to avoid the interpolation process in time, and at the same time, the whole calculation area can use a larger time step, some mixed sub-grid techniques are proposed, that is, weak condition stable and unconditional stable FDTD methods are adopted in the fine grid area, such as: the HIE method, the AID-FDTD method, the CN-FDTD method and the spatial filtering FDTD method are used for calculating the fine grid region, and because the methods have the characteristic of expandable stability conditions or unconditional stability, the fine grid region can adopt time steps meeting the CFL conditions of the coarse grid, so that the whole calculation can adopt larger time steps to carry out numerical calculation, the calculation efficiency is improved, the interpolation processing in time is avoided, and the numerical realization is easier. However, the prior art only shows the application of ADI-FDTD in the two-dimensional sub-grid technology, and the ADI-FDTD has low calculation efficiency. In the prior art, the mixed sub-grid technology adopts a CN-FDTD method, but the method also has the problem of higher calculation complexity. The hybrid subgrid technique described in the prior art adopts a spatial filtering-FDTD method, and although the spatial filtering FDTD method retains the explicit iteration characteristics of the conventional FDTD method, the calculation accuracy is greatly affected by the dielectric properties of the materials in the calculation region. In summary, the GPR forward modeling method based on the conventional FDTD method has the problems of low calculation efficiency, large calculation complexity and low calculation accuracy.
Disclosure of Invention
The invention aims to provide a ground penetrating radar forward modeling method and system based on LOD-FDTD, which can solve the problems of low calculation efficiency, high calculation complexity and low calculation accuracy.
In order to achieve the above object, the present invention provides the following solutions:
a ground penetrating radar forward modeling method based on LOD-FDTD comprises the following steps:
acquiring an iterative formula of electric flux density and magnetic field strength based on an LOD-FDTD method in a ground penetrating radar scene;
acquiring constitutive relation of electric flux density and electric field strength in a frequency domain when a target is a dispersive medium in a complex ground penetrating radar scene;
according to the iterative formula and the constitutive relation, a Z transformation technology is adopted to process the multistage Debye dispersion medium, so that an iterative formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium is obtained;
acquiring a calculation region in a ground penetrating radar scene, wherein the calculation region comprises a coarse grid region and a fine grid region;
performing numerical calculation on the fine grid area in a fine grid mode by adopting an LOD-FDTD method to obtain first detection target information;
performing numerical calculation on the coarse grid region in a coarse grid form by adopting an FDTD method to obtain second detection target information;
Determining detection target information according to the first detection target information and the second detection target information;
and carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation to obtain the detection target information after the calculation domain is reduced.
Optionally, the iterative formula for acquiring the electric flux density and the magnetic field strength based on the LOD-FDTD method in the ground penetrating radar scene specifically includes:
applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field intensity based on the LOD-FDTD method:
wherein I is 6 6×6 identity matrix, A 1 、A 2 And u can be expressed as:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T mu is the magnetic permeability coefficient, A 1 And A 2 Is an auxiliary variable, +.>Is a spatial differential operator in x, y and z directions, D x 、D y And D z The electric flux densities in the x, y and z directions, H x 、H y And H z The magnetic field strength in the x, y and z directions respectively, Δt is the mathematical variable of the differential equation, and n is the time step in the iterative process.
Optionally, the constitutive relation of the electric flux density and the electric field strength is: d (ω) =ε 0 ε r (ω)E(ω);
Wherein D is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r For the dielectric constant as a function of frequency, ω represents D (ω) =ε 0 ε r Constitutive relation in (ω) E (ω) frequency domain.
Optionally, the performing numerical calculation on the fine grid area in the form of fine grid by using an LOD-FDTD method to obtain first detection target information specifically includes:
interpolation is carried out on the interface fine grid area by a spatial interpolation method, and a fine grid area boundary and an electric field value in the fine grid area are obtained;
calculating the magnetic field value of the fine grid region by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in a multistage Debye dispersion medium;
correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
and determining first detection target information according to the corrected electric field value and the corrected magnetic field value.
A ground penetrating radar forward modeling system based on LOD-FDTD comprises:
the first iterative formula acquisition module is used for acquiring an iterative formula of electric flux density and magnetic field strength based on an LOD-FDTD method in a ground penetrating radar scene;
the constitutive relation acquisition module is used for acquiring constitutive relation of electric flux density and electric field strength in a frequency domain when a target is a dispersive medium in a complex ground penetrating radar scene;
according to the iterative formula and the constitutive relation, a Z transformation technology is adopted to process the multistage Debye dispersion medium, so that an iterative formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium is obtained;
The second iteration formula determining module is used for processing the multistage Debye dispersion medium by adopting a Z transformation technology according to the iteration formula and the constitutive relation to obtain an iteration formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium;
the computing area acquisition module is used for acquiring a computing area in the ground penetrating radar scene, wherein the computing area comprises a coarse grid area and a fine grid area;
the first detection target information determining module is used for carrying out numerical calculation on the fine grid area in a fine grid mode by adopting an LOD-FDTD method to obtain first detection target information;
the second detection target information determining module is used for carrying out numerical calculation on the coarse grid area in a coarse grid mode by adopting an FDTD method to obtain second detection target information;
the detection target information determining module is used for determining detection target information according to the first detection target information and the second detection target information;
and the detection target information determining module is used for carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation after the calculation domain is reduced to obtain the detection target information after the calculation domain is reduced.
Optionally, the first iterative formula obtaining module specifically includes:
The first iterative formula obtaining unit is used for applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field strength based on the LOD-FDTD method:
wherein I is 6 6×6 identity matrix, A 1 、A 2 And u can be expressed as:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T mu is the magnetic permeability coefficient, A 1 And A 2 Is an auxiliary variable, +.>Is a spatial differential operator in x, y and z directions, D x 、D y And D z The electric flux densities in the x, y and z directions, H x 、H y And H z The magnetic field strength in the x, y and z directions respectively, Δt is the mathematical variable of the differential equation, and n is the time step in the iterative process.
Optionally, the constitutive relation of the electric flux density and the electric field strength is: d (ω) =ε 0 ε r (ω)E(ω);
Wherein D is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r For the dielectric constant as a function of frequency, ω represents D (ω) =ε 0 ε r Constitutive relation in (ω) E (ω) frequency domain.
Optionally, the first detection target information determining module specifically includes:
the electric field value determining unit is used for interpolating the interface fine grid area by adopting a spatial interpolation method to obtain the boundary of the fine grid area and the electric field value in the fine grid area;
the magnetic field value determining unit is used for obtaining the magnetic field value of the fine grid area by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in the multistage Debye dispersion medium;
The correction unit is used for correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
and the first detection target information determining unit is used for determining the first detection target information according to the corrected electric field value and the corrected magnetic field value.
According to the specific embodiment provided by the invention, the invention discloses the following technical effects:
the invention provides a ground penetrating radar forward modeling method based on LOD-FDTD, which is characterized in that the LOD-FDTD method is adopted in a fine grid area to carry out numerical calculation, and the traditional FDTD is adopted in a coarse grid area to carry out calculation, so that the fine grid subdivision of any proportion of a local coarse grid area can be realized, the universality is higher, and the problems of low calculation efficiency and low calculation accuracy can be solved. In addition, a ZT-LOD-FDTD method based on Z-transform technology is proposed for processing a multistage Debye-type dispersive medium with respect to the dispersion properties of soil, a detection target, and the like. Compared with an ADE method adopted in the traditional technology, the Z change method does not need to reserve a plurality of values at the first several moments, the memory occupation is less, the numerical implementation is direct, and the calculation complexity is reduced. Meanwhile, for the constructed back-shadow soil dispersion region, the CPML technology developed based on the Z transformation technology is adopted to conduct calculation space phase truncation, so that the calculated amount can be greatly reduced, and the calculation efficiency is improved.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings that are needed in the embodiments will be briefly described below, and it is obvious that the drawings in the following description are only some embodiments of the present invention, and other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of a forward modeling method of a ground penetrating radar based on LOD-FDTD;
FIG. 2 is a schematic view of a three-dimensional coarse and fine grid interface;
FIG. 3 is a schematic illustration of interpolation of the electric field ex at the interface;
FIG. 4 is a computational flow diagram of a hybrid subgrid LOD-FDTD method;
FIG. 5 is a schematic diagram of a typical three-dimensional GPR system with a detection target being a dielectric sphere;
FIG. 6 is a schematic diagram of time domain waveforms of electric field Ez at receiving antenna calculated by three methods, respectively;
FIG. 7 is a graph of relative calculation errors for a hybrid subgrid LOD-FDTD method;
fig. 8 is a schematic diagram of time domain waveforms comparing receiving points at different transmit-receive antenna distances;
FIG. 9 is a schematic diagram of a typical three-dimensional GPR system for detecting multiple targets;
FIG. 10 is a schematic diagram of time domain waveforms of electric field Ez at a receiving antenna calculated by three methods, respectively;
FIG. 11 is a graph of relative calculation errors for the hybrid subgrid LOD-FDTD method of different scales;
FIG. 12 is a schematic diagram of a time snapshot of a multi-objective scene for Ez amplitude calculation using a hybrid subgrid LOD-FDTD method;
FIG. 13 is a schematic diagram of time domain waves of the electric field Ez at the receiving antenna calculated by three methods, respectively;
FIG. 14 is a graph of relative calculation errors for the hybrid subgrid LOD-FDTD method of different scales;
FIG. 15 is a block diagram of a forward modeling system for a ground penetrating radar based on LOD-FDTD.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and completely with reference to the accompanying drawings, in which it is apparent that the embodiments described are only some embodiments of the present invention, but not all embodiments. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to be within the scope of the invention.
The invention aims to provide a ground penetrating radar forward modeling method and system based on LOD-FDTD, which can solve the problems of low calculation efficiency, high calculation complexity and low calculation accuracy.
In order that the above-recited objects, features and advantages of the present invention will become more readily apparent, a more particular description of the invention will be rendered by reference to the appended drawings and appended detailed description.
The prior art hybrid HIE subgrid method does not address the analysis of dispersive media. The hybrid LOD-FDTD subgrid technique in the prior art is used to simulate a periodic metal nanoparticle array, but only a two-dimensional case of computational scheme is presented. GPR forward modeling method based on traditional FDTD method often needs to adopt finer uniform Yee grid to ensure accuracy of numerical simulation result. While a large number of calculation units and a small time step resulting from a fine grid will greatly reduce the numerical calculation advantage of the FDTD method. In order to solve the problem, the invention adopts a mixed LOD-FDTD subgrid technology to numerically simulate several classical GPR scenes. The fine grid area is calculated by adopting an unconditionally stable LOD-FDTD method, and the residual coarse grid area is calculated by adopting a traditional FDTD method, so that the problem of local area oversampling caused by adopting uniform fine grids in the whole calculation is avoided. Meanwhile, due to unconditional stable characteristics of the LOD-FDTD method, the whole calculation area can adopt time steps meeting the CFL of the coarse grid to carry out numerical calculation, and long-time interpolation processing is avoided. In addition, the subgrid technology provided by the invention can realize any proportion subdivision of the coarse and fine grids, and has higher universality. For a dispersion medium with multi-stage Debye type dispersion attribute, the FDTD method and the LOD-FDTD based on the Z change method are adopted for processing, compared with the ADE method, the calculation steps are fewer, and the numerical implementation is more direct. Finally, through numerical simulation of three classical GPR scenes and comparison with the traditional FDTD, the correctness and the effectiveness of the method provided by the invention are verified.
FIG. 1 is a flow chart of a forward modeling method of a ground penetrating radar based on LOD-FDTD. As shown in fig. 1, a forward modeling method of a ground penetrating radar based on LOD-FDTD includes:
step 101: and acquiring an iterative formula of the electric flux density and the magnetic field strength based on the LOD-FDTD method in the ground penetrating radar scene.
Step 102: and acquiring constitutive relation of electric flux density and electric field strength in a frequency domain when the target is a dispersive medium in a complex scene of the ground penetrating radar.
Step 103: and processing the multistage Debye dispersion medium by adopting a Z transformation technology according to the iteration formula and the constitutive relation to obtain an iteration formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium.
Step 104: a calculation region in a ground penetrating radar scene is acquired, wherein the calculation region comprises a coarse grid region and a fine grid region.
Step 105: and carrying out numerical calculation on the fine grid area in a fine grid form by adopting an LOD-FDTD method to obtain first detection target information.
The step 105 specifically includes:
interpolation is carried out on the interface fine grid area by a spatial interpolation method, and a fine grid area boundary and an electric field value in the fine grid area are obtained;
calculating the magnetic field value of the fine grid region by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in a multistage Debye dispersion medium;
Correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
and determining first detection target information according to the corrected electric field value and the corrected magnetic field value.
Step 106: and carrying out numerical calculation on the coarse grid region in the form of a coarse grid by adopting an FDTD method to obtain second detection target information.
Step 107: and determining the detection target information according to the first detection target information and the second detection target information.
Step 108: and carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation to obtain the detection target information after the calculation domain is reduced.
The step 101 specifically includes:
applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field intensity based on the LOD-FDTD method:
the calculation formula is written in the form of a matrix:
wherein:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T (3)
wherein I is 6 6×6 identity matrix, A 1 And A 2 Is an auxiliary variable, mu is the magnetic permeability coefficient,is a spatial differential operator in x, y and z directions, D x 、D y And D z The electric flux densities in the x, y and z directions, H x 、H y And H z The magnetic field strengths in the x, y and z directions, respectively, are obtained by differential discretization of the Max Wei Xuandu equation over time, wherein Δt has no specific meaning, but is commonly used in discretization equations N is the number of time steps in the iterative process.
Meanwhile, the constitutive relation of D and E is:
D(ω)=ε 0 ε r (ω)E(ω) (5)
d is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r ω represents the constitutive relation in the frequency domain of the formula for the frequency dependent dielectric coefficient. Equation (5) represents the processing of the constitutive relation in the frequency domain.
For a multi-stage Debye dispersive medium, the constitutive relationship in equation (5) can be expressed as:
wherein ε 0 Is the free space dielectric constant, σ is the conductivity, Δε p =ε s,p∞,p Is the difference between the relative permittivity at zero frequency and the relative permittivity at infinite frequency. τ p Is the pole relaxation time, N p Is the number of poles in the sensitivity response.
Transforming equation (6) from the frequency domain to the z-domain form:
wherein F is p R is an auxiliary variable:
converting (7), (8), (9) from z-domain to time domain, i.e. using P n Instead of P (z), P n-1/2 Instead of z -1/2 P(z)(P=(D,E,F p ,R)),
N+1/2 to n+1 sub-time steps:
equations 10-13 are iterative equations for electric field strength and magnetic field strength.
Wherein:
substituting (14) and (15) into (10) to obtain an iterative formula (taking the x direction as an example) of E at the time of n+1/2:
wherein:
also, E can be obtained by the method described above x Iterative formula at time n+1
Wherein:
The invention combines the traditional FDTD method and the LOD-FDTD method to provide a mixed multi-scale grid computing scheme. The electromagnetic components of the coarse grid area are calculated by adopting a traditional FDTD method, and the electromagnetic components of the fine grid area are calculated by adopting an LOD-FDTD method (such as formulas (17) and (18)). By virtue of the unconditional stability of the LOD-FDTD method, the fine grid region can adopt a time step which meets the CFL stability condition of the coarse grid, so that the coarse and fine grid regions only need to be subjected to spatial interpolation processing, and do not need to be subjected to temporal interpolation processing. The purpose of the spatial interpolation is that the fine-grid electric field value e at the interface of the coarse and fine grids cannot be solved iteratively by the LOD-FDTD method, such as the electric field value e of the fine-grid region at the interface in FIG. 2 x The coarse grid electric field E on the interface is needed to be interpolated to the fine grid electric field value E on the interface, so that the electromagnetic field component in the fine grid area can be solved iteratively. In the three-dimensional computation space, the coarse and fine grid regions share six interfaces. For descriptive convenience, the invention uses the electric field value e at the interface of the thick and thin grids shown in FIG. 2 x For example, the related interpolation is described, and the interpolation processing method on the rest interfaces is the same.
FIG. 3 shows interface j of FIG. 2 F =j A (j f Schematic diagram of electric field distribution of coarse and fine grids on =1), subscript F represents coarse grid coordinates, and F represents fine grid coordinates. To be used forFig. 3 in-frame fine-grid electric fieldTo the point ofFor example, wherein: s is(s) 1 Is a positive integer, s is more than or equal to 1 1 R is less than or equal to R, R is the proportion of coarse and fine meshes, i F Grid coordinates (i) representing the coarse grid at the interface in the x-direction A ≤i F <i B )。i f Grid coordinate i representing fine grid at interface in x-direction f =R×(i F -i A ),(i A ,i B ),(j A ,j B ),(k A ,k B ) The front and rear interfaces in the x, y, z directions of the coarse mesh region are shown, respectively. />
In-frame fine grid electric fieldTo->The method can be calculated by adopting the following spatial interpolation formula:
when (when)
When (when)
Wherein the method comprises the steps of
/>
It should be noted that: fig. 3 shows that the interpolation calculation formula is the same in the case that the coarse and fine mesh division ratio is odd and the coarse and fine mesh division ratio is even.
After the above spatial interpolation process, j is in the fine grid region f =2 field on planeThe iterative formula of (c) needs to be modified as follows, taking the n+1/2 time instant calculation formula as an example (formula (17)).
When (when)Time->
When (when)
Wherein the method comprises the steps ofObtained from the above equations (19) (20). s is(s) 2 Is a positive integer, s is more than or equal to 1 2 ≤R,k F Represents grid coordinates (k) of the coarse grid at the interface in the z-direction A ≤k F <k B )。k f Grid coordinate k representing the z-direction of a fine grid at an interface f =R×(k F -k A )。
After the iteration of the electromagnetic components in the fine grid area is completed, interpolation processing is needed to be carried out on the electromagnetic components of the coarse grid in the fine grid area, and the electromagnetic information exchange between the coarse grid and the fine grid is completed (taking the x direction as an example, other directions have the same interpolation correction formula):
if R is an odd number
/>
If R is an even number
The final mixed subgrid LOD-FDTD method comprises the following steps:
1) The electromagnetic field values (E, H) of the coarse mesh region are calculated using the FDTD method.
2) Interpolation is carried out on the fine mesh area of the interface by using the formulas (19) and (20) to obtain an electric field value e of the boundary of the fine mesh area.
3) Method (e) using equations (17), (18), (26), (27) y And e z The direction calculation formula and formulas (17) (18) (26) (27)) calculate the electric field e in the fine cell region.
4) The magnetic field value h of the fine mesh region is calculated using the formulas (11), (13).
5) The electromagnetic field values of the coarse mesh region are corrected using equations (28) - (31).
In order to test the correctness and effectiveness of the proposed hybrid LOD-FDTD subgrid method, three classical ground penetrating radar scenes are simulated by using the method, and a simulation platform of a numerical example of the invention is as follows: 32.0GB RAM, intel (R) Core (TM) i7-11700 CPU@3.80GHz. And cutting off the adopted numerical simulation area by utilizing a CPML technology constructed based on Z transformation, wherein an underground detection target comprises two electromagnetic structures of a medium ball and a medium cylinder. The transmitting antenna is used as an excitation source point in the free space, and the receiving antenna is used as a detection point for respectively transmitting and receiving detection signals. The time domain form of the excitation source is:
Wherein the center frequency is 300MHz, and the time period τ=1.7/f c Coefficient a 1 =-0.493,a 2 =0.144,a 3 =-0.01。
Scene one: medium ball
In practical situations, the density, humidity and other attributes of the soil can be different due to different depths. FIG. 5 is a schematic diagram of a typical three-dimensional GPR system with a detection target being a dielectric sphere; tx is the transmitting end, rx is the receiving end, and as shown in FIG. 5, the simulation area in this example contains three layers of media, namely air, soil 1 and soil 2. The detection target is a medium ball, and the upper part and the lower part of the medium ball are respectively positioned in two soil media. The dispersion attribute of the soil is described by adopting a multistage Debye model, and the soil parameter of the first layer is epsilon =3.0,Δε 1 =0.10,Δε 2 =0.25,τ 1 =2.78ns,τ 2 =0.113 ns, conductivity σ=0.415×10 -3 S/m. Second layerSoil parameters epsilon =4.5,Δε 1 =2.10,Δε 2 =0.70,τ 1 =4.08ns,τ 2 =0.261 ns, conductivity σ=1.11×10 -3 S/m. The length of the calculation space in the x, y and z directions is 1.8m,1.8m and 2.4m respectively, and the CPML layer number is 20. Where the coarse mesh size is Δx=Δy=Δz=0.03 m. The transmitting and receiving antennas are placed in free space, 0.18m from the ground, and the transmitting and receiving antennas (Tx-Rx) are 0.45m apart in the y-direction. The size of the subgrid area is 0.3m multiplied by 0.3m, the detection target is a dielectric sphere buried 0.3m below the ground, the relative dielectric constant is 30, and the radius is 0.12m. In order to verify the correctness and the high efficiency of the calculation result of the mixed subgrid LOD-FDTD method provided by the invention, a uniform fine grid FDTD method (delta x is adopted f =Δy f =Δz f =0.01m) as a reference solution.
Fig. 6 shows the results of the calculation using the hybrid LOD-FDTD subgrid method with two grid split ratios (split ratios 1:3 and 1:4, respectively) and the calculation using the uniform fine grid FDTD method. From the waveform consistency, the proposed subgrid method is applicable to odd and even subdivision. Meanwhile, in order to quantitatively analyze numerical calculation errors, the relative calculation errors of the mixed subgrid LOD-FDTD method are calculated by adopting a formula (33), as shown in fig. 7, the relative errors generated by odd or even times are smaller, and as the grid section ratio becomes larger, the numerical errors are smaller, so that the universality of the method is further verified, and the fact that the fine grid area adopts a larger time step and cannot cause larger influence on numerical dispersion difference is also shown.
Fig. 8 shows a received signal when the distance between the receiving antenna and the transmitting antenna is changed from 0 to 0.39m, and it can be found that the information of the electromagnetic radiation of the detection target received by the receiving end increases with the increase of the Tx-Rx distance by the normalized time domain waveform.
Table 1 shows the memory occupation and the CPU running time of the two methods, and it can be seen from table 1 that the proposed hybrid subgrid LOD-FDTD method can greatly reduce the computational memory and the CPU running time.
TABLE 1 memory occupancy (MB) and computation time(s) for fine-grid (FDTD) and hybrid subgrid methods
Scene II: multi-target detection including a media column and a media sphere
Considering the diversity of subsurface detection targets, in a second example, a GPR scene containing a plurality of detection targets is simulated. In the embodiment, the medium sphere and the medium column are placed in the same soil coal quality, and the correctness and the effectiveness of the proposed mixed subgrid LOD-FDTD method in the simulation practical application are further verified. Fig. 9 is a schematic diagram of a typical three-dimensional GPR system for detecting multiple targets, where Tx is a transmitting end and Rx is a receiving end, and as shown in fig. 9, the size of the whole computation space is the same as that of the computation example a, and the number of CPML layers is 10. The soil parameter is epsilon =3.2,Δε 1 =0.75,Δε 2 =0.3,τ 1 =2.71ns,τ 2 =0.108 ns, conductivity σ=1.11×10-3S/m. The two sub-grid areas are 0.3m x 0.3m, the two detected targets are buried under the ground by 0.3m, and the centers of the two detected targets are 0.15m apart in the x direction. Dielectric constant epsilon of dielectric sphere r =30, radius 0.12m; dielectric constant of dielectric cylinder epsilon r =50, radius 0.09m, length 0.24m. In order to verify the correctness and effectiveness of the hybrid subgrid LOD-FDTD method, the calculation result of the FDTD method of uniform fine mesh subdivision is used as a comparison.
Fig. 10 shows time domain waveforms at the receiving antenna calculated by the mixed LOD-FDTD method with different split ratios and the FDTD method with uniform and precise split, and the calculated results are better matched. Meanwhile, fig. 11 shows the calculation error of the mixed LOD-FDTD and the reference solution under the multi-target detection scene, and the correctness of the proposed method is verified through the results. Fig. 12 is a time snapshot diagram of a multi-target scene for calculating Ez amplitude using a hybrid subgrid LOD-FDTD method, where (a) is a snapshot 1 diagram, (b) is a snapshot 2 diagram, (c) is a snapshot 3 diagram, and (d) is a snapshot 4 diagram. Fig. 12 extracts the instantaneous changes of the electric field at different moments by means of electric field snapshot, and the transition from the coarse grid region to the fine grid region can be smoothed by the electric field, reflecting the stability and correctness of the proposed coarse-fine field value exchange of the hybrid subgrid LOD-FDTD.
Table 2 shows the memory occupation and CPU running time of the two methods, and it can be seen from table ii that the proposed hybrid subgrid LOD-FDTD method can greatly reduce the computational memory and CPU running time.
TABLE 2 memory occupancy (MB) and computation time(s) for fine-grid (FDTD) and hybrid subgrid methods
Scene III: dispersion ball
In the case of actual GPR detection, a detection target whose dielectric coefficient is frequency dependent is often encountered. The example simulates the GPR scene, changes the medium balls in the example A into dispersion medium balls, adopts the same calculation space, and only considers the condition of only one layer of soil. The dispersion medium sphere is constructed by a multi-stage Debye model, and the parameter is epsilon =168.0151,Δε 1 =34.2007,Δε 2 =12.0813,τ 1 =13.56ms,τ 2 =0.108ms。
Fig. 13 calculates the time domain waveform of the electric field of the receiving antenna by using a mixed subgrid LOD-FDTD subdivision scheme with a coarse/fine grid ratio of 1:3 and 1:4 and a uniform fine grid FDTD method. With the calculation result of the uniform fine mesh FDTD as a reference solution, fig. 14 quantitatively analyzes errors of the two split modes. As shown in the figure, the time domain waveform is better in fit and the relative error of the calculation result is small, so that the correctness of the proposed hybrid subgrid LOD-FDTD method is verified.
Table 3 shows the memory occupation and CPU running time of the three methods, and it can be seen from table iii that the proposed hybrid subgrid LOD-FDTD method can greatly reduce the computational memory and CPU running time.
TABLE 3 memory occupancy (MB) and computation time(s) for fine-grid (FDTD) and hybrid subgrid methods
The invention also provides a ground penetrating radar forward modeling system based on LOD-FDTD, and FIG. 14 is a structure diagram of the ground penetrating radar forward modeling system based on LOD-FDTD, as shown in FIG. 14, the system comprises:
The first iterative formula obtaining module 201 is configured to obtain an iterative formula of electric flux density and magnetic field strength based on the LOD-FDTD method in a ground penetrating radar scene;
the constitutive relation obtaining module 202 is configured to obtain constitutive relation of electric flux density and electric field strength in a frequency domain when the target is a dispersive medium in a complex scene of the ground penetrating radar. The constitutive relation of the electric flux density and the electric field intensity is as follows: d (ω) =ε 0 ε r (ω)E(ω);
Wherein D is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r ω represents the constitutive relation in the frequency domain of the formula for the frequency dependent dielectric coefficient.
The second iterative formula determining module 203 is configured to process the multi-stage Debye dispersion medium by using a Z-transform technique according to the iterative formula and the constitutive relation, so as to obtain an iterative formula of electric field strength and magnetic field strength of ZT-LOD-FDTD in the multi-stage Debye dispersion medium;
a calculation region acquisition module 204, configured to acquire a calculation region in a ground penetrating radar scene, where the calculation region includes a coarse grid region and a fine grid region;
the first detection target information determining module 205 is configured to perform numerical calculation on the fine grid area in a fine grid form by adopting an LOD-FDTD method, so as to obtain first detection target information;
The second detection target information determining module 206 is configured to perform numerical calculation on the coarse mesh area in the form of a coarse mesh by using an FDTD method to obtain second detection target information;
a detection target information determining module 207, configured to determine detection target information according to the first detection target information and the second detection target information;
and the detection target information determining module 208 after the calculation domain is reduced is configured to spatially truncate the detection target information by adopting a CPML method based on Z transformation, so as to obtain detection target information after the calculation domain is reduced.
The first iterative formula obtaining module 201 specifically includes:
the first iterative formula obtaining unit is used for applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field strength based on the LOD-FDTD method:
wherein I is 6 6×6 identity matrix, A 1 、A 2 And u can be expressed as:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T/>mu is the magnetic permeability coefficient, A 1 And A 2 Is an auxiliary variable, +.>Is a spatial differential operator in x, y and z directions, D x 、D y And D z The electric flux densities in the x, y and z directions, H z 、H y And H z The magnetic field strength in the x, y and z directions respectively, Δt is the mathematical variable of the differential equation, and n is the time step in the iterative process.
The first detection target information determining module 205 specifically includes:
the electric field value determining unit is used for interpolating the interface fine grid area by adopting a spatial interpolation method to obtain the boundary of the fine grid area and the electric field value in the fine grid area;
the magnetic field value determining unit is used for obtaining the magnetic field value of the fine grid area by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in the multistage Debye dispersion medium;
the correction unit is used for correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
and the first detection target information determining unit is used for determining the first detection target information according to the corrected electric field value and the corrected magnetic field value.
The explicit differential characteristic of the traditional FDTD algorithm makes the calculation accuracy and efficiency of the traditional FDTD algorithm difficult to balance in GPR scene simulation calculation. In order to ensure accurate numerical calculation, smaller space step length is generally adopted to conduct grid subdivision, so that the application advantage of the FDTD method in the GPR system is greatly reduced. Based on the method, the invention provides a mixed subgrid LOD-FDTD method, an unconditionally stable LOD-FDTD method is adopted in a fine grid area, an FDTD method is adopted in other coarse grid areas, and the multi-scale mixed subgrid technology can effectively simulate a three-dimensional GPR scene. In general, the main contributions of the work of the invention are as follows: 1) The ZT-LOD-FDTD method based on the Z transformation method simulates the multi-stage Debye dispersion coal quality, and compared with a cyclic convolution method and an ADE method, the numerical implementation is more direct; 2) The mixed subgrid technology combining the LOD-FDTD method and the traditional FDTD method is provided, the unconditional stability of the LOD-FDTD method enables the whole calculation area to adopt a uniform coarse grid time step to carry out numerical calculation, so that only spatial interpolation processing is needed, and no temporal interpolation processing is needed, and the application difficulty of the subgrid technology is reduced; 3) The proposed subgrid technology can carry out fine grid subdivision of any proportion on the coarse grid in the calculation area, and has higher universality. Finally, through simulating three GPR scenes, the correctness and the effectiveness of the hybrid LOD-FDTD subgrid technology are verified.
In the present specification, each embodiment is described in a progressive manner, and each embodiment is mainly described in a different point from other embodiments, and identical and similar parts between the embodiments are all enough to refer to each other. For the system disclosed in the embodiment, since it corresponds to the method disclosed in the embodiment, the description is relatively simple, and the relevant points refer to the description of the method section.
The principles and embodiments of the present invention have been described in detail with reference to specific examples, which are provided to facilitate understanding of the method and core ideas of the present invention; also, it is within the scope of the present invention to be modified by those of ordinary skill in the art in light of the present teachings. In view of the foregoing, this description should not be construed as limiting the invention.

Claims (6)

1. The ground penetrating radar forward modeling method based on LOD-FDTD is characterized by comprising the following steps of:
acquiring an iterative formula of electric flux density and magnetic field strength based on an LOD-FDTD method in a ground penetrating radar scene;
acquiring constitutive relation of electric flux density and electric field strength in a frequency domain when a target is a dispersive medium in a complex ground penetrating radar scene;
According to the iterative formula and the constitutive relation, a Z transformation technology is adopted to process the multistage Debye dispersion medium, so that an iterative formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium is obtained;
acquiring a calculation region in a ground penetrating radar scene, wherein the calculation region comprises a coarse grid region and a fine grid region;
performing numerical calculation on the fine grid area in a fine grid mode by adopting an LOD-FDTD method to obtain first detection target information;
performing numerical calculation on the coarse grid region in a coarse grid form by adopting an FDTD method to obtain second detection target information;
determining detection target information according to the first detection target information and the second detection target information;
carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation to obtain detection target information after a calculation domain is reduced;
the method for carrying out numerical calculation on the fine grid area in a fine grid form by adopting an LOD-FDTD method to obtain first detection target information specifically comprises the following steps:
1) Interpolation is carried out on the interface fine grid area by a spatial interpolation method, and a fine grid area boundary and an electric field value in the fine grid area are obtained;
Interpolation is carried out on the fine grid area of the interface by using formulas (19) and (20), and an electric field value of the boundary of the fine grid area is obtained;
equation 19 is:
when (when)
Equation 20 is:
when (when)
Wherein the method comprises the steps of
Calculating the electric field in the fine cell region using the methods of formulas (17), (18), (26), (27); equation (17) is:
wherein:
equation (18) is:
equation (26) is:
when (when)Time of day
Equation (27) is:
when (when)
Wherein the method comprises the steps ofObtained from the above equations (19) (20); s is(s) 2 Is a positive integer, s is more than or equal to 1 2 ≤R,k F Represents grid coordinates, k, of the coarse grid at the interface in the z direction A ≤k F <k B ;k f Represents grid coordinates, k, of the fine grid at the interface in the z direction f =R×(k F -k A ) The method comprises the steps of carrying out a first treatment on the surface of the Subscript F represents coarse grid coordinates and F represents fine grid coordinates;
2) Calculating the magnetic field value of the fine grid region by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in a multistage Debye dispersion medium;
calculating a magnetic field value h of the fine mesh region using the formulas (11), (13);
the formula (11) is:
equation (13) is:
3) Correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
correcting the electromagnetic field value of the coarse mesh region using equations (28) - (31);
equation (28) is:
if R is an odd number of times,
equation (29) is:
if R is an odd number of times,
Equation (30) is:
if R is an even number and if R is an even number,
equation (31) is:
if R is an even number
Wherein s is 1 Is a positive integer, s is more than or equal to 1 1 R is not more than R, R is the proportion of coarse and fine grids; i.e F Grid coordinates of the coarse grid at the interface in the x direction are represented; i.e f Grid coordinate i representing fine grid at interface in x-direction f =R×(i F -i A ) The method comprises the steps of carrying out a first treatment on the surface of the Ex represents the electric field strength in the x direction; e, e x An electric field value representing a fine-meshed region at the interface; (i) A ,i B ),(j A ,j B ),(k A ,k B ) Front and rear interfaces in x, y and z directions of the coarse grid region are respectively represented; epsilon 0 Is the free space dielectric constant; sigma is the conductivity; τ p Is the pole relaxation time; delta epsilon p =ε s,p∞,p Is the difference between the relative permittivity of zero frequency and the relative permittivity at infinite frequency; mu is the magnetic permeability coefficient;
4) And determining first detection target information according to the corrected electric field value and the corrected magnetic field value.
2. The forward modeling method of ground penetrating radar based on LOD-FDTD according to claim 1, wherein the iterative formula for obtaining the electric flux density and the magnetic field strength based on LOD-FDTD in the ground penetrating radar scene specifically comprises:
applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field intensity based on the LOD-FDTD method:
Wherein I is 6 6×6 identity matrix, A 1 、A 2 And u can be expressed as:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T mu is the magnetic permeability coefficient, A 1 And A 2 Is an auxiliary variable, +.>Is a spatial differential operator in x, y and z directions, D x 、D y And D z In the x, y and z directions respectivelyElectric flux density of H x 、H y And H z The magnetic field strength in the x, y and z directions respectively, Δt is the mathematical variable of the differential equation, and n is the time step in the iterative process.
3. The forward modeling method of ground penetrating radar based on LOD-FDTD according to claim 1, wherein constitutive relation of electric flux density and electric field strength is: d (ω) =ε 0 ε r (ω)E(ω);
Wherein D is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r For the dielectric constant as a function of frequency, ω represents D (ω) =ε 0 ε r Constitutive relation in (ω) E (ω) frequency domain.
4. Ground penetrating radar forward modeling system based on LOD-FDTD, characterized by comprising:
the first iterative formula acquisition module is used for acquiring an iterative formula of electric flux density and magnetic field strength based on an LOD-FDTD method in a ground penetrating radar scene;
the constitutive relation acquisition module is used for acquiring constitutive relation of electric flux density and electric field strength in a frequency domain when a target is a dispersive medium in a complex ground penetrating radar scene;
The second iteration formula determining module is used for processing the multistage Debye dispersion medium by adopting a Z transformation technology according to the iteration formula and the constitutive relation to obtain an iteration formula of the electric field intensity and the magnetic field intensity of ZT-LOD-FDTD in the multistage Debye dispersion medium;
the computing area acquisition module is used for acquiring a computing area in the ground penetrating radar scene, wherein the computing area comprises a coarse grid area and a fine grid area;
the first detection target information determining module is used for carrying out numerical calculation on the fine grid area in a fine grid mode by adopting an LOD-FDTD method to obtain first detection target information;
the second detection target information determining module is used for carrying out numerical calculation on the coarse grid area in a coarse grid mode by adopting an FDTD method to obtain second detection target information;
the detection target information determining module is used for determining detection target information according to the first detection target information and the second detection target information;
the detection target information determining module is used for carrying out space truncation on the detection target information by adopting a CPML method based on Z transformation after the calculation domain is reduced to obtain the detection target information after the calculation domain is reduced;
the first detection target information determining module specifically includes:
1) The electric field value determining unit is used for interpolating the interface fine grid area by adopting a spatial interpolation method to obtain the boundary of the fine grid area and the electric field value in the fine grid area;
interpolation is carried out on the fine grid area of the interface by using formulas (19) and (20), and an electric field value of the boundary of the fine grid area is obtained;
equation (19) is:
when (when)
Equation (20) is:
when (when)
Wherein the method comprises the steps of
Calculating the electric field in the fine cell region using the methods of formulas (17), (18), (26), (27);
equation (17) is:
wherein:
equation (18) is:
wherein:
equation (26) is:
when (when)Time of day
Equation (27) is:
when (when)
Wherein the method comprises the steps ofObtained from the above equations (19) (20); s is(s) 2 Is a positive integer, s is more than or equal to 1 2 ≤R,k F Represents grid coordinates (k) of the coarse grid at the interface in the z-direction A ≤k F <k B );k f Grid coordinate k representing the z-direction of a fine grid at an interface f =R×(k F -k A );
2) The magnetic field value determining unit is used for obtaining the magnetic field value of the fine grid area by adopting an iterative formula of the electric flux density and the magnetic field strength of ZT-LOD-FDTD in the multistage Debye dispersion medium;
calculating a magnetic field value h of the fine mesh region using the formulas (11), (13);
the formula (11) is:
equation (13) is:
3) The correction unit is used for correcting the electric field value and the magnetic field value to obtain a corrected electric field value and a corrected magnetic field value;
Correcting the electromagnetic field value of the coarse mesh region using equations (28) - (31);
equation (28) is:
if R is an odd number of times,
equation (29) is:
if R is an odd number of times,
equation (30) is:
if R is an even number and if R is an even number,
equation (31) is:
if R is an even number
Wherein s is 1 Is a positive integer, s is more than or equal to 1 1 R is not more than R, R is the proportion of coarse and fine grids; i.e F Grid coordinates of the coarse grid at the interface in the x direction are represented; i.e f Grid coordinate i representing fine grid at interface in x-direction f =R×(i F -i A ) The method comprises the steps of carrying out a first treatment on the surface of the Ex represents the electric field strength in the x direction; e, e x An electric field value representing a fine-meshed region at the interface; (i) A ,i B ),(j A ,j B ),(k A ,k B ) Front and rear interfaces in x, y and z directions of the coarse grid region are respectively represented; epsilon 0 Is the free space dielectric constant; sigma is the conductivity; τ p Is the pole relaxation time; delta epsilon p =ε s,p∞,p Is the difference between the relative permittivity of zero frequency and the relative permittivity at infinite frequency; mu is the magnetic permeability coefficient;
4) And the first detection target information determining unit is used for determining the first detection target information according to the corrected electric field value and the corrected magnetic field value.
5. The forward modeling system of ground penetrating radar based on LOD-FDTD of claim 4, wherein said first iterative formula obtaining module specifically comprises:
The first iterative formula obtaining unit is used for applying the local one-dimension based on the split idea to the traditional Maxwell equation set to obtain an iterative formula of the electric flux density and the magnetic field strength based on the LOD-FDTD method:
wherein I is 6 6×6 identity matrix, A 1 、A 2 And u can be expressed as:
u=[D x ,D y ,D z ,H x ,H y ,H z ] T mu is the magnetic permeability coefficient, A 1 And A 2 Is an auxiliary variable, +.>Is a spatial differential operator in x, y and z directions, D x 、D y And D z The electric flux densities in the x, y and z directions, H x 、H y And H z The magnetic field strength in the x, y and z directions respectively, Δt is the mathematical variable of the differential equation, and n is the time step in the iterative process.
6. The LOD-FDTD based ground penetrating radar forward modeling system of claim 4, wherein the constitutive relation of electric flux density and electric field strength is: d (ω) =ε 0 ε r (ω)E(ω);
Wherein D is the electric flux density, E is the electric field strength, ε 0 Is the dielectric constant of free space epsilon r For the dielectric constant as a function of frequency, ω represents D (ω) =ε 0 ε r Constitutive relation in (ω) E (ω) frequency domain.
CN202211164835.1A 2022-09-23 2022-09-23 Ground penetrating radar forward modeling method and system based on LOD-FDTD Active CN115902875B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211164835.1A CN115902875B (en) 2022-09-23 2022-09-23 Ground penetrating radar forward modeling method and system based on LOD-FDTD

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211164835.1A CN115902875B (en) 2022-09-23 2022-09-23 Ground penetrating radar forward modeling method and system based on LOD-FDTD

Publications (2)

Publication Number Publication Date
CN115902875A CN115902875A (en) 2023-04-04
CN115902875B true CN115902875B (en) 2023-09-22

Family

ID=86492329

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211164835.1A Active CN115902875B (en) 2022-09-23 2022-09-23 Ground penetrating radar forward modeling method and system based on LOD-FDTD

Country Status (1)

Country Link
CN (1) CN115902875B (en)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104408256A (en) * 2014-12-01 2015-03-11 天津工业大学 Implementation algorithm for truncating one dimensional Debye medium Crank-Nicolson perfectly matched layer
CN110133644A (en) * 2019-06-03 2019-08-16 中南大学 Ground Penetrating Radar D integral pin-fin tube method based on Interpolating scaling functions method

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104408256A (en) * 2014-12-01 2015-03-11 天津工业大学 Implementation algorithm for truncating one dimensional Debye medium Crank-Nicolson perfectly matched layer
CN110133644A (en) * 2019-06-03 2019-08-16 中南大学 Ground Penetrating Radar D integral pin-fin tube method based on Interpolating scaling functions method

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
FDTD及其混合算法在微带天线计算中的理论与应用研究;李磊;《中国优秀硕博士(博士)学位论文全文数据库 信息科技辑》;正文第11-107页 *
Modeling Plasmonic Structures Using LOD-FDTD Methods With Accurate Dispersion Models of Metals at Optical Wavelengths;Konstantinos P. Prokopidis等;《 Journal of Lightwave Technology》;第193-200页 *
Three-Dimensional Unconditionally Stable LOD-FDTD Methods With Low Numerical Dispersion in the Desired Directions;Alok Kumar Saxena等;《IEEE Transactions on Antennas and Propagation》;第3055-3067页 *
基于LOD-FDTD的微带线边缘奇异性处理技术研究;李磊;张昕;孙亚秀;;电子与信息学报(03);D第746-752页 *
基于UPML边界条件的LOD-FDTD方法GPR波的数值模拟;韩晓冰;魏海亮;李奇;;西安科技大学学报(05);第703-708页 *

Also Published As

Publication number Publication date
CN115902875A (en) 2023-04-04

Similar Documents

Publication Publication Date Title
Warren et al. Creating finite-difference time-domain models of commercial ground-penetrating radar antennas using Taguchi’s optimization method
CN112949134B (en) Earth-well transient electromagnetic inversion method based on non-structural finite element method
Chaturvedi et al. Electromagnetic imaging of underground targets using constrained optimization
Lambot et al. Electromagnetic inversion of GPR signals and subsequent hydrodynamic inversion to estimate effective vadose zone hydraulic properties
Cassidy et al. The application of finite-difference time-domain modelling for the assessment of GPR in magnetically lossy materials
CN104992029B (en) DISCRETE RANDOM MEDIUM modeling method in a kind of multiple dimensioned non-homogeneous lunar soil layer
CN111239827B (en) Three-dimensional seismic data multiple suppression method based on local similarity coefficient
US20110155389A1 (en) System and Method For Providing A Physical Property Model
CN104408021A (en) Electric dipole source three-dimensional time domain finite difference direct interpretation imaging method
CN103777248A (en) TEM one-dimensional forward modeling method applicable to irregular transmitting loop
Mozaffari et al. 2.5 D crosshole GPR full-waveform inversion with synthetic and measured data
CN110414182B (en) Ground penetrating radar FRTM algorithm introducing antenna directional diagram
CN115238550A (en) Self-adaptive unstructured grid landslide rainfall geoelectric field numerical simulation calculation method
CN110119586B (en) Axial conductivity anisotropy transient electromagnetic three-component three-dimensional FDTD forward modeling method
CN104267443B (en) Magnetotelluric field static displacement correction method based on inversion model
Song et al. Finite element method for modeling 3D resistivity sounding on anisotropic geoelectric media
CN112733364B (en) Foil cloud scattering rapid calculation method based on impedance matrix partitioning
CN112327374B (en) DGTD forward modeling method for GPU ground penetrating radar complex medium
CN115902875B (en) Ground penetrating radar forward modeling method and system based on LOD-FDTD
CN117406190A (en) Non-excavation pull rod corrosion detection method, device and equipment based on radar signals
Qi et al. Full waveform modeling of transient electromagnetic response based on temporal interpolation and convolution method
CN109655910A (en) The two-parameter full waveform inversion method of Ground Penetrating Radar based on phasing
CN112882022B (en) Time-frequency combined full waveform inversion method for ground penetrating radar data
CN115542315B (en) ADI-FDTD-based ground penetrating radar forward modeling method and system
Stephan et al. Adding realistic noise models to synthetic ground‐penetrating radar data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant