我从 1988 年的报告中复制了一个 FORTRAN 66/IV 程序，我正在尝试用 gfortran（windows 的 mingw）编译它。我已将一长串错误减少到 3 个（加上 2 个警告），尽管我尽了最大努力，但我无法再进一步。我将不胜感激任何帮助和建议。
错误：
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)
这是源代码：
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