SUBROUTINE TYPE75 (TIME,XIN,OUT,T,DTDT,PAR,INFO,ICNTRL,*) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C* Revised for TRNSYS by N.Blair and R.Schwarz C************************************************************************ C* SUBROUTINE: BYPASS C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate apparatus dew point and bypass C* factor given entering and leaving moist C* air conditions of cooling coil. C*********************************************************************** C* INPUT VARIABLES SAMPLE VALUES C* XIN(1) TEnt Entering air temperature (C) 25.0 C* XIN(2) WEnt Entering air humidity ratio (-) 0.08 C* XIN(3) TLvg Leaving air temperature (C) 30.0 C* XIN(4) WLvg Leaving air humidity ratio (-) 0.1 C* C* OUTPUT VARIABLES C* OUT(1) TAdP Apparatus dewpoint temperature (C) 28.9 C* OUT(2) WAdp Apparatus dewpoint humidity ratio (-) 0.02546 C* OUT(3) BF Bypass factor (-) -0.2844 C* OUT(4) ErrStat Error status indicator, 0 = ok, 1 = error (-) 0.0 C*********************************************************************** C MAJOR RESTRICTIONS: None C C DEVELOPER: Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C Hugh Henderson & Kannan Rengarajan C Florida Solar Energy Center C C DATE: January 1, 1992 C C INCLUDE FILES: prop.inc C SUBROUTINES CALLED: None C FUNCTIONS CALLED: DEWPOINT C ENTHALPY C HUMRATIO C SATPRESS C XITERATE C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C*********************************************************************** C INTERNAL VARIABLES: C hEnt Entering air enthalpy C hLvg Leaving air enthalpy C hAdp Air enthalpy at apparatus dew point C slope Ratio temperature difference to humidity difference C between entering and leaving air states C tAdpEst Estimate of TAdp from slope C error Deviation of dependent variable in iteration C iter Iteration counter C icvg Iteration convergence flag C F1,F2 Previous values of dependent variable in XITERATE C X1,X2 Previous values of independent variable in XITERATE C************************************************************************ DOUBLE PRECISION XIN,OUT INTEGER Patm,CpAir,CpWat,CpLiq,CpVap,DViscAir, & DViscLiq,KAir,KLiq,RhoLiq,Hfg,RAir, & TKelMult,TAbsAdd,PaMult,PAbsAdd REAL Prop(16) DIMENSION XIN(4),OUT(4),INFO(15) INTEGER ErrStat,INFO,IOPT,NI,NP,ND CHARACTER*3 YCHECK(4),OCHECK(4) COMMON /LUNITS/LUR,LUW,IFORM,LUK DATA itmax/20/ DATA YCHECK/'TE1','DM1','TE1','DM1'/ DATA OCHECK/'TE1','DM1','DM1','DM1'/ PARAMETER (Patm = 1) PARAMETER (CpAir = 2) PARAMETER (CpWat = 3) PARAMETER (CpVap = 4) PARAMETER (CpLiq = 5) PARAMETER (DViscAir = 6) PARAMETER (DViscLiq = 7) PARAMETER (KAir = 8) PARAMETER (KLiq = 9) PARAMETER (RhoLiq = 10) PARAMETER (Hfg = 11) PARAMETER (RAir = 12) PARAMETER (TKelMult = 13) PARAMETER (TAbsAdd = 14) PARAMETER (PaMult = 15) PARAMETER (PAbsAdd = 16) prop(patm) = 101325.0 prop(cpair) = 1006.0 prop(cpwat) = 4186.0 prop(cpvap) = 1805.0 prop(cpliq) = 4186.0 prop(dviscair) = .0000182 prop(dviscliq) = .00144 prop(kair) = .026 prop(kliq) = .604 prop(rholiq) = 998.0 prop(hfg) = 2501000.0 prop(rair) = 287.055 prop(tkelmult) = 1.0 prop(tabsadd) = 273.15 prop(pamult) = 1.0 prop(pabsadd) = 0.0 IOPT = -1 NI = 4 !CORRECT NUMBER OF INPUTS NP = 0 !CORRECT NUMBER OF PARAMETERS ND = 0 !CORRECT NUMBER OF DERIVATIVES tent = xin(1) went = xin(2) tlvg = xin(3) wlvg = xin(4) IF (INFO(7).EQ.-1) THEN CALL TYPECK(IOPT,INFO,NI,NP,ND) C CHECKS TO SEE IF USER'S INFO MATCHES CORRECT NUMBERS CALL RCHECK(INFO,YCHECK,OCHECK) C CHECKS TO SEE IF INPUT AND OUTPUT UNITS MATCH ENDIF C1*** Iterate to determine apparatus dewpoint at which the ADP C1*** equals the temperature calculated by extending the line between C1*** entering and leaving conditions to the saturation curve C1*** Calculate "slope" of temperature vs. humidity ratio between C1*** entering and leaving states slope = (TEnt-TLvg)/(WEnt-WLvg) C1*** Initialize iteration parameters TAdp = DEWPOINT(Prop,WLvg) DO 100 iter=1,itmax C1*** Calculate apparatus dewpoint and compare with predicted value C1*** using entering conditions and slope WAdp = HUMRATIO(Prop(Patm),SATPRESS(Prop,TAdp)) TAdpEst = TEnt - slope*(WEnt-WAdp) error = TAdp-TAdpEst TAdp = XITERATE(TAdp,error,X1,F1,X2,F2,iter,icvg) C1*** If converged, exit loop IF (icvg .EQ. 1) GO TO 110 100 CONTINUE C1*** Apparatus dewpoint has not converged after maximum iterations. C1*** Print error message, set return error flag, and RETURN WRITE(LUW,1001) itmax 1001 FORMAT(/1X,'*** ERROR IN SUBOUTINE BYPASS ***'/ & 1X,' Apparatus dewpoint has not ' & 'converged after ',I2,' iterations'/) ErrStat = 1 110 CONTINUE C1*** Calculate bypass factor from enthalpies hLvg = ENTHALPY(Prop,TLvg,WLvg) hEnt = ENTHALPY(Prop,TEnt,WEnt) hAdp = ENTHALPY(Prop,TAdp,WAdp) BF = (hLvg-hAdp)/(hEnt-hAdp) out(1) = tadp out(2) = wadp out(3) = bf out(4) = errstat 999 RETURN 1 END REAL FUNCTION DEWPOINT (Prop,W) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* FUNCTION: DEWPOINT C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate the dewpoint temperature given C* humidity ratio C*********************************************************************** C* INPUT VARIABLES C* W Humidity ratio (-) C* C* OUTPUT VARIABLES C* DewPoint Dew point temperature of air (C) C* C* PROPERTIES C* Patm Atmospheric pressure (Pa) C*********************************************************************** C MAJOR RESTRICTIONS: None C C DEVELOPER: Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: None C SUBROUTINES CALLED: None C FUNCTIONS CALLED: None C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C*********************************************************************** C INTERNAL VARIABLES: C pw Partial water vapor pressure (Pa) C small Small number C*********************************************************************** INTEGER Patm,CpAir,CpWat,CpLiq,CpVap,DViscAir, & DViscLiq,KAir,KLiq,RhoLiq,Hfg,RAir, & TKelMult,TAbsAdd,PaMult,PAbsAdd REAL Prop(16) PARAMETER (Patm = 1) PARAMETER (CpAir = 2) PARAMETER (CpWat = 3) PARAMETER (CpVap = 4) PARAMETER (CpLiq = 5) PARAMETER (DViscAir = 6) PARAMETER (DViscLiq = 7) PARAMETER (KAir = 8) PARAMETER (KLiq = 9) PARAMETER (RhoLiq = 10) PARAMETER (Hfg = 11) PARAMETER (RAir = 12) PARAMETER (TKelMult = 13) PARAMETER (TAbsAdd = 14) PARAMETER (PaMult = 15) PARAMETER (PAbsAdd = 16) DATA small/1.E-9/ C1*** Test for "dry" air IF (W .LT. small) THEN DewPoint = -999 ELSE C1*** Calculate the partial water vapor pressure as a function of C1*** humidity ratio. pw= Prop(Patm)*W/(.62198+W) C1*** Calculate dewpoint as saturation temperature at water vapor C1*** partial pressure DewPoint = SATTEMP(Prop,pw) ENDIF 999 RETURN END REAL FUNCTION ENTHALPY (Prop,TDB,W) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* FUNCTION: ENTHALPY C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate the enthalpy of moist air. C*********************************************************************** C* INPUT VARIABLES: C* TDB Dry bulb temperature (C) C* W Humidity ratio (-) C* C* OUTPUT VARIABLES: C* Enthalpy Enthalpy of moist air (J/kg) C* C* PROPERTIES: C* CpAir Specific heat of air (J/kg C) C* CpVap Specific heat of water vapor (J/kg C) C* Hfg Reference heat of vaporization of water (J/kg) C*********************************************************************** C MAJOR RESTRICTIONS Uses perfect gas relationships C Fit for enthalpy of saturated water vapor C C DEVELOPER: Shauna Gabel C Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: PROP.INC C SUBROUTINES CALLED: None C FUNCTIONS CALLED: None C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C*********************************************************************** INTEGER Patm,CpAir,CpWat,CpLiq,CpVap,DViscAir, & DViscLiq,KAir,KLiq,RhoLiq,Hfg,RAir, & TKelMult,TAbsAdd,PaMult,PAbsAdd REAL Prop(16) PARAMETER (Patm = 1) PARAMETER (CpAir = 2) PARAMETER (CpWat = 3) PARAMETER (CpVap = 4) PARAMETER (CpLiq = 5) PARAMETER (DViscAir = 6) PARAMETER (DViscLiq = 7) PARAMETER (KAir = 8) PARAMETER (KLiq = 9) PARAMETER (RhoLiq = 10) PARAMETER (Hfg = 11) PARAMETER (RAir = 12) PARAMETER (TKelMult = 13) PARAMETER (TAbsAdd = 14) PARAMETER (PaMult = 15) PARAMETER (PAbsAdd = 16) C1*** Calculate the enthalpy as a function of dry bulb temperature and C1*** humidity ratio. hDryAir = Prop(CpAir)*TDB hSatVap = Prop(Hfg) + Prop(CpVap)*TDB Enthalpy = hDryAir + W*hSatVap RETURN END REAL FUNCTION HUMRATIO (Patm,Pw) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* FUNCTION: HUMRATIO C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate the humidity ratio from water C* vapor pressure and atmospheric pressure C*********************************************************************** C* INPUT VARIABLES C* Patm Atmospheric pressure (Pa) C* Pw Partial water vapor pressure (Pa) C* C* OUTPUT VARIABLES C* HumRatio Humidity ratio (-) C*********************************************************************** C MAJOR RESRICTIONS: None C C DEVELOPER: Shauna Gabel C Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: None C SUBROUTINES CALLED: None C FUNCTIONS CALLED: None C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C************************************************************************ C1*** Calculate the humidity ratio. HumRatio = 0.62198*Pw/(Patm-Pw) RETURN END REAL FUNCTION SATPRESS (Prop,T) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* SUBROUTINE: SATPRESS C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate saturation pressure of water C* vapor as a function of temperature C*********************************************************************** C* INPUT VARIABLES C* T Temperature (C) C* C* OUTPUT VARIABLES C* SatPress Saturation pressure (Pa) C* C* PROPERTIES C* TKelMult Multiplying factor to convert user T to Kelvin C* TAbsAdd Additive factor to convert user T to absolute T C* tKel = Prop(TKelMult) * (T + Prop(TAbsAdd)) C* PaMult Multiplying factor to convert user P to Pascals C* PAbsAdd Additive factor to convert user P to absolute P C* Pa = Prop(PaMult) * (P + Prop(PAbsAdd)) C*********************************************************************** C MAJOR RESTRICTIONS: 173.16 K <= Temp <= 473.15 K C C DEVELOPER: Shauna Gabel C Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: prop.inc C SUBROUTINES CALLED: None C FUNCTIONS CALLED: None C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C C Hyland, R.W., and A. Wexler. 1983. C Formulations for the thermodynamic C properties of the saturated phases of H2O C from 173.15 K to 473.15 K. ASHRAE C Transactions, Vol. 89, No. 2A, pp. 500-519 C*********************************************************************** C INTERNAL VARIABLES: C tKel Temperature in Kelvin (K) C pascals Saturation pressure (Pa) C*********************************************************************** INTEGER Patm,CpAir,CpWat,CpLiq,CpVap,DViscAir, & DViscLiq,KAir,KLiq,RhoLiq,Hfg,RAir, & TKelMult,TAbsAdd,PaMult,PAbsAdd REAL Prop(16) PARAMETER (Patm = 1) PARAMETER (CpAir = 2) PARAMETER (CpWat = 3) PARAMETER (CpVap = 4) PARAMETER (CpLiq = 5) PARAMETER (DViscAir = 6) PARAMETER (DViscLiq = 7) PARAMETER (KAir = 8) PARAMETER (KLiq = 9) PARAMETER (RhoLiq = 10) PARAMETER (Hfg = 11) PARAMETER (RAir = 12) PARAMETER (TKelMult = 13) PARAMETER (TAbsAdd = 14) PARAMETER (PaMult = 15) PARAMETER (PAbsAdd = 16) DATA C1/-5674.5359/,C2/6.3925247/,C3/-0.9677843E-2/ DATA C4/0.62215701E-6/,C5/0.20747825E-8/,C6/-0.9484024E-12/ DATA C7/4.1635019/,C8/-5800.2206/,C9/1.3914993/,C10/-0.048640239/ DATA C11/0.41764768E-4/,C12/-0.14452093E-7/,C13/6.5459673/ C1*** Convert temperature from user units to Kelvin. tKel = Prop(TKelMult)*(T+Prop(TAbsAdd)) C1*** If below freezing, calculate saturation pressure over ice. IF (tKel .LT. 273.15) THEN pascals = EXP(C1/tKel+C2+C3*tKel+C4*tKel**2+C5*tKel**3+C6* & tKel**4+C7*ALOG(tKel)) C1*** If above freezing, calculate saturation pressure over liquid water. ELSE IF (tKel .GE. 273.15) THEN pascals = EXP(C8/tKel+C9+C10*tKel+C11*tKel**2+C12*tKel**3+C13 & *ALOG(tKel)) ENDIF C1*** Convert pressure from Pascals to user units SatPress = pascals/Prop(PaMult) - Prop(PAbsAdd) RETURN END REAL FUNCTION XITERATE (X0,F0,X1,F1,X2,F2,ICount,ICvg) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* C* SUBROUTINE: XITERATE C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Iterately solves for the value of X which C* satisfies F(X)=0. Given Xi,F(Xi) pairs, C* the subroutine tests for convergence and C* provides a new guess for the value of the C* independent variable X. C*********************************************************************** C* INPUT VARIABLES C* F0 Current value of the function F(X) C* X0 Current value of X C* F1,F2 Two previous values of F(Xi) C* X1,X2 Two previous values of X C* C* NOTE: F1,X1,F2,X2 MUST BE STORED AND SAVED IN CALLING C* ROUTINE. THEY NEED NO INITIALIZATION C* C* ICount Number of iterations C* C* OUTPUT VARIABLES C* XIterate New estimate of X for F(X)=0 C* ICvg Convergence flag ICvg = 0: Not converged C* ICvg = 1: Converged C*********************************************************************** C DEVELOPER: Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: None C SUBROUTINES CALLED: None C FUNCTIONS CALLED: None C C REVISION HISTORY: None C C REFERENCE: None C*********************************************************************** C INTERNAL VARIABLES C small Small number used in place of zero C mode Number of points used in fit C mode = 1: Use XPerburb to get new X C mode = 2: Linear equation to get new X C mode > 2: Quadratic equation to get new X C coef(i) Coefficients for quadratic fit C F(X) = coef(1) + coef(2)*X + coef(3)*X*X C check Term under radical in quadratic solution C FiQ,XiQ Double precision values of Fi,Xi C slope Slope for linear fit C tolRel Relative error tolerance C xPerturb Perturbation applied to X to initialize iteration C*********************************************************************** DOUBLE PRECISION coef(3),check,F0Q,F1Q,F2Q,X0Q,X1Q,X2Q DATA tolRel/1.E-5/,xPerturb/0.1/,small/1.E-9/ C1*** Check for convergence by comparing change in X IF ((ABS(X0-X1) .LT. tolRel*MAX(ABS(X0),small) .AND. & ICount .NE. 1) .OR. F0 .EQ. 0.) THEN XIterate = X0 ICvg=1 RETURN ENDIF C1*** Not converged. C2*** If after the second iteration there are enough previous points to C2 fit a quadratic for the new X. If the quadratic fit is not C2 applicable, mode will be set to 1 or 2 and a new X will be C2 determined by incrementing X from xPerturb or from a linear fit. ICvg=0 mode=ICount 10 IF (mode .EQ. 1) THEN C1*** New guess is specified by xPerturb IF (ABS(X0) .GT. small) THEN XIterate = X0*(1.+xPerturb) ELSE XIterate = xPerturb ENDIF ELSEIF (mode .EQ. 2) THEN C1*** New guess calculated from LINEAR FIT of most recent two points SLOPE=(F1-F0)/(X1-X0) IF(slope.EQ.0) THEN mode=1 GO TO 10 ENDIF XIterate=X0-F0/SLOPE ELSE C1*** New guess calculated from QUADRATIC FIT C1*** If two Xi are equal, set mode for linear fit and return to top IF (X0 .EQ. X1) THEN X1=X2 F1=F2 mode=2 GO TO 10 ELSEIF (X0 .EQ. X2) THEN mode=2 GO TO 10 ENDIF C1*** Determine quadratic coefficients from the three data points C1*** using double precision. F2Q=F2 F1Q=F1 F0Q=F0 X2Q=X2 X1Q=X1 X0Q=X0 coef(3)=((F2Q-F0Q)/(X2Q-X0Q)-(F1Q-F0Q)/(X1Q-X0Q))/(X2Q-X1Q) coef(2)=(F1Q-F0Q)/(X1Q-X0Q)-(X1Q+X0Q)*coef(3) coef(1)=F0-(coef(2)+coef(3)*X0Q)*X0Q C1*** If points are colinear, set mode for linear fit and return to top IF (ABS(coef(3)) .LT. 1.D-10) THEN mode=2 GO TO 10 ENDIF C1*** Check for precision. If the coefficients do not accurately C1*** predict the given data points due to round-off errors, set C1*** mode for a linear fit and return to top. IF (ABS((coef(1)+(coef(2)+coef(3)*X1Q)*X1Q-F1Q)/F1Q) .GT. & 1.D-4) THEN mode=2 GO TO 10 ENDIF C1*** Check for imaginary roots. If no real roots, set mode to C1*** estimate new X by simply incrementing by xPerturb check=coef(2)**2-4*coef(1)*coef(3) IF (check .LT. 0) THEN C1*** Imaginary roots -- go back to linear fit mode=2 GO TO 10 ELSEIF (check .GT. 0) THEN C1*** Real unequal roots -- determine root nearest to most recent guess XIterate=(-coef(2)+SQRT(check))/coef(3)/2 xOther=-XIterate-coef(2)/coef(3) IF (ABS(XIterate-X0) .GT. ABS(xOther-X0)) XIterate=xOther ELSE C1*** Real Equal Roots -- one solution XIterate=-coef(2)/coef(3)/2 ENDIF ENDIF C1*** Set previous variable values for the next iteration IF (mode .LT. 3) THEN C1*** No valid previous points to eliminate. X2=X1 F2=F1 X1=X0 F1=F0 ELSE C1*** Eliminate one previous point based on sign and magnitude of F(X) C2*** Keep the current point and eliminate one of the previous ones. IF (F1*F0 .GT. 0 .AND. F2*F0 .GT. 0) THEN C2*** All previous points of same sign. Eliminate one with biggest F(X) IF (ABS(F2) .GT. ABS(F1)) THEN X2=X1 F2=F1 ENDIF ELSE C1*** Points of different sign. C1*** Eliminate the previous one with the same sign as current F(X). IF (F2*F0 .GT. 0) THEN X2=X1 F2=F1 ENDIF ENDIF X1=X0 F1=F0 ENDIF RETURN END C*********************************************************************** C* FILE: PROP.INC C* C* This file assigns a numbers to air and water property names to be C* used in the "Prop" array. C*********************************************************************** C DEVELOPER: Inger Andresen C Michael J. Brandemuehl, PhD, PE C C DATE: July 1, 1991 C C FILES REQUIRED: None C*********************************************************************** C INTERNAL VARIABLES: C Patm Atmospheric pressure (Pa) C CpAir Specific heat of dry air (J/kg C) C CpLiq Specific heat of liquid water (J/kg C) C CpVap Specific heat of saturated water vapor (J/kg C) C DViscAir Air dynamic viscosity (kg/m s) C DViscLiq Liquid dynamic viscosity (kg/m s) C KAir Air thermal conductivity (W/m C) C KLiq Liquid thermal conductivity (W/m C) C RhoLiq Liquid density (kg/m3) C Hfg Latent heat of vaporization of water (J/kg) C RAir Gas constant for air (J/kg C) C TKelMult Multiplying factor to convert user T to Kelvin C TAbsAdd Additive factor to convert user P to Kelvin C tKel = Prop(TKelMult)*T + Prop(TKelAdd) C PaMult Multiplying factor to convert user P to Pascals C PAbsAdd Additive factor to convert user P to Pascals C Pa = Prop(PaMult)*P + Prop(PaAdd) C*********************************************************************** REAL FUNCTION SATTEMP (Prop,P) C*********************************************************************** C* Copyright ASHRAE. Toolkit for HVAC System Energy Calculations C*********************************************************************** C* FUNCTION: SATTEMP C* C* LANGUAGE: FORTRAN 77 C* C* PURPOSE: Calculate the saturation (boiling) C* temperature of water given pressure C*********************************************************************** C* INPUT VARIABLES C* P Pressure (Pa) C* C* OUTPUT VARIABLES C* SatTemp Saturation temperature of water vapor (C) C*********************************************************************** C MAJOR RESTRICTIONS: None C C DEVELOPER: Shauna Gabel C Michael J. Brandemuehl, PhD, PE C University of Colorado at Boulder C C DATE: January 1, 1992 C C INCLUDE FILES: prop.inc C SUBROUTINES CALLED: None C FUNCTIONS CALLED: SATPRESS C XITERATE C C REVISION HISTORY: None C C REFERENCE: 1989 ASHRAE Handbook - Fundamentals C*********************************************************************** C INTERNAL VARIABLES: C tSat Water temperature guess (C) C pSat Pressure corresponding to temp. guess (Pa) C error Deviation of dependent variable in iteration C iter Iteration counter C icvg Iteration convergence flag C F1,F2 Previous values of dependent variable in XITERATE C X1,X2 Previous values of independent variable in XITERATE C*********************************************************************** DATA itmax/50/ C1*** Use an iterative process to determine the saturation temperature C1*** at a given pressure using a correlation of saturated water vapor C1*** pressure as a function of temperature C1*** Initial guess of boiling temperature tSat = 100. C1*** Iterate to find the saturation temperature C1*** of water given the total pressure C2*** Set iteration loop parameters DO 100 iter = 1,itmax C1*** Calculate saturation pressure for estimated boiling temperature pSat = SATPRESS(Prop,tSat) C1*** Compare with specified pressure and update estimate of temperature error = P - pSat tSat = XITERATE (tSat,error,X1,F1,X2,F2,iter,icvg) C2*** If converged leave loop iteration IF (icvg .EQ. 1) GO TO 110 C2*** Water temperature not converged, repeat calculations with new C2*** estimate of water temperature 100 CONTINUE C1*** Saturation temperature has not converged after maximum specified C1*** iterations. Print error message, set return error flag, and RETURN WRITE(LUW,1001) itmax 1001 FORMAT(/1X,'*** ERROR IN FUNCTION SatTemp ***'/ & 1X,' Saturation temperature has not ' & 'converged after ',I2,' iterations'/) 110 SatTemp = tSat RETURN END