1 : C CARWRK SOURCE PV 11/03/07 21:15:15 6885 2 : SUBROUTINE CARWRK(XF,YF,KFLAG,U,UX,UY,UXX,UXY,UYY,ITYPE, 3 : $ XORG,YORG,NPN,EDICON,ICOO,ITRAVA) 4 : C 5 : C EVALUATE POTENTIAL AND ITS FIRST AND SECOND DERIVATIVES 6 : C----------------------------------------------------------------------- 7 : C 8 : IMPLICIT INTEGER(I-N) 9 : IMPLICIT REAL*8(A-H,O-Z) 10 : C----------------------------------------------------------------------- 11 : C 12 : COMPLEX*16 Z, ZE, ZN, ZP, ZPP 13 : C 14 : C----------------------------------------------------------------------- 15 : SEGMENT EDICON 16 : INTEGER KSTRT, KSTEP, NMIR, IS 17 : REAL*8 CROT, SROT, SYMFCT 18 : LOGICAL LREAL, LIMAG 19 : ENDSEGMENT 20 : SEGMENT ICOO 21 : REAL*8 X(MV),Y(MV),P(MV),WNODE(MV) 22 : INTEGER LISVO(MV) 23 : ENDSEGMENT 24 : SEGMENT ITRAVA 25 : REAL*8 KENN(M42,2),SIGMA(M42),DELRHO(M42),C(M50,M50) 26 : REAL*8 AK(M50),UM(M50),RM(M50) 27 : INTEGER IL(M50) 28 : ENDSEGMENT 29 : C 30 : C COMMON // NEU, NPN, NPU 31 : C 32 : C TABLES AND WORKING STORE FOR POTENTIAL AND FIELD EDIT 33 : C 34 : C 35 : C SET UP EQUATIONS FOR LEAST-SQUARES FIT 36 : C 37 : M50=RM(/1) 38 : NPU = NPN 39 : NTERM =MIN(( NPU - 4),15) 40 : IF(LREAL .AND. LIMAG) NTERM = MIN( ((NPU - 5) / 2 ),7) 41 : KEND = KSTRT + NTERM * KSTEP 42 : XE = XF - XORG 43 : YE = YF - YORG 44 : ZE = CMPLX(XE,YE) 45 : IFLAG = 0 46 : J = 0 47 : C 48 : DO 90 JJ = 1, NPN 49 : C 50 : IF(WNODE(JJ) .EQ. 0.0) GO TO 90 51 : XN = X(JJ) - XORG 52 : YN = Y(JJ) - YORG 53 : ZN = CMPLX(XN,YN) 54 : U = P(JJ) 55 : C 56 : C LOGARITHMIC TERMS DUE TO FILAMENTS 57 : C TERMS DUE TO OTHER SOURCES 58 : C 59 : C 60 : C COEFFICIENT MATRIX 61 : C 62 : J = J + 1 63 : C(J,1) = WNODE(JJ) 64 : I = 1 65 : C 66 : IF(KEND .LT. KSTRT) GO TO 80 67 : C 68 : DO 70 KI = KSTRT, KEND, KSTEP 69 : Z = ZN**KI 70 : C 71 : IF(.NOT. LREAL) GO TO 60 72 : I = I + 1 73 : C(J,I) = WNODE(JJ) * REAL(Z) 74 : C 75 : 60 IF(.NOT. LIMAG) GO TO 70 76 : I = I + 1 77 : C(J,I) = WNODE(JJ) * IMAG(Z) 78 : C 79 : 70 CONTINUE 80 : C 81 : 80 C(J,I+1) = WNODE(JJ) * U 82 : C 83 : 90 CONTINUE 84 : C 85 : C----------------------------------------------------------------------- 86 : C PERFORM THE LEAST SQUARES FIT 87 : C 87 : CALL GLSQ(C,AK,IL,J,I,ALFA,1.0D-10,1.0D-10,M50) 88 : C 89 : C----------------------------------------------------------------------- 90 : C PERFORM THE ANALYTIC DIFFERENTIATION 91 : C 92 : IFLAG = KFLAG 93 : U = AK(1) 94 : UX = 0.0 95 : UY = 0.0 96 : UXX = 0.0 97 : UXY = 0.0 98 : UYY = 0.0 99 : C 100 : C HARMONIC TERMS 101 : C 102 : IF(KEND .LT. KSTRT) GO TO 150 103 : I = 1 104 : C 105 : DO 140 KI = KSTRT, KEND, KSTEP 106 : Z = ZE**KI 107 : ZP = (1.0,0.0) 108 : IF(KI .GT. 1) ZP = FLOAT(KI) * ZE ** (KI-1) 109 : ZPP = (0.0,0.0) 110 : IF(KI .EQ. 2) ZPP = (2.0,0.0) 111 : IF(KI .GT. 2) ZPP = FLOAT(KI*(KI-1)) * ZE ** (KI-2) 112 : C 113 : IF(.NOT. LREAL) GO TO 120 114 : I = I + 1 115 : U = U + AK(I) * REAL(Z) 116 : UX = UX + AK(I) * REAL(ZP) 117 : UY = UY - AK(I) * IMAG(ZP) 118 : UXX = UXX + AK(I) * REAL(ZPP) 119 : UXY = UXY - AK(I) * IMAG(ZPP) 120 : C 121 : 120 IF(.NOT. LIMAG) GO TO 140 122 : I = I + 1 123 : U = U + AK(I) * IMAG(Z) 124 : UX = UX + AK(I) * IMAG(ZP) 125 : UY = UY + AK(I) * REAL(ZP) 126 : UXX = UXX + AK(I) * IMAG(ZPP) 127 : UXY = UXY + AK(I) * REAL(ZPP) 128 : C 129 : 140 CONTINUE 130 : C 131 : 150 UYY = - UXX 132 : C 133 : C 134 : C TERMS DUE TO OTHER SOURCES 135 : C 136 : 200 KFLAG = IFLAG 137 : RETURN 138 : C*********************************************************************** 139 : END 140 : 141 : 142 :
© Cast3M 2003 - All rights reserved.
Disclaimer