## Which subroutine does complex Hermitian tridiagonal matrix?

Open discussion regarding features, bugs, issues, vendors, etc.

### Which subroutine does complex Hermitian tridiagonal matrix?

which subroutine computes all eigenvalues and, optionally, eigenvectors of complex symmetric tridiagonal matrix A?

I used dstev.f for real symmetric tridiagonal matrix, but now I need to use matrix in complex number, I am wondering which subroutine I can call for?

Thank you!
Last edited by carolasu on Wed Aug 15, 2012 1:09 am, edited 1 time in total.
carolasu

Posts: 9
Joined: Wed Mar 21, 2012 12:50 am

### Re: Which subroutine does complex symmetric tridiagonal matr

Do you mean complex symmetric tridiagonal or complex Hermitian tridiagonal?
Julien Langou

Posts: 835
Joined: Thu Dec 09, 2004 12:32 pm
Location: Denver, CO, USA

### Re: Which subroutine does complex symmetric tridiagonal matr

Thanks. I mean complex Hermitian tridiagonal.
carolasu

Posts: 9
Joined: Wed Mar 21, 2012 12:50 am

### Re: Which subroutine does complex symmetric tridiagonal matr

This is easy then. For any complex Hermitian tridiagonal matrix, there exists a symmetric scaling (so scaling of the rows and the columns) such that the matrix becomes a real symmetric tridiagonal matrix. In another word, given A, complex Hermitian tridiagonal matrix, there exists S, unitary complex diagonal, such that T=S*A*conj(S) is real symmetric tridiagonal matrix. Note that S diagonal unitary simply means that the diagonal entries have modulus 1. So actually T = S*A*conj(S) is a similarity transformation so A and T have the same eigenvalues.

Therefore, you quickly convert your input A to S, then you use LAPACK real symmetric tridiagonal eigensolver (dstev, dstemr, dstedc, dsteqr, ... ) to get T = V * W * V^H (where V is the matrix of eigenvectors, it is real orthogonal, and W is the matrix of eigenvalues, it is real diagonal) and you do not forget to scale back the V to get your true eigenvectors: A = ( conj(S)*V ) * W * ( conj(S)*V )^H. W stays real, the eigenvectors conj(S)*V becomes complex. conj(S)*V is a unitary matrix.

(Note: this is the reason why there is no complex Hermitian tridiagonal eigensolver. The real symmetric tridiagonal eigensolvers are all that are needed.)

Cannot find any reference for the algo online so here is a quick demo in matlab.

Code: Select all
`   clear;    n=10;%  A is complex Hermitian tridiagonal %  create D diagonal of A and E sub-diagonal of A   E(1:n-1,1)=randn(n-1,1)+i*randn(n-1,1);   D(1:n,1)=randn(n,1); %  create AA n-by-n from D and E (needed for check only)   AA=zeros(n); for i=1:n, AA(i,i)=D(i); end; for i=1:n-1, AA(i,i+1)=E(i); AA(i+1,i)=conj(E(i)); end;%  compute S a unitary complex diagonal matrix such that T = S*A*S' complex symmetric tridiagonal%  note that S diagonal unitary simply means that the diagonal entries have modulus 1%  this also means that S'*S = I so actually T = S*A*S' is a similarity transformation so A and T %  have the same eigenvalues%  compute at the same time the matrix T%  the matrix T is stored on D and E and replace the initial matrix A   S=zeros(n,1);   S(1) = 1;   if ( abs(E(1)) == 0 )     S(2) = 1;   else     S(2) = E(1)/abs(E(1));     E(1) = abs(E(1));   end   for i=2:n-1,      if ( abs(E(i)) == 0 )         S(i+1) = 1;      else         S(i+1) = S(i)*E(i)/abs(E(i));         E(i) = abs(E(i)); % ( note that E is real now )      end;   end;%  create TT n-by-n symmetric real from D and E %  this is needed in Matlab to call the eigensolver%  in Fortran: not needed: call directly dstemr, dstedc, dsteqr etc. directly with D and E   TT=zeros(n); for i=1:n, TT(i,i)=D(i); end; for i=1:n-1, TT(i,i+1)=E(i); TT(i+1,i)=E(i); end;%  calling the eigensolver to obtain eigenvectors V and associated eigenvalues W of T%  Note: T is real symmetric tridiagonal   [V,W] = eig(TT);%  backtransform the eigenvectors V of T to obtain the eigenvectors of A%  so perform: V <- S' * V   for i=2:n,      V(i,1:n) = conj(S(i))*V(i,1:n);   end%  check   fprintf('|| A - V W V^H || / || A || = %2.2e \n', norm( AA - V*W*V' )/norm( AA ));   fprintf('|| I - V V^H ||             = %2.2e \n', norm( eye(n) - V*V' ));`

( In a real code, you work with array of size n for D, array of size n-1 for E and array of size n for S. No need to store n-by-n matrices. )

Cheers,
Julien.
Last edited by Julien Langou on Wed Aug 15, 2012 3:05 am, edited 1 time in total.
Julien Langou

Posts: 835
Joined: Thu Dec 09, 2004 12:32 pm
Location: Denver, CO, USA

### Re: Which subroutine does complex Hermitian tridiagonal matr

I see. It's very clear.

(1)complex A stored as real D, complex E.

(2)Use S to scale complex E to real E.

(3)Call dstev.f with real E and real D to get V, W.

(4)Use S to scale real V to complex V.

AA and TT are not needed, but we can use those to check.

Thank you!
Last edited by carolasu on Wed Aug 15, 2012 3:51 am, edited 1 time in total.
carolasu

Posts: 9
Joined: Wed Mar 21, 2012 12:50 am

### Re: Which subroutine does complex Hermitian tridiagonal matr

Does anyone know if there is a theorem, which says that you can transform a hermitian tridiagonalmatrix in a real symmetric tridiagonalmatrix by similarity transform with a unitary diagonalmatrix? And where i can find the theorem?

Thanks for help!!!
mattcimd

Posts: 1
Joined: Sun Nov 11, 2012 5:49 pm