1. OVERVIEW OF THE PACKAGE Below is a list giving a brief description of the files found in the package. ============================================================================ Readme : This description. Main routines ------------- lanbpro.m : Lanczos bidiagonalization with partial reorthogonalization (BPRO). lansvd.m : Sparse SVD routine built on BPRO. lanpro.m : Lanczos tridiagonalization with partial reorthogonalization (PRO). laneig.m : Sparse symmetric eigenvalue routine built on PRO. Similar to the "lanso" subroutine in the LANSO package. Computational routines ---------------------- refinebounds.m: Refines the error bounds based on the gap structure of the Ritz values. compute_int.m : Routine for determining which Lanczos vectors to purge in a reorthogonalization. reorth.m : Iterated (modified or classical) Gram-Schmidt reorthogonalization. reorth.mex<...> : MEX file using a FORTRAN 77 version of reorth.m. tqlb.m : Computes the eigenvalues and the top and bottom elements of the eigenvectors for a symmetric tridiagonal matrix. tqlb.mex<...> : The corresponding MEX files call a modified version of the EISPACK routine TQL1 (see below). bdsqr.m : Computes the singular values and the bottom elements of the left singular vectors for a lower bidiagonal matrix. bdsqr.mex<...>: The corresponding MEX files call the LAPACK routine DBDSQR. Test scripts and test problems ------------------------------ test.m : Main script for testing sparse SVD calculations. testtqlb.m : Script for comparing different versions of TQLB. Afunc.m : Function that computes Y <- A*X. Atfunc.m : Function that computes Y <- A'*X. AtAfunc.m : Function that computes Y <- A'*A*X. Cfunc.m : Function that computes Y <- [A*X(m+1:m+n,:); A'*X(1:m,:)] i.e. Y = [ 0 A ] * X [ A' 0 ]. helio.mat : Testproblem from helioseismology. It is a full 212 x 100 matrix with condition number 2.9*10^12. Harwell-Boeing/ : Directory with various testmatrices from the Harwell-Boeing collection (in Matrix Market format). The directory also contains Matlab function for querying, reading and writing matrices stored in Matrix Market file format. ============================================================================ The routine TQLB is used for computing the eigenvalues and top and bottom elements of the eigenvectors of the Lanczos tridiagonal, which enter into the error bounds. TQLB is based on an optimized version of the EISPACK routine TQL1 (see below), and for the sake of speed FORTRAN versions of TQLB, the reorthogonalization routine REORTH and the LAPACK routine DBDSQR have been included as MEX files. The current versions of the MEX files have been compiled for SGI64, SUN and Linux x86 machines. Purely Matlab based versions are included for use on other platforms, but these can be rather slow in comparison. 2. IMPLEMENTATION DETAILS 2.1 Modifications to tqlb relative to the version used in the lanso package: In TQLB I have added a little extra code to compute the bottom elements of the normalized eigenvectors. In addition I have replaced the calls to PYTHAG in the inner loop of the QL iteration with in-line code according to the suggestions in A. Cline and J. Meyering, "Converting EISPACK to run efficiently on a vector processor", Tech. Memo., Pleasant Valley Software, Austin TX, 1989. (see their code in http://www.netlib.org/eispack/3090vf/double/tql1.f). Apart from decreasing the execution time by a factor of 3 compared to the original version, it also causes the small eigenvalues to be computed more accurately. As a test I computed the eigenvalues of the standard testmatrix: ( 2 -1 ) ( -1 2 -1 ) T = ( -1 2 . ) , ( . . -1 ) ( -1 2 ) where the true eigenvalues are lambda_i = 4 * cos(pi/2 * i / n+1)^2, i = 1,2,...,n . Below are the results (T=time in seconds, E=max relative error) obtained on an SGI Origin 200 with a 180MHz MIPS R10000 CPU, 32Kbytes primary, 1Mb secondary cache: Original TQLB Modified TQLB MATLAB EIG (eigenvalues only) ----------------------------------------------------------------------------- n T E T E T E 100 0.020 1.4e-12 0.007 3.1e-13 0.024 6.2e-14 1000 1.68 1.2e-9 0.55 2.9e-10 1.98 1.5e-10 10000 154.1 1.9e-6 52.2 8.2e-8 181.2 3.0e-8 ============================================================================== Rasmus Munk Larsen, Stanford University, 2000