C**** VERSION:10/28/88 C By F.D. Drake ! This component has been assigned Type Number 201. If that number conflicts with ! another user Type number, you will need to change it and recompile the appropriate ! dll. SUBROUTINE TYPE201(TIME,XIN,OUT,T,DTDT,PAR,INFO,ICNTRL,*) !DEC$ATTRIBUTES DLLEXPORT :: TYPE201 C******************************************************************** C**** C**** THIS SUBROUTINE CALCULATES THE HEAT LOSSES AND PRESSURE DROPS C**** IN A STEAM PIPELINE C**** C**** SAMPLE INPUT AND OUTPUT FOLLOW THE NOMENCLATURE C**** C******************************************************************** C******************************************************************** C NOMENCLATURE: (IN ALPHABETICAL ORDER) C C ALPHA : THERMAL DIFFUSITIVY C CF : FRICTION COEFFICIENT C CP : SPECIFIC HEAT OF STEAM IN J/KG/K C DE : OUTER DIAMETER OF PIPE IN M OR FT C DELP : PRESSURE DROP IN MPA OR PSI C DI : INNER DIAMETER OF PIPE IN M OR FT C DINS : OUTER DIAMETER OF INSULATION IN M OR FT C FAK : FACTOR C HEXT : EXTERIOR HEAT TRANSFER COEFF. IN W/M**2*K OR BTU/HR/FT^2/F C K : THERMAL CONDUCTIVITY OF STEAM IN W/MK OR BTU/FT/F C KINS : THERMAL " OF INSULATION IN W/MK OR BTU/FT/F C KP : THERMAL CONDUCTIVITY OF PIPE IN W/MK OR BTU/FT/F C L : LENGTH OF PIPE IN M OR FT C MDOT : MASS FLOW RATE IN KG/S OR LB/HR C MU : STEAM VISCOSITY IN PA*S C MVP : MEAN VALUE OF PRESSURE IN PIPE C MVT : MEAN VALUE OF TEMPERATURE IN PIPE C NU : NUSSELT-NUMBER C PEN : PECLET-NUMBER C PINP : INPUT PRESSURE IN MPA OR PSIA C POUTP : OUTPUT PRESSURE IN MPA OR PSIA C PRN : PRANDTL-NUMBER C PRO : PIPE ROUGHNESS IN MM OR INCH C QLOSS : TOTAL HEAT LOSS IN W OR BTU/HR C RE : REYNOLDSNUMBER C ROAV : STEAM DENSITY IN KG/M**3 C RTOT : TOTAL THERMAL RESISTANCE IN W/K C TA : AMBIENT TEMP. IN C OR F C TINP : INPUT TEMP. IN C OR F C TOUTP : OUTPUT TEMP. IN C OR F C UNITS : 1: INPUT/OUTPUT IN SI UNITS C : 2: IN/OUT IN BRITISH UNITS C WST : STEAM VELOCITY IN M/S C C************************************************************************** C C SAMPLE INPUT AND OUTPUT VALUES C C PARAMETERS SAMPLE VALUES C PAR(1) MODE 2 C PAR(2) UNIT 0 C PAR(3) LENGTH 13000.0 C PAR(4) DI 2.0 C PAR(5) DE 2.04 C PAR(6) DINS 2.44 C PAR(7) KP 26.0 C PAR(8) KINS .033 C PAR(9) HEXT .93 C PAR(10) PROUGH .0036 C C INPUTS C XIN(1) MSTEAM 200000.0 C XIN(2) TOUT 440.0 C XIN(3) POUT 190.0 C XIN(4) TA 60.0 C C OUTPUTS C OUT(1) MSTEAM 200000.0 C OUT(2) TIN 446.0 C OUT(3) PIN 196.7 C OUT(4) QLOSS 699900.0 C OUT(5) DELTAP 6.726 C C******************************************************************** DOUBLE PRECISION XIN,OUT REAL NU,L,KP,KINS,MU,MDOT,MODE,K1,MU1,MVT,MVTCPH,MVTCPL,K,MVP DIMENSION XIN(4),OUT(5),PAR(10),INFO(15),A(8),B(7) COMMON /LUNITS/LUR,LUW,LUK C**** DATA FOR CORRELATIONS TO DETERMINE KP AND MU DATA(A(I),I=1,8)/1.76E01,5.87E-02,1.04E-04,-4.51E-08, +1.0351E02,4.198E-01,-2.771E-05,2.1482E14/ DATA(B(I),I=1,7)/4.07E-01,8.04E01,1.858E03,-5.90, +3.53E02,6.765E02,1.021E02/ DATA PI /3.14159/ ! Set the version information for TRNSYS IF (INFO(7).EQ.-2) THEN INFO(12) = 15 RETURN 1 ENDIF C******************************************************************** C**** SET PARAMETERS AND INPUTS C**** MODE=1:INPUT PRESSURE AND TEMP. ARE INLET CONDITIONS C**** MODE=2:INPUT PRESSURE AND TEMP. ARE DESIRED OUTLET COND. C******************************************************************** MODE=PAR(1) UNITS=PAR(2) L =PAR(3) DI =PAR(4) DE =PAR(5) DINS=PAR(6) KP =PAR(7) KINS=PAR(8) HEXT=PAR(9) PRO =PAR(10) MDOT=XIN(1) TINP=XIN(2) PINP=XIN(3) TA =XIN(4) C**** UNIT CONVERSION IF(UNITS.EQ.0.) THEN L=L*0.3048 DI=DI*0.3048 DE=DE*0.3048 DINS=DINS*0.3048 KP=KP*0.57782 KINS=KINS*0.57782 HEXT=HEXT*0.17612 PRO=PRO*25.4 MDOT=MDOT/3600*0.454 TINP=(TINP-32)/1.8 PINP=PINP*0.00689 TA=(TA-32)/1.8 END IF C**** CHECK RANGE OF PARAMETERS AND INPUTS IF(MODE.NE.1.AND.MODE.NE.2) GOTO 60 IF(UNITS.NE.1.AND.UNITS.NE.0) GOTO 60 GOTO 62 60 WRITE(LUW,61) TIME, INFO(1), INFO(2) 61 FORMAT(T2,'********** WARNING *********',/, +T2,' TIME : ',E12.5,' UNIT : ',I3,' TYPE : ',I3,/, +'CHECK RANGE OF PARAMTERS AND INPUTS !',//) C**** SET INFO 62 INFO(6)=5. INFO(9)=0. C**** INITIAL VALUES FOR ITERATION POUTP=PINP TOUTP=TINP C******************************************************************** C**** PROPERTY CALCULATIONS C******************************************************************** C**** CALCULATION OF SPECIFIC HEAT CP MVT=TINP MVP=PINP MVTCPH=MVT+10 MVTCPL=MVT-10 CALL STEAM('SI',MVTCPH,MVP,HST1,SST,QST,VSTAV,U,12) CALL STEAM('SI',MVTCPL,MVP,HST2,SST,QST,VSTAV,U,12) CP=(HST1-HST2)/20.*1000. C**** DENSITY CALCULATION CALL STEAM('SI',MVT,MVP,HST,SST,QST,VSTAV,UST,12) ROAV=1./VSTAV C**** CHECK IF(HST.LE.2107.) THEN WRITE(LUW,50) TIME, INFO(1), INFO(2) 50 FORMAT(T2,'********** WARNING *********',/, +T2,' TIME : ',E12.5,' UNIT : ',I3,' TYPE : ',I3,/, + ' WARNING !!! STEAM CONDENSATION IN PIPE !!!',//) CALL TYPECK(2,INFO,0,0,0) END IF C**** CALCULATION OF KP AND MU FROM CORRELATIONS K1=A(1)+A(2)*MVT+A(3)*MVT**2.+A(4)*MVT**3. K=(K1+(A(5)+A(6)*MVT+A(7)*MVT**2.)*ROAV/1000.+A(8)* *(ROAV/1000.)**2./MVT**(4.2))/1000. MU1=B(1)*MVT+B(2) IF(MVT.GT.300.) GOTO 3 MU=(MU1-ROAV/1000.*(B(3)+B(4)*MVT))*1.E-07 GOTO 4 3 MU=(MU1+B(5)*ROAV/1000.+B(6)*(ROAV/1000.)**2.+B(7)* *(ROAV/1000.)**3.)*1.E-07 C******************************************************************** C**** CALCULATION OF DIMENSIONLESS NUMBERS AND OTHER VALUES C******************************************************************** 4 RE=4.*MDOT/(PI*DI*MU) WST=4.*MDOT/(DI**2.*PI*ROAV) ALPHA=K/(CP*ROAV) PRN=MU/(ROAV*ALPHA) PEN=WST*DI/ALPHA C******************************************************************** C******************************************************************** C**** CALCULATION OF HEATLOSSES C******************************************************************** C******************************************************************** C**** CHECK: LAMINAR OR TURBULENT IF (RE.GT.2300.) GOTO 5 C**** CALC. OF NUSSELT NUMBER FROM EMPIRICAL CORRELATION C**** LAMINAR FLOW NU=3.65+(0.0668*PEN*DI/L)/(1.+0.045*(PEN*DI/L)**(2./3.)) GOTO 10 C**** TURBULENT FLOW 5 NU=0.024*RE**0.786*PRN**0.45*(1.+(DI/L)**(2./3.)) C**** HEATLOSS CALCULATIONS 10 HINT=NU*K/DI RTOT=1./(HINT*2.*PI*DI/2.*L)+ALOG(DE/DI)/(2.*PI*KP*L)+ +ALOG(DINS/DE)/(2.*PI*KINS*L)+1./(HEXT*2.*PI*DINS/2*L) IF(MODE.EQ.1) FAK=-1. IF(MODE.EQ.2) FAK=1. TOUTP=TA+(TINP-TA)*EXP(FAK/(RTOT*MDOT*CP)) QLOSS=MDOT*CP*(TINP-TOUTP)*(-FAK) C******************************************************************** C******************************************************************** C**** DETERMINATION OF PRESSURE DROP C******************************************************************** C******************************************************************** IF (RE.LE.2300.) THEN CF=16./RE ELSE IF(RE.LE.(DI*1000/PRO*LOG10(0.1*DI*1000/PRO))) THEN C******************************************************************** C**** CALCUALTION OF FRICTION COEFFICIENTS FROM EMPIRICAL CORRELATIONS C******************************************************************** C**** IN CASE OF SMOOTH PIPE CF=0.25*0.309/(LOG10(RE/7.))**2. END IF IF(RE.GT.(DI*1000./PRO*LOG10(0.1*DI*1000./PRO)).AND.RE.LT. + (400.*DI*1000./PRO*LOG10(3715*DI/PRO))) THEN C**** TRANSITION ZONE CF=0.25*0.0055*(1.+(2.*10.*PRO/DI+10.**6./RE)**(1./3.)) END IF IF(RE.GE.(400.*DI*1000./PRO*LOG10(3715*DI/PRO))) THEN C**** ROUGH PIPE CF=0.25/(2.*LOG10(3715.*DI/PRO))**2. END IF END IF C**** CALC. OF PRESSURE DROP IF (MODE.EQ.1) THEN IF ((4.*CF*L/DI*WST**2./PINP/10.**6.*ROAV).GE.1.) GOTO 30 DELP=PINP*(1.-SQRT(1.-4.*CF*L/DI*WST**2./PINP/10.**6.*ROAV)) POUTP=PINP-DELP ELSE POUTP=SQRT(PINP**2.+(2.*CF*ROAV*L/DI*WST**2./10.**6.)**2.)+ + 2.*CF*ROAV*L/DI*WST**2./10.**6. J=0 25 CALL STEAM ('SI',TOUTP,POUTP,HOUT,SOUT,QOUT,VO,UO,12) RO=1./VO WST=4.*MDOT/(DI**2*PI*RO) DELP=POUTP*(1.-SQRT(1.-4.*CF*L/DI*WST**2./POUTP/10.**6.*RO)) PEXNEW=POUTP-DELP IF(J.GT.2000) GOTO 30 IF(PEXNEW.LE.PINP) GOTO 28 POUTP=POUTP-0.01 J=J+1 GOTO 25 28 DELP=POUTP-PINP END IF C**** CHECK IF FLOW RATE POSSIBLE IF (POUTP.LE.0) THEN 30 WRITE(LUW,52) TIME, INFO(1), INFO(2) 52 FORMAT(T2,'********** WARNING *********',/, +T2,' TIME : ',E12.5,' UNIT : ',I3,' TYPE : ',I3,/, + 'MASSFLOW*PIPEROUGHNESS/DIAMETER RATIO IS TOO HIGH ;',/ +,'OPERATION OF PIPE IMPOSSIBLE: CHOSE NEW VALUES !',//) CALL TYPECK(2,INFO,0,0,0) END IF C**** CHECK IF CONDENSATION IN PIPE CALL STEAM ('SI',TOUTP,POUTP,HOUT,SOUT,QOUT,VOUT,UOUT,12) IF(HOUT.LE.2107.) THEN WRITE(LUW,54) TIME, INFO(1), INFO(2) 54 FORMAT(T2,'********** WARNING *********',/, +T2,' TIME : ',E12.5,' UNIT : ',I3,' TYPE : ',I3,/, + ' WARNING !!! STEAM CONDENSATION IN PIPE !!!',//) CALL TYPECK(2,INFO,0,0,0) END IF C******************************************************************** C**** SET OUTPUT C******************************************************************** 40 IF(UNITS.EQ.1.) THEN OUT(1)=MDOT OUT(2)=TOUTP OUT(3)=POUTP OUT(4)=QLOSS OUT(5)=DELP ELSE OUT(1)=MDOT*3600/0.454 OUT(2)=TOUTP*1.8+32 OUT(3)=POUTP/0.00689 OUT(4)=QLOSS*3.4123 OUT(5)=DELP/0.00689 END IF RETURN 1 END CSTEAM CGiven any two steam properties, this routine evaluates the remaining Cproperties associated with the state. A call to this routine has Cthe following form: C C CALL STEAM(UNITS,T,P,H,S,Q,V,U,ITYPE) C Cwhere UNITS is a character variable (declared as CHARACTER*2) and/or set Cequal to "SI" or "EN". T, P, H, S, Q, V, and U are the steam temperature C(deg. C, deg. F), pressure (MPa, psia), enthalpy (kJ/kg, Btu/lbm), Centropy (kJ/kg-deg.C, Btu/lbm-deg.F), quality, specific volume C(m**3/kg, ft**3/lbm), and internal energy (kJ/kg, Btu/lbm), and ITYPE Crefers to which of the 7 arguments are specified as known proper- Cties. For example, in the following call to STEAM, properties for steam Cwould be evaluated at a temperature of 1000 degrees fahrenheit and pressure Cof 900 psi absolute. C C CALL STEAM('EN',1000.,900.,H,S,Q,V,U,12) C CA quality larger than 1 is a flag that the steam is a superheated Cvapor. In order to specify a saturated vapor or liquid condition, use Ca quality of 1. or 0. along with any additional property. C SUBROUTINE STEAM(UNITS,T,P,H,S,Q,V,U,ITYPE) CHARACTER*2 UNITS C C UNITS = 'EN' OR 'SI' C T : TEMPERATURE (DEG. C OR F) C P : PRESSURE (MPA OR PSIA) C H : ENTHALPY (KJ/KG OR BTU/LBM) C S : ENTROPY (KJ/KG K OR BTU/LBM R) C Q : QUALITY (01 FOR SUPERHEATED VAPOR C KI=0 20100 I=ITYPE/10 J=ITYPE-I*10 K=I IF(I.LT.J) GO TO 200 I=J J=K 200 IF (I.GT.6.OR.J.GT.7) GO TO 202 M=1 IF (UNITS.NE.'SI') GO TO 11110 GO TO (11001,11002,11003,11004,11100,11006,11007),I 11009 GO TO (11001,11002,11003,11004,11100,11006,11007),J 11110 GO TO (210,220,230,240,250,260),I 202 WRITE (LUW,201) ITYPE 201 FORMAT(' ITYPE OUT OF RANGE. ITYPE=',I10) CALL MYSTOP(1001) RETURN 210 J=J-1 GO TO (11,12,13,51,5000,8500),J 220 J=J-2 GO TO (14,15,52,4000,202), J 230 J=J-3 GO TO (16,5400,7400,202),J 240 J=J-4 GO TO (4400,6400,202),J 250 J=J-5 GO TO (3400,8000),J 260 GO TO 202 C**************** GIVEN TEMPERATURE AND PRESSURE ******************************* 11 IF ( T - TSL(P) ) 64, 1, 2 64 IF ( ABS( T-TSL(P) ) ) 1, 1, 65 1 Q = 1. H = HG(T) S = SGT(T) V=VG(T,P) U=H-(P*V)*0.18509 GO TO 100 2 Q=1.5E+38 H = HTP(T,P) S = STP(T,P) V=VG(T,P) U=H-(P*V) *0.18509 GO TO 100 C SAT. LIQUID- T,P 65 H = HF(T) S = SFT(T) Q = 0. V=VF(T) U=H-(P*V) *0.18509 GO TO 100 C***************** GIVEN TEMPERATURE AND ENTHALPY ****************************** 12 IF ( H - HG(T) ) 4, 4, 3 C SUPER HEAT - GUESS P 3 P1 = 100. DO 70 I=1, 10 P2 = P1* 1.01 HTPP1 = HTP(T, P1) IF( ABS (H-HTPP1)- .01) 5,5,6 6 P1 = (( H - HTPP1)/(HTP(T, P2) - HTPP1))* ( P2- P1) + P1 70 CONTINUE 5 P = P1 S = STP(T,P) V=VG(T,P) U=H-(P*V) *0.18509 Q=1.5E+38 GO TO 100 C SATURATED VAPOR OR MIXTURE 4 IF ( H - HG(T) ) 10, 9, 9 10 HFT1 = HF(T) IF ( H- HFT1) 7, 7, 8 C LIQUID 7 S = SFT(T) P = PSL(T) Q = 0. V=VF(T) U=H-(P*V) *0.18509 GO TO 100 8 Q = ( H - HFT1) /(HG(T) - HFT1) S = SFT(T) + Q*( SGT(T)-SFT(T)) P = PSL(T) V= VF(T)+(Q*(VG(T,P)-VF(T))) U=H-(P*V) *0.18509 GO TO 100 9 S = SGT(T) P = PSL(T) V=VG(T,P) U=H-(P*V) *0.18509 Q = 1. GO TO 100 C******************* GIVEN TEMPERATURE AND ENTROPY ***************************** C GUESS P GIVEN (T,S)- SUP. HEAT 13 IF ( S - SGT(T) ) 24, 24, 23 23 P1 = 100. DO 71 I=1, 10 P2 = P1*1.01 STPP1 = STP(T, P1) IF (ABS(S-STPP1)-.0000001) 21,21,22 22 P1 = P1 + (( S -STPP1)/( STP(T,P2) - STPP1)) * ( P2 - P1) 71 CONTINUE 21 P=P2 H = HTP(T,P) V=VG(T,P) U=H-(P*V) *0.18509 Q=1.5E+38 GO TO 100 C SATURATED REGION 24 IF( S - SGT(T) ) 26, 255, 255 26 IF( S- SFT(T)) 265, 265, 25 265 Q = 0. P = PSL(T) H = HF(T) IF (KI.NE.1) V=VF(T) U=H-(P*V) *0.18509 GO TO 100 255 H = HG(T) Q = 1. P = PSL(T) V=VG(T,P) U=H-(P*V) *0.18509 GO TO 100 25 Q = ( S - SFT(T))/( SGT(T)-SFT(T)) H = HF(T) + Q * (HG(T)- HF(T)) P = PSL(T) V= VF(T)+Q*(VG(T,P)-VF(T)) U=H-(P*V) *0.18509 GO TO 100 14 IF ( H - HG2(P) ) 28, 28, 27 C****************** GIVEN PRESSURE AND ENTHALPY ******************************** C SUP. HEAT - P, H. 27 T = THP(H,P) S = STP(T,P) Q=1.5E+38 V=VG(T,P) U=H-(P*V) *0.18509 GO TO 100 28 IF ( ABS( H - HG2(P) ) ) 29, 29, 30 29 T = TSL(P) S = SGP(P) V=VG(T,P) U=H-(P*V) *0.18509 Q = 1. GO TO 100 30 T = TSL(P) HFT1 = HF(T) IF ( H-HFT1) 31, 31, 32 32 SFP1 = SFP(P) Q = ( H- HFT1) / ( HG(T)- HFT1) S = SFP1 + Q * ( SGP(P) - SFP1) V= VF(T)+Q*(VG(T,P)-VF(T)) U=H-(P*V) *0.18509 GO TO 100 31 S = SFP(P) Q = 0. V=VF(T) U=H-(P*V) *0.18509 GO TO 100 C SUP. HEAT GIVEN P AND S. 34 H = HPS(P,S) T = THP(H,P) V=VG(T,P) U=H-(P*V) *0.18509 Q=1.5E+38 GO TO 100 C VAPOR OR MIXTURE GIVEN P AND S 33 IF( ABS( S-SGP(P)) ) 37, 37, 38 C VAPOR, GIVEN P AND S. 37 T = TSL(P) H = HG2(P) V=VG(T,P) U=H-(P*V) *0.18509 Q = 1. GO TO 100 38 IF ( S-SFP(P)) 135,35,36 C LIQUID 35 T = TSL(P) H = HF(T) V=VF(T) U=H-(P*V) *0.18509 Q = 0. GO TO 100 C****************** GIVEN PRESSURE AND ENTROPY ********************************* 15 IF ( S - SGP(P) ) 33, 33, 34 36 Q = ( S - SFP(P))/ ( SGP(P)- SFP(P)) T = TSL(P) H = HF(T) + Q*(HG(T)-HF(T)) V= VF(T)+ Q*(VG(T,P)-VF(T)) U=H-(P*V) *0.18509 GO TO 100 C COMPRESSED LIQUID 135 T1=100 DO1000 I=1,10 T2=T1*1.01 IF ((ABS(S-SFT(T1))).LT.0.00001) GO TO 1001 1000 T1=T1+(((S-SFT(T1))/ (SFT(T2)-SFT(T1)))*(T2-T1)) 1001 T=T1 IF (KI.EQ.1) GO TO 6050 H=HF(T) IF (KI.EQ.1) GO TO 1007 V=VF(T) 1007 IF (KI.EQ.1) P=PSL(T) U=H-(P*V) *0.18509 Q=0 GO TO 100 C****************** GIVEN ENTHALPY AND ENTROPY ********************************* C GUESS P GIVE H AND S , SUP. HEAT 16 P1 = 100. DO 72 I =1,10 P2 = P1* 1.01 IF ( ABS( H - HPS( P1,S)) - .01) 39, 39, 40 40 P1= P1 + (( H-HPS(P1, S))/( HPS(P2, S) -HPS( P1,S)))*(P2- P1) IF ( P1.LE. 0.) GO TO 107 72 CONTINUE GO TO 39 107 WRITE (LUW,56) 56 FORMAT( '0', ' GIVEN ENTHALPY AND ENTROPY , TEMP. AND PRESS. 1 CAN BE SPECIFIED FOR SUP. HEATED STEAM ONLY') GO TO 75 39 P = P1 T = THP(H,P) V=VG(T,P) U=H-(P*V) *0.18509 Q=1.5E+38 GO TO 100 C****************** GIVEN TEMPERATURE AND QUALITY ****************************** C T AND Q 51 P = PSL(T) 54 H1 = HF(T) H = H1 + Q* ( HG(T) - H1) S1 = SFT(T) S = S1 + Q*(SGT(T) - S1) V=VF(T)+Q*(VG(T,P)-VF(T)) U=H-(P*V) *0.18509 GO TO 100 C****************** GIVEN PRESSURE AND QUALITY ********************************* C P AND Q 52 T = TSL(P) GO TO 54 C****************** GIVEN QUALITY AND SPECIFIC VOLUME ************************** C********************** GUESS TEMPERATURE ************************************** 3400 T=100.0 DO 3500 I=1,10 T2= T*1.01 VFT=VF(T) V1= VFT+ (VG(T,PSL(T))-VFT)*Q IF ((ABS(V-V1))-0.0001 ) 3510,3510,3501 3501 VFT2=VF(T2) V2= VFT2 + (VG(T2,PSL(T2))-VFT2)*Q 3500 T=T+ ((V-V1)/(V2-V1))*(T2-T) WRITE (LUW,3525) 3525 FORMAT (10X, '***** EQUATIONS DID NOT CONVERGE *****' ) 3510 P=PSL(T) S=SFT(T)+ (SGT(T)-SFT(T))*Q H=HF(T)+(HG(T)-HF(T))*Q U=H-(P*V)*0.18509 GO TO 100 C************ GIVEN PRESSURE AND VOLUME *************************************** 4000 T=TSL(P) IF (VG(T,P)-V) 2115,2114,2116 2115 IF ((ABS(V-VG(T,P))-0.00100)) 2114,2114,2111 2116 IF (V-VF(T)) 2117,2117,2110 2117 IF (ABS(V-VF(T))-0.0001) 2113,2113,2112 2110 Q=(V-VF(T))/(VG(T,P)-VF(T)) H= HF(T)+Q*(HG2(P)-HF(T)) S= SFP(P) + Q*(SGP(P)-SFP(P)) U=H-(P*V)*0.18509 GO TO 100 2113 Q=0 H=HF(T) S=SFT(T) U=H-(P*V)*0.18509 GO TO 100 2114 Q=1 H=HG(T) S=SGT(T) U=H-(P*V)*0.18509 GO TO 100 2112 WRITE (LUW,4009) 4009 FORMAT (10X,'***** STATE IN COMPRESSED LIQUID REGION. STEAM2 CANNO 1T SPECIFY STATE. *******') GO TO 100 C*********** GUESS TEMPERATURE GIVEN PRESSURE AND VOLUME *********************** 2111 T1=T DO 2121 IVP=1,10 T2=T1*1.01 IF (ABS(V-VG(T1,P)-0.00001))2122,2122,2121 2121 T1=T1+ ((V-VG(T1,P))/(VG(T2,P)-VG(T1,P)))*(T2-T1) 2122 T=T1 H= HTP(T,P) S= STP (T,P) Q= 1.5E+38 U=H-(P*V)*0.18509 GO TO 100 C************ GIVEN QUALITY AND ENTROPY **************************************** C***************** GUESS TEMPERATURE ******************************************* 4400 T=100.0 DO 4500 I=1,10 T2=T*1.01 SFTT=SFT(T) S1= SFTT+(SGT(T)-SFTT)*Q IF ((ABS(S-S1))-0.00001) 4510 ,4510 , 4501 4501 SFTT2= SFT(T2) S2=SFTT2 + (SGT(T2)-SFTT2)*Q 4500 T=T+ ((S-S1)/(S2-S1))*(T2-T) WRITE (LUW,3525) 4510 P=PSL(T) H=HF(T)+(HG(T)-HF(T))*Q V=VF(T)+(VG(T,P)-VF(T))*Q U=H-(P*V)*0.18509 GO TO 100 C************** GIVEN TEMPERATURE AND VOLUME *********************************** 5000 IF (V-VF(T)) 1002,1002,1020 1002 P=PSL(T) H=HF(T) S=SFT(T) U=H-(P*V)*0.18509 Q=0 GO TO 100 1020 P=PSL(T) IF (V-VG(T,P)) 1010 , 1011 , 1012 1010 Q= (V-VF(T))/(VG(T,P)-VF(T)) IF (ABS(Q-1).LT.0.0005) GO TO 1011 H= HF(T)+Q*(HG(T)-HF(T)) S= SFT(T)+ Q* (SGT(T)-SFT(T)) U=H-(P*V)*0.18509 GO TO 100 1011 S=SGT(T) H=HG(T) Q=1 U=H-(P*V)*0.18509 GO TO 100 1012 P=PG(T,V) S=STP(T,P) H=HTP(T,P) Q=1.5E+38 U=H-(P*V)*0.18509 GO TO 100 C*************** GIVEN QUALITY AND ENTHALPY ************************************ C******************** GUESS TEMPERATURE **************************************** 5400 T=100.0 DO 5500 I=1,10 T2= T*1.01 HFT=HF(T) H1= HFT+(HG(T)-HFT)*Q IF ((ABS(H-H1))-0.00001) 5510,5510,5501 5501 HFT2= HF(T2) H2= HFT2+ (HG(T2)-HFT2)*Q 5500 T=T+((H-H1)/(H2-H1))*(T2-T) WRITE (LUW,3525) 5510 P=PSL(T) S=SFT(T)+(SGT(T)-SFT(T))*Q V=VF(T)+(VG(T,P)-VF(T))*Q U=H-(P*V)*0.18509 GO TO 100 C*********** GIVEN ENTROPY AND SPECIFIC VOLUME ********************************* 6400 T=100 NTHRU=0 KI=1 6401 P=PSL(T) NTHRU=NTHRU+1 IF (NTHRU.EQ.26) GO TO 6600 SG1=SGT(T) VG1=VG(T,P) IF ((ABS(SG1-S)).LE.0.0100) GO TO 6000 6402 IF ((SG1-S).LE.0.0) GO TO 6199 VF1=VF(T) SF1=SFT(T) IF ((ABS(SF1-S)).LE.0.0150) GO TO 6300 6403 IF ((SF1-S).LE.0.0) GO TO 6500 GO TO 135 6050 VFF=VF(T) IF ((V-VFF).GE.0.0) GO TO 6100 H=HF(T) Q=0.0 P=PSL(T) U=H-(P*V)*0.18509 GO TO 100 6000 IF ((ABS(VG1-V)).GT.0.0100) GOTO 6402 H=HG(T) U=H-(P*V)*0.18509 Q=1 GO TO 100 6199 P1=PSL(T) 6195 DO 6198 I=1,10 STPP1=STP(T,P1) IF ((ABS( S-STPP1))-0.00100) 6197,6197,6196 6196 P2=P1*0.9 PDIF= ((S-STPP1)/( STP(T,P2) -STPP1))*(P2-P1) IF (P1+PDIF) 6194,6198,6198 6198 P1=P1+PDIF 6194 P1=P1-(0.5*P1) GO TO 6195 6197 P=P1 6200 V1=VG(T,P) IF ((ABS(V-V1)).LE.0.0100) GO TO 6210 T2=T*1.01 P=PSL(T2) 6202 DO 6204 I=1,10 P2=P*0.99 STPP=STP(T2,P) IF((ABS (S-STPP)).LE.0.001) GO TO 6205 PDFF=((S-STPP)/(STP(T2,P2)-STPP))*(P2-P) IF (P+PDFF) 6203,6204,6204 6204 P=P+PDFF 6203 P=P-P/2 GO TO 6202 6205 V2=VG(T2,P) T=T+((V-V1)/(V2-V1))*(T2-T) GO TO 6401 6210 H=HTP(T,P) U=H-(P*V)*0.18509 Q=1.5E+38 GO TO 100 6300 IF ((ABS(VF1-V)).GT.0.0100) GO TO 6403 P=PSL(T) H=HF(T) U=H-(P*V)*0.18509 Q=0.0 GO TO 100 6500 Q= (S-SF1)/(SG1-SF1) V1=Q*(VG1-VF1)+VF1 IF ((ABS(V1-V)).LE.0.0100) GO TO 6510 T2=T*1.01 P2=PSL(T2) V2=Q*(VG(T2,P2)-VF(T2))+VF(T2) T=T+((V-V1)/(V2-V1))*(T2-T) GO TO 6401 6510 H=Q*(HG(T)-HF(T)) +HF(T) P=PSL(T) U=H-(P*V)*0.18509 GO TO 100 6100 T2=T*1.01 VF2=VF(T2) T=T+((V-VFF)/(VF2-VFF))*(T2-T) GO TO 6401 6600 WRITE (LUW,3525) GO TO 100 7400 GO TO 202 8000 GO TO 202 8500 GO TO 202 C*********************** UNIT CONVERSIONS ************************************** 11001 T=T*1.8+32.0 GO TO 11100 11002 P=P/6894.757E-6 GO TO 11100 11003 H=H/2.326 GO TO 11100 11004 S=S/4.1868 GO TO 11100 11006 V=V/0.06242797 GO TO 11100 11007 U=U/2.326 GO TO 11100 11100 M=M+1 IF (M.EQ.3) GO TO 11110 GO TO 11009 100 CONTINUE IF (UNITS.NE.'SI') GO TO 75 P=P*6894.757E-6 S=S*4.1868 V=V*0.06242797 T=(T - 32.0)/1.8 H=H*2.326 U=U*2.326 75 RETURN END C********************* STEAM2 FUNCTIONS **************************************** C************* ENTROPY OF A FLUID AS A FUNCTION OF PRESSURE ******************** FUNCTION SFP(P) TB = ( TSL(P)- 360.0)/3100. SFP = 0.515755 + (3.96796-(4.59799-(34.2517-(60.7233+(367.036- 1 (12035.9+123466.0*TB)*TB)*TB)*TB)*TB)*TB)*TB RETURN END C************* ENTROPY OF A FLUID AS A FUNCTION OF TEMPERATURE ***************** FUNCTION SFT(T) TB = ( T-360.0)/3100.0 SFT = 0.515755+(3.96796-(4.59799-(34.2517-(60.7233+(367.036- 1 (12035.9+123466.0*TB)*TB)*TB)*TB)*TB)*TB)*TB RETURN END C************* ENTROPY OF A SUP. HEATED GAS AS FUNCTION OF TEMP. AND PRESS.***** FUNCTION STP(T,P) T1= T P1 =P T =255.38+T/1.8 P = P/14.6959 B1 = (2641.62 * 10. **(80870.0/(T*T)))/T B0= 1.89-B1 B2= 82.546 B3= 162460.0/T B4= 0.21828*T B5= 126970.0/T F0= 1.89-B1*(372420.0/(T*T)+2.0) B6 = B0*B3- 2.0*F0*(B2-B3) B7= 2.0*F0*(B4-B5) -B0*B5 B8= 0.43429448*ALOG(T) B9= B0*P*P/(2.0* (T*T) ) BETA = ((B0-F0)*P+B9*(B6+B9*B0*(B0*(B4-B5)-2.0*B7)))/T STP = 0.809691*B8-0.253801*0.43429448* ALOG(P)+0.00018052*T 1 -11.4276/T- 0.355579-0.0241983*BETA T=T1 P = P1 RETURN END C******** ENTHALPY OF SUP. HEATED GAS AS FUNCTION OF TEMP. AND PRESS.********** FUNCTION HTP(T,P) T1 =T P1 = P T = 255.38 + T/1.8 P = P/14.6959 B1 = ( 2641.62 *10.**( 80870.0/ ( T*T)))/T B0 = 1.89 - B1 B2 = 82.546 B3 = 162460.0 /T B4 = 0.21828 * T B5 = 126970.0 / T F0 = 1.89 - B1 * ( 372420.0 / (T*T) + 2.0) B6 = B0 *B3 - 2.0 * F0 * ( B2 - B3) B7 = 2.0 * F0 * ( B4- B5) - B0* B5 B8 = 0.43429448 * ALOG(T) F = 775.596 + ( 0.63296 + 0.000162467 *T) *T + 47.3635 *B8 B9 = B0 * P*P/ (2.0 *T*T) HTP = F+ 0.043557 *(F0*P+B9*(B0*(B2-B3 +2.0 *B7*B9)-B6)) T=T1 P = P1 RETURN END C******** ENTHALPY OF A GAS AS FUNCTION OF PRESSURE **************************** FUNCTION HG2(P) X = ALOG10(P) HG2= 1105.9387 + ( 32.756807 +( 4.6198474 + ( 0.20672996 1 +(-0.54116930+ (0.49241362-0.17884885*X)*X)*X)*X)*X)*X RETURN END C********* PRESSURE OF A SATURATED LIQUID AS FUNCTION OF TEMPERATURE *********** FUNCTION PSL(T) TK = (T-32.0)/ 1.8 + 273.16 X = 647.27 - TK Y = X*(3.2437814 + (5.86826E - 03 + 1.1702379E - 08*X*X)*X)/ 1 (TK *(1.0+ 2.1878462E - 03*X) ) PSL = 14.696 *218.167 /(10.0**Y) RETURN END C************ TEMPERATURE OF A SATURATED LIQUID AS FUNCTION OF PRESSURE ******** FUNCTION TSL(P) X = 0.43429448* ALOG(P) X4= X*X*X*X TSL = 101.74419 + (77.052576 + (11.951549+ 2.0562054*X)*X)*X 1 +( 0.42070502+ ( -0.068410987 + 0.0625368*X)*X) *X4 2 - 0.0065948781*X*X*X*X4 RETURN END C******* ENTROPY OF A SATURATED GAS AS FUNCTION OF PRESSURE ******************** FUNCTION SGP(P) IF (P-100) 1, 1, 2 2 IF( P-2000) 3,3,4 3 SGP = 2.2411/ (P** 0.070069) RETURN 4 SGP = 5.66316/ (P**0.194936) RETURN 1 SGP = 1.98473/ (P** 0.0458907) RETURN END C********** ENTROPY OF A SATURATED GAS AS FUNCTION OF TEMPERATURE ************** FUNCTION SGT(T) P = PSL(T) IF(P-100)1,1,2 2 IF (P-2000)3,3,4 3 SGT = 2.2411/(P** 0.070069) RETURN 4 SGT = 5.66316/ ( P** 0.194936) RETURN 1 SGT = 1.98473/ ( P** 0.0458907) RETURN END C************* ENTHALPY OF A SATURATED FLUID AS FUNCTION OF TEMPERATURE ******* FUNCTION HF (T) T4 = T*T*T*T HF = -32.459904+(1.0249349+(-4.1497723E-04+3.0776818E -06*T)*T 1 )*T +(-1.2602857E -08+(3.065808E -11-3.8340495E -14*T)*T)*T4 2 + 1.9906834E -17*T*T*T*T4 RETURN END C*********** ENTHALPY OF SATURATED GAS AS FUNCTION OF TEMPERATURE ************** FUNCTION HG(T) TK= (T-32.0)/1.8+273.16 X= 647.27-TK Y= X*(3.2437814+(5.86826E-03+1.1702379E-08*X*X)*X)/(TK*(1.0+ 1 2.1878462E-03*X)) P = 14.696*218.167/(10.0**Y) X = ALOG10(P) HG = 1105.9387+( 32.756807+(4.6198474+(0.20672996+(-0.5411693 1 +(0.49241362-0.17884885*X)*X)*X)*X)*X)*X RETURN END C*********** ENTHALPY OF SUP. HEATED GAS AS FUNCTION OF PRES. AND ENTROPY******* FUNCTION HPS(P,S) IF(P-100.0) 1,1,2 1 SG= 1.9847258/(P** 0.04589071) GO TO 5 2 IF (P-2000.0)3,4,4 3 SG = 2.2410983/ ( P**0.070069) GO TO 5 4 SG=5.6631597/(P**0.1949355) 5 IF (SG-S) 6,12,12 6 X= 0.43429448* ALOG(P) X4 = X*X*X*X IF(P-10.0) 7, 7, 8 7 S0 = 2.150098+(-0.25438439+(2.17448E -04-9.3986E-04*X)*X)*X A0 = 1223.2933+( -0.57781294+( 0.2303143-1.0434265*X)*X)*X A1 = 820.09617 + (-1.9634176+(2.6069465-0.76847051*X)*X)*X A2 = 895.12074+(-10.468214+( 7.0858389- 10.321004*X)*X)*X A3 = 547.70336+(195.11068+(-313.48831+166.94769*X)*X)*X A4 = 0.0 GO TO 11 8 IF (P-450.0)9,9,10 9 S0 = 2.3335568+(-0.3201163+(0.086914914-0.056646879*X)*X)*X S0 = S0 +(0.018338725- 0.002447789*X)*X4 A0 = 1357.2272+ (73.791077+(-75.924681+34.275666*X)*X)*X 1 -6.0270154*X4 A1= 1144.6178+(33.297322+(-26.451758+8.9579684*X)*X)*X 1 -1.0968016*X4 A2= 993.78383+(521.1334+(-506.58014+220.41684*X)*X)*X 1 -37.982498*X4 A3= 1424.0878+(-1663.6047+(1345.659-489.18341*X)*X)*X 1 +73.075686*X4 A4= 3431.7851+(-7341.2575+(5997.1054-2208.4202*X)*X)*X 1 +297.74553*X4 GO TO 11 10 S0 =1.7066779+(0.54400879+(-0.37780533+0.077093291*X)*X)*X 1 - 0.0054871968*X4 A0= 1400.0 A1= 742.2428+(661.0354+(-321.27928+53.456926*X)*X)*X A2= -3491.438+(4615.4327+(-1470.6537+145.94655*X)*X)*X A3= 34807.748+(-35596.564+(12288.438-1388.0814*X)*X)*X A4= 0.0 11 A5 = S-S0 HPS = A0 +(A1+(A2+ (A3+ A4*A5)*A5)*A5)*A5 RETURN 12 X=0.43429448*ALOG(P) X4 = X*X*X*X A0 = -4.7169141+(-10.049146+(-7.0532835-1.9473822*X)*X)*X 1 + (0.1175487-0.25473452*X)*X4 A1 = 561.46162+(76.93328+(12.117678+2.1291364*X)*X)*X 1 +( 0.12850077+ 0.14437713*X)*X4 HPS = A0 +A1*S RETURN END C******** TEMP. OF SUP. HEATED GAS AS FUNCTION OF ENTHALPY AND PRES.************ FUNCTION THP(H,P) ITER = 0 MAX = 10 EPS = 0.01 DT = 2.0 T = 1.68 *H - 1110.0 H1 = HTP(T,P) 1 ITER = ITER + 1 TDT = T + DT T = T + DT*(H-H1)/(HTP(TDT,P) - H1) H1 = HTP (T,P) IF (ABS(H-H1) - EPS)4, 4, 2 2 IF( ITER - MAX) 1, 1, 3 3 WRITE (LUW,5) CALL MYSTOP(1001) RETURN 4 THP = T RETURN 5 FORMAT ( 10X, ' DID NOT CONVERGE IN FUNCTION THP') END FUNCTION VF(T) T1= T/1000.0 VF= (1.5852859 + (0.2603053 +( -0.7238563 +(10.972684 + (-25.3429 1 69 + 23.071255*T1)*T1)*T1)*T1)*T1)/100.0 RETURN END FUNCTION VG(T,P) ATM =P/14.6959 TDEGK=255.38 +0.5555556*T BE1= (2641.62 * 10.0 **(80870.0/(TDEGK*TDEGK)))/TDEGK BEO= 1.89-BE1 BE2=82.546 BE3= 162460.0/TDEGK BE4= 0.21828 * TDEGK BE5= 126970.0/TDEGK PC1= BEO*ATM PC2= PC1*PC1 TDEGK2= TDEGK*TDEGK BFX= BEO*(1.0+PC1*(BE2-BE3+PC2*(BE4-BE5)/ TDEGK2 )/TDEGK2) VG= 0.0160185 * (4.55504*TDEGK/ATM+ BFX) RETURN END C************* FUNCTION TO CALCULATE PRESSURE GIVEN TEMP. AND VOL. ************* FUNCTION PG(T,V) IMPLICIT REAL(K) KO=5.5 K1=3.62 K2=2.97 B=0.8061588 BPRIME= 0.23217 TSMALL=(T+459.67)/1165.45 RHO= 0.05078/V PSMALL= -(KO+(K1/TSMALL))*(RHO*RHO) + (K2*(-TSMALL+(1/TSMALL)))* 1 ((RHO*RHO)*RHO) + ((RHO* TSMALL)/0.232)/(1-(B*RHO)+ 2 (BPRIME*RHO)*RHO) PG=PSMALL*3208.2 RETURN END