Hi, I need help to calculate Eigenvectors using MKL routines.
I have some codes that were compiled using Compaq Visual Fortran and they had IMSL libraries within. Now I'm using Visual Studio 2010+Visual Fortran XE 2013. And I'm trying to convert those old codes without IMSL, but using MKL (Lapack) functions. In this case is the function DEVCRG ( http://www.roguewave.com/Portals/0/products/imsl-numerical-libraries/for... ). As you can check they have an example in that page. That's the same matrix i'm trying in the code bellow. I've checked with MATLAB and I've found the same results.
So, here we go with MKL, following the decision tree ( http://software.intel.com/sites/products/documentation/hpc/mkl/mklman/GU... ) in the manual I'm trying in this sequence those routines: ?gebal , ?gehrd , ?unghr , ?hseqr , ?trevc , ?gebak
I've arrived to find the eigenvalues, but I can't find the eigenvectors!!
program eigenvalvec use MKL95_LAPACK implicit none double precision, allocatable :: x(:,:),scalev(:) double complex, allocatable :: xc(:,:),tauc(:),w(:),uph(:,:),q(:,:),vr(:,:),vl(:,:),qz(:,:) logical, allocatable :: select(:) integer :: n = 3.0 integer :: nouput,i,ilo,ihi,info character*14 file10 allocate(x(n,n),scalev(n),xc(n,n),tauc(n-1),w(n),uph(n,n),q(n,n),select(n),vr(n,n),vl(n,n),qz(n,n)) nouput=1 x(1,1)=8.0 ; x(1,2)=-1.0 ; x(1,3)=-5.0 x(2,1)=-4.0 ; x(2,2)=4.0 ; x(2,3)=-2.0 x(3,1)=18.0 ; x(3,2)=-5.0 ; x(3,3)=-7.0 xc = CMPLX(x) call gebal(xc,scalev,ilo,ihi,'P',info) !Balances a general matrix to improve the accuracy of computed eigenvalues and eigenvectors uph = xc call gehrd(uph,tauc,ilo,ihi,info) !Reduces a general matrix to upper Hessenberg form q = uph call unghr(q,tauc,ilo,ihi,info) !Generates the complex unitary matrix Q determined by ?gehrd. call hseqr(uph,w,ilo,ihi,q,'S','I',info) !Computes all eigenvalues and (optionally) the Schur factorization of a matrix reduced to Hessenberg form. !w: EIGENVALUES!!! vl = q vr = q call trevc(uph,'B',select,vl,vr,n,info) !Computes selected eigenvectors of an upper (quasi-)triangular matrix computed by ?hseqr. call gebak(vr,scalev,ilo,ihi,'B','R',info) !Transforms eigenvectors of a balanced matrix to those of the original nonsymmetric matrix. end program eigenvalvec
Please, I really need some light on that... specially with "?trevc". Thank you very much!!!