Home > Explicit Time-Domain Finite-Element Method Stabilized for an Arbitrarily Large Time Step

Page 1 |

Explicit Time-Domain Finite-Element Method Stabilized for an Arbitrarily Large Time Step

Qing He, Houle Gan

I. INTRODUCTION

netics can be categorized into two classes. One is the ex- plicit time-domain method; the other being the implicit time-do- main method. In an explicit time-domain method, the field solu- tion at each time step is evaluated from the field solutions at pre- vious time steps; whereas in an implicit time-domain method, the field solution at each time step involves the field solution that is unknown. Explicit methods can avoid solving a matrix equation, while implicit methods generally require a matrix so- lution. Despite its advantage of avoiding a matrix solution, an explicit method requires the time step to be restricted by the smallest space step for ensuring stability. For problems that have fine features relative to working wavelength like on-chip integrated circuits, explicit methods require a large number of time steps to finish one simulation, which is computationally expensive. In an unconditionally stable method, there is no restriction between space step and time step for obtaining a stable solution. It also permits the use of an arbitrarily large time step without becoming unstable. Existing unconditionally

Manuscript received August 10, 2011; revised November 02, 2011; accepted May 30, 2012. Date of publication July 10, 2012; date of current version October 26, 2012. This work was supported in part by a grant from Intel Corporation, in part by a grant from the Office of Naval Research under award N00014-10-1- 0482, and in part by grants from the NSF under Awards 0747578, 0802178, and 1065318. The authors are with the School of Electrical and Computer Engineering, Purdue University, West Lafayette, IN 47907 USA (e-mail: djiao@purdue.edu). Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org. Digital Object Identifier 10.1109/TAP.2012.2207666

stable methods are all implicit methods. In the FDTD based methods, a family of implicit schemes [1]–[9] such as the ADI (alternating direction implicit)-FDTD [1]–[3], the CN (Crank–Nicolson)-FDTD [4], the LOD (locally 1-D)-FDTD [5], the Laguerre-FDTD [6], and the split-step FDTD [8] methods have been developed to achieve unconditional sta- bility. In [7], it is shown that these unconditionally stable implicit FDTD methods can be derived from a general formula- tion based on generalized matrix operator equations pertaining to some classical splitting formulae, from which a variety of other unconditionally stable implicit schemes can further be deduced. Similarly, a group of unconditionally stable time-do- main finite-element methods (TDFEM) such as the Newmark method [10], ADI-FETD [11], the CN-FETD [12], and the backward difference have also been developed. These methods require the solution of a weighted sum of the mass matrix and the stiffness matrix, whereas the explicit TDFEM only requires the solution of the mass matrix, which is either diagonal in nature or can be diagonalized by the orthogonal vector basis functions [13]–[15], hence eliminating the need for solving a matrix equation. In summary, the large time step provided by the existing unconditionally stable scheme is achieved by resorting to implicit time-domain methods that sacrifice com- putational efficiency. Moreover, late-time instability has also been observed from implicit methods. In [16], a new FDTD method using the alternating-direction explicit (ADE) method was developed for efficient electromagnetic field simulation. Though the method is explicit, yet it is not unconditionally stable. In fluid dynamics, techniques have been developed to filter unstable high frequency waves to increase the time step of an explicit scheme [28]. However, the time step is only extended by a limited amount. The filtering technique has also been developed to extend the stability limit of the explicit FDTD method [17] in electromagnetics. As to-date, no explicit methods have been made unconditionally stable. The research question considered in this work is: can an ex- plicit method be made unconditionally stable so that its strength, which lies in avoiding a matrix solution, is retained while elim- inating its shortcoming in time step? The contribution of this work is the successful development of an explicit time-domain method that is unconditionally stable, a capability that previously did not exist. It is stable for any large time step without being restricted by the space step. If the time step is chosen based on accuracy, the proposed method produces not only stable but also accurate results no matter how large the time step is. We have already conducted preliminary research on the proposed method in [18], [19]. In this paper, we complete it from both theoretical and numerical perspectives. The paper is organized as follows. In Section II, we present

0018-926X/$31.00 © 2012 IEEE

Page 2 |

the background of a time-domain finite-element method. In Section III, we describe the proposed theory for making an explicit time-domain method unconditionally stable. In Section IV, we propose an explicit time-domain finite-element method that is unconditionally stable. The linear computational complexity of the explicit method at each time step is preserved by the proposed method. Section V demonstrates the stability, accuracy, and efficiency of the proposed method. It is also shown that the proposed unconditionally stable explicit method outperforms both the conditionally stable explicit method and the unconditionally stable implicit method in computational efficiency as well as stability. Although the proposed method is presented in the framework of a time-domain finite-element method, the essential idea can be applied to other time domain methods, hence, contributing to the removal of one major com- putational bottleneck in time-domain electromagnetic analysis. II. BACKGROUND OF A TIME-DOMAIN FINITE-ELEMENT METHOD Consider the second-order vector wave equation (1) where , , , , and are the electric field, free-space permeability, relative permeability, permittivity, and current density respectively. A time-domain finite-element based solu- tion of (1) and its boundary conditions results in the following system of linear equations [20]: (2) in which is called a mass matrix, is called a stiffness ma- trix, is the unknown field vector, and is a current excitation vector. The and are sparse and symmetric. Typically they have only tens of nonzero elements in each row regardless of the matrix size . These matrices can be assembled in linear time and stored in linear computational resources from their el- emental contributions as follows: (3) where and are the vector basis functions used to expand while denotes a volume integration in each element. Compared to other time-domain methods, a time-domain finite- element method possesses flexibility in both geometrical and material modeling. III. PROPOSED THEORY FOR MAKING AN EXPLICIT TIME-DOMAIN METHOD UNCONDITIONALLY STABLE

Page 3 |

A central-difference based explicit solution to the above yields (9) where . To analyze the stability of (9), we set the excitation to be zero and perform a -transform of (9), obtaining (10) For an explicit time marching like (9) to be stable, of (10) should be bounded by 1. As a result, the following condition must be satisfied: (11) where eigenvalue is the -th entry of the diagonal matrix , and is matrix size. Due to the property of (positive definite) and (semi-pos- itive definite), the eigenvalues of (5), , are nonnegative. The smallest eigenvalue of (5) is zero, which is due to the null space of the stiffness matrix . These zero eigenvalues always exist. The corresponding eigenvector is called a DC mode. Despite its zero eigenvalue, a DC mode can also have a complicated field distribution such as the DC mode of an integrated circuit made of multiple metallic wires immersed in dielectric mate- rials. As for the largest eigenvalue of (5), although theoretically speaking, the resonance frequency of a structure, and hence the eigenvalues of (5), can be infinitely large, the largest eigen- value that can be numerically found is limited by the smallest space resolution. To be specific, the square root of the maximum eigenvalue of (5), (whose unit is angular frequency), is in- versely proportional to the smallest space resolution as the following: (12) where is speed of light. This is because given the smallest space step , the minimum wavelength that can be captured by the space discretization is . The angular frequency cor- responding to such a wavelength is . From (12), it is clear that the smaller the space step, the larger the maximum eigenvalue that can be numerically identified from (5). The meaning of (11) is significant. It demonstrates that when an explicit method becomes unstable, among all the eigenvec- tors (modes) that are contained in the field solution shown in (7), not every mode becomes unstable. Only a subset of the modes, whose eigenvalues are so large that (11) is violated, is unstable. The rest of the modes are stable. For example, the DC modes, the eigenvalues of which are zero, are always stable irre- spective of the choice of time step. Therefore, we conclude that the set of modes that violate (11) in the field solution are the root cause of the instability associated with an explicit time-do- main method when a large time step is used. The in (11) is an angular resonance frequency of the system. When (11) is vi- olated, , where is the frequency cor- responding to . Therefore, given a time step , the unstable modes are also those modes which vary with space at such a high frequency that it cannot be accurately captured by the given time step based on sampling theorem. The number of these unstable modes is also finite. The remaining question is why these unstable modes exist? They exist because of the fine discretization as can be seen from (12). A fine discretization cannot be avoided in problems having fine feature sizes relative to working wavelength. The finer the space discretization, the larger the maximum eigenvalue of (5). Once these unstable modes are inherent in (5), even though the right hand side does not have a projection onto them, the nu- merical roundoff error will have a projection onto them. This can be seen clearly from (9), which is a diagonal system of equa- tions with the -th entry in vector , , representing the coeffi- cient of the -th mode of (5). It is clear that even though is zero, the roundoff error can make nonzero. As a result, in the field expansion shown in (7), the coefficients of the unstable modes would not be zero, and hence the unstable modes exist in the field solution at each time instant. Meanwhile their eigen- values are so large, i.e., these modes vary with space at such a high frequency that they cannot be accurately simulated by the given time step, and hence instability occurs. The above analysis also clearly shows why in cases, where fine features do not exist and hence the space step can be solely determined by accuracy, the time step suggested by the sta- bility criterion has a good correlation with that required by accu- racy. In this case, the frequency corresponding to the maximum eigenvalue of (5) agrees well with the physically important max- imum frequency to be captured. Therefore, , where is the physically important maximum frequency to be captured. As a result, the time step suggested by stability cri- terion (4) has a good correlation with that dictated by accuracy, for sampling an -based system.

Page 4 |

the given time step , denoted by . We expand the field so- lution in the space of stable modes, , as the following: (14) This is equivalent to setting in (13) to be zero at each time instant. Meanwhile, we project the numerical system (2) onto the space of , which is achieved by multiplying both sides of (2) by . We thereby obtain (15) By using the property shown in (6), we have and hence (16) where is a diagonal matrix comprising the eigenvalues cor- responding to stable modes, which are no greater than . The resultant explicit marching of (16) is stable for the given time step regardless of how large the time step is because (11) is always satisfied. It is worth mentioning that by removing unstable modes from the field solution

Page 5 |

that have a large weight are important, whereas the modes, whose eigenvalues are so far away from the working frequency, can be truncated based on prescribed accuracy. Thus, (21) can be computed as (22) with controlled accuracy , where is composed of eigenmodes whose relative weights in the field solution are greater than for the given frequency or the given spectrum. The frequency corresponding to the maximum eigenvalue of these modes represents the maximum frequency of space variation that is physically important in a system response, which is . It is clear that by removing beyond- modes, the accuracy of the field solution is not affected. This is true for the field solution at any point in the computational domain, no matter how far or close the point is from the source since the weights of beyond- modes are negligible due to the large gap between their eigenvalues and working frequency square . In addition, as can be seen from (12), if the space step is chosen based on the accuracy requirement, i.e., capturing a space variation up to frequency, then the maximum eigenvalue of (5) does not go beyond that associated with . Then, no beyond- modes exist in the numerical system. Once beyond- modes exist, that is because the space step is finer than that required by accuracy for discretizing an -based system. Therefore, they are not necessary for the required accuracy, and hence can be eliminated without affecting accuracy. IV. PROPOSED EXPLICIT TIME DOMAIN FINITE ELEMENT METHOD THAT IS UNCONDITIONALLY STABLE The proposed explicit time-domain finite-element method that is stable for any large time step has two steps. The first step is a preprocessing step for building a complete and also accurate space that spans all the stable modes for a given time step . The second step is to perform a march-on in time with the given time step in an explicit time-domain method without violating stability. Both steps retain the strength of an explicit method in avoiding a matrix solution, and hence achieving sta- bility for any time step without sacrificing the linear (optimal) complexity of an explicit method based simulation.

Page 6 |

solving the eigenvalue problem shown in (5), we only need to solve a reduced eigenvalue problem as (27) Denoting the union of the eigenvectors of (27) by , and the eigenvalue matrix by . The solution of (25) can be expanded in the space of . We therefore have (28) where contains the time-variant weights of the eigenvectors , which is similar to in (7) but with a significantly reduced size . Since is symmetric positive definite and is sym- metric, we have [25] (29) Substituting (28) into (25), multiplying both sides of (25) by , and using the property (29), we obtain (30) where can be solved via a central-difference based scheme like (9).

Page 7 |

we select those modes whose eigenvalues are no greater than to form . Thus where are the physically important eigenvectors of (27) whose eigenvalues are no greater than . When is chosen based on accuracy requirement, is equal to , i.e., all the physically important modes will be included in . Hence, the simulation is not only stable but also accurate. It is also worth mentioning that we store by separately storing and , which has a linear cost. We do not need to multiply them together since their direct product is not required in the computation, which will be seen very clearly from Section IV-B.

Page 8 |

Fig. 1. Illustration of the unconditional stability of the proposed method.

mode whose eigenvalue is zero from the reduced eigenvalue so- lution (27), we allow for the use of any large time step without making the simulation unstable. The same applies to other ex- amples if one would like to use a time step that is infinitely large to examine stability. In addition, in this example, since the time step of 0.0001 s, 0.001 s, and 0.01 s also satisfies accuracy re- quirements for the given input source, the results shown in Fig. 1 for these three time steps are not only stable but also accurate, which reveal that the voltages at the near and far end of the par- allel plate operating in the input spectrum are nothing but the in- tegration of the input current divided by a capacitive parameter of the structure. With 0.1 s time step, the proposed method still generates a stable result. But the result is not accurate because the time step of 0.1 s is greater than that required by accuracy.

Fig. 2. Voltage waveforms of an on-chip interconnect. Fig. 3. Simulation of a millimeter-scale waveguide with thin films. (a) Struc- ture. (b). Voltage waveforms.

need 14,800 steps in preprocessing and the cost after prepro- cessing is negligible. Thus, the speedup is 330. Fig. 2 shows an excellent agreement between the proposed method and the conventional central-difference based method.

Page 9 |

To simulate this example, the conventional explicit scheme re- quires a time step of s to maintain stability. Because of this small time step, over ten thousand steps were needed to finish the simulation. With the proposed explicit method, we were able to use a time step of s solely determined by accuracy (because ) to generate accurate and stable results within 1,000 steps. We collected the field solution every 50 steps during the preprocessing. The solver automatically simulated for 3,000 time steps in preprocessing from which 60 vectors were sampled, among which 47 orthog- onal vectors were constructed to form . Using the space , the original large eigenvalue problem was transformed to a small eigenvalue problem of size 47, from which 33 physically important modes were identified with the first two eigenvalues found to be and rad s, with error being less than 1%. Since the conventional explicit method needs steps and the proposed method requires only 3,000 steps in preprocessing, the speedup of the proposed method is 6. In this millimeter-scale example, because the space discretization is not significantly smaller than that required by accuracy, the gap between the time step allowed by accuracy and that permitted by the stability criterion is not that large. Therefore, the speedup of the proposed scheme is not as signifi- cant as that observed in previous example that involves a space resolution much smaller than that required by accuracy due to the existence of fine features relative to working wavelength. In Fig. 3(b), we plot the voltages sampled at the near and far ends of the waveguide simulated by the proposed method in comparison with those generated by the conventional explicit method as well as implicit Newmark-based unconditionally stable scheme. Excellent agreement is observed.

Fig. 4. Simulation of an on-chip bus. (a) Illustration of the structure and mate- rial. (b) Total CPU time comparison between 3 methods. Fig. 5. Comparison of two methods in late-time stability.

explicit method and the explicit method in [15] both employ orthogonal prism vector bases and hence exhibit linear com- plexity. However, the proposed method was able to use 3000 steps to finish the entire simulation with a time step of s whereas the conditionally stable explicit method in [15] required 300,000 steps to finish the simulation. As for the im- plicit Newmark method, although it permitted the same large time step as that used by the proposed method, it failed to fac- torize the matrix for large problem sizes.

Page 10 |

modes only. In contrast, we do observe late-time instability of implicit methods in some of the examples we simulated. For example, when simulating the millimeter-scale waveguide ex- ample described in Section V-B-2), we observed late-time in- stability from the Newmark-based implicit scheme, as shown in Fig. 5, whereas the proposed method is stable. VI. CONCLUSION In this work, we propose an unconditionally stable explicit time-domain finite-element method that removes the depen- dence of the time step on the space step. It is stable for an arbitrarily large time step. If the time step is chosen based on accuracy, the proposed method is not only stable but also accurate no matter how large the time step is. The method is applicable to general 3-D problems with inhomogeneous mate- rials and irregular structures. Different from previous methods for achieving such a capability, in which the source that is responsible for the instability still remains in the numerical system, and hence one has to resort to implicit methods to bound the error amplification factor by 1 to ensure stability, the proposed method completely eliminates the root cause that is responsible for the instability associated with an explicit time-domain method. As a result, an explicit method can also be made unconditionally stable, allowing the time step to be solely determined by accuracy requirements. The complete removal of the unstable modes for a given time step is achieved by strictly expanding the field solution in the space of stable modes, and also projecting the numerical system onto the space of the stable modes. By doing so, the original system matrix that spans both stable and unstable modes is transformed to a system of stable modes only. The stable modes are found in the preprocessing step by a time-domain solution based fast eigenvalue solution of , with . This also allows the explicit marching to be efficiently performed in a reduced system of . In the proposed method, the strength of an explicit method in avoiding computationally intensive matrix solutions is re- tained, while its shortcoming of requiring a small time step is overcome for problems having fine features relative to working wavelength. Numerical experiments have demonstrated that the proposed unconditionally stable explicit method outperforms both the conditionally stable explicit method and the uncondi- tionally stable implicit method in efficiency. The method is also shown to be stable at late time, while late-time instability is ob- served in unconditionally stable implicit methods. Although the proposed method is presented in a time-domain finite-element method, the essential idea can be applied to other time domain methods. APPENDIX Here, we prove that the eigenvalues of the large system (5) can be found from the reduced system (27) as long as the space used to reduce (5) to (27) contains the information of the eigenvectors corresponding to these eigenvalues. Consider an eigenpair of (5). It satisfies (A-1) If the space contains the information of , the can be expanded in the space as the following: (A-2) where is a coefficient vector. Substituting (A-2) into (A-1) and testing both sides of (A-1) by , we obtain (A-3) From (26), the above can be further written as (A-4) As a result, the eigenpair is the solution of (27). There- fore, the eigenvalues that satisfy (5) also satisfy (27) as long as the used to reduce (5) to (27) contains the information of the eigenvectors corresponding to these eigenvalues. In ad- dition, by front multiplying the eigenvector obtained from (27) corresponding to by , one can obtain the eigenvector of (5) as can be seen from (A-2). REFERENCES

[1] T. Namiki, ��A new FDTD algorithm based on alternating-direction im- plicit method,��

Page 11 |

[16] S. Aono, M. Unno, and H. Asai, ��Alternating-direction explicit FDTD method for three-dimensional full-wave simulation,�� in

- 2-Mode Concepts in Social Network Analysis
- WP Admin Protection We congratulate you on installation of our plug-in that will not only help you to protect your website but w
- INTERIOR-POINT ALGORITHMS FOR CONVEX OPTIMIZATION BASED ON PRIMAL-DUAL METRICS
- W001-0045 Critical Half.pmd
- 62 1 Ministério da Fazenda
- The Coupling of Alternative Splicing and Nonsense‑Mediated mRNA Decay
- Automating Cognitive Model Improvement by A*Search and Logistic Regression
- LE PEP
- Gabriel Iglesias's official app is not fat... It's fluffy in a good way.
- Untitled
- ыWhat is survival analysis? ыWhen do we use survival analysis? ыOne-sample survival analysis ыTwo-sample survival analysis
- +1/2 -1/2 low energy spin state oriented with (parallel to) magnetic field E Two Dimensional NMR FT (t1) 1D “pulse and
- WLWLSA hopes that all of our members are having a wonderful summer
- FACUNDO CABRAL
- Historical Roots of Police Behavior: Chicago, 1890-1925 Mark H. Haller Law & Society Review, Vol. 10, No. 2, Essays in Honor
- DAFTAR RIWAYAT HIDUP
- forestplot: Advanced Forest Plot Using 'grid' Graphics
- DEDICATORIA A mi madre con mucho amor y cariño le dedico todo mi esfuerzo y trabajo puesto para la realización de esta tesis. Verónica P
- Sample Paper
- MANAGEMENT of INFORMATION SECURITY Third Edition
- Research in Media Effects

- Wen Xin Diao Long
- philosophy basis
- literary aesthetics thought
- the theory of literature writing
- text aesthetics
- non profit organizations
- public services
- participation
- OTN protection
- TMF
- Analysis of Protection Switching
- Ring Protection Group
- Ecological Civilization
- Environmental Education
- School of Ecological Civilization Education
- Ideological and Political Education
- Private Banking
- Risk Management
- Mode
- Maximum Likelihood
- VHF Radar

All Rights Reserved Powered by Free Document Search and Download

Copyright © 2011This site does not host pdf,doc,ppt,xls,rtf,txt files all document are the property of their respective owners. complaint#nuokui.com