BLSVD-2: SPARSE SVD VIA HYBRID BLOCK LANCZOS PROCEDURE (VER. 1.0) FOR EIGENSYSTEMS OF THE FORM A'A BLSVD IS A PROGRAM WRITTEN IN FORTRAN 77 DESIGNED TO COMPUTE SINGULAR VALUES AND SINGULAR VECTORS OF A LARGE SPARSE MATRIX A. THE SINGULAR VALUES (AND SINGULAR VECTORS) OF A ARE DETERMINED BY THE EIGENVALUES (AND EIGENVECTORS) OF THE MATRIX B, WHERE B = A'A. THE EIGENVALUES OF B ARE THE SQUARES OF THE SINGULAR VALUES OF A, THE EIGENVECTORS CORRESPOND TO THE RIGHT SINGULAR VECTORS ONLY. THE LEFT SINGULAR VECTORS OF A ARE THEN DETERMINED BY U = 1/SIGMA A*V, WHERE {U,SIGMA,V} IS A SINGULAR TRIPLET OF A. THIS PARTICULAR IMPLEMENTATION IS DISCUSSED IN "MULTIPROCESSOR SPARSE SVD ALGORITHMS AND APPLICATIONS", PH.D. THESIS BY M. BERRY, UNIVERSITY OF ILLINOIS AT URBANA-CHAMPAIGN, OCTOBER 1990. THIS "HYBRID" BLOCK LANCZOS PROCEDURE CONSISTS OF FIVE PHASES: PHASE 1: BLOCK LANCZOS OUTER ITERATION TO YIELD A SYMMETRIC BLOCK TRIDIAGONAL MATRIX S WHICH SHARES THE SAME EIGENVALUES OF THE MATRIX B=A'A. TOTAL OR COMPLETE RE-ORTHOGONALIZATION IS USED HERE. PHASE 2: LANCZOS METHOD FOR TRIDIAGONALIZING THE S MATRIX FROM PHASE 1 TO YIELD THE TRIDIAGONAL MATRIX T WHICH PRESERVES THE SAME EIGENVALUES. COMPLETE OR TOTAL RE-ORTHOGONAL- IZATION IS USED FOR THIS LANCZOS RECURSION. A POINT LANCZOS METHOD IS USED IF A BLOCKSIZE (NB) OF 1 IS ENCOUNTERED. PHASE 3: APPLY THE STANDARD QL ITERATION TO DIAGONALIZE T AND HENCE PRODUCE APPROXIMATE EIGENVALUES (ARRAY ALPHA) OF MATRIX B=A'A AND HENCE THE SQUARES OF THE SINGULAR VALUES OF THE ORIGINAL SPARSE MATRIX A. PHASE 4: CONVERGENCE TEST USING A USER-SUPPLIED RESIDUAL TOLERANCE (TOL) PHASE 5: ITERATION RESTART WITH ORTHOGONAL PROJECTION WITH RESPECT TO ANY (ALL) CONVERGED EIGENVECTORS OF B=A'A. - HOW TO USE SUBROUTINE BLKLAN FOR THE SVD (CALLED BY MAIN PROGRAM) THE CALLING SEQUENCE FOR BLKLAN IS . CALL BLKLAN(N,IK,LDV,V,EIG,IC,IB,ITER,TOL,RES,MAXIT, IK0,IC0,IB0,IMEM) . THE USER SPECIFIES AS PART OF THE PARAMETER LIST: . N ... COLUMN DIMENSION OF THE SPARSE MATRIX A WHOSE SVD IS SOUGHT {INTEGER}. THIS IS THE ORDER OF THE EIGENSYSTEM FOR A'A. IK ... NUMBER OF SINGULAR TRIPLETS (SINGULAR VALUES AND THEIR CORRESPONDING SINGULAR VECTORS) DESIRED {INTEGER}. IB ... INITIAL BLOCK SIZE FOR OUTER ITERATION {INTEGER}. IC ... UPPER BOUND FOR DIMENSION OF KRYLOV SUBSPACE GENERATED VIA OUTER ITERATION {INTEGER}. IC IS THE MAXIMUM DIMENSION FOR THE BLOCK UPPER BIDIAGONAL MATRIX S GENERATED IN PHASE 1 ABOVE. LDV ... LEADING DIMENSION OF ARRAY V WHICH CONTAINS APPROXIMATE RIGHT SINGULAR VECTORS ON OUTPUT {INTEGER}. LDV MUST BE GREATER THAN OR EQUAL TO N. TOL ... USER-SPECIFIED TOLERANCE FOR APPROXIMATE SINGULAR TRIPLETS {REAL*8}. MAXIT ... MAXIMUM NUMBER OF OUTER ITERATIONS ALLOWED {INTEGER}. . BLKLAN RETURNS VIA ITS PARAMETER LIST THE FOLLOWING ITEMS: . ITER ... NUMBER OF OUTER ITERATIONS REQUIRED {INTEGER}. IMEM ... ESTIMATE OF STORAGE NEEDED (IN BYTES) {INTEGER}. IB0 ... FINAL BLOCK SIZE USED IN OUTER ITERATION {INTEGER}. IC0 ... LAST BOUND FOR DIMENSION OF KRYLOV SUBSPACE USED WITHIN OUTER ITERATION {INTEGER}. IK0 ... NUMBER OF SINGULAR TRIPLETS APPROXIMATED {INTEGER}. EIG ... ONE-DIMENSIONAL ARRAY CONTAINING THE IK0 APPROXIMATE EIGENVALUES OF A'A {REAL*8}, AND HENCE THE SQUARES OF THE SINGULAR VALUES OF THE ORIGINAL SPARSE MATRIX A. V ... TWO-DIMENSIONAL ARRAY CONTAINING THE IK0 APPROXIMATE EIGENVECTORS CORRESPONDING TO THE APPROXIMATE EIGENVALUES IN ARRAY ALPHA {REAL*8}. RES ... ONE-DIMENSIONAL ARRAY CONTAINING THE IK0 RESIDUALS OF THE APPROXIMATE EIGENPAIRS OF A'A {REAL*8}. . SPARSE MATRIX-VECTOR MULTIPLICATION SUBROUTINES OPA, OPM, AND OPN MUST BE SUPPLIED BY THE USER. . THE CALLING SEQUENCE FOR OPA (ONLY USED IN MAIN PROGRAM) IS . CALL OPA(X,Y) . OP TAKES AN N-VECTOR X AND RETURNS THE M-VECTOR Y = A*X, WHERE A IS AN M BY N SPARSE MATRIX. IN OUR TEST PROGRAM, BLS2, WE EMPLOY THE HARWELL-BOEING SPARSE MATRIX FORMAT FOR ACCESSING ELEMENTS OF THE M BY N SPARSE MATRIX A AND ITS TRANSPOSE. OTHER SPARSE MATRIX FORMATS CAN BE USED, OF COURSE. . SUBROUTINE OPA IS SPECIFICALLY USED TO DERIVE APPROXIMATE LEFT SINGULAR VECTORS (MATRIX U ABOVE) FROM THE RIGHT SINGULAR VECTORS {MATRIX V} SOLELY OBTAINED BY SUBROUTINE BLKLAN (SEE MAIN PROGRAM). . THE CALLING SEQUENCE FOR OPN IS . CALL OPN(N,X,Y) OPT TAKES AN N-VECTOR X AND RETURNS THE N-VECTOR Y = A'A*X, WHERE A' DENOTES THE TRANSPOSE OF THE SPARSE MATRIX A. . THE CALLING SEQUENCE FOR OPM IS . CALL OPM(N,NC,X,LDX,Y,LDY) OPM TAKES A BLOCK OF N-VECTORS STORED IN THE TWO-DIMENSIONAL ARRAY X, AND RETURNS THE BLOCK OF N-VECTORS Y = A'A*X, WHERE A' DENOTES THE TRANSPOSE OF THE SPARSE MATRIX A. HERE, NC IS THE NUMBER OF COLUMNS OF THE TWO-DIMENSIONAL ARRAYS X AND Y WHOSE LEADING DIMENSIONS ARE LDX AND LDY, RESPECTIVELY. . IN ADDITION, BLKLAN REQUIRES THE FOLLOWING USER-SUPPLIED PARAMETERS FOR TEMPORARY ARRAY ALLOCATIONS. . USER-SUPPLIED PARAMETERS FOR BLKLAN: . USER SHOULD SET IB2 = INITIAL BLOCK SIZE (IB) IK2 = NUMBER OF DESIRED SINGULAR VALUES (IK) IC2 = INITIAL KRYLOV SUBSPACE DIMENSION (IC) N2 = COL DIMENSION OF SPARSE MATRIX A ( N) . SAMPLE DECLARATION FROM BLS2 FILE FOR LAI7 DATASET: . PARAMETER(IB2=10, IK2=10, IC2=50, N2=82) . THIS VERSION OF BLKLAN IS DESIGNED TO APPROXIMATE THE IK-LARGEST SINGULAR TRIPLETS OF A. USERS INTERESTED IN THE IK-SMALLEST SINGULAR TRIPLETS NEED ONLY REMOVE LOOPS {740, 741, 742, 743} FROM THE BLS2 SOURCE. FOLLOWING THE LINE CALL TQL2(IC2,NN,ALPHA,BETA,QQ,INFO) THE COLUMNS OF THE TWO-DIMENSIONAL ARRAY QQ WHICH ARE USED TO OBTAIN APPROXIMATE EIGENPAIRS OF THE MATRIX A'A CORRESPOND TO THE ELEMENTS OF ALPHA (IN ASCENDING ORDER) WHICH ARE APPROXIMATE EIGENVALUES OF A'A (AND HENCE SQUARES OF THE SINGULAR VALUES OF THE MATRIX A). . EXPLICIT CALLS TO THE UNIX TIMER "ETIME" IS MADE IN THE BLS2 SOURCE. A CALL TO ANOTHER TIMING ROUTINE (ELAPSED CPU TIME) MADE BE NEEDED. . IF THE PARAMETER "VECTORS" (READ FROM BLI1) IS SET TO "TRUE", THE UNFORMATTED OUTPUT FILE BLO3 WILL CONTAIN THE APPROXIMATE SINGULAR VECTORS WRITTEN IN THE ORDER U[1], V[1], U[2], V[2], ..., U[IK0], V[IK0]. HERE U[I] AND V[I] DENOTE THE LEFT AND RIGHT SINGULAR VECTORS, RESPECTIVELY, CORRESPONDING TO THE I-TH APPROXIMATE SINGULAR VALUE, SING[I]. - INFORMATION PLEASE ADDRESS ALL QUESTIONS, COMMENTS, OR CORRECTIONS TO: M. W. BERRY DEPARMENT OF COMPUTER SCIENCE UNIVERSITY OF TENNESSEE 107 AYRES HALL KNOXVILLE, TN 37996-1301 EMAIL: BERRY@CS.UTK.EDU PHONE: (615) 974-5067