|
"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in
message <fuas3d$du7$1@fred.mathworks.com>...
> Aviral, if you substitute z = exp(j*w) in your definite integrals, they
become
> contour integrals taken with respect to z and moving counterclockwise
> around the unit circle in the complex plane for the matrix-valued
integrands
>
> -j*inv(z*I-A)*B*C*inv(1/z*I-A)/z.
>
> Cramer's rule shows us that both these two inverses are rational functions
of
> z for each of the elements of the matrices. Thus each of the elements of
the
> above integrands is some rational and therefore meromorphic function of
z.
> We can find their integrals around the closed circle by simply evaluating all
> their single poles which lie within the unit circle, by virtue of the theory of
> residues.
>
> In principle, matlab's symbolic toolbox ought to be able to determine
these
> residues and therefore evaluate the integrals without resorting to numerical
> quadrature. I am not sure how well the toolbox can actually perform at this
> for large size matrices because the complexity of expressions increases
> rapidly as their size increases, but you can try putting it through its paces
> along these lines and see how well it performs.
>
> Roger Stafford
------------
It has occurred to me on thinking of the matter at greater length, that the
poles of the meromorphic function I referred to earlier are precisely the
eigenvalues of the matrix A and their reciprocals, a fact which should
considerably ease the burden of symbolic toolbox manipulation.
Again referring to Cramer's rule, the inverse of z*I-A can be seen to be the
ratio of two polynomials in z. The denominator is just det(z*I-A), which in
turn equals (z-L1)*(z-L2)*...*(z-Ln) where L1, L2, ..., Ln are the eigenvalues of
A. The numerator will be some polynomial of no more than the n-1 degree.
A similar situation exists for the inverse of 1/z*I-A in respect to the
reciprocals of A's eigenvalues. If the numerator and denominator from
Cramer's rule are multiplied by (-z)^n, and divided by the product of the
eigenvalues, then the denominator becomes (z-1/L1)*(z-1/L2)*...*(z-1/Ln)
and the numerator will become z times a polynomial of degree no more than
n-1. This first z factor will cancel out the z in the denominator of the original
integrand, namely at the right end of
-j*inv(z*I-A)*B*C*inv(1/z*I-A)/z.
Provided no eigenvalue lies right on the unit circle itself, (which would create
a singularity smack on the path of integration,) for each eigenvalue, either it
lies in the interior of the unit circle or its reciprocal does, but not both. This
gives us precisely n poles to be found in the interior of our meromorphic
integrand. The definite integral you seek is then equal to 2*pi*j times the
sum of the residues at these n poles.
One way of determining these residues is this. Suppose eigenvalue L lies
inside the unit circle. If we take the symbolic derivative of the above
denominator, det(z*I-A), with respect to z, and then assign z in it to the
numerical value z = L, that is the equivalent of having divided det(z*I-A) by
z-L symbolically and then set z numerically equal to L in the remaining
factor. However, doing it this way avoids the problem of not possessing a
symbolic expression for L for sizes of A greater than 4x4. The entire
integrand, but corrected by having the z-L factor in the denominator
removed in this fashion, can then be evaluated numerically at z = L and the
result would be just the residue at this pole. An analogous method would be
used to adjust det(1/z*I-A) in a similar way for an eigenvalue occurring
outside the circle.
This is admittedly something I haven't tried out on my computer as yet, but
the theory appears valid to me. The main practical difficulty I foresee is that
the toolbox may have trouble obtaining the above symbolic inverses for large
size matrices, A. Some of the expressions may not bear looking at by us
mere mortals for fear we will turn to stone, though we probably don't really
need to look at them to obtain numerical results from them. (On my system
no expression string is permitted to exceed 8192 characters, so I am
protected to a degree from excessive shock from such mind-boggling
visions. On the other hand, I can't get the symbolic inverse of z*I-A for
anything past a 3 x 3 matrix, A.)
As to whether this method is superior to doing the integration by numerical
quadrature I cannot be sure. Numerical quadrature could encounter severe
accuracy problems if some of the eigenvalues happen to lie very near the unit
circle, as might be anticipated. Using the residue method would avoid such
problems for the most part, except that the residue obtained for the pole just
inside the circle, even though accurate, would presumably be very large. On
the other hand, it is probably much less laborious to prepare the code for
doing this problem numerically. In any case, you could regard this little
dissertation as strengthening your theoretical understanding of the nature of
the task facing you.
Roger Stafford
|
Feed for this Thread
