Quantcast Testing a Lanczos (tridiagonalization) algorithm Text - Physics Forums Library

PDA

View Full Version : Testing a Lanczos (tridiagonalization) algorithm


senorbum
Jul9-08, 10:56 AM
So I implemented a tridiagonalization algorithm, however I don't know if the result that I am receiving is accurate. Does anyone know of an easy way of testing this? I suppose it is possible to scour around for code, but it would be difficult to find matching code. So if anyone knows anything that could help me, that would be much appreciated!

Edit: One other question. If you tridiagonalize a matrix, then compute the Singular Value Decomposition, then re-compute the matrix with say 'k' singular values, isn't the resulting matrix going to also be tridiagonalized? In which case, how can you apply the SVD of a tridiagonal matrix back to the original matrix?

vkroom
Jul19-08, 08:13 PM
what u can do is test the code for a smaller matrix of the same type. and if u r trying to find the extremal eigenvalue then u can do the following :

check for convergence of the extremal eigenvalue of ur tridiagonal matrix after it has grown by 10 dimensions, i.e. find the lowest/heighest eig val when ur tridiagonal matrix is 10 dim, then chk again when its 20 dim, then 30 dim, and so forth. once the value has reached a fixed point (it generally does so by 70-80 dim) u have ur extremal eigen value.
this saves a lot of resource, especially when u r generating the matrix with a lanczos.

senorbum
Jul21-08, 08:48 AM
what u can do is test the code for a smaller matrix of the same type. and if u r trying to find the extremal eigenvalue then u can do the following :

check for convergence of the extremal eigenvalue of ur tridiagonal matrix after it has grown by 10 dimensions, i.e. find the lowest/heighest eig val when ur tridiagonal matrix is 10 dim, then chk again when its 20 dim, then 30 dim, and so forth. once the value has reached a fixed point (it generally does so by 70-80 dim) u have ur extremal eigen value.
this saves a lot of resource, especially when u r generating the matrix with a lanczos.

Thanks. I have been doing that recently. Is simple lanczos even good for more than the extremal eigenvalue? It seems that I'd have to reorthogonalize if I wanted more than one, but I'm not 100% sure. Basically I am wondering if it will preserve say the top x eigenvalues where x << matrix dimension (say like 20 in a 200x200 matrix). It seems like it may generate some spurious values, but I was wondering if you could confirm this.

vkroom
Jul27-08, 01:57 AM
the lanczos vectors lose their orthonormality rapidly. so some kind of re-orthonormalization scheme is required. but thus far i didn't have any problem with the naive lanczos iterations.
You might want to consult the book "Matrix Computations" by Golub. Its pretty good.

the bottomline is that for practical and reliable calculations its better to use some ready made routine, rather than code one for urself. there are a few fine tuning required to produce a robust lanczos routine, and u mite wanna use some readymade ones for the purpose.

one such is ARPACK. this one is coded in fortran. there's a c++ inrterface : ARPACK++


good luck

senorbum
Jul28-08, 08:49 AM
the lanczos vectors lose their orthonormality rapidly. so some kind of re-orthonormalization scheme is required. but thus far i didn't have any problem with the naive lanczos iterations.
You might want to consult the book "Matrix Computations" by Golub. Its pretty good.

the bottomline is that for practical and reliable calculations its better to use some ready made routine, rather than code one for urself. there are a few fine tuning required to produce a robust lanczos routine, and u mite wanna use some readymade ones for the purpose.

one such is ARPACK. this one is coded in fortran. there's a c++ inrterface : ARPACK++


good luck

Well, my portion of this project is done, but thanks for the response. I used a naive version of Lanczos, but will need somebody who has a better understanding of linear algebra to implement a reorthogonalization technique. The reason I couldn't use a pre-written routine is because I was coding it using CUDA, which is NVidia's GPU programming language. If I had more time I might have been able to figure out how ARPACK has theirs coded and replicate it, but my time is up within a week and a half, so I have more pressing issues to attend to (The project is for a summer internship, hence the random cutoff time).

Anyways, thanks again for the response.