I have copied a FORTRAN 66/IV program from a 1988 report, and I am trying to compile it with gfortran (mingw for windows). I have reduced a long list of errors down to 3 (plus 2 warnings), and I cannot get any further despite my best efforts. I would be grateful for any help and advice.
The errors:
green.f:298.16: RDE=(EXPR(J)/REYLOC(J)-EXPR(J-1)/REYLOC(J-1))/ZDIFF 1 Error: PROCEDURE attribute conflicts with COMMON attribute in 'expr' at (1) green.f:412.7: 1,5X,F10.4) 1 Error: Nonnegative width required in format string at (1) green.f:390.19: WRITE(OUT,11)(Z(J),CUR(J),CP(J),PH(J),RMSQ(J),U(J),Q(J),J=1,P) 1 Error: FORMAT label 11 at (1) not defined green.f:249.61: CALL OPTION(MSQLOC,RAD,LOCR,RDASH,X,Y,FEQ,HALFCO,H,H1,IMBAL, 1 Warning: Type mismatch in argument 'imbalcrxn' at (1); passed REAL(4) to INTEGER(4) green.f:122.72: 1OUT) 1 Warning: Missing actual argument for argument 'out' at (1) This is the source code:
00000001 PROGRAM GREEN 00000002 REAL TRAD(100),CUR(100),V(100),EXPR(100),Z(100),MSQ(100),U(100) 00000003 1,XX(100),REYLOC(100),Y(3),YDASH(3),AA(21),NU,INCR,LSCALE,LAT, 00000004 2IMB,LIMIT,MINF,MINFSQ,Z1(100),LOCR,MSQLOC,Q(100),PH(100), 00000005 3CP(100),RI(100),RSTAR(100),RK(100),FACT,TW(100),DR(100),DRAG, 00000006 4XEXP(20),DEXP(20),THEXP(20),HEXP(20),FCEXP(20),DSTP(100), 00000007 5THPLT(100),HPLT(100),CFPLT(100) 00000008 INTEGER AXIS,COND,CURV,FSTART,STRET,DIL,P,HQ,TQ,DEV,SURF,OUT 00000009 1,P1,NSTA 00000010 CHARACTER*8 LABEL 00000011 COMMON/CB1/AXIS,CURV,COND,FSTART,STRET,DIL,P,DEV,JTE,RATIO 00000012 COMMON/CB2/REC,TRAD,UDASH,Z,J,U,RK 00000013 COMMON/CB3/CUR,EXPR,REYLOC 00000014 COMMON/CB4/MINFSQ,TINF,RC 00000015 COMMON/CB5/RTHETA,THETA,FEQ,HALFC,CRXN,H,H1,RICH,FACT 00000016 COMMON/CB6/MINF,TSTAG,M,TQ,HQ,IRPT,IDENT,KK,HDASH 00000017 COMMON/CB7/NSTA,XEXP,DEXP,THEXP,HEXP,CFEXP,P1,LABEL 00000018 CALL DIG 00000019 CALL SIMBEG 00000020 DEV=15 00000021 OUT=16 00000022 IRPT=1 00000023 35 CALL INPUT(TRAD,CUR,EXPR,V,Z,Z1,XX,Y1,Y2,RC,REC,SURF,LS,RK) 00000024 00000025 00000026 P1=P 00000027 Y(1)=Y1 00000028 Y(2)=Y2 00000029 WRITE(OUT,20)IDENT,MINF,RC 00000030 CALL EVALFP(PINF,HQINF,REC,Q,V,MSQ,U,P,REYLOC,PH,CP,RK) 00000031 00000032 00000033 IF(IRPT.GT.2)GOTO 21 00000034 CALL INDATA(Z,Z1,XX,CUR,PH,CP,MSQ,U,Q,SURF,AXIS,P,CURV,OUT,RK) 00000035 00000036 00000037 21 IF(HQ.NE.0)GOTO 14 00000038 H=(Y(1)+1)*(1+0.2*REC*MSQ(LS))-1 00000039 GOTO 15 00000040 14 H=Y(1) 00000041 Y(1)=(Y(1)+1)/(1+0.2*REC*MSQ(LS))-1 00000042 15 IF(TQ.NE.1) GOTO 16 00000043 Y(2)=Y(2)*TRAD(LS) 00000044 RTHETA=Y(2)*REYLOC(LS)/TRAD(LS) 00000045 GOTO 17 00000046 16 RTHETA=Y(2) 00000047 Y(2)=TRAD(LS)*Y(2)/REYLOC(LS) 00000048 17 J=LS+1 00000049 X=Z(LS) 00000050 CALL FUNC(5,YDASH,X,Y,LS) 00000051 00000052 00000053 RI(LS)=RICH 00000054 RSTAR(LS)=Y(2) 00000055 SLOPE=UDASH/U(LS) 00000056 DH1DHB=-1.72/(Y(1)-1)**2.0-0.2*(Y(1)-1) 00000057 IF(FSTART.EQ.0)GOTO 9 00000058 Y(3)=Y(2)*DH1DHB*HDASH/TRAD(LS)+H1*(HALFC-(H+1)*Y(2)*SLOP 00000059 1E/TRAD(LS)) 00000060 Y(3)=Y(3)-CRXN/TRAD(LS) 00000061 9 CF=2*HALFC 00000062 TW(LS)=HALFC*1.2*REYLOC(LS)*REYLOC(LS)*1.51E-5*1.51E-5 00000063 RW3=THETA*H 00000064 G=SQRT(1/HALFC)*(1-1/Y(1)) 00000065 PI=-H*THETA*UDASH/(U(LS)*HALFC) 00000066 CFS=CF*Q(LS) 00000067 IF(COND.EQ.0.AND.CURV.EQ.0.AND.STRET.EQ.0.AND.DIL.EQ.0)GOTO4 00000068 WRITE(OUT,5) 00000069 IF(COND.GT.0)WRITE(OUT,6) 00000070 IF(CURV.GT.0)WRITE(OUT,8) 00000071 IF(DIL.GT.0)WRITE(OUT,18) 00000072 IF(STRET.GT.0.AND.(AXIS.EQ.1.OR.COND.GT.0))WRITE(OUT,12) 00000073 4 IF(FSTART.GT.0)WRITE(OUT,11)HDASH 00000074 WRITE(OUT,2) 00000075 WRITE(OUT,3)Z1(LS),Y(1),CF,RTHETA,THETA,RW3,H,CFR,G,PI,Y(3),FEQ 00000076 00000077 00000078 DSTP(LS)=RW3 00000079 THPLT(LS)=THETA 00000080 HPLT(LS)=H 00000081 CFPLT(LS)=CF 00000082 00000083 DEL=Z(LS+1)-Z(LS) 00000084 CALL STEP(YDASH,Y,DEL) 00000085 00000086 00000087 LL=LS+1 00000088 DRAG=0.0 00000089 DO 10 J=LL,P 00000090 00000091 00000092 SW1=Z(J)-Z(J-1) 00000093 LIMIT=Z(J) 00000094 IF(SQ1.LT.DEL)DEL=SW1 00000095 13 CALL VINT(3,DEL,X,Y,LIMIT,0.5E-5,1000,YDASH,AA,21,OUT,J) 00000096 00000097 00000098 IF(Y(3).LE.-0.009)Y(3)=-0.009 00000099 X=Z(J) 00000100 CALL FUNC(0,YDASH,X,Y,J) 00000101 CF=2.0*HALFC 00000102 TW(J)=HALFC*1.2*REYLOC(J)*REYLOC(J)*1.51E-5*1.51E-5 00000103 RSTAR(J)=Y(2) 00000104 RI(J)=RICH 00000105 RW3=THETA*H 00000106 CFR=CF*Q(J) 00000107 G=SQRT(1.0/HALFC)*(1.0-1/Y(1)) 00000108 PI=-H*THETA*UDASH/(U(J)*HALFC) 00000109 DR(J)=((TW(J)+TW(J-1))/2)*(Z1(J)-Z1(J-1))*0.0325 00000110 DRAG=DRAG+DR(J) 00000111 WRITE(OUT,3)Z1(J),Y(1),CF,RTHETA,THETA,RW3,H,CFR,G,PI,Y(3),FEQ 00000112 00000113 00000114 DSTP(J)=RW3 00000115 THPLT(J)=THETA 00000116 HPLT(J)=H 00000117 CFPLT(J)=CF 00000118 00000119 10 CONTINUE 00000120 WRITE(OUT,40)DRAG 00000121 IF(AXIS.EQ.1.OR.CURV.GT.0) CALL CUORAX(AXIS,CURV,P,RSTAR,RI,Z1,LS, 00000122 1OUT) 00000123 00000124 00000125 CALL ENDPLT 00000126 00000127 READ(DEV,*)IRPT 00000128 IF(IRPT.GT.5)STOP 00000129 GOTO 35 00000130 2 FORMAT(3(/),1H ,5X,1HX,7X,2HHT,5X,5HCFLOC,4X,6HRTHETA,6X,5HTHETA 00000131 1,5X,7HDELSTAR,5X,1HH,6X,5HCFREF,8X,1HG,9X,2HPI,8X,1HF,7X,3HFEQ) 00000132 3 FORMAT(1H ,1PE10.3,1X,0PF6.3,1X,F8.5,1X,3(1PE10.3,1X),0PF7.3,1X, 00000133 1F8.5,1X,2(1PE10.3,1X),0PF8.5,1X,F8.5) 00000134 40 FORMAT(1H,20X,26HTOTAL SKIN FRICTION DRAG= ,F12.6) 00000135 5 FORMAT(1H0,20HALLOWANCES MADE FOR:) 00000136 6 FORMAT(1H ,20X,26HCONVERGENCE AND DIVERGENCE) 00000137 18 FORMAT(1H ,20X,10HDILATATION) 00000138 8 FORMAT(1H ,20X,22HLONGITUDINAL CURVATURE) 00000139 11 FORMAT(1H0,23HINITIAL VALUE OF DH/DX=,F9.4) 00000140 12 FORMAT(1H ,20X,18HLATERAL STRETCHING) 00000141 20 FORMAT(1H1,3X,73HLAG ENTRAINMENT B.L CALCULATION FOR TWO DIMENSION 00000142 1AL AND AXISYMMETRIC FLOW,5(/),1H,10HIDENT. NO.,1X,I5,5X,6HMINF= , 00000143 2E11.4,5X,4HRC= ,E11.4) 00000144 STOP 00000145 END 00000146 00000147 SUBROUTINE INPUT(TRAD,CUR,EXPR,V,Z,Z1,XX,Y1,Y2,RC,REC,SURF,LS,RK) 00000148 REAL TRAD(100),CUR(100),V(100),EXPR(100),Z(100),Z1(100),XX(100), 00000149 1MINF,RK(100),XEXP(20),DEXP(20),THEXP(20),HEXP(20),CFEXP(20) 00000150 INTEGER AXIS,COND,CURV,STRET,FSTART,DIL,P,SURF,DEV,TQ,HQ,NSTA,P1 00000151 CHARACTER*8 LABEL 00000152 COMMON/CB1/AXIS,CURV,COND,FSTART,STRET,DIL,P,DEV,JTE,RATIO 00000153 COMMON/CB6/MINF,TSTAG,M,TQ,HQ,IRPT,IDENT,KK,HDASH 00000154 COMMON/CB7/NSTA,XEXP,DEXP,THEXP,HEXP,CFEXP,P1,LABEL 00000155 GOTO(20,5,6,8,9),IRPT 00000156 20 READ(DEV,*)P,AXIS,SURF 00000157 READ(DEV,*)JTE 00000158 DO 10 J=1,P 00000159 TRAD(J)=1.0 00000160 10 CONTINUE 00000161 IF(SURF.EQ.0.AND.AXIS.EQ.0)GOTO 3 00000162 READ(DEV,*)(Z1(J),XX(J),J=1,P) 00000163 IF(AXIS.EQ.0) GOTO 4 00000164 DO 11 J=1,P 00000165 TRAD(J)=XX(J) 00000166 11 CONTINUE 00000167 4 IF(SURF.EQ.0)GOTO 1 00000168 Z(1)=0.0 00000169 DO 13 J=2,P 00000170 DELZ=Z1(J)-Z1(J-1) 00000171 DELX=XX(J)-XX(J-1) 00000172 Z(J)=SQRT(DELZ*DELZ*DELX*DELX)+Z(J-1) 00000173 13 CONTINUE 00000174 GOTO 5 00000175 3 READ(DEV,*)(Z1(J),J=1,P) 00000176 1 DO 14 J=1,P 00000177 Z(J)=Z1(J) 00000178 14 CONTINUE 00000179 5 READ(DEV,*)M 00000180 READ(DEV,*)(V(J),J=1,P) 00000181 READ(DEV,*)(RK(J),J=1,P) 00000182 6 READ(DEV,*)COND,CURV,STRET,FSTART,DIL 00000183 IF(COND.NE.1) GOTO7 00000184 READ(DEV,*)(EXPR(J),J=1,P) 00000185 7 IF(CURV.NE.1) GOTO 8 00000186 READ(DEV,*)(CUR(J),J=1,P) 00000187 8 READ(DEV,*)Y2,Y1,RC,MINF,TSTAG,REC 00000188 READ(DEV,*),TQ,HQ,LS 00000189 9 IF(FSTART.EQ.1) READ(DEV,*)HDASH 00000190 READ(DEV,*)IDENT 00000191 00000192 00000193 READ(DEV,*)LABEL 00000194 READ(DEV,*)NSTA 00000195 READ(DEV,*)(XEXP(I),I=1,NSTA) 00000196 READ(DEV,*)(DEXP(I),I=1,NSTA) 00000197 READ(DEV,*)(THEXP(I),I=1,NSTA) 00000198 READ(DEV,*)(HEXP(I),I=1,NSTA) 00000199 READ(DEV,*)(CFEXP(I),I=1,NSTA) 00000200 RETURN 00000201 END 00000202 00000203 SUBROUTINE FUNC(N,YDASH,X,Y,LL) 00000204 REAL TRAD(100),CUR(100),Z(100),U(100),MSQLOC,LOCR,Y(3),YDASH(3), 00000205 1LAMDA,IMBAL,LAT,LSQ,NEWK,RK(100),FACT 00000206 INTEGER AX,CU,CO,FS,ST,DIL,P,DEV,LL 00000207 COMMON/CB1/AX,CU,CO,FS,ST,DIL,P,DEV,JTE,RATIO 00000208 COMMON/CB2/REC,TRAD,UDASH,Z,J,U,RK 00000209 COMMON/CB5/RTHETA,THETA,FEQ,HALFC,CRXN,H,H1,RICH,FACT 00000210 CRXN=0.0 00000211 LAMDA=1 00000212 IMBAL=0.0 00000213 RAD=1.0 00000214 SW1=Z(J)-Z(J-1) 00000215 IF(AX.GT.0)RAD=(TRAD(J-1)+(X-Z(J-1))*(TRAD(J)-TRAD(J-1))/SW1) 00000216 UDASH=(U(J)-U(J-1))/SW1 00000217 ULOC=U(J-1)+UDASH*(X-Z(J-1)) 00000218 CALL VELMR(MSQLOC,TLOC,LOCR,ULOC) 00000219 00000220 00000221 IF(Y(1).GT.2.65)Y(1)=2.65 00000222 H=(1+Y(1))*(1+0.2*REC*MSQLOC)-1.0 00000223 THETA=Y(2)/RAD 00000224 DELTA=0.0 00000225 10 RTHETA=LOCR*THETA 00000226 00000227 00000228 FACT=1+(2000-RTHETA)*0.00003734 00000229 IF(FACT.LT.1.0)THEN 00000230 FACT=1.0 00000231 END IF 00000232 HALFCO=FACT*(0.005065/(0.4342945*ALOG(RTHETA*(1+.056*MSQLOC) 00000233 1)-1.02)-0.000375)/SQRT(1+0.2*MSQLOC) 00000234 00000235 IF(J.GT.JTE)HALFCO=0 00000236 HBO=1/(1-6.55*SQRT(HALFCO*(1+0.04*MSQLOC))) 00000237 HALFC=HALFCO*(0.9/(Y(1)/HBO-0.4)-0.5) 00000238 CALL SFCR(HALFC,RK(LL),AAA) 00000239 00000240 00000241 HALFC=HALFC*AAA 00000242 IF(HALFC.LT.0.000001) HALFC=0.000001 00000243 H1=3.15+1.72/(Y(1)-1)-0.01*(Y(1)-1)**2.0 00000244 DH1DHB=-1.72/(Y(1)-1)**2.0-0.02*(Y(1)-1) 00000245 Q=((Y(1)-1.0)/(6.432*Y(1)))**2.0/(1.0+0.04*MSQLOC) 00000246 FEQ=H1*(HALFC-(H+1.0)*(HALFC-Q)/(0.8*H)) 00000247 RDASH=RAD*HALFC-(H+2-MSQLOC)*Y(2)*UDASH/ULOC 00000248 DU=THETA*(H1+H)*UDASH/ULOC 00000249 CALL OPTION(MSQLOC,RAD,LOCR,RDASH,X,Y,FEQ,HALFCO,H,H1,IMBAL, 00000250 1CRXN,RICH,LAMDA,DU) 00000251 00000252 00000253 IF(FEQ.LT.0.0)FEQ=0.0 00000254 IF(N.EQ.0)RETURN 00000255 IF(N.EQ.5) Y(3)=FEQ 00000256 DUEQ=(H1+H)*(HALFC-FEQ/H1)/(H+1) 00000257 FMA=SQRT(1+0.1*MSQLOC) 00000258 FMB=1+0.075*MSQLOC*(1+0.2*MSQLOC)/(1+0.1*MSQLOC) 00000259 YDASH(1)=((RAD+DELTA)*Y(3)-RAD*(H1*HALFC)+H1*(H+1)*Y(2)*UDASH/ 00000260 1ULOC+CRXN)/Y(2)/DH1DHB 00000261 YDASH(2)=RDASH+IMBAL 00000262 IF(Y(3).LE.-0.009)Y(3)=0.009 00000263 YDASH(3)=(Y(3)*(Y(3)+0.02)+0.5333*HALFCO)/(Y(3)+0.01)*(2.8* 00000264 1LAMDA*FMA*(SQRT(0.64*HALFCO+0.024*FEQ+1.2*FEQ*FEQ)-SQRT(0.64* 00000265 2HALFCO+0.024*Y(3)+1.2*Y(3)*Y(3)))+DUEQ-DU*FMB)/(Y(2)*(H1+H)) 00000266 3 RETURN 00000267 END 00000268 00000269 SUBROUTINE VELMR(MSQLOC,TLOC,LOCR,ULOC) 00000270 REAL MSQLOC,MINFSQ,LOCR 00000271 COMMON/CB4/MINFSQ,TINF,RC 00000272 MSQLOC=ULOC*ULOC*MINFSQ/(1+0.2*MINFSQ*(1-ULOC*ULOC)) 00000273 TLOC=TINF*(1+0.2*MINFSQ)/(1+0.2*MSQLOC) 00000274 LOCR=RC*ULOC*TLOC/TINF*(TLOC+114)/(TINF+114) 00000275 RETURN 00000276 END 00000277 00000278 SUBROUTINE OPTION(MSQLOC,RAD,LOCR,RDASH,X,Y,FEQ,HALFCO,H,H1,IMBAL, 00000279 1 CRXN,RICH,LAMBDA,DU) 00000280 REAL REYLOC(100),EXP(100),Z(100),Y(3),LOCR,CUR(100),IMBAL,LDASH, 00000281 1LAT,LSCALE,U(100),TRAD(100),MSQLOC,LAMDA,LSQ,L1,L2,L3 00000282 INTEGER AX,CU,CO,FS,ST,DIL,P,DEV,JTE,RATIO 00000283 COMMON/CB1/AX,CU,CO,FS,ST,DIL,P,DEV,JTE,RATIO 00000284 COMMON/CB3/CUR,EXPR,REYLOC 00000285 COMMON/CB2/REC,TRAD,UDASH,Z,J,U,RK 00000286 IMBAL=0.0000 00000287 RSLOPE=0.0000 00000288 RICH=0.0 00000289 DPHIDZ=0.0 00000290 CRXN=0.0000 00000291 L1=1.0 00000292 L2=1.0 00000293 L3=1.0 00000294 SW4=1+0.2*REC*MSQLOC 00000295 ZDIFF=Z(J)-Z(J-1) 00000296 IF(CO.EQ.0)GOTO 2 00000297 LDASH=(REYLOC(J)-REYLOC(J-1))/ZDIFF 00000298 RDE=(EXPR(J)/REYLOC(J)-EXPR(J-1)/REYLOC(J-1))/ZDIFF 00000299 IMBAL=RDE-RDASH 00000300 CRXN=-IMBAL*2*(H1*(Y(1)-1)-Y(1))/(2*Y(1)-1) 00000301 IF(J.GT.JTE)GOTO 4 00000302 2 IF(CU.EQ.0.AND.ST.EQ.0.AND.DIL.EQ.0)GOTO 5 00000303 IF(CU.EQ.0)GOTO 3 00000304 CURV=CUR(J-1)+(CUR(J)-CUR(J-1))*(X-Z(J-1))/ZDIFF 00000305 RICH=0.6667*(1+0.2*MSQLOC)*(0.3+H1/Y(1))*Y(2)*(H1+H)*CURV/RAD 00000306 3 IF(ST.EQ.0)GOTO 4 00000307 DPHIDZ=-IMBAL/Y(2)/(2*Y(1)-1) 00000308 RSLOPE=(TRAD(J)-TRAD(J-1))/(ZDIFF*RAD) 00000309 4 ALPHA=4.5 00000310 IF(CURV.GT.0.0)ALPHA=7.0 00000311 L1=1+ALPHA*RICH 00000312 IF(ST.EQ.1)L2=1-2.33*(H1/Y(1)+0.3)*(H+H1)*Y(2)/RAD*(DPHIDZ+RSLOPE 00000313 1) 00000314 IF(DIL.EQ.1)L3=1+2.33*MSQLOC*(1+H1/Y(1))*DU/RAD 00000315 LAMDA=L1*L2*L3 00000316 IF(J.GT.JTE)LAMDA=0.5 00000317 IF(LAMDA.LT.0.499)LAMDA=0.499 00000318 LSQ=LAMDA*LAMDA 00000319 C=0.5333*HALFCO*(1-1/LSQ)-0.02*FEQ/LSQ-FEQ*FEQ/LSQ 00000320 IF(C.GT.0.0000999)C=0.0000999 00000321 FEQ=SQRT(0.0001-C)-0.01 00000322 5 RETURN 00000323 END 00000324 00000325 SUBROUTINE EVALFP(PINF,HQINF,REC,Q,V,MSQ,U,P,REYLOC,PH,CP,RK) 00000326 REAL MINF,V(100),MINFSQ,U(100),MSQ(100),REYLOC(100),Q(100), 00000327 1PH(100),CP(100),RK(100) 00000328 INTEGER P,TQ,HQ 00000329 COMMON/CB4/MINFSQ,TINF,RC 00000330 COMMON/CB6/MINF,TSTAG,M,TQ,HQ,IRPT,IDENT,KK,HDASH 00000331 MINFSQ=MINF*MINF 00000332 TR=1.0+0.2*MINFSQ 00000333 PINF=TR**(-3.5) 00000334 HQINF=TR**3.5/(0.7*MINFSQ) 00000335 TINF=TSTAG/TR 00000336 GOTO(1,2,3,4,5),M 00000337 1 DO 10 JJ=1,P 00000338 MSQ(JJ)=MINFSQ*V(JJ)*V(JJ)/(1.0+0.2*MINFSQ*(1.0-V(JJ)*V(JJ))) 00000339 CALL COMPTU(M,MSQ,U,V,TR,REYLOC,JJ,Q,CP,PH,PINF,HQINF) 00000340 10 CONTINUE 00000341 GOTO 5 00000342 2 DO 11 JJ=1,P 00000343 QQ=1.0+0.7*MINFSQ*V(JJ) 00000344 MSQ(JJ)=5.0*(TR*QQ**(-0.2857143)-1.0) 00000345 CALL COMPTU(M,MSQ,U,V,TR,REYLOC,JJ,Q,CP,PH,PINF,HQINF) 00000346 11 CONTINUE 00000347 GOTO 5 00000348 3 DO 12 JJ=1,P 00000349 MSQ(JJ)=5.0*(V(JJ)**(-0.2857143)-1.0) 00000350 CALL COMPTU(M,MSQ,U,V,TR,REYLOC,JJ,Q,CP,PH,PINF,HQINF) 00000351 12 CONTINUE 00000352 GOTO 5 00000353 4 DO 13 JJ=1,P 00000354 MSQ(JJ)=V(JJ)*V(JJ) 00000355 CALL COMPTU(M,MSQ,U,V,TR,REYLOC,JJ,Q,CP,PH,PINF,HQINF) 00000356 13 CONTINUE 00000357 5 RETURN 00000358 END 00000359 00000360 SUBROUTINE COMPTU(M,MSQ,U,V,TR,REYLOC,JJ,Q,CP,PH,PINF,HQINF) 00000361 COMMON/CB4/MINFSQ,TINF,RC 00000362 REAL MSQ(100),MINFSQ,U(100),V(100),REYLOC(100),Q(100),PH(100), 00000363 1CP(100) 00000364 IF(M-1)3,1,3 00000365 3 U(JJ)=SQRT(TR*MSQ(JJ)/(MINFSQ*(1.0+0.2*MSQ(11)))) 00000366 GOTO 2 00000367 1 U(JJ)=V(JJ) 00000368 2 TLOC=TINF*(1+0.2*MINFSQ)/(1.0+0.2*MSQ(JJ)) 00000369 REYLOC(JJ)=RC*U(JJ)*TLOC/TINF*(TLOC+114)/(TINF+114) 00000370 Q(JJ)=U(JJ)*U(JJ)*(TLOC/TINF)**2.5 00000371 PH(JJ)=(1+0.2*MSQ(JJ))**(-3.5) 00000372 CP(JJ)=(PH(JJ)-PINF)*HQINF 00000373 RETURN 00000374 END 00000375 00000376 SUBROUTINE INDATA(Z,Z1,XX,CUR,PH,CP,MSQ,U,Q,SURF,AXIS,P,CURV, 00000377 1OUT,RK) 00000378 REAL Z(100),Z1(100),XX(100),CUR(100),PH(100),CP(100),RMSQ(100), 00000379 1MSQ(100),U(100),Q(100),RK(100) 00000380 INTEGER SURF,CURV,AXIS,P,OUT 00000381 DO 15 J=1,P 00000382 RMSQ(J)=SQRT(MSQ(J)) 00000383 15 CONTINUE 00000384 IF(SURF.GT.0.OR.AXIS.GT.0)GOTO 1 00000385 IF(CURV.GT.0)GOTO 2 00000386 WRITE(OUT,20) 00000387 WRITE(OUT,10)(Z(J),CP(J),PH(J),RMSQ(J),U(J),Q(J),RK(J),J=1,P) 00000388 RETURN 00000389 2 WRITE(OUT,21) 00000390 WRITE(OUT,11)(Z(J),CUR(J),CP(J),PH(J),RMSQ(J),U(J),Q(J),J=1,P) 00000391 RETURN 00000392 1 IF(AXIS.GT.0)GOTO 3 00000393 IF(CURV.GT.0)GOTO 4 00000394 WRITE(OUT,22) 00000395 GOTO 6 00000396 4 WRITE(OUT,23) 00000397 GOTO 7 00000398 3 IF(CURV.GT.9)GOTO 5 00000399 WRITE(OUT,24) 00000400 6 WRITE(OUT,12)(Z1(J),XX(J),Z(J),CP(J),PH(J),RMSQ(J),U(J),Q(J), 00000401 1RK(J),J=1,P) 00000402 RETURN 00000403 5 WRITE(OUT,25) 00000404 7 DO 16 J=1,P 00000405 WRITE(OUT,13)Z1(J),XX(J),Z(J),CUR(J),CP(J),PH(J),RMSQ(J), 00000406 1U(J),Q(J) 00000407 16 CONTINUE 00000408 RETURN 00000409 10 FORMAT(1H ,20X,1PE11.4,3X,0PF10.4,3X,F7.5,5X,F6.4,5X,F7.4,3X, 00000410 1F10.4,3X,1PE11.4) 00000411 11 FORMAT(1H ,20X,1PE11.4,5X,E11.4,3X,0PF10.4,3X,F7.5,5X,F6.4,5X,F7.4 00000412 1,5X,F10.4) 00000413 12 FORMAT(1H ,10X,1PE11.4,5X,E11.4,5X,E11.4,3X,0PF10.4,3X,F7.5,5X, 00000414 1F6.4,5X,F7.4,3X,F10.4,3X,F7.4) 00000415 13 FORMAT(1H ,4X,1PE11.4,4X,E11.4,4X,E11.4,4X,E11.4,2X,0PF10.4,2X, 00000416 1F7.5,4X,F6.4,4X,F7.4,2X,F10.4) 00000417 20 FORMAT(1H0,25X,1HX,13X,2HCP,9X,3HP/H,9X,1HM,8X,6HU/UREF,6X, 00000418 16HQ/QREF,6X,2HRK) 00000419 21 FORMAT(1H0,25X,1HX,14X,5HLCURV,10X,2HCP,9X,3HP/H,9X,1HM,8X, 00000420 16HU/UREF,6X,6HQ/QREF) 00000421 22 FORMAT(1H0,15X,1HX,15X,1HZ,15X,1HS,13X,2HCP,9X,3HP/H,9X,1HM,8X, 00000422 16HU/UREF,6X,6HQ/QREF) 00000423 23 FORMAT(1H0,10X,1HX,14X,1HZ,14X,1HS,13X,5HLCURV,9X,2HCP,8X,3HP/H, 00000424 18X,1HM,7X,6HU/UREF,5X,6HQ/QREF) 00000425 24 FORMAT(1H0,15X,1HX,14X,3HRAD,14X,1HS,13X,2HCP,9X,3HP/H,9X,1HM, 00000426 18X,6HU/UREF,6X,6HQ/QREF,6X,2HRK) 00000427 25 FORMAT(1H0,10X,1HX,13X,3HRAD,13X,1HS,13X,5HLCURV,9X,2HCP,8X, 00000428 13HP/H,8X,1HM,7X,6HU/UREF,5X,6HQ/QREF) 00000429 END 00000430 00000431 SUBROUTINE STEP(YDASH,Y,DEL) 00000432 REAL YDASH(3),Y(3),S(3) 00000433 DO 1 I=1,3 00000434 IF(ABS(YDASH(I)).LT.0.1E-05)YDASH(I)=0.1E-05 00000435 1 CONTINUE 00000436 VAR=0.01 00000437 S(1)=ABS(VAR/YDASH(1)) 00000438 S(2)=ABS(VAR*Y(2)/YDASH(2)) 00000439 S(3)=DEL*VAR 00000440 DEL=S(1) 00000441 DO 2 I=1,3 00000442 IF(S(I).LT.DEL)DEL=S(I) 00000443 2 CONTINUE 00000444 RETURN 00000445 END 00000446 00000447 SUBROUTINE CUORAX(AXIS,CURV,P,RSTAR,RI,Z1,LS,OUT) 00000448 REAL RSTAR(100),RI(100),Z1(100) 00000449 INTEGER AXIS,CURV,P,OUT 00000450 IF(AXIS.EQ.1.AND.CURV.GT.0)GOTO 1 00000451 IF(AXIS.EQ.0)GOTO 2 00000452 WRITE(OUT,20) 00000453 WRITE(OUT,10)(Z1(J),RSTAR(J),J=LS,P) 00000454 RETURN 00000455 2 WRITE(OUT,21) 00000456 WRITE(OUT,10)(Z1(J),RI(J),J=LS,P) 00000457 RETURN 00000458 1 WRITE(OUT,22) 00000459 WRITE(OUT,12)(Z1(J),RSTAR(J),RI(J),J=LS,P) 00000460 RETURN 00000461 10 FORMAT(1H ,1PE11.4,5X,E13.6) 00000462 12 FORMAT(1H ,1PE11.4,5X,E13.6,5X,E13.6) 00000463 20 FORMAT(///1H0,6X,1HX,13X,7HR*THETA) 00000464 21 FORMAT(///1H0,6X,1HX,9X,14HRICHARDSON NO.) 00000465 22 FORMAT(///1HO,6X,1HX,13X,7HR*THETA,7X,14HRICHARDSON NO.) 00000466 END 00000467 00000468 SUBROUTINE VINT(N,H,X,Y,XD,E,NS,DY,RK,N7,OUT,J) 00000469 INTEGER QD,QDP,TD,TDP,V0,VE,U0,RD,OUT,LL 00000470 DIMENSION Y(N),DY(N),RK(N7) 00000471 1 IF(ABS(XD-X).LT.1E-20)GOTO 100 00000472 IF(ABS(H).LT.1E-20)GOTO 100 00000473 G0=5*E 00000474 ED=0.03125*G0 00000475 2 V0=0 00000476 GOTO 23 00000477 3 QD=0 00000478 H0=0 00000479 9 X0=X 00000480 F0=XD-X 00000481 Y0=F0-H 00000482 IF(H.GT.0.)GOTO 10 00000483 Y0=-Y0 00000484 10 IF(Y0.GT.0.)GOTO 11 00000485 HD=F0 00000486 U0=-1 00000487 GOTO 12 00000488 11 HD=H 00000489 U0=0 00000490 12 V0=V0+1 00000491 IF(V0.GT.NS) GOTO 100 00000492 QDP=QD+1 00000493 13 DO 22 TDP=QDP,7,1 00000494 TD=TDP-1 00000495 X=X0+H0 00000496 IF(TD.EQ.QD) GOTO 15 00000497 CALL FUNC(N,DY,X,Y,J) 00000498 15 DO 21 RD=1,N,1 00000499 GOTO (120,121,122,123,124,125,126),TDP 00000500 120 RK(5*N*RD)=Y(RD) 00000501 GOTO 21 00000502 121 RK(RD)=HD*DY(RD) 00000503 H0=0.5*HD 00000504 F0=0.5*RK(RD) 00000505 GOTO 20 00000506 122 RK(N+RD)=HD*DY(RD) 00000507 F0=0.25*(RK(RD)+RK(N+RD)) 00000508 GOTO 20 00000509 123 RK(2*N+RD)=HD*DY(RD) 00000510 H0=HD 00000511 F0=-RK(N+RD)+2.*RK(2*N+RD) 00000512 GOTO 20 00000513 124 RK(3*N+RD)=HD*DY(RD) 00000514 H0=0.66666666667*HD 00000515 F0=(7.*RK(RD)+10.*RK(N+RD)+RK(3*N+RD))/27. 00000516 GOTO 20 00000517 125 RK(4*N+RD)=HD*DY(RD) 00000518 H0=0.2*HD 00000519 F0=(28.*RK(RD)-125.*RK(N+RD)+546.*RK(2*N+RD)+54.*RK(3*N+RD)- 00000520 1378.*RK(4*N+RD))/625. 00000521 GOTO 20 00000522 126 RK(6*N+RD)=HD*DY(RD) 00000523 F0=0.1666666667*(RK(RD)+4.*RK(2*N+RD)+RK(3*N+RD)) 00000524 X=X0+HD 00000525 ER=(-42.*RK(RD)-224.*RK(2*N+RD)-21.*RK(3*N+RD)+162.*RK(4*N+RD) 00000526 1+125.*RK(6*N+RD))/67.2 00000527 YN=RK(5*N+RD)+F0 00000528 IF(ABS(YN).LT.1E-8) YN=1 00000529 ER=ABS(ER/YN) 00000530 IF(ER.GT.G0) GOTO 115 00000531 IF(ED.GT.ER) GOTO 20 00000532 QD=-1 00000533 20 Y(RD)=RK(5*N+RD)+F0 00000534 21 CONTINUE 00000535 22 CONTINUE 00000536 IF(QD.LT.0)GOTO 23 00000537 IF(U0.LT.0)GOTO 23 00000538 H=2.*H 00000539 23 F0=XD-X 00000540 IF(H.GT.0) GOTO 25 00000541 F0=-F0 00000542 25 IF(F0.GT.0.)GOTO 3 00000543 RETURN 00000544 115 DO 24 RD=1,N 00000545 DY(RD)=RK(RD)/HD 00000546 24 CONTINUE 00000547 H=0.5*HD 00000548 QD=1 00000549 GOTO 11 00000550 100 WRITE(OUT,101)H,XD,X,VO 00000551 101 FORMAT(19H VINT HAS FAILE H=,E11.4,3HXD=,E11.4,2HX=,E11.4,3HV0=, 00000552 1I4) 00000553 STOP 00000554 102 RETURN 00000555 END 00000556 00000557 SUBROUTINE SFCR(HALFC,RK,AAA) 00000558 REAL Y1,Y2,K1,K2,AAA,NUM,DEN 00000559 Y1=HALFC 00000560 K1=2.439*ALOG(RK)-3.0-(1.0/SQRT(HALFC)) 00000561 K2=1.2195 00000562 5 NUM=SQRT(Y1)+Y1*(K1+(K2*ALOG(Y1))) 00000563 DEN=K1/2.0+K2*(1.0+(ALOG(Y1))/2.0) 00000564 Y2=Y1-NUM/DEN 00000565 IF (ABS((Y2-Y1)/Y2).LT.0.000001) GOTO 10 00000566 Y1=Y2 00000567 GOTO 5 00000568 10 AAA=Y2/HALFC 00000569 RETURN 00000570 END 
