I thought that the function eigen(A) will return a matrix with eigenvectors that are independent of each other (thus forming a base and the matrix being invertible). This seems not to be the case in the following example A=matrix(c(1,2,0,1),nrow=2,byrow=T) eigen(A) ->ev solve(ev$vectors) note that I try to get the upper triangular form with eigenvalues on the diagonal and (possibly) 1 just atop the eigenvalues to be used to solve a linear differential equation x' = Ax, x(0)=x0 x(t) = P exp(D t) P^-1 x0 where D is this upper triangular form and P is the "passage matrix" (not sure about the correct english name) given by a base of eigenvectors. So the test would be solve(ev$vectors) %*% A %*% ev$vectors - D should be 0 Thanks for any help, Christian. ps: please copy reply also to my address, my subscription to the R-help list seems to have delays -- *********************************************************** http://cognition.ups-tlse.fr/vas-y.php?id=chj jost at cict.fr Christian Jost (PhD, MdC) Centre de Recherches sur la Cognition Animale Universite Paul Sabatier, Bat IV R3 118 route de Narbonne 31062 Toulouse cedex 4, France Tel: +33 5 61 55 64 37 Fax: +33 5 61 55 61 54
Christian Jost wrote:> I thought that the function > eigen(A) > will return a matrix with eigenvectors that are independent of each > other (thus forming a base and the matrix being invertible). This seems > not to be the case in the following example > A=matrix(c(1,2,0,1),nrow=2,byrow=T) > eigen(A) ->ev > solve(ev$vectors) > > note that I try to get the upper triangular form with eigenvalues on the > diagonal and (possibly) 1 just atop the eigenvalues to be used to solve > a linear differential equation > x' = Ax, x(0)=x0 > x(t) = P exp(D t) P^-1 x0 > where D is this upper triangular form and P is the "passage matrix" (not > sure about the correct english name) given by a base of eigenvectors. So > the test would be > solve(ev$vectors) %*% A %*% ev$vectors - D > should be 0 > > Thanks for any help, Christian. > > ps: please copy reply also to my address, my subscription to the R-help > list seems to have delaysThat particular matrix has repeated eigenvalues and a degenerate eigenspace. > A <- matrix(c(1,0,2,1),nc=2) > A [,1] [,2] [1,] 1 2 [2,] 0 1 > eigen(A) $values [1] 1 1 $vectors [,1] [,2] [1,] 1 -1.000000e+00 [2,] 0 1.110223e-16
Christian Jost wrote:> I thought that the function > eigen(A) > will return a matrix with eigenvectors that are independent of each > other (thus forming a base and the matrix being invertible). This > seems not to be the case in the following example > A=matrix(c(1,2,0,1),nrow=2,byrow=T) > eigen(A) ->ev > solve(ev$vectors) >I guess eigen tries to get independent eigenvectors, butr that is not always possible, and your matrix is a case of that. Note that all eigenvectors of A are a multiple of (0,1)^T, so there cannot be two independent ones. Kjetil> note that I try to get the upper triangular form with eigenvalues on > the diagonal and (possibly) 1 just atop the eigenvalues to be used to > solve a linear differential equation > x' = Ax, x(0)=x0 > x(t) = P exp(D t) P^-1 x0 > where D is this upper triangular form and P is the "passage matrix" > (not sure about the correct english name) given by a base of > eigenvectors. So the test would be > solve(ev$vectors) %*% A %*% ev$vectors - D > should be 0 > > Thanks for any help, Christian. > > ps: please copy reply also to my address, my subscription to the > R-help list seems to have delays-- Kjetil Halvorsen. Peace is the most effective weapon of mass construction. -- Mahdi Elmandjra
Hi, Doug:
What about returning list(T., Z) or list(quT, Z), where T. or quT
is quasi-upper triangular and Z is orthogonal, consistent with the
documentation? This would make it easier for people like me to build an
intuition for its properties by playing with it. If it finds
substantial use, someone may later suggest a special storage format for
the quasi-upper-triangular part.
Thank for your work on this.
Spencer Graves
Douglas Bates wrote:
> Spencer Graves wrote:
>
>> Does R have a function for the Schur decomposition? The
>> documentation for library(Matrix) describes a function
"Schur", but
>> it seems to be missing from the Windows version 0.8-14 (2004-09-14)
>> and 0.8-15 (2004-10-02).
>> The R 2.0.0 pat documentation for "eigen" refers to
>> "http://www.netlib.org/lapack/lug/lapack_lug.html", and the
>> description there for eigen analysis of a non-symmetric matrix says,
>> "This problem can be solved via the Schur factorization of A,
defined
>> in the real case as
>> A = ZTZT,
>>
>> where Z is an orthogonal matrix and T is an upper quasi-triangular
>> matrix with 1-by-1 and 2-by-2 diagonal blocks, the 2-by-2 blocks
>> corresponding to complex conjugate pairs of eigenvalues of A."
>> Thanks,
>> Spencer Graves
>
>
> That documentation entry was a note to myself to add the Schur
> factorization to the Matrix package but I have not yet done so. If
> anyone can suggest a reasonable data structure for the result, it
> would be fairly easy to add the Schur decomposition.
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
>
--
Spencer Graves, PhD, Senior Development Engineer
O: (408)938-4420; mobile: (408)655-4567