 General CULA Dense (LAPACK & BLAS) support and troubleshooting. Use this forum if you are having a general problem or have encountered a bug.

Hello,
my problem is an inaccuracy of eigenvectors computed by culaDsteqr. In the file "eivectors.pdf" is the 660'th eigenvector, computed with three different algorithms (matlab, Numerical Recipies(Press), CULA) plotted. Whereby matlab and the NR delivers the same eigenvector, some components of the CULA eigenvector differ more than one decade! Why?

The file "test.c" was my code, to make data to check CULA's accuracy. "data.txt" contains the diagonals of the input tridiagonal matrix. I use SUSE 11.3, CULA premium R10, GeForce GTX 275 and CUDA 3.2 .

Yours sincerely
Malte Sommer
Attachments test.c
code to test CULA data.txt
first row: diagonal
second row: off diagonal eivectors.pdf
plot of one eigenvector, computed with CULA, matlab and NR
msommer

Posts: 2
Joined: Mon Apr 19, 2010 12:30 am

msommer,

I checked your data using the properties of eigenvalues for tridiagonal (T) matrices: T = V * E * V'

Doing this, I calculated a different answer than MATLAB, but an equally (if not more) valid solution. The function culaEig is a MEX wrapper to culaSteqr.

Code: Select all
`>> Ein = [ ... your data ... ];>> Din = [ ... your data ... ];>> T = diag(Din) + diag(Ein,-1) + diag(Ein,+1);>> [Vec,Eig] = eig(T);>> T_recon = Vec * Eig * Vec';>> norm(T) / norm(T_recon)ans =    1.0000`

This shows that MATLAB is finding a valid Eigen solution.

Code: Select all
`>> [culaVec,culaEig] = culaEig(T);>> culaT_recon = culaVec * culaEig * culaVec';>> norm(T) / norm(culaT_recon);ans =    1.0000`

This shows the answer produced by CULA is also a valid Eigen solution.

Code: Select all
`>> diffVec = culaVec - Vec; >> norm( diffVec )ans =    2.0000>> norm( diffVec )`

This shows that the answer produced by MATLAB and CULA are different, yet acceptable.

Hope this helps with your accuracy issues! Let me know if you have any other questions.

Kyle
kyle

Posts: 301
Joined: Fri Jun 12, 2009 7:47 pm

Dear kyle,

thank you, for your immediate reply, but I repeated your computations in matlab with "format long". The norm ratio of matlab was "norm(T) / norm(T_recon)=1.000000000000003" whereby cula's norm was "norm(T) / norm(culaT_recon)= 1.000000000017061". Ok, the differency is small, but in my physical problem, I compute definitely different results for particular eigenvalues. Comparisons with analytic solutions have shown, that the code is stable.

Yours sincerely
Malte Sommer
msommer

Posts: 2
Joined: Mon Apr 19, 2010 12:30 am 