An Unconditionally Stable Cylindrical FDTD Method to Analyze the EM Ground Wave Propagation

We extended the unconditionally stable associated Hermite (AH) FDTDmethod to the cylindrical coordinate system for analyzing the electromagnetic (EM) ground wave propagation. With AH domain differentiate operator technology and paralleling-in-order solution scheme, the cylindrical time-domain Maxwell equations with complex frequency shifted perfectly matched layer (CFS-PML) are transformed to AH domain generating a five-point banded equation, which also unified the formula from the central axis of cylinder. The numerical results represent the accuracy and efficiency of the proposed method when comparing for the conventional FDTD method.


Introduction
In the finite-difference time-domain (FDTD) method [1], the time step size is limited by the Courant-Friedrich-Levy (CFL) stability condition [2].So for the problems with fine structures, it should choose a small time step and then a large number of marching-on-in-time steps for simulation, which is a time-consuming process.To eliminate the stability condition for fine structure, many unconditionally stable (US) methods were proposed, such as alternating direction implicit (ADI) method [3] and orthogonal function-based methods [4][5][6][7][8][9].The orthogonal function-based methods can be classified into two main categories: the marchingon-in-order scheme and the paralleling-in-order scheme.The former is based on weighted Laguerre polynomials (WLP), and the latter is based on associated Hermite (AH) functions.For the latter, all of the works are based on in Cartesian grid previously.
In many applications, such as EM ground wave propagation from a lightning channel, a cylindrical coordinate is preferred.Here, we extend the AH FDTD method to cylindrical coordinates system to analyze this case.The time-domain Maxwell equations in cylindrical coordinates are transformed to AH domain using a AH differential operator technology [6].Then a five-point banded equation is derived for calculation of expanding coefficients.And a paralleling-inorder scheme [6] is used for indirectly but efficiently solution.For the CFS-PML is wildly used to effectively attenuate evanescent waves [10,11], we apply it to AH FDTD method and make a unified formula including a treatment of central axis of cylinder when following the Ampere law.One should note that although CFS-PML is used in [6,9], it is implemented in Cartesian coordinate system.In this paper, we firstly extended it to cylindrical coordinate system, which might give AH FDTD method much more convenient to analyze the scattering problems for the cylindrical shapebased objective.Finally, a numerical example is used to validate the accuracy and efficiency of the proposed method.
Then, a temporal Galerkin testing procedure is used to eliminate the time variables.By applying central difference scheme, we can transform (1), (2), and (3) to AH domain as The matrix a above is AH domain differential operator [6], which can be considered as the similar time-domain operator d/dt or frequency-domain operator jω and respective identity matrix I is for 1.By assembling (4), (5), and ( 6) and eliminating the electric field components, we can derive a five-point banded matrix equation with H φ as where Δz k+1 /Δz k , International Journal of Antennas and Propagation When for the boundary condition of central axis, it needs to be treated separately owing to singularities.According to Ampere's law [12], Then, the differential form for (14) can be obtained as follows:

z, t 15
Converting it to the AH domain by using the AH domain differential operator a and performing discretization for i = 1 By applying ( 16) to ( 5), the modified magnetic field (12) in central axis can be derived as follows: After the central axis formula to modify ( 12) is updated, we can apply the lower-upper (LU) decomposition procedure to calculate this equation.Also, we can perform parallelingin-order solution scheme [6] for (12).If every magnetic field variable H φ | i,j in the entire computational space is calculated, the results of electronic field variables E r | i,j and E z | i,j can be gotten from ( 4) and ( 5).In addition, only the lossy medium is considered above, while for the case of dispersive medium, the formulation derivation is almost the same as before, except for a slight modification of ( 8) with ε ξ being replaced by its AH domain form from [7].

Numerical Verification
To validate the performance of the proposed AH FDTD method, a numerical example with an EM ground wave propagation problem is presented, as shown in Figure 1.The computational domain is discretized into 50 × 50 cells with Δr = Δz = 0 1 m, including an axial symmetric boundary condition for the left boundary and the other three boundaries terminated by 10 PML layers.The background material is set to be earth with medium parameters: ε ξ = 12, σ eξ = 0 05, μ φ = 1, and σ mz = 0.The PML parameters are scaled following the expressions in [11].
where ξ 0 represents the interface between FDTD and PML grids, d is the thickness of the PML, and m = 4 is the order of the polynomial.A Gaussian pulse source J z modulated by sinusoidal signal is located at the grid (0, 27).
To verify the performance of the proposed method, the result from conventional FDTD method is used as reference.The time step for conventional FDTD method is chosen as Δt = 0 21 ns to meet the CFL condition, and the time duration is chosen as 127.4 μs, which leads to the total timemarching number of 600,000.While for AH FDTD method, Δt can be chosen as 21 ns to evaluate expanding coefficients of source.AH parameters are selected as the number of AH functions Q = 48, the time-translating parameter T f = 63 7 μs, and time-scaling parameter l = 6 8 × 10 −6 , respectively.PML parameters are selected to meet the best performances in this example as η ξ = 0 005, κ max = 7, and σ max = 1 2σ opt .
The magnetic fields H φ at p 1 and p 2 are shown in Figure 2, where we can see the good agreement between the proposed method and the conventional FDTD method.
The comparison of computational resources between the AH FDTD method and the conventional FDTD method is listed in Table 1.Although the proposed method has much 3 International Journal of Antennas and Propagation more memory consumption, it gives much higher efficiency in CPU time and there is about 95% reduction of computing time in contrast with conventional FDTD method.

Conclusion
The AH FDTD method is developed to analyze the ground wave propagation problem under the cylindrical coordinates in this paper.With the CFS-PML absorbing boundary, a unified banded equation, including the central axis formula, is derived in AH domain by using the paralleling-in-order solution scheme and AH differential operator technology.The numerical example with a cylindrical ground wave propagation case is given, which also makes a good verification of the performance for the proposed method.Therefore, for its good performance of numerical dispersion, the AH FDTD method has a potential to analyze the propagation of waves on long distances such as ground penetrating radar (GPR) applications, where segmentation technology should be explored to deal with a probable long-time response (e.g., signals as biexponentials) in the future.

Figure 1 :
Figure 1: The configuration of the computational domain.

2 Figure 2 :
Figure 2: Comparison of the magnetic fields calculated by proposed method and conventional FDTD at two observation points.

Table 1 :
Comparison of the computational resources.