Modeling effectiveness and identification of multi-scale objects in farmland soils with improved Yee-FDTD methods

Finite-Difference Time-Domain (FDTD) is the most popular time-domain approach in computational electromagnetics. Due to the Courant-Friedrich-Levy (CFL) condition and the perfect match layer (PML) boundary precision, FDTD cannot simulate soil medium whose surface is connected by multiple straight lines or curves (multi-scale) accurately and efficiently, which greatly limits the application of FDTD method to simulate buried objects in soils. Firstly, this study proposed the absorption boundary and adopted two typical perfect matching layers (UPML and CPML) to compare their absorption effects, and then using the three forms of improved Yee-FDTD algorithm, alternating-direction implicit (ADI-FDTD), unconditionally stable (US-FDTD) and hybrid implicit explicit finite time-domain (HIE-FDTD) to divide and contrast the boundary model effects. It showed that the HIE-FDTD was suitable for inversion of multi-scale structure object modeling, while ADI-FDTD and US-FDTD were ideal for single-boundary objects in both uniaxial perfectly matched layer (UMPL) and convolution perfectly matched layer (CPML) finite element space. After that, all the models were tested by computer performance for their simulated efficiency. When simulating single boundary objects, UPML-US-FDTD and ADI-FDTD could achieve the ideal results, and in the boundary inversion of multi-scale objects, HIE-FDTD modeling results and efficiency were the best. Test modeling speeds of CPML-HIE-FDTD were compared with three kinds of waveform sources, Ricker, Blackman-Harris and Gaussian. Finally, under the computer condition in which the CPU was i5-8250, the HIE-FDTD model still had better performance than the traditional Yee-FDTD forward modeling algorithm. For modeling multi-scale objects in farmland soils, the methods used CPML combined with the HIE-FDTD were the most efficient and accurate ways. This study can solve the problem that the traditional FDTD algorithm cannot construct non-mesh objects by utilizing the diversity characteristics of Yee cell elements.


Introduction 
The Yee-FDTD-based GPR forward modeling algorithm has become a hot research direction in geophysics and computer science gradually. Yee [1] applied a spatial meshing method based on Maxwell's equation and the set of finite difference equations was applicable for those boundary conditions or electromagnetic pulse. This method has been pioneered the FDTD method to simulate the electromagnetic field. Panagiotis et al. [2] proposed a GPR simulation system of three-dimensional based on FDTD models and consisted of a parabolic reflector transmitter and a multi-static receiver array. Their results showed that the approach had made the simulated model powerful and more accurate [3] . However, due to a large amount of calculation when simulating the antenna and the 3D module, the portability was not high. Lyu et al. [4,5] accorded to derivations of formulas and FDTD algorithm to establish a short axis of hyperbolic curves, and the radius of buried objects in the hidden fractures can be accurately detected. Ramadan [6] presented a new unsplit-field PML formulation for truncating frequency dispersive media based on the ADI-FDTD method. Garcia et al. [7] found that it can use the ADI-FDTD method to improve the Crank-Nicolson (CN) and Yee-FDTD schemes and simulate a simple parallel-plate structure excited by a low-frequency voltage source. Arimal et al. [8] used a multi-angle synthetic aperture (SA-GPR) method to obtain clear FDTD forward images, and proposed this effective modeling method but only for targets with a certain shape or cell size. Niu et al. [9] proposed a 3-D artificial anisotropy hybrid implicit-explicit finite-difference time-domain (AA-HIE-FDTD) with the nearly unchanged CFL condition and has higher accuracy. Iraklis Giannakis et al. [10,11] were fitting antennas model by FDTD and tested through a comparison between numerical and laboratory measurements. The experimental results showed that background features can affect the GPR frequencies performance of landmine detection but it must be validated using diverse scenarios and did not provide a reliable test outcome and advanced signal processing FDTD methods. In view of the shortcomings of the current FDTD simulation research, this paper will compare the performance of several improvement methods, US-FDTD, HIE-FDTD, and ADI-FDTD, and analyze those electromagnetic characteristics of objects in farmland soil.

Yee-FDTD's cell and TM-z polarization
The Yee-FDTD (FDTD) can simulate rich time-domain information of electromagnetic field problem directly, and simplify the medium physical process into a differential format electromagnetic parameter [12,13] . The relationship between the field quantities obeys the six rotation equations of the Maxwell equation.
After decomposing the current density into the conductor current density and applied current density, Maxwell's form can be denoted as: where, E is the electric field intensity; H is the magnetic-field intensity; σ is the conductivity of the medium, S/m; σ e denoted as the conductivity in the electric field and σ m denoted as conductivity in a magnetic field; t denoted as time, s; ε represents the medium dielectric constant, F/m; μ is the magnetic permeability; J is the current density, A/m 2 . M is the magnetic flux density;  is the vector differential operator. In this way, the geometric space problem can be discretized into orthogonal spatial grid points, and the excitation source can classify the electric and magnetic fields in discrete positions. FDTD relies on Maxwell's theoretical basis to decompose discrete 3D geometric problems in the medium into Yee units (Figure 1). i, j and k represent a certain three-dimensional coordinate position in Figure 1. Yee's cell constitutes an interlaced grid and uses the cell size as the minimum model resolution [14,15] . shows that the electric field component place in the middle of each edge of the Yee unit, and parallel to each edge; the magnetic field component is placed at the center of each face of the Yee unit, and parallel to the normal of each face. Each magnetic field vector was surrounded by four electric fields to form the curvature of the magnetic field. Similarly, each electric field vector is surrounded by four magnetic field vectors. Due to the properties of Faraday and Ampere's law, the field strength only displays the response amplitudes in the E z , H x , and H y vector directions, this can be denoted as the TM-z polarization mode [16,17] .

Figure 2 TMz model with magnetic and electric field components
If there was no numerical change in the direction of a certain coordinate (for example, the model constructed in Section 2.2 has a fixed value on the Z-axis in the coordinate system), the partial derivative related to this coordinate must exist [18][19][20] . Obviously, the electromagnetic classification in the two-dimensional case can be divided into two independent groups. E x , E y , and H z belong to a group which can be denoted as TE electromagnetic field group (Equation (3)). Similarly, E z , H x , and H y can be denoted as TM electromagnetic field group (Equation (4)). The two-dimensional FDTD equation can be simplified as: where, σ represents the medium's conductivity, S/m; ε represents the dielectric constant of the medium; represents the permeability of the medium, F/m; J is the current density, A/m 2 ; J z denoted as the current density along z-direction.

Selection of excitation source and modeling parameters
The gprMax is an open-source code and control platform that simulates electromagnetic wave propagation based on python code [21,22] . This experiment used this platform to for Yee-FDTD modeling. The depth of detected objects which can retard agricultural machinery operations is generally within 2.20 m, so the antenna model takes the bowtie shape of MALA 1.2G as a spatial solution Reference [23][24][25]. All GPR center frequency was set as 5 GHz. The time step preset as 3e-12 and the time-window preset as 5e−9. As for the source in Figure 3, there were two different waveforms' time domain and the change of the fast Fourier-transform. Obviously, the Ricker wave was a zero-phase wavelet, and the waveform had only one positive peak, which had a short duration and fast convergence [26,27] . There exist two distinct side lobes among the waveform, and this feature makes it easier for Yee-FDTD to be more accurate when dividing a finite grid. GPR detecting targets by means of ground coupling, a horizontal direct wave was formed which was reflected in the signal grayscale and produced a strong energy reflection of a direct wave in the early period or the soil surface. The Gaussian-dot waveform had multiple peaks of amplitude and troughs during its spread process. Those peaks and troughs values will increase within a certain range while the spread time is increasing (Figure 3a). The experiment set the excitation source type for different demand scenarios to ensure that the characteristics of these sources will perform the best analysis in the FDTD simulation boundary.
The total size of the soil layer 3D model set to 240 mm× 210 mm×120 mm, and the 2D model set to 240 mm×210 mm×2 mm which represents directions X, Y, and Z, respectively. Obviously, the smaller the step size was, the higher the model accuracy [28,29] . The time window must be the minimum to satisfy the electromagnetic wave that can pass through the target body and received in the simulated soil layer. It means that the ratio of spatial resolution should be 1/10 of the smallest wavelength because the wavelength and minimum spatial resolution can be calculated according to the following equation:

UPML and CPML boundary condition
The virtual absorbing boundary can simulate the non-physical reflection processing need to be constructed in the truncation boundary of the calculation domain. Traditionally, when we use the PML method to load the absorption boundary, it is set a special dielectric layer with a finite thickness at the cutoff boundary [30,31] .
The wave impedance of the dielectric layer is matched with the impedance of the adjacent dielectric wave completely, and the incident wave can be rapidly attenuated after entering the PML layer [32] . However, this PML method has its limitations. The target body of the simulation comparison must reach high physical mechanism accuracy. It cannot model or simulate Yee's cell of grids with a large incident angle and will difficult to control the low-frequency electromagnetic wave absorption effect.
The UPML algorithm can select a suitable uniaxial anisotropic medium constitutive parameter matching layer [33,34] . It avoids effectively the splitting of the electric and the magnetic field, and this also can be used in a Non-orthogonal and unstructured grid [35] . According to the principle of Ampere's and Faraday's laws to satisfy Maxwell's function, it made the UPML algorithm more efficient than PML when there is no need to take special care of the interface plane between the boundary and the interior regions [36] . Also, the uniaxial media can be matched perfectly in the internal loss media without any modification. The starting point of UPML theory is based on the modified form of Maxwell's lossless media equations. For the existence of the model's coordinate stretched, the Maxwell equation of the coordinate expansion medium with equivalent meaning can be expressed as follows: where, S x , S y , S z are coordinate scaling factor; x, y, z are the coordinates corresponding to the component forward coordinate system.
PML and UPML for telescopic coordinates are equivalent, So Equations (7) and (8) can also be expressed as: where, 1 1 1 and  denoted as diagonal dielectric constant;  denoted as magnetic permeability and has the uniaxial anisotropic medium characteristics. In summary, after the FDTD elements equation connection to the boundary condition, the function can be updated as: ( 1)/ 2 , 1 ( where, point (i, j) represents the projected coordinates at any point in space. The PML is a special medium composed of fictitious constitutive eigenvalue, which realizes electromagnetic wave non-reflective propagation through matching conditions, independent of frequency and incident angle. UPML proposes the use of anisotropic materials without changing the form of the electromagnetic field. And the CPML has high efficiency in absorbing litter waves and long-term time-domain calculations, but the storage space is very expensive. In the case that the appropriate parameter distribution can ensure the calculation accuracy, the calculation space can be effectively reduced to achieve the precision effect. The CPML frequency-domain equation in the coordinate expansion can be expressed as Equations (13)-(15): 11 11 where, S ew = k ew +σ pew /(α ew + jωε 0 ), S mw = k mw + σ pmw /(α ew + jωμ 0 ); and (w = x, y, z) is the coordinate scaling factor. The k ew, σ pew, α ew, κ mw, σ pmw and α mw denoted as the newly introduced parameters. The σ w and ε w denoted as conductivity and dielectric constant of the background medium. The σ pew, σ pmw are denoted as conductivity and magnetic permeability artificially that was added for the CPML region. The TE-z polarized plane wave (field's value is E x , E y and H z ) within any direction in two dimensions can be expressed as: In the air medium (σ =0) of the CPML region, the absorption efficiency of CPML is related to frequency. λ=ωlm was used to calculate the attenuation coefficient of the electromagnetic wave in the X-direction. When θ = 0, the attenuation coefficient can be denoted as: The larger of the σ p is, the stronger the absorption effect is. In the Yee differential meshing, the electric field and the magnetic field are spatially different from each other by half a grid, and the parameters of the proportional scaling factors S eia and S mi were also different by half a grid.

Multi-scale structuring
Each object composed of Yee cells has two planes that were identical and perpendicular to the coordinate system. The "Slack Variable" must be set when processing the data [37,38] . The model was constructed using a grid form with a side length of 10 mm/cube. For example, the model d belongs to a two-layer type, where the first layer is 40 mm+40 mm, the second layer is 60 mm. Here it takes the convex and ladder model as an example, the transmitter and receiver were in the middle of the X-axis, The outer rectangular box is the UPML boundary layer. The unit of scale in the X-axis and Y-axis direction is meter. The origin of the three-dimensional coordinate axis (0, 0, 0) is in the lower-left corner, and the height of the upper surface of the simulated soil layer is 0.17 m. The spatial solution of 2 mm was chosen based on the previous section. There are 5 Yee cells/mm. Firstly, the data set was chosen which was generated randomly 200 GPR images with thicknesses between 10 mm and 50 mm using the 2-D module in gprMax, where this dataset was known angles. The number of horizontal grids in the sample data was selected randomly within 1200 cells (240 mm). It represents the object width (X-axis direction). The number of longitudinal grids (Y-axis direction) is selected randomly from 1050 cells (210 mm) to be constructed. It should be noted that those grids have a corresponding relationship with the YEE cells in the FDTD. There are 300 GPR images of data in 1050 Yee cells unit. Because the FDTD can construct regular cell boundaries only, it must refer to Newton-Leibnitz's idea of fitting the curve to meet the accuracy range of the demand when we try to construct a curve dataset. Each grid (square) in this paper contains 5 cells. The simulated soil parameter refers to typical tillage soil in southern China. The relative permittivity in this experiment was set ϵ soil as 3, 5, 7, 9, respectively, conductivity was set as σ soil =0.001 s/m to 0.01 s/m. In addition, it can be set the relative permittivity ϵ air as 1 and conductivity σ air =0.00001 s/m to simulate air (ideal non-magnetic material). As showing in Figure 5, the domain boundary is 0.240 m in the X-axis direction, 0.210 m in the Y-axis direction and the model thickness was set to 0.002 m. The outermost rectangular box has the exact matching layer boundary.
The coordinate positions of the transmitting antenna and receiver have been fixed in the coordinate system already. The lower-left corner in the figure is the origin of the 3D coordinate system (0, 0, 0). The simulated soil layer can be set by the 'BOX' command in the gprMax material model and decided its size. Taking one model as an example, firstly, it should be set the contour coordinates values of the point A, B, and C in the coordinate system, so the next step can guarantee the size of the object not to exceed the finite boundary. Secondly, it can mesh irregular boundaries in the soil layer using the three steps described above (marking: boundary 1, 2 and 3), but the discretization step must ensure the accuracy of the research problem.

E z field components analysis
The signal between 0-0.40 ns was a direct wave between transmitter and receiver. Secondly, the electromagnetic signal reached the object upper surface and passes through a thickness of 10 mm at 0.46 ns. According to the time node, the average thickness of the reflected wavelet will be calculated to be approximately 9.982 mm and the absolute error was 0.018 mm. Similarly, when the thickness of the model increased to 30 mm and 50 mm, the absolute error in the Y direction can be calculated to be 0.021 mm and 0.024 mm, respectively. As shown in Figure 6, the time response of the E z field recorded 50 cells under the Ricker source point. The electromagnetic wave intensity of Yee-FDTD was more effective than HIE-FDTD and ADI-FDTD, but US-FDTD is the most significant because US-FDTD uses a second-order iteration equation to solve the boundary space problems when in the set soil multi-scale parameters [39,40] . As time goes by, the "distance" of US-FDTD and Yee-FDTD's E z in the cell unit becomes larger, e.g., D 3 in the figure is greater than D 5 , and D 6 is greater than D 4 . Those effects were relative to ADI-FDTD and HIE-FDTD. Because there is only a need to find the most significant electromagnetic field forwards in order to avoid the re-projection in the computational efficiency of processing UMPL or CMPL.
The US-FDTD statistical results are shown in Table 1, conclusions can be verified that the reflection coefficient of the object interface was reduced gradually as the thickness of the model was becoming larger. The more layers between the object models, the more obvious the A-scan interference, and when the energy loss of the electromagnetic wave was increasing gradually, the intensity of the reflected wave will be reducing. The mean absolute error (MAE) and root mean squared error (RMSE) were two of the most common metrics used to measure accuracy for continuous variables. In this article, there will use these two indicators to judge and evaluate the object recognition indicators of the model.   Figure 7a shows the model A，there were the radar surface direct wave between 0.03-0.05 ns. When time changed from 1.4 ns to 1.7 ns, a strong reflection arose, the waveform was flat and the curve opening was wide. There was also a noticeable electromagnetic wave signal between 2.7 ns and 3.0 ns, which corresponds to be the upper and lower sides of model A (cuboid). The flatter the object boundary was, the larger the response signal curve open. Figure 7b shows model B, there were obvious parabolic reflections between the two cubes. It was different from model A. There existed two relatively high and weak reflected waves between the trace numbers 5 to 15 and 55 to 65. Figure 7c is the model C and this model had the less diffraction wave (interactive interference signal). The closer the object boundary was to the parabolic dispersion direction of the excitation source waveform reflected on the object surface, the less diffraction occurs after the electromagnetic signal passes through the object. Figure  7d shows that the response wave of model D and this concave was constructed by adding a connection "bridge" in the middle of model b. It is obvious that at a position where the object boundary distance was similar, it was easier to generate a strong reflection of the signal. If the surface of the object was a flat wave reflection and there is a strong reflection on both wings, then the object is likely to have a double reflection boundary. Figure  7e is the model E, the upper surface width of model A was 80 mm, and this model had a width of 40 mm. The model E's parallel state of the reflected wave was shorter than model A. There were also strong reflected waves formed on both sides near the 2.0 ns due to the cubes distributed on the left and right sides. Figure  value of 0.83 mm. The conclusion proves that the width and depth features of the boundary can be predicted and using the US-FDTD forward modeling to estimate can meet the accuracy requirements. As shown in Figure 9, after deal with absorbs surface waves by CPML, The magnetic field signal of ADI-FDTD was not much different from that of US-FDTD. It can be focused on the CN condition for the solution involving current sources and lossy media, so the time iteration was longer than other algorithms.

Modeling accuracy and analysis
To discuss the results when the number of the feature patterns was 200, 160, 120, 80, and 40, respectively. It should be noted that among the 200 testing samples, 100 GPR images will meet the single-scale features, and another 100 are curved or multi-scale features. Each object composed by the Yee-FDTD cells has two planes that are identical and perpendicular to the ZOX coordinate system, and that the object has freedom degrees in XOY. The best performances can achieve when the testing sample was 200 and the values of RMSE and MAE were obtained optimally throughout this experimental group. As the number of testing patterns decreases by 40 each time, the value of RMSE does not increase. When the number of testing patterns was between 200 and 80, the value of MAE does not change significantly, which indicated that the performance of the model with a sample size of 80 or more was not significant. In addition, it can be seen from the figure that the larger the training sample capacity, the smaller the value of RMSE.

Figure 10 RMSE versus surface depth
This section tested how the depth and width of the object's UPML and CPML can affect feature extraction. The height in the Y-axis direction was chosen to be 40 mm (200 cells). First, the test experimental group was divided into ten types in surface depth. Those simulated objects all have apex features on the hyperbola of the electromagnetic wave signals, and their waveforms correspond to objects in different depths. According to the above section, we know that when the testing samples are between 80 and 200, the sensitivity of the MAE is not high. As shown in Figure 10 below, the value of RMSE decreases as the depth of the training sample consisting of simplified boundaries increase in soil. It indicated that the deeper the object's depth in soil, the worse PML scattering influence accuracy obtained. The multi-scale type is a rule that exhibits a larger value of RMSE as the depth was deeper both under the CPML and UMPL conditions. This shows that the boundary influence of the object's depth on the single-scales was not obvious, and the influence on the multi-scales was significant. To divide the model boundary from 20 mm (100 cells) to 100 mm (500 cells) into 5 kinds for analysis. The boundary size of the simulated soil target and the number of multi-scale boundaries were generated by gpr-Max in advance randomly. In addition to controlling the width of the object, the rest of the parameters were unchanged. As can be seen from Table 2, the US-FDTD's RMSE had affected by size on both single boundary (SB) and multiple boundary (MB) types. As the size changes, the apparent single boundary sample achieves the best performances in 100 mm. This shows that the larger the boundary width was, the better the model results were.
The SB and MB only have three classification errors caused by the calculation error (marking "miss value" in Table 2). According to the testing results, the boundary quantity feature is ideal. This is because the electromagnetic wave signals will interfere with each other. Although there had added the UMPL setting (10 cells, 15 cells, 20 cells) at the soil layer boundary, the distance between the object and the UMPL boundary is still affected by the excitation source pulse or the object size. The ADI-FDTD's RMSE data of the multi-boundary indicates that the larger the object boundary size extension was, the more significant the effect was. This shows that although the settings between boundaries also affect each other, in the HIE-FDTD cells, there have already been divided unite of the object according to the spatial solution before the experiment, so the multi-boundary effect will offset some of the signal interference. In summary, it is feasible to use multiple FDTD methods to simulate the maximum class interval variance between boundary types. The HIE-FDTD is suitable for inversion of multi-scale structure object modeling in farmland, while ADI-FDTD and US-FDTD are ideal for single-boundary objects in both UMPL and CPML finite element space.

Modeling efficiency and computer performance
FDTD-based method's accuracy has been verified, then, this chapter will compare their computer performance. The modeling computer has 4 CPU threads, and its related parameter is 1×Inter (R) Core (TM) i5-8250U CPU@1.60GHz (4 cores, 8 cores with Hyper-Threading). Used three kinds of waveform sources, Ricker, Blackman-Harris and Gaussian, to compare their CPU running time. The spatial discretization of the main compare model was set 0.002 m×0.002 m×0.002 m. The source frequency was 5GHz and the time-window was 5 ns. The simulation time step was set as 4.717e−12 s, and the domain size is 0.24 m×0.21 m×0.002 m (12600 cells). Test mode is 2D TM-z and B-scan computer performance. Figure 11 shows that all three incentive sources have achieved good performance or results, When the number of 2D slices is more than 75, Ricker's CPU runtime is rising significantly. When the number of slices is more than 85, the Blackman-Harris and Gaussian CPU running time are similar. As the number of slices increases, the Ricker model is stable relatively. As shown in the preliminary conclusions of Section 2.2.1, the Ricker wavelet is a zero-phase wavelet with a short duration and fast convergence. In addition, the speed of structure cells is also a very important indicator of FDTD simulation computer performance. This study uses iMac15, Inter (R)-i7-4790KCPU module in gprMax as the benchmark, and the models feature different cubic domain sizes (from (100×3) cells to (400×3) cells). Simulation waveform used built-in dipole source which is included in the benchmarking python module in gprMax. Figure 12a shows results and it was found that the model proposed in this paper has achieved good performance in the same model of computer operating environment. Even at the frequency of CPU 1.60GHz, when the cubic domain size is 400 cells, the modeling speed reaches 86 M-cell/s. Comparing the benchmark of the i7-4790k, in the cubic domain size between 100-200 cells, the efficiency of the i5-8250 is also ideal. Figure 12b shows the comparison of modeling speeds among ADI-FDTD, Yee-FDTD and HIE-FDTD. Under the computer condition that the CPU is i5-8250, the performance of HIE-FDTD is above 80 M-cell/s when the domain length is between 50 and 100. As the modeling cell's space increases, the modeling speeds of ADI-FDTD and Yee-FDTD decrease. When the domain length is 250, their speeds of the three models were reduced to zero because of the CPU frequency limit. Summed up and compared with the traditional Yee-FDTD forward modeling algorithm, the HIE-FDTD model is superior. It can meet the expected standards for computer performance and fulfills the need for efficient use of FDTD to simulate different shape boundaries in the soil. a. Side length of cubic domain/cell b. Length of the multi-scale structure Figure 12 Computer performance of three kinds of FDTD model

Conclusions
This study proposed an efficient modeling method that used many improved FDTD methods to simulate the complex boundaries or shapes of objects in the soil, and solved the problem that traditional FDTD algorithm cannot construct non-mesh objects by utilizing the diversity characteristics of YEE cell elements, and proposes the idea of simplifying multiple boundaries. Then compares the simulation results with ADI-FDTD, HIE-FDTD and US-FDTD computer performance indicators. After discussion and research, the conclusions of this paper were divided into the following three points: 1) The HIE-FDTD is suitable for inversion of multi-scale structure object modeling in farmland, while ADI-FDTD and US-FDTD are ideal for single-boundary objects in both UMPL and CPML finite element space. This provides important theoretical support for the non-destructive detection technology, especially for agricultural machinery; 2) The HIE-FDTD Modeling results and efficiency are best in the boundary inversion of multi-scale objects. Even if the computer condition is not high, the HIE-FDTD model still has better than the traditional Yee-FDTD forward modeling algorithm. Its modeling speed can reach 86 M-cell/s; 3) For the multi-scale objects in farmland soils, the modeling methods using CPML combine with the HIE-FDTD are the most efficient and accurate ways. The FDTD modeling method proposed can efficiently and accurately simulate the electromagnetic effects of multi-boundary objects. This will provide an efficient modeling method for GPR ground-penetrating radar and an efficient reference for computer performance of 3D-FDTD based modeling.