Abstract:
A method for determining electrical and magnetic field effects determines Lanczos matrices by performing a preconditioned conjugate gradient method using a nested multi-grid, vector and scalar potential preconditioner so that Pade via Lanczos frequency expansion may be used to determine the electrical and magnetic field effects over a frequency range without having to perform computationally slow and memory intensive matrix decomposition.

Description:
CROSS REFERENCE TO RELATED APPLICATIONS  
       [0001]     This application claims priority to U.S. provisional application Ser. No. 60/551,294 filed Mar. 8, 2004, which is incorporated herein by reference. 
     
    
     FIELD OF THE INVENTION  
       [0002]     The present invention generally relates to computer implemented methods for solving Maxwell&#39;s equations and in particular, to a computer implemented method for determining electrical and magnetic field effects.  
       BACKGROUND OF THE INVENTION  
       [0003]     Electronic circuit designs are getting smaller and more complex in geometry, while their operating frequencies are getting higher. As a result, electronic circuit designers rely on numerical simulation to predict electrical and magnetic fields of their designs before manufacturing. Currently available simulation tools, however, fail to solve many practical problems due to large computer memory requirements and long computation times.  
         [0004]     For a given electronic circuit design, electrical and magnetic fields may be calculated using Maxwell&#39;s equations: 
 
∇× {right arrow over (E)}=−jω{right arrow over (B)} 
 
∇× {right arrow over (H)}=jω{right arrow over (D)}+{right arrow over (J)} 
 
∇· {right arrow over (D)}=ρ 
 
∇· {right arrow over (B)}= 0   (1) 
 
 where E is the electrical field, H is the magnetic field, D is the electric flux density, B is the magnetic flux density, J is the current density, ω is the angular frequency, ρ is the electric charge density, and j={square root}{square root over (−1)}. 
 
         [0005]     Numerical methods such as the finite element method or the integral equation method can be used to solve the Maxwell equations. For details on the finite element method, see, e.g., Silvester, P. P., and Ferrari, R. L, “Finite Elements for Electrical Engineers”, 3 rd  ed., Cambridge University Press, Chapter 9, 1996. For details on the integral equation method, see, e.g., Sun, D. K., Lee, J. F., and Cendes, Z. J., “ALPS-A New fast Frequency-Sweep Procedure for Microwave Devices,” IEEE Trans. Microwave Theory Tech., pp. 398-402, 2001.  
         [0006]     Using these numerical techniques, the Maxwell equations may be transformed into a set of algebraic simultaneous equations: 
 
([ K   o   ]−ω[K   1   ]−ω   2   [M   o ]){ x}=ω{b}   (2) 
 
 where [K o ], [K 1 ], and [M 0 ] are known matrices, {x} is the unknown electrical field, and {b} is a known electrical excitation. 
 
         [0007]     Computer simulation may then be performed to solve equation (2) for the electrical field, {x}, at all frequencies of interest. Upon obtaining the electrical field, {x}, the corresponding magnetic field can be calculated using Maxwell&#39;s equations (1).  
         [0000]     The Pade Via Lanczos Method  
         [0008]     When analyzing an electronic circuit design, it may be necessary to calculate the electrical and magnetic fields at, for example, about 1,000 different frequencies to cover the application frequency range. One way to do this is to first calculate the electrical field, {x} in equation (2), at just one frequency point and then, expand the solution at this frequency to other frequencies. This frequency expansion technique is known as the “Pade via the Lanczos” method, or alternatively, the Krylov subspace method.  
         [0009]      FIG. 1  illustrates, as an example, various steps performed in the Pade via the Lanczos method. In  101 , Maxwell&#39;s equations (1) are transformed into the set of algebraic simultaneous equations (2). In  102 , an LU decomposition of [K 0 ] is performed so that matrices [L] and [U] are found such that [K 0 ]=[L][U], where [L] is a lower triangular matrix (i.e., zero for all elements above the diagonal) and [U] is an upper triangular matrix (i.e., zero for all elements below the diagonal). For an example of how the LU decomposition may be performed, see, e.g., Press, W. H., et al, “Numerical Recipes,” Cambridge Press, London, pp. 31-38, 1986.  
         [0010]     In  103 , Lanczos Matrices of [K 0 ] −1 [K 1 ] are then determined, where [K 0 ] −1  is the inverse matrix of [K 0 ] and [K 0 ] −1 =[U] −1 [L] −1 . For an example of how the Lanczos Matrices may be determined, see, e.g., Z. Bai, J. Demmel, J. Dongarra, A. Ruhe and H. van der Vorst, editors, “Templates for the solution of Algebraic Eigenvalue Problems: A Practical Guide,” SIAM, Philadelphia, pp. 189-196, 2000.  
         [0011]     In  104 , the Pade via Lanczos frequency expansion is then performed. For examples of how the Pade via Lanczos frequency expansion may be performed, see, e.g., Feldmann, P., and Freund, R. W., “efficient Linear Circuit Analysis by Pade Approximation via the Lanczos Process,” IEEE trans. Computer-Aided Design, Vol. 34, pp. 639-649, 1995; Bracken, J. E., Sun, D. K., and Cendes, Z. J., “S-domain Methods for Simultaneous Time and Frequency Characterization of Electromagnetic Devices,” IEEE Trans. Microwave Theory Tech., pp. 1277-1290, 1998; and Sun, D. K., Lee, J. F., and Cendes, Z. J., “ALPS-A New fast Frequency-Sweep Procedure for Microwave Devices,” IEEE Trans. Microwave Theory Tech., pp. 398-402, 2001.  
         [0012]     The LU decomposition of matrix [K 0 ] performed in  102 , however, is not only very slow but also very memory intensive. As a result, it is generally not practical to simulate a typical electronic circuit design using the Pade via the Lanczos method.  
         [0000]     The NMGAV Method  
         [0013]     Another method for solving Maxwell&#39;s equations is the Nested Multi-Grid, Vector (A) and Scalar (V) Potential Preconditioner (NMGAV) method. Details on this method may be found in Zhu, Yu, “Methodologies and Methods for Fast Finite-Element Analysis of Electromagnetic Devices and Systems”, Ph.D. Thesis, University of Illinois at Urbana-Champaign, Chapter 7, 2002 (the “Zhu reference.  
         [0014]     The NMGAV method uses an iterative procedure to find solutions instead of using a direct matrix solver. In an iterative solution method, the solution is approximated iteratively until the residual (i.e., difference between the exact and approximated solutions) is very small. A good iterative solver can converge sooner than others.  
         [0015]     In the computation of Maxwell&#39;s equations, there may be spurious modes that significantly slow down the convergence of an iterative solver, however. For details and examples of such spurious modes, see, e.g., Dyczij-Edlinger, R., Peng, G., and Lee, J. F., “Stability Conditions for Using TVFEMs to Solve Maxwell Equations in the Frequency Domain,” International Journal of Numerical Modeling: Electronic Networks, Devices and Fields, Vol. 13, Issue 2-3, pp.245-260, 2000.  
         [0016]     The NMGAV method uses the AV preconditioner and the nested multi-grid method to screen out the spurious modes. As a result, a convergent solution can be obtained in less than, for example, 30 iterations, versus thousands of iterations by other iterative methods. As a result of its quick convergence, the NMGAV method can calculate electrical and magnetic solutions at a single frequency two orders of magnitude faster than the corresponding portion of the Pades via Lanczos method described above. In addition to its speed advantage, computer memory requirements for this method are relatively quite small. This is because no matrix inversion or decomposition is needed in this method. Also, because of its low memory requirement, the NMGAV method is practical for solving large problems.  
         [0017]     An example of the NMGAV method follows: 
        1. If the coarsest grid, then solve A 0   h z E   h =r E   h  with an LU direct solver or an iterative solver.     2. else 
            a. pre-smooth A o   h z E   h =r E   h  using forward GS v 1  times with z E =0     b. r v   h =G(r E   h −A 0   h z h )     c. pre-smooth A vv   h z v   h =r v   h  using forward GS v 1  times with z v =0     d. z E   h =z E   h +G T z v   h       e. r E   2h =I h   2h (r h −A 0   h z E   h )     f. z E   2h =NMGAV (r E   2h )     g. z E   h =z E   h +I 2h   h z E   2h       h. r v   h =G(r E   h −A 0   h z h )     i. post-smooth A vv   h z v   h =r v   h  using backward GS v 1  times with z v =0     j. r E   h =(r E   h −A 0   h G T z v   h )     k. post-smooth A 0   h z E   h =r E   h  using backward GS v 1  times with z v =0     l. z E   h =z E   h +G T z v   h      
            3. end if        
 
         [0033]      FIG. 2  illustrates an example of an application of the NMGAV method according to the Zhu reference. In 201, a frequency is initially set to a lower limit frequency of an application frequency range. In  202 , the algebraic simultaneous equations (2) are solved at the current frequency using the NMGAV method. In  203 , a determination is made whether the current frequency is greater than or equal to an upper limit frequency of the application frequency range. If the determination in  203  is YES, then the method is completed and an END condition results. On the other hand, if the determination in  203  is NO, then in  204 , the current frequency is incremented, and the method loops back to  202 , so that the algebraic simultaneous equations (2) may be solved again at the incremented frequency.  
         [0034]     Although the NMGAV method can speed up the solution of the simultaneous equations (2) at a single frequency by, for example, two orders of magnitude with a relatively much smaller computer memory requirement than the corresponding portion of the Pade via Lanczos method, solutions at all the frequencies of interest have to be calculated one by one. Repeating the solution process many times (such as, for example, 1,000 times) for the entire frequency range more than gives back the gain in the single frequency solution speed.  
         [0035]     On the other hand, although the Pade via Lanczos method can obtain solutions at all frequencies by solving the simultaneous equations (2) only once, it requires an LU decomposition of the matrix [K 0 ] in the set of algebraic simultaneous equations (2), which is not only computationally slow but also computer memory intensive, making the Pade via Lanczos method impractical for large problems.  
       OBJECTS AND SUMMARY OF THE INVENTION  
       [0036]     Accordingly, it is an object of one or more aspects of the present invention is to provide a computer implemented method for determining electrical and magnetic fields that is computationally faster than the Pade via Lanczos and NMGAV methods.  
         [0037]     Another object of one or more aspects of the present invention is to provide a computer implemented method for determining electrical and magnetic fields that requires less memory than the Pade via Lanczos method.  
         [0038]     Still another object of one or more aspects of the present invention is to provide a computer implemented method for determining electrical and magnetic fields that provides results substantially the same as or of better accuracy than the Pade via Lanczos and NMGAV methods.  
         [0039]     These and additional objects are accomplished by the various aspects of the present invention, wherein briefly stated, one aspect is a method for determining electrical and magnetic field effects, comprising: determining Lanczos matrices by performing a preconditioned conjugate gradient using a nested multi-grid, vector and scalar potential preconditioner; and performing a Pade via Lanczos frequency expansion to determine electrical and magnetic field effects over a frequency range.  
         [0040]     Another aspect is an apparatus for determining electrical and magnetic field effects, comprising a processor configured to: determine Lanczos matrices by performing a preconditioned conjugate gradient using a nested multi-grid, vector and scalar potential preconditioner; and perform a Pade via Lanczos frequency expansion to determine electrical and magnetic field effects over a frequency range.  
         [0041]     Still another aspect is a computer readable medium encoded with a method executable on a computer for determining electrical and magnetic field effects, comprising: determining Lanczos matrices by performing a preconditioned conjugate gradient using a nested multi-grid, vector and scalar potential preconditioner; and performing a Pade via Lanczos frequency expansion to determine electrical and magnetic field effects over a frequency range.  
         [0042]     Additional objects, features and advantages of the various aspects of the present invention will become apparent from the following description of its preferred embodiment, which description should be taken in conjunction with the accompanying drawings.  
     
    
     BRIEF DESCRIPTION OF THE DRAWINGS  
       [0043]      FIG. 1  illustrates a flow diagram of a prior art Pade via Lanczos method for solving Maxwell&#39;s equations over an application frequency range.  
         [0044]      FIG. 2  illustrates a flow diagram of a prior art NMGAV method for solving Maxwell&#39;s equations over an application frequency range.  
         [0045]      FIG. 3  illustrates a flow diagram of a method for determining electrical and magnetic fields utilizing aspects of the present invention.  
         [0046]      FIG. 4  illustrates a top plan view of a three-pole, low-pass filter interconnect structure.  
         [0047]      FIG. 5  illustrates comparative frequency responses for S-parameters calculated by a point-by-point procedure and S-Parameters calculated by the method of the present invention for the three-pole, low-pass filter interconnect structure of  FIG. 4 .  
         [0048]      FIG. 6  illustrates a computer system useful for practicing the method for determining electrical and magnetic fields utilizing aspects of the present invention.  
     
    
     DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT  
       [0049]      FIG. 3  illustrates, as an example, a computer implemented method for determining electrical and magnetic field effects so as to be useful for simulating integrated circuit, package, and board designs. The method combines certain aspects of the NMGAV and Pade via Lanczos methods in order to achieve superior performance while minimizing memory resource requirements.  
         [0050]     In  301 , the set of algebraic simultaneous equations (2) is formed from Maxwell&#39;s equations (1) using a conventional technique, such as the finite element method, the boundary element method, or the integral equation method.  
         [0051]     In  302 , an initial expansion frequency ω o  is selected and the set of algebraic simultaneous equations (2) are re-written as follows: 
 
([ K   o ]−ω o   [K   1 ]−ω o   2   [M   o ]){ x }−(ω−ω o )([ K   1 ]+2ω o   [M   o ]){ x }−(ω−ω o ) 2   [M   o   ]{x}=ω{b} 
 
or 
 
([ A   o   ]−s[A   1   ]−s   2   [A   2 ]){ x }=(1 +s ){ b   o }  (3) 
 
 where  
             s   =       ω   -     ω   o         ω   o               (   4   )             
 
[ A   o   ]=[M   EE   ]=[K   o   ]−ω   o   [K   1 ]−ω o   2   [M   o ]  (5) 
 
[ A   1 ]=ω o   [K   1 ]+2ω o   2   [M   o ]  (6) 
 
[A 2 ]=ω o   2 [M o ]  (7) 
 
{b o }=ω o {b}  (8) 
 
         [0052]     Preferably, the initial expansion frequency ω o  is approximately in the middle of the application frequency range of interest.  
         [0053]     By using the HZR mass lumping method (as described, for example, in Cook, R. D., “Concept and Applications of Finite Element Analysis”, John Wiley &amp; Sons, pp. 370-373, 1976), the matrix [M o ] can be approximated by a diagonal matrix [m o ]. The HZR mass lumping method is useful in this case to get accurate solutions, especially when prism elements are used (as described, for example, in Sacks, Z. S. and Lee, J. F., “A Finite Element Time Domain Method Using Prism Elements for Microwave Cavities,” IEEE Trans. Electromagnetic Compatibilities, Nov. 1995).  
         [0054]     Solving the set of equations (3) with a Lanczos solver has been found to be very slow. This is because it is desirable to capture the smallest eigenvalues of equations (3) to reach a solution convergence, but the Lanczos method captures far more largest eigenvalues than smallest eigenvalues.  
         [0055]     Therefore, the set of equations (3) are rearranged as follows to accelerate the solution convergence. 
 
(└ A′   o   ┘−s└A′   1   ┘−s   2   └A′   2 ┘){ y }=(1 +s ){ b′   o }  (9) 
 
 where  
             s   =       ω   -     ω   o         ω   o               (   10   )             
 
[ A′   o   ]=[A   2 ] −1    (11) 
 
[ A′   1   ]=[A   2 ] −1   [A   1   ][A   0 ] −1    (12) 
 
[ A′   2   ]=[A   0 ] −1    (13) 
 
{b′ o   }=[A   2 ] −1   {b   o }  (14) 
 
{ x}=[A   0 ] −1   {y}   (15) 
 
         [0056]     By solving equation (9) for electrical fields {y}, the electrical fields {x} which are to be determined for equation (3) may be solved through equation (15).  
         [0057]     In  303 , the Lanczos matrices are constructed by performing a preconditioned conjugate gradient (PCG) method using the NMGAV method for the preconditioner. For details on the PCG method, see, e.g., Barrett, R. et al., “Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods,” SIAM, Philadephia, p. 15, 1994.  
         [0058]     [A′ 0 ]=[A 2 ] −1  is the preconditioner in the following PCG method. The relationship of [A′ 0 ] −1 =[A 2 ] is used in the PCG method.  
           1.   ⁢           ⁢     y   0       =   0     ;           ⁢       r   0     =     b   o   ′       ;           ⁢       z   0     =       A   2     ⁢     r   o         ;           ⁢       p   0     =     z   0       ;       
           q   0     =       z   0         (       r   0     ,     z   0       )       1   /   2           ;       r   0   ″     =       r   0         (       r   0     ,     z   0       )       1   /   2               
        2. For j  32  0, 1, 2, 3, . . . , until convergence Do 
            a. solve A o t j =p j  for t j  with the NMGAV method, t j =NMGAV (p j ) (see below)  
         b   .           ⁢     α   j       =       (       r   j     ,     z   j       )       (       t   j     ,     p   j       )           
    c. x j+1 =x j +α j p j       d. r j+1 =r j −α j t j       e. z j+1 =A 2 r j+1   
         f   .           ⁢     β   j       =       (       r     j   +   1       ,     z     j   +   1         )       (       r   j     ,     z   j       )           
    g. p j+1 =z j+1 +β j p j   
           h   .           ⁢     q     j   +   1         =       z     j   +   1           (       r     j   +   1       ,     z     j   +   1         )       1   /   2           ,       r     j   +   1     ″     ⁢       r     j   +   1           (       r     j   +   1       ,     z     j   +   1         )       1   /   2               
   
            3. Enddo 
 
 where [Q]=(q 0 , q 1 , q 2 , . . . ) is the matrix containing all Lanczos vectors for use later in the frequency expansion, and the matrix [R]=(r′ 0 , r′ 0 , r′ 0 , . . . ) is useful in determining the signs of the Lanczos tri-diagonal matrix. 
       
 
         [0066]     The Lanczos Matrices are then constructed, for example, as follows: 
        1. Form the following tri-diagonal matrix with the PCG constants α k  and β k  obtained in  302 .  
               ⁢         (   16   )     ⁢     
     [     T   m     ]     =     [           (     1     α   0       )             ±       β   0           α   0           0       0       0               ±       β   0           α   0             (       1     α   1       +       β   0       α   0         )             ±       β   1           α   1           0       0           0       .       .       .       0           0       0       .       .           ±       β     m   -   2             α     m   -   2                 0       0       0           ±       β     m   -   2             α     m   -   2               (       1     α     m   -   1         +       β     m   -   2         α     m   -   2           )           ]           
    2. (Optional) Obtain eigenvalues (λ 1 , λ 2 , . . . , λ m ) of [T m ] using the TQLI method [11].     3. (Optimal) Obtain eigenvectors [U] m×m =(u 1 , u 2 , . . . , u m ) by substituting (λ 1 , λ 2 , . . . , λ m ) into [T m ].        
 
         [0070]     Equation (16) is an important contribution and aspect of the present invention that serves to bridge the PCG (using NMGAV) and Pade via Lanczos methods so that the NMGAV method can be employed to quickly determine a solution at an expansion frequency and the Pade via Lanczos method may be employed to perform frequency expansion around the expansion frequency for the application frequency range without repeatedly performing the NMGAV method as described in reference to  FIG. 2 . A proof of equation (16) is given in Appendix A.  
         [0071]     The sign of each of the off diagonal terms in step 1 is determined by the famous Lanczos relationship of [A′ 2 ][Q]=[A 0 ] −1 [Q]=[R][T], where the Lanczos matrices [Q] and [R] are obtained above during the PCG method.  
                     ±       β   0           α   0       ⁡     [     A   0     ]       ⁢     {     r   1     }       =       {     q   0     }     -         1     α   0       ⁡     [     A   0     ]       ⁢     {     r   0     }                 (   18   )                         ±       β   1           α   1       ⁡     [     A   0     ]       ⁢     {     r   2     }       =       {     q   1     }     -         (       1     α   1       +       β   0       α   0         )     ⁡     [     A   0     ]       ⁢     {     r   1     }       -           ±       β   0           α   0       ⁡     [     A   0     ]       ⁢     {     r   0     }           ⁢     
     ⁢     …   ⁢           .     
     ⁢   …   ⁢           .             (   19   )                       ±       β   j           α   j       ⁡     [     A   0     ]       ⁢     {     r     j   +   1       }       =       {     q   j     }     -         (       1     α   j       +       β     j   -   1         α     j   -   1           )     ⁡     [     A   0     ]       ⁢     {     r   j     }       -           ±       β     j   -   1             α     j   -   1         ⁡     [     A   0     ]       ⁢     {     r     j   -   1       }                 (   20   )             
 
         [0072]     In the above calculation for the signs of the off diagonal terms, not all elements of the vectors are required. Only one non-zero element of the vectors, such as the element associated with a port, will be sufficient. Determining the signs of the off diagonal terms is an important step in the entire process. Without doing this, the final solutions may not be right.  
         [0073]     In  304 , a Pade via Lanczos frequency expansion is performed to calculate solutions at frequencies in the application frequency range. Applying the Lanczos congruent transformation to equation (9) and using the following Lanczos relationships: 
 
[ A   o ] −1   [Q]−[R][T   m   ]=[A   2 ] −1   [Q][T   m ]  (21) 
 
[ Q]   T   [A   2 ] −1   [Q]=[I]   (22) 
 
 we have 
 
[ Q]   T ([ A′   0   ]−s[A′   1   ]−s   2   [A′   2 ])[ Q]{z }=(1 +s )[ Q]   T   {b′   0 }  (22) 
 
or 
 
([ T   0   ]−s[T   1   ]−s   2   [T   2 ]){ z }−(1 +s )[ Q]   T   {b′   o }  (23) 
 
 where 
 
{y}=[Q]{z}
 
[T 0 ]=[Q] T [A′ 0 ][Q]=[I]
 
[ T   1   ]=[Q]   T   [A′   1   ][Q]=[Q]   T   [A   2 ] −1   [A   1   ][A   2 ] −1   [Q][T   m ]
 
[T 2 ]=[Q] T [A′ 2 ][Q]=[T m ]  (24) 
 
         [0074]     In the above equations, [T 1 ] can be obtained very quickly because [A 2 ]=ω o   2 [m o ] is a diagonal matrix.  
         [0075]     Therefore,  
                     {   x   }     =       ⁢           [     A   0     ]       -   1       ⁢     {   y   }       =           [     A   0     ]       -   1       ⁡     [   Q   ]       ⁢     {   z   }                     =       ⁢           (     1   +   s     )     ⁡     [     A   0     ]         -   1       ⁡     [   Q   ]                       ⁢             (       [   I   ]     -     s   ⁡     [     T   1     ]       -       s   2     ⁡     [     T   m     ]         )       -   1       ⁡     [   Q   ]       T     ⁢     {     b   o   ′     }                   =       ⁢       (     1   +   s     )     ⁢               (       b   o   ′     ,     b   o   ′       )       1   /   2       ⁡     [     A   2     ]         -   1       ⁡     [   Q   ]       ⁡     [     T   m     ]                         ⁢         (       [   I   ]     -     s   ⁡     [     T   1     ]       -       s   2     ⁡     [     T   m     ]         )       -   1       ⁢     {     e   1     }                     (   25   )             
 
         [0076]     Since the PCG method converges very quickly (e.g., about 50 iterations), dimensions of [I], [T 1 ], and [T m ] are very small (e.g., less than 50). As a result, inversion of the matrix ([I]−s[T 1 ]−s 2 [T m ]) at any frequency s can be computed very quickly with any matrix inversion numerical method. Inversion of [A 2 ]=ω o   2 [m] is trivial, since it is a diagonal matrix.  
         [0077]     One way to ensure that equation (25) is correct is to check the solution at the expansion frequency ω o  (or s=0). At s=0,  
                     {   x   }     =       ⁢               [     A   0     ]       -   1       ⁡     [   Q   ]       ⁡     [   Q   ]       T     ⁢     {     b   o   ′     }                   =       ⁢                     [     A   0     ]       -   1       ⁡     [   Q   ]       ⁡     [   Q   ]       T     ⁡     [     A   2     ]         -   1       ⁢     {     b   o     }       =         [     A   0     ]       -   1       ⁢     {     b   o     }                         ⁢     (               [   Q   ]     ⁡     [   Q   ]       T     ⁡     [     A   2     ]         -   1       =               [   Q   ]     T     ⁡     [     A   2     ]         -   1       ⁡     [   Q   ]       =     [   I   ]         )                   (   26   )             
 
 which is indeed the solution of equation (3) at s=0. 
 
         [0078]      FIGS. 4 and 5  are provided, as an example, to illustrate the effectiveness of the method for determining electrical and magnetic field effects as described in reference to  FIG. 3  compared to the prior art NMGAV method described in reference to  FIG. 2 .  
         [0079]     Referring to  FIG. 4 , a top plan view of a microstrip, three-pole, low-pass filter interconnect structure is illustrated, having the following dimensions and characteristics:  
         [0080]     L 1 =1.651 mm, L 2 =1.143 mm, L 3 =0.635 mm 
        W 1 =0.635 mm, W 2 =1.524 mm, W 3 =3.175 mm     Metal Thickness=0.06 mm     Substrate Thickness=0.635 mm     Substrate Dielectric Constant=9.6        
 
         [0085]     Now referring to  FIG. 5 , comparative frequency responses for S-parameters S 12   511  and S 11   512  calculated by a point-by-point procedure using the NMGAV method described in reference to  FIG. 2 , and S-Parameters S 12   501  and S 11   502  calculated by the method of the present invention described in reference to  FIG. 3  are illustrated for the interconnect structure of  FIG. 4 .  
         [0086]     For both methods, the finite element model used for the interconnect structure consisted of the following: 
        Number of Triangular Prism Elements=16,644     Number of Nodes=9,912     Number of Edges=37,819        
 
         [0090]     As can be readily seen from the figure, S-parameter S 12  and S 11  solutions,  501  and  502 , obtained by using the proposed frequency expansion method are very close to the s-parameter solutions S 12  and S 11 ,  511  and  512 , calculated by the point-by-point NMGAV method (taken at 0.2 GHz increments). Agreement between the two solutions is especially excellent near the expansion point of 10.5 GHz, which happens to be approximately in the middle of the application frequency range of 1 GHz (lower limit frequency) to 20 GHz (upper limit frequency).  
         [0091]     However, the computer time used by the proposed frequency expansion method was found to be considerably less than that of the point-by-point NMGAV method. For example, using a Pentium 4 (3.0 GHz) computer with 2.0 GB of DRAM, the point-by-point NMGAV method took 2,302.7 seconds whereas the proposed frequency expansion method only took 218.5 seconds for a speedup ratio of 10.53.  
         [0092]     As is apparent from the various references cited and incorporated herein by reference, the proposed frequency expansion method as described in reference to  FIG. 3  has many applications, especially in the simulation of integrated circuits and other electronic devices, electronic device packaging, and printed wiring or circuit boards.  
         [0093]      FIG. 6  illustrates, as an example, a computer system  600  useful for performing such simulations and practicing the method for determining electrical and magnetic fields described in reference to  FIG. 3 . The computer system  600  includes a central processing unit (CPU)  601  which may perform the method according to a computer program stored in mass storage  660  and executed using system memory  602  in a conventional manner. A CD-ROM drive  612  (or other computer readable medium drive) may be used to install the program from a CD-ROM disc (or other computer readable medium).  
         [0094]     Although the various aspects of the present invention have been described with respect to a preferred embodiment, it will be understood that the invention is entitled to full protection within the full scope of the appended claims.  
         [0095]     Appendix A—Derivation of the Relationship between Lanczos and PCG Methods.  
         [0096]     We first review the Lanczos method for solving Ax−λMx=0.  
                   1.   ⁢           ⁢     x   0       =   0     ,       r   0     =       M     -   1       ⁢   b       ,       v   0     =   0     ,     γ   =       (       r   0     ,     r   0       )       1   /   2         ,       v   1     =       r   0     γ       ,       η   1     =   0       ⁢     
     ⁢       2.   ⁢           ⁢   For   ⁢           ⁢   j     =   1     ,   2   ,   3   ,   …   ⁢           ,     m   ⁢           ⁢   Do             (   A1   )                 a   .           ⁢     w   j       =         M     -   1       ⁢     Av   j       -       η   j     ⁢     v     j   -   1                   (   A2   )                 b   .           ⁢     δ   j       =       (       w   j     ,     v   j       )     M             (   A3   )                 c   .           ⁢     w     j   +   1         =       w   j     -       δ   j     ⁢     v   j                 (   A4   )                 d   .           ⁢     η     j   +   1         =       (       w     j   +   1       ,     w     j   +   1         )     M     1   /   2               (   A5   )                 e   .           ⁢     v     j   +   1         =       w     j   +   1         η     j   +   1                 (   A6   )             
        3. EndDo  
                 4.   ⁢           [     T   m     ]     =     [           δ   1           η   2         0       …             η   2           δ   2           η   3         …           0                               …           …       …         η   m           δ   m           ]       ,       [     V   m     ]     =     [           v   1           v   2         …         v   m           ]               (   A7   )             
       
 
         [0098]     The residual vector of the approximation solution can be written as 
 
(A−Λ m M)[V m ]∝v m+1    (A8) 
 
         [0099]     In the Lanczos method 
 
 V   m   T   MV   m   =I   m  or ( v   i   , v   j ) M =δ ij    (A9) 
 
 AV   m   =MV   m   T   m +residual  T   m   ≈V   m   T   AV   m    (A10) 
 
         [0100]     We then outline the PCG (preconditioned conjugate gradient) method for Ax=b with M as the preconditioner 
        1. x 0 =0; r o =b o ; z 0 =(M) −1 r 0 ; p o =z o (A 11)     2. For j=0, 1, 2, 3, . . . , until convergence Do  
               a   .           ⁢     α   j       =       (       r   j     ,     z   j       )       (       Ap   j     ,     p   j       )               (   A12   )                 b   .           ⁢     x     j   +   1         =       x   j     +       α   j     ⁢     p   j                 (   A13   )                 c   .           ⁢     r     j   +   1         =       r   j     -       α   j     ⁢     Ap   j                 (   A14   )                 d   .           ⁢     z     j   +   1         =         (   M   )       -   1       ⁢     r     j   +   1                 (   A15   )                 e   .           ⁢     β   j       =       (       r     j   +   1       ,     z     j   +   1         )       (       r   j     ,     z   j       )               (   A16   )                 f   .           ⁢     p     j   +   1         =       z     j   +   1       +       β   j     ⁢     p   j                 (   A17   )             
    3. EndDo        
 
         [0104]     In the above method, r j  is residual of Ax=b, while z j  is residual of the preconditioned system M −1 Ax=M −1 b. From Eq (A8), it can be deduced that z j  in the PCG method is proportional to v j+1 , and thus w j+1 , in the Lanczos method. From Eqs (A2) and (A3) we have  
               δ     j   +   1       =       (       w     j   +   1       ,     v     j   +   1         )     M                         ⁢     =       (         M     -   1       ⁢     Av     j   +   1         ,     v     j   +   1         )     M                           ⁢     =         (         M     -   1       ⁢     Aw     j   +   1         ,     w     j   +   1         )     M         (       w     j   +   1       ,     w     j   +   1         )     M                             ⁢     =         (         M     -   1       ⁢     Az   j       ,     z   j       )     M         (       z   j     ,     z   j       )     M                             ⁢     =       (       Az   j     ,     z   j       )       (       r   j     ,     z   j       )                               (       Az   j     ,     z   j       )     =     (       A   ⁡     (       p   j     -       β     j   -   1       ⁢     p     j   -   1           )       ,     (       p   j     -       β     j   -   1       ⁢     p     j   -   1           )       )                         ⁢     =       (       Ap   j     ,     p   j       )     +       β     j   -   1     2     ⁡     (       Ap     j   -   1       ,     p     j   -   1         )                   ⁢           ⁢           (     p   ⁢           ⁢   are   ⁢           ⁢   A                         ⁢   orthogonal   )                         δ     j   +   1       =         (       Ap   j     ,     p   j       )       (       z   j     ,     r   j       )       +       β     j   -   1     2     ⁢       (       Ap     j   -   1       ,     p     j   -   1         )       (       z     j   -   1       ,     r     j   -   1         )       ⁢       (       z     j   -   1       ,     r     j   -   1         )       (       z   j     ,     r   j       )                               ⁢     =       1     α   j       +       β     j   -   1         α     j   -   1                       
 
         [0105]     From Eq (A2), we have 
 
( w   j+1   , v   j ) M =0=( M   −1   Av   j+1   , v   j ) M −η j+1 ( v   j   , v   j ) M =( M   −1   Av   j+1   , v   j ) M −η j+1  
 
         [0106]     Therefore,  
               η     j   +   1       =       (         M     -   1       ⁢     Av     j   +   1         ,     v   j       )     M                 =         (         M     -   1       ⁢     Aw     j   +   1         ,     w   j       )     M                w   j          M     ⁢            w     j   -   1            M                     =         (         M     -   1       ⁢     Az     j   +   1         ,     z   j       )     M                z   j          M     ⁢            z     j   -   1            M                     =       (       Az   j     ,     z     j   -   1         )                z   j          M     ⁢            z     j   -   1            M                     =       (       A   ⁡     (       p   j     -       β     j   -   1       ⁢     p     j   -   1           )       ,     (       p     j   -   1       -       β     j   -   2       ⁢     p     j   -   2           )       )           (       r   j     ,     z   j       )       1   /   2       ⁢       (       r     j   -   1       ,     z     j   -   1         )       1   /   2                       =       (       Ap     j   -   1       ,     p     j   -   1         )           (       r   j     ,     z   j       )       1   /   2       ⁢       (       r     j   -   1       ,     z     j   -   1         )       1   /   2                       =     ±         β     j   -   1           α     j   -   1