C$TEST LBAN C TO RUN AS A MAIN PROGRAM REMOVE NEXT LINE SUBROUTINE LBAN C*********************************************************************** C C EXAMPLE OF USE OF THE PORT PROGRAM BABS C C*********************************************************************** INTEGER N, I, IWRITE, I1MACH REAL G(3, 200), EVEC(100) N=10 DO 10 I=1,N G(1,I)=-1.0 G(2,I)=1.0 G(3,I)=-1.0 10 CONTINUE G(2,1)=-.75 G(2,N)=-.75 G(3,1)=-.5 G(1,2)=-.5 G(1,N)=-.5 G(3,N-1)=-.5 IWRITE=I1MACH(2) CALL EIGVEC(N,3,2,G,3,-1.0,EVEC,2) DO 20 I=1,N WRITE(IWRITE,21)EVEC(I) 20 CONTINUE 21 FORMAT(12H EIGENVECTOR,F16.8) STOP END SUBROUTINE EIGVEC(N,M,ML,G,IG,EVAL,EVEC,LIMIT) C C GIVEN A BANDED MATRIX PACKED INTO G WITH C N ROWS, M NONZERO DIAGONALS AND ML NONZERO DIAGONALS C ON AND BELOW THE DIAGONAL AND GIVEN AN EIGENVALUE OF THE C MATRIX IN EVAL, THIS SUBROUTINE USES INVERSE ITERATION TO C DETERMINE THE CORRESPONDING EIGENVECTOR AND RETURNS IT C IN EVEC. C LIMIT IS A BOUND ON THE NUMBER OF ITERATIONS C INTEGER N, M, ML, IG, LIMIT INTEGER I, JAL, ISTKGT, JINTER, JX, MU, IERR, NERROR INTEGER LIM, JJ, ISAMAX, JXI, IST(1000) REAL G(IG, N), EVEC(N), EVAL REAL SIZE, R1MACH, EPS, SC, BET, D1, SC2, ABS REAL R(1000) DOUBLE PRECISION D(500) COMMON /CSTAK/ D EQUIVALENCE (D(1),IST(1)),(R(1),D(1)) CALL ENTER(1) C DETERMINE ITERATION TOLERANCE SIZE = BANM(N,ML,M,G,IG) EPS=SIZE*R1MACH(4) C SUBTRACT EIGENVALUE FROM DIAGONAL OF G DO 10 I=1,N G(ML,I)=G(ML,I) - EVAL 10 CONTINUE C GET SPACE FROM STACK FOR AL,INTER, AND SCRATCH VECTOR JAL =ISTKGT(N*(ML-1),3) JINTER=ISTKGT(N,2) JX=ISTKGT(N,3) C GET LU DECOMPOSITION OF MATRIX CALL BALU(N,ML,M,G,IG,R(JAL),ML-1,IST(JINTER),MU,EPS) C OBTAIN INITIAL RIGHT HAND SIDE IF (NERROR(IERR).NE.0) CALL ERROFF DO 20 I=1,N EVEC(I)=1.0 20 CONTINUE CALL BABS(N,G,IG,EVEC,N,1,MU) LIM=0 JJ=ISAMAX(N,EVEC,1) SC=1.0/EVEC(JJ) C SCALE FIRST RHS TO HAVE INFINITY NORM OF 1 CALL SSCAL(N,SC,EVEC,1) C ITERATIVE PHASE BEGINS HERE 30 LIM=LIM+1 C MAKE A COPY OF OLD APPROXIMATION CALL MOVEFR(N,EVEC,R(JX)) C GET NEW APPROXIMATION OF EIGNVECTOR CALL BAFS(N,ML,R(JAL),ML-1,IST(JINTER),EVEC,N,1) CALL BABS(N,G,IG,EVEC,N,1,MU) BET=1.0/EVEC(JJ) JJ=ISAMAX(N,EVEC,1) SC2=1.0/EVEC(JJ) C COMPUTE CONVERGENCE CRITERIA D1=0.0 DO 40 I=1,N JXI=JX-1+I D1=AMAX1(D1,ABS((R(JXI)-BET*EVEC(I))*SC2)) 40 CONTINUE SC=SC2 CALL SSCAL(N,SC,EVEC,1) C TEST FOR CONVERGENCE AND IF ITERATION LIMIT EXCEEDED IF (D1.GT.EPS.AND.LIM.LT.LIMIT) GO TO 30 CALL LEAVE RETURN END .