Apparatus and method for analyzing circuits using reduced-order modeling of large linear subcircuits

Abstract

A method and apparatus for generating and analyzing reduced-order models of linear circuits. The method and apparatus generates a characteristic tridiagonal matrix using a look-ahead Lanczos procedure and then checks the matrix for stability and passivity. If the matrix is not stable and passive, the matrix is modified via a partial Pade via Lanczos approximation algorithm and then rechecked for passivity and stability. The process is iterated until a passive and stable model is generated. When a satisfactory matrix is generated, a Pade approximant of a transfer function of the circuit is computed and the frequency response of this circuit is generated and displayed.

Description

FIELD OF THE INVENTION

The invention pertains to a method and apparatus for analyzing electrical circuits. More particularly, the invention pertains to a method and apparatus for analyzing linear subcircuits of a larger circuit by means of a Lanczos-type procedure to determine a scalar Pade approximant of the transfer function of the linear subcircuit.

BACKGROUND OF THE INVENTION

Electrical circuits typically are comprised of linear, passive circuit elements, such as resistors, inductors, and capacitors as well as non-linear and active circuit elements, such as transistors and power sources. In the process of designing circuits, particularly large-scale integrated circuits, it is normal practice to mathematically model the circuit to simulate performance. Particularly, the outputs of the circuit are modeled as a function of the inputs of the circuit. The mathematical model is used to determine various response characteristics of the circuit.

In our prior U.S. patent application Ser. No. 08/269,230, filed Jun. 30, 1994, now U.S. Pat. No. 5,537,329, Ser. No. 08/489,270, filed Jun. 9, 1995, now U.S. Pat. No. 5,689,685, and Ser. No. 08/904,233, filed Jul. 31, 1997, all of which are incorporated herein by reference, we disclosed new methods and apparatus for modeling linear subcircuits within a larger electrical circuit.

An electrical circuit, such as an integrated circuit, frequently contains large portions thereof comprised entirely of linear circuit elements. For instance, a wire or other connector usually can be modeled as a network of resistances and capacitances. Accordingly, the entire interconnect structure between a last non-linear circuit element in an integrated circuit and the end of a package pin can be mathematically modeled as a linear subcircuit.

The modeling of the larger overall circuit can be extremely complex in large scale integrated (LSI) circuits and very large scale integrated (VLSI) circuits. Our aforementioned patent applications disclose techniques for modeling the linear subcircuit portions within the larger circuit by a greatly reduced model compared to the equations describing the linear subcircuit.

Referring to our aforementioned patent applications, in U.S. patent application Ser. No. 08/269,230, we disclosed a Pade via Lanczos method for determining a scalar Pade approximant of a transfer function of a linear circuit for use in determining, with substantial accuracy, the frequency response behavior of the linear circuit.

In U.S. patent application Ser. No. 08/489,270, we disclosed an extension to the previous application in which the Pade via Lanczos method was extended to matrix transfer functions. More specifically, U.S. patent application Ser. No. 08/489,270 discloses a numerically stable method to determine a matrix Pade approximation of the frequency response of a linear circuit or subcircuit with little numerical degradation, and providing an acceptable computational cost per order of approximation. The method reduces the very large matrices used to represent the p by p matrix transfer function to a pair of much smaller matrices so that the resulting approximate matrix transfer function has basically the same characteristics as the original matrix transfer function. Using the method and apparatus disclosed therein, the original matrix can be reduced to any particular size desired. The accuracy of the approximation depends on the order of reduction in size of the matrices. Generally, the matrices can be reduced to a very large extent without significant decrease in accuracy.

In U.S. patent application Ser. No. 08/904,223, we disclosed a specialized adaptation of the method and apparatus of U.S. patent application Ser. No. 08/489,270, discussed above, for use in connection with passive, linear circuits or subcircuits in which the matrices are guaranteed to be symmetric. That method can reduce circuit simulation processing load by as much as one half for passive linear circuits. Further, when the additional condition is met that the circuit contains only two types of linear circuit elements, the matrices are guaranteed to be not only symmetric, but also positive definite, which guarantees stability and passivity of the approximation.

In the method and apparatus of application Ser. No. 08/904,223, however, linear circuits that contain R, L and C circuit elements as well as linear circuits with different input and output vectors are not guaranteed to result in a symmetric Lanczos matrix. Therefore, they are not guaranteed to result in a stable and passive Pade approximate transfer function. However, if the original circuit is stable and possibly passive, then, in some applications, it is crucial that the reduced-order models also have these properties.

Therefore, it is an object of the present invention to provide an improved method and apparatus for generating a Pade approximant of a transfer function of a linear circuit.

It is a further object of the present invention to provide an improved method and apparatus for generating a partial Pade approximant of a transfer function of a linear circuit that is passive and stable.

It is another object of the present invention to provide an improved method and apparatus for modifying a Lanczos matrix corresponding to a Pade approximant of a transfer function of a linear circuit that is not passive and stable so as to make it passive and stable with as minimal a sacrifice in the accuracy of the approximant as possible.

SUMMARY OF THE INVENTION

The present invention is a method and apparatus for analyzing linear circuits using reduced-order mathematical approximations thereof generated by a partial Pade via Lanczos type procedure. It includes a processing unit having a processor, associated memory, and stored programs. The memory stores circuit characteristic data representing the circuit. The stored programs include a matrix reduction program for causing the processor to process the circuit characteristic data and to generate block tri-diagonal matrix data representing a block tri-diagonal matrix representative of the circuit being modeled. Specifically, the circuit characteristic data representing the circuit is input through an input device and stored in the memory. A characteristic matrix is generated from the circuit characteristic data and that matrix is iteratively reduced to the aforementioned block tri-diagonal matrix by a look-ahead Lanczos procedure. This matrix is herein termed the original Lanczos matrix.

The stored programs then check the tri-diagonal matrix to determine if the associated reduced-order model is stable and passive. In order to determine stability, the matrix is checked to determine if there are any poles of the rational transfer function of the circuit in the right half of the complex plane. Passivity is checked by determining if there are any zeros of the rational function in the right half of the complex plane and whether additional conditions are met.

If the matrix is passive and stable, the programs generate frequency response data for the circuit from the matrix. However, if the matrix is not stable and passive, the stored programs determine alternate locations for these poles and/or zeros and compute a corresponding partial Pade approximation based thereupon. Particularly, a vector x.sub.n is generated that is added to the last column in the Lanczos matrix. The partial Pade approximation alters the location of all poles and zeros, not just the poles and zeros in the right half of the complex plane. The new resultant matrix is then checked again for stability and passivity using eigenvalue based criteria. If the new matrix is not passive and stable, the process of altering the matrix is repeated until passivity and stability are achieved.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a block diagram of the components of an exemplary apparatus in accordance with the present invention.

FIG. 2 is a flow chart illustrating steps for generating a Pade approximant transfer function of a linear circuit in accordance with the present invention.

FIG. 3 is a flow chart illustrating in more detail step 95 of FIG. 2.

FIG. 4 is a diagram of the actual frequency response and the corresponding original Pade approximant frequency response of a first exemplary circuit.

FIG. 5 is a diagram of the actual frequency response and the corresponding partial Pade approximant frequency response in accordance with the present invention of the first exemplary circuit.

FIG. 6 is a diagram illustrating the poles of the first exemplary circuit both as determined based upon the original Pade approximant and based upon the partial Pade approximant in accordance with the present invention.

FIG. 7 is a diagram of the actual frequency response and the corresponding original Pade approximant frequency response of a second exemplary circuit.

FIG. 8 is a diagram of the actual frequency response and the corresponding partial Pade approximant frequency response in accordance with the present invention of the second exemplary circuit.

DETAILED DESCRIPTION OF THE INVENTION

I. Introduction

Referring to the drawings (in which like reference numerals identify similar or identical elements) and particularly FIG. 1, the present disclosure describes an apparatus and method implementing a circuit analyzer 10, which includes a processing unit 12 having a processor 14, memory 16, and stored programs 18, including a transfer function program 20; an input device 22; and an output device 24.

The present invention is a further improvement on our work disclosed in the aforementioned patent application Ser. Nos. 08/269,230 (now U.S. Pat. No. 5,537,329), 08/489,270, now U.S. Pat. No. 5,689,685, and 08/904,233. A reduced order tri-diagonal Lanczos matrix (hereinafter termed the "original" Lanczos matrix) is generated to represent the transfer function of a linear circuit generally in accordance with the method and apparatus disclosed in our application Ser. No. 08/904,233. That matrix is then checked to determine if it will generate a stable and passive Pade approximant transfer function. If so, the transfer function is generated and a frequency response of the transfer function is generated and displayed to the user. However, if the original Lanczos matrix is not stable, the matrix is modified by a partial Pade approximation algorithm in which a vector, x, is generated and then added to the last column of the original Lanczos matrix to generate a modified tri-diagonal Lanczos matrix that is, hopefully, stable and passive and that has a corresponding transfer function that is very close to the transfer function of the original Lanczos matrix. The new Lanczos matrix is then checked to assure that it is stable and passive and, if so, the transfer function corresponding to that matrix is generated and the frequency response of the transfer function is generated and displayed to the user. If this matrix also is not stable and passive, the process is repeated until a stable and passive partial Pade approximant transfer function can be generated.

In an exemplary embodiment, the processing unit 12 is a SPARC workstation, available from SUN MICROSYSTEMS, INC., having about 32 MB associated RAM memory and a hard or fixed drive memory 15. The processor 14 operates using the UNIX operating system to run application software 18 implementing the disclosed circuit analyzer 10 and method. It should be understood that the above described processor is merely exemplary and that any computer or processing system having sufficient processing capability and memory may be employed.

A description of the circuit to be analyzed is provided from a circuit characteristic data source 26 through an input device 22 to the processor 14. The data source 26 may be any source for providing information to a computer, such as a floppy disk, CD-ROM or DAT containing a circuit characteristic file. It also may be a human operator that enters the data manually. The input device 22 may be a keyboard, a floppy disk drive, a CD-ROM drive or another data reading device.

The processor 14 receives the circuit characteristic data from the circuit characteristic data source 26 through the input device 22. The received circuit characteristic data is stored in memory 16 for processing in accordance with the present invention to generate a frequency response signal representing the approximated or predicted frequency response of the circuit. The generated frequency response data of a linear subcircuit may be a sub-goal in the overall analysis of a larger, non-linear circuit. Alternately, the circuit analyzed may be an independent linear circuit.

The generated frequency response data may be converted by the processor 14 into a graphic representation for plotting the approximated frequency response in volts, amperes, or other units, as a function of frequency. The graphic representation is sent to an output device 24 such as a display for displaying the predicted frequency response of the circuit under analysis. Alternately, the output device 24 may include specialized graphics programs to convert the generated frequency response signal to a displayed graphic.

In even further possible embodiments, the generated frequency response data may include determined poles, zeroes, and/or residues listed in a file for output as columns or tables of text by the output device 24. The output device 24 may be a CRT display or a hard copy printer.

The circuit analyzer 10 performs the application programs and subroutines described hereinbelow which are implemented from compiled source code in the C++ and/or FORTRAN programming languages. The application programs and subroutines are described herein below in conjunction with FIGS. 2 and 3.

Referring to FIG. 2, the circuit analyzer 10 includes software for generating a frequency response signal of a circuit. In step 40, the program is started. In step 45, the circuit characteristic data is received and stored in memory 16 for use with the program.

In step 50, the program processes the circuit characteristic data to generate a characteristic matrix and a characteristic vector that mathematically describe the circuit.

In step 55, the matrix and vector are operated upon to reduce the characteristic matrix using the matrix reduction program of the stored programs 22 to produce a block tri-diagonal matrix T.sub.n using the look-ahead Lanczos procedure. In step 60, the program computes the eigendecompositions of the Lanczos vector T.sub.n and obtains the approximate poles and zeros of the transfer function corresponding to matrix T.sub.n.

In step 65, the program generates the matrix T.sub.n and the corresponding poles and zeros. In step 70, the program checks the matrix T.sub.n, including the locations of it poles and zeros, to determine if the transfer function corresponding to the original Lanczos matrix T.sub.n will be passive (if the original circuit is passive) and stable. If the original circuit is not passive, then it is not necessary that the transfer function be passive, rather, it only need be stable. If so, processing flows to step 75, in which the processor calculates residues from the eigendecompositions to generate a partial Pade approximant of a transfer function of the circuit.

Finally, in step 80, the partial Pade approximant is applied to generate a graphic representation of the circuit response predicted by the approximant. The representation may be a graph of the signal (e.g., voltage or current magnitude) at a particular terminal of the circuit plotted against frequency, such as illustrated in FIGS. 4, 5, 7, and 9, to be discussed in further detail below.

If, in step 70, the matrix T.sub.n is found not to be both stable and passive (if the original circuit was passive), processing flows to step 85. In step 85, the processor selects partial or complete prescribed poles and/or zeros of T.sub.n in accordance with the present invention as described more fully below. In step 90, a vector x is computed using the partial Pade via Lanczos procedure. In step 95, the vector x is added to the last column of the matrix T.sub.n to create a matrix T.sub.n.

In step 100, the processor computes the eigendecompositions of T.sub.n and obtains the approximate poles and zeros of the modified matrix T.sub.n. At this point, processing flows back to step 60 where it is determined if the modified matrix T.sub.n is stable and passive. If so, processing flows through to step 75 through 80 as previously described. If not, steps 85 through 100 are repeated to generate a further T.sub.n vector.

The details of step 55 can be found in aforementioned U.S. patent application Ser. No. 08/489,270 and is also reviewed and summarized below in the section entitled Pade Approximation and Krylov Subspaces.

The details of step 90 are presented in more detail in FIG. 3.

Referring now to FIG. 3, in step 300, the row vectors .phi.(.theta..sub.j) for j=1,2, . . . k are computed. In step 302, the processor computes the row vectors .psi.(.tau..sub.j) for j=1,2, . . . ,l. In step 304, the processor calculates the last column R.sub.k of the k by n matrix T.sub.n. In step 306, the processor calculates the last column S.sub.l of the l by (n-1) matrix T'.sub.n.

In step 308, the matrix equation: ##EQU1## is solved for x.sub.m or the equation: ##EQU2## is solved for x.

Note that X(:,j) and X(:,i:j) denote the j-th column and the i-th to the j-th columns of the matrix X, respectively.

In step 310, the value X.sub.n or x is returned and provided for step 85 of FIG. 2.

The steps in FIG. 3 are described in more detail in Section IV.

Section II below is a brief review of reduced-order modeling techniques based on Krylov subspaces. That section describes the details of steps 50, 55 and 60 in FIG. 2. For an even more detailed discussion of Krylov subspace-based reduced-order modeling and steps 50, 55 and 60, reference should be made to our aforementioned prior U.S. patent application Ser. No. 08/489,270.

II. Krylov's Subspaced-based Reduced-order Modeling

A. Circuit Transfer Functions

A lumped, linear(ized), time-invariant circuit can be described by the following system of first-order algebraic differential equations:

C x(t)=-Gx(t)+bu(t), y(t)=c.sup.T x(t). (Eq. 1)

Here, C and G are N.times.N matrices and represent the contributions of memory elements (such as capacitors and inductors) and memoryless elements (such as resistors), respectively, and x(t) is a N-dimensional state vector representing the circuit variables. The functions u(t) and y(t) represent the excitation and the output of interest, respectively. By applying the Laplace transform to the system (1), assuming zero initial condition, we obtain

sCX(s)=-GX(s)+bU(s), Y(s)=c.sup.T X(s)

Here X(s), U(s) and Y(s) denote the Laplace transforms of x(t), u(t) and y(t), respectively. It follows that the input U(s) and output Y(s) are related by the transfer function

H(s)=c.sup.T (G+sC).sup.-1 b (Eq. 2) Let s.sub.0 .epsilon. be a real expansion points such that the matrix G+s.sub.0 C is nonsingular. Using any factorization of G+s.sub.0 C of the form G+s.sub.0 C=L.multidot.U.sup.-1,L, U.epsilon..sup.N.times.N,(Eq. 3) we now rewrite the transfer function H(s) as follows: H(s)=1.sup.T (I-(s-s.sub.0)A).sup.-1 r (Eq. 4) where A:=-L.sup.-1 CU.sup.-1, r:=L.sup.-1 b 1:=U.sup.-T c.

For example, in (Eq. 3), we can choose an LU factorization of G+s.sub.0 C, or even the trivial factorization L=G+s.sub.0 C and U=I.

By expanding the transfer function (Eq. 4) about s.sub.0 and setting .sigma.=s-s.sub.0, we obtain the representation ##EQU3##

B. Pade Approximation and Krylov Subspaces

Note that H is a rational function. In general, H is of type (N-1,N), i.e., its numerator and denominator polynomials are of degree N-1 and N, respectively. This suggests basing reduced-order modeling on approximate transfer functions H.sub.n of the form ##EQU4## where .psi..sub.n-1 and .phi..sub.b are polynomials of degree at most n-1 and n, respectively, and n<<N. The function (Eq. 6) is called an n-th Pade approximant if the polynomials .psi..sub.n-1 and .phi..sub.n are chosen such that the expansion of H.sub.n about s.sub.0 agrees with (Eq. 5) in its first 2n terms, i.e.,

H(s.sub.0 +.sigma.)=H.sub.n (s.sub.0 +.sigma.)+O(.sigma..sup.2n),(Eq. 7)

Fortunately, the explicit use of moments can be avoided by computing H.sub.n by means of suitable basis vectors for the Krylov subspaces

K.sub.n (A,r)=span{r,A,r,A.sup.2 r, . . . ,A.sup.N-1 r} and K.sub.n (A.sup.T,l)=span{1,A.sup.T1,.sup.2 1, . . . ,(A.sup.T).sup.n-1 l}(Eq. 8) More precisely, one generates two sequences of vectors {j.sub.h }.sub.j=1.sup.n and {w.sub.j }.sub.j=1.sup.n (Eq. 9) such that span{v.sub.1, v.sub.2, . . . ,v.sub.n }=K.sub.n (A, r) and span{w.sub.1, w.sub.2, . . . ,w.sub.n }=K-n(A.sup.T, l). (Eq. 10) Note that m.sub.2j =l.sup.T A.sup.2j r=((A.sup.T).sup.j l).sup.T (A.sup.j r) and m.sub.2j=1 =l.sup.T A.sup.2j+1 r=((A.sup.T).sup.j l).sup.T A(A.sup.j r).(Eq. 11) By (Eq. 11) and (Eq. 10), the information contained in the moments m.sub.0, m.sub.1, . . . ,m.sub.2n-1 is also contained in the quantities w.sub.j.sup.T v.sub.j and w.sub.j.sup.T Av.sub.j, j=1,2, . . . ,n.(Eq. 12) Hence, the point of using the basis vectors (Eq. 10) is that they allow one to obtain all necessary "moment information" for constructing H.sub.n from (Eq. 12), without ever explicitly forming the moments.

C. The Lanczos-Based Approach

The Lanczos process constructs the basis vectors (Eq. 9) so that they are biorthogonal, i.e.,

w.sub.j.sup.T v.sub.k =0 for all j.noteq.k, j,k=1,2, . . . ,n.(Eq. 13) This can be done by means of simple three-term recurrences. The coefficients of these recurrences form a tri-diagonal matrix T.sub.n from which one can directly obtain the n-th Pade approximant H.sub.n. More precisely, we have H.sub.n (s.sub.0 +.sigma.)=(l.sup.T r)e.sub.1.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.1 (Eq. 14) The use of the Lanczos process to compute H.sub.n can also be expressed in terms of the vectors (Eq. 9). Let v.sub.n =.left brkt-bot.v.sub.1 v.sub.2 . . . v.sub.n .right brkt-bot. and W.sub.n =.left brkt-bot.w.sub.1 w.sub.2 . . . w.sub.n .right brkt-bot.(Eq. 15) denote the N.times.n matrices whose columns are just the vectors (Eq. 9). Then, the matrix T.sub.n is given by T.sub.n =(W.sub.n.sup.T V.sub.n).sup.-1 W.sub.n.sup.T A V.sub.n.(Eq. 16) Note that, by (Eq. 13) and (Eq. 15), the matrix W.sub.n.sup.T V.sub.n in (16) is diagonal. III. Partial Pade Approximation

Turning now to the key of the present invention, let H.sub.n (s) be an (n-1, n) Pade approximation to the original transfer function H(s), matching the maximal number of moments: 2n. The partial Pade approximation, proposed in C. Brezinski, Partial Pade Approximation, J. Approx. Theory, vol. 54, pp. 210-233, 1988 (hereinafter Brezinski), is a function of the same form, ##EQU5## where some of the poles and zeros are prescribed. As in the classical Pade approximation, the rational function should maximize the number of matched coefficients in a Taylor expansion about s.sub.0. Here the numerator polynomial .psi..sub.n-1 is of degree n-1 with l prescribed roots and denominator polynomial .phi..sup.n degree n with k prescribed roots. The remaining degrees of freedom, 2n-l-k, could be used to ensure the matching of an equal number of Taylor coefficients, i.e.,

H(s.sub.0 +.sigma.)=H.sub.n (s.sub.0 +.sigma.)+O(.sigma..sup.2n-l-k)(Eq. 24) In Brezinski, the Taylor expansion of the function to be approximated is explicitly applied to construct a partial Pade approximation. The more restricted problem of prescribing just some of the poles in Pade approximations is widely studied in literature. The common terminology used is Pade-type approximations.

It has been shown that the Pade approximation H.sub.n (s) of an RLC circuit transfer function will occasionally be unstable and/or non-passive despite the fact that RLC circuits are always stable and passive. The instability or the non-passivity of the transfer function H(n) is, in general, caused by a small number of misplaced poles and zeros. The effects of these poles and zeros on the approximate transfer function in the domain of interest must be minimal, once the accuracy of the Pade approximation becomes satisfactory.

We propose using the partial Pade approximation to generate stable and/or passive reduced-order models of circuit transfer functions through post-processing of the Pade approximation originally computed with the PVL algorithm.

Our method comprises: (1) computation of the Pade approximation of the transfer function via a Lanczos process (the PVL algorithm), (2) identification of the misplaced poles and/or zeros that are responsible for the instability or the non-passivity of the reduced-order model, if any, (3) prescription of alternative locations for these poles and/or zeros, (4) computation of the corresponding partial Pade approximation using a novel algorithm, and (5) a check of the stability and/or passivity of the result using the eigenvalue based criterion described in Z. Bai and R. W. Freund "An Eigenvalue-Based Test for Positive Realness of Transfer Functions," Numerical Analysis Manuscript, Bell Laboratories, Murray Hill, N.J. (Appendix A attached hereto). Since we argue that the effect of these poles and zeros on the transfer function in the domain of interest must be small, we expect little or no perturbation to the already converged poles and/or zeros.

Below, we propose a novel method for computing a partial Pade approximation and also demonstrate via numerical experiments its usefulness for enforcing stability and passivity in reduced-order models. The significant novelty in our approach for constructing a partial Pade approximation comprises the use of a matrix update technique on the tridiagonal matrix obtained by a Lanczos process instead of using moments explicitly as proposed by Brezinski. The resulting partial Pade approximation obtained via the Lanczos process and post-processing, which we coin PVL.pi., is significantly superior in terms of both numerical stability and computational efficiency.

IV. Post Processing of PVL: The PVL.pi. Algorithm

As remarked in the previous section, instead of using moments explicitly for the computation of the partial Pade approximation, we prefer to post-process the matrix T.sub.n generated by the PVL algorithm. The partial Pade approximation problem is therefore recast as the problem of finding an updating matrix F.sub.n such that some of the eigenvalues of the matrices T.sub.n =T.sub.n +F.sub.n and T'.sub.n =T'.sub.n +F'.sub.n are the prescribed values. Here, for any given matrix X, X' denotes the matrix obtained by deleting the first row and column of X. Recall that the zeros of the reduced-order transfer function H.sub.n are given by the eigenvalues of T'.sub.n. For a properly chosen updating matrix F.sub.n, we will have modified reduced-order transfer function associated with the updated matrix T.sub.n such that

H.sub.n (s.sub.0 +.sigma.)=(l.sup.T r)e.sub.1.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.1 =H(s.sub.0 +.sigma.)+O(.sigma..sup.2n-l-k) (Eq. 25) Using the formulation of the Lanczos process in terms of polynomials instead of vectors, it can be shown that the updating matrix F.sub.n must be a rank-one matrix of the form F.sub.n =xe.sub.n.sup.T Furthermore, if there are a total of m prescribed poles and zeros, where m=k+l.ltoreq.n, then the vector x can be chosen as ##EQU6## To verify the order of approximation of the modified reduced-order transfer function defined by T.sub.n =T.sub.n +xe.sub.n.sup.T by the Sherman-Morrison formula, we have e.sub.n.sup.T (I-.sigma.T.sub.n).sup.-1 =e.sub.1.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.1 -.tau..sigma.e.sub.1.sup.T (I-.sigma.T.sub.n).sup.-1 xe.sub.n.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.1 Here, the scalar, .tau., is given by .tau..sup.-1 =1-.sigma.e.sub.n.sup.T (I-.sigma.T.sub.n).sup.-1 x. By using the formula ##EQU7## where adj (I-.sigma.T.sub.n) stands for the adjoint matrix made up of (n-1).times.(n-1) cofactors of I-.sigma.T.sub.n, one can verify that ##EQU8## Therefore, the modified reduced-order transfer function defined by (Eq. 25) has the approximation order 2n-k-l, as defined in (Eq. 24). One also can verify directly by straightforward calculation that l.sup.T A.sup.j r=(l.sup.T r)e.sub.1.sup.T T.sub.n.sup.j e.sub.1 =(l.sup.T r)e.sub.1.sup.T T.sub.n.sup.j e.sub.1 for j=0,1,2, . . . ,2n-m-1. This again justifies that the modified reduce-order transfer function H.sub.n (s.sub.0 +.sigma.) is a partial Pade approximation to the original transfer function H(s.sub.0 +.sigma.).

The following algorithm computes the desired vector x and corresponds to the algorithm outlined in the flowchart of FIG. 3.

Algorithm 1 (Computation of the vector x)

INPUT:

{.theta..sub.j }.sub.j=1.sup.k and {.tau..sub.j }.sub.j=1.sup.l, where m=k+l.ltoreq.n.

OUTPUT: Vector x such that the matrices

T.sub.n =T.sub.n +xe.sub.n.sup.T and T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T have the prescribed eigenvalues {.theta..sub.j }.sub.j=1.sup.k and {.tau..sub.j }.sub.j=1.sup.l, respectively. 1) Compute the row vector .phi.(.theta..sub.j) for j=1,2, . . . ,k: if .theta..sub.j .epsilon..lambda.(T.sub.n), then .phi.(.theta..sub.j) is the corresponding left eigenvector, else .phi.(.theta..sub.j)=e.sub.n.sup.T (.tau..sub.j I-T.sub.n).sup.-1. 2) Compute the row vector .psi.(.tau..sub.j) for j=1,2, . . . ,l: if .tau..sub.j .epsilon..lambda.(T'.sub.n), then .psi.(.tau..sub.j) is the corresponding left eigenvector, else .psi.(.tau..sub.j)=e.sub.n-1.sup.T (.theta..sub.j I-T'.sub.n).sup.-1. 3) Compute the last column R.sub.k (:,n) of the k.times.n matrix ##EQU9## 4) Compute the last column S.sub.l (:,n-1) of the l.times.(n-1) matrix S.sub.l =.UPSILON..sub.l .PSI..sub.l -.PSI..sub.l T'.sub.n =diag(.tau..sub.j)[.psi.(.tau..sub.j)]-[.psi.(.tau..sub.j)]T'.sub.n. 5) Solve for x.sub.m : ##EQU10##

Or solve for x: ##EQU11##

Note that X(:,j) and X(:,i:j) denote the j-th column and the i-th to the j-th columns of the matrix X, respectively. Moreover, .lambda.(X) denotes the set of eigenvalues of X. If all matrices and vectors are real, and the prescribed complex eigenvalues .theta..sub.j and/or .tau..sub.j occur in complex conjugate pairs, if any, then the above algorithm can be modified slightly such that the vector x (and therefore the modified matrix T.sub.n) stays real.

Two special assignments of poles and zeros are of particular interest. The first case is when all poles are prescribed, i.e., k=n. Then we have

.PHI..sub.n xe.sub.n.sup.T =R.sub.n, where .PHI..sub.n is an n.times.n matrix. If .PHI..sub.n is nonsingular, then the vector x is given by x=.PHI..sub.n.sup.-1 R.sub.n (:,n). In this case, all eigenvalues of the matrix T.sub.n =T.sub.n +xe.sub.n.sup.T are the prescribed values, which implies that all poles of the modified reduced-order transfer function defined by T.sub.n are the prescribed values. Let p.sub.j =s.sub.0 +1/.lambda..sub.j for j=1, . . . n be the poles of the reduced-order model. Without loss of generality, assume that all poles are stable, except the last two poles p.sub.n-1 and p.sub.n, which are unstable. Let .theta..sub.j =1/(p.sub.j -s.sub.0) for j=1,2, . . . ,n-2. Let the unstable poles p.sub.n-1, and p.sub.n be assigned to the stable poles p.sub.n-1 and p.sub.n, and let .theta..sub.n-1 =1/(p.sub.n-1 -s.sub.0) and .theta..sub.n =1/(p.sub.n -s.sub.0). Then with all these prescribed poles, the modified reduced-order transfer function H.sub.n (s.sub.0 +.sigma.) is guaranteed to be stable.

The second special case is when all zeros are prescribed. In this case

.PSI..sub.n-1 x'.multidot.e.sub.n-1.sup.T =S.sub.n-1, where .PSI..sub.n-1 is an (n-1).times.(n-1) matrix. If it is nonsingular, then we have x'=.PSI..sub.n-1.sup.-1 S.sub.n-1 (:,n-1). As a result, all eigenvalues of the matrix T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T are the prescribed values, which means that all zeros of the modified reduced-order transfer function defined by T.sub.n are the prescribed values. With a similar discussion as for assigning poles, it can be guaranteed that the modified reduced-order model has no zeros in the right half-plane.

With the above algorithm to update partial or complete poles and/or zeros of the reduced-order transfer function, we are ready to present a new algorithm, which generates a partial Pade approximation via the Lanczos process. The new algorithm is based on the PVL algorithm and does post-processing in the presence of unstable poles and right-half-plane zeros of the reduced-order transfer function. The resulting computational procedure is herein coined the PVL.pi. algorithm. Algorithm 2 below corresponds to the algorithm outlined in the flowchart of FIG. 2.

Algorithm 2 (Sketch of the PVL.pi. algorithm) 1) Run n steps of the Lanczos process to obtain the tridiagonal matrix T.sub.n. This step corresponds to step 55 in FIG. 2. 2) Compute eigendecompositions of T.sub.n =S.sub.n diag(.lambda..sub.j)S.sup.-1.sub.n and T'.sub.n diag(.lambda.'.sub.j) (S'.sub.n).sup.-1 and obtain the approximate poles and zeros p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/'.sub.j, respectively. This step corresponds to step 60 in FIG. 2. 3) Check stability and/or passivity. If the reduced-order model defined by T.sub.n is stable and/or passive, then set T.sub.n =T.sub.n and go to Step 8). 4) Select partial or complete prescribed poles ##EQU12## Set .theta..sub.j =1/(p.sub.j -s.sub.0) and .tau..sub.j =1/(z.sub.j -s.sub.0). This step corresponds to step 85 in FIG. 2. 5) Use Algorithm 1 to compute the vector x. This step corresponds to step 90 in FIG. 2. 6) Compute eigendecompositions of T.sub.n =T.sub.n +xe.sub.n.sup.T =S.sub.n diag(.lambda..sub.j)S.sub.n.sup.-1 and T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T =S'.sub.n diag(.lambda.'.sub.j)(S'.sub.n).sup.-1 and obtain the modified approximate poles and zeros p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/.lambda.'.sub.j. This step corresponds to step 100 in FIG. 2. 7) Check stability and/or passivity. Go to Step 4) if necessary. 8) Compute the poles (zeros) and residues of the modified reduced order model H.sub.n (s.sub.0 +.sigma.)=(1.sup.T r)e.sub.1.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.1. This step corresponds to step 75 in FIG. 2.

Recall from Section 2 that all zeros being in the left half-plane is only a necessary condition for passivity, and that an additional positivity condition needs to be checked. In Algorithm 2, the "passivity check" refers to this necessary condition only.

V. EXAMPLES

In this section, we present two circuit examples to demonstrate the effectiveness of the PVL.pi. method.

A. The PEEC Circuit

Our first example is a circuit resulting from the PEEC (partial element equivalent circuit) discretization of an electromagnetic problem. This circuit has been used often as a test problem in the model-reduction literature. The circuit consists of only inductors, capacitors, and inductive couplings, and it is driven by a finite resistance source. We are interested in analyzing the response of this circuit, in this case the current flowing through one of the inductors. The circuit being mainly LC, its poles would lie exactly on the imaginary axis in the absence of the 100.OMEGA. source resistance. The source resistance causes the poles to move into the left half of the complex plane, but leaves most of them in the proximity of the imaginary axis. Due to this proximity, during the computation of a reduced-order model, numerical and approximation errors can easily produce unstable poles, i.e., poles lying in the right half of the complex plane. Indeed, running the PVL algorithm (with shift s.sub.0 -2.pi.10.sup.9) for n=60 iterations produces an almost exact transfer function in the 0-5 GHz range, as shown in FIG. 4, but results in 8 unstable poles and 7 nearly stable poles, see FIG. 6.

In order to produce a stable reduced-order model, we use PVL.pi. to force all 15 unstable poles, p.sub.j =Re(p.sub.j)+i IM(p.sub.j), into the left half-plane. This is done by setting, in Algorithm 1, ##EQU13## for all 15 unstable poles, while the stable poles are unchanged. FIG. 6 shows the poles of the reduced-order transfer function H.sub.60 (i.omega.) produced by the PVL algorithm, and the modified reduced-order transfer function

H.sub.60 (i.omega.) generated by the PVL.pi. algorithm. The modified, reduced-order transfer function H.sub.60 (s) has all its poles in the left half plane and is stable. FIG. 5 shows the frequency responses of the exact transfer function H.sub.60 (i.omega.) compared to H(i.omega.) produced by the PVL.pi. algorithm with 60 iterations, in the 0-5 GHz range. The match remains as good as before.

B. A Package Model

The next application is the analysis of a 64-pin package model used for an RF integrated circuit. Only eight of the package pins carry signal, the rest being either unused or carrying supply voltages. The package model is described by approximately 4000 circuit elements, resistors, capacitors, inductors, and inductive couplings. The size of the nodal circuit matrices for the package is about 2000. Here, we compute the input impedance of one of the pins. After executing the PVL algorithm for n=80 iterations, the resulting reduced-order model has a transfer function H.sub.80 defined by T.sub.80. The approximation H.sub.80 is sufficiently accurate up to 10 GHz, but has 2 poles and 4 zeros violating stability and passivity conditions.

We use PVL.pi. to force the 2 unstable poles, p.sub.j =Re(p.sub.j)+i Im(p.sub.j), and the 4 unstable zeros z.sub.j =Re(z.sub.j)+iIm(z.sub.j) into the left half-plane. This is done by setting, in Algorithm 1, ##EQU14## for the 2 unstable poles, and ##EQU15## for the 4 unstable zeros. The reduced-order model produced by the PVL.pi. has all poles and zeros of the modified reduced order transfer function given by

T.sub.80 =T.sub.80 +xe.sup.T in the left complex half-plane, where the vector x has only 2+4=6 non-zero entries at the bottom. By construction, the modified, reduced-order model is stable, and it also satisfies the necessary condition that a passive system cannot have zeros in the right half-plane. To verify passivity, an additional positivity condition has to be satisfied; see, e.g. B. D. O. Anderson and S. Vongpanitlerd, "Network Analysis and Synthesis," Englewood Cliffs, N.J.: Prentice-Hall, 1973. In Z. Bai and R. W. Freund, "An Eigenvalue-Based Test for Positive Realness of Transfer Functions," Numerical Analysis Manuscript, Bell Laboratories, Murray Hill, N.J. (Appendix A attached hereto), we derived a simple eigenvalue-based test for this positivity condition. Using this test, one easily verifies that the reduced-order model is indeed passive. FIG. 7 shows the frequency response of the exact transfer function and the reduced-order transfer function generated by the PVL algorithm. FIG. 7 shows the frequency response of the exact transfer function and the reduced-order transfer function generated by the PVL.pi. algorithm. FIGS. 7 and 8 clearly illustrate that no accuracy has been lost by enforcing stability and passivity by means of our post-processing procedure.

We term the method and apparatus of the present invention partial Pade via Lanczos modeling, or PVL .pi. for short. The scalar method and apparatus disclosed herein should be extensible to matrix Pade approximants. Several algorithms for checking the stability and passivity of the original Lanczos matrix T.sub.n are known. One particularly elegant solution is found in appendix A hereto which is a paper entitled "Eigenvalue-Based Characterization and Test for Positive Realness of Transfer Functions" by Zhaojoun Bai and Roland Freund.

Having thus described a particular embodiment of the invention, various alterations, modifications, and improvements will readily occur to those skilled in the art. Such alterations, modifications and improvements as are made obvious by this disclosure are intended to be part of this description though not expressly stated herein, and are intended to be within the spirit and scope of the invention. Accordingly, the foregoing description is by way of example only, and not limiting. The invention is limited only as defined in the following claims and equivalents thereto.

Claims

1. An apparatus for generating a frequency response signal of a circuit, comprising:

a processing unit having:

a first program for generating, from circuit characteristic data, a block tridiagonal matrix relating to said circuit characteristic data;

a second program for determining if said block tridiagonal matrix is stable and passive;

a third program for modifying said block tridiagonal matrix to make it stable and passive, if said block tridiagonal matrix is not stable and passive;

a fourth program for generating a frequency response signal as a function of said stable and passive block tridiagonal matrix or modified block tridiagonal matrix; and

an output device for displaying said circuit frequency response signal.

2. An apparatus as set forth in claim 1 wherein said second program comprises a program for generating eigenvalues and eigenvectors of said tridiagonal matrix.

3. An apparatus as set forth in claim 2 wherein said third program comprises a program for determining a vector that can be added to a last column of said tridiagonal matrix that will make said matrix passive and stable and adds said vector to said last column to generate said modified block tridiagonal matrix.

4. An apparatus as set forth in claim 3 wherein said third program further comprises a program for determining eigenvalues and eigenvectors of said modified block tridiagonal matrix and determining if said modified block tridiagonal matrix is passive and stable.

5. An apparatus as set forth in claim 4 further comprising a program for running said third program repeatedly to further modify said modified block tridiagonal matrix, if said modified block tridiagonal matrix is not passive and stable.

6. An apparatus as set forth in claim 5 wherein said fourth program determines a transfer function of said block tridiagonal matrix or modified block tridiagonal matrix by calculating poles, zeros and residues from said eigenvalues and eigenvectors to generate a Pade approximant of a transfer function of said circuit.

7. An apparatus as set forth in claim 6 wherein said first program comprises a program for processing said circuit characteristic data to obtain a characteristic matrix and a characteristic vector of said circuit, and reducing the characteristic matrix to produce said block tridiagonal matrix by a look-ahead Lanczos procedure.

8. An apparatus as set forth in claim 1 wherein said third program comprises:

a program for computing a Pade approximation of a transfer function of said block tridiagonal matrix;

a program for identifying any misplaced poles and zeros that contribute to a lack of stability or passivity of said block tridiagonal matrix;

a program for prescribing alternative locations for said identified poles and zeros;

a program for computing a partial Pade approximation of a transfer function of said block tridiagonal matrix; and

a program for checking stability and passivity of said partial Pade approximation.

9. An apparatus as set forth in claim 8 wherein said third program further comprises a program for repeating said identifying, prescribing, computing and checking programs of claim 8 upon said modified block tridiagonal matrix to further modify it, if said partial Pad6 approximation is not passive and stable.

10. An apparatus as set forth in claim 1 wherein said third program comprises:

a program for computing a row vector .phi.(.theta..sub.j) for j=1,2, . . . ,k;

a program for computing a row vector .psi.(.tau..sub.j) for j=1,2, . . . ,l

a program for computing a last column R.sub.k (:,n) of a matrix: ##EQU16## a program for computing a last column Sl(:,n-1) of a matrix:

S.sub.l =.UPSILON..sub.l .psi..sub.l -.psi..sub.l T'.sub.n :=diag(.tau..sub.j)[.psi.(.tau..sub.j)]-[.psi.(.tau..sub.j)]T'.sub.n ; and

a program for solving for x.sub.m the equation ##EQU17## or solving for x: ##EQU18##

11. An apparatus as set forth in claim 10 wherein said second program comprises: a program for computing eigendecompositions of T.sub.n =S.sub.n diag (.lambda..sub.j)S.sup.-1.sub.n and T'.sub.n diag(.lambda.'.sub.j) (S'.sub.n).sup.-1 and obtaining approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/'.sub.j,

respectively; and

a program for checking stability and/or passivity. If the reduced-order model defined by T.sub.n is stable and/or passive.

12. An apparatus as set forth in claim 11 wherein said third program comprises:

a program for computing eigendecompositions of

T.sub.n =T.sub.n +xe.sub.n.sup.T =S.sub.n diag(.lambda..sub.j)S.sub.n.sup.-1

and

T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T =S'.sub.n diag(.lambda.'.sub.j)(S'.sub.n).sup.-1

and obtaining the modified approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/.lambda.'.sub.j.

a program for checking stability and passivity; and

a program for computing the poles and zeros and residues of the modified reduced order model

H.sub.n (s.sub.0 +.sigma.)=(l.sup.T r)e.sub.l.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.l.

13. An apparatus for generating a frequency response signal indicative of a frequency response of a circuit, the apparatus comprising:

a data source providing circuit characteristic data representing the circuit; and

a processing unit having:

a memory for receiving and storing the circuit characteristic data from the data source;

a first program for `generating`, from circuit characteristic data, a block tridiagonal matrix relating to said circuit characteristic data;

a second program for determining if said block tridiagonal matrix is stable and passive;

a third program for modifying said block tridiagonal matrix to make it stable and passive, if said block tridiagonal matrix is not stable and passive; and

a fourth program for generating a frequency response signal as a function of said stable and passive block tridiagonal matrix or said modified block tridiagonal matrix.

14. An apparatus as set forth in claim 13 wherein said second program comprises a program for generating eigenvalues and eigenvectors of said tridiagonal matrix.

15. An apparatus as set forth in claim 14 wherein said third program comprises a program for determining a vector that can be added to a last column of said tridiagonal matrix that will make said matrix passive and stable and a program for adding said vector to said last column to generate said modified block tridiagonal matrix.

16. An apparatus as set forth in claim 15 wherein said third program further comprises a program for determining eigenvalues and eigenvectors of said modified block tridiagonal matrix and determining if said modified block tridiagonal matrix is passive and stable.

17. An apparatus as set forth in claim 16 further comprising a program for running said third program repeatedly to further modify said modified block tridiagonal matrix, if said modified block tridiagonal matrix is not passive and stable.

18. An apparatus as set forth in claim 17 wherein said fourth program determines a transfer function of said stable and passive block tridiagonal matrix or said modified block tridiagonal matrix by calculating poles, zeros and residues from said eigenvalues and eigenvectors to generate a Pade approximant of a transfer function of said circuit.

19. An apparatus as set forth in claim 13 wherein said third program comprises:

a program for computing a Pade approximation of a transfer function of said block tridiagonal matrix;

a program for identifying any misplaced poles and zeros that contribute to a lack of stability or passivity of said block tridiagonal matrix;

a program for prescribing alternative locations for said identified poles and zeros;

a program for computing a partial Pade approximation of a transfer function of said matrix; and

a program for checking stability and passivity of said partial Pade approximation.

20. An apparatus as set forth in claim 19 wherein said third program comprises a program for repeating said identifying, prescribing, computing and checking programs set forth in claim 19 upon said modified block tridiagonal matrix to further modify it, if said partial Pade approximation is not passive and stable.

21. An apparatus as set forth in claim 13 wherein said third program comprises:

a program for computing a row vector .phi.(.theta..sub.j) for j=1,2, . . . ,k;

a program for computing a row vector .psi.(.tau..sub.j) for j=1,2, . . . ,l

a program for computing a last column R.sub.k (:,n) of a matrix: ##EQU19## a program for computing a last column Sl(:,n-1) of a matrix:

S.sub.l =.UPSILON..sub.l .psi..sub.l -.psi..sub.l T'.sub.n :=diag(.tau..sub.j)[.psi.(.tau..sub.j)]-[.psi.(.tau..sub.j)]T'.sub.n ; and

a program for solving for x.sub.m the equation ##EQU20## or solving for x: ##EQU21##

22. An apparatus as set forth in claim 21 wherein said second program comprises: a program for computing eigendecompositions of

T.sub.n =S.sub.n diag(.lambda..sub.j)S.sup.-1.sub.n and T'.sub.n diag(.lambda.'.sub.j)(S'.sub.n).sup.-1 and obtaining approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/'.sub.j,

respectively; and

a program for checking stability and/or passivity. If the reduced-order model defined by T.sub.n is stable and/or passive.

23. An apparatus as set forth in claim 22 wherein said third program comprises:

a program for computing eigendecompositions of

T.sub.n =T.sub.n +xe.sub.n.sup.T =S.sub.n diag(.lambda..sub.j)S.sub.n.sup.-1

and

T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T =S'.sub.n diag(.lambda.'.sub.j)(S'.sub.n).sup.-1

and obtaining the modified approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/.lambda.'.sub.j.

a program for checking stability and passivity; and

a program for computing the poles and zeros and residues of the modified reduced order model

H.sub.n (s.sub.0 +.sigma.)=(l.sup.T r)e.sub.l.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.l.

24. A method for generating a frequency response signal of a circuit, the method comprising the steps of:

providing circuit characteristic data associated with the circuit;

generating, from circuit characteristic data, a block tridiagonal matrix relating to the circuit characteristic data;

determining if said block tridiagonal matrix is stable and passive;

modifying said block tridiagonal matrix to make it stable and passive, if said block tridiagonal matrix is not stable and passive; and

producing a frequency response signal as a function of said stable and passive block tridiagonal matrix or modified block tridiagonal matrix.

25. A method as set forth in claim 24 wherein said determining step comprises generating eigenvalues and eigenvectors of said tridiagonal matrix.

26. A method as set forth in claim 25 wherein said modifying step comprises the steps of:

determining a vector that can be added to a last column of said tridiagonal matrix that will make said matrix passive and stable; and

adding said vector to said last column to generate said modified block tridiagonal matrix.

27. A method as set forth in claim 26 wherein said modifying step further comprises;

determining eigenvalues and eigenvectors of said modified block tridiagonal matrix; and

determining if said modified block tridiagonal matrix is passive and stable.

28. A method as set forth in claim 27 further comprising the step of:

repeating said modifying step to further modify said modified block tridiagonal matrix, if said modified block tridiagonal matrix is not passive and stable.

29. A method as set forth in claim 28 wherein said producing step comprises the step of determining said transfer function by calculating poles, zeros and residues from said eigenvalues and eigenvectors to generate a Pade approximant of a transfer function of said circuit.

30. A method as set forth in claim 29 wherein said generating step comprises the steps of:

processing said circuit characteristic data to obtain a characteristic matrix and a characteristic vector of said circuit; and

reducing the characteristic matrix to produce said block tridiagonal matrix by a look-ahead Lanczos procedure.

31. A method as set forth in claim 24 wherein said modifying step comprises the steps of:

identifying any misplaced poles and zeros that contribute to a lack of stability or passivity of said block tridiagonal matrix;

prescribing alternative locations for said identified poles and zeros to generate said modified block tridiagonal matrix;

computing a partial Pade approximation of a transfer function of said modified block tridiagonal matrix; and

checking stability and passivity of said partial Pade approximation.

32. An apparatus as set forth in claim 31 wherein said modifying step further comprises repeating said identifying, prescribing, computing and checking steps set forth in claim 31 upon said modified block tridiagonal matrix to further modify it, if said partial Pade approximation is not passive and stable.

33. An apparatus as set forth in claim 24 wherein said modifying step comprises the steps of:

computing a row vector .phi.(.theta..sub.j) for j=1,2, . . . ,k;

computing a row vector .psi.(.tau..sub.j) for j=1,2, . . . ,l

computing a last column R.sub.k (:,n) of a matrix: ##EQU22## computing a last column Sl(:,n-1) of a matrix:

S.sub.l =.UPSILON..sub.l .psi..sub.l -.psi..sub.l T'.sub.n :=diag(.tau..sub.j)[.psi.(.tau..sub.j)]-[.psi.(.tau..sub.j)]T'.sub.n ; and

solving for x.sub.m the equation ##EQU23## Or solve for x: ##EQU24##

34. An apparatus as set forth in claim 33 wherein said determining step comprises the steps of: computing eigendecompositions of T.sub.n =S.sub.n diag(.lambda..sub.j)S.sup.-1.sub.n and T'.sub.n diag(.lambda.'.sub.j)(S'.sub.n) and obtaining approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/'.sub.j,

respectively; and

checking stability and/or passivity; and

if the reduced-order model defined by T.sub.n is stable and/or passive, then setting

T.sub.n =T.sub.n.

35. An apparatus as set forth in claim 34 wherein said modifying step comprises the steps of:

computing eigendecompositions of

T.sub.n =T.sub.n +xe.sub.n.sup.T =S.sub.n diag(.lambda..sub.j)S.sub.n.sup.-1

and

T'.sub.n =T'.sub.n +x'e.sub.n-1.sup.T =S'.sub.n diag(.lambda.'.sub.j)(S'.sub.n).sup.-1

and obtaining modified approximate poles and zeros

p.sub.j =s.sub.0 +1/.lambda..sub.j and z.sub.j =s.sub.0 +1/.lambda.'.sub.j.

checking stability and passivity; and

computing poles and zeros and residues of said modified block tridiagonal matrix

H.sub.n (s.sub.0 +.sigma.)=(l.sup.T r)e.sub.l.sup.T (I-.sigma.T.sub.n).sup.-1 e.sub.l.

36. An apparatus comprising:

a processing unit having:

a memory for storing circuit characteristic data representing a circuit;

a first program for generating, from circuit characteristic data, a block tridiagonal matrix relating to the circuit characteristic data;

a second program for determining if said block tridiagonal matrix is stable and passive;

a third program for modifying said block tridiagonal matrix to make it stable and passive, if said block tridiagonal matrix is not stable and passive; and

a fourth program for generating a frequency response signal as a function of said stable and passive block tridiagonal matrix or modified block tridiagonal matrix; and

an output device for displaying said circuit frequency response signal.

37. An apparatus comprising:

a memory for storing circuit characteristic data representing a circuit;

a first program for generating, from circuit characteristic data, a block tridiagonal matrix relating to the circuit characteristic data;

a second program for determining if said block tridiagonal matrix is stable and passive; and

a third program for modifying said block tridiagonal matrix to make it stable and passive, if said block tridiagonal matrix is not stable and passive;

a processor for executing the programs; and

means for generating a frequency response signal in response to the matrix transfer function data.