Eigenvalues of a tridiagonal matrix

[petr1243]

Could some one please explain the logic behind the subroutine tqli(d,e,z) found on pg.1228 of numerical recipes in fortran. For this subroutine, d and e are one dimensional arrays and z is an optional multi-dimensional array(only used if also seeking for eigenvectors of matrix). The 3 modules, and the whole jacobi rotation technique has me confused. Thank you very much.

 

[mda]

It looks like you're looking at the F90 version...
Have you looked at the online F77 version that provides some explanation?

 

[petr1243]

Are you referring to the explanation from indiana university? If it isn't, then could you please send me the link.

 

[Maxwell]

Can you post the contents of the subroutine?

I've actually worked on this problem before - in C, however - but maybe I can help provide some input.

 

[petr1243]

I'm actually working on a tridiagonal matrix in which the diagonal term is some formula in terms of integer n, and the super and sub diagonal terms are equivalent to one another(constant, independent of n). I would like to implement this subroutine into my program. I would prefer an explanation in f90, but if you can explain it in f77 then that's fine. Here is the f77 program: http://www.physics.louisville.edu/help/nr/bookfpdf/f11-3.pdf (on pg.474)

I am not too skilled in f77, so could you please be explicit in your explanation. Thank you for your time.
Here is the subroutine in f90, where the 3 modules are easily accessible on the net:

SUBROUTINE tqli(d,e,z)

USE nrtype; USE nrutil, ONLY : assert_eq, nrerror
USE nr, ONLY : pythag

IMPLICIT NONE

REAL(sp), DIMENSION(:), INTENT(inout) :: d, e
!REAL(sp), DIMENSION(:,:), OPTIONAL, INTENT(inout) :: z

INTEGER(i4b) :: i, iter, l, m, n, ndum
REAL(sp) :: b, c, dd, f, g, p, r, s
REAL(sp), DIMENSION(SIZE(e)) :: ff

n = assert_eq(SIZE(d), SIZE(e), 'tqli: n')
!IF(PRESENT(z)) ndum=assert_eq(n, SIZE(z, 1), SIZE(z, 2), 'tqli: ndum')
e(:)=EOSHIFT(e(:),1)
DO l=1, n
iter=0
iterate: DO
DO m=l, n-1
dd=ABS(d(m)) + ABS(d(m+1))
IF(ABS(e(m))+dd == dd) EXIT
END DO
IF (m==l) EXIT iterate
IF(iter == 30) CALL nrerror('too many iterations in tqli')
iter = iter+1
g=(d(l+1)-d(l))/(2.0_sp*e(l))
r = pythag(g, 1.0_sp)
g = d(m) - d(l) + e(l) / (g + SIGN(r, g))
s = 1.0
c = 1.0
p = 0.0
DO i = m-1, l, -1
f = s*e(i)
b = c*e(i)
r = pythag(f, g)
e(i+1) = r
IF(r == 0.0) THEN
d(i+1) = d(i+1) - p
e(m) = 0.0
CYCLE iterate
END IF
s = f/r
c = g/r
g = d(i + 1) - p
r = (d(i) - g) * s + 2.0_sp * c * b
p = s * r
d(i + 1) = g + p
g = c * r - b
!IF(PRESENT(z)) THEN
!ff(1:n) = z(1:n,i+1)
!z(1:n, i+1) = s*z(1:n,i) + c*ff(1:n)
!z(1:n,i) = c*z(1:n, i) - s*ff(1:n)
!END IF
END DO
d(l) = d(l) - p
e(l) = g
e(m) = 0.0
END DO iterate
END DO
END SUBROUTINE tqli

!I have commented z out, since I'm only interested in the eigenvalues of this matrix.