An implementation of the dqds algorithm (positive case)

Abstract The dqds algorithm was introduced in 1994 to compute singular values of bidiagonal matrices to high relative accuracy but it may also be used to compute eigenvalues of tridiagonal matrices. This paper discusses in detail the issues that have to be faced when the algorithm is to be realized on a computer: criteria for accepting a value, for splitting the matrix, and for choosing a shift to reduce the number of iterations, as well as the relative advantages of using IEEE arithmetic when available. Ways to avoid unnecessary over/underflows are described. In addition some new formulae are developed to approximate the smallest eigenvalue from a twisted factorization of a matrix. The results of extensive testing are presented at the end. The list of contents is a valuable guide to the reader interested in specific features of the algorithm.