LASVD-1: SPARSE SVD VIA LANSO (VER. 1.0) USING 2-CYCLIC MATRICES LANSO IS A SUBROUTINE PACKAGE WRITTEN IN FORTRAN 77 DESIGNED TO FIND SOME EIGENVALUES AND EIGENVECTORS OF A LINEAR OPERATOR OPB THAT IS REAL SYMMETRIC. THIS IS A MODIFIED VERSION OF THE LANSO CODE DISTRBUTED BY B. PARLETT AND HIS COLLEAGUES AT BERKELEY. THE OPERATOR B IS ASSUMED TO BE OF THE FORM B = [ O A ], WHERE A IS M BY N (M>>N) AND SPARSE. [ A' O ] HENCE, THE SINGULAR TRIPLETS OF A ARE COMPUTED AS THE EIGENPAIRS OF B. THE EIGENVALUES OF B ARE + OR - THE SINGULAR VALUES OF A, THE FIRST M COMPONENTS OF THE EIGENVECTORS CORRESPOND TO THE LEFT SINGULAR VECTORS OF A, AND THE LAST N COMPONENTS OF THE EIGENVECTORS OF B CORRESPOND TO THE RIGHT SINGULAR VECTORS OF A. LANSO IMPLEMENTS THE SIMPLE LANCZOS ALGORITHM WITH SIMON'S SELECTIVE ORTHOGONALIZATION TO ACTIVELY MAINTAIN EXTENDED SEMI- ORTHOGONALITY AMONGST THE COMPUTED LANCZOS VECTORS. LANSO IS MOST EFFICIENT WHEN ONLY A MODEST FRACTION (< 1/4) OF B-EIGENPAIRS ARE WANTED. ONE POTENTIALLY USEFUL FEATURE OF LANSO IS ITS DELIVERY OF CERTAIN EIGENVALUES NOT YET GOOD TO FULL PRECISION, FOR SOME TASKS THEY ARE OFTEN GOOD ENOUGH. THE LARGEST EIGENVALUES OF B TEND TO BE CAPTURED BEFORE THE SMALL ONES BUT CONVERGENCE DEPENDS ON SEPARATION RATIOS AND THESE VARY WITH THE APPLICATION. THE PROGRAM WILL FIND EIGENVALUES OF MULTIPLICITY GREATER THAN ONE BUT THE COPIES STABILIZE ONE AT A TIME, THE INTERVAL BETWEEN COPIES DEPENDS ON THE PROBLEM AND WORDLENGTH. EXPERIENCE SUGGESTS THE GAP IS BETWEEN 4 AND 20 STEPS. - HOW TO USE LANSO FOR THE SVD ONE INVOKES LANSO THROUGH THE TOP LEVEL SUBROUTINE NAMED LANDR. THE CALLING SEQUENCE FOR LANDR IS . CALL LANDR(N,MAXPRS,LANMAX,ENDL,ENDR,EV,KAPPA,NW,W,J,NEIG,RITZ,BND,ERR) . THE USER SPECIFIES AS PART OF THE PARAMETER LIST: . N ... DIMENSION OF THE EIGENPROBLEM (NROW+NCOL), WHERE A IS NROW BY NCOL. MAXPRS ... MAXIMUM NUMBER OF EIGENPAIRS OF B WANTED, (MAXPRS SHOULD BE < N/4), LANMAX ... MAXIMUM NUMBER OF LANCZOS STEPS ALLOWED, (MAXPRS <= LANMAX <= N) [ENDL,ENDR] ... 2 END-POINTS OF AN INTERVAL WITHIN WHICH ALL UNWANTED EIGENVALUES LIE, (ENDL < ENDR) (ONE CHOICE IS ENDL=-MACHEPS,ENDR=MACHEPS) EV ... {INTEGER} IF EV.LE.0 THEN ONLY EIGENVALUES AND NO EIGENVECTORS ARE WANTED, OTHERWISE EV IS TAKEN TO BE 8, AND THE OUTPUT FILE LAO8 (UNFORMATTED) WILL CONTAIN THE EIGENVECTORS OF B (SINGULAR VECTORS OF A). THE EIGENVECTORS OF B ARE WRITTEN TO LAO8 CORRESPONDING TO THE ASCENDING ORDER OF THE RITZ VALUES IN THE LAO9 OUTPUT FILE. THE FIRST NROW ELEMENTS OF EACH EIGENVECTOR IS AN EQUIVALENT LEFT SINGULAR VECTOR OF A, AND THE LAST NCOL ELEMENTS OF EACH EIGENVECTOR IS AN EQUIVALENT RIGHT SINGULAR VECTOR OF A. KAPPA ... {REAL*8} RELATIVE ACCURACY OF RITZ VALUES ACCEPTABLE AS EIGENVALUES, (KAPPA WILL BE RESET TO MAX(KAPPA,MACHEPS**(3/4)) IF EV.GT.0, OTHERWISE KAPPA IS UNTOUCHED) NW ... SIZE OF THE WORKING ARRAY, (NW >= 6*N+1+4*LANMAX IF EV.LE.0, OTHERWISE NW >= 6*N+1+4*LANMAX+LANMAX*LANMAX) AND W ... THE WORKING ARRAY ITSELF. . THE FIRST N ENTRIES OF W MUST BE INITIALIZED; IF INITIALIZED TO ALL ZEROES, A RANDOMLY CHOSEN STARTING VECTOR WILL BE USED, OTHERWISE IT WILL BE TAKEN AS A USER-SPECIFIED STARTING VECTOR FOR THE LANCZOS RUN. . LANDR RETURNS VIA ITS PARAMETER LIST THE FOLLOWING ITEMS: . J ... NUMBER OF LANCZOS STEPS ACTUALLY TAKEN, NEIG ... NUMBER OF RITZ VALUES STABILIZED(SEE KAPPA), RITZ ... ARRAY THAT HOLDS THE RITZ VALUES, BND ... ARRAY THAT HOLDS THE ERROR BOUNDS AND IERR ... AN ADVISORY ERROR FLAG. . IF EV.GT.0 THEN . IN ADDITION, LANDR RETURNS VIA FORTRAN CHANNEL EV=8 ALL THE EIGENVECTORS OF B WHICH CORRESPOND TO LEFT AND RIGHT SINGULAR VECTORS OF A. THE I-TH LEFT SINGULAR VECTOR IS WRITTEN FIRST, AND THE I-TH RIGHT SINGULAR VECTOR FOLLOWS: (U(K,I),K=1,NROW), (V(K,I),K=1,NCOL). NOTE THAT (OUTPUT) CHANNEL EV=8 IS OPEN 'UNFORMATTED' TO SAVE DISK SPACE AND PRESERVE ACCURACY. 1 SUBROUTINE NAMED OPB MUST BE SUPPLIED BY THE USER. . THE CALLING SEQUENCE FOR OPB IS . CALL OPB(N,X,Y) . OPB TAKES AN N-VECTOR X AND RETURNS Y = B*X AS DESCRIBED ABOVE. OPB WILL BE CALLED WITH AN N ALWAYS EQUAL TO THE DIMENSION OF THE EIGENPROBLEM (NROW+NCOL). IN OUR TEST PROGRAM, LAS1, WE EMPLOY THE HARWELL-BOEING SPARSE MATRIX FORMAT FOR ACCESSING ELEMENTS OF THE NROW BY NCOL SPARSE MATRIX A AND ITS TRANSPOSE. OTHER SPARSE MATRIX FORMATS CAN BE USED, OF COURSE. IF SECONDARY STORAGE (SUCH AS A DISK DRIVE) MUST BE USED TO STORE THE COMPUTED LANCZOS VECTORS, THE SUBROUTINE NAMED STORE MUST BE SUPPLIED BY THE USER IN PLACE OF OUR TRIVIALLY IMPLEMENTED STORE WHICH ASSUMES THE EXISTENCE OF VIRTUAL MEMORY. THE CALLING SEQUENCE FOR STORE IS . CALL STORE(N,ISW,J,S) . WHERE . N ... THE DIMENSION OF THE EIGENPROBLEM (NROW+NCOL), ISW ... ACTION TYPE SWITCH, CAN ONLY BE 1, 2, 3 OR 4. ISW = 1 REQUESTS STORING OF J-TH LANCZOS VECTOR Q(J), ISW = 2 REQUESTS RETRIEVAL OF J-TH LANCZOS VECTOR Q(J), ISW = 3 REQUESTS STORING OF M*Q(J) FOR J = 1 OR 2, ISW = 4 REQUESTS RETRIEVAL OF M*Q(J) FOR J = 1 OR 2. J ... INDEX OF THE J-TH LANCZOS VECTOR, S ... THE N-VECTOR TO BE STORED/RETRIEVED. LANDR WILL END ITS RUN IF EITHER . 1) MAXPRS EIGENPAIRS HAVE BEEN FOUND, OR 2) LANMAX LANCZOS STEPS HAVE BEEN TAKEN, OR 3) TWO RITZ VALUES INSIDE [ENDL,ENDR] HAVE STABILIZED. . WHICHEVER OCCURS FIRST. NOTE THAT AN EIGENVALUE OF MULTIPLICITY HIGHER THAN ONE MAY NOT HAVE ALL ITS COPIES DELIVERED. - ACKNOWLEDGEMENTS BERESFORD PARLETT AND HIS COLLEAGUES (DAVID SCOTT, HORST SIMON, BAHRAM NOUR-OMID, JOHN LEWIS, AND JEREMY DU CROZ) ARE TO BE CREDITED FOR THE LANSO STRATEGY AND INITIAL FORTRAN-77 CODE. - 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