home *** CD-ROM | disk | FTP | other *** search
/ Programmer 7500 / MAX_PROGRAMMERS.iso / INFO / FORTRAN / D2SOURCE.ZIP / PLR3.FOR < prev    next >
Encoding:
Text File  |  1985-11-18  |  4.6 KB  |  151 lines
  1. C  STACK PROGRAM/ PLUME RISE/ GROUND FIRE/ CGW
  2.       SUBROUTINE PLRS (UR,TMP,PMM,IL,IM,IR,XI,HT,HML,IPC,IRTP)
  3.       COMMON PDM(60),HS,DS,TSC,VS,RDE,P,HR,CR,PD2(42)
  4.       DIMENSION DELF(6),DTDZ(6),PT(6,3),ITT(12) 
  5.       DATA DXT /10./
  6.       DATA G /9.8/  
  7.       DATA ZR /2./  
  8.       DATA CP /.24/ 
  9.       DATA GI /.64/ 
  10.       DATA GC /.5/  
  11.       DATA DELF/1.2,1.1,1.,1.,.9,.8/
  12.       DATA DTDZ/-.01,-.008,-.006,0.,.01,.037/   
  13.       DATA PT/.05,.05,.1,.1,.15,.15,.05,.1,.15,.2,.25,.3,.1,.15,.2,.25, 
  14.      1.3,.35/   
  15.       DATA ITT/3,2,3,1,3,2,3,2,2,3,3,0/ 
  16.       IF (XI.EQ.0.) GO TO 30
  17.       IF (IPC.GT.1) GO TO 10
  18.       IRTP=1
  19.       RETURN
  20.   10  X=XI  
  21.       IF (XI.LT.XMX) GO TO 20   
  22.       X=XMX 
  23.       IRTP=IRTP+1   
  24.   20  IF (IR-6) 280,290,320
  25.   30  IRTP=-1   
  26.       TA=TMP+273.   
  27.       S=G*DTDZ(IM)/TA   
  28.       IF (IL.NE.12) GO TO 60
  29. C     16. ATMOSPHERIC PRESSURE  
  30.   40  CALL READA (16,IRT,IA,PMM)
  31.       IF (IRT.LT.0) GO TO 40
  32.   60  PA=1013.*PMM/760. 
  33. C     23. OUTPUT CONTROL CODE   
  34.   50  CALL DEF (23,IRT) 
  35.       IF (IRT.EQ.0) READ(*,*,ERR=70) IPC
  36.       GO TO 80  
  37.   70  CALL DEF (63,IRT) 
  38.       GO TO 50  
  39.   80  IF (IR.GT.6) GO TO 300
  40. C     24. HEIGHT OF STACK   
  41.   90  CALL READA (24,IRT,IA,HS) 
  42.       IF (IRT.LT.0) GO TO 90
  43. C     25. DIAMETER OF STACK 
  44.  100  CALL READA (25,IRT,IA,DS) 
  45.       IF (IRT.LT.0) GO TO 100   
  46. C     26. TEMPERATURE OF STACK  
  47.  110  CALL READA (26,IRT,IA,TSC)
  48.       IF (IRT.LT.0) GO TO 110   
  49. C     27. VELOCITY OF EFFLUENT  
  50.  120  CALL READA (27,IRT,IA,VS) 
  51.       IF (IRT.LT.0) GO TO 120   
  52. C     28. RELATIVE DENSITY OF EFFLUENT  
  53.  130  CALL DEF (28,IRT) 
  54.       IF (IRT.EQ.0) READ(*,*) RDE   
  55.  140  TS=TSC+273.   
  56.       F=0.  
  57.       IF (TS.LT.TA) WRITE(*,150)
  58.  150  FORMAT(' DHH/DHB/DHBT NOTE: STK TMP LESS THAN AIR TMP')   
  59.       IF (TS.LT.TA) GO TO 160   
  60.       F=(TS-TA)/TS*G*VS*DS**2/4.
  61.       IF (F.LE.55.) XA=14.*F**.625  
  62.       IF (F.GT.55.) XA=34.*F**.4
  63.       XMX=3.5*XA
  64.  160  IT=ITT(IL)
  65.       IF (IT.NE.0) P=PT(IM,IT)  
  66.       IF (IT.NE.0) GO TO 170
  67. C     29. FROST PROFILE EXPONENT
  68.       CALL DEF (29,IRT) 
  69.       IF (IRT.EQ.0) READ(*,*) P 
  70.  170  UZ=UR*(HS/ZR)**P
  71.       IF (S.GT.0.) XMX=2.4*UZ/S**.5 
  72.       FM=RDE*VS*VS*DS*DS/4. 
  73.       IF (UZ.LT.1..AND.S.GE.0.) GO TO 200   
  74.       VR=VS/UZ  
  75.       IF (VR.LT.4) WRITE(*,180) 
  76.  180  FORMAT (' DHJ NOTE: VS/UZ LT 4')  
  77.       IF (S.LT.0.) WRITE(*,190) 
  78.  190  FORMAT (' DHJ NOTE: UNSTABLE MET CONDITIONS') 
  79.       DHJ=3.*VR*DS  
  80.       GO TO 210 
  81.  200  DHJ=4.*(FM/S)**.25
  82.  210  IF (S.LE.0.) GO TO 220
  83.       DHJB=1.5*(FM/UZ)**.333/S**.167
  84.       IF (DHJB.LT.DHJ) DHJ=DHJB 
  85.  220  X=1.  
  86.       IF (IR.EQ.5.AND.IPC.EQ.0) X=XMX   
  87.       DELH=(VS*DS/UR)*(1.5+(2.68E-3*PA*((TS-TA)/TS)*DS))
  88.       DHH=DELH*DELF(IM) 
  89.       WRITE(*,230)  
  90.  230  FORMAT (/8X,'X',8X,'DHH',7X,'DHB',6X,'DHBT',5X,'DHJ') 
  91.       DEL=1.6*(F**.333)/UZ  
  92.  240  IF (IR.EQ.5.AND.X.GT.XMX) X=XMX   
  93.       DHB=DEL*X**.667
  94.       IF (S.LE.0.) GO TO 250
  95.       DHB=2.5*(F/(UZ*S))**.333  
  96.       IF (UR.GE.1.) GO TO 250   
  97.       DHMTT=5.*(F**.25)/(S**.375)   
  98.       IF (DHMTT.LT.DHB) DHB=DHMTT   
  99.  250  DHJX=1.44*DS*(VS/UZ)**.667*(X/DS)**.333   
  100.       IF (DHJX.GT.DHJ) DHJX=DHJ 
  101.       DHBT=F**.333*X**.667/UZ*(1.065-6.25*DTDZ(IM)) 
  102.       IF (IPC.EQ.1.OR.IPC.EQ.3) WRITE(*,260) X,DHH,DHB,DHBT,DHJX
  103.       IF (IR.EQ.5.AND.X.GE.XMX) GO TO 270   
  104.       IF (IR.EQ.6.AND.DHJX.GE.DHJ) GO TO 270
  105.  260  FORMAT (6F10.2)   
  106.       IF (X.EQ.1.) X=0. 
  107.       X=X+DXT   
  108.       GO TO 240 
  109.  270  WRITE(*,260) X,DHH,DHB,DHBT,DHJ,P 
  110.       IF (IR.EQ.5) HT=HS+DHBT   
  111.       IF (IR.EQ.6) XMX=X
  112.       IF (IR.EQ.6) HT=HS+DHJ
  113.       RETURN
  114.  280  DHBT=F**.333*X**.667/UZ*(1.065-6.25*DTDZ(IM)) 
  115.       HT=HS+DHBT
  116.       GO TO 360
  117.  290  DHJX=1.44*DS*(VS/UZ)**.667*(X/DS)**.333   
  118.       IF (DHJX.GT.DHJ) DHJX=DHJ 
  119.       HT=HS+DHJX
  120.       GO TO 360 
  121. C     30. HEAT RELEASED 
  122.  300  CALL READA (30,IRT,IA,HR) 
  123.       IF (IRT.LT.0) GO TO 300   
  124. C     31. CLOUD RADIUS  
  125.  305  CALL READA (31,IRT,IA,CR) 
  126.       IF (IRT.LT.0) GO TO 305   
  127.       HT=0. 
  128.       IF (IM.LT.4) GO TO 330
  129.       IF (IM.EQ.4) S=G*3.322E-4/TA  
  130.       ROA=352320.*PMM/760./TA   
  131.       RTS=S**.5 
  132.       XMX=3.14159*UR/RTS
  133.       X=XMX 
  134.       WRITE(*,310) XMX  
  135.  310  FORMAT (' XMX=',F7.0) 
  136.  320  FT=1.-COS(RTS*X/UR)   
  137.       FB=G*HR/(3.14159*ROA*CP*TA)   
  138.       IF (IR.EQ.8) GO TO 350
  139.       HT=(3.*FB*FT/(S*GI**3)+(CR/GI)**4)**.25-(CR/GI)
  140.       GO TO 360
  141.  330  WRITE(*,340)
  142.  340  FORMAT (' HEIGHT DEFINED FOR STABLE CONDITIONS ONLY')
  143.       RETURN
  144.  350  HT=(3.*FB*FT/(UR*GC*GC*S)+(CR/GC)**3)**.333-(CR/GC)
  145.  360  IF (HT.LT.HML) GO TO 370
  146.       HT=HML
  147.  370  WRITE(*,380) HT
  148.  380  FORMAT(' HTS=',F7.2)
  149.       RETURN
  150.       END
  151.