Eigen decomposition of a matrix of form W * diag(S) * W' for matrix exponential in MATLAB

W is a tall and skinny real valued matrix, and diag(S) is a diagonal matrix consists of +1 or -1 on the diagonal. I want the eigen decomposition of A = W * diag(S) * W' where single quote denotes transposition. The main problem is that A is pretty big. Since A is symmetric, rank deficient, and I actually know the maximum rank of A (from W), I think I should be able to do this efficiently. Any idea how to approach this?

My eventual goal is to compute the matrix exponential of A without using MATLAB's expm which is pretty slow for big matrices and does not take advantage of rank deficiency. If A = U * diag(Z) * U' is the eigen decomposition, exp(A) = U * diag(exp(Z)) * U'.

While finding an orthogonal U such that W * diag(S) * W' = U' * diag(Z) * U' looks promising to have an easy algorithm, I need some linear algebr开发者_如何学Pythona help here.

I'd first perform the so called 'thin' QR factorization of W, then compute the eigenvalue decomposition of R*diag(S)*R', then use this to compute the eig decomposition of A.

N = 10;
n=3;
S = 2*(rand(1,n)>0.5)-1;
W = rand(N,n);

[Q,R] = qr(W,0);
[V,D] = eig(R*diag(S)*R');

%this is the non rank-deficient part of eig(W*diag(S)*W')
D_A = D;
V_A = Q*V;

%compare with
[V_full,D_full] = eig(W*diag(S)*W');

Hope this helps.

A.

MATLAB actually has an implementation for retrieving the largest (or smallest) eigen values and vectors. Use eigs(A,k) to get the k largest.

To get the largest only, one can use the Power iteration method, or better yet Rayleigh quotient iteration.