C FIGURE 5.4.1
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      EXTERNAL FUNA                                                     
C                              N_ = NUMBER OF SUBINTERVALS           
C                              M_ = NUMBER OF UNKNOWNS
      PARAMETER (N_=4,M_=2*N_+1)                                       
C                              BREAK POINTS ARE XPTS(K), K=0,N          
      COMMON /BLKX/ N,XPTS(0:N_)                                     
C                              UNKNOWNS ARE A(I), I=1,M                 
      COMMON /BLKA/ M,A(M_)                                      
      DIMENSION AMAT(M_,-3:3),B(M_)
      N = N_                                                         
      M = M_                                                     
      ALPHA = 3.636                                                     
      DO 5 K=0,N                                                        
         XPTS(K) = (DBLE(K)/N)**ALPHA                                   
    5 CONTINUE                                                          
C                              CALCULATE RIGHT HAND SIDE VECTOR AND     
C                              COEFFICIENT MATRIX                       
      L = 3
      DO 15 K=1,M                                                       
         B(K) = 1.1D0*PHI(K,0,1.D0)                                     
         DO 10 I=MAX(K-L,1),MIN(K+L,M)
            AMAT(K,I-K) = GAUSS(FUNA,K,I)
   10    CONTINUE                                                       
   15 CONTINUE                                                          
C                              SOLVE LINEAR SYSTEM USING BAND SOLVER    
      CALL LBAND(AMAT,A,B,M,L)                                        
C                              CALCULATE MAXIMUM ERROR                  
      ERMAX = 0.0                                                       
      DO 20 J=0,100                                                     
         X = J/100.D0                                                   
         ERR = ABS(USOL(0,X) - X**1.1D0)
         ERMAX = MAX(ERMAX,ERR)                                         
   20 CONTINUE                                                          
      PRINT 25, ERMAX                                                   
   25 FORMAT (E15.5)                                                    
      STOP                                                              
      END                                                               
      FUNCTION FUNA(X,K,I)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
C                              INTEGRAND FOR (K,I)-TH ELEMENT OF        
C                              COEFFICIENT MATRIX
      PK  = PHI(K,0,X)
      PKx = PHI(K,1,X)
      PI  = PHI(I,0,X)
      PIx = PHI(I,1,X)
      FUNA = PKx*PIx + 0.11D0/X**2*PK*PI
      RETURN                                                            
      END                                                               
      SUBROUTINE LOCATE(I,HORS,K)                                       
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      CHARACTER*1 HORS                                                  
C                              DETERMINE THE TYPE ('H' OR 'S') OF BASIS 
C                              FUNCTION I, AND THE POINT (XPTS(K)) ON   
C                              WHICH IT IS CENTERED                     
      IF (MOD(I,2).EQ.0) HORS = 'H'                                     
      IF (MOD(I,2).EQ.1) HORS = 'S'                                     
      K = I/2                                                           
      RETURN                                                            
      END                                                               
      FUNCTION OMEGA(IDER,X)                                            
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
C                              EVALUATE OMEGA (IDER=0) OR ITS IDER-TH   
C                              DERIVATIVE AT X                          
      OMEGA = 0.0                                                       
      RETURN                                                            
      END