C CYLWRK    SOURCE    PV        11/03/07    21:16:15     6885           
      SUBROUTINE CYLWRK(RF,ZF,KFLAG,U,UR,UZ,URR,URZ,UZZ,ITYPE,                  
     $NPN,EDICON,ICOO,ITRAVA)                                                   
C                                                                                
C      EVALUATE POTENTIAL AND ITS DERIVATIVES FOR CYLINDER COORDINATES           
C-----------------------------------------------------------------------        
C                                                                               
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)                                                  
C      POTENTIAL VALUES                                                          
      SEGMENT EDICON                                                            
        INTEGER  KSTRT, KSTEP, NMIR, IS                                         
        REAL*8  CROT, SROT, SYMFCT                                              
        LOGICAL LREAL, LIMAG                                                    
      ENDSEGMENT                                                                
      SEGMENT ICOO                                                              
        REAL*8 X(MV),Y(MV),P(MV),WNODE(MV)                                      
        INTEGER LISVO(MV)                                                       
      ENDSEGMENT                                                                
      SEGMENT ITRAVA                                                            
        REAL*8 KENN(M42,2),SIGMA(M42),DELRHO(M42),C(M50,M50)                    
        REAL*8 AK(M50),UM(M50),RM(M50)                                          
       INTEGER IL(M50)                                                          
      ENDSEGMENT                                                                
C       COMMON // NEU, NPN, NPU                                                  
C                                                                                
C      CONSTANTS DESCRIBING SYMMETRY CONDITIONS                                  
C                                                                                
C                                                                                
C      SOURCE TERMS                                                              
C                                                                                
C      TABLES AND WORKING STORE FOR POTENTIAL AND FIELD EDIT                     
C                                                                                
C-----------------------------------------------------------------------        
C                                                                                
      M50=RM(/1)                                                                
      NPU = MIN(NPN,19) 
      KEND = KSTRT + KSTEP * (NPU - IS - 2)                                     
      IFLAG = 0                                                                 
      J = 0                                                                     
C                                                                               
      DO 90 JJ = 1, NPN                                                         
C                                                                               
      IF(WNODE(JJ) .EQ. 0.0) GO TO 90                                          
      RN  = X(JJ)                                                               
      ZN  = Y(JJ)                                                               
      UN  = P(JJ)                                                               
C                                                                                
C      SET UP HARMONIC FUNCTIONS FOR LOCATION (RN,ZN)                            
C      NOTE.         U(SUBSCRIPT M) IS STORED IN UM(M+3)                         
C                    R(SUBSCRIPT M) IS STORED IN RM(M+3)                         
C                    THE DO LOOP INDEX K IS (M+1)                                
C                                                                                
      UM(3) = 1.0                                                               
      IF(IS .NE. 0) UM(3) = RN                                                  
      RM(3) = FLOAT(IS+IS) * ZN                                                 
C                                                                                
      DO 40 K = 2, KEND                                                         
      UM(K+2) = (RN*RM(K+1) + ZN*UM(K+1))                                       
      RM(K+2) = (ZN*RM(K+1) - RN*UM(K+1)) * FLOAT(K+IS+IS-1) / FLOAT(K)         
   40 CONTINUE                                                                  
C                                                                                
C      SET UP COEFFICIENT MATRIX FOR FIT                                         
C                                                                                
      J = J + 1                                                                 
      I = IS                                                                    
      C(J,1) = WNODE(JJ)                                                        
C                                                                                
      IF(KEND .LT. KSTRT) GO TO 80                                              
C                                                                                
      DO 70 K = KSTRT, KEND, KSTEP                                              
      I = I + 1                                                                 
      C(J,I) = WNODE(JJ) * UM(3) * UM(K+2)                                      
   70 CONTINUE                                                                  
C                                                                                
   80 C(J,I+1) = WNODE(JJ) * UN                                                 
C                                                                                
   90 CONTINUE                                                                  
C      WRITE(6,FMT= ' ('' I ET J MAX '',/,(2I6))') I,J                           
C      WRITE(6,FMT= ' ('' KSTRT KSTEP KEND ITYPE  '',/,(4I5))') KSTRT,           
C     $KSTEP,KEND,ITYPE                                                          
C                                                                                
C-----------------------------------------------------------------------        
C      PERFORM LEAST-SQUARES FIT                                                 
C                                                                                
      CALL GLSQ(C,AK,IL,J,I,ALFA,1.0D-10,1.0D-10,M50)                           
C                                                                                
C-----------------------------------------------------------------------        
C      COMPUTE POTENTIAL AND DERIVATIVES                                         
C                                                                                
      UM(1) = 0.0                                                               
      UM(2) = 0.0                                                               
      UM(3) = 1.0                                                               
      IF(IS .NE. 0) UM(3) = RF                                                  
      RM(1) = 0.0                                                               
      RM(2) = FLOAT(IS)                                                         
      RM(3) = 2.0 * RM(2) * ZF                                                  
C                                                                               
      DO 110 K = 2, KEND                                                        
      UM(K+2) = (RF*RM(K+1) + ZF*UM(K+1))                                       
      RM(K+2) = (ZF*RM(K+1) - RF*UM(K+1)) * FLOAT(K+IS+IS-1) / FLOAT(K)         
  110 CONTINUE                                                                  
C                                                                               
      IFLAG = KFLAG                                                             
      U   = AK(1) * FLOAT(IS)                                                   
      UR  = 0.0                                                                 
      UZ  = 0.0                                                                 
      URZ = 0.0                                                                 
      UZZ = 0.0                                                                 
      I   = IS                                                                  
C                                                                                
      IF(KEND .LT. KSTRT) GO TO 150                                             
C                                                                                
      DO 120 K = KSTRT, KEND, KSTEP                                             
      I = I + 1                                                                 
      U   = U   + AK(I) * UM(3) * UM(K+2)                                       
      XK  = AK(I) * FLOAT(K+IS+IS-1)                                            
      UR  = UR  + XK * RM(K+1)                                                  
      UZ  = UZ  + XK * UM(K+1)                                                  
      XK  = XK * FLOAT(K+IS+IS-2)                                               
      URZ = URZ + XK * RM(K)                                                    
      UZZ = UZZ + XK * UM(K)                                                    
  120 CONTINUE                                                                  
C                                                                                
  150 KFLAG = IFLAG                                                             
      URR = - UZZ                                                               
      RETURN                                                                    
C***********************************************************************        
      END