Sources Esope | Cast3M

carwrk
C CARWRK    SOURCE    CB215821  17/07/21    21:15:01     9513           
      SUBROUTINE CARWRK(XF,YF,KFLAG,U,UX,UY,UXX,UXY,UYY,ITYPE,
     $ XORG,YORG,NPN,EDICON,ICOO,ITRAVA)
C
C     EVALUATE POTENTIAL AND ITS FIRST AND SECOND DERIVATIVES
C-----------------------------------------------------------------------
C
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
C-----------------------------------------------------------------------
C
      COMPLEX*16 Z, ZE, ZN, ZP, ZPP
C
C-----------------------------------------------------------------------
      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
C      COMMON // NEU, NPN, NPU
C
C     TABLES AND WORKING STORE FOR POTENTIAL AND FIELD EDIT
C
C
C     SET UP EQUATIONS FOR LEAST-SQUARES FIT
C
      M50=RM(/1)
      NPU = NPN
      NTERM =MIN(( NPU - 4),15)
      IF(LREAL .AND. LIMAG) NTERM = MIN( ((NPU - 5) / 2 ),7)
      KEND = KSTRT + NTERM * KSTEP
      XE = XF - XORG
      YE = YF - YORG
      ZE = CMPLX(XE,YE)
      IFLAG = 0
      J = 0
C
      DO 90 JJ = 1, NPN
C
      IF(WNODE(JJ) .EQ. 0.0) GO TO 90
      XN = X(JJ) - XORG
      YN = Y(JJ) - YORG
      ZN = CMPLX(XN,YN)
      U = P(JJ)
C
C     LOGARITHMIC TERMS DUE TO FILAMENTS
C     TERMS DUE TO OTHER SOURCES
C
C
C     COEFFICIENT MATRIX
C
      J = J + 1
      C(J,1) = WNODE(JJ)
      I = 1
C
      IF(KEND .LT. KSTRT) GO TO 80
C
      DO 70 KI = KSTRT, KEND, KSTEP
      Z = ZN**KI
C
      IF(.NOT. LREAL) GO TO 60
      I = I + 1
      C(J,I) = WNODE(JJ) * REAL(Z)
C
   60 IF(.NOT. LIMAG) GO TO 70
      I = I + 1
      C(J,I) = WNODE(JJ) * AIMAG(Z)
C
   70 CONTINUE
C
   80 C(J,I+1) = WNODE(JJ) * U
C
   90 CONTINUE
C
C-----------------------------------------------------------------------
C     PERFORM THE LEAST SQUARES FIT
C
      CALL GLSQ(C,AK,IL,J,I,ALFA,1.0D-10,1.0D-10,M50)
C
C-----------------------------------------------------------------------
C      PERFORM THE ANALYTIC DIFFERENTIATION
C
      IFLAG = KFLAG
      U   = AK(1)
      UX  = 0.0
      UY  = 0.0
      UXX = 0.0
      UXY = 0.0
      UYY = 0.0
C
C      HARMONIC TERMS
C
      IF(KEND .LT. KSTRT) GO TO 150
      I = 1
C
      DO 140 KI = KSTRT, KEND, KSTEP
      Z = ZE**KI
      ZP = (1.0,0.0)
      IF(KI .GT. 1) ZP = FLOAT(KI) * ZE ** (KI-1)
      ZPP = (0.0,0.0)
      IF(KI .EQ. 2) ZPP = (2.0,0.0)
      IF(KI .GT. 2) ZPP = FLOAT(KI*(KI-1)) * ZE ** (KI-2)
C
      IF(.NOT. LREAL) GO TO 120
      I = I + 1
      U   = U   + AK(I) * REAL(Z)
      UX  = UX  + AK(I) * REAL(ZP)
      UY  = UY  - AK(I) * AIMAG(ZP)
      UXX = UXX + AK(I) * REAL(ZPP)
      UXY = UXY - AK(I) * AIMAG(ZPP)
C
  120 IF(.NOT. LIMAG) GO TO 140
      I = I + 1
      U   = U   + AK(I) * AIMAG(Z)
      UX  = UX  + AK(I) * AIMAG(ZP)
      UY  = UY  + AK(I) * REAL(ZP)
      UXX = UXX + AK(I) * AIMAG(ZPP)
      UXY = UXY + AK(I) * REAL(ZPP)
C
  140 CONTINUE
C
  150 UYY = - UXX
C
C
C      TERMS DUE TO OTHER SOURCES
C
  200 KFLAG = IFLAG
      RETURN
C***********************************************************************
      END