Abstract
According to the characteristics of the polarizability in frequency domain of three common models of dispersive media, the relation between the polarization vector and electric field intensity is converted into a time domain differential equation of second order with the polarization vector by using the conversion from frequency to time domain. Newmark βγ difference method is employed to solve this equation. The electric field intensity to polarizability recursion is derived, and the electric flux to electric field intensity recursion is obtained by constitutive relation. Then FDTD iterative computation in time domain of electric and magnetic field components in dispersive medium is completed. By analyzing the solution stability of the above differential equation using central difference method, it is proved that this method has more advantages in the selection of time step. Theoretical analyses and numerical results demonstrate that this method is a general algorithm and it has advantages of higher accuracy and stability over the algorithms based on central difference method.
1. Introduction
In 1966, the finite difference time domain (FDTD) method [1] is put forward by Yee. With the further development, the FDTD method probably becomes the most popular numerical technique for solving complex electromagnetic problems [1–8]. FDTD is a time domain algorithm but the constitutive relation of dispersive medium is given by a frequency domain expression. Hence, it is difficult to use FDTD method directly in analyzing electromagnetic characteristics of dispersive medium. According to the different treatments of constitutive relation, the FDTD methods used in dispersive medium can be classified as RC method [5, 6], JEC method [7], PLJERC method [8], ADE method [9], transform method [10–12], SO method [13–15], SARC method [16], and so on.
In the case of nondispersive medium, the spatial discretization of Maxwell curl equation can be determined by the requirement from numerical dispersion to stability, and the time discrete can be determined by Courant stability condition. So the calculation involves just the first order partial differential equation. The time domain expression of constitutive relation often can be written as a second order partial differential equation in the case of dispersive medium. To obtain an accurate solution of the partial differential equation, time step must be very small. In this paper, the general differential equation of dispersive medium is given and the Newmark method [17–19] is introduced to the treatment of the constitutive relation in dispersive medium. This new method has a higher accuracy and stability than the central difference method. The detailed process is as follows.
According to the polarizability characteristics in frequency domain of three common models (Debye, Drude, and Lorentz models) of dispersive medium, the relation in frequency domain between and is converted into a second order differential equation (has a general form to the three common models) by using the conversion from frequency to time domain. Then the Newmark method which is used in FETD is employed to solve the equation. Hence the recursion is derived, and the recursion is obtained by constitutive relation.
2. The General Newmark-FDTD Method
2.1. Maxwell Equation in Dispersive Medium and Constitutive Relation
Maxwell curl equation in linear isotropic dispersive medium is The electric constitutive relation is where and are permittivity coefficient and polarizability in frequency domain, respectively. is the relative permittivity when frequency is infinite.
The FDTD iterative computation of can be accomplished by the discretization of (1) according to the standing Yee cell. And it is needed to deal with constitutive relation in the iterative computation of .
2.2. The United Expression of Polarizability Function of Common Dispersive Medium Models
Three common dispersive medium models are Debye model, Drude model, and Lorentz model [2]. Their polarizabilities can be written as [20] where is numbers of poles, represents the th pole, , is the relative permittivity at zero frequency, is the relaxation time of pole, is collision frequency, is Drude frequency, and is the natural frequency of oscillator.
Let of the three dispersive medium models in (3) can be written as a united form: where
2.3. Constitutive Relations in Time Domain and Its Newmark’s Solution
Polarization vector is introduced: The constitutive relation equation (2) can be rewritten as Let Then from (3) and (7) we know that From (4) and (9), we have After substituting (5) into (11) and using conversion from frequency to time domain, we have in the above differential equation with second order time derivative is equivalent to a excitation source. In the year of 1959, a time-stepping algorithm for solving the differential equations like (12) is given by Newmark according to the time derivative difference approximate solution of Taylor series and the mean value theorem. This algorithm provides higher accuracy and stability [17]. By the deriving method of Zienkiewicz, solution of (12) is [18] where Equation (13) is called two-step or Newmark method [19]. The calculation of involves the values of and . The coefficients in (13) include two variables and . The stability and accuracy of (12) when selecting different and are also discussed in [19] in detail. Usually, and . To assure the accuracy, it is needed to choose an appropriate and . According to [19], unconditional convergence can be obtained by taking , . The time-stepping calculation from to can be finished by (13).
After converting the constitutive relation equation (8) into discrete time domain, we have Substituting (13) into (15), we know that Then the time-stepping calculation from , to can be accomplished.
In conclusion, the general Newmark-FDTD method in dispersive medium can be reduced as follows.(1)Use the difference discrete expressions in curl (1) to get .(2)Use (16) to get .(3)Use (13) to get .(4)Back to step 1.
3. The Advantages of Algorithm in This Paper
For the second order differential equations like After using central difference method, we have
The time domain solution stability of the above differential equation can be analyzed by using the transform. To a time series of , its transform is defined as is the transform of the sampling sequence .
transform is proposed to (18); we know that
According to the final value theorem, time domain function has time stable solutions when the poles of located in a unit circle on complex plane in domain. An analysis of the solution stability of second order homogeneous differential equation without first order derivative is made in [21]. In this paper, the case of including first order derivative is analyzed. Equation (20) is considered as a quadratic equation with one unknown using as dependent variable. The equation can be solved and let ; we know that
The corresponding cubic equation with one unknown to the above inequality has three solutions:
Let be all positive quantities. The derivative of the cubic equation with one unknown is at , which is an increasing function satisfying the requirements of the inequality. Hence, when central difference method is used to solve (17), the time-stepping is limited to According to (23), the requirements of time-stepping when central difference method is used to solve the second order time domain differential equations corresponding to the three models of dispersive medium in (6) are That is to say, for Debye medium, corrected results cannot be obtained whatever the value of the time-stepping takes when central difference method is used to solve (17). For Drude medium, corrected results can be obtained whatever the value of the time-stepping takes. But for Lorentz medium, the solution stability is relevant to its resonance frequency at poles.
However, the solution of (12) is unconditionally stable when using Newmark difference algorithm. Hence, the selection of time-stepping is more flexible. Meanwhile, a higher accuracy is obtained because the Newmark method has a combination of Taylor expansions and means value theorem is employed.
To analyze the difference between using Newmark method and central difference method, the results obtained from these two methods and Runge-Kutta solutions [20] are given. Figure 1 is the absolute error comparing the classical Runge-Kutta solutions of using these two methods, in which the solid line and the circle are the result from Newmark and central difference method, respectively. As it is shown in the figure, the time it takes to get steady state for Newmark method is shorter than central difference method. Let in (17).
The theoretical analysis and numerical results show that the algorithm presented in this paper has higher steady state and calculation accuracy than central difference method.
4. Numerical Results
To demonstrate the correctness of our algorithm, the backscatter of three different dispersive media sphere is given and compared with Mie series solution. Let , , , where is spatial discrete grid-scale of Yee cell. Consider where and is the pulse peak time.
Example 1 (the backscatter of Debye sphere). Let the radius of the sphere be 0.25 m and filled with single pole Debye medium. The complex relative permittivity is where , , s. The backscatter RCS is shown in Figure 2. The solid line and circle are the results from this method and Mie series solution, respectively. It is shown that these two results are in good agreement with each other. Consider min FDTD calculation.
Example 2 (the backscatter of plasma sphere). The radius of the sphere is 3.75 mm and filled with single pole Drude medium. The complex relative permittivity is where rad/s, Hz. The backscatter RCS is shown in Figure 3. The solid line and circle are the result from this method and Mie series solution, respectively. It is shown that these two results are in good agreement with each other. Consider mm in FDTD calculation.
Example 3 (the backscatter of Lorentz sphere). The radius of the sphere is m and filled with single pole Lorentz medium. The complex relative permittivity is where , , Hz, Hz. The backscatter is shown in Figure 4. The solid line and circle are the result from this method and Mie series solution, respectively. It is shown that these two results are in good agreement with each other. Consider m in FDTD calculation.
5. Conclusions
In this paper, a general Newmark-FDTD algorithm is given to deal with the electromagnetic problems in dispersive medium. This new method combines the Newmark difference method which is widely used in FETD calculation in dispersive medium. Theoretical analyses and numerical results demonstrate that this algorithm has advantages of higher accuracy and stability over the algorithms based on central difference method.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This project is supported by the National High Technology Research and Development Program of China (2012AA01308), the National Basic Research Program, and the National Natural Scientific Foundation of China (61231003 and 61401344).