C PROPATH VER. 10.1 (HELIUM 4) C C==================================================================== C May 2, 2001: modified by Ryo Akasaka for ver.12.1 C The following 18 functions were added but not implemented. C AKPD(8A), AKPDD(8B), AKTD(8C), AKTDD(8D), CVPD(7A), CVTD(7B), C EPSPD(2A), EPSPDD(2B), EPSTD(2C), EPSTDD(2D), GAMPD(9A), C GAMTD(9B), TPH2(6H), TPS2(6S), WPD(8E), WPDD(8F), WTD(8G), C WTDD(8H) C C------------------------------------------------- F8A = AKPD REAL FUNCTION AKPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='AKPD' IF (MESS.NE.0) CALL S18HE(FUN) AKPD=-1.0E+30 RETURN END C------------------------------------------------- F8B = AKPDD REAL FUNCTION AKPDD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='AKPDD' IF (MESS.NE.0) CALL S18HE(FUN) AKPDD=-1.0E+30 RETURN END C------------------------------------------------- F8C = AKTD REAL FUNCTION AKTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='AKTD' IF (MESS.NE.0) CALL S18HE(FUN) AKTD=-1.0E+30 RETURN END C------------------------------------------------- F8D = AKTDD REAL FUNCTION AKTDD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='AKTDD' IF (MESS.NE.0) CALL S18HE(FUN) AKTDD=-1.0E+30 RETURN END C------------------------------------------------- F7A = CVPD REAL FUNCTION CVPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='CVPD' IF (MESS.NE.0) CALL S18HE(FUN) CVPD=-1.0E+30 RETURN END C------------------------------------------------- F7B = CVTD REAL FUNCTION CVTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='CVTD' IF (MESS.NE.0) CALL S18HE(FUN) CVTD=-1.0E+30 RETURN END C------------------------------------------------- F2A = EPSPD REAL FUNCTION EPSPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='EPSPD' IF (MESS.NE.0) CALL S18HE(FUN) EPSPD=-1.0E+30 RETURN END C------------------------------------------------- F2B = EPSPDD REAL FUNCTION EPSPDD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='EPSPDD' IF (MESS.NE.0) CALL S18HE(FUN) EPSPDD=-1.0E+30 RETURN END C------------------------------------------------- F2C = EPSTD REAL FUNCTION EPSTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='EPSTD' IF (MESS.NE.0) CALL S18HE(FUN) EPSTD=-1.0E+30 RETURN END C------------------------------------------------- F2D = EPSTDD REAL FUNCTION EPSTDD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='EPSTDD' IF (MESS.NE.0) CALL S18HE(FUN) EPSTDD=-1.0E+30 RETURN END C------------------------------------------------- F9A = GAMPD REAL FUNCTION GAMPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='GAMPD' IF (MESS.NE.0) CALL S18HE(FUN) GAMPD=-1.0E+30 RETURN END C------------------------------------------------- F9B = GAMTD REAL FUNCTION GAMTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='GAMTD' IF (MESS.NE.0) CALL S18HE(FUN) GAMTD=-1.0E+30 RETURN END C------------------------------------------------- F6H = TPH2 REAL FUNCTION TPH2(P,H) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='TPH2' IF (MESS.NE.0) CALL S18HE(FUN) TPH2=-1.0E+30 RETURN END C------------------------------------------------- F6S = TPS2 REAL FUNCTION TPS2(P,S) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='TPS2' IF (MESS.NE.0) CALL S18HE(FUN) TPS2=-1.0E+30 RETURN END C------------------------------------------------- F8E = WPD REAL FUNCTION WPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='WPD' IF (MESS.NE.0) CALL S18HE(FUN) WPD=-1.0E+30 RETURN END C------------------------------------------------- F8F = WPDD REAL FUNCTION WPDD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='WPDD' IF (MESS.NE.0) CALL S18HE(FUN) WPDD=-1.0E+30 RETURN END C------------------------------------------------- F8G = WTD REAL FUNCTION WTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='WTD' IF (MESS.NE.0) CALL S18HE(FUN) WTD=-1.0E+30 RETURN END C------------------------------------------------- F8H = WTDD REAL FUNCTION WTDD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='WTDD' IF (MESS.NE.0) CALL S18HE(FUN) WTDD=-1.0E+30 RETURN END C C==================================================================== C------------------------------------------------- F1 = AIPPT REAL FUNCTION AIPPT(P,T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='AIPPT' A=P+T IF (MESS.NE.0) CALL S18HE(FUN) AIPPT=-1.0E+30 RETURN END C------------------------------------------------- F94 = AJTPT [K/Pa] REAL FUNCTION AJTPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F94HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AJTPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F94HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') AJTPT=FF RETURN END C------------------------------------------------- F82 = AKPT REAL FUNCTION AKPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F82HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AKPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F82HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') AKPT=FF RETURN END C------------------------------------------------- F2 = ALAPP REAL FUNCTION ALAPP(P) CHARACTER FUN*6 DOUBLE PRECISION F2HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALAPP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F2HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') ALAPP=FF RETURN END C------------------------------------------------- F3 = ALAPT REAL FUNCTION ALAPT(T) CHARACTER FUN*6 DOUBLE PRECISION F3HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALAPT' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F3HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') ALAPT=FF RETURN END C------------------------------------------------- F4 = ALHP REAL FUNCTION ALHP(P) CHARACTER FUN*6 DOUBLE PRECISION F4HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALHP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F4HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') ALHP=FF RETURN END C------------------------------------------------- F5 = ALHT REAL FUNCTION ALHT(T) CHARACTER FUN*6 DOUBLE PRECISION F5HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALHT' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F5HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') ALHT=FF RETURN END C------------------------------------------------- F6 = ALMPD REAL FUNCTION ALMPD(P) CHARACTER FUN*6 DOUBLE PRECISION F6HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALMPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F6HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') ALMPD=FF RETURN END C------------------------------------------------- F7 = ALMPDD REAL FUNCTION ALMPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F7HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALMPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F7HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') ALMPDD=FF RETURN END C------------------------------------------------- F8 = ALMPT REAL FUNCTION ALMPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F8HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALMPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F8HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') ALMPT=FF RETURN END C------------------------------------------------- F9 = ALMTD REAL FUNCTION ALMTD(T) CHARACTER FUN*6 DOUBLE PRECISION F9HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALMTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F9HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') ALMTD=FF RETURN END C------------------------------------------------- F10 = ALMTDD REAL FUNCTION ALMTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F10HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='ALMTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F10HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') ALMTDD=FF RETURN END C------------------------------------------------- F11 = AMUPD REAL FUNCTION AMUPD(P) CHARACTER FUN*6 DOUBLE PRECISION F11HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AMUPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F11HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') AMUPD=FF RETURN END C------------------------------------------------- F12 = AMUPDD REAL FUNCTION AMUPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F12HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AMUPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F12HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') AMUPDD=FF RETURN END C------------------------------------------------- F13 = AMUPT REAL FUNCTION AMUPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F13HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AMUPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F13HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') AMUPT=FF RETURN END C------------------------------------------------- F14 = AMUTD REAL FUNCTION AMUTD(T) CHARACTER FUN*6 DOUBLE PRECISION F14HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AMUTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F14HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') AMUTD=FF RETURN END C------------------------------------------------- F15 = AMUTDD REAL FUNCTION AMUTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F15HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='AMUTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F15HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') AMUTDD=FF RETURN END C------------------------------------------------- F92 = BPPT [1/K] REAL FUNCTION BPPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F92HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='BPPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F92HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') BPPT=FF RETURN END C------------------------------------------------- F90 = BSPT [1/Pa] REAL FUNCTION BSPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F90HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='BSPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F90HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') BSPT=FF RETURN END C------------------------------------------------- F91 = BTPT [1/Pa] REAL FUNCTION BTPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F91HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='BTPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F91HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') BTPT=FF RETURN END C------------------------------------------------- F93 = BVPT [1/K] REAL FUNCTION BVPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F93HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='BVPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F93HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') BVPT=FF RETURN END C------------------------------------------------- F16 = CPPD REAL FUNCTION CPPD(P) CHARACTER FUN*6 DOUBLE PRECISION F16HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CPPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F16HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') CPPD=FF RETURN END C------------------------------------------------- F17 = CPPDD REAL FUNCTION CPPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F17HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CPPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F17HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') CPPDD=FF RETURN END C------------------------------------------------- F18 = CPPT REAL FUNCTION CPPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F18HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CPPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F18HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') CPPT=FF RETURN END C------------------------------------------------- F19 = CPTD REAL FUNCTION CPTD(T) CHARACTER FUN*6 DOUBLE PRECISION F19HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CPTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F19HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') CPTD=FF RETURN END C------------------------------------------------- F20 = CPTDD REAL FUNCTION CPTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F20HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CPTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F20HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') CPTDD=FF RETURN END C------------------------------------------------- F21 = CRP REAL FUNCTION CRP(A) CHARACTER A*1 DOUBLE PRECISION F21HE,PBAR,T0K COMMON/UNIT/KPA,MESS CALL S15HE(KPA,PBAR,T0K) FF = F21HE(A) IF (MESS.EQ.0) GO TO 50 IF(FF.EQ.-1.0E+20) THEN WRITE(6,6010) A 6010 FORMAT(1H ,5X,'***** OUT OF RANGE AT CRP FOR HELIUM WHEN ', & 'A = ',A1,' *****') FF=-1.0E+20 END IF 50 IF(A.EQ.'T') THEN IF(FF.EQ.-1.0E+20) T0K=0.0 FF=FF-REAL(T0K) ELSE IF(A.EQ.'P') THEN IF(FF.EQ.-1.0E+20) PBAR=1.0 FF=FF/REAL(PBAR) END IF CRP=FF RETURN END C------------------------------------------------- F76 = CVPDD REAL FUNCTION CVPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F76HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CVPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F76HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') CVPDD=FF RETURN END C------------------------------------------------- F77 = CVPT REAL FUNCTION CVPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F77HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CVPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F77HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') CVPT=FF RETURN END C------------------------------------------------- F78 = CVTDD REAL FUNCTION CVTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F78HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='CVTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F78HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') CVTDD=FF RETURN END C------------------------------------------------- F22 = EPSPT REAL FUNCTION EPSPT(P,T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='EPSPT' A=P+T IF (MESS.NE.0) CALL S18HE(FUN) EPSPT=-1.0E+30 RETURN END C------------------------------------------------- F89 = FC C************************************************ C FUNCTION FOR FUNDDAMENTAL CONSTANTS C PROPATH VER.7.1, MAY 8, 1990 C USAGE: B=FC(A) C A, B : CHARACTER TYPE VALIABLES C B='4.0026' WHEN A='M' C B='2077.2' WHEN A='R' C************************************************ REAL FUNCTION FC(A) CHARACTER A*1,MSG*120 COMMON/UNIT/KPA,MESS IF (A.EQ.'M') THEN FC=4.0026 ELSE IF (A.EQ.'R') THEN FC=2077.2 ELSE FC=-1.E+20 IF (MESS.NE.0) THEN MSG='**** OUT OF RANGE AT FC FOR HELIUM4 WHEN A=''' & //A//''' ****' WRITE(6,'(1H ,A)') MSG END IF END IF RETURN END C------------------------------------------------- F95= GAMPT REAL FUNCTION GAMPT(P,T) CHARACTER FUN*6 REAL P,T,FF INTEGER KPA DOUBLE PRECISION F95HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMPT'/ CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F95HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') GAMPT=FF RETURN END C------------------------------------------------- F96 = GAMPDD REAL FUNCTION GAMPDD(P) CHARACTER FUN*6 REAL P,FF INTEGER KPA DOUBLE PRECISION F96HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMPDD'/ CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F96HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') GAMPDD=FF RETURN END C------------------------------------------------- F97 = GAMTDD REAL FUNCTION GAMTDD(T) CHARACTER FUN*6 REAL T,FF INTEGER KPA DOUBLE PRECISION F97HE,PBAR,DBT,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMTDD'/ CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F97HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') GAMTDD=FF RETURN END C------------------------------------------------- F23 = HPD REAL FUNCTION HPD(P) CHARACTER FUN*6 DOUBLE PRECISION F23HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F23HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') HPD=FF RETURN END C------------------------------------------------- F24 = HPDD REAL FUNCTION HPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F24HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F24HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') HPDD=FF RETURN END C------------------------------------------------- F71 = HPS REAL FUNCTION HPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F71HE,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HPS' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBS=DBLE(S) FF = F71HE(DBP,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,S,'P','S') HPS=FF RETURN END C------------------------------------------------- F25 = HPT REAL FUNCTION HPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F25HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F25HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') HPT=FF RETURN END C------------------------------------------------- F26 = HPX REAL FUNCTION HPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F26HE,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HPX' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBX=DBLE(X) FF = F26HE(DBP,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,X,'P','X') HPX=FF RETURN END C------------------------------------------------- F27 = HTD REAL FUNCTION HTD(T) CHARACTER FUN*6 DOUBLE PRECISION F27HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F27HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') HTD=FF RETURN END C------------------------------------------------- F28 = HTDD REAL FUNCTION HTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F28HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F28HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') HTDD=FF RETURN END C------------------------------------------------- F29 = HTX REAL FUNCTION HTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F29HE,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='HTX' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBX=DBLE(X) FF = F29HE(DBT,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,X,'T','X') HTX=FF RETURN END C------------------------------------------------- F84 = IDENTF C************************************************ C FUNCTION FOR IDENTIFICATION OF SUBSTANCE C PROPATH VER.12.1, MAY 2, 2001 C USAGE: B=IDENTF(A) C A, B : CHARACTER TYPE VALIABLES C B='HELIUM4' WHEN A='S' C B='HE4' WHEN A='C' C B='12.1' WHEN A='V' C************************************************ CHARACTER*20 FUNCTION IDENTF(A) CHARACTER A*1,MSG*120 COMMON/UNIT/KPA,MESS IF (A.EQ.'S') THEN IDENTF='HELIUM 4(IPTS 1968)' ELSE IF (A.EQ.'C') THEN IDENTF='HE4' ELSE IF (A.EQ.'V') THEN IDENTF='12.1' ELSE IDENTF='????????????????????' IF (MESS.NE.0) THEN MSG='**** OUT OF RANGE AT IDENTF FOR HELIUM4 WHEN A=''' & //A//''' ****' WRITE(6,'(1H ,A)') MSG END IF END IF RETURN END C------------------------------------------------- F66 = PLDT REAL FUNCTION PLDT(T) CHARACTER FUN*6 DOUBLE PRECISION F66HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PLDT' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F66HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) PBAR=1.0 PLDT=FF/REAL(PBAR) RETURN END C------------------------------------------------- F68 = PMLT REAL FUNCTION PMLT(T) CHARACTER FUN*6 DOUBLE PRECISION F68HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PMLT' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F68HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) PBAR=1.0 PMLT=FF/REAL(PBAR) RETURN END C------------------------------------------------- F85 = PRPD REAL FUNCTION PRPD(P) CHARACTER FUN*6 DOUBLE PRECISION F85HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PRPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F85HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') PRPD=FF RETURN END C------------------------------------------------- F86 = PRPDD REAL FUNCTION PRPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F86HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PRPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F86HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') PRPDD=FF RETURN END C------------------------------------------------- F81 = PRPT REAL FUNCTION PRPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F81HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PRPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F81HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') PRPT=FF RETURN END C------------------------------------------------- F87 = PRTD REAL FUNCTION PRTD(T) CHARACTER FUN*6 DOUBLE PRECISION F87HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PRTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F87HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') PRTD=FF RETURN END C------------------------------------------------- F88 = PRTDD REAL FUNCTION PRTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F88HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PRTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F88HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') PRTDD=FF RETURN END C------------------------------------------------- F99 = PSBT REAL FUNCTION PSBT(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='PSBT' A=T IF (MESS.NE.0) CALL S18HE(FUN) PSBT=-1.0E+30 RETURN END C------------------------------------------------- F30 = PST REAL FUNCTION PST(T) CHARACTER FUN*6 DOUBLE PRECISION F30HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='PST' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F30HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) PBAR=1.0 PST=FF/REAL(PBAR) RETURN END C------------------------------------------------- F72 = PSTD REAL FUNCTION PSTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='PSTD' A=T IF (MESS.NE.0) CALL S18HE(FUN) PSTD=-1.0E+30 RETURN END C------------------------------------------------- F73 = PSTDD REAL FUNCTION PSTDD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='PSTDD' A=T IF (MESS.NE.0) CALL S18HE(FUN) PSTDD=-1.0E+30 RETURN END C------------------------------------------------- F31 = SIGP REAL FUNCTION SIGP(P) CHARACTER FUN*6 DOUBLE PRECISION F31HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SIGP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F31HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') SIGP=FF RETURN END C------------------------------------------------- F32 = SIGT REAL FUNCTION SIGT(T) CHARACTER FUN*6 DOUBLE PRECISION F32HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SIGT' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F32HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') SIGT=FF RETURN END C------------------------------------------------- F33 = SPD REAL FUNCTION SPD(P) CHARACTER FUN*6 DOUBLE PRECISION F33HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F33HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') SPD=FF RETURN END C------------------------------------------------- F34 = SPDD REAL FUNCTION SPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F34HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F34HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') SPDD=FF RETURN END C------------------------------------------------- F35 = SPT REAL FUNCTION SPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F35HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F35HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') SPT=FF RETURN END C------------------------------------------------- F36 = SPX REAL FUNCTION SPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F36HE,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='SPX' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBX=DBLE(X) FF = F36HE(DBP,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,X,'P','X') SPX=FF RETURN END C------------------------------------------------- F37 = STD REAL FUNCTION STD(T) CHARACTER FUN*6 DOUBLE PRECISION F37HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='STD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F37HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') STD=FF RETURN END C------------------------------------------------- F38 = STDD REAL FUNCTION STDD(T) CHARACTER FUN*6 DOUBLE PRECISION F38HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='STDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F38HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') STDD=FF RETURN END C------------------------------------------------- F39 = STX REAL FUNCTION STX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F39HE,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='STX' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBX=DBLE(X) FF = F39HE(DBT,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,X,'T','X') STX=FF RETURN END C------------------------------------------------- F67 = TLDP REAL FUNCTION TLDP(P) CHARACTER FUN*6 DOUBLE PRECISION F67HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TLDP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F67HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TLDP=FF-REAL(T0K) RETURN END C------------------------------------------------- F69 = TMLP REAL FUNCTION TMLP(P) CHARACTER FUN*6 DOUBLE PRECISION F69HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TMLP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F69HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TMLP=FF-REAL(T0K) RETURN END C------------------------------------------------- F64 = TPH REAL FUNCTION TPH(P,H) CHARACTER FUN*6 DOUBLE PRECISION F64HE,DBP,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TPH' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBH=DBLE(H) FF = F64HE(DBP,DBH) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,H,'P','H') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TPH=FF-REAL(T0K) RETURN END C------------------------------------------------- F65 = TPS REAL FUNCTION TPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F65HE,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TPS' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBS=DBLE(S) FF = F65HE(DBP,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,S,'P','S') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TPS=FF-REAL(T0K) RETURN END C------------------------------------------------- F98 = TPSEUP C--------------------------------------------------------------- C----- FUNCTION SUBRROGRM TPSEUP(P) TO FIND PSEUDO BOILING POINT C----- AS A FUNCTION OF PRESSSURE TPSEUP.FOR C--------------------------------------------------------------- FUNCTION TPSEUP(P) DOUBLE PRECISION DBP,PBAR,F98HE,T0K CHARACTER FUN*6 COMMON /UNIT/KPA,MESS FUN='TPSEUP' C--- SET OF UNIT --- CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR C--- FUNCTION CALL --- FF=F98HE(DBP) C--- ERROR CHECK & MESSAGE --- IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') C--- SUBSDTTUDTON OF THE VALUE INTO THE FUNCTION --- IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TPSEUP=FF-REAL(T0K) RETURN END C------------------------------------------------- F70 = TPV REAL FUNCTION TPV(P,V) CHARACTER FUN*6 DOUBLE PRECISION F70HE,DBP,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TPV' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBV=DBLE(V) FF = F70HE(DBP,DBV) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,V,'P','V') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TPV=FF-REAL(T0K) RETURN END C------------------------------------------------- F41 = TRPL REAL FUNCTION TRPL(A) CHARACTER A*1 DOUBLE PRECISION F41HE,PBAR,T0K COMMON/UNIT/KPA,MESS CALL S15HE(KPA,PBAR,T0K) FF = F41HE(A) IF (MESS.EQ.0) GO TO 50 IF(FF.EQ.-1.0E+20) THEN WRITE(6,5000) A 5000 FORMAT(1H ,5X,'***** OUT OF RANGE AT TRPL FOR HELIUM WHEN', & ' A = ',A1,' *****') FF=-1.0E+20 END IF 50 IF(A.EQ.'T') THEN IF(FF.EQ.-1.0E+20) T0K=0.0 FF=FF-REAL(T0K) ELSE IF(A.EQ.'P') THEN IF(FF.EQ.-1.0E+20) PBAR=1.0 FF=FF/REAL(PBAR) END IF TRPL=FF RETURN END C------------------------------------------------- F100 = TSBP REAL FUNCTION TSBP(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='TSBP' A=P IF (MESS.NE.0) CALL S18HE(FUN) TSBP=-1.0E+30 RETURN END C------------------------------------------------- F40 = TSP REAL FUNCTION TSP(P) CHARACTER FUN*6 DOUBLE PRECISION F40HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='TSP' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F40HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') IF((FF.EQ.-1.0E+10).OR.(FF.EQ.-1.0E+20)) T0K=0.0 TSP=FF-REAL(T0K) RETURN END C------------------------------------------------- F74 = TSPD REAL FUNCTION TSPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='TSPD' A=P IF (MESS.NE.0) CALL S18HE(FUN) TSPD=-1.0E+30 RETURN END C------------------------------------------------- F75 = TSPDD REAL FUNCTION TSPDD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN='TSPDD' A=P IF (MESS.NE.0) CALL S18HE(FUN) TSPDD=-1.0E+30 RETURN END C------------------------------------------------- F42 = UPD REAL FUNCTION UPD(P) CHARACTER FUN*6 DOUBLE PRECISION F42HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F42HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') UPD=FF RETURN END C------------------------------------------------- F43 = UPDD REAL FUNCTION UPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F43HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F43HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') UPDD=FF RETURN END C------------------------------------------------- F79 = UPS REAL FUNCTION UPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F79HE,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UPS' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBS=DBLE(S) FF = F79HE(DBP,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,S,'P','S') UPS=FF RETURN END C------------------------------------------------- F44 = UPT REAL FUNCTION UPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F44HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F44HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') UPT=FF RETURN END C------------------------------------------------- F45 = UPX REAL FUNCTION UPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F45HE,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UPX' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBX=DBLE(X) FF = F45HE(DBP,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,X,'P','X') UPX=FF RETURN END C------------------------------------------------- F46 = UTD REAL FUNCTION UTD(T) CHARACTER FUN*6 DOUBLE PRECISION F46HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F46HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') UTD=FF RETURN END C------------------------------------------------- F47 = UTDD REAL FUNCTION UTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F47HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F47HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') UTDD=FF RETURN END C------------------------------------------------- F48 = UTX REAL FUNCTION UTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F48HE,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='UTX' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBX=DBLE(X) FF = F48HE(DBT,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,X,'T','X') UTX=FF RETURN END C------------------------------------------------- F49 = VPD REAL FUNCTION VPD(P) CHARACTER FUN*6 DOUBLE PRECISION F49HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VPD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F49HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') VPD=FF RETURN END C------------------------------------------------- F50 = VPDD REAL FUNCTION VPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F50HE,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VPDD' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR FF = F50HE(DBP) IF (MESS.NE.0) CALL S17HE(FF,FUN,P,'P') VPDD=FF RETURN END C------------------------------------------------- F80 = VPS REAL FUNCTION VPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F80HE,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VPS' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBS=DBLE(S) FF = F80HE(DBP,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,S,'P','S') VPS=FF RETURN END C------------------------------------------------- F51 = VPT REAL FUNCTION VPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F51HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F51HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') VPT=FF RETURN END C------------------------------------------------- F52 = VPX REAL FUNCTION VPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F52HE,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VPX' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBX=DBLE(X) FF=F52HE(DBP,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,X,'P','X') VPX=FF RETURN END C------------------------------------------------- F53 = VTD REAL FUNCTION VTD(T) CHARACTER FUN*6 DOUBLE PRECISION F53HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VTD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F53HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') VTD=FF RETURN END C------------------------------------------------- F54 = VTDD REAL FUNCTION VTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F54HE,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VTDD' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K FF = F54HE(DBT) IF (MESS.NE.0) CALL S17HE(FF,FUN,T,'T') VTDD=FF RETURN END C------------------------------------------------- F55 = VTX REAL FUNCTION VTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F55HE,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='VTX' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBX=DBLE(X) FF = F55HE(DBT,DBX) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,X,'T','X') VTX=FF RETURN END C------------------------------------------------- F83 = WPT REAL FUNCTION WPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F83HE,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='WPT' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBT=DBLE(T)+T0K FF = F83HE(DBP,DBT) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,T,'P','T') WPT=FF RETURN END C------------------------------------------------- F56 = XPH REAL FUNCTION XPH(P,H) CHARACTER FUN*6 DOUBLE PRECISION F56HE,DBP,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XPH' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBH=DBLE(H) FF = F56HE(DBP,DBH) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,H,'P','H') XPH=FF RETURN END C------------------------------------------------- F57 = XPS REAL FUNCTION XPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F57HE,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XPS' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBS=DBLE(S) FF = F57HE(DBP,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,S,'P','S') XPS=FF RETURN END C------------------------------------------------- F58 = XPU REAL FUNCTION XPU(P,U) CHARACTER FUN*6 DOUBLE PRECISION F58HE,DBP,DBU,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XPU' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBU=DBLE(U) FF = F58HE(DBP,DBU) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,U,'P','U') XPU=FF RETURN END C------------------------------------------------- F59 = XPV REAL FUNCTION XPV(P,V) CHARACTER FUN*6 DOUBLE PRECISION F59HE,DBP,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XPV' CALL S15HE(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR DBV=DBLE(V) FF = F59HE(DBP,DBV) IF (MESS.NE.0) CALL S16HE(FF,FUN,P,V,'P','V') XPV=FF RETURN END C------------------------------------------------- F60 = XTH REAL FUNCTION XTH(T,H) CHARACTER FUN*6 DOUBLE PRECISION F60HE,DBT,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XTH' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBH=DBLE(H) FF = F60HE(DBT,DBH) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,H,'T','H') XTH=FF RETURN END C------------------------------------------------- F61 = XTS REAL FUNCTION XTS(T,S) CHARACTER FUN*6 DOUBLE PRECISION F61HE,DBT,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XTS' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBS=DBLE(S) FF = F61HE(DBT,DBS) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,S,'T','S') XTS=FF RETURN END C------------------------------------------------- F62 = XTU REAL FUNCTION XTU(T,U) CHARACTER FUN*6 DOUBLE PRECISION F62HE,DBT,DBU,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XTU' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBU=DBLE(U) FF = F62HE(DBT,DBU) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,U,'T','U') XTU=FF RETURN END C------------------------------------------------- F63 = XTV REAL FUNCTION XTV(T,V) CHARACTER FUN*6 DOUBLE PRECISION F63HE,DBT,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN='XTV' CALL S15HE(KPA,PBAR,T0K) DBT=DBLE(T)+T0K DBV=DBLE(V) FF = F63HE(DBT,DBV) IF (MESS.NE.0) CALL S16HE(FF,FUN,T,V,'T','V') XTV=FF RETURN END C----------------------------------------------------S15HE SUBROUTINE S15HE(KPA,PBAR,T0K) DOUBLE PRECISION PBAR,T0K IF (KPA.EQ.0) THEN PBAR=1.0D-5 T0K=0.0D00 ELSE IF (KPA.EQ.1) THEN PBAR=1.0D00 T0K=273.15D00 ELSE IF (KPA.EQ.2) THEN PBAR=1.0D00 T0K=0.0D00 ELSE IF (KPA.EQ.3) THEN PBAR=1.0D-5 T0K=273.15D00 ELSE PBAR=1.0D-5 T0K=0.0D00 ENDIF RETURN END C--- ERROR MESSAGE (P,T)-------------------------------S16HE SUBROUTINE S16HE(FF,FUN,P,T,C1,C2) CHARACTER FUN*6,FLUID*6,C1,C2 FLUID='HELIUM' C--- LEVEL 1 ERROR CHECK & MESSAGE --- IF (FF.EQ.-1.0E10) THEN WRITE(6,910) FUN 910 FORMAT(1H ,5X,'***** NO CONVERGENCE AT ',A6,' FOR ',A6,' *****') ELSE IF (FF.EQ.-1.0E20) THEN C--- LEVEL 2 ERROR CHECK & MESSAGE --- WRITE(6,920) FUN, FLUID, C1, P, C2, T 920 FORMAT(1H ,5X,'***** OUT OF RANGE AT ',A6,' FOR ',A6,' WHEN ',A1, 1 ' =', G14.6,' AND ',A1,' =', G14.6,' *****') END IF RETURN END C--- ERROR MESSAGE P OR T------------------------------S17HE SUBROUTINE S17HE(FF,FUN,T,C1) CHARACTER FUN*6,FLUID*6,C1 FLUID='HELIUM' C--- LEVEL 1 ERROR CHECK & MESSAGE --- IF (FF.EQ.-1.0E10) THEN WRITE(6,910) FUN, FLUID 910 FORMAT(1H ,5X,'***** NO CONVERGENCE AT ',A6,' FOR ',A6,' *****') ELSE IF (FF.EQ.-1.0E20) THEN C--- LEVEL 2 ERROR CHECK & MESSAGE --- WRITE(6,920) FUN, FLUID, C1, T 920 FORMAT(1H ,5X,'***** OUT OF RANGE AT ',A6,' FOR ',A6,' WHEN ',A1, 1 ' =',G14.6,' *****') END IF RETURN END C--- ERROR MESSAGE ---------------------------------S18HE SUBROUTINE S18HE(FUN) CHARACTER FUN*6 C--- LEVEL 3 ERROR MESSAGE WRITE(6,100) FUN 100 FORMAT(1H ,5X,'***** FUNCTION ',A6,' UNAVAILABLE FOR HELIUM', 1 ' *****') RETURN END C***F2HE * ALAPP(P) LAPLACE COEFFICIENT DOUBLE PRECISION FUNCTION F2HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA GA/9.80665D00/, PLL/0.504D-01/, PC/2.18797D00/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F2HE=0. RETURN ENDIF XE=F49HE(P) IF (XE.LT.0) GO TO 910 RHOL=1.D00/XE XE=F50HE(P) IF (XE.LT.0) GO TO 910 RHOV=1.D00/XE XE=F31HE(P) IF (XE.LT.0) GO TO 910 SIG=XE F2HE=DSQRT(SIG/(GA*(RHOL-RHOV))) RETURN 900 XE=-1.E+20 910 F2HE=XE RETURN END C***F3HE * ALAPT(T) LAPLACE COEFFICIENT DOUBLE PRECISION FUNCTION F3HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.15001D0/,TMIN/2.1773D0/,GA/9.80665D0/,TC/5.15D0/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TC)/TC).LE.1.D-05) THEN F3HE=0. RETURN ENDIF XE=F53HE(TK) IF (XE.LT.0) GO TO 910 RHOL=1.D0/XE XE=F54HE(TK) IF (XE.LT.0) GO TO 910 RHOV=1.D0/XE XE=F32HE(TK) IF (XE.LT.0) GO TO 910 SIG=XE F3HE=DSQRT(SIG/(GA*(RHOL-RHOV))) RETURN 900 XE=-1.E+20 910 F3HE=XE RETURN END C***F4HE * ALHP(P) LATENT HEAT OF VAPORIZATION DOUBLE PRECISION FUNCTION F4HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, PLL/0.504D-01/, PC/2.2746D00/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F4HE=0. RETURN ENDIF XE=F23HE(P) IF (XE.LT.EV) GO TO 910 HL=XE XE=F24HE(P) IF (XE.LT.EV) GO TO 910 F4HE=XE-HL RETURN 900 XE=-1.E+20 910 F4HE=XE RETURN END C***F5HE * ALHT(T) LATENT HEAT OF VAPORIZATION DOUBLE PRECISION FUNCTION F5HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/2.1773D00/,EV/-1.D05/,TC/5.2014D0/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LE.1.D-05) THEN F5HE=0. RETURN ENDIF XE=F27HE(TK) IF (XE.LT.EV) GO TO 910 HL=XE XE=F28HE(TK) IF (XE.EQ.EV) GO TO 910 F5HE=XE-HL RETURN 900 XE=-1.E+20 910 F5HE=XE RETURN END C***F6HE * ALMPD(P) THERMAL CONDUCTIVITY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F6HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.4699D00/, DM/4.0026D00/ IF (P.LT.PLL.OR.P.GT.2.19001D00) GO TO 900 XE=F40HE(P) IF (XE.LT.0.) GO TO 910 TK=XE XE=F49HE(P) IF (XE.LT.0.) GO TO 910 RHO=1.D0/XE DD=RHO/DM CALL S04HE(TK,DD,CC,CT1,CR1,1) DPDT=101325.D00*CT1 DPDR=25314.69D00*CR1 IF (DPDR.LT.0) GO TO 900 CALL S08HE(RHO,TK,DPDT,DPDR,RLAM) F6HE=RLAM RETURN 900 XE=-1.E+20 910 F6HE=XE RETURN END C***F7HE * ALMPDD(P) THERMAL CONDUCTIVITY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F7HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.4699D00/, DM/4.0026D00/ IF (P.LT.PLL.OR.P.GT.2.19001D00) GO TO 900 XE=F40HE(P) IF (XE.LT.0.) GO TO 910 TK=XE XE=F50HE(P) IF (XE.LT.0.) GO TO 910 RHO=1.D0/XE DD=RHO/DM CALL S04HE(TK,DD,CC,CT1,CR1,1) DPDT=101325.D00*CT1 DPDR=25314.69D00*CR1 IF (DPDR.LT.0) GO TO 900 CALL S08HE(RHO,TK,DPDT,DPDR,RLAM) F7HE=RLAM RETURN 900 XE=-1.E+20 910 F7HE=XE RETURN END C***F8HE * ALMPT(P,T) THERMAL CONDUCTIVITY DOUBLE PRECISION FUNCTION F8HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D00/, RHOC/17.3987D00/,TC/5.2014D00/ IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 IF (P.LT.0.504D-01.OR.P.GT.700.D00) GO TO 900 IF (TK.LT.3.5D00.OR.TK.GT.300.D00) GO TO 900 XE=F51HE(P,TK) IF (XE.LT.0.) GO TO 910 RHO=1.D0/XE DD=RHO/DM C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) DPDT=101325.D00*((1.D0-XX)*CT3+XX*(CT1+DPT1-DPT2)) DPDR=25314.69D00*((1.D0-XX)*CR3+XX*(CR1+DPR1-DPR2)) GO TO 300 30 DL=RHOC IF (TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) DPDT=101325.D00*(CT1+DPT1-DPT2) DPDR=25314.69D00*(CR1+DPR1-DPR2) 300 IF (DPDR.LT.0) GO TO 900 CALL S08HE(RHO,TK,DPDT,DPDR,RLAM) F8HE=RLAM RETURN 900 XE=-1.E+20 910 F8HE=XE RETURN END C***F9HE * ALMTD(T) THERMAL CONDUCTIVITY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F9HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/3.5D00/, DM/4.0026D00/ IF (TK.LT.TMIN.OR.TK.GT.5.15D00) GO TO 900 XE=F53HE(TK) IF (XE.LT.0.) GO TO 910 RHO=1.D0/XE DD=RHO/DM CALL S04HE(TK,DD,CC,CT1,CR1,1) DPDT=101325.D00*CT1 DPDR=25314.69D00*CR1 IF (DPDR.LT.0) GO TO 900 CALL S08HE(RHO,TK,DPDT,DPDR,RLAM) F9HE=RLAM RETURN 900 XE=-1.E+20 910 F9HE=XE RETURN END C***F10HE * ALMTDD(T) THERMAL CONDUCTIVITY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F10HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/3.5D00/, DM/4.0026D00/ IF (TK.LT.TMIN.OR.TK.GT.5.15D00) GO TO 900 XE=F54HE(TK) IF (XE.LT.0.) GO TO 910 RHO=1.D0/XE DD=RHO/DM CALL S04HE(TK,DD,CC,CT1,CR1,1) DPDT=101325.D00*CT1 DPDR=25314.69D00*CR1 IF (DPDR.LT.0) GO TO 900 CALL S08HE(RHO,TK,DPDT,DPDR,RLAM) F10HE=RLAM RETURN 900 XE=-1.E+20 910 F10HE=XE RETURN END C***F11HE * AMUPD(P) COEFFICIENT OF VISCOSITY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F11HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.4699D00/, PC/2.2746D00/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 XE=F40HE(P) IF(XE.LT.EV) GO TO 910 TK=XE XE=F53HE(XE) IF(XE.LT.0.) GO TO 910 RHO=1.D0/XE CALL S09HE(RHO,TK,RMU) F11HE=RMU RETURN 900 XE=-1.E+20 910 F11HE=XE RETURN END C***F12HE * AMUPDD(P) COEFFICIENT OF VISCOSITY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F12HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.4699D00/, PC/2.2746D00/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 XE=F40HE(P) IF(XE.LT.EV) GO TO 910 TK=XE XE=F54HE(XE) IF(XE.LT.0.) GO TO 910 RHO=1.D0/XE CALL S09HE(RHO,TK,RMU) F12HE=RMU RETURN 900 XE=-1.E+20 910 F12HE=XE RETURN END C***F13HE * AMUPT(P,T) COEFFICIENT OF VISCOSITY DOUBLE PRECISION FUNCTION F13HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PP/700.D00/, TP/300.D00/,TL/3.5D00/ IF (P.LT.PLL.OR.P.GT.PP.OR.TK.LT.TL.OR.TK.GT.TP) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 RHO=1.D0/XE CALL S09HE(RHO,TK,RMU) F13HE=RMU RETURN 900 XE=-1.E+20 910 F13HE=XE RETURN END C***F14HE * AMUTD(T) COEFFICIENT OF VISCOSITY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F14HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.2014D00/, TMIN/3.5D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 XE=F53HE(TK) IF(XE.LT.0.) GO TO 910 RHO=1.D0/XE CALL S09HE(RHO,TK,RMU) F14HE=RMU RETURN 900 XE=-1.E+20 910 F14HE=XE RETURN END C***F15HE * AMUTDD(T) COEFFICIENT OF VISCOSITY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F15HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.2014D00/, TMIN/3.5D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 XE=F54HE(TK) IF(XE.LT.0.) GO TO 910 RHO=1.D0/XE CALL S09HE(RHO,TK,RMU) F15HE=RMU RETURN 900 XE=-1.E+20 910 F15HE=XE RETURN END C***F16HE * CPPD(P) ISOBARIC SPECIFIC HEAT OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F16HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,C3/405565.13D0/,DM/4.0026D0/, 1 PLL/0.504D-01/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.2.19001D00) GO TO 900 XE=F40HE(P) IF(XE.LT.EV) GO TO 910 TK=XE XE=F49HE(P) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) IF(DPDR.LT.0) GO TO 900 RHO=DD*DM F16HE=C1-C2*CC+C3*TK*DPDT*DPDT/(RHO*RHO*DPDR) RETURN 900 XE=-1.E+20 910 F16HE=XE RETURN END C***F17HE * CPPDD(P) ISOBARIC SPECIFIC HEAT OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F17HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,C3/405565.13D0/,DM/4.0026D0/, 1 PLL/0.504D-01/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.2.19001D00) GO TO 900 XE=F40HE(P) IF(XE.LT.EV) GO TO 910 TK=XE XE=F50HE(P) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) IF(DPDR.LT.0) GO TO 900 RHO=DD*DM F17HE=C1-C2*CC+C3*TK*DPDT*DPDT/(RHO*RHO*DPDR) RETURN 900 XE=-1.E+20 910 F17HE=XE RETURN END C***F18HE * CPPT(P,T) ISOBARIC SPECIFIC HEAT CAPACITY DOUBLE PRECISION FUNCTION F18HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, RHOC/17.3987D0/,TC/5.2014D0/, 2 CS1/3115.8D0/, CS2/2.531469D04/, CS3/1.01325D05/ IF(P.GT.2.19001D00.AND.P.LT.2.275D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.203D00) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) CV=CS2*((1.D0-XX)*C3+XX*(C1+C1DC-C2DC)) DPDT=CS3*((1.D0-XX)*CT3+XX*(CT1+DPT1-DPT2)) DPDR=CS2*((1.D0-XX)*CR3+XX*(CR1+DPR1-DPR2)) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) CV=CS2*(C1+C1DC-C2DC) DPDT=CS3*(CT1+DPT1-DPT2) DPDR=CS2*(CR1+DPR1-DPR2) 300 RHO=DD*DM IF(DPDR.LT.0) GO TO 900 F18HE=CS1-CV+TK*DPDT*DPDT/(RHO*RHO*DPDR) RETURN 900 XE=-1.E+20 910 F18HE=XE RETURN END C***F19HE * CPTD(T) ISOBARIC SPECIFIC HEAT OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F19HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,C3/405565.13D0/,DM/4.0026D0/, 1 TMIN/2.1773D00/ IF (TK.LT.TMIN.OR.TK.GT.5.15D00) GO TO 900 XE=F53HE(TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) IF(DPDR.LT.0) GO TO 900 RHO=DD*DM F19HE=C1-C2*CC+C3*TK*DPDT*DPDT/(RHO*RHO*DPDR) RETURN 900 XE=-1.E+20 910 F19HE=XE RETURN END C***F20HE * CPTDD(T) ISOBARIC SPECIFIC HEAT OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F20HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,C3/405565.13D0/,DM/4.0026D0/, 1 TMIN/2.1773D00/ IF (TK.LT.TMIN.OR.TK.GT.5.15D00) GO TO 900 XE=F54HE(TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) IF(DPDR.LT.0) GO TO 900 RHO=DD*DM F20HE=C1-C2*CC+C3*TK*DPDT*DPDT/(RHO*RHO*DPDR) RETURN 900 XE=-1.E+20 910 F20HE=XE RETURN END C***F21HE * CRP(A) QUANTITIES AT THE CRITICAL POINT DOUBLE PRECISION FUNCTION F21HE(A) CHARACTER*1 A,B(5) DATA B/'H','P','S','T','V'/ IF (A.EQ.B(1)) THEN F21HE=6740.62D00 ELSEIF (A.EQ.B(2)) THEN F21HE=2.2746D00 ELSEIF (A.EQ.B(3)) THEN F21HE=5698.79D00 ELSEIF (A.EQ.B(4)) THEN F21HE=5.2014D00 ELSEIF (A.EQ.B(5)) THEN F21HE=1.43596D-02 ELSE F21HE=-1.E+20 ENDIF RETURN END C***F23HE * HPD(P) SPECIFIC ENTHALPY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F23HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D00/, DLT/1.D-05/, PLL/0.504D-01/, EV/-1.D05/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 250 XE=F49HE(P) IF(XE.LT.0) GO TO 910 RHO=1.D0/XE XE=F42HE(P) IF(XE.LT.EV) GO TO 910 F23HE=XE+P/RHO*1.D05 RETURN 250 F23HE=6740.62D00 RETURN 900 XE=-1.E+20 910 F23HE=XE RETURN END C***F24HE * HPDD(P) SPECIFIC ENTHALPY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F24HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D00/, DLT/1.D-05/, PLL/0.504D-01/, EV/-1.D05/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 250 XE=F50HE(P) IF(XE.LT.0) GO TO 910 RHO=1.D0/XE XE=F43HE(P) IF(XE.LT.EV) GO TO 910 F24HE=XE+P/RHO*1.D05 RETURN 250 F24HE=6740.62D00 RETURN 900 XE=-1.E+20 910 F24HE=XE RETURN END C***F25HE * HPT(P,T) SPECIFIC ENTHALPY DOUBLE PRECISION FUNCTION F25HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D0/, TC/5.2014D0/, DLT/1.D-05/, EV/-1.D05/, 1 TMAX/1400.0D00/, PP/7.00D02/, PLL/0.504D-01/ IF(DABS((P-PC)/PC).LT.DLT.AND.DABS((TK-TC)/TC).LT.DLT) GO TO 250 IF(P.GT.PP.OR.TK.GT.TMAX.OR.P.LT.PLL) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.0) GO TO 910 RHO=1.D0/XE XE=F44HE(P,TK) IF(XE.LT.EV) GO TO 910 F25HE=XE+P/RHO*1.D05 RETURN 250 F25HE=6740.62D00 RETURN 900 XE=-1.E+20 910 F25HE=XE RETURN END C***F26HE * HPX(P,X) SPECIFIC ENTHALPY OF MIXTURE DOUBLE PRECISION FUNCTION F26HE(P,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/,PLL/0.504D-01/,PC/2.2746D0/ IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 XE=F23HE(P) IF(XE.LT.EV) GO TO 910 HL=XE XE=F24HE(P) IF(XE.LT.EV) GO TO 910 F26HE=HL+X*(XE-HL) RETURN 900 XE=-1.E+20 910 F26HE=XE RETURN END C***F27HE * HTD(T) SPECIFIC ENTHALPY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F27HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.2014D00/, TMIN/2.1773D00/, EV/-1.D05/, 1 DLT/1.D-05/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 250 XE=F53HE(TK) IF(XE.LT.0) GO TO 910 RHO=1.D0/XE XE=F46HE(TK) IF(XE.LT.EV) GO TO 910 F27HE=XE+F30HE(TK)/RHO*1.D05 RETURN 250 F27HE=6740.62D00 RETURN 900 XE=-1.E+20 910 F27HE=XE RETURN END C***F28HE * HTDD(T) SPECIFIC ENTHALPY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F28HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.2014D00/, TMIN/2.1773D00/, EV/-1.D05/, 1 DLT/1.D-05/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 250 XE=F54HE(TK) IF(XE.LT.0) GO TO 910 RHO=1.D0/XE XE=F47HE(TK) IF(XE.LT.EV) GO TO 910 F28HE=XE+F30HE(TK)/RHO*1.D05 RETURN 250 F28HE=6740.62D00 RETURN 900 XE=-1.E+20 910 F28HE=XE RETURN END C***F29HE * HTX(T,X) SPECIFIC ENTHALPY OF MIXTURE DOUBLE PRECISION FUNCTION F29HE(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.2014D00/, TMIN/2.1773D00/, EV/-1.D05/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F27HE(TK) IF(XE.LT.EV) GO TO 910 HL=XE XE=F28HE(TK) IF(XE.LT.EV) GO TO 910 F29HE=HL+X*(XE-HL) RETURN 900 XE=-1.E+20 910 F29HE=XE RETURN END C***F30HE * PST(T) SATURATION PRESSURE DOUBLE PRECISION FUNCTION F30HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION C(10),A(7) DATA C/ -3.9394635287D00, 141.27497598D00, -1640.7741565D00, 111974.557102D00,-55283.309818D00,166219.56504D00,-325212.82840D00, 2 398843.22750D00, -277718.06992D00, 83395.204183D00/, A/ 3 2.8815021673423D06,-3.4741773255001D06,1.7450414426205D06, 4-4.674013495016D05,7.040900480648D04,-5.655830372955D03, 5 1.892711006669D02/, PC/2.2746D00/, TC/5.2014D00/, PL/5.04D-02/, 6 TL/2.1773D00/, DLT/1.D-05/ IF(TK.LT.TL.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 300 IF(DABS((TK-TL)/TL).LT.DLT) GO TO 350 IF(TK.GT.4.8) GO TO 100 T1=(TK-1.D-03)/1.002D00 P=C(1)*T1+C(2) T2=1./T1 DO 5 I=3,10 P=P+C(I)*T2 5 T2=T2/T1 F30HE=DEXP(P)*1.3332237D-06 RETURN 100 P=0. DO 200 K=1,7 200 P=P+A(K)*TK**(K-1) F30HE=P RETURN 300 F30HE=PC RETURN 350 F30HE=PL RETURN 900 F30HE=-1.E+20 RETURN END C***F31HE * SIGP(P) SURFACE TENSION DOUBLE PRECISION FUNCTION F31HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.15D00/, PLL/0.504D-01/, PC/2.18797D00/,TM/3.15D00/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F31HE=0. RETURN ENDIF TA=F40HE(P) CN=0.5428D00 IF (TA.GT.TM) CN=1.0638D00 F31HE=0.239D-3*((TC-TA)/(TC-TM))**CN RETURN 900 XE=-1.E+20 910 F31HE=XE RETURN END C***F32HE * SIGT(T) SURFACE TENSION DOUBLE PRECISION FUNCTION F32HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.1501D0/,TMIN/2.1773D0/,TM/3.15D0/,TC/5.15D0/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TC)/TC).LE.1.D-05) THEN F32HE=0. RETURN ENDIF CN=0.5428D00 IF (TK.GT.TM) CN=1.0638D00 F32HE=0.239D-3*((TC-TK)/(TC-TM))**CN RETURN 900 XE=-1.E+20 910 F32HE=XE RETURN END C***F33HE * SPD(P) SPECIFIC ENTROPY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F33HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, PC/2.2746D0/, DLT/1.D-05/, PLL/0.504D-01/, 1 C1/1.901264178D03/, C2/2.531479538D04/, C3/5.193049518D03/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 250 XE=F49HE(P) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) XE=F40HE(P) IF(XE.LT.0) GO TO 910 TK=XE CALL S03HE(TK,DD,SS,1) F33HE=C1+C2*SS+C3*DLOG(TK) RETURN 250 F33HE=5698.79D00 RETURN 900 XE=-1.E+20 910 F33HE=XE RETURN END C***F34HE * SPDD(P) SPECIFIC ENTROPY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F34HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, PC/2.2746D0/, DLT/1.D-05/, PLL/0.504D-01/, 1 C1/1.901264178D03/, C2/2.531479538D04/, C3/5.193049518D03/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 250 XE=F50HE(P) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) XE=F40HE(P) IF(XE.LT.0) GO TO 910 TK=XE CALL S03HE(TK,DD,SS,1) F34HE=C1+C2*SS+C3*DLOG(TK) RETURN 250 F34HE=5698.79D00 RETURN 900 XE=-1.E+20 910 F34HE=XE RETURN END C***F35HE * SPT(P,T) SPECIFIC ENTROPY DOUBLE PRECISION FUNCTION F35HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, TC/5.2014D0/, RHOC/17.3987D0/, 2 C1/1.901264178D03/, C2/2.531479538D04/, C3/5.193049518D03/ XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) S1DC=0 S2DC=0 CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,50,60), INDX 10 DL=RHOC CALL S03HE(TK,DL,S1DC,1) CALL S03HE(TK,DL,S2DC,2) 20 CALL S03HE(TK,DD,S1,M) CALL S03HE(TK,DD,S3,3) SS=(1.D0-XX)*S3+XX*(S1+S1DC-S2DC) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S03HE(TK,DL,S1DC,1) CALL S03HE(TK,DL,S2DC,2) 40 CALL S03HE(TK,DD,S1,M) SS=S1+S1DC-S2DC 300 F35HE=C1+C2*SS+C3*DLOG(TK) RETURN 50 F35HE=5698.79D00 RETURN 60 XE=-1.E+20 910 F35HE=XE RETURN END C***F36HE * SPX(P,X) SPECIFIC ENTROPY OF MIXTURE DOUBLE PRECISION FUNCTION F36HE(P,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D00/, PLL/0.504D-01/, EV/-1.D05/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F33HE(P) IF(XE.LT.EV) GO TO 910 SL=XE XE=F34HE(P) IF(XE.LT.EV) GO TO 910 F36HE=SL+X*(XE-SL) RETURN 900 XE=-1.E+20 910 F36HE=XE RETURN END C***F37HE * STD(T) SPECIFIC ENTROPY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F37HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, TC/5.2014D0/, DLT/1.D-05/, TMIN/2.1773D00/, 1 C1/1.901264178D03/, C2/2.531479538D04/, C3/5.193049518D03/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 250 XE=F53HE(TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S03HE(TK,DD,SS,1) F37HE=C1+C2*SS+C3*DLOG(TK) RETURN 250 F37HE=5698.79D00 RETURN 900 XE=-1.E+20 910 F37HE=XE RETURN END C***F38HE * STDD(T) SPECIFIC ENTROPY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F38HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, TC/5.2014D0/, DLT/1.D-05/, TMIN/2.1773D00/, 1 C1/1.901264178D03/, C2/2.531479538D04/, C3/5.193049518D03/ IF (TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 250 XE=F54HE(TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S03HE(TK,DD,SS,1) F38HE=C1+C2*SS+C3*DLOG(TK) RETURN 250 F38HE=5698.79D00 RETURN 900 XE=-1.E+20 910 F38HE=XE RETURN END C***F39HE * STX(T,X) SPECIFIC ENTROPY OF MIXTURE DOUBLE PRECISION FUNCTION F39HE(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F37HE(TK) IF(XE.LT.EV) GO TO 910 SL=XE XE=F38HE(TK) IF(XE.LT.EV) GO TO 900 F39HE=SL+X*(XE-SL) RETURN 900 XE=-1.E+20 910 F39HE=XE RETURN END C***F40HE * TSP(P) SATURATION TEMPERATURE DOUBLE PRECISION FUNCTION F40HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION C(10),A(7) DATA C/ -3.9394635287D00, 141.27497598D00, -1640.7741565D00, 111974.557102D00,-55283.309818D00,166219.56504D00,-325212.82840D00, 2 398843.22750D00, -277718.06992D00, 83395.204183D00/, A/ 3 2.8815021673423D06,-3.4741773255001D06,1.7450414426205D06, 4-4.674013495016D05,7.040900480648D04,-5.655830372955D03, 5 1.892711006669D02/, PC/2.2746D00/, TC/5.2014D00/, PL/5.04D-02/, 6 TL/2.1773D00/, DLT/1.D-05/ IF(P.LT.PL.OR.P.GT.PC) GO TO 900 IC=0 IF(DABS((P-PC)/PC).LT.DLT) GO TO 300 IF(DABS((P-PL)/PL).LT.DLT) GO TO 350 IF(P.GT.1.738) GO TO 200 P0=DLOG(P/1.3332237D-06) T1=2.1773 3 F=-P0 IC=IC+1 DF=0. DO 5 K=1,10 F=F+C(K)*T1**(2-K) 5 DF=DF+DBLE(2-K)*C(K)*T1**(1-K) T0=T1-F/DF DELT=DABS((T0-T1)/T1) IF(DELT.LE.1.D-07) GO TO 10 T1=T0 IF(IC.GT.5000) GO TO 910 GO TO 3 10 F40HE=1.002D00*T0+1.D-03 RETURN 200 T1=5.3 IF(P.GT.2.04) T1=5.05 205 F=-P IC=IC+1 DF=0. DO 210 K=1,7 DF=DF+DBLE(K-1)*A(K)*T1**(K-2) 210 F=F+A(K)*T1**(K-1) T0=T1-F/DF IF(DABS((T0-T1)/T1).LE.1.D-07) GO TO 20 IF(IC.GT.5000) GO TO 910 T1=T0 GO TO 205 20 F40HE=T0 RETURN 300 F40HE=TC RETURN 350 F40HE=TL RETURN 900 F40HE=-1.E+20 RETURN 910 F40HE=-1.E+10 RETURN END C***F41HE * TRPL(A) QUANTITIES AT THE TRIPLE POINT DOUBLE PRECISION FUNCTION F41HE(A) CHARACTER*1 A,B(2) DATA B/'P', 'T'/ IF (A.EQ.B(1)) THEN F41HE=0.0504D0 ELSEIF (A.EQ.B(2)) THEN F41HE=2.1773D0 ELSE F41HE=-1.E+20 ENDIF RETURN END C***F42HE * UPD(P) SPECIFIC INTERNAL ENERGY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F42HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, CR/0.820558D-01/, PC/2.2746D00/, PL/0.504D-01/, 1 C1/4.821865787D01/, C2/2.531479538D04/, C3/5.193049518D03/ IF(P.LT.PL.OR.P.GT.PC) GO TO 900 XE=F49HE(P) IF(XE.LT.0) GO TO 910 DD=1.D0/(XE*DM) TK=F40HE(P) CALL S02HE(TK,DD,HH,1) F42HE=C1+C2*(HH-CR*TK)+C3*TK RETURN 900 XE=-1.E+20 910 F42HE=XE RETURN END C***F43HE * UPDD(P) SPECIFIC INTERNAL ENERGY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F43HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, CR/0.820558D-01/, PC/2.2746D00/, PL/0.504D-01/, 1 C1/4.821865787D01/, C2/2.531479538D04/, C3/5.193049518D03/ IF(P.LT.PL.OR.P.GT.PC) GO TO 900 XE=F50HE(P) IF(XE.LT.0) GO TO 910 DD=1.D0/(XE*DM) TK=F40HE(P) CALL S02HE(TK,DD,HH,1) F43HE=C1+C2*(HH-CR*TK)+C3*TK RETURN 900 XE=-1.E+20 910 F43HE=XE RETURN END C***F44HE * UPT(P,T) SPECIFIC INTERNAL ENERGY DOUBLE PRECISION FUNCTION F44HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/,RHOC/17.3987D0/,CR/0.820558D-01/, 1 C1/4.821865787D01/, C2/2.531479538D04/, C3/5.193049518D03/, 2 TC/5.2014D0/ XE=F51HE(P,TK) IF(XE.LT.0) GO TO 910 DD=1.D0/(XE*DM) H1DC=0 H2DC=0 CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,60), INDX 10 DL=RHOC CALL S02HE(TK,DL,H1DC,1) CALL S02HE(TK,DL,H2DC,2) 20 CALL S02HE(TK,DD,H1,M) CALL S02HE(TK,DD,H3,3) HH=(1.D0-XX)*H3+XX*(H1+H1DC-H2DC) GO TO 300 30 DL=RHOC IF(T.LT.TC) DL=G02HE(T)/DM CALL S02HE(TK,DL,H1DC,1) CALL S02HE(TK,DL,H2DC,2) 40 CALL S02HE(TK,DD,H1,M) HH=H1+H1DC-H2DC 300 F44HE=C1+C2*(HH-CR*TK)+C3*TK RETURN 60 XE=-1.E+20 910 F44HE=XE RETURN END C***F45HE * UPX(P,X) SPECIFIC INTERNAL ENERGY OF MIXTURE DOUBLE PRECISION FUNCTION F45HE(P,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PC/2.2746D00/, EV/-1.D05/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F42HE(P) IF(XE.LT.EV) GO TO 910 UL=XE XE=F43HE(P) IF(XE.LT.EV) GO TO 910 F45HE=UL+X*(XE-UL) RETURN 900 XE=-1.E+20 910 F45HE=XE RETURN END C***F46HE * UTD(T) SPECIFIC INTERNEL ENERGY OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F46HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/4.821865787D01/, C2/2.531479538D04/, C3/5.193049518D03/, 2 TMIN/2.1773D00/, TMAX/5.2014D00/, DM/4.0026D00/, 2 CR/0.820558D-01/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 XE=F53HE(TK) IF(XE.LT.0) GO TO 910 DD=1.D0/(XE*DM) CALL S02HE(TK,DD,HH,1) F46HE=C1+C2*(HH-CR*TK)+C3*TK RETURN 900 XE=-1.E+20 910 F46HE=XE RETURN END C***F47HE * UTDD(T) SPECIFIC INTERNEL ENERGY OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F47HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/4.821865787D01/, C2/2.531479538D04/, C3/5.193049518D03/, 1 TMIN/2.1773D00/, TMAX/5.2014D00/, DM/4.0026D00/, 2 CR/0.820558D-01/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 XE=F54HE(TK) IF(XE.LT.0) GO TO 910 DD=1.D0/(XE*DM) CALL S02HE(TK,DD,HH,1) F47HE=C1+C2*(HH-CR*TK)+C3*TK RETURN 900 XE=-1.E+20 910 F47HE=XE RETURN END C***F48HE * UTX(T,X) SPECIFIC INTERNAL ENERGY OF MIXTURE DOUBLE PRECISION FUNCTION F48HE(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F46HE(TK) IF(XE.LT.EV) GO TO 910 UL=XE XE=F47HE(TK) IF(XE.LT.EV) GO TO 910 F48HE=UL+X*(XE-UL) RETURN 900 XE=-1.E+20 910 F48HE=XE RETURN END C***F49HE * VPD(P) SPECIFIC VOLUME OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F49HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA DM/4.0026D0/, PC/2.2746D0/, DLT/1.D-05/, PLL/0.504D-01/, 1 C0/358.56031D0/, C1/1846.5856D0/ CALL S12HE(A,AG,DA) IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 500 TK=F40HE(P) IF(TK.GT.5.15D0) GO TO 350 RX0=G02HE(TK)/DM CALL S05HE(TK,P,RX0,2) IF(RX0.LT.0) GO TO 910 F49HE=1.D00/(RX0*DM) RETURN 350 F49HE=1.D00/(88.07D00-C0*TK+C1) RETURN 500 F49HE=1.43596D-02 RETURN 900 F49HE=-1.E+20 RETURN 910 F49HE=-1.E+10 RETURN END C***F50HE * VPDD(P) SPECIFIC VOLUME OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F50HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA DM/4.0026D0/, PC/2.2746D0/, DLT/1.D-05/, PLL/0.504D-01/, 1 C0/338.91051D0/, C1/1745.3891D0/ CALL S12HE(A,AG,DA) IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LT.DLT) GO TO 500 TK=F40HE(P) IF(TK.GT.5.15D0) GO TO 350 RX0=G03HE(TK)/DM CALL S05HE(TK,P,RX0,1) IF(RX0.LT.0) GO TO 910 F50HE=1.D00/(RX0*DM) RETURN 350 F50HE=1.D00/(52.22D00+C0*TK-C1) RETURN 500 F50HE=1.43596D-02 RETURN 900 F50HE=-1.E+20 RETURN 910 F50HE=-1.E+10 RETURN END C***F51HE * VPT(P,T) SPECIFIC VOLUME DOUBLE PRECISION FUNCTION F51HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, PATM/1.01325D0/, TC/5.2014D0/, RHOC/17.3987D0/, 1 PC/2.2746D0/, DLT/1.D-05/ FNC(Y1,Y2)=DABS((Y1-Y2)/Y2) IF(P.LE.PC) THEN PS=F30HE(TK) IF(FNC(P,PS).LT.DLT) GO TO 500 ENDIF IC=0 PIDC=0 PIIDC=0 RX0=50. XP=P/PATM CALL S13HE(P,TK,M,INDX,XX) IF(INDX.LE.2) RX0=35. IF(INDX.EQ.4.AND.M.EQ.1) RX0=0.01 GO TO (10,20,30,40,50,60), INDX 10 DD=RHOC CALL S01HE(TK,DD,PIDC,X,1,0) CALL S01HE(TK,DD,PIIDC,X,2,0) 20 CALL S01HE(TK,RX0,FR1,DFR1,M,1) IC=IC+1 CALL S01HE(TK,RX0,FR3,DFR3,3,1) F1=(1.D00-XX)*FR3+XX*(FR1+PIDC-PIIDC)-XP RX1=RX0-F1/((1.D00-XX)*DFR3+XX*DFR1) IF(FNC(RX1,RX0).LE.1.D-07) GO TO 300 RX0=RX1 IF(IC.LT.5000) GO TO 20 GO TO 910 30 DD=RHOC IF(TK.LT.TC) DD=G02HE(TK)/DM CALL S01HE(TK,DD,PIDC,X,1,0) CALL S01HE(TK,DD,PIIDC,X,2,0) 40 CALL S01HE(TK,RX0,FR0,DFR,M,1) IC=IC+1 FR=FR0-XP+PIDC-PIIDC RX1=RX0-FR/DFR IF(FNC(RX1,RX0).LE.1.D-07) GO TO 300 RX0=RX1 IF(IC.LE.5000) GO TO 40 GO TO 910 500 F51HE=F54HE(TK) RETURN 300 F51HE=1.D00/(RX1*DM) RETURN 50 F51HE=1.43596D-02 RETURN 60 F51HE=-1.E+20 RETURN 910 F51HE=-1.E+10 RETURN END C***F52HE * VPX(P,X) SPECIFIC VOLUME OF MIXTURE DOUBLE PRECISION FUNCTION F52HE(P,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PC/2.2746D00/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F49HE(P) IF(XE.LT.0) GO TO 910 VL=XE XE=F50HE(P) IF(XE.LT.0) GO TO 910 F52HE=VL+X*(XE-VL) RETURN 900 XE=-1.E+20 910 F52HE=XE RETURN END C***F53HE * VTD(T) SPECIFIC VOLUME OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F53HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA DM/4.0026D0/, TC/5.2014D0/, DLT/1.D-05/, 1 C0/358.56031D0/, C1/1846.5856D0/, TMIN/2.1773D00/ CALL S12HE(A,AG,DA) IF(TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 500 IF(TK.GT.5.15D0) GO TO 350 P=F30HE(TK) RX0=G02HE(TK)/DM CALL S05HE(TK,P,RX0,2) IF(RX0.LT.0) GO TO 910 F53HE=1.D00/(RX0*DM) RETURN 350 F53HE=1.D00/(88.07D00-C0*TK+C1) RETURN 500 F53HE=1.43596D-02 RETURN 900 F53HE=-1.E+20 RETURN 910 F53HE=-1.E+10 RETURN END C***F54HE * VTDD(T) SPECIFIC VOLUME OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F54HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA DM/4.0026D0/, TC/5.2014D0/, DLT/1.D-05/, 1 C0/338.91051D0/, C1/1745.3891D0/, TMIN/2.1773D00/ CALL S12HE(A,AG,DA) IF(TK.LT.TMIN.OR.TK.GT.TC) GO TO 900 IF(DABS((TK-TC)/TC).LT.DLT) GO TO 500 IF(TK.GT.5.15D0) GO TO 350 P=F30HE(TK) RX0=G03HE(TK)/DM CALL S05HE(TK,P,RX0,1) IF(RX0.LT.0) GO TO 910 F54HE=1.D00/(RX0*DM) RETURN 350 F54HE=1.D00/(52.22D00+C0*TK-C1) RETURN 500 F54HE=1.43596D-02 RETURN 900 F54HE=-1.E+20 RETURN 910 F54HE=-1.E+10 RETURN END C***F55HE * VTX(T,X) SPECIFIC VOLUME OF MIXTURE DOUBLE PRECISION FUNCTION F55HE(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(X.LT.0.OR.X.GT.1.D00) GO TO 900 XE=F53HE(TK) IF(XE.LT.0) GO TO 910 VL=XE XE=F54HE(TK) IF(XE.LT.0) GO TO 910 F55HE=VL+X*(XE-VL) RETURN 900 F55HE=-1.E+20 RETURN 910 F55HE=-1.E+10 RETURN END C***F56HE * XPH(P,H) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F56HE(P,H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, PLL/0.504D-01/, PC/2.2746D00/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F56HE=0. RETURN ENDIF XE=F23HE(P) IF(XE.LT.EV) GO TO 910 HL=XE XE=F24HE(P) IF(XE.LT.EV) GO TO 910 IF(H.LT.HL.OR.H.GT.XE) GO TO 900 F56HE=(H-HL)/(XE-HL) RETURN 900 XE=-1.E+20 910 F56HE=XE RETURN END C***F57HE * XPS(P,S) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F57HE(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, PLL/0.504D-01/, PC/2.2746D00/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F57HE=0. RETURN ENDIF XE=F33HE(P) IF(XE.LT.EV) GO TO 910 SL=XE XE=F34HE(P) IF(XE.LT.EV) GO TO 910 IF(S.LT.SL.OR.S.GT.XE) GO TO 900 F57HE=(S-SL)/(XE-SL) RETURN 900 XE=-1.E+20 910 F57HE=XE RETURN END C***F58HE * XPU(P,U) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F58HE(P,U) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, PLL/0.504D-01/, PC/2.2746D00/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F58HE=0. RETURN ENDIF XE=F42HE(P) IF(XE.LT.EV) GO TO 910 UL=XE XE=F43HE(P) IF(XE.LT.EV) GO TO 910 IF(U.LT.UL.OR.U.GT.XE) GO TO 900 F58HE=(U-UL)/(XE-UL) RETURN 900 XE=-1.E+20 910 F58HE=XE RETURN END C***F59HE * XPV(P,V) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F59HE(P,V) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PC/2.2746D00/ IF(P.LT.PLL.OR.P.GT.PC) GO TO 900 IF(DABS((P-PC)/PC).LE.1.D-05) THEN F59HE=0. RETURN ENDIF XE=F49HE(P) IF(XE.LT.0) GO TO 910 VL=XE XE=F50HE(P) IF(XE.LT.0) GO TO 910 IF(V.LT.VL.OR.V.GT.XE) GO TO 900 F59HE=(V-VL)/(XE-VL) RETURN 900 XE=-1.E+20 910 F59HE=XE RETURN END C***F60HE * XTH(T,H) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F60HE(TK,H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TMAX)/TMAX).LE.1.D-05) THEN F60HE=0. RETURN ENDIF XE=F27HE(TK) IF(XE.LT.EV) GO TO 910 HL=XE XE=F28HE(TK) IF(XE.LT.EV) GO TO 910 IF(H.LT.HL.OR.H.GT.XE) GO TO 900 F60HE=(H-HL)/(XE-HL) RETURN 900 XE=-1.E+20 910 F60HE=XE RETURN END C***F61HE * XTS(T,S) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F61HE(TK,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TMAX)/TMAX).LE.1.D-05) THEN F61HE=0. RETURN ENDIF XE=F37HE(TK) IF(XE.LT.EV) GO TO 910 SL=XE XE=F38HE(TK) IF(XE.LT.EV) GO TO 910 IF(S.LT.SL.OR.S.GT.XE) GO TO 900 F61HE=(S-SL)/(XE-SL) RETURN 900 XE=-1.E+20 910 F61HE=XE RETURN END C***F62HE * XTU(T,U) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F62HE(TK,U) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA EV/-1.D05/, TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TMAX)/TMAX).LE.1.D-05) THEN F62HE=0. RETURN ENDIF XE=F46HE(TK) IF(XE.LT.EV) GO TO 910 UL=XE XE=F47HE(TK) IF(XE.LT.EV) GO TO 910 IF(U.LT.UL.OR.U.GT.XE) GO TO 900 F62HE=(U-UL)/(XE-UL) RETURN 900 XE=-1.E+20 910 F62HE=XE RETURN END C***F63HE * XTV(P,V) DRYNESS FRACTION DOUBLE PRECISION FUNCTION F63HE(TK,V) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/2.1773D00/, TMAX/5.2014D00/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TMAX)/TMAX).LE.1.D-05) THEN F63HE=0. RETURN ENDIF XE=F53HE(TK) IF(XE.LT.0) GO TO 910 VL=XE XE=F54HE(TK) IF(XE.LT.0) GO TO 910 IF(V.LT.VL.OR.V.GT.XE) GO TO 900 F63HE=(V-VL)/(XE-VL) RETURN 900 XE=-1.E+20 910 F63HE=XE RETURN END C***F64HE * TPH(P,H) TEMPERATURE DOUBLE PRECISION FUNCTION F64HE(P,H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D00/, TC/5.2014D00/, PLL/0.504D-01/, PP/700.D00/, 1 TMAX/1400.D00/,DLT/1.D-05/ FNC(H1,H2)=DABS((H1-H2)/H2) IF(P.LT.PLL.OR.P.GT.PP) GO TO 900 HMAX=F25HE(P,TMAX) IF(FNC(H,HMAX).LE.DLT) GO TO 500 IF(P.GT.30.13D0) THEN T1=F69HE(P) TMIN=T1 ELSE T1=F67HE(P) TMIN=T1 ENDIF HMIN=F25HE(P,T1) IF(P.LT.0.06) THEN HL=F23HE(P) HMIN=DMIN1(HL,HMIN) ENDIF IF(FNC(H,HMIN).LE.DLT) GO TO 500 IF(H.LT.HMIN.OR.H.GT.HMAX) GO TO 900 IF(DABS((P-PC)/PC).LE.DLT) GO TO 10 IF(P.LE.PC) GO TO 10 T0=TC IF(P.LT.70.) GO TO 8 TM=F69HE(P) T0=TM 8 IF(H.GT.2.5D06) T0=TMAX H0=F25HE(P,T0) T1=T0+(H-H0)/F18HE(P,T0) IF(T1.LT.TMIN) T1=TMIN H1=F25HE(P,T1) IF(H1.LT.-1.D05) GO TO 910 GO TO 100 10 HD=F23HE(P) HDD=F24HE(P) TS=F40HE(P) IF(DABS((H-HD)/HD).LE.1.D-4.OR.DABS((H-HDD)/HDD).LE.1.D-4) 1 GO TO 30 IF(H.GE.HD) GO TO 20 T0=TS H0=HD T1=TMIN H1=HMIN IF(H1.LT.-1.D05) GO TO 910 GO TO 100 20 IF(H.LE.HDD) GO TO 30 T0=TMAX H0=F25HE(P,T0) T1=T0+(H-H0)/F18HE(P,T0) IF(T1.LT.TMIN) T1=TMIN H1=F25HE(P,T1) IF(H1.LT.-1.D05) GO TO 910 GO TO 100 30 F64HE=TS RETURN 100 CALL S10HE(P,H,T0,T1,H0,H1,TMIN) IF(T1.LT.-1.D05) GO TO 910 500 F64HE=T1 RETURN 900 F64HE=-1.E+20 RETURN 910 F64HE=-1.E+10 RETURN END C***F65HE * TPS(P,S) TEMPERATURE DOUBLE PRECISION FUNCTION F65HE(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D00/, TC/5.2014D00/, PLL/0.504D-01/, PP/700.D00/, 1 TMAX/1400.D00/,DLT/1.D-05/ FNC(S1,S2)=DABS((S1-S2)/S2) IF(P.LT.PLL.OR.P.GT.PP) GO TO 900 SMAX=F35HE(P,TMAX) IF(FNC(S,SMAX).LE.DLT) GO TO 500 IF(P.GT.30.13D00) THEN T1=F69HE(P) TMIN=T1 ELSE T1=F67HE(P) TMIN=T1 ENDIF SMIN=F35HE(P,T1) IF(P.LT.0.06) THEN SL=F33HE(P) SMIN=DMIN1(SL,SMIN) ENDIF IF(FNC(S,SMIN).LE.DLT) GO TO 500 IF(S.LT.SMIN.OR.S.GT.SMAX) GO TO 900 IF(P.LT.PC) GO TO 10 T0=TC IF(P.LT.70.) GO TO 8 TM=F69HE(P) T0=TM 8 S0=F35HE(P,T0) T1=T0+(S-S0)*T0/F18HE(P,T0) IF(T1.LT.TMIN) T1=TMIN S1=F35HE(P,T1) IF(S1.LT.-1.D05) GO TO 910 GO TO 100 10 SD=F33HE(P) SDD=F34HE(P) TS=F40HE(P) IF(DABS((S-SD)/SD).LE.1.D-4.OR.DABS((S-SDD)/SDD).LE.1.D-4) 1 GO TO 30 IF(S.GE.SD) GO TO 20 T0=TS S0=SD T1=TMIN S1=SMIN IF(S1.LT.-1.D05) GO TO 910 GO TO 100 20 IF(S.LE.SDD) GO TO 30 T0=TS S0=SDD T1=TMAX S1=SMAX GO TO 100 30 F65HE=TS RETURN 100 CALL S11HE(P,S,T0,T1,S0,S1,TMIN) IF(T1.LT.-1.D05) GO TO 910 500 F65HE=T1 RETURN 900 F65HE=-1.E+20 RETURN 910 F65HE=-1.E+10 RETURN END C***F66HE * PLDT(T) PRESSURE ON LAMBDA LINE DOUBLE PRECISION FUNCTION F66HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(6) DATA A/ -4151.5D00, -8196.2D00, -21289.D00, -34686.D00, 7.605D00, 1 91.769D00/, PL/0.504D-01/, TL/2.1773D0/, DLT/1.D-05/, TU/1.7678 2 D0/, PU/30.13D0/ IF(TK.LT.TU.OR.TK.GT.TL) GO TO 900 IF(DABS((TK-TL)/TL).LT.DLT) GO TO 30 IF(DABS((TK-TU)/TU).LT.DLT) GO TO 50 X=TK/TL-1.D00 SM=0. DO 10 K=1,4 10 SM=SM+A(K)*X**K F66HE=(SM+A(5)*(1.D00-DEXP(A(6)*X))+1.D00)*PL RETURN 30 F66HE=PL RETURN 50 F66HE=PU RETURN 900 F66HE=-1.E+20 RETURN END C***F67HE * TLDP(P) TEMPERATURE ON LAMBDA LINE DOUBLE PRECISION FUNCTION F67HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(6) DATA A/ -4151.5D00, -8196.2D00, -21289.D00, -34686.D00, 7.605D00, 1 91.769D0/, PL/0.504D-01/, TL/2.1773D0/, DLT7/1.D-07/,DLT/1.D-5/, 2 PU/30.13D0/,TU/1.7678D0/ IF(P.LT.PL.OR.P.GT.PU) GO TO 910 IC=0 IF(DABS((P-PL)/PL).LT.DLT) GO TO 30 IF(DABS((P-PU)/PU).LT.DLT) GO TO 50 T1=1. 3 F=-(P/PL-1.D00) DF=0. DO 10 K=1,4 F=F+A(K)*T1**K 10 DF=DF+DBLE(K)*A(K)*T1**(K-1) F=F+A(5)*(1.D00-DEXP(A(6)*T1)) DF=DF-A(6)*DEXP(A(6)*T1) T0=T1-F/DF IF(DABS((T0-T1)/T1).LE.DLT7) GO TO 20 T1=T0 IC=IC+1 IF(IC.GT.5000) GO TO 900 GO TO 3 20 F67HE=(T0+1.D00)*TL RETURN 30 F67HE=TL RETURN 50 F67HE=TU RETURN 900 F67HE=-1.E+10 RETURN 910 F67HE=-1.E+20 RETURN END C***F68HE * PMLT(TK) PRESSURE ON MELTING LINE DOUBLE PRECISION FUNCTION F68HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(3) DATA A/ 2.0017D00, 0.52725D00, -0.018944D00/, PT/30.43D00/, 1 TT/1.7778D0/, DLT/1.D-05/, C1/30.D0/, C2/22.904D0/,PU/30.13D0/, 2 TU/1.7678D0/, TMAX/11.023D00/ IF(TK.LT.TU.OR.TK.GT.TMAX) GO TO 900 IF(DABS((TK-TT)/TT).LT.DLT) GO TO 30 IF(DABS((TK-TU)/TU).LT.DLT) GO TO 50 IF(TK.LT.TT) GO TO 100 TS=TK/TT-1.D00 SM=0. DO 10 K=1,3 10 SM=SM+A(K)*TS**K F68HE=(SM+1.D00)*PT RETURN 30 F68HE=PT RETURN 50 F68HE=PU RETURN 100 F68HE=C1*TK-C2 RETURN 900 F68HE=-1.E+20 RETURN END C***F69HE * TMLP(P) TEMPERATURE ON MELTING LINE DOUBLE PRECISION FUNCTION F69HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(3) DATA A/ 2.0017D0, 0.52725D0, -0.018944D0/, PT/30.43D0/, DLT/1.D-5/ 1 , TT/1.7778D0/, DLT7/1.D-07/, C1/30.D0/, C2/0.7634666667D0/, 2 TU/1.7678D0/, PU/30.13D0/ IF(P.LT.30.13D00.OR.P.GT.700.D00) GO TO 910 IC=0 IF(DABS((P-PT)/PT).LT.DLT) GO TO 30 IF(DABS((P-PU)/PU).LT.DLT) GO TO 50 IF(P.LT.PT) GO TO 100 T1=1.0 3 F=-(P/PT-1.D00) DF=0. DO 10 K=1,3 F=F+A(K)*T1**K 10 DF=DF+DBLE(K)*A(K)*T1**(K-1) T0=T1-F/DF IF(DABS((T0-T1)/T1).LE.DLT7) GO TO 20 T1=T0 IC=IC+1 IF(IC.GT.5000) GO TO 900 GO TO 3 20 F69HE=(T0+1.D00)*TT RETURN 30 F69HE=TT RETURN 50 F69HE=TU RETURN 100 F69HE=P/C1+C2 RETURN 900 F69HE=-1.E+10 RETURN 910 F69HE=-1.E+20 RETURN END C***F70HE * TPV(P,V) TEMPERATURE DOUBLE PRECISION FUNCTION F70HE(P,V) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, PC/2.2746D0/, TC/5.2014D0/, PUL/30.13D0/, 1 DC/17.3987D0/, DLT/1.D-05/, PP/700.D00/, PLL/0.504D-01/ FNC(Y1,Y2)=DABS((Y1-Y2)/Y2) XP=P/1.01325D0 IRG=0 IC=0 DPT1=0 DPT2=0 PIDC=0 PIIDC=0 IF(P.GT.PP.OR.P.LT.PLL) GO TO 900 DD=1.D0/(V*DM) IF(FNC(DD,DC).LT.DLT.AND.FNC(P,PC).LT.DLT) GO TO 850 D14=1.D0/(F51HE(P,1126.85D0)*DM) IF(DD.LT.D14) GO TO 900 IF(P.GT.PUL) THEN TML=F69HE(P) DML=1.D0/(F51HE(P,TML)*DM) IF(DD.LT.DML) THEN GO TO 5 ELSE IF(FNC(DD,DML).LE.DLT) THEN GO TO 5 ELSE GO TO 900 ENDIF ELSE TLD=F67HE(P) DLD=1.D0/(F51HE(P,TLD)*DM) IF(P.LT.0.06) THEN DL=1.D0/(F49HE(P)*DM) DLD=DMAX1(DLD,DL) ENDIF IF(DD.LT.DLD) THEN GO TO 5 ELSE IF(FNC(DD,DLD).LE.DLT) THEN GO TO 5 ELSE GO TO 900 ENDIF ENDIF 5 D15=1.D0/(F51HE(P,15.0D0)*DM) IF(DD.LE.D15) GO TO 30 D10=1.D0/(F51HE(P,10.0D0)*DM) P3=F68HE(10.0D0) IF(P.GT.P3) D10=DML IF(DD.LT.D10.AND.DD.LE.DC) GO TO 40 IF(DD.LT.D10.AND.DD.GT.DC) GO TO 50 IF(P.GT.PC.AND.DD.LT.DC) GO TO 10 IF(P.GT.PC.AND.DD.GT.DC) GO TO 20 D1=1.D0/(F49HE(P)*DM) D2=1.D0/(F50HE(P)*DM) EPS1=DABS((DD-D1)/D1) EPS2=DABS((DD-D2)/D2) IF (EPS1.LT.DLT.OR.EPS2.LT.DLT) GO TO 70 IF(DD.GT.D1) GO TO 20 IF(DD.LT.D2) GO TO 10 70 F70HE=F40HE(P) RETURN 10 M=1 T1=10 GO TO 100 20 M=2 T1=20 GO TO 200 30 M=3 T1=50 GO TO 100 40 M=1 T1=50 GO TO 400 50 M=2 T1=50 510 CALL S04HE(T1,DC,Z1,DPT1,Z2,1) CALL S04HE(T1,DC,Z1,DPT2,Z2,2) CALL S01HE(T1,DC,PIDC,Z1,1,0) CALL S01HE(T1,DC,PIIDC,Z1,2,0) 400 X5=(15.D0-T1)/5.D0 X1=1.D0-X5 CALL S04HE(T1,DD,Z1,DT1,Z2,M) CALL S04HE(T1,DD,Z1,DT3,Z2,3) CALL S01HE(T1,DD,FR1,Z1,M,0) CALL S01HE(T1,DD,FR3,Z1,3,0) FR=X1*FR3+X5*(FR1+PIDC-PIIDC)-XP DPDT=X1*DT3+X5*(DT1+DPT1-DPT2) T2=T1-FR/DPDT IF(DABS((T2-T1)/T1).LT.1.D-07) GO TO 800 IC=IC+1 IF(IC.GT.5000) GO TO 910 T1=T2 IF(M.EQ.2) GO TO 510 GO TO 400 200 DL=DC 220 IF(T1.LT.TC) DL=G02HE(T1)/DM 250 CALL S04HE(T1,DL,Z2,DPT1,Z2,1) CALL S04HE(T1,DL,Z2,DPT2,Z2,2) CALL S01HE(T1,DL,PIDC,Z2,1,0) CALL S01HE(T1,DL,PIIDC,Z2,2,0) 100 CALL S04HE(T1,DD,Z1,DT1,Z2,M) CALL S01HE(T1,DD,FR1,Z2,M,0) DPDT=DT1+DPT1-DPT2 FF1=FR1+PIDC FF2=PIIDC+XP FR=FR1+PIDC-PIIDC-XP T2=T1-FR/DPDT IF(DABS((T2-T1)/T1).LT.1.D-07) GO TO 800 IC=IC+1 IF(IC.GT.5000) GO TO 910 T1=T2 IF(M.EQ.2) GO TO 220 GO TO 100 800 F70HE=T2 RETURN 850 F70HE=TC RETURN 900 F70HE=-1.E+20 RETURN 910 F70HE=-1.E+10 RETURN END C***F71HE * HPS(P,S) SPECIFIC ENTHALPY DOUBLE PRECISION FUNCTION F71HE(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PP/700.D00/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.PP) GO TO 900 XE=F65HE(P,S) IF(XE.LT.EV) GO TO 910 IF (P.GE.2.2746D00) GO TO 200 SL=F33HE(P) SV=F34HE(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 HL=F23HE(P) HV=F24HE(P) XX=F57HE(P,S) IF (XX.LT.EV) GO TO 900 F71HE=HL+XX*(HV-HL) RETURN 200 F71HE=F25HE(P,XE) RETURN 900 XE=-1.E+20 910 F71HE=XE RETURN END C***F76HE * CVPDD(P) ISOCHORIC SPECIFIC HEAT OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F76HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,DM/4.0026D0/, 1 PLL/0.504D-01/,PC/2.2746D00/,EV/-1.0D05/ IF (P.LT.PLL.OR.P.GT.PC) GO TO 900 XE=F40HE(P) IF (XE.LT.EV) GO TO 910 TK=XE XE=F50HE(P) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) F76HE=C1-C2*CC RETURN 900 XE=-1.E+20 910 F76HE=XE RETURN END C***F77HE * CVPT(P,T) ISOCHORIC SPECIFIC HEAT CAPACITY DOUBLE PRECISION FUNCTION F77HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, RHOC/17.3987D0/,TC/5.2014D0/, 2 CS1/3115.8D0/, CS2/2.531469D04/ XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) CV=CS2*((1.D0-XX)*C3+XX*(C1+C1DC-C2DC)) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) CV=CS2*(C1+C1DC-C2DC) 300 F77HE=CS1-CV RETURN 900 XE=-1.E+20 910 F77HE=XE RETURN END C***F78HE * CVTDD(T) ISOCHORIC SPECIFIC HEAT OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F78HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA C1/3115.8D0/,C2/25314.69D0/,DM/4.0026D0/, 1 TMIN/2.1773D00/, TMAX/5.2014D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 XE=F54HE(TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) CALL S04HE(TK,DD,CC,DPDT,DPDR,1) F78HE=C1-C2*CC RETURN 900 XE=-1.E+20 910 F78HE=XE RETURN END C***F79HE * UPS(P,S) SPECIFIC INTERNAL ENERGY DOUBLE PRECISION FUNCTION F79HE(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PP/700.D00/, EV/-1.D05/ IF (P.LT.PLL.OR.P.GT.PP) GO TO 900 XE=F65HE(P,S) IF (XE.LE.EV) GO TO 910 IF (P.GE.2.2746D00) GO TO 200 SL=F33HE(P) SV=F34HE(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 XE=F57HE(P,S) IF (XE.LT.EV) GO TO 910 UL=F42HE(P) UV=F43HE(P) F79HE=UL+XE*(UV-UL) RETURN 200 F79HE=F44HE(P,XE) RETURN 900 XE=-1.0E+20 910 F79HE=XE RETURN END C***F80HE * VPS(P,S) SPECIFIC VOLUME DOUBLE PRECISION FUNCTION F80HE(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PLL/0.504D-01/, PP/700.D00/, EV/-1.0E05/ IF (P.LT.PLL.OR.P.GT.PP) GO TO 900 XE=F65HE(P,S) IF (XE.LE.EV) GO TO 910 IF (P.GE.2.2746D00) GO TO 200 SL=F33HE(P) SV=F34HE(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 XE=F57HE(P,S) IF (XE.LT.EV) GO TO 910 VL=F49HE(P) VV=F50HE(P) F80HE=VL+XE*(VV-VL) RETURN 200 F80HE=F51HE(P,XE) RETURN 900 XE=-1.0E+20 910 F80HE=XE RETURN END C***F81HE * PRPT(P,T) PRANDTL NUMBER DOUBLE PRECISION FUNCTION F81HE(P,T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) AMU=F13HE(P,T) IF(AMU.LT.0) GO TO 910 ALM=F8HE(P,T) IF(ALM.LT.0) GO TO 910 CP=F18HE(P,T) IF(CP.LT.0) GO TO 910 F81HE=AMU*CP/ALM RETURN 910 F81HE=-1.0E20 RETURN END C***F82HE * AKPT(P,T) ADIABATIC EXPONENT DOUBLE PRECISION FUNCTION F82HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, RHOC/17.3987D0/,TC/5.2014D0/, 2 CS1/3115.8D0/, CS2/2.531469D04/, CS3/1.01325D05/ XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) CV=CS2*((1.D0-XX)*C3+XX*(C1+C1DC-C2DC)) DPDT=CS3*((1.D0-XX)*CT3+XX*(CT1+DPT1-DPT2)) DPDR=CS2*((1.D0-XX)*CR3+XX*(CR1+DPR1-DPR2)) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) CV=CS2*(C1+C1DC-C2DC) DPDT=CS3*(CT1+DPT1-DPT2) DPDR=CS2*(CR1+DPR1-DPR2) 300 RHO=DD*DM IF(DPDR.LT.0) GO TO 900 CPP=CS1-CV+TK*DPDT*DPDT/(RHO*RHO*DPDR) CVV=CS1-CV F82HE=CPP*RHO/(CVV*P*1.0D05)*DPDR RETURN 900 XE=-1.E+20 910 F82HE=XE RETURN END C***F83HE * WPT(P,T) VELOCITY OF SOUND DOUBLE PRECISION FUNCTION F83HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) VV=F51HE(P,TK) IF(VV.LT.0) GO TO 910 AK=F82HE(P,TK) IF(AK.LT.0) GO TO 910 F83HE=DSQRT(P*1.0D05*VV*AK) RETURN 910 F83HE=-1.0E20 RETURN END C***F85HE * PRPD(P) PRANDTL NUMBER OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F85HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) AMU=F11HE(P) IF(AMU.LT.0) GO TO 910 ALM=F6HE(P) IF(ALM.LT.0) GO TO 910 CP=F16HE(P) IF(CP.LT.0) GO TO 910 F85HE=AMU*CP/ALM RETURN 910 F85HE=-1.0E+20 RETURN END C***F86HE * PRPDD(P) PRANDTL NUMBER OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F86HE(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) AMU=F12HE(P) IF(AMU.LT.0) GO TO 910 ALM=F7HE(P) IF(ALM.LT.0) GO TO 910 CP=F17HE(P) IF(CP.LT.0) GO TO 910 F86HE=AMU*CP/ALM RETURN 910 F86HE=-1.0E+20 RETURN END C***F87HE * PRTD(T) PRANDTL NUMBER OF SATURATED LIQUID DOUBLE PRECISION FUNCTION F87HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) AMU=F14HE(TK) IF(AMU.LT.0) GO TO 910 ALM=F9HE(TK) IF(ALM.LT.0) GO TO 910 CP=F19HE(TK) IF(CP.LT.0) GO TO 910 F87HE=AMU*CP/ALM RETURN 910 F87HE=-1.0E+20 RETURN END C***F88HE * PRTDD(T) PRANDTL NUMBER OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F88HE(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) AMU=F15HE(TK) IF(AMU.LT.0) GO TO 910 ALM=F10HE(TK) IF(ALM.LT.0) GO TO 910 CP=F20HE(TK) IF(CP.LT.0) GO TO 910 F88HE=AMU*CP/ALM RETURN 910 F88HE=-1.0E+20 RETURN END C***F90HE * BSPT(P,T) ADIABATIC COMPRESSIBILITY DOUBLE PRECISION FUNCTION F90HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 VV=XE XE=F83HE(P,TK) IF (XE.LT.0) GO TO 910 WC=XE F90HE=VV/(WC*WC) RETURN 900 XE=-1.E+20 910 F90HE=XE RETURN END C***F91HE * BTPT(P,T) ISOTHERMAL COMPRESSIBILITY DOUBLE PRECISION FUNCTION F91HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2014D00) GO TO 900 XE=F18HE(P,TK) IF(XE.LT.0.) GO TO 910 CP=XE XE=F77HE(P,TK) IF(XE.LT.0.) GO TO 910 CV=XE XE=F90HE(P,TK) IF(XE.LT.0.) GO TO 910 BS=XE F91HE=CP/CV*BS RETURN 900 XE=-1.E+20 910 F91HE=XE RETURN END C***F92HE * BPPT(P,T) VOLUMETRIC COEFFICIENT OF EXPANSION DOUBLE PRECISION FUNCTION F92HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, RHOC/17.3987D0/,TC/5.2014D0/, 2 CS2/2.531469D04/, CS3/1.01325D05/ IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 FVV=XE DD=1.D0/(XE*DM) C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) DPDT=CS3*((1.D0-XX)*CT3+XX*(CT1+DPT1-DPT2)) DPDR=CS2*((1.D0-XX)*CR3+XX*(CR1+DPR1-DPR2)) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) DPDT=CS3*(CT1+DPT1-DPT2) DPDR=CS2*(CR1+DPR1-DPR2) 300 IF(DPDR.LT.0) GO TO 900 F92HE=FVV*(DPDT/DPDR) RETURN 900 XE=-1.E+20 910 F92HE=XE RETURN END C***F93HE * BVPT(P,T) PRESSURE COEFFICIENT DOUBLE PRECISION FUNCTION F93HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA DM/4.0026D0/, RHOC/17.3987D0/,TC/5.2014D0/, 2 CS3/1.01325D05/ IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 PP=1.D05*P XE=F51HE(P,TK) IF(XE.LT.0.) GO TO 910 DD=1.D0/(XE*DM) C1DC=0. C2DC=0. DPT1=0. DPT2=0. DPR1=0. DPR2=0. CALL S13HE(P,TK,M,INDX,XX) GO TO (10,20,30,40,40,900), INDX 10 DL=RHOC CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 20 CALL S04HE(TK,DD,C1,CT1,CR1,M) CALL S04HE(TK,DD,C3,CT3,CR3,3) DPDT=CS3*((1.D0-XX)*CT3+XX*(CT1+DPT1-DPT2)) GO TO 300 30 DL=RHOC IF(TK.LT.TC) DL=G02HE(TK)/DM CALL S04HE(TK,DL,C1DC,DPT1,DPR1,1) CALL S04HE(TK,DL,C2DC,DPT2,DPR2,2) 40 CALL S04HE(TK,DD,C1,CT1,CR1,M) DPDT=CS3*(CT1+DPT1-DPT2) 300 F93HE=1.D00/PP*DPDT RETURN 900 XE=-1.E+20 910 F93HE=XE RETURN END C***F94HE * AJTPT(P,T) JOULE-THOMSON COEFFICIENT DOUBLE PRECISION FUNCTION F94HE(P,TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 XE=F51HE(P,TK) IF(XE.LT.-1.D08) GO TO 910 VV=XE XE=F92HE(P,TK) IF(XE.LT.-1.D08) GO TO 910 BP=XE XE=F18HE(P,TK) IF(XE.LT.-1.D08) GO TO 910 CP=XE F94HE=VV/CP*(TK*BP-1.0D0) RETURN 900 XE=-1.E+20 910 F94HE=XE RETURN END C***F95HE * GAMPT(P,T) RATIO OF SPECIFIC HEATS DOUBLE PRECISION FUNCTION F95HE(P,TK) IMPLICIT DOUBLE PRECISION (A-H,L-Z) IF(P.GT.2.19001D00.AND.P.LT.2.2799D00.AND.TK.GT.5.15D00.AND. 1 TK.LT.5.2019D00) GO TO 900 XE=F18HE(P,TK) IF(XE.LT.0.) GO TO 910 CP=XE XE=F77HE(P,TK) IF(XE.LT.0.) GO TO 910 CV=XE F95HE=CP/CV RETURN 900 XE=-1.E+20 910 F95HE=XE RETURN END C***F96HE * GAMPDD(P) RATIO OF SPECIFIC HEATS OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F96HE(P) IMPLICIT DOUBLE PRECISION (A-H,L-Z) DATA PLL/0.504D-01/ IF (P.LT.PLL.OR.P.GT.2.19001D00) GO TO 900 XE=F17HE(P) IF(XE.LT.0) GO TO 910 CPDD=XE XE=F76HE(P) IF(XE.LT.0) GO TO 910 CVDD=XE F96HE=CPDD/CVDD RETURN 900 XE=-1.0E+20 910 F96HE=XE RETURN END C***F97HE * GAMTDD(T) RATIO OF SPECIFIC HEATS OF SATURATED VAPOR DOUBLE PRECISION FUNCTION F97HE(TK) IMPLICIT DOUBLE PRECISION (A-H,L-Z) DATA TMIN/2.1773D00/ IF (TK.LT.TMIN.OR.TK.GT.5.15D00) GO TO 900 XE=F20HE(TK) IF(XE.LT.0) GO TO 910 CPDD=XE XE=F78HE(TK) IF(XE.LT.0) GO TO 910 CVDD=XE F97HE=CPDD/CVDD RETURN 900 XE=-1.0E+20 910 F97HE=XE RETURN END C*** F98HE * TPSEUP(P) PSEUDO BOILING POINT DOUBLE PRECISION FUNCTION F98HE(P) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION T(2),C(2),TL(3),TR(3),CL(3),CR(3) P1=F21HE('P') PP=DABS((P-P1)/P1) T1=F21HE('T') IF (P.LT.P1.OR.P.GT.500.001D00) THEN F98HE=-1.0E+20 RETURN ENDIF IF (PP.LT.1.0D-5) THEN F98HE=T1 RETURN ENDIF T2=0.85*T1 P2=F30HE(T2) TM0=T1+(T1-T2)*(P-P1)/(P1-P2) C-----TM0 KINJICHI TC=T1 T(1)=TM0 DEL=T1*0.001D00 IF (P.GT.5.0.AND.P.LE.50.0) THEN T(1)=0.362*P+4.5 DEL=T(1)*0.005D00 ELSE IF (P.GT.50.0.AND.P.LE.120.0) THEN T(1)=0.257*P+10.0 DEL=T(1)*0.01D00 ELSE T(1)=0.181*P+19.0 DEL=T(1)*0.01D00 ENDIF 150 EPS=1.0D-7 IREP=0 IREM=5000 KCONT=0 ICONT=0 C(1)=F18HE(P,T(1)) T(2)=T(1)-DEL C(2)=F18HE(P,T(2)) 1000 KCONT=KCONT+1 RINC=C(2)-C(1) IF(RINC.GT.0.5)THEN GOTO 1500 ELSE T(2)=T(1) C(2)=C(1) T(1)=T(1)+DEL C(1)=F18HE(P,T(1)) IF (KCONT.GT.2000) GO TO 8000 GOTO 1000 ENDIF 1500 KCONT=0 C(1)=-C(1) C(2)=-C(2) 2000 IREP=IREP+1 IF(IREP.GT.IREM) GO TO 8000 TT=T(2)+1.3*(T(2)-T(1)) CC=-F18HE(P,TT) 3000 CONV=DABS((CC-C(2))/CC) IF(CONV.LT.EPS) THEN F98HE=TT RETURN ENDIF 4000 DEC=C(1)-CC ICONT=ICONT+1 IF(DEC.GE.0.0) THEN T(1)=TT C(1)=CC ELSE IF (ICONT.GT.5) GO TO 6000 TT=(TT+T(2))*0.5 CC=-F18HE(P,TT) GOTO 4000 ENDIF 5000 DEC=C(1)-C(2) IF(DEC.LT.0.0) THEN TT=T(1) CC=C(1) T(1)=T(2) C(1)=C(2) T(2)=TT C(2)=CC ENDIF GOTO 2000 6000 TA=T(1) TB=TT IF (TB.LT.TA) THEN TA=TT TB=T(1) ENDIF TC=TA+0.5*(TB-TA) CA=F18HE(P,TA) CB=F18HE(P,TB) 6050 KCONT=KCONT+1 IF (KCONT.GT.IREM) GO TO 8000 DELT=DABS((TA-TB)/TA) IF (DELT.LT.EPS) GO TO 7000 CC=F18HE(P,TC) DTA=(TC-TA)*0.3 TL(1)=TA TL(2)=TC-DTA TL(3)=TC DTB=(TB-TC)*0.3 TR(1)=TC TR(2)=TC+DTB TR(3)=TB CL(1)=CA CL(3)=CC CR(1)=CC CR(3)=CB CL(2)=F18HE(P,TL(2)) CR(2)=F18HE(P,TR(2)) CMXL=CL(1) ML=1 DO 6120 I=2,3 IF(CL(I).GT.CMXL) THEN CMXL=CL(I) ML=I ENDIF 6120 CONTINUE CMXR=CR(1) MR=1 DO 6130 I=2,3 IF(CR(I).GT.CMXR) THEN CMXR=CR(I) MR=I ENDIF 6130 CONTINUE IF(CMXL.GT.CMXR) THEN IF(ML.EQ.1) THEN TA=TL(1)-DTA CA=F18HE(P,TA) ELSE TA=TL(ML-1) CA=CL(ML-1) ENDIF IF(ML.EQ.3) THEN TB=TR(2) CB=CR(2) ELSE TB=TL(ML+1) CB=CL(ML+1) ENDIF IF (ML.NE.2) GO TO 6133 DELT=DABS((TL(ML)-TL(ML-1))/TL(ML)) IF(DELT.LT.1.0D-7) GO TO 7000 6133 TC=TL(ML) CC=CL(ML) ELSE IF(MR.EQ.1) THEN TA=TL(2) CA=CL(2) ELSE TA=TR(MR-1) CA=CR(MR-1) ENDIF IF(MR.EQ.3) THEN TB=TR(3)+DTB CB=F18HE(P,TB) ELSE TB=TR(MR+1) CB=CR(MR+1) ENDIF IF (MR.NE.2) GO TO 6135 DELT=DABS((TR(MR)-TR(MR-1))/TR(MR)) IF(DELT.LT.1.0D-7) GO TO 7000 6135 TC=TR(MR) CC=CR(MR) ENDIF GO TO 6050 7000 F98HE=TC RETURN 8000 F98HE=-1.0E+10 RETURN END C *** SPECIFIC VOLUME *** SUBROUTINE S01HE(TT,RR,FR0,DFR,M,KC) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA CR/8.20558D-02/, RHOC/17.3987D00/, TC/5.2014D00/, X0/0.0/, 1 X1/1.D00/ CALL S12HE(A,AG,DA) TR=TC/TT GAM=AG(M) Y=RR/RHOC SUM1=0. DO 110 K=1,9 X=DBLE(K-1) 110 SUM1=SUM1+A(K)*TR**((X-2.D0)/2.D0) T2=DSQRT(TR) T3=DSQRT(T2) CALL S06HE(M,1,8,X1,X0,X1,T2,SUM2) CALL S06HE(M,9,12,X1,X0,T2,TR,SUM3) CALL S06HE(M,13,18,X1,X0,T2,T3,SUM4) T4=TR*TR T5=T4*TR SUM5=DA(19,M)*TR+DA(20,M)*T4 SUM6=DA(21,M)*TR+DA(22,M)*T4+DA(23,M)*T5 SUM7=DA(24,M)*TR+DA(25,M)*T4+DA(26,M)*T5 Y2=Y*Y Z1=Y2*DEXP(GAM*Y2) Z2=3.D00+2.D00*GAM*Y2 Z3=5.D00+2.D00*GAM*Y2 FR0=RR*CR*TT*(1.D00+Y*(SUM1+Y*(SUM2+Y*(SUM3+Y*(SUM4+Y*SUM5)))) 1 +Z1*(SUM6+Y2*SUM7)) IF(KC.EQ.0) GO TO 300 DFR=CR*TT*(1.D00+Y*(2.D00*SUM1+Y*(3.D00*SUM2+Y*(4.D00*SUM3+Y* 1 (5.D00*SUM4+6.D00*Y*SUM5))))+Z1*(Z2*SUM6+Z3*Y2*SUM7)) RETURN 300 DFR=0. RETURN END C *** SPECIFIC INTERNAL ENERGY *** SUBROUTINE S02HE(TT,RR,HH,M) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA CR/8.20558D-2/, RHOC/17.3987D00/, TC/5.2014D0/, X0/0.0/, 1 X1/1.D00/, X2/0.5D00/, X4/0.25D00/ TR=TC/TT GAM=AG(M) Y=RR/RHOC SM1=0. T1=1.0/TR T2=DSQRT(TR) X5=-X1 X3=X2 DO 110 K=1,9 Z=A(K)*T1 SM1=SM1+X5*Z T1=T1*T2 110 X5=X5+X3 CALL S06HE(M,1,8,X0,X2,X1,T2,SM2) CALL S06HE(M,9,12,X2,X1,T2,TR,SM3) T3=DSQRT(T2) CALL S06HE(M,13,18,X2,X4,T2,T3,SM4) SM5=TR*(DA(19,M)+2.D0*DA(20,M)*TR) SM6=TR*(DA(21,M)+TR*(2.D00*DA(22,M)+3.D00*DA(23,M)*TR)) SM7=TR*(DA(24,M)+TR*(2.D00*DA(25,M)+3.D00*DA(26,M)*TR)) YY=Y*Y Z1=GAM*YY Z2=DEXP(Z1) H1=Y*(SM1+Y*(SM2/2.D0+Y*(SM3/3.D0+Y*(SM4/4.D0+Y*SM5/5.D0)))) 1+(Z2-1.)/(2.D0*GAM)*SM6+(Z2*(Z1-1.)+1.)/(2.D0*GAM*GAM)*SM7 HH=CR*TT*H1 RETURN END C *** SPECIFIC ENTROPY *** SUBROUTINE S03HE(TT,RR,S1,M) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA CR/8.20558D-2/, RHOC/17.3987D00/, TC/5.2014D0/, X1/1.D00/, 1 X2/0.5D00/, X4/0.25D00/ TR=TC/TT GAM=AG(M) Y=RR/RHOC SM1=0. T1=1.0/TR T2=DSQRT(TR) X5=-2.0D0 X3=X2 DO 110 K=1,9 SM1=SM1+X5*A(K)*T1 T1=T1*T2 110 X5=X5+X3 SSR=CR*DLOG(CR*RR*TT) CALL S06HE(M,1,8,-X1,X2,X1,T2,SM2) CALL S06HE(M,9,12,-X2,X1,T2,TR,SM3) T3=DSQRT(T2) CALL S06HE(M,13,18,-X2,X4,T2,T3,SM4) T4=TR*TR SM5=DA(20,M)*T4 SM6=DA(22,M)*T4+2.D00*DA(23,M)*TR*T4 SM7=DA(25,M)*T4+2.D00*DA(26,M)*TR*T4 Z1=GAM*Y*Y Z2=DEXP(Z1) S1=CR*(Y*(SM1+Y*(SM2/2.D0+Y*(SM3/3.D0+Y*(SM4/4.D0+Y*SM5/5.D0)))) 1+(Z2-1.)/(2.D0*GAM)*SM6+(Z2*(Z1-1.)+1.)/(2.D0*GAM*GAM)*SM7)-SSR RETURN END C *** SPECIFIC HEAT CAPACITY *** SUBROUTINE S04HE(T,RR,CV,DPT,DPR,M) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/A(9),AG(4),DA(26,3) DATA CR/8.20558D-2/, RHOC/17.3987D00/, TC/5.2014D0/, X1/1.D00/, 1 X2/0.5D00/, X4/0.25D00/, X0/0.0/ TR=TC/T GAM=AG(M) Y=RR/RHOC SR1=0. ST1=0. CV1=0. T1=1.0/TR T2=DSQRT(TR) XI=-2.0D0 XJ=-X1 XS=X2 DO 110 K=1,9 Z1=A(K)*T1 Z2=Z1*XI SR1=SR1+Z1 ST1=ST1+Z2 CV1=CV1+Z2*XJ T1=T1*T2 XJ=XJ+XS 110 XI=XI+XS CALL S07HE(M,1,8,-X1,X0,X2,X1,T2,SR2,ST2,CV2) CALL S07HE(M,9,12,-X2,X2,X1,T2,TR,SR3,ST3,CV3) T3=DSQRT(T2) CALL S07HE(M,13,18,-X2,X2,X4,T2,T3,SR4,ST4,CV4) ST5=DA(20,M)*TR*TR SR5=ST5+DA(19,M)*TR CV5=2.D00*ST5 CALL S07HE(M,21,23,X0,X1,X1,TR,TR,SR6,ST6,CV6) CALL S07HE(M,24,26,X0,X1,X1,TR,TR,SR7,ST7,CV7) Y2=Y*Y Z1=GAM*Y2 Z2=DEXP(Z1) CV=CR*(Y*(CV1+Y*(CV2/2.D0+Y*(CV3/3.D0+Y*(CV4/4.D0+Y*CV5/5.D0)))) 1+(Z2-1.)/(2.D0*GAM)*CV6+(Z2*(Z1-1.)+1.)/(2.D0*GAM*GAM)*CV7) DPT=-CR*RR*(-1.D0+Y*(ST1+Y*(ST2+Y*(ST3+Y*(ST4+Y*ST5))))+Y2* 1 Z2*(ST6+Y2*ST7)) DPR=CR*T*(1.D0+Y*(2.D0*SR1+Y*(3.D0*SR2+Y*(4.D0*SR3+Y*(5.D0*SR4+ 1 6.D0*Y*SR5))))+Y2*Z2*((3.D0+2.D0*Z1)*SR6+Y2*(5.D0+2.D0*Z1)*SR7)) RETURN END C *** SPECIFIC VOLUME OF SATURATED *** SUBROUTINE S05HE(T,P,RX0,M) IMPLICIT DOUBLE PRECISION(A-H,O-Z) PATM=1.01325D00 PIDC=0 PIIDC=0 IC=0 XP=P/PATM IF(M.EQ.1) GO TO 120 CALL S01HE(T,RX0,PIDC,X,1,0) CALL S01HE(T,RX0,PIIDC,X,2,0) 120 CALL S01HE(T,RX0,FR0,DFR,M,1) IC=IC+1 FR=FR0-XP+PIDC-PIIDC RX1=RX0-FR/DFR IF(DABS((RX0-RX1)/RX0).LT.1.D-07) GO TO 300 RX0=RX1 IF(IC.LT.5000) GO TO 120 RX0=-1.D04 RETURN 300 RX0=RX1 RETURN END C SUBROUTINE S06HE(M,K1,K2,X1,XS,T1,TS,SS) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/AA(9),AG(4),DA(26,3) SS=0. XI=X1 TI=T1 DO 100 L=K1,K2 SS=SS+XI*DA(L,M)*TI TI=TI*TS 100 XI=XI+XS RETURN END C SUBROUTINE S07HE(M,K1,K2,X1,X2,X3,T1,T2,S1,S2,S3) IMPLICIT DOUBLE PRECISION(A-H,O-Z) COMMON /L01HE/AA(9),AG(4),DA(26,3) S1=0. S2=0. S3=0. XI=X1 XJ=X2 XS=X3 TI=T1 DO 100 L=K1,K2 Z1=DA(L,M)*TI Z2=Z1*XI S1=S1+Z1 S2=S2+Z2 S3=S3+Z2*XJ TI=TI*T2 XJ=XJ+XS 100 XI=XI+XS RETURN END C *** THERMAL CONDUCTIVITY *** SUBROUTINE S08HE(RHO,T,DPDT,DPDD,RLAM) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION F(12),C(5) DATA X0/0.392D0/, F/3.726229668D0,0.186297053D-3,-0.7275964435D-6, 1-0.1427549651D-3,0.3290833592D-4,-0.5213335363D-7,0.4492659933D-7, 2-0.5924416513D-8,0.7087321137D-5,-0.6013335678D-5,0.8067145814D-6, 3 0.3995125013D-6/, E1/2.8461D0/, E2/0.27156D0/, BETA/0.3554D0/, 4 DELTA/4.304D0/, DC/69.158D0/, TC/5.18992D0/, PC/227460.D0/, 5 GAMMA/0.1743D0/, CONST/-5.882788298D00/, CON/3.4685233D-17/, 6 C/0.7034007057D0, 3.739232544D0, -26.20316969D0, 59.82252246D0, 7 -49.26397634D0/, GAMBT/0.245216657288D0/, RBT/2.81373100732D0/ T1=T**(1.D0/3.D0) T2=T1*T1 D2=RHO*RHO D3=RHO*D2 IF(T.GT.12..OR.T.LT.3.5) GO TO 5 RKT=1./DPDD/RHO DELD=DABS((RHO-DC)/DC) DELT=DABS((T-TC)/TC) R2=(DELT/0.2D0)**2+(DELD/0.25D0)**2 IF(R2.GE.1.) GO TO 20 W=DELT/DELD**RBT X1=(W+X0)/X0 XX2B=X1**0.7108D0 XX2BE=(1.+E2*XX2B)**GAMBT H=E1*X1*XX2BE DHDX=E1*XX2BE/X0+E1*E2/X0*XX2B*XX2BE/(1.+E2*XX2B)*GAMMA D2KT=(DELTA*H-W*DHDX/BETA)*(DELD**3.304D0) RKT1=DC*DC/D2/D2KT/PC RKT=R2*RKT+(1.-R2)*RKT1 IF(RKT.GT.0) GO TO 20 RLAM=-1.E+20 RETURN 20 PDT=DPDT*DPDT CALL S09HE(RHO,T,RMU) RKRT=CON*DSQRT(RKT)*T*T/RHO/RMU*PDT*DEXP(-18.66D0 1 *DELT**2-4.25D0*DELD**4) GO TO 6 5 RKRT=0. 6 A=0 TT=T DO 100 I=2,5 A=C(I)/TT+A 100 TT=TT*T RK0=T**(C(1))*DEXP(A+CONST) DL=D2*DLOG(RHO/68.D0) RLAM=RK0+F(1)*RKRT+RHO*(F(2)+F(3)*T+F(4)*T1+F(5)*T2)+D3*(F(6)+ 1 F(7)*T1+F(8)*T2)+DL*(F(9)+F(10)*T1+F(11)*T2+F(12)/T) RETURN END C *** COEFICIENT OF VISCOSITY *** SUBROUTINE S09HE(DD,T,RMU) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(5,4) DATA A/ -1.35311743D-01, 1.00347841D00, 1.20654649D00, 1 -1.49564551D-01, 1.25208416D-02, -4.75295259D01, 8.76799309D01, 2 -4.20741589D01, 8.33128289D00, -5.89252385D-01, 5.47309267D02, 3 -9.04870586D02, 4.31404928D02, -8.14504854D01, 5.37008433D00, 4 -1.68439324D03, 3.33108630D03, -1.63219172D03, 3.08804413D02, 5 -2.02936367D01/ RHO=DD*1.D-03 SM=0. DO 100 I=1,4 DO 100 J=1,5 DTLG=DLOG(T) 100 SM=SM+A(J,I)*RHO**(I-1)*DTLG**(J-2) IF(T.GT.100) GO TO 200 RMU=DEXP(SM)*1.D-07 RETURN 200 S5=0 DO 250 J=1,5 250 S5=S5+A(J,1)*DTLG**(J-2) C1=12.451D00/T-295.67D00/T**2-4.1249D00 ETA0=196.D00*T**0.71938D00*DEXP(C1) ETAE=DEXP(SM)-DEXP(S5) RMU=(ETA0+ETAE)*1.D-07 RETURN END C *** SUB. TPH(P,H,....) *** SUBROUTINE S10HE(P,H,T0,T1,H0,H1,TMIN) IMPLICIT DOUBLE PRECISION(A-H,O-Z) IC=0 100 DELT=(H-H1)*(T1-T0)/(H1-H0) RDELT=DABS(DELT/T1) IF(RDELT.GT.1.D-08) GO TO 300 IF(DABS((T0-T1)/T1).LE.1.D-07) RETURN 300 IC=IC+1 IF(IC.GT.20000) GO TO 900 T0=T1 H0=H1 T1=T1+DELT IF(T1.LT.TMIN) T1=TMIN H1=F25HE(P,T1) IF(T1.GT.0.OR.H1.GT.-1.D05) GO TO 100 900 T1=-1.E+10 RETURN END C *** SUB. TPS(P,S,...) *** SUBROUTINE S11HE(P,S,T0,T1,S0,S1,TMIN) IMPLICIT DOUBLE PRECISION(A-H,O-Z) IC=0 100 DELT=(S-S1)*(T1-T0)/(S1-S0) RDELT=DABS(DELT/T1) IF(RDELT.GT.1.D-08) GO TO 300 IF(DABS((T0-T1)/T1).LT.1.D-07) RETURN 300 IC=IC+1 IF(IC.GT.20000) GO TO 900 T0=T1 S0=S1 T1=T1+DELT IF(T1.LT.TMIN) T1=TMIN S1=F35HE(P,T1) IF(T1.GT.0.OR.S1.GT.-1.D05) GO TO 100 900 T1=-1.E+10 RETURN END C *** P,T RANGE CHECK *** SUBROUTINE S13HE(P,TK,M,IR,XX) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC/2.2746D0/, TC/5.2014D0/, DLT/1.D-05/, DELT7/1.D-07/, 1 PUL/30.13D0/, PLL/0.504D-01/, PP/700.D00/,TMAX/1400.00D00/, 2 TLL/2.1773D0/ FNC(Y1,Y2)=DABS((Y1-Y2)/Y2) XX=0 M=1 IF(P.LT.PLL.OR.P.GT.PP.OR.TK.GT.TMAX) GO TO 900 IF(FNC(P,PC).LT.DLT.AND.FNC(TK,TC).LT.DLT) GO TO 250 IF(P.LE.PUL) THEN TLD=F67HE(P) IF(TK.GT.TLD) THEN GO TO 8 ELSEIF(FNC(TK,TLD).LE.DELT7) THEN GO TO 8 ELSE GO TO 900 ENDIF ELSE TML=F69HE(P) IF(TK.GT.TML) THEN GO TO 8 ELSEIF(FNC(TK,TML).LE.DELT7) THEN GO TO 8 ELSE GO TO 900 ENDIF ENDIF 8 IF(TK.LE.TC) GO TO 18 IF(TK.GE.15.D00) GO TO 30 PL=G01HE(TK) IF(TK.LE.10.D00) GO TO 15 45 XX=(15.D00-TK)/5.D00 IF(PL-P) 50,40,40 15 IF(PL-P) 20,10,10 18 IF(TK.LT.TLL) GO TO 20 PS=F30HE(TK) IF(FNC(P,PS).LE.DLT) GO TO 10 IF(P.GT.PS) GO TO 20 10 M=1 IR=4 RETURN 30 M=3 IR=4 RETURN 40 M=1 IR=2 RETURN 50 M=2 IR=1 RETURN 20 M=2 IR=3 RETURN 250 IR=5 RETURN 900 IR=6 RETURN END C *** T,V RANGE CHECK *** SUBROUTINE S20HE(T,RHO,M,IR,XX) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.2014D0/,DLT/1.D-05/,RHOC/17.3987D00/,TMAX/1400.0D00/ FNC(Y1,Y2)=DABS((Y1-Y2)/Y2) XX=0 M=1 IF(T.GT.TMAX) GO TO 900 IF(FNC(RHO,RHOC).LT.DLT.AND.FNC(T,TC).LT.DLT) GO TO 250 IF(T.GT.15) GO TO 30 IF(T.LT.10) GO TO 10 XX=(15.D00-T)/5.D00 IF(RHO.GT.RHOC) GO TO 50 GO TO 40 10 IF(RHO.GT.RHOC) GO TO 20 M=1 IR=4 RETURN 20 M=2 IR=3 RETURN 30 M=3 IR=4 RETURN 40 M=1 IR=2 RETURN 50 M=2 IR=1 RETURN 250 IR=5 RETURN 900 IR=6 RETURN END C *** EXTENTED VAPOR PRESSURE CURVE (ISOCHORIC)(T=5.2---15.0) *** DOUBLE PRECISION FUNCTION G01HE(T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(10) DATA A/ -0.575138215565D02, 0.578150540028D02,-0.267911399538D02, 1 0.727441153705D01, -0.124019907062D01, 0.137975208666D00, 2 -0.100252994192D-1, 0.459031903684D-3, -0.120265993233D-4, 3 0.137488629073D-6/ SS=0. DO 100 K=1,10 100 SS=SS+A(K)*T**(K-1) G01HE=SS RETURN END C *** IUPAC RHOL(T) *** DOUBLE PRECISION FUNCTION G02HE(T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(6) DATA A/ 1.84869708D00, -6.19302374D00, 25.6345661D00, 1 -48.1183580D00, 43.5729630D00, -15.7686525D00/, TC/5.2014D00/, 2 RC/69.64D00/, DLT/1.D-05/ SS=0. IF(DABS((T-TC)/TC).LT.DLT) GO TO 200 DO 100 K=1,6 X=DBLE(K) 100 SS=SS+A(K)*(1.D0-T/TC)**(X/3.D0) 200 G02HE=(SS+1.D0)*RC RETURN END C *** IUPAC RHOV(T) *** DOUBLE PRECISION FUNCTION G03HE(T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(6) DATA A/ -0.994650995D00, -1.85602032D00, 4.21416868D00, 1 -5.85965799D00, 5.14208869D00, -1.62486795D00/, TC/5.2014D00/, 2 RC/69.64D00/, DLT/1.D-05/ SS=0. IF(DABS((T-TC)/TC).LT.DLT) GO TO 200 DO 100 K=1,6 X=DBLE(K) 100 SS=SS+A(K)*(1.D0-T/TC)**(X/3.D0) 200 G03HE=(SS+1.D0)*RC RETURN END C *** COMMON DATA *** SUBROUTINE S12HE(A,G,DA) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DIMENSION A(10),G(4),C1(26),C2(26),C3(26),DA(26,3),C5(9),C6(3) DATA C1/ 0.545717516825D00, -0.526218785883D01, 1 0.138438162246D02, -0.328771057024D02, 0.452233867349D02, 2-0.305702975807D02, 0.104706024245D02, -0.178805251528D01, 3 0.105274082899D01, -0.346849108795D01, 0.193273083846D01, 4 0.148012644478D00, 0.202536392881D02, -0.122920322279D03, 5 0.296428402559D03, -0.356639003545D03, 0.223090549893D03, 6-0.596732203082D02, -0.563950731866D00, 0.464238971532D00, 7 0.347980556526D01, -0.356551573938D01, 0.897149773720D00, 8 0.116282299693D01, -0.441170030307D00, -0.710142832414D-1/ DATA C2/ 0.754851995161D-1, 0.206192112563D01, 1-0.143787315035D02, 0.232146598728D02, -0.234782470378D02, 2 0.305245092459D02, -0.105330479443D02, 0.308405355051D01, 3-0.904001973585D00, 0.766114882473D00, 0.153184346013D01, 4-0.523659211353D00, -0.501716678458D01, 0.276430266139D02, 5-0.540258084825D02, 0.467422000131D02, -0.998075965096D01, 6-0.666243942985D01, -0.157966867600D00, 0.319801581443D00, 7 0.560781782190D01, -0.121896063466D02, -0.395245008176D01, 8 0.577174618615D00, -0.139495783311D01, -0.194397852603D00/ DATA C3/ -0.132911108089D00, 0.260105053391D01, 1-0.194632822744D02, 0.458381300052D02, -0.593970872032D02, 2 0.804044275465D02, -0.446727519233D02, 0.103755618662D02, 3 0.553343739698D-1, 0.926486049266D00, 0.239789555903D00, 4-0.193125134907D00, -0.114221631904D00, 0.132468957538D00, 5 0.234038221334D01, -0.711715967549D01, 0.106499189447D02, 6-0.781699044091D01, -0.552943290960D-1, 0.299638181518D00, 7 0.403394887953D01, -0.196877124146D02, 0.885278967090D-1, 8 0.230895686650D00, -0.186842775250D01, 0.306305701360D00/ DATA C5/ -0.459869970852D-4, -0.443178627200D-2, 1 0.202738009935D00, 0.568152317410D00, -0.177764093289D01, 2-0.140448471303D01, 0.253215569538D01, -0.144853086875D01, 3 0.257224547009D00/ DATA C6/-0.756786904225D00,-0.151357380845D00,-0.151357380845D00/ DO 10 K=1,26 10 DA(K,1)=C1(K) DO 15 K=1,26 15 DA(K,2)=C2(K) DO 20 K=1,26 20 DA(K,3)=C3(K) DO 25 K=1,9 25 A(K)=C5(K) DO 30 K=1,3 30 G(K)=C6(K) G(4)=1000.D00 RETURN END REAL FUNCTION T90(T) * Conversion of temperature scal from IPTS-68 to ITS-90 * Range -200C<=T68<=630C anf T68>=1064 REAL T68,T0K REAL*8 CA(8) INTEGER KPA,MESS COMMON/UNIT/KPA,MESS DATA CA/-0.14542, -0.26722, 1.06471, 1.131286, -4.25835, + -1.59924, 7.29176, -3.52573/ IF((KPA.EQ.1).OR.(KPA.EQ.3)) THEN T0K=0.0 ELSE T0K=273.15 END IF T68=T-T0K * Check of argument range IF((-200.0.LE.T68).AND.(T68.LE.630.0)) THEN RT68=T68/630.0 TEMP=T68+RT68*(CA(1)+RT68*(CA(2)+RT68*(CA(3)+RT68*(CA(4) + +RT68*(CA(5)+RT68*(CA(6)+RT68*(CA(7)+RT68*CA(8)))))))) ELSE IF((1064.0.LE.T68).AND.(T68.LE.3000.0)) THEN TEMP=T68-0.25*(T68+273.15)*(T68+273.15)/1337.33/1337.33 ELSE IF(MESS.NE.0) WRITE(6,6020) T 6020 FORMAT(1H ,5X,'**** OUT OF RANGE AT T90', & ' WHEN T68 =',1PE14.7,' ****') T90=-1.0E+20 RETURN END IF T90=TEMP+T0K RETURN END REAL FUNCTION T68(T) * Conversion of temperature scal from ITS-90 to IPTS-68 * Range -200C<=T68<=630C anf T68>=1064 T1=T90(T) T2=T ITER=0 DT=T-T1 1000 CONTINUE ITER=ITER+1 IF(ITER.GE.10000) THEN WRITE(6,*) ' ITERATION EXCEEDED MAXIMUM ' T68=-1.0E+10 RETURN END IF T2=T2+DT T1=T90(T2) DT=T-T1 IF(ABS(DT).GE.1.0E-4) GOTO 1000 T68=T2 RETURN END * Subroutine Program Specifying KPA and MESS * MS-FORTRAN and MS-C Mixed Langage Programing * SUBROUTINE KPAMES(KPAC,MESSC) COMMON/UNIT/KPA,MESS INTEGER KPAC INTEGER MESSC KPA=KPAC MESS=MESSC RETURN END