C FIGURE 5.3.3b
      FUNCTION USOL(IDER,X)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      COMMON /BLKA/ M,A(1)                                              
C                              EVALUATE THE APPROXIMATE SOLUTION (IDER=0)
C                              OR ITS IDER-TH DERIVATIVE AT X           
      USOL = OMEGA(IDER,X)
      DO 10 I=1,M
         USOL = USOL + A(I)*PHI(I,IDER,X)                               
   10 CONTINUE                                                          
      RETURN                                                            
      END                                                               
      FUNCTION PHI(I,IDER,X)                                            
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      CHARACTER*1 HORS                                                  
C                              EVALUATE THE I-TH BASIS FUNCTION (IDER=0)
C                              OR ITS IDER-TH DERIVATIVE AT X           
      CALL LOCATE(I,HORS,K)                                             
      PHI = HERM(HORS,K,IDER,X)                                         
      RETURN                                                            
      END                                                               
      FUNCTION HERM(HORS,K,IDER,X)                                      
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      CHARACTER*1 HORS                                                  
      COMMON /BLKX/ N,XPTS(0:1)                                         
C                              EVALUATE AT X THE IDER-TH DERIVATIVE OF  
C                              THE CUBIC HERMITE "H" OR "S" FUNCTION    
C                              CENTERED AT XPTS(K).                     
      XK = XPTS(K)                                                      
      HH = XPTS(N)-XPTS(0)                                              
      IF (K.GT.0) XKM1 = XPTS(K-1)                                      
      IF (K.LE.0) XKM1 = XPTS(0)-HH                                     
      IF (K.LT.N) XKP1 = XPTS(K+1)                                      
      IF (K.GE.N) XKP1 = XPTS(N)+HH                                     
      HERM = 0.0                                                        
      IF (X.LT.XKM1) RETURN                                             
      IF (X.GT.XKP1) RETURN                                             
      IF (X.LE.XK) THEN                                                 
         H = XK-XKM1                                                    
         S = (X-XKM1)/H                                                 
         IF (HORS.EQ.'H') THEN                                          
            IF (IDER.EQ.0) HERM = 3.0*S**2 - 2.0*S**3                   
            IF (IDER.EQ.1) HERM = 6.0*S/H - 6.0*S**2/H                  
            IF (IDER.EQ.2) HERM = 6.0/H**2 - 12.0*S/H**2                
         ELSE                                                           
            IF (IDER.EQ.0) HERM = -H*S**2 + H*S**3                      
            IF (IDER.EQ.1) HERM = -2.0*S + 3.0*S**2                     
            IF (IDER.EQ.2) HERM = -2.0/H + 6.0*S/H                      
         ENDIF                                                          
      ELSE                                                              
         H = XKP1-XK                                                    
         S = (XKP1-X)/H                                                 
         IF (HORS.EQ.'H') THEN                                          
            IF (IDER.EQ.0) HERM = 3.0*S**2 - 2.0*S**3                   
            IF (IDER.EQ.1) HERM = -6.0*S/H + 6.0*S**2/H                 
            IF (IDER.EQ.2) HERM = 6.0/H**2 - 12.0*S/H**2                
         ELSE                                                           
            IF (IDER.EQ.0) HERM = H*S**2 - H*S**3                       
            IF (IDER.EQ.1) HERM = -2.0*S + 3.0*S**2                     
            IF (IDER.EQ.2) HERM = 2.0/H - 6.0*S/H                       
         ENDIF                                                          
      ENDIF                                                             
      RETURN                                                            
      END                                                               
      FUNCTION GAUSS(F,K,I)                                             
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)                                
      CHARACTER*1 HORS                                                  
      DIMENSION S(3),W(3)                                               
      COMMON /BLKX/ N,XPTS(0:1)                                         
C                             USE THE GAUSS 3 POINT FORMULA TO          
C                             INTEGRATE F(X,K,I).                           
      DATA S                                                            
     & /0.1127016653792583D0,0.5000000000000000D0,0.8872983346207417D0/ 
      DATA W                                                            
     & /0.2777777777777778D0,0.4444444444444444D0,0.2777777777777778D0/ 
      GAUSS = 0.0                                                       
C                             THE FUNCTION F INVOLVES PHI(K,...), WHICH 
C                             IS NONZERO ONLY FROM XPTS(IK-1) TO        
C                             XPTS(IK+1)                                
      CALL LOCATE(K,HORS,IK)                                            
      DO 10 J=MAX(IK,1),MIN(IK+1,N)                                     
         H = XPTS(J)-XPTS(J-1)                                          
         DO 5 JX=1,3                                                    
            X = XPTS(J-1)+S(JX)*H                                       
            GAUSS = GAUSS + H*W(JX)*F(X,K,I)
    5    CONTINUE                                                       
   10 CONTINUE                                                          
      RETURN                                                            
      END