C FIGURE 0.4.4
      SUBROUTINE LBAND(A,X,B,N,L)                                    
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
C                                                                       
C     ARGUMENT DESCRIPTIONS                                             
C                                                                       
C     A - (INPUT) A IS AN N BY 2*L+1 ARRAY CONTAINING THE BAND MATRIX
C         IN BAND STORAGE MODE.  A(I,J-I) CONTAINS THE ELEMENT IN 
C         ROW I, COLUMN J OF THE BAND MATRIX.
C     X - (OUTPUT) X IS THE SOLUTION VECTOR OF LENGTH N.                
C     B - (INPUT) B IS THE RIGHT HAND SIDE VECTOR OF LENGTH N.          
C     N - (INPUT) N IS THE NUMBER OF EQUATIONS AND NUMBER OF UNKNOWNS        
C         IN THE LINEAR SYSTEM.                                         
C     L - (INPUT) L IS THE HALF-BANDWIDTH, DEFINED AS THE MAXIMUM VALUE 
C         OF ABS(I-J) SUCH THAT AIJ IS NONZERO.                         
C                                                                       
      DIMENSION A(N,-L:L),X(N),B(N),A2(N,-L:2*L)
C                              COPY B TO X
      DO 5 I=1,N
      X(I) = B(I)
      DO 5 J=-L,2*L
C                              COPY A TO A2, WHICH HAS L EXTRA
C                              DIAGONALS ABOVE, TO HOLD FILL-IN
         IF (J.LE.L) THEN
            A2(I,J) = A(I,J)
         ELSE
            A2(I,J) = 0.0
         ENDIF
    5 CONTINUE
C                              BEGIN FORWARD ELIMINATION                       
      DO 40 K=1,N-1
         IBIG = K
         BIG = ABS(A2(K,0))
C                              FIND THE LARGEST POTENTIAL PIVOT
         DO 10 I=K,MIN(K+L,N)
            IF (ABS(A2(I,K-I)).GT.BIG) THEN
               BIG = ABS(A2(I,K-I))
               IBIG = I
            ENDIF                                                              
   10    CONTINUE
         IF (BIG.EQ.0.0) GO TO 70
C                              SWITCH ROW K WITH THE ROW (IBIG) 
C                              CONTAINING THE LARGEST POTENTIAL PIVOT
         DO 15 J=K,MIN(K+2*L,N)
            TEMP = A2(IBIG,J-IBIG)
            A2(IBIG,J-IBIG) = A2(K,J-K)
            A2(K,J-K) = TEMP
   15    CONTINUE
         TEMP = X(IBIG)
         X(IBIG) = X(K)
         X(K) = TEMP
         DO 30 I=K+1,MIN(K+L,N)
            AMUL = -A2(I,K-I)/A2(K,0)
            IF (AMUL.EQ.0.0) GO TO 30
C                              ADD AMUL TIMES ROW K TO ROW I
            DO 20 J=K,MIN(K+2*L,N)
               A2(I,J-I) = A2(I,J-I) + AMUL*A2(K,J-K)
   20       CONTINUE
            X(I) = X(I) + AMUL*X(K)
   30    CONTINUE                                                           
   40 CONTINUE
      IF (A2(N,0).EQ.0.0) GO TO 70
C                              BEGIN BACK SUBSTITUTION
      X(N) = X(N)/A2(N,0)
      DO 60 K=N-1,1,-1
         DO 50 J=K+1,MIN(K+2*L,N)
            X(K) = X(K) - A2(K,J-K)*X(J)
   50    CONTINUE
         X(K) = X(K)/A2(K,0)
   60 CONTINUE
      RETURN                                                                   
C                              MATRIX IS SINGULAR
   70 PRINT 80
   80 FORMAT (' THE MATRIX IS SINGULAR')
      RETURN                                                                   
      END