Abstract

A novel hybrid algorithm is proposed to reduce the computation cost of the finite-difference time-domain (FDTD) method in calculating the transient near-field scattering from rough surface. The scattering problem is split into the FDTD calculation of equivalent sources on the contour enclosing rough surface and the calculation of the near-field radiation with the two-dimensional (2-D) time-domain Huygens’ principle. The radiation fields are found from a surface integral of the temporal convolution for which the direct numerical integration of the convolution is computationally expensive. In this paper, the 2-D time-domain Green’s function as the convolution kernel is approximated with a sum of exponential terms by using the Prony’s method. Then, the semianalytical recursive convolution (SARC) approach is applied to complete the update of the near-field radiation. Compared with the traditional FDTD, this hybrid algorithm can significantly reduce the memory usage and run time, especially for the large distance between the rough surface and observation point.

1. Introduction

With the extensive application of the wide-band radar, the transient electromagnetic (EM) scattering from the rough surface has gained increasing attention, which plays important role in the field of radar surveillance, EM imagining, target identification, and stealth, etc. [14]. The previous researches are mainly focused on the simulation of the EM scattering from rough surface in the far-field zone [57]. In fact, the near-field scattering problem is the widespread presence of the target working above background, such as the low-altitude cable over ground or sea. Only by mastering the near-field scattering characteristics from rough surface, we can systematically analyze the EM interference from background and obtain the relatively complete coupling result of target. Therefore, it is of critical importance to study an efficient time-domain algorithm on the near-field scattering for the coupling analysis of the low-altitude target and the protection of the sensitive electronic device.

In the past decade, the numerical method has been widely used in the EM analysis owning to the rapid development of the high performance computer [811]. Among them, the finite‐difference time-domain (FDTD) method is a powerful and popular technique for numerically solving the transient field problem, and it has achieved successful applications in the far-field scattering from rough surface [1214]. However, when the FDTD is used to simulate the near-field scattering, a brute-force way is to enlarge the computational domain to enclose the rough surface and observation point, which will lead a heavy burden to the computation resource available. Moreover, the farther the near-zone observation point is from rough surface, the heavier the computation burden will be. By now, some hybrid algorithms combining the FDTD and analytical-approximate method have been proposed to reduce the FDTD space lattice and improve the computational efficiency. In [15], the FDTD associated with the Kirchhoff surface integral was proposed to calculate the near-field response of the dipole antenna in the three-dimensional (3-D) case, where the dipole antenna was assigned in the FDTD region, and the near-field response outside the FDTD space was calculated by the Kirchhoff surface integral. For the two-dimensional (2-D) problem, by means of the 2-D time-domain Huygens’ principle, the temporal convolution integral needs to be solved to obtain the near-field response exterior to the FDTD region, which complicates the computation even more. The direct integration of the temporal convolution is computationally expensive, and hence it is applicable for the EM analysis of the small model such as the infinitely long line source in [16, 17]. Xu et al. in [18] employed the trapezoidal approximation to solve the temporal convolution, and then the near-field scattering exterior to the FDTD region was computed by applying the scheme of treating the radiation field as the equivalent cylindrical wave, which can be accurately used for the target with cylinder shape. However, for the large rough surface model, the characteristics of the structure and scattering wave make the simulation of the near-field scattering with these approaches difficult. In recent years, a fruitful semianalytical recursive convolution (SARC) approach was introduced in the dispersive medium problem, realizing the high-efficient computation of the convolution integral [19, 20]. But when the 2-D time-domain Green’s function is acted as the convolution kernel, the recursive operation is difficult to realize due to its specific structure. In this paper, this difficulty is overcome by applying the Prony’s method to approximate the convolution kernel with a sum of exponential terms, which avoids the complex analysis of the convergence of the series expansion scheme in [21]. The hybrid algorithm combining the FDTD with the 2-D time-domain Huygens’ principle is proposed to calculate the transient near-field scattering from the rough surface, in which the normal part of the convolution integral is calculated by the SARC approach and the singular part is done using the linear interpolation. Because a great deal of FDTD meshes between the rough surface and near-zone observation point are removed, the memory usage and run time in the proposed algorithm are dramatically reduced than those in the traditional FDTD method. Meanwhile, this hybrid algorithm retains the unique advantages of the FDTD method.

The rest of this paper is organized as follows: in Section 2, the hybrid algorithm combining the FDTD with the 2-D time-domain Huygens’ principle is elaborated, including the Prony approximation of the convolution kernel and the SARC approach for solving the near-field radiation of the equivalent sources. In Section 3, the accuracy and efficiency of the proposed hybrid algorithm are verified through the comparison of the numerical results by the hybrid algorithm and by the traditional FDTD method. Section 4 draws the conclusions of this paper and offers the recommendations for future work.

2. Hybrid Algorithm for Near-Field Scattering from Rough Surface

For the FDTD calculation of the transient near-field scattering, a common way is to expand the FDTD region to enclose the rough surface and the near-zone observation point, as shown in Figure 1(a). With the increase in the height of observation point, the memory usage and run time will sharply increase, thus often making it difficult to realize the FDTD calculation of the near-field scattering. For the proposed hybrid algorithm, however, the FDTD region only contains the rough surface as indicated in Figure 1(b); the near-field response at the observation point is regarded as the radiation of the equivalent sources on the output boundary of the FDTD region. Comparing Figures 1(a) and 1(b), it can be seen that the computational domain in the hybrid algorithm is much smaller than that in the FDTD method, especially in the case with high observation point.

2.1. Theory of Near-Field Transformation in 2-D Case

Assume that and are respectively the electric and magnetic fields on the output boundary surface of the FDTD domain. According to the equivalence principle, the surface electric and magnetic currents and charges can be obtained bywhere denotes the unit normal vector on the output boundary. For the 2-D problem, the vector and scalar potentials can be expressed with the 2-D time-domain Green’s function as

Here, and represent the permeability and permittivity, is the distance from the observation point to source field , and is the speed of the EM wave. In accordance with the vector and scalar potentials, the 2-D radiation field exterior to the FDTD region can be written as

Substituting equation (2) into equation (3) and shifting the time reference backward by , the formals for calculating the 2-D radiation field of the equivalent sources can be obtained from [21], in which, for the TM case, are

The corresponding formulas in the TE mode can be derived with the principle of duality, which are not described here. It is obvious that the convolution kernel in equation (4) has a square-root singularity at . An efficient way in addressing the singularity is to approximately calculate the convolution integral in the interval using the linear interpolation. Assuming , the time derivative of the equivalent source, e.g., , in the time-reversed and offset form can be given by

For instance, substituting (5) into the second term of (4a), we have

A similar procedure can be implemented for other terms in equation (4) to obtain the corresponding approximate results in the interval . The convolution integral in the interval is computed by using the high-efficient SARC approach, as will be elaborated later.

2.2. Prony Approximation of Convolution Kernel

Although equation (4) is applicable for the computer execution, the direct integration of it is such cumbersome as to be impractical, due to the need for storing and processing all past values of the time derivative of the equivalent source. The computational complexity of this scheme scales as , where is the total time step. In this paper, the SARC approach will be applied to decrease the complexity of the convolution integral. It can be found from equation (4) that the convolution kernel mainly consists of two forms: and , which make it difficult to execute the recursive operation of the convolution. Based on the idea of the signal estimation in the high-speed circuit [22], the Prony’s method is used to approximate the convolution kernel with a sum of the exponential terms:

Here, is the discrete sample of the convolution kernel at , where denotes the equally spaced time step, represents the total number of poles, and and represent the coefficient and amplitude of the th pole. The and are evaluated in the Prony’s method by solving two sets of linear equations with an intermediate nonlinear equation with th degree. Firstly, equally spaced time samples of the convolution kernel are taken to set up a set of linear equation aswhere . The roots can be determined by solving this oversampled set of linear equations with the least-squares approach. By virtue of equation (8), the can be found as the roots of the following th degree equation:

If the complex conjugate roots arise in the solution to equation (9), use their magnitude to substitute for the complex pair of , and cancel any duplicate. Substituting to equation (7) and introducing discrete samples of the convolution kernel, a new set of the linear equations with the values of the can be established as follows:where . The coefficients can be obtained by using the least-squares approach to equation (10).

Figure 2 shows the comparison of the fitting results by the Prony’s approximation with the analytical results, in which the red short dash denotes the analytical result and blue line denotes the fitting result with Prony’s method. For the fitting calculation, the space increment is , the time increment is , and the observation point and source point are located at and , respectively. The sampling of the convolution kernel begins at , and the interval is 20. In Figures 2(a) and 2(b), the and are respectively expressed as convolution kernels 1 and 2 for convenience. The total numbers of the samples and poles are selected to be 30 and 10, respectively. It can be observed from Figure 2 that the nearly identical results are produced with the two methods, implying the reliability of the Prony’s method.

2.3. Semianalytical Recursive Convolution Approach

The RASC approach is a fruitful scheme for the calculation of the convolution integral, which is based on the idea of the efficient analysis of the transfer response in the digital signal processing (DSP) [23]. For a linear time-invariance and dynamical system, as shown in Figure 3, the system response can be expressed as the convolution of the input signal with the impulse response . We have

In general, the impulse response in an asymptotically stable system is of the form,where stands for the number of poles, and is the switch function:

Substituting equation (12) into equation (11), the system response at can be obtained as

The system response for the th pole is

Dividing the integral of equation (15) into two intervals, and , we get to the recursive formula of system response as

In this paper, the near-field radiated by the equivalent source is analogous to the system response in the DSP. The convolution kernel can be regard as the impulse response after it is approximated with a sum of exponential terms. The time derivative of the equivalent source corresponds to the input signal. Thus, the near-field radiation can be updated by applying the SARC approach. Next, focus on the SARC approach for evaluating the radiation field in equation (4), where the representative convolutions are of the forms,

For equation (17a) as example, using the Prony’s method mentioned above can approximate the convolution kernel as

According to equation (16), the temporal convolution in the interval can be expressed aswhere can be approximated in accordance with the linear interpolation approach as follows:

Substituting equation (20) into equation (19) and analytically solving the integral over the interval , we find the following field-update equation:where

From equation (21), it is clear that the SARC approach has a high efficiency in calculating the near-field radiation. This is because this approach avoids the requirement for storing and processing of the complete history of the equivalent source on the output boundary of the FDTD region, and its computational complexity is reduced as only . Beyond that, the SARC approach has a better accuracy than the general recursive approach, due to the analytical calculation of the convolution integral in the discrete time interval.

3. Results and Discussion

Figure 4 shows the geometry of the near-field scattering from rough surface. The observation point is located above the rough surface, and denotes the incident angle. Without loss of generality, the one-dimensional Weierestrass-Mandelbrot function is used here to generate the rough surface [24]. In the numerical calculation, the length of the rough surface is taken to be and the spatial increment is taken to be to allow for reasonable resolution of the rough surface model discretized with the staircase approximation, and the time increment is set to be in accordance with Courant stability condition [22]. The 2-D FDTD space lattice is terminated by the uniaxial perfectly matched layer (UPML) [25], and its thickness is chosen to be . Consider the rough surface is exposed to the Gaussian pulse, whose expression is . The parameters of the incident wave are adopted as follows: and .

To illustrate the accuracy and efficiency of the proposed hybrid algorithm, Figures 5 and 6 present the comparison of the near-field scattering results by the hybrid algorithm with those by the traditional FDTD method for the TE and TM polarized conditions, respectively. In Figures 5(a) and 5(b), the near-zone observation points are respectively located at and , the incident angle is selected to be , and the root-mean-square (RMS) height is , the fractal dimension is , and the permittivity of the rough surface is . In Figures 6(a) and 6(b), the near-zone observation points are positioned at and ; the rest of the parameters are shown as follows: , , , and . It can be observed from Figures 5 and 6 that the time-domain scattered waveforms for the two methods are in good agreement for different conditions. This result demonstrates the accuracy of the hybrid algorithm proposed in this paper. Beyond that, hybrid algorithm can dramatically lessen the computation cost compared with the traditional FDTD method. For the TE polarized wave incidence and 4000 time-steps updating as example, the computation costs required for two methods are provided in Table 1. When the height of the near-zone observation point is , the number of FDTD meshes, the memory usage, and the run time in the hybrid algorithm are reduced to about 72%, 76%, and 74% of those in the traditional FDTD method. For , the key factors above using the hybrid algorithm are only 39%, 42%, and 33% of those using the traditional method, respectively. It is obvious that the hybrid algorithm becomes more efficient with the increase in the height of the near-zone observation point. This is because that the FDTD region in the traditional method has to be greatly enlarged to enclose the rough surface and observation point, whereas that in the hybrid algorithm remains unchanged when the observation point rises. In addition, we can find from Table 1 that the proposed algorithm differs little from the traditional FDTD method in the run time for . This is owing to the fact that, in this case of very low observation point, the advantage of the hybrid algorithm in reducing the FDTD meshes is almost offset by its calculation of the near-field radiation. As the height of the near-zone observation point increases, the proposed hybrid algorithm becomes more favorable in the computation efficiency.

4. Conclusions

In this paper, a novel hybrid algorithm has been proposed to calculate the transient near-field scattering from rough surface for both TE and TM polarization. This hybrid algorithm is based on the combination of the FDTD and the 2-D time-domain Huygens’ principle, and it is the development of the traditional FDTD method. In the hybrid algorithm, the FDTD is only performed in the computational domain containing the rough surface, which makes a large number of FDTD meshes eliminated in the space between the rough surface and observation point. The proposed hybrid algorithm is validated by comparing the near-field scattered waveforms of the rough surface for both the hybrid algorithm and the traditional FDTD method. The numerical results show that the hybrid algorithm is an accurate and efficient technique in the calculation of the transient near-field scattering from the rough surface. Moreover, the relative efficiency of the hybrid algorithm becomes higher with the increase in the height of the observation point. In the future, the hybrid algorithm based on SARC approach will be developed to study the transient composite scattering and coupling of the dispersive target and background.

Data Availability

The simulation data used to support the finding of this study are included within the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China (Grant nos. 62001345, 61901324, and 61961041), in part by the China Postdoctoral Science Foundation (Grant nos. 2019M653548 and 2019M663928XB), and in part by the Fundamental Research Funds for the Central Universities (Grant nos. XJS200501, XJS200507, JB200501, and YJS2105).