Semi-analytical recursive convolution finite-element time-domain method for electromagnetic analysis of dispersive media

Journal Pre-proof Semi-analytical recursive convolution finite-element time-domain method for electromagnetic analysis of dispersive media Linqian Li, Bing Wei, Qian Yang, Debiao Ge

PII:

S0030-4026(19)31652-3

DOI:

https://doi.org/10.1016/j.ijleo.2019.163754

Reference:

IJLEO 163754

To appear in:

Optik

Received Date:

20 August 2019

Accepted Date:

7 November 2019

Please cite this article as: Li L, Wei B, Yang Q, Ge D, Semi-analytical recursive convolution finite-element time-domain method for electromagnetic analysis of dispersive media, Optik (2019), doi: https://doi.org/10.1016/j.ijleo.2019.163754

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier.

Semi-analytical recursive convolution finite-element time-domain method for electromagnetic analysis of dispersive media Linqian Li1, 2, Bing Wei1, 2,* , Qian Yang1, 2, Debiao Ge1, 2

1 2

School of Physics and Optoelectronic Engineering, Xidian University, Xi’an 710071, China.

Collaborative Innovation Center of Information Sensing and Understanding, Xidian University.

[email protected]

ro



of

Xi’an 710071, China

ABSTRACT: Based on the idea of semi-analytical convolution in digital signal processing (DSP), a new technique of finite-element time-domain (FETD) for dealing with dispersive media is

-p

presented. By comparison with semi-analytical convolution in DSP, a unified recursive formulation of semi-analytical convolution for three kinds of dispersive media models i.e. Drude

re

model, Debye model and Lorentz model is described. This formulation includes electric field E and complex polarization vector ψ . On the other hand, the weak form solution of the FETD

lP

equation on account of the idea of DSP and the iteration equation including E and ψ are obtained. Then the achievement of semi-analytical recursive convolution finite-element time-domain (SARC-FETD) method is developed by combining the above two equations. Finally,

ur na

the feasibility of this algorithm is validated with three-dimension numerical examples.

Keywords: Semi-analytical recursive-convolution, Finite-element time-domain, Dispersive media. I.

Introduction

Jo

Most materials in nature are dispersive media. Requirement of electromagnetic properties

analysis in bioelectromagnetics, aircraft stealth etc. is urgent. In dealing with dispersive media, with Fourier transformation, electromagnetic algorithm in time domain can obtain wideband information by once calculation.

Finite element time domain (FETD) method paid more

attention and studied by scholars because of its advantages of unstructured grids whose discrete elements fit well with the surface of the complex geometry [1]-[11]. In dispersive media, due to its dielectric relationship is frequency-dependent; the analysis and

discussion of electromagnetic property directly is complicated in time domain. The previous work on dispersion media processing in FETD method is relatively few. In 1999, the hybrid method of finite difference time domain (FDTD) and FETD is used to handle the dispersive materials through recursive convolution [12]. From 2001 to 2002, the FETD method for three kinds of common dispersive media as well as its stability are studied by Jiao Dan etc by recursively evaluated convolution integral [13] [14]. In 2012, the explicit finite-element time-domain and modeling of dispersive media is investigated by Li Xiaolei [15]. But it's actually a different approach called DGTD. Semi-adaptive time step scheme is adopted in edge-based FETD to

of

simulate the electromagnetic action by Cai etc.[16]. A coordinate stretching-based PML formulation is implemented for the truncation of computational domains containing general linear

ro

and nonlinear dispersive media by David etc.[17]. Most time the relation is dealt with convolution integral, and its stability and generality are remain to enhance [18].

In this paper, a new technique of Finite-element time-domain (FETD) for dealing with

-p

dispersive media is presented. First of all, corresponding to the shock response of system, the time domain polarizability function of three kinds of dispersive media is expressed as a sum of

re

exponent form. By the idea of semi-analytical recursive convolution (SARC) in digital signal processing (DSP) [19]-[21], a unified recursive formulation including electric field E and complex polarization vector ψ is obtained. On the other hand, from Maxwell’s equations, the

lP

wave equation including E and ψ is obtained, then the residual equation is got by Galerkin method [22]. Then, the calculation region is partitioned into many small tetrahedrons elements. Afterwards E and ψ is expanded by basic function, respectively, the matrix equation in the entire

ur na

region is received. The iteration equation including E and ψ is obtained after the above equation is discrete by Newmark’s method [23]. Combining the two equations obtained above, the field can be solved gradually along time axis. Finally, three-dimensional numerical examples are given to demonstrate the validity of our technique.

Jo

II. SARC Techniques in DSP

Based on DSP theory, for a time-invariant and linear system, the system response y  t  is equal

to a convolution of impulse response h  t  and input signal x  t  y t   h t   x t    h  t   x   d t

0

The shock response function is expressed as a sum of exponent form [21]

(1)

Q

h  t    H q exp   q t U  t 

(2)

q 1

where



represents the sum of poles, Q indicates the total number of pole, H q and  q are

coefficients, U  t  is switching function 0, t  0 U t    1, t  0 n n Substituting(2) into(1), and introducing discrete time step t  nt , y  y  nt  , x  x  nt  ,

Q nt

H

q 1 0

exp   q  nt    x   d

q

For single pole, there is [21] yqn  H q exp   q nt 

nt

(3)

ro

yn  

of

the system response can be written in the discretized form as follows

 exp     x   d q

(4)

-p

0

n where yq represents the q -th output of h  t  . Decomposing the integral over the time interval

re

0, nt  into integral over 0,  n  1 t  and  n  1 t , nt  , we can obtain nt



 n 1 t

exp  q   x   d

(5)

lP

yqn  exp   q t  yqn 1  H q exp   q nt 

From the formulation(5), the output signal yqn corresponding to the pole can be calculated by a form of recursive iteration. The formulation described above is deduced based on the exponent

ur na

function (2) of shock response, so it is called semi-analytical recursive convolution (SARC) method.

III. SARC Techniques in Dispersive Media The time domain expression of polarizability functions of three kinds of dispersive media

Jo

models are described as [24]

 Drude  t  

 p2 vc

exp  0  t  

 p2 vc

exp  vc t 

 Debye  t    c  exp   c t   Lorentz  t  

  2 0

  c 2 0

2





Im exp    c  j 02  c 2 

 t 

(6)

The formulations (6) can be described as a general sum of complex exponent function form [21]

Q

  t    Im Gq exp   q t 

(7)

q 1

where Q indicates the total number of pole, Im is the imaginary operator, Gq and  q are coefficients, their values are shown in Table 1 [21] .

The constitutive relation in frequency domain is given by D     0  E     0 P  

(8)

where D , E , P ,  0 and   are displacement vector, electric field, polarization vector, permittivity

of

in vacuum and relative permittivity in infinite frequency state, respectively. And P is described as P         E  

where  is polarizability.

ro

(9)

After Fourier transform, the constitutive relation in time domain is given by

-p

D  t    0  E  t    0 P  t 

where P is

(10)

re

P t    t   E t 

  E  t      d t

(11)

0

lP

From Eq.(11), the conclusion is obtained that the polarization function is similar to the transfer function, E and P are corresponding to the output and input of system as in Fig. 1 shows. So the idea of SARC can be used into the constitutive relation of E and P in time domain.

ur na

E t 

Fig. 1

P t 

 t 

The constitutive relation of E and P in time domain

n n The signal is discrete by intervals t , and D  D  nt  , E  E  nt  . Based on Eq.(7), P can

Jo

be written as

P n   Im  ψqn  Q

(12)

q1

where ψ is complex polarization vector. Eq.(10) is described as [21] D n   0  E n   0  Im ψqn  Q

q 1

According to(5), ψ can be written as

(13)

ψqn = exp   q t  ψqn 1  Gq exp   q nt 

nt



 n 1 t

exp  q  E   d

(14)

E   is replace over the time interval    n  1 t , nt  with the Newton polynomial as

follows E    Ene 

E n  E n 1    n  1 t    t

(15)

where Ene  E n 1 , here we take the first order approximation for(15), substituting it into (14), and transforming integral variable  into      n  1 t , we obtain

of

ψqn 1  exp   q t  ψqn  c0, q E n 1  c1, q E n

(16)

where

(17)

q= 0

q  0

lP

re

-p

ro

 Gq t ， q= 0   2 c0, q =   1  exp   q t    Gq 1  ，  q  0 q   q t     Gq t ，   2 c1, q =    1  exp   q t   Gq   exp   q t  ， q   q t   

Thus, Eq. (16) is the unified recursive formulation for general three kind of dispersive media models. It not only inherits the stability of the SARC method in liner system, but also has a

ur na

favorable versatility. For different models, if only the pole of polarizability  q and the coefficient G are given, the calculation can be completed. q

IV. Matrix Equation of FETD in Dispersive Media

Jo

The second-order wave equation of electric field is given from Maxwell’s equations

where J denotes current density.

1  2 D J   E   2   t   t

(18)

tˆ  zˆ  n

 0 , 0

nt ns

 , ,

t  ext

s

Fig. 2

of

 ABC

Computational region

ro

Substituting (13) into(18), then

-p

Q  2 Im ψ  q    J 1  2 E      E    0  2   0  t t t 2 q 1  

Suppose computation region (shown in Fig. 2) contains dispersive medium, and

(19) 

denotes the

re

surface of region. The first vector Absorbing Boundary Condition (ABC) is used in this paper [25]

(20)

lP

1   n     E   Y n   n  E  t  

where n , Y  1 Z    and n   n  E    Et represents out normal unit vector, admittance of

ur na

medium near to truncate boundary and tangential electric field component, respectively. Eq. (19) with vector boundary conditions (20) is multiplied by weighting functions, and taken integration over computation region as well as along its boundary; we can get the matrix equation. The matrix equation is discrete by Newmark-β’s method, the iteration equation in time domain is obtained

Jo

P  E n 1  Q  E n   R E n -1  S Im  ψqn  exp   q t   2 ψqn   ψqn 1  hn  Q

(21)

q 1

where  P  , Q  ,  R  and  S  are coefficient matrixes,

h  n

is incentive source vector.

Combining Eqs.(16) and(21), then the spatial electric field of three kind of dispersive models

can be solved. V. Numerical results One of the advantages of SARC-FETD method is that it can calculate the mix model of general

dispersive media. The cross section of mixed spheres which consists of Debye, Lorentz and Debye model is shown in Fig. 3. From inside to outside, the Debye model, the Lorentz model and the Debye model are followed. They have the radius of 2×10-3 m, 4×10-3 m and 6×10-3 m in turn. The parameters of

Debye

sphere are:

 s  1.16 ,

   1.01 ,

  2.95 104 Ω

and

 0  6.5 1010 s . The parameters of Lorentz sphere are:  s  2.25 ,  c  4.2 107 Hz and

0  6.0 108 Hz . The parameters of Drude sphere are:  p  6.75 108 rad / s and c  7.5 108 Hz .

The result of the monostatic RCS of mixed sphere of three general dispersive media is shown in

of

Fig. 4. In Figure, the solid line represents results of SARC-FETD method, and the “o” the analytical. Figure shows that the result of our method is in good agreement with the analytic

Cross section of mixed dispersive media sphere

ur na

Fig. 3

lP

re

-p

ro

method.

-30 -40

Jo

RCS/dBsm

-50 -60 -70

SARC-FETD Mie

-80 -90

0

20

40

60

f /GHz

80

100

Fig. 4

Monostatic RCS of mixed sphere consists of Debye, Lorentz and Debye model

The Drude plate is 2 m×2 m×1 m, and geometry is shown inFig. 5. The plane wave is illuminated along the –z axis. The back scattering of Drude plate is calculated when the electron concentrations N e is equal to 3.0 × 1014 m-3, 1.4 × 1014 m-3 and 2.4 × 1013 m-3, respectively. The cyclotron frequency is  c  7.5 107 Hz . The result of monostatic RCS is shown in Fig. 6. It is indicated by the curves in Fig. 6 that different electron concentrations corresponding to different scattering. When the electron concentration is large, the backscattering

of

is relatively large, and when the electron concentration is decreased, the backscattering is

lP

re

-p

ro

decreased.

Fig. 5

ur na

20

The geometry of the Drude plate

10

0

Jo

RCS/dBsm

-10 -20 -30 -40

13

-3

Ne=2.410 m 14 -3 Ne=1.410 m 14 -3 Ne=3.010 m

-50 -60 -70 0.00

Fig. 6

0.05

0.10

0.15 0.20 0.25 0.30 f/GHz Monostatic RCS of Drude plate in different electron concentration

VI. Conclusion The idea of SARC in DSP is introduced into the FETD method. The algorithm of SARC-FETD dealing with dispersive media is presented. This technique not only keeps the flexibility of FETD for complicated Geometry, but also absorbs the advantages of stability and high accuracy of SARC method in linear system. The calculation can be finished by the SARC-FETD iteration formulations and unified program only the polarizability of pole and the corresponding coefficient are given. Finally, the feasibility and versatility of this technique is validated with

of

three- dimension numerical results. Acknowledgments

ro

This research is supported by the National Natural Science Foundation of China under Grant No. 61571348; Pre-research Field Foundation under Grant No.6140518020206; State key

-p

Laboratory Foundation of National Defense Science and Technology under Grant No. 201702007; State Key Laboratory Open Project of Simulation and Effects of Intense Pulse Radiation

re

Environment under Grant No. SKLIPR1705.

lP

References

[1] J. F. Lee and Z. Sacks, "Whitney elements time domain (WETD) methods," IEEE. T. Magn., vol. 31, no.3, pp. 1325-1329, 1995.

[2] S. D. Gedney and U. Navsariwala, "An unconditionally stable finite element time-domain

ur na

solution of the vector wave equation," IEEE. Microw. Guided. W., vol. 5, no.10, pp. 332-334, 1995.

[3] D. J. Riley, J. M. Jin, L. Zheng and L. E. R. Petersson, "Total-and Scattered-Field Decomposition Technique for the Finite-Element Time-Domain Method," IEEE T. Antenn. Propag., vol. 54, no.1, pp. 35-41, 2006.

Jo

[4] X. Wu, Y. Jin and L. Z. Zhou, "Efficient TDFEM schemes with second-order PML based on different temporal basis functions," Waves. Random. Complex., vol. 21, no.2, pp. 199-219, 2011.

[5] Q. He and D. Jiao, "Explicit and Unconditionally Stable Time-Domain Finite-Element Method with a More Than “Optimal” Speedup," Electromagnetics., vol. 34, no. 3-4, pp. 199-209, 2014. [6] J. M. Jin, The finite element method in electromagnetic, New York: Jorn Wiley &Sons. 2014.

[7] D. Y. Na, F. L. Teixeira, H. Moon and Y. A. Omelchenko, “Full-wave FETD-based PIC algorithm with local explicit update,” IEEE International Symposium on Antennas & Propagation. IEEE, 2016. [8] K. Zhang, C. F. Wang and J. M. Jin, “A Hybrid FETD-GSM Algorithm for Broadband Full-Wave Modeling of Resonant Waveguide Devices,” IEEE Transactions on Microwave Theory & Techniques, vol. 99, pp. 1-12, 2017. [9] D. Feng, X. Wang and B. Zhang, “Specific evaluation of tunnel lining multi-defects by all-refined GPR simulation method using hybrid algorithm of FETD and FDTD,”

of

Construction and Building Materials, vol. 185, pp. 220-229, 2018. [10] P. Li, L. J. Jiang, Y. J. Zhang, S. Xu, and H. Bağci, “An Efficient Mode-Based Domain

ro

Decomposition Hybrid 2-D/Q-2D Finite-Element Time-Domain Method for Power/Ground Plate-Pair Analysis,” IEEE T Microw Theory, vol. 66, no.10, pp. 4357-4366, 2018. [11] L. Q. Li, B. Wei, Q. Yang, D. B. Ge, “Shift-Operator Finite-Element Time-Domain Dealing

-p

with Dispersive Media,” Optik, vol. 182, no. 1, pp. 832-838, 2019.

[12] M. S. Yeung, "Application of the hybrid FDTD-FETD method to dispersive materials,"

re

Microw. Opt. Techn. Let., vol. 23, no.04, pp. 238-242, 1999.

[13] D. Jiao and J. M. Jin, "Time-domain finite-element modeling of dispersive media," IEEE Microw. Wirel. Co., vol. 16, no.5, pp. 220-222, 2001.

lP

[14] D. Jiao and J. M. Jin, “A General Approach for the Stability Analysis of the Time-Domain Finite-Element Method for Electromagnetic Simulations,” IEEE T. Antenn. Propag., vol. 50, no.11, pp. 1624-1632, 2002.

ur na

[15] X. L. Li, “Investigation of explicit finite-element time-domain methods and modeling of dispersive media and 3D high-speed circuits,” Dissertations & Theses - Gradworks, 2012. [16] H. Cai, X. Hu, B. Xiong, M. S. Zhdanov, “Finite-element time-domain modeling of electromagnetic data in general dispersive medium using adaptive Padé series,” Computers and Geosciences, vol. 109, pp. 194-205, 2017.

Jo

[17] D. S. Abraham, D. D. Giannacopoulos, “A Perfectly Matched Layer for the Nonlinear Dispersive Finite-Element Time-Domain Method,” IEEE Transactions on Magnetics, 2019.

[18] X. Ai, Y. Tian, Z. W. Cui, Y. P. Han and X. W. Shi, “A Dispersive Conformal FDTD Technique for Accurate Modeling Electromagnetic Scattering of THz Waves by Inhomogeneous Plasma Cylinder Array,” Prog. Electromagn. Res., vol. 142, pp. 353-368, 2013.

[19] W. Janke and G. Blakiewicz, “Semi-analytical recursive algorithms for convolution calculations,” IEE Proceedings-Circuits, Devices and Systems, vol. 142, no.2, pp. 125-130, 1995. [20] W. Pietrenko, W. Janke and M. K. Kazimierczuk, "Application of semianalytical recursive convolution algorithms for large-signal time-domain simulation of switch-mode power converters," IEEE. T. Circuitsir. -I., vol. 48, no.10, pp. 1246-1252, 2001. [21] Y. Q. Zhang and D. B. Ge, "A unified FDTD approach for electromagnetic analysis of dispersive objects," Prog. Electromagn. Res., vol. 96, pp. 155-172, 2009.

of

[22] O. C. Zienkiewicz, and R. L. Taylor, The Finite Element Method: Its Basis and Fundamentals, Amsterdam: Elsevier. 2005.

ro

[23] N, M. Newmark, "A Method of Computation for Structural Dynamics," J. Een. Mech. Div., vol. 85, pp. 67-94, 1959.

[24] A. Taflove and S. C. Hagness, Computational Electrodynamics-he Finite Difference Time

-p

Domain Method, third edition. London: Artech House. 2005.

[25] C. Salvatore and C. Gaia, "Assessment of the performances of first-and second-order

re

time-domain ABC's for the truncation of finite element grids," Microw. Opt. Techn. Let., vol.

Jo

ur na

lP

38, no.1, pp. 11-16, 2003.

Table 1

Coefficients values of polarizability function

Dispersive media models

Gq

q

Debye

j   0

1 0

G1  j p2 vc

0

G2   j p2 vc

vc

Drude

02 

02  c 2

 c  j 02  c 2

Jo

ur na

lP

re

-p

ro

of

Lorentz