C PROPATH VER. 10.1 HELIUM 4(NIST-ITS 1990) 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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(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 S18A04(FUN) WTDD=-1.0E+30 RETURN END C C==================================================================== C------------------------------------------- F1 : AIPPT [(mol/kg)^2] REAL FUNCTION AIPPT(P,T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'AIPPT' A = P+T IF (MESS.NE.0) CALL S18A04(FUN) AIPPT=-1.0E30 RETURN END C------------------------------------------- F94 : AJTPT [K/Pa] REAL FUNCTION AJTPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F94A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AJTPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F94A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') AJTPT = FF RETURN END C------------------------------------------- F82 : AKPT [-] REAL FUNCTION AKPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F82A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AKPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F82A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') AKPT = FF RETURN END C------------------------------------------- F2 : ALAPP [m] REAL FUNCTION ALAPP(P) CHARACTER FUN*6 DOUBLE PRECISION F2A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALAPP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F2A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') ALAPP = FF RETURN END C------------------------------------------- F3 : ALAPT [m] REAL FUNCTION ALAPT(T) CHARACTER FUN*6 DOUBLE PRECISION F3A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALAPT' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F3A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') ALAPT = FF RETURN END C------------------------------------------- F4 : ALHP [J/kg] REAL FUNCTION ALHP(P) CHARACTER FUN*6 DOUBLE PRECISION F4A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALHP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F4A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') ALHP = FF RETURN END C------------------------------------------- F5 : ALHT [J/kg] REAL FUNCTION ALHT(T) CHARACTER FUN*6 DOUBLE PRECISION F5A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALHT' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F5A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') ALHT = FF RETURN END C------------------------------------------- F6 : ALMPD [W/(m*K)] REAL FUNCTION ALMPD(P) CHARACTER FUN*6 DOUBLE PRECISION F6A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALMPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F6A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') ALMPD = FF RETURN END C------------------------------------------- F7 : ALMPDD [W/(m*K)] REAL FUNCTION ALMPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F7A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALMPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F7A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') ALMPDD = FF RETURN END C------------------------------------------- F8 : ALMPT [W/(m*K)] REAL FUNCTION ALMPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F8A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALMPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F8A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') ALMPT = FF RETURN END C------------------------------------------- F9 : ALMTD [W/(m*K)] REAL FUNCTION ALMTD(T) CHARACTER FUN*6 DOUBLE PRECISION F9A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALMTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F9A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') ALMTD = FF RETURN END C------------------------------------------- F10 : ALMTDD [W/(m*K)] REAL FUNCTION ALMTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F10A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'ALMTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F10A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') ALMTDD = FF RETURN END C------------------------------------------- F11 : AMUPD [Ps*s] REAL FUNCTION AMUPD(P) CHARACTER FUN*6 DOUBLE PRECISION F11A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AMUPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F11A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') AMUPD = FF RETURN END C------------------------------------------- F12 : AMUPDD [Ps*s] REAL FUNCTION AMUPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F12A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AMUPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F12A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') AMUPDD = FF RETURN END C------------------------------------------- F13 : AMUPT [Ps*s] REAL FUNCTION AMUPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F13A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AMUPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F13A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') AMUPT = FF RETURN END C------------------------------------------- F14 : AMUTD [Ps*s] REAL FUNCTION AMUTD(T) CHARACTER FUN*6 DOUBLE PRECISION F14A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AMUTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F14A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') AMUTD = FF RETURN END C------------------------------------------- F15 : AMUTDD [Ps*s] REAL FUNCTION AMUTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F15A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'AMUTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F15A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') AMUTDD = FF RETURN END C------------------------------------------- F92 : BPPT [1/K] REAL FUNCTION BPPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F92A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'BPPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F92A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(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 F90A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'BSPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F90A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(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 F91A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'BTPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F91A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(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 F93A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'BVPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F93A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') BVPT = FF RETURN END C------------------------------------------- F16 : CPPD [J/(kg*K)] REAL FUNCTION CPPD(P) CHARACTER FUN*6 DOUBLE PRECISION F16A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CPPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F16A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') CPPD = FF RETURN END C------------------------------------------- F17 : CPPDD [J/(kg*K)] REAL FUNCTION CPPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F17A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CPPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F17A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') CPPDD = FF RETURN END C------------------------------------------- F18 : CPPT [J/(kg*K)] REAL FUNCTION CPPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F18A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CPPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F18A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') CPPT = FF RETURN END C------------------------------------------- F19 : CPTD [J/(kg*K)] REAL FUNCTION CPTD(T) CHARACTER FUN*6 DOUBLE PRECISION F19A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CPTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F19A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') CPTD = FF RETURN END C------------------------------------------- F20 : CPTDD [J/(kg*K)] REAL FUNCTION CPTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F20A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CPTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F20A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') CPTDD = FF RETURN END C------------------------------------------- F21 : CRP REAL FUNCTION CRP(A) CHARACTER A*1 DOUBLE PRECISION F21A04,PBAR,T0K COMMON/UNIT/KPA,MESS CALL S15A04(KPA,PBAR,T0K) FF = F21A04(A) IF (MESS.EQ.0) GO TO 50 IF(FF.EQ.-1.0E20) THEN WRITE(6,6010) A 6010 FORMAT(1H ,5X,'***** OUT OF RANGE AT CRP FOR HELIUM WHEN ', & 'A = ',A1,' *****') FF=-1.0E20 END IF 50 IF(A.EQ.'T') THEN IF(FF.EQ.-1.0E20) T0K=0.0 FF = FF-REAL(T0K) ELSE IF(A.EQ.'P') THEN IF(FF.EQ.-1.0E20) PBAR=1.0 FF = FF/REAL(PBAR) END IF CRP = FF RETURN END C------------------------------------------- F76 : CVPDD [J/(kg*K)] REAL FUNCTION CVPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F76A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CVPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F76A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') CVPDD = FF RETURN END C------------------------------------------- F77 : CVPT [J/(kg*K)] REAL FUNCTION CVPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F77A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CVPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F77A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') CVPT = FF RETURN END C------------------------------------------- F78 : CVTDD [J/(kg*K)] REAL FUNCTION CVTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F78A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'CVTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F78A04(DBT)) IF (MESS.NE.0) CALL S17A04(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 S18A04(FUN) EPSPT = -1.0E30 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.E20 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 F95A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMPT'/ CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F95A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(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 F96A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMPDD'/ CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F96A04(DBP)) IF (MESS.NE.0) CALL S17A04(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 F97A04,PBAR,DBT,T0K COMMON/UNIT/KPA,MESS DATA FUN/'GAMTDD'/ CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F97A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') GAMTDD = FF RETURN END C------------------------------------------- F23 : HPD [J/kg] REAL FUNCTION HPD(P) CHARACTER FUN*6 DOUBLE PRECISION F23A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F23A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') HPD = FF RETURN END C------------------------------------------- F24 : HPDD [J/kg] REAL FUNCTION HPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F24A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F24A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') HPDD = FF RETURN END C------------------------------------------- F71 : HPS [J/kg] REAL FUNCTION HPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F71A04,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HPS' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBS = DBLE(S) FF = REAL(F71A04(DBP,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,S,'P','S') HPS = FF RETURN END C------------------------------------------- F25 : HPT [J/kg] REAL FUNCTION HPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F25A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F25A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') HPT = FF RETURN END C------------------------------------------- F26 : HPX [J/kg] REAL FUNCTION HPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F26A04,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HPX' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBX = DBLE(X) FF = REAL(F26A04(DBP,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,X,'P','X') HPX = FF RETURN END C------------------------------------------- F27 : HTD [J/kg] REAL FUNCTION HTD(T) CHARACTER FUN*6 DOUBLE PRECISION F27A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F27A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') HTD = FF RETURN END C------------------------------------------- F28 : HTDD [J/kg] REAL FUNCTION HTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F28A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F28A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') HTDD = FF RETURN END C------------------------------------------- F29 : HTX [J/kg] REAL FUNCTION HTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F29A04,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'HTX' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBX = DBLE(X) FF = REAL(F29A04(DBT,DBX)) IF (MESS.NE.0) CALL S16A04(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='HE' WHEN A = 'C' C B='12.1' WHEN A = 'V' C************************************************ CHARACTER*25 FUNCTION IDENTF(A) CHARACTER A*1,MSG*120 COMMON/UNIT/KPA,MESS IF (A.EQ.'S') THEN IDENTF='HELIUM 4(ITS 1990)' ELSE IF (A.EQ.'C') THEN IDENTF='HE' 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 [Pa] REAL FUNCTION PLDT(T) CHARACTER FUN*6 DOUBLE PRECISION F66A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PLDT' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F66A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) PBAR=1.0 PLDT = FF/REAL(PBAR) RETURN END C------------------------------------------- F68 : PMLT [Pa] REAL FUNCTION PMLT(T) CHARACTER FUN*6 DOUBLE PRECISION F68A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PMLT' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F68A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) PBAR=1.0 PMLT = FF/REAL(PBAR) RETURN END C------------------------------------------- F85 : PRPD [-] REAL FUNCTION PRPD(P) CHARACTER FUN*6 DOUBLE PRECISION F85A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PRPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F85A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') PRPD = FF RETURN END C------------------------------------------- F86 : PRPDD [-] REAL FUNCTION PRPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F86A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PRPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F86A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') PRPDD = FF RETURN END C------------------------------------------- F81 : PRPT [-] REAL FUNCTION PRPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F81A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PRPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F81A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') PRPT = FF RETURN END C------------------------------------------- F87 : PRTD [-] REAL FUNCTION PRTD(T) CHARACTER FUN*6 DOUBLE PRECISION F87A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PRTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F87A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') PRTD = FF RETURN END C------------------------------------------- F88 : PRTDD [-] REAL FUNCTION PRTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F88A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PRTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F88A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') PRTDD = FF RETURN END C------------------------------------------- F99 : PSBT [Pa] REAL FUNCTION PSBT(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'PSBT' A = T IF (MESS.NE.0) CALL S18A04(FUN) PSBT=-1.0E30 RETURN END C------------------------------------------- F30 : PST [Pa] REAL FUNCTION PST(T) CHARACTER FUN*6 DOUBLE PRECISION F30A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'PST' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F30A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) PBAR=1.0 PST = FF/REAL(PBAR) RETURN END C------------------------------------------- F72 : PSTD [Pa] REAL FUNCTION PSTD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'PSTD' A = T IF (MESS.NE.0) CALL S18A04(FUN) PSTD=-1.0E30 RETURN END C------------------------------------------- F73 : PSTDD [Pa] REAL FUNCTION PSTDD(T) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'PSTDD' A = T IF (MESS.NE.0) CALL S18A04(FUN) PSTDD=-1.0E30 RETURN END C------------------------------------------- F31 : SIGP [N/m] REAL FUNCTION SIGP(P) CHARACTER FUN*6 DOUBLE PRECISION F31A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SIGP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F31A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') SIGP = FF RETURN END C------------------------------------------- F32 : SIGT [N/m] REAL FUNCTION SIGT(T) CHARACTER FUN*6 DOUBLE PRECISION F32A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SIGT' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F32A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') SIGT = FF RETURN END C------------------------------------------- F33 : SPD [J/(kg*K)] REAL FUNCTION SPD(P) CHARACTER FUN*6 DOUBLE PRECISION F33A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F33A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') SPD = FF RETURN END C------------------------------------------- F34 : SPDD [J/(kg*K)] REAL FUNCTION SPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F34A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F34A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') SPDD = FF RETURN END C------------------------------------------- F35 : SPT [J/(kg*K)} REAL FUNCTION SPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F35A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F35A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') SPT = FF RETURN END C------------------------------------------- F36 : SPX [J/(kg*K)] REAL FUNCTION SPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F36A04,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'SPX' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBX = DBLE(X) FF = REAL(F36A04(DBP,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,X,'P','X') SPX = FF RETURN END C------------------------------------------- F37 : STD [J/(kg*K)] REAL FUNCTION STD(T) CHARACTER FUN*6 DOUBLE PRECISION F37A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'STD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F37A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') STD = FF RETURN END C------------------------------------------- F38 : STDD [J/(kg*K)] REAL FUNCTION STDD(T) CHARACTER FUN*6 DOUBLE PRECISION F38A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'STDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F38A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') STDD = FF RETURN END C------------------------------------------- F39 : STX [J/(kg*K)] REAL FUNCTION STX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F39A04,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'STX' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBX = DBLE(X) FF = REAL(F39A04(DBT,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,X,'T','X') STX = FF RETURN END C------------------------------------------- F67 : TLDP [K] REAL FUNCTION TLDP(P) CHARACTER FUN*6 DOUBLE PRECISION F67A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TLDP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F67A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TLDP = FF-REAL(T0K) RETURN END C------------------------------------------- F69 : TMLP [K] REAL FUNCTION TMLP(P) CHARACTER FUN*6 DOUBLE PRECISION F69A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TMLP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F69A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TMLP = FF-REAL(T0K) RETURN END C------------------------------------------- F64 : TPH [K] REAL FUNCTION TPH(P,H) CHARACTER FUN*6 DOUBLE PRECISION F64A04,DBP,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TPH' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBH = DBLE(H) FF = REAL(F64A04(DBP,DBH)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,H,'P','H') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TPH = FF-REAL(T0K) RETURN END C------------------------------------------- F65 : TPS [K] REAL FUNCTION TPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F65A04,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TPS' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBS = DBLE(S) FF = REAL(F65A04(DBP,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,S,'P','S') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TPS = FF-REAL(T0K) RETURN END C------------------------------------------- F98 : TPSEUP [K] 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,F98A04,T0K CHARACTER FUN*6 COMMON /UNIT/KPA,MESS FUN = 'TPSEUP' C--- SET OF UNIT --- CALL S15A04(KPA,PBAR,T0K) DBP=DBLE(P)*PBAR C--- FUNCTION CALL --- FF=REAL(F98A04(DBP)) C--- ERROR CHECK & MESSAGE --- IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') C--- SUBSDTTUDTON OF THE VALUE INTO THE FUNCTION --- IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TPSEUP = FF-REAL(T0K) RETURN END C------------------------------------------- F70 : TPV [K] REAL FUNCTION TPV(P,V) CHARACTER FUN*6 DOUBLE PRECISION F70A04,DBP,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TPV' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBV = DBLE(V) FF = REAL(F70A04(DBP,DBV)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,V,'P','V') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TPV = FF-REAL(T0K) RETURN END C------------------------------------------- F41 : TRPL REAL FUNCTION TRPL(A) CHARACTER A*1 DOUBLE PRECISION F41A04,PBAR,T0K COMMON/UNIT/KPA,MESS CALL S15A04(KPA,PBAR,T0K) FF = REAL(F41A04(A)) IF (MESS.EQ.0) GO TO 50 IF(FF.EQ.-1.0E20) THEN WRITE(6,5000) A 5000 FORMAT(1H ,5X,'***** OUT OF RANGE AT TRPL FOR HELIUM WHEN', & ' A = ',A1,' *****') FF = -1.0E20 END IF 50 IF(A.EQ.'T') THEN IF(FF.EQ.-1.0E20) 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 [K] REAL FUNCTION TSBP(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'TSBP' A = P IF (MESS.NE.0) CALL S18A04(FUN) TSBP = -1.0E30 RETURN END C------------------------------------------- F40 : TSP [K] REAL FUNCTION TSP(P) CHARACTER FUN*6 DOUBLE PRECISION F40A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'TSP' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F40A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') IF((FF.EQ.-1.0E10).OR.(FF.EQ.-1.0E20)) T0K = 0.0 TSP = FF-REAL(T0K) RETURN END C------------------------------------------- F74 : TSPD [K] REAL FUNCTION TSPD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'TSPD' A = P IF (MESS.NE.0) CALL S18A04(FUN) TSPD = -1.0E30 RETURN END C------------------------------------------- F75 : TSPDD [K] REAL FUNCTION TSPDD(P) CHARACTER FUN*6 COMMON/UNIT/KPA,MESS FUN = 'TSPDD' A = P IF (MESS.NE.0) CALL S18A04(FUN) TSPDD = -1.0E30 RETURN END C------------------------------------------- F42 : UPD [J/kg] REAL FUNCTION UPD(P) CHARACTER FUN*6 DOUBLE PRECISION F42A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F42A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') UPD = FF RETURN END C------------------------------------------- F43 : UPDD [J/kg] REAL FUNCTION UPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F43A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F43A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') UPDD = FF RETURN END C------------------------------------------- F79 : UPS [J/kg] REAL FUNCTION UPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F79A04,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UPS' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBS = DBLE(S) FF = REAL(F79A04(DBP,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,S,'P','S') UPS = FF RETURN END C------------------------------------------- F44 : UPT [J/kg] REAL FUNCTION UPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F44A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F44A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') UPT = FF RETURN END C------------------------------------------- F45 : UPX [J/kg] REAL FUNCTION UPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F45A04,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UPX' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBX = DBLE(X) FF = REAL(F45A04(DBP,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,X,'P','X') UPX = FF RETURN END C------------------------------------------- F46 : UTD [J/kg] REAL FUNCTION UTD(T) CHARACTER FUN*6 DOUBLE PRECISION F46A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F46A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') UTD = FF RETURN END C------------------------------------------- F47 : UTDD [J/kg] REAL FUNCTION UTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F47A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F47A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') UTDD = FF RETURN END C------------------------------------------- F48 : UTX [J/kg] REAL FUNCTION UTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F48A04,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'UTX' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBX = DBLE(X) FF = REAL(F48A04(DBT,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,X,'T','X') UTX = FF RETURN END C------------------------------------------- F49 : VPD [m^3/kg] REAL FUNCTION VPD(P) CHARACTER FUN*6 DOUBLE PRECISION F49A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VPD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F49A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') VPD = FF RETURN END C------------------------------------------- F50 : VPDD [m^3/kg] REAL FUNCTION VPDD(P) CHARACTER FUN*6 DOUBLE PRECISION F50A04,DBP,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VPDD' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR FF = REAL(F50A04(DBP)) IF (MESS.NE.0) CALL S17A04(FF,FUN,P,'P') VPDD = FF RETURN END C------------------------------------------- F80 : VPS [m^3/kg] REAL FUNCTION VPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F80A04,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VPS' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBS = DBLE(S) FF = REAL(F80A04(DBP,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,S,'P','S') VPS = FF RETURN END C------------------------------------------- F51 : VPT [m^3/kg] REAL FUNCTION VPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F51A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F51A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') VPT = FF RETURN END C------------------------------------------- F52 : VPX [m^3/kg] REAL FUNCTION VPX(P,X) CHARACTER FUN*6 DOUBLE PRECISION F52A04,DBP,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VPX' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBX = DBLE(X) FF = REAL(F52A04(DBP,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,X,'P','X') VPX = FF RETURN END C------------------------------------------- F53 : VTD [m^3/kg] REAL FUNCTION VTD(T) CHARACTER FUN*6 DOUBLE PRECISION F53A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VTD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F53A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') VTD = FF RETURN END C------------------------------------------- F54 : VTDD [m^3/kg] REAL FUNCTION VTDD(T) CHARACTER FUN*6 DOUBLE PRECISION F54A04,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VTDD' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K FF = REAL(F54A04(DBT)) IF (MESS.NE.0) CALL S17A04(FF,FUN,T,'T') VTDD = FF RETURN END C------------------------------------------- F55 : VTX [m^3/kg] REAL FUNCTION VTX(T,X) CHARACTER FUN*6 DOUBLE PRECISION F55A04,DBT,DBX,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'VTX' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBX = DBLE(X) FF = REAL(F55A04(DBT,DBX)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,X,'T','X') VTX = FF RETURN END C------------------------------------------- F83 : WPT [m/s] REAL FUNCTION WPT(P,T) CHARACTER FUN*6 DOUBLE PRECISION F83A04,DBP,DBT,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'WPT' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBT = DBLE(T)+T0K FF = REAL(F83A04(DBP,DBT)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,T,'P','T') WPT = FF RETURN END C------------------------------------------- F56 : XPH [-] REAL FUNCTION XPH(P,H) CHARACTER FUN*6 DOUBLE PRECISION F56A04,DBP,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XPH' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBH = DBLE(H) FF = REAL(F56A04(DBP,DBH)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,H,'P','H') XPH = FF RETURN END C------------------------------------------- F57 : XPS [-] REAL FUNCTION XPS(P,S) CHARACTER FUN*6 DOUBLE PRECISION F57A04,DBP,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XPS' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBS = DBLE(S) FF = REAL(F57A04(DBP,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,S,'P','S') XPS = FF RETURN END C------------------------------------------- F58 : XPU [-] REAL FUNCTION XPU(P,U) CHARACTER FUN*6 DOUBLE PRECISION F58A04,DBP,DBU,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XPU' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBU = DBLE(U) FF = REAL(F58A04(DBP,DBU)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,U,'P','U') XPU = FF RETURN END C------------------------------------------- F59 : XPV [-] REAL FUNCTION XPV(P,V) CHARACTER FUN*6 DOUBLE PRECISION F59A04,DBP,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XPV' CALL S15A04(KPA,PBAR,T0K) DBP = DBLE(P)*PBAR DBV = DBLE(V) FF = REAL(F59A04(DBP,DBV)) IF (MESS.NE.0) CALL S16A04(FF,FUN,P,V,'P','V') XPV = FF RETURN END C------------------------------------------- F60 : XTH [-] REAL FUNCTION XTH(T,H) CHARACTER FUN*6 DOUBLE PRECISION F60A04,DBT,DBH,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XTH' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBH = DBLE(H) FF = REAL(F60A04(DBT,DBH)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,H,'T','H') XTH = FF RETURN END C------------------------------------------- F61 : XTS [-] REAL FUNCTION XTS(T,S) CHARACTER FUN*6 DOUBLE PRECISION F61A04,DBT,DBS,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XTS' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBS = DBLE(S) FF = REAL(F61A04(DBT,DBS)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,S,'T','S') XTS = FF RETURN END C------------------------------------------- F62 : XTU [-] REAL FUNCTION XTU(T,U) CHARACTER FUN*6 DOUBLE PRECISION F62A04,DBT,DBU,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XTU' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBU = DBLE(U) FF = REAL(F62A04(DBT,DBU)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,U,'T','U') XTU = FF RETURN END C------------------------------------------- F63 : XTV [-] REAL FUNCTION XTV(T,V) CHARACTER FUN*6 DOUBLE PRECISION F63A04,DBT,DBV,PBAR,T0K COMMON/UNIT/KPA,MESS FUN = 'XTV' CALL S15A04(KPA,PBAR,T0K) DBT = DBLE(T)+T0K DBV = DBLE(V) FF = REAL(F63A04(DBT,DBV)) IF (MESS.NE.0) CALL S16A04(FF,FUN,T,V,'T','V') XTV = FF RETURN END C----------------------------------------------------S15A04 SUBROUTINE S15A04(KPA,PBAR,T0K) DOUBLE PRECISION PBAR,T0K IF (KPA.EQ.0) THEN PBAR=1.0D00 T0K=0.0D00 ELSE IF (KPA.EQ.1) THEN PBAR=1.0D05 T0K=273.15D00 ELSE IF (KPA.EQ.2) THEN PBAR=1.0D05 T0K=0.0D00 ELSE IF (KPA.EQ.3) THEN PBAR=1.0D00 T0K=273.15D00 ELSE PBAR=1.0D-5 T0K=0.0D00 ENDIF RETURN END C--- ERROR MESSAGE (P,T)-------------------------------S16A04 SUBROUTINE S16A04(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 S17A04(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 S18A04(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 THIS FILE IS EXTRACTED FROM THE APRIL 1990 VERSION OF NIST'S HEPROP C CODE. REFERENCE: STATE EQUATION OF LIQUID HELIUM FROM 0.8 TO 2.5 C JOURNAL OF LOW TEMPERATURE PHYSICS 79, 93, 1990, V. ARP C C THE COMPLETE HEPROP CODE, VALID FROM 0.8 TO 1500 K, IS SOLD BY NIS C FOR INFORMATION, CONTACT DR. DAN FRIEND, THERMOPHYSICS DIVISION, C NIST, 325 BROADWAY, BOULDER CO 80303. C C UPDATED VERSIONS, UNDER THE NAME HEPAK, ARE SOLD BY CRYODATA, C PO BOX 558, NIWOT CO 80544. IMPROVEMENTS HAVE BEEN MADE IN SPEED C OF EXECUTION AND IN "USER-FRIENDLINESS", BUT WITH NO CHANGE IN C NUMERICAL OUTPUT. C C V. ARP, APRIL 5, 1994. C C****************************************************** ALAPP(PA) DOUBLE PRECISION FUNCTION F2A04(PA) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA GA/9.80665D00/, PMIN/1.47514D00/, PMAX/2.164543D05/ IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 RHOL = 1.0D0/F49A04(PA) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 RHOV = 1.0D0/F50A04(PA) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 SIG = F31A04(PA) IF (SIG.LT.-1.E18) GO TO 920 IF (SIG.LT.-1.E08) GO TO 910 F2A04 = DSQRT(SIG/(GA*(RHOL-RHOV))) RETURN 910 F2A04 = -1.E10 RETURN 920 F2A04 = -1.E20 RETURN END C****************************************************** ALAPT(T) DOUBLE PRECISION FUNCTION F3A04(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.13001D0/, TMIN/0.8D0/, GA/9.80665D0/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 RHOL = 1.0D0/F53A04(TK) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 RHOV = 1.0D0/F54A04(TK) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 SIG = F32A04(TK) IF (SIG.LT.-1.E18) GO TO 920 IF (SIG.LT.-1.E08) GO TO 910 F3A04 = DSQRT(SIG/(GA*(RHOL-RHOV))) RETURN 910 F3A04 = -1.0E10 RETURN 920 F3A04 = -1.0E20 RETURN END C****************************************************** ALHP(PA) DOUBLE PRECISION FUNCTION F4A04(PA) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D0/, PMAX/2.164543D05/ IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 HL = F23A04(PA) IF (HL.LT.-1.E18) GO TO 920 IF (HL.LT.-1.E08) GO TO 910 HV = F24A04(PA) IF (HV.LT.-1.E18) GO TO 920 IF (HV.LT.-1.E08) GO TO 910 F4A04 = HV-HL RETURN 910 F4A04 = -1.0E10 RETURN 920 F4A04 = -1.0E20 RETURN END C****************************************************** ALHT(T) DOUBLE PRECISION FUNCTION F5A04(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D0/, TMAX/5.13001D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 HL = F27A04(TK) IF (HL.LT.-1.E18) GO TO 920 IF (HL.LT.-1.E08) GO TO 910 HV = F28A04(TK) IF (HV.LT.-1.E18) GO TO 920 IF (HV.LT.-1.E08) GO TO 910 F5A04=HV-HL RETURN 910 F5A04 = -1.0E10 RETURN 920 F5A04 = -1.0E20 RETURN END C ***************************************************** ALMPD(PA) DOUBLE PRECISION FUNCTION F6A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 4.70436D04, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) RHOL = 1.0D0/F53A04(T) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 DPT = G22A04(RHOL, T) IF (DPT.LT.-1.E08) GO TO 920 DPR = G23A04(RHOL, T) IF (DPR.LT.-1.E08) GO TO 920 CALL S6A04(RHOL, T, DPT, DPR, RLAML) F6A04 = RLAML RETURN 910 F6A04 = -1.E10 RETURN 920 F6A04 = -1.E20 RETURN END C ***************************************************** ALMPDD(PA) DOUBLE PRECISION FUNCTION F7A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 4.70436D04, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) RHOV = 1.0D0/F54A04(T) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 DPT = G22A04(RHOV, T) IF (DPT.LT.-1.E08) GO TO 920 DPR = G23A04(RHOV, T) IF (DPR.LT.-1.E08) GO TO 920 CALL S6A04(RHOV, T, DPT, DPR, RLAMV) F7A04 = RLAMV RETURN 910 F7A04 = -1.E10 RETURN 920 F7A04 = -1.E20 RETURN END C ***************************************************** ALMPT(PA,T) DOUBLE PRECISION FUNCTION F8A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 3.5D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 RHOG = 1.0D0/F51A04(PA, T) IF (RHOG.LT.-1.E18) GO TO 920 IF (RHOG.LT.-1.E08) GO TO 910 DPT = G22A04(RHOG, T) IF (DPT.LT.-1.E08) GO TO 920 DPR = G23A04(RHOG, T) IF (DPR.LT.-1.E08) GO TO 920 CALL S6A04(RHOG, T, DPT, DPR, RLAMG) F8A04 = RLAMG RETURN 910 F8A04 = -1.E10 RETURN 920 F8A04 = -1.E20 RETURN END C ***************************************************** ALMTD(T) DOUBLE PRECISION FUNCTION F9A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 3.5D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 RHOL = 1.0D0/F53A04(T) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 DPT = G22A04(RHOL, T) IF (DPT.LT.-1.E08) GO TO 920 DPR = G23A04(RHOL, T) IF (DPR.LT.-1.E08) GO TO 920 CALL S6A04(RHOL, T, DPT, DPR, RLAML) F9A04 = RLAML RETURN 910 F9A04 = -1.E10 RETURN 920 F9A04 = -1.E20 RETURN END C ***************************************************** ALMTDD(T) DOUBLE PRECISION FUNCTION F10A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 3.5D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 RHOV = 1.0D0/F54A04(T) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 DPT = G22A04(RHOV, T) IF (DPT.LT.-1.E08) GO TO 920 DPR = G23A04(RHOV, T) IF (DPR.LT.-1.E08) GO TO 920 CALL S6A04(RHOV, T, DPT, DPR, RLAMV) F10A04 = RLAMV RETURN 910 F10A04 = -1.E10 RETURN 920 F10A04 = -1.E20 RETURN END C ***************************************************** AMUPD(PA) DOUBLE PRECISION FUNCTION F11A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX/ 4.70436D04, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) RHOL = 1.0D0/F53A04(T) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 CALL S7A04(RHOL, T, RMUL) F11A04 = RMUL RETURN 910 F11A04 = -1.0E10 RETURN 920 F11A04 = -1.0E20 RETURN END C **************************************************** AMUPDD(PA) DOUBLE PRECISION FUNCTION F12A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX/ 4.70436D04, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) RHOV = 1.0D0/F54A04(T) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 CALL S7A04(RHOV, T, RMUV) F12A04 = RMUV RETURN 910 F12A04 = -1.0E10 RETURN 920 F12A04 = -1.0E20 RETURN END C****************************************************** AMUPT(PA,T) DOUBLE PRECISION FUNCTION F13A04(PA, T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 3.5D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 RHOG = 1.0D0/F51A04(PA, T) IF (RHOG.LT.-1.E18) GO TO 920 IF (RHOG.LT.-1.E08) GO TO 910 CALL S7A04(RHOG, T, RMUG) F13A04 = RMUG RETURN 910 F13A04 = -1.0E10 RETURN 920 F13A04 = -1.0E20 RETURN END C****************************************************** AMUTD(PA) DOUBLE PRECISION FUNCTION F14A04(T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN, TMAX / 3.5D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 RHOL = 1.0D0/F53A04(T) IF (RHOL.LT.-1.E18) GO TO 920 IF (RHOL.LT.-1.E08) GO TO 910 CALL S7A04(RHOL, T, RMUL) F14A04 = RMUL RETURN 910 F14A04 = -1.0E10 RETURN 920 F14A04 = -1.0E20 RETURN END C****************************************************** AMUTDD(T) DOUBLE PRECISION FUNCTION F15A04(T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN, TMAX / 3.5D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 RHOV = 1.0D0/F54A04(T) IF (RHOV.LT.-1.E18) GO TO 920 IF (RHOV.LT.-1.E08) GO TO 910 CALL S7A04(RHOV, T, RMUV) F15A04 = RMUV RETURN 910 F15A04 = -1.0E10 RETURN 920 F15A04 = -1.0E20 RETURN END C ***************************************************** CPPD(PA) DOUBLE PRECISION FUNCTION F16A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / DD2 = 5000. TA = 2.0 TM = 3.0 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) IF (T.GT.2.3) DD2=G25A04(PA,T) IF (DD2.LT.-1.0E08) GO TO 910 CP2 = G14A04(DD2, T) IF (T.LT.2.98.AND.DD2.GT.140.) THEN DD1=G10A04(PA,T) IF (DD1.LT.-1.0E08) GO TO 910 CP1 = G17A04(DD1, T) TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD1.GE.140.AND.DD1.LE.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F16A04=CP1*((TA-T)/(TA-TM))+CP2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F16A04 = CP2 RETURN ELSEIF (T.LT.TM) THEN F16A04 = CP1 RETURN ENDIF ENDIF ELSE F16A04 = CP2 RETURN ENDIF 910 F16A04 = -1.E10 RETURN 920 F16A04 = -1.E20 RETURN END C ***************************************************** CPPDD(PA) DOUBLE PRECISION FUNCTION F17A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) DD2=G29A04(PA) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 F17A04 = G14A04(DD2, T) RETURN 910 F17A04 = -1.E10 RETURN 920 F17A04 = -1.E20 RETURN END C ***************************************************** CPPT(PA,T) DOUBLE PRECISION FUNCTION F18A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 0.79999D00 / DATA PC/ 2.2746D05/, DLT/ 1.0D-5/ DD1 = 10. DD2 = 5000. TA = 2.0 TM = 3.0 IF (PA.LE.PC) THEN PS = F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 250 ENDIF IF (PA.GT.PS.AND.T.LT.2.35) GO TO 55 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (T.GT.13.89429D0) GO TO 50 IF (PA.GE.5041.8D00.AND.PA.LE.30.134D05) THEN TLMD = F67A04(PA) IF (DABS((T-TLMD)/TLMD).LT.4.8D-4) T = TLMD+0.001D00 ENDIF IF (PA.GE.25.328D05) THEN PML = F68A04(T) IF (DABS((PA-PML)/PML).LE.DLT) GO TO 50 IF (PA.GT.PML) GO TO 920 ENDIF 50 IF (PA.GT.PC) PS = 1.474 IF (T.GT.2.3) DD2 = G11A04(PA,T) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 CP2 = G14A04(DD2,T) IF (T.LT.2.98.AND.DD2.GT.140.AND.DD2.LT.189.) THEN IF (PA.GT.PS) THEN DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 CP1 = G17A04(DD1,T) ENDIF TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD2.GT.140.AND.DD2.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F18A04 = CP1*((TA-T)/(TA-TM))+CP2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F18A04 = CP2 RETURN ELSEIF (T.LT.TM) THEN F18A04 = CP1 RETURN ENDIF ENDIF ELSE F18A04 = CP2 RETURN ENDIF 55 DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 F18A04 = G17A04(DD1,T) RETURN 250 F18A04 = F17A04(PA) RETURN 910 F18A04 = -1.E10 RETURN 920 F18A04 = -1.E20 RETURN END C ***************************************************** CPTD(T) DOUBLE PRECISION FUNCTION F19A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) CP = F16A04(PA) IF (CP.LT.-1.0E18) GO TO 920 IF (CP.LT.-1.0E08) GO TO 910 F19A04 = CP RETURN 910 F19A04 = -1.E10 RETURN 920 F19A04 = -1.E20 RETURN END C ***************************************************** CPTDD(T) DOUBLE PRECISION FUNCTION F20A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) CPV = F17A04(PA) IF (CPV.LT.-1.0E18) GO TO 920 IF (CPV.LT.-1.0E08) GO TO 910 F20A04 = CPV RETURN 910 F20A04 = -1.E10 RETURN 920 F20A04 = -1.E20 RETURN END C ***************************************************** HPD(PA) DOUBLE PRECISION FUNCTION F23A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) VV = F49A04(PA) IF (VV.LT.-1.0E18) GO TO 920 IF (VV.LT.-1.0E08) GO TO 910 UU = F42A04(PA) F23A04 = UU+PA*VV RETURN 910 F23A04 = -1.E10 RETURN 920 F23A04 = -1.E20 RETURN END C ***************************************************** HPDD(PA) DOUBLE PRECISION FUNCTION F24A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 VV=F50A04(PA) IF (VV.LT.-1.0E18) GO TO 920 IF (VV.LT.-1.0E08) GO TO 910 UU = F43A04(PA) F24A04 = UU+PA*VV RETURN 910 F24A04 = -1.E10 RETURN 920 F24A04 = -1.E20 RETURN END C ***************************************************** HPT(PA,T) DOUBLE PRECISION FUNCTION F25A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 0.79999D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 UU = F44A04(PA, T) IF (UU.LT.-1.0E18) GO TO 920 IF (UU.LT.-1.0E08) GO TO 910 VV = F51A04(PA, T) IF (VV.LT.-1.0E18) GO TO 920 IF (VV.LT.-1.0E08) GO TO 910 F25A04 = UU+PA*VV RETURN 910 F25A04 = -1.E10 RETURN 920 F25A04 = -1.E20 RETURN END C****************************************************** HPX(P,X) DOUBLE PRECISION FUNCTION F26A04(PA,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / DATA ER10/-1.0D08/, ER20/-1.0D18/ IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 HL = F23A04(PA) IF(HL.LT.ER20) GO TO 920 IF(HL.LT.ER10) GO TO 910 HV = F24A04(PA) IF(HV.LT.ER20) GO TO 920 IF(HV.LT.ER10) GO TO 910 F26A04 = HL+X*(HV-HL) RETURN 910 F26A04 = -1.E+10 RETURN 920 F26A04 = -1.E+20 END C ***************************************************** HTD(T) DOUBLE PRECISION FUNCTION F27A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) HHL = F23A04(PA) IF (HHL.LT.-1.0E18) GO TO 920 IF (HHL.LT.-1.0E08) GO TO 910 F27A04 = HHL RETURN 910 F27A04 = -1.E10 RETURN 920 F27A04 = -1.E20 RETURN END C ***************************************************** HTDD(T) DOUBLE PRECISION FUNCTION F28A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) IF (PA.LT.-1.0E08) GO TO 910 HHV=F24A04(PA) IF (HHV.LT.-1.0E18) GO TO 920 IF (HHV.LT.-1.0E08) GO TO 910 F28A04 = HHV RETURN 910 F28A04 = -1.E10 RETURN 920 F28A04 = -1.E20 RETURN END C****************************************************** HTX(T,X) DOUBLE PRECISION FUNCTION F29A04(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMAX/5.13001D00/, TMIN/0.79999D00/, ER10/-1.D08/ 1 ER20/-1.0D18/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 HL = F27A04(TK) IF(HL.LT.ER20) GO TO 920 IF(HL.LT.ER10) GO TO 910 HV = F28A04(TK) IF(HV.LT.ER20) GO TO 920 IF(HV.LT.ER10) GO TO 910 F29A04 = HL+X*(HV-HL) RETURN 910 F29A04 = -1.E+10 RETURN 920 F29A04 = -1.E+20 END C ***************************************************** SIGP(PA) DOUBLE PRECISION FUNCTION F31A04(PA) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TC/5.15D00/, PMIN/1.47514D00/, PMAX/2.164543D05/, & TM/3.15D00/ IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 TA = F40A04(PA) CN = 0.5428D00 IF (TA.GT.TM) CN = 1.0638D00 F31A04 = 0.239D-3*((TC-TA)/(TC-TM))**CN RETURN 920 F31A04 = -1.E20 RETURN END C ***************************************************** SIGT(T) DOUBLE PRECISION FUNCTION F32A04(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.8D00/,TM/3.15D00/,TC/5.15D00/,TMAX/5.13001D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 900 IF (DABS((TK-TC)/TC).LE.1.D-05) THEN F32A04 = 0. RETURN ENDIF CN = 0.5428D00 IF (TK.GT.TM) CN = 1.0638D00 F32A04 = 0.239D-3*((TC-TK)/(TC-TM))**CN RETURN 900 F32A04 = -1.E+20 RETURN END C ***************************************************** SPD(PA) DOUBLE PRECISION FUNCTION F33A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / DD2 = 5000. TA = 2.0 TM = 3.0 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) IF (T.GT.2.3) DD2=G25A04(PA,T) IF (DD2.LT.-1.0E08) GO TO 910 SS2 = G13A04(DD2, T) IF (T.LT.2.98.AND.DD2.GT.140.) THEN DD1=G10A04(PA,T) IF (DD1.LT.-1.0E08) GO TO 910 SS1 = G6A04(DD1, T) TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD1.GE.140.AND.DD1.LE.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F33A04=SS1*((TA-T)/(TA-TM))+SS2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F33A04 = SS2 RETURN ELSEIF (T.LT.TM) THEN F33A04 = SS1 RETURN ENDIF ENDIF ELSE F33A04=SS2 RETURN ENDIF 910 F33A04 = -1.E10 RETURN 920 F33A04 = -1.E20 RETURN END C ***************************************************** SPDD(PA) DOUBLE PRECISION FUNCTION F34A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) DD2 = G29A04(PA) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 F34A04 = G13A04(DD2, T) RETURN 910 F34A04 = -1.E10 RETURN 920 F34A04 = -1.E20 RETURN END C ***************************************************** SPT(PA,T) DOUBLE PRECISION FUNCTION F35A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 0.79999D00 / DATA PC/ 2.2746D05/, TC/ 5.1953D00/, DLT/ 1.0D-5/ DD1 = 10. DD2 = 5000. DPC = DABS((PA-PC)/PC) DTC = DABS((T-TC)/TC) IF (DPC.LE.DLT.AND.DTC.LE.DLT) THEN F35A04 = 5.7685D03 RETURN ENDIF TA = 2.0 TM = 3.0 IF (PA.LE.PC) THEN PS = F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 250 ENDIF IF (PA.GT.PS.AND.T.LT.2.35) GO TO 55 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (T.GT.13.89429D0) GO TO 50 IF (PA.GE.5041.8D00.AND.PA.LE.30.134D05) THEN TLMD = F67A04(PA) IF (DABS((T-TLMD)/TLMD).LT.4.8D-4) T = TLMD+0.001D00 ENDIF IF (PA.GE.25.328D05) THEN PML = F68A04(T) IF (DABS((PA-PML)/PML).LE.DLT) GO TO 50 IF (PA.GT.PML) GO TO 920 ENDIF 50 IF (PA.GT.PC) PS = 1.474 IF (T.GT.2.3) DD2 = G11A04(PA,T) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 SS2 = G13A04(DD2,T) IF (T.LT.2.98.AND.DD2.GT.140.AND.DD2.LT.189.) THEN IF (PA.GT.PS) THEN DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 SS1 = G6A04(DD1,T) ENDIF TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD2.GT.140.AND.DD2.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F35A04 = SS1*((TA-T)/(TA-TM))+SS2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F35A04 = SS2 RETURN ELSEIF (T.LT.TM) THEN F35A04 = SS1 RETURN ENDIF ENDIF ELSE F35A04 = SS2 RETURN ENDIF 55 DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 F35A04 = G6A04(DD1,T) RETURN 250 F35A04 = F34A04(PA) RETURN 910 F35A04 = -1.E10 RETURN 920 F35A04 = -1.E20 RETURN END C****************************************************** SPX(P,X) DOUBLE PRECISION FUNCTION F36A04(PA,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMAX/2.164543D05/, PMIN/1.47514D00/, ER10/-1.D08/, & ER20/-1.D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 SL = F33A04(PA) IF(SL.LT.ER20) GO TO 920 IF(SL.LT.ER10) GO TO 910 SV = F34A04(PA) IF(SV.LT.ER20) GO TO 920 IF(SV.LT.ER10) GO TO 910 F36A04 = SL+X*(SV-SL) RETURN 910 F36A04 = -1.E+10 RETURN 920 F36A04 = -1.E+20 END C ***************************************************** STD(T) DOUBLE PRECISION FUNCTION F37A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) SS = F33A04(PA) IF (SS.LT.-1.0E18) GO TO 920 IF (SS.LT.-1.0E08) GO TO 910 F37A04 = SS RETURN 910 F37A04 = -1.E10 RETURN 920 F37A04 = -1.E20 RETURN END C ***************************************************** STDD(T) DOUBLE PRECISION FUNCTION F38A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) SPV = F34A04(PA) IF (SPV.LT.-1.0E18) GO TO 920 IF (SPV.LT.-1.0E08) GO TO 910 F38A04 = SPV RETURN 910 F38A04 = -1.E10 RETURN 920 F38A04 = -1.E20 RETURN END C****************************************************** STX(T,X) DOUBLE PRECISION FUNCTION F39A04(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.13001D00/, ER10/-1.D08/, 1 ER20/-1.D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 SL = F37A04(TK) IF(SL.LT.ER20) GO TO 920 IF(SL.LT.ER10) GO TO 910 SV = F38A04(TK) IF(SV.LT.ER20) GO TO 920 IF(SV.LT.ER10) GO TO 910 F39A04 = SL+X*(SV-SL) RETURN 910 F39A04 = -1.E10 RETURN 920 F39A04 = -1.E20 END C ***************************************************** UPD(PA) DOUBLE PRECISION FUNCTION F42A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / DD2 = 5000. TA = 2.0 TM = 3.0 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 50 PS=1.474 T = F40A04(PA) IF (T.GT.2.3) DD2=G25A04(PA,T) IF (DD2.LT.-1.0E08) GO TO 910 UU2 = G12A04(DD2, T) IF (T.LT.2.98.AND.DD2.GT.140.) THEN DD1=G10A04(PA,T) IF (DD1.LT.-1.0E08) GO TO 910 UU1 = G15A04(DD1, T) TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD1.GE.140.AND.DD1.LE.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F42A04=UU1*((TA-T)/(TA-TM))+UU2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F42A04 = UU2 RETURN ELSEIF (T.LT.TM) THEN F42A04 = UU1 RETURN ENDIF ENDIF ELSE F42A04=UU2 RETURN ENDIF 910 F42A04 = -1.E10 RETURN 920 F42A04 = -1.E20 RETURN END C ***************************************************** UPDD(PA) DOUBLE PRECISION FUNCTION F43A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) DD2=G29A04(PA) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 F43A04 = G12A04(DD2, T) RETURN 910 F43A04 = -1.E10 RETURN 920 F43A04 = -1.E20 RETURN END C ***************************************************** UPT(PA,T) DOUBLE PRECISION FUNCTION F44A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0001D00, 0.79999D00 / DATA PC/ 2.2746D05/, TC/ 5.1953D00/, DLT/ 1.0D-5/ DD1 = 10. DD2 = 5000. DPC = DABS((PA-PC)/PC) DTC = DABS((T-TC)/TC) IF (DPC.LE.DLT.AND.DTC.LE.DLT) THEN F44A04 = 18.682D03 RETURN ENDIF TA = 2.0 TM = 3.0 IF (PA.LE.PC) THEN PS = F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 250 ENDIF IF (PA.GT.PS.AND.T.LT.2.35) GO TO 55 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (T.GT.13.89429D0) GO TO 50 IF (PA.GE.5041.8D00.AND.PA.LE.30.134D05) THEN TLMD = F67A04(PA) IF (DABS((T-TLMD)/TLMD).LT.4.8D-4) T = TLMD+0.001D00 ENDIF IF (PA.GE.25.328D05) THEN PML = F68A04(T) IF (DABS((PA-PML)/PML).LE.DLT) GO TO 50 IF (PA.GT.PML) GO TO 920 ENDIF 50 IF (PA.GT.PC) PS = 1.474 IF (T.GT.2.3) DD2 = G11A04(PA,T) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 UU2 = G12A04(DD2, T) IF (T.LT.2.98.AND.DD2.GT.140.AND.DD2.LT.189.) THEN IF (PA.GT.PS) THEN DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 UU1 = G15A04(DD1, T) ENDIF TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD2.GT.140.AND.DD2.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F44A04 = UU1*((TA-T)/(TA-TM))+UU2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F44A04 = UU2 RETURN ELSEIF (T.LT.TM) THEN F44A04 = UU1 RETURN ENDIF ENDIF ELSE F44A04 = UU2 RETURN ENDIF 55 DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 F44A04 = G15A04(DD1, T) RETURN 250 F44A04 = F43A04(PA) RETURN 910 F44A04 = -1.E10 RETURN 920 F44A04 = -1.E20 RETURN END C****************************************************** UPX(PA,X) DOUBLE PRECISION FUNCTION F45A04(PA,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PMAX/2.164543D05/, ER10/-1.D08/, & ER20/-1.D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 UL = F42A04(PA) IF(UL.LT.ER20) GO TO 920 IF(UL.LT.ER10) GO TO 910 UV = F43A04(PA) IF(UV.LT.ER20) GO TO 920 IF(UV.LT.ER10) GO TO 910 F45A04 = UL+X*(UV-UL) RETURN 910 F45A04 = -1.E+10 RETURN 920 F45A04 = -1.E+20 END C ***************************************************** UTD(T) DOUBLE PRECISION FUNCTION F46A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) UUD = F42A04(PA) IF (UUD.LT.-1.0E18) GO TO 920 IF (UUD.LT.-1.0E08) GO TO 910 F46A04 = UUD RETURN 910 F46A04 = -1.E10 RETURN 920 F46A04 = -1.E20 RETURN END C ***************************************************** UTDD(T) DOUBLE PRECISION FUNCTION F47A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) IF (PA.LT.-1.0E18) GO TO 920 IF (PA.LT.-1.0E08) GO TO 910 UUV = F43A04(PA) IF (UUV.LT.-1.0E18) GO TO 920 IF (UUV.LT.-1.0E08) GO TO 910 F47A04 = UUV RETURN 910 F47A04 = -1.E10 RETURN 920 F47A04 = -1.E20 RETURN END C****************************************************** UTX(T,X) DOUBLE PRECISION FUNCTION F48A04(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.13001D00/, ER10/-1.D08/, 1 ER20/-1.D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 UL = F46A04(TK) IF(UL.LT.ER20) GO TO 920 IF(UL.LT.ER10) GO TO 910 UV = F47A04(TK) IF(UV.LT.ER20) GO TO 920 IF(UV.LT.ER10) GO TO 910 F48A04 = UL+X*(UV-UL) RETURN 910 F48A04 = -1.E10 RETURN 920 F48A04 = -1.E20 END C ***************************************************** VPD(PA) DOUBLE PRECISION FUNCTION F49A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX, PS51/ 1.47514D00, 2.2746D05, 2.11583D05 / DD2 = 5000. TA = 2.0 TM = 3.0 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (PA.LT.PS51) GO TO 50 T = F40A04(PA) IF (T.GT.5.17D00) GO TO 20 RHOL = -0.250905513690478D05 & +T*( 0.994417261904769D04 -0.981547619047626D03*T) GO TO 30 20 RHOL = -0.835862524846021D06 & +T*( 0.323205980838379D06 -0.312406145324988D05*T) 30 F49A04 = 1.0D0/RHOL RETURN 50 PS=1.474 T = F40A04(PA) IF (T.GT.2.3) DD2=G25A04(PA,T) IF (DD2.LT.-1.0E08) GO TO 910 IF (T.LT.2.98.AND.DD2.GT.140.) THEN DD1=G10A04(PA,T) IF (DD1.LT.-1.0E08) GO TO 910 TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD1.GE.140.AND.DD1.LE.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F49A04=1.0D00/(DD1*((TA-T)/(TA-TM))+DD2*((T-TM)/(TA-TM))) RETURN ELSEIF (T.GT.TA) THEN F49A04 = 1.0D00/DD2 RETURN ELSEIF (T.LT.TM) THEN F49A04 = 1.0D00/DD1 RETURN ENDIF ENDIF ELSE F49A04 = 1.0D00/DD2 RETURN ENDIF 910 F49A04 = -1.E10 RETURN 920 F49A04 = -1.E20 RETURN END C ***************************************************** VPDD(PA) DOUBLE PRECISION FUNCTION F50A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX, PS51/ 1.47514D00, 2.2746D05, 2.11583D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (PA.LT.PS51) GO TO 50 T = F40A04(PA) IF (T.GT.5.17D0) GO TO 20 RHOV = 0.169123666666669D05 1 +T*( -0.670666666666677D04 +0.666666666666676D03*T) GO TO 30 20 RHOV = 0.515417476016971D06 1 +T*( -0.199331153501323D06 +0.192743717408819D05*T) 30 F50A04 = 1.0D0/RHOV RETURN 50 DD2=G29A04(PA) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 F50A04 = 1.0D00/DD2 RETURN 910 F50A04 = -1.E10 RETURN 920 F50A04 = -1.E20 RETURN END C ***************************************************** VPT(PA,T) DOUBLE PRECISION FUNCTION F51A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.79999D00, 1500.0001D00/ DATA PMAX, PMIN/ 100.0001D06, 1.47514D00/ DATA PC/ 2.2746D05/, TC/ 5.1953D00/, DLT/ 1.0D-5/ DD1 = 10. DD2 = 5000. DPC = DABS((PA-PC)/PC) DTC = DABS((T-TC)/TC) IF (DPC.LE.DLT.AND.DTC.LE.DLT) THEN F51A04 = 14.360D-3 RETURN ENDIF TA = 2.4 TM = 3.0 IF (PA.LE.PC) THEN PS=F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 250 ENDIF IF (PA.GT.PS.AND.T.LT.2.35) GO TO 55 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (T.GT.13.89429D0) GO TO 50 IF (PA.GE.5041.8D00.AND.PA.LE.30.134D05) THEN TLMD = F67A04(PA) IF (DABS((T-TLMD)/TLMD).LT.4.8D-4) T = TLMD+0.001D00 ENDIF IF (PA.GE.25.328D05) THEN PML=F68A04(T) IF (DABS((PA-PML)/PML).LE.DLT) GO TO 50 IF (PA.GT.PML) GO TO 920 ENDIF 50 IF (PA.GT.PC) PS=1.474 IF (T.GT.2.3) DD2=G11A04(PA,T) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 IF (T.LT.2.98.AND.DD2.GT.140.AND.DD2.LT.189.) THEN IF (PA.GT.PS) THEN DD1=G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 ENDIF TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD2.GT.140.AND.DD2.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F51A04 =1.0D00/(DD1*((TA-T)/(TA-TM))+DD2*((T-TM)/(TA-TM))) RETURN ELSEIF (T.GT.TA) THEN F51A04 = 1.0D00/DD2 RETURN ELSEIF (T.LT.TM) THEN F51A04 = 1.0D00/DD1 RETURN ENDIF ENDIF ELSE F51A04 =1.0D00/DD2 RETURN ENDIF 55 DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 F51A04 = 1.D00/DD1 RETURN 250 F51A04 = 1.0D0/64.64D00 RETURN 910 F51A04 = -1.E10 RETURN 920 F51A04 = -1.E20 RETURN END C****************************************************** VPX(P,X) DOUBLE PRECISION FUNCTION F52A04(PA,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PC/2.2746D05/, ER10/-1.D08/, ER20/-1.D18/ IF(PA.LT.PMIN.OR.PA.GT.PC) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 VL = F49A04(PA) IF(VL.LT.ER20) GO TO 920 IF(VL.LT.ER10) GO TO 910 VV = F50A04(PA) IF(VV.LT.ER20) GO TO 920 IF(VV.LT.ER10) GO TO 910 F52A04 = VL+X*(VV-VL) RETURN 910 F52A04 = -1.E10 RETURN 920 F52A04 = -1.E20 END C ***************************************************** VTD(T) DOUBLE PRECISION FUNCTION F53A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.1953D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) VV = F49A04(PA) IF (VV.LT.-1.0E18) GO TO 920 IF (VV.LT.-1.0E08) GO TO 910 F53A04 = VV RETURN 910 F53A04 = -1.E10 RETURN 920 F53A04 = -1.E20 RETURN END C ***************************************************** VTDD(T) DOUBLE PRECISION FUNCTION F54A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.1953D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) IF (PA.LT.-1.0E08) GO TO 910 VV2=F50A04(PA) IF (VV2.LT.-1.0E18) GO TO 920 IF (VV2.LT.-1.0E08) GO TO 910 F54A04 = VV2 RETURN 910 F54A04 = -1.E10 RETURN 920 F54A04 = -1.E20 RETURN END C****************************************************** VTX(T,X) DOUBLE PRECISION FUNCTION F55A04(TK,X) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.1953D00/, ER10/-1.D08/, 1 ER20/-1.D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 IF(X.LT.0.OR.X.GT.1.D00) GO TO 920 VL = F53A04(TK) IF(VL.LT.ER20) GO TO 920 IF(VL.LT.ER10) GO TO 910 VV = F54A04(TK) IF(VV.LT.ER20) GO TO 920 IF(VV.LT.ER10) GO TO 910 F55A04 = VL+X*(VV-VL) RETURN 910 F55A04 = -1.E10 RETURN 920 F55A04 = -1.E20 RETURN END C****************************************************** XPH(P,H) DOUBLE PRECISION FUNCTION F56A04(PA,H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PMAX/2.164543D05/, ER10/-1.0D08/, & ER20/-1.0D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 HL = F23A04(PA) IF(HL.LT.ER20) GO TO 920 IF(HL.LT.ER10) GO TO 910 HV = F24A04(PA) IF(HV.LT.ER20) GO TO 920 IF(HV.LT.ER10) GO TO 910 IF(H.LT.HL.OR.H.GT.HV) GO TO 920 F56A04 = (H-HL)/(HV-HL) RETURN 910 F56A04 = -1.E10 RETURN 920 F56A04 = -1.E20 RETURN END C****************************************************** XPS(P,S) DOUBLE PRECISION FUNCTION F57A04(PA,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PMAX/2.164543D05/, ER10/-1.0D08/, & ER20/-1.0D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 SL = F33A04(PA) IF(SL.LT.ER20) GO TO 920 IF(SL.LT.ER10) GO TO 910 SV = F34A04(PA) IF(SV.LT.ER20) GO TO 920 IF(SV.LT.ER10) GO TO 910 IF(S.LT.SL.OR.S.GT.SV) GO TO 920 F57A04 = (S-SL)/(SV-SL) RETURN 910 F57A04 = -1.E10 RETURN 920 F57A04 = -1.E20 END C****************************************************** XPU(P,U) DOUBLE PRECISION FUNCTION F58A04(PA,U) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PMAX/2.164543D05/, ER10/-1.0D08/, & ER20/-1.0D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 UL = F42A04(PA) IF(UL.LT.ER20) GO TO 920 IF(UL.LT.ER10) GO TO 910 UV = F43A04(PA) IF(UV.LT.ER20) GO TO 920 IF(UV.LT.ER10) GO TO 910 IF(U.LT.UL.OR.U.GT.UV) GO TO 920 F58A04 = (U-UL)/(UV-UL) RETURN 910 F58A04 = -1.E10 RETURN 920 F58A04 = -1.E20 RETURN END C****************************************************** XPV(P,V) DOUBLE PRECISION FUNCTION F59A04(PA,V) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PMAX/2.2746D05/, ER10/-1.0D08/, & ER20/-1.0D18/ IF(PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 VL = F49A04(PA) IF(VL.LT.ER20) GO TO 920 IF(VL.LT.ER10) GO TO 910 VV = F50A04(PA) IF(VV.LT.ER20) GO TO 920 IF(VV.LT.ER10) GO TO 910 IF(V.LT.VL.OR.V.GT.VV) GO TO 920 F59A04 = (V-VL)/(VV-VL) RETURN 910 F59A04 = -1.E10 RETURN 920 F59A04 = -1.E20 RETURN END C****************************************************** XTH(T,H) DOUBLE PRECISION FUNCTION F60A04(TK,H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.13001D00/, ER10/-1.0D08/, 1 ER20/-1.0D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 HL = F27A04(TK) IF(HL.LT.ER20) GO TO 920 IF(HL.LT.ER10) GO TO 910 HV = F28A04(TK) IF(HV.LT.ER20) GO TO 920 IF(HV.LT.ER10) GO TO 910 IF(H.LT.HL.OR.H.GT.HV) GO TO 920 F60A04 = (H-HL)/(HV-HL) RETURN 910 F60A04 = -1.E10 RETURN 920 F60A04 = -1.E20 RETURN END C****************************************************** XTS(T,S) DOUBLE PRECISION FUNCTION F61A04(TK,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.13001D00/, ER10/-1.0D08/, 1 ER20/-1.0D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 SL = F37A04(TK) IF(SL.LT.ER20) GO TO 920 IF(SL.LT.ER10) GO TO 910 SV = F38A04(TK) IF(SV.LT.ER20) GO TO 920 IF(SV.LT.ER10) GO TO 910 IF(S.LT.SL.OR.S.GT.SV) GO TO 920 F61A04 = (S-SL)/(SV-SL) RETURN 910 F61A04 = -1.E10 RETURN 920 F61A04 = -1.E20 RETURN END C****************************************************** XTU(T,U) DOUBLE PRECISION FUNCTION F62A04(TK,U) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.13001D00/, ER10/-1.0D08/, 1 ER20/-1.0D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 UL = F46A04(TK) IF(UL.LT.ER20) GO TO 920 IF(UL.LT.ER10) GO TO 910 UV = F47A04(TK) IF(UV.LT.ER20) GO TO 920 IF(UV.LT.ER10) GO TO 910 IF(U.LT.UL.OR.U.GT.UV) GO TO 920 F62A04 = (U-UL)/(UV-UL) RETURN 910 F62A04 = -1.E10 RETURN 920 F62A04 = -1.E20 END C****************************************************** XTV(P,V) DOUBLE PRECISION FUNCTION F63A04(TK,V) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA TMIN/0.79999D00/, TMAX/5.1953D00/, ER10/-1.0D08/, 1 ER20/-1.0D18/ IF(TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 VL = F53A04(TK) IF(VL.LT.ER20) GO TO 920 IF(VL.LT.ER10) GO TO 910 VV = F54A04(TK) IF(VV.LT.ER20) GO TO 920 IF(VV.LT.ER10) GO TO 910 IF(V.LT.VL.OR.V.GT.VV) GO TO 920 F63A04 = (V-VL)/(VV-VL) RETURN 910 F63A04 = -1.E10 RETURN 920 F63A04 = -1.E20 RETURN END C****************************************************** TPH(P,H) DOUBLE PRECISION FUNCTION F64A04(PA, H) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC1/2.164543D05/, TC/5.1953D00/, PMIN/1.47514D0/, 1 PP/100.D06/, TMAX/1500.D00/, DLT/1.D-05/ FNC(H1,H2) = DABS((H1-H2)/H2) IF (PA.LT.PMIN.OR.PA.GT.PP) GO TO 920 HMAX = F25A04(PA,TMAX) IF (FNC(H,HMAX).LE.DLT) GO TO 500 IF (PA.GT.25.328D05) THEN T1 = F69A04(PA) TMIN = T1 ELSE T1 = 0.8D00 TMIN = T1 ENDIF HMIN = F25A04(PA,T1) IF (PA.LT.0.06D05) THEN HL = F23A04(PA) HMIN = DMIN1(HL,HMIN) ENDIF IF (FNC(H,HMIN).LE.DLT) GO TO 500 IF (H.LT.HMIN.OR.H.GT.HMAX) GO TO 920 IF (PA.LE.PC1) GO TO 10 T0 = TC IF (PA.LT.70.D05) GO TO 8 TM = F69A04(PA) T0 = TM 8 IF (H.GT.2.5D06) T0 = TMAX H0 = F25A04(PA,T0) T1 = T0+(H-H0)/F18A04(PA,T0) IF (T1.LT.TMIN) T1 = TMIN H1 = F25A04(PA,T1) IF (H1.LT.-1.E18) GO TO 920 IF (H1.LT.-1.E08) GO TO 910 GO TO 100 10 HD = F23A04(PA) HDD = F24A04(PA) TS = F40A04(PA) IF (H.GE.HD.AND.H.LE.HDD) GO TO 30 IF (DABS((H-HD)/HD).LE.1.D-4.OR.DABS((H-HDD)/HDD).LE.1.D-4) & GO TO 30 IF (H.GE.HD) GO TO 20 T0 = TS H0 = HD T1 = TMIN H1 = HMIN IF (H1.LT.-1.E08) GO TO 910 GO TO 100 20 IF (H.LE.HDD) GO TO 30 T0 = TMAX H0 = F25A04(PA,T0) T1 = T0+(H-H0)/F18A04(PA,T0) IF (T1.LT.TMIN) T1 = TMIN H1 = F25A04(PA,T1) IF (H1.LT.-1.E18) GO TO 920 IF (H1.LT.-1.E08) GO TO 910 GO TO 100 30 F64A04 = TS RETURN 100 CALL S4A04 (PA,H,T0,T1,H0,H1,TMIN) IF (T1.LT.-1.0E08) GO TO 910 500 F64A04 = T1 RETURN 910 F64A04 = -1.E10 RETURN 920 F64A04 = -1.E20 RETURN END C****************************************************** TPS(P,S) DOUBLE PRECISION FUNCTION F65A04 (PA,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PC1/2.164543D05/, TC/5.1953D00/, PMIN/1.47514D00/, & PP/100.D06/, TMAX/1500.D00/, DLT/1.D-05/ FNC(S1,S2) = DABS((S1-S2)/S2) IF (PA.LT.PMIN.OR.PA.GT.PP) GO TO 920 SMAX = F35A04(PA,TMAX) IF (FNC(S,SMAX).LE.DLT) GO TO 500 IF (PA.GT.25.328D05) THEN T1 = F69A04(PA) TMIN = T1 ELSE T1 = 0.8D00 TMIN = T1 ENDIF SMIN = F35A04(PA,T1) IF (PA.LT.0.06D05) THEN SL = F33A04(PA) SMIN = DMIN1(SL,SMIN) ENDIF IF (FNC(S,SMIN).LE.DLT) GO TO 500 IF (S.LT.SMIN.OR.S.GT.SMAX) GO TO 920 IF (PA.LT.PC1) GO TO 10 T0 = TC IF (PA.LT.70.D05) GO TO 8 TM = F69A04(PA) T0 = TM 8 S0 = F35A04(PA,T0) T1 = T0+(S-S0)*T0/F18A04(PA,T0) IF (T1.LT.TMIN) T1 = TMIN S1 = F35A04(PA,T1) IF (S1.LT.-1.E18) GO TO 920 IF (S1.LT.-1.E08) GO TO 910 GO TO 100 10 SD = F33A04(PA) SDD = F34A04(PA) TS = F40A04(PA) IF (S.GE.SD.AND.S.LE.SDD) GO TO 30 IF (DABS((S-SD)/SD).LE.1.D-4.OR.DABS((S-SDD)/SDD).LE.1.D-4) & GO TO 30 IF (S.GE.SD) GO TO 20 T0 = TS S0 = SD T1 = TMIN S1 = SMIN IF (S1.LT.-1.E18) GO TO 920 IF (S1.LT.-1.E08) 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 F65A04 = TS RETURN 100 CALL S5A04 (PA,S,T0,T1,S0,S1,TMIN) IF (T1.LT.-1.E08) GO TO 910 500 F65A04 = T1 RETURN 910 F65HE = -1.E10 RETURN 920 F65HE = -1.E20 RETURN END C ***************************************************** TPV(PA,V) DOUBLE PRECISION FUNCTION F70A04(PA, V) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0D06, 5.0418D03/ DATA TMAX, DELT8 / 1500.0D00, 1.0D-10 / DATA PC, PML, DLT/ 2.2746D05, 30.134D05, 1.0D-5 / P0 = PA/1.0D06 DD = 1.0D00/V RHO = DD/4.0026D00 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (PA.GT.PML) THEN T1 = F69A04(PA) DDMAX = 1.0D00/F51A04(PA,T1) ELSE T1 = G19A04(PA) DDMAX = 1.0D00/F51A04(PA,T1) ENDIF DDMIN = 1.0D00/F51A04(PA,TMAX) IF (DABS((DD-DDMAX)/DD).LT.DLT) THEN T2 = T1 GO TO 200 ELSE IF (DABS((DD-DDMIN)/DD).LT.DLT) THEN T2 = TMAX GO TO 200 ENDIF IF (DD.GT.DDMAX.OR.DD.LT.DDMIN) GO TO 920 IF (PA.LT.PC) THEN DDL = 1.0D00/F49A04(PA) DDV = 1.0D00/F50A04(PA) IF (DABS((DD-DDL)/DDL).LE.DLT.OR.DABS((DD-DDV)/DDV).LE.DLT) & GO TO 300 IF (DD.LE.DDL.AND.DD.GE.DDV) GO TO 300 IF (DD.GT.DDL) THEN T2 = F40A04(PA) GO TO 55 ENDIF ENDIF IC = 0 TDLT = 100.0D00 T2 = 0.0 DO 50 I=1,14 T2 = T2+TDLT DDM = 1.0D00/F51A04(PA,T2) IF (DABS((DD-DDM)/DD).LE.DELT8) GO TO 200 IF (DDM.LE.DD) GO TO 55 T1 = T2 50 CONTINUE T2 = 1500.D00 55 F1 = G26A04(PA,DD,T1) 60 TM = (T1+T2)/2.0 F2 = G26A04(PA,DD,TM) IF (F1*F2.GT.0) THEN T1 = TM ELSE T2 = TM ENDIF IC = IC+1 IF (IC.GT.1000) GO TO 910 DLTT2 = DABS((T1-T2)/T2) IF (DLTT2.GT.DELT8) GO TO 60 200 F70A04 = T2 RETURN 300 F70A04 = F40A04(PA) RETURN 910 F70A04 = -1.E10 RETURN 920 F70A04 = -1.E20 RETURN END C****************************************************** HPS(P,S) DOUBLE PRECISION FUNCTION F71A04(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PP/100.0D06/, ER10/-1.D08/, ER20/-1.D18/ DATA PC1/2.164543D05/ IF (P.LT.PMIN.OR.P.GT.PP) GO TO 920 T = F65A04(P,S) IF (T.LT.ER20) GO TO 920 IF (T.LT.ER10) GO TO 910 IF (P.GE.PC1) GO TO 200 SL = F33A04(P) SV = F34A04(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 HL = F23A04(P) HV = F24A04(P) XX = F57A04(P,S) IF (XX.LT.ER20) GO TO 920 IF (XX.LT.ER10) GO TO 910 F71A04=HL+XX*(HV-HL) RETURN 200 F71A04 = F25A04(P, T) RETURN 910 F71A04 = -1.E10 RETURN 920 F71A04 = -1.E20 RETURN END C ***************************************************** CVPDD(PA) DOUBLE PRECISION FUNCTION F76A04(PA) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMIN, PMAX / 1.47514D00, 2.164543D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 T = F40A04(PA) DD2=G29A04(PA) IF (DD2.LT.-1.0E08) GO TO 910 F76A04 = G18A04(DD2, T) RETURN 910 F76A04 = -1.E10 RETURN 920 F76A04 = -1.E20 RETURN END C ***************************************************** CVPT(PA,T) DOUBLE PRECISION FUNCTION F77A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA PMAX, PMIN/ 100.0D06, 1.47514D00/ DATA TMAX, TMIN/ 1500.0D00, 0.79999D00 / DATA PC/ 2.2746D05/, TC/ 5.1953D00/, DLT/ 1.0D-5/ DD1 = 10. DD2 = 5000. DPC = DABS((PA-PC)/PC) DTC = DABS((T-TC)/TC) IF (DPC.LE.DLT.AND.DTC.LE.DLT) THEN F77A04 = 2.8882D03 RETURN ENDIF TA = 2.0 TM = 3.0 IF (PA.LE.PC) THEN PS = F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 250 ENDIF IF (PA.GT.PS.AND.T.LT.2.35) GO TO 55 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 IF (T.GT.13.89429D0) GO TO 50 IF (PA.GE.5041.8D00.AND.PA.LE.30.134D05) THEN TLMD = F67A04(PA) IF (DABS((T-TLMD)/TLMD).LT.4.8D-4) T = TLMD+0.001D00 ENDIF IF (PA.GE.25.328D05) THEN PML = F68A04(T) IF (DABS((PA-PML)/PML).LE.DLT) GO TO 50 IF (PA.GT.PML) GO TO 920 ENDIF 50 IF (PA.GT.PC) PS = 1.474 IF (T.GT.2.3) DD2 = G11A04(PA,T) IF (DD2.LT.-1.0E18) GO TO 920 IF (DD2.LT.-1.0E08) GO TO 910 CV2 = G18A04(DD2,T) IF (T.LT.2.98.AND.DD2.GT.140.AND.DD2.LT.189.) THEN IF (PA.GT.PS) THEN DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 CV1 = G5A04(DD1,T) ENDIF TA=2.98-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) TM=2.53-5.6D-3*(DD1-140.)-3.5D-2*(DMAX1(0.D00,DD1-180.D00)) IF (DD2.GT.140.AND.DD2.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN F77A04 = CV1*((TA-T)/(TA-TM))+CV2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN F77A04 = CV2 RETURN ELSEIF (T.LT.TM) THEN F77A04 = CV1 RETURN ENDIF ENDIF ELSE F77A04 = CV2 RETURN ENDIF 55 DD1 = G10A04(PA,T) IF (DD1.LT.-1.0E18) GO TO 920 IF (DD1.LT.-1.0E08) GO TO 910 F77A04 = G5A04(DD1,T) RETURN 250 F77A04 = F76A04(PA) RETURN 910 F77A04 = -1.E10 RETURN 920 F77A04 = -1.E20 RETURN END C ***************************************************** CVTDD(T) DOUBLE PRECISION FUNCTION F78A04(T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX / 0.79999D00, 5.13001D00 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PA = F30A04(T) CVV = F76A04(PA) IF (CVV.LT.-1.0E18) GO TO 920 IF (CVV.LT.-1.0E08) GO TO 910 F78A04 = CVV RETURN 910 F78A04 = -1.E10 RETURN 920 F78A04 = -1.E20 RETURN END C****************************************************** UPS(P,S) DOUBLE PRECISION FUNCTION F79A04(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D0/, PP/100.0D06/, ER10/-1.D08/, ER20/-1.D18/ DATA PC1/2.164543D05/ IF (P.LT.PMIN.OR.P.GT.PP) GO TO 920 T = F65A04(P,S) IF (T.LE.ER20) GO TO 920 IF (T.LE.ER10) GO TO 910 IF (P.GE.PC1) GO TO 200 SL = F33A04(P) SV = F34A04(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 XE = F57A04(P,S) IF (XE.LT.ER20) GO TO 920 IF (XE.LT.ER10) GO TO 910 UL = F42A04(P) UV = F43A04(P) F79A04 = UL+XE*(UV-UL) RETURN 200 F79A04 = F44A04(P,T) RETURN 910 F79A04 = -1.E10 RETURN 920 F79A04 = -1.E20 RETURN END C****************************************************** VPS(P,S) DOUBLE PRECISION FUNCTION F80A04(P,S) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA PMIN/1.47514D00/, PP/100.0D06/, ER10/-1.0D08/, ER20/-1.0D18/ DATA PC1/2.164543D05/ IF (P.LT.PMIN.OR.P.GT.PP) GO TO 920 T = F65A04(P,S) IF (T.LE.ER20) GO TO 920 IF (T.LE.ER10) GO TO 910 IF (P.GE.PC1) GO TO 200 SL = F33A04(P) SV = F34A04(P) IF (S.LT.SL.OR.S.GT.SV) GO TO 200 XE = F57A04(P,S) IF (XE.LT.ER20) GO TO 920 IF (XE.LT.ER10) GO TO 910 VL = F49A04(P) VV = F50A04(P) F80A04 = VL+XE*(VV-VL) RETURN 200 F80A04 = F51A04(P,T) RETURN 910 F80A04 = -1.0E10 RETURN 920 F80A04 = -1.0E20 RETURN END C****************************************************** PRPT(PA,T) DOUBLE PRECISION FUNCTION F81A04(PA,T) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ AMUG = F13A04(PA,T) IF (AMUG.LT.ER20) GO TO 920 IF (AMUG.LT.ER10) GO TO 910 ALMG = F8A04(PA,T) IF (ALMG.LT.ER20) GO TO 920 IF (ALMG.LT.ER10) GO TO 910 CPG = F18A04(PA,T) IF (CPG.LT.ER20) GO TO 920 IF (CPG.LT.ER10) GO TO 910 F81A04 = AMUG*CPG/ALMG RETURN 910 F81A04 = -1.0E10 RETURN 920 F81A04 = -1.0E20 RETURN END C ***************************************************** AKPT(PA,T) DOUBLE PRECISION FUNCTION F82A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ CPG = F18A04(PA,T) IF (CPG.LT.ER20) GO TO 920 IF (CPG.LT.ER10) GO TO 910 CVG = F77A04(PA,T) IF (CVG.LT.ER20) GO TO 920 IF (CVG.LT.ER10) GO TO 910 VV = F51A04(PA,T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 RHOG = 1.0D00/VV DPR = G23A04(RHOG, T) IF (DPR.LT.ER20) GO TO 920 IF (DPR.LT.ER10) GO TO 910 F82A04 = CPG*RHOG*DPR/(CVG*PA) RETURN 910 F82A04 = -1.0E10 RETURN 920 F82A04 = -1.0E20 RETURN END C ***************************************************** WPT(PA,T) DOUBLE PRECISION FUNCTION F83A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ VV = F51A04(PA, T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 AKP = F82A04(PA, T) IF (AKP.LT.ER20) GO TO 920 IF (AKP.LT.ER10) GO TO 910 F83A04 = DSQRT(PA*VV*AKP) RETURN 910 F83A04 = -1.0E10 RETURN 920 F83A04 = -1.0E20 RETURN END C****************************************************** PRPD(P) DOUBLE PRECISION FUNCTION F85A04(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ AMUL = F11A04(P) IF (AMUL.LT.ER20) GO TO 920 IF (AMUL.LT.ER10) GO TO 910 ALML = F6A04(P) IF (ALML.LT.ER20) GO TO 920 IF (ALML.LT.ER10) GO TO 910 CPL = F16A04(P) IF (CPL.LT.ER20) GO TO 920 IF (CPL.LT.ER10) GO TO 910 F85A04 = AMUL*CPL/ALML RETURN 910 F85A04 = -1.0E+10 RETURN 920 F85A04 = -1.0E+20 RETURN END C****************************************************** PRPDD(P) DOUBLE PRECISION FUNCTION F86A04(P) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ AMUV = F12A04(P) IF (AMUV.LT.ER20) GO TO 920 IF (AMUV.LT.ER10) GO TO 910 ALMV = F7A04(P) IF (ALMV.LT.ER20) GO TO 920 IF (ALMV.LT.ER10) GO TO 910 CPV = F17A04(P) IF (CPV.LT.ER20) GO TO 920 IF (CPV.LT.ER10) GO TO 910 F86A04 = AMUV*CPV/ALMV RETURN 910 F86A04 = -1.0E10 RETURN 920 F86A04 = -1.0E20 RETURN END C****************************************************** PRTD(T) DOUBLE PRECISION FUNCTION F87A04(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ AMUL = F14A04(TK) IF (AMUL.LT.ER20) GO TO 920 IF (AMUL.LT.ER10) GO TO 910 ALML = F9A04(TK) IF (ALML.LT.ER20) GO TO 920 IF (ALML.LT.ER10) GO TO 910 CPL = F19A04(TK) IF (CPL.LT.ER20) GO TO 920 IF (CPL.LT.ER10) GO TO 910 F87A04 = AMUL*CPL/ALML RETURN 910 F87A04 = -1.0E10 RETURN 920 F87A04 = -1.0E20 RETURN END C****************************************************** PRTDD(T) DOUBLE PRECISION FUNCTION F88A04(TK) IMPLICIT DOUBLE PRECISION(A-H,O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ AMUV = F15A04(TK) IF (AMUV.LT.ER20) GO TO 920 IF (AMUV.LT.ER10) GO TO 910 ALMV = F10A04(TK) IF (ALMV.LT.ER20) GO TO 920 IF (ALMV.LT.ER10) GO TO 910 CPV = F20A04(TK) IF (CPV.LT.ER20) GO TO 920 IF (CPV.LT.ER10) GO TO 910 F88A04 = AMUV*CPV/ALMV RETURN 910 F88A04 = -1.0E10 RETURN 920 F88A04 = -1.0E20 RETURN END C ***************************************************** BSPT(PA,T) DOUBLE PRECISION FUNCTION F90A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ VV = F51A04(PA, T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 WP = F83A04(PA, T) IF (WP.LT.ER20) GO TO 920 IF (WP.LT.ER10) GO TO 910 F90A04 = VV/(WP*WP) RETURN 910 F90A04 = -1.0E10 RETURN 920 F90A04 = -1.0E20 RETURN END C ***************************************************** BTPT(PA,T) DOUBLE PRECISION FUNCTION F91A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ BS = F90A04(PA, T) IF (BS.LT.ER20) GO TO 920 IF (BS.LT.ER10) GO TO 910 CPG = F18A04(PA, T) IF (CPG.LT.ER20) GO TO 920 IF (CPG.LT.ER10) GO TO 910 CVG = F77A04(PA, T) IF (CVG.LT.ER20) GO TO 920 IF (CVG.LT.ER10) GO TO 910 F91A04 = BS*CPG/CVG RETURN 910 F91A04 = -1.0E10 RETURN 920 F91A04 = -1.0E20 RETURN END C ***************************************************** BPPT(PA,T) DOUBLE PRECISION FUNCTION F92A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ CPG = F18A04(PA, T) IF (CPG.LT.ER20) GO TO 920 IF (CPG.LT.ER10) GO TO 910 CVG = F77A04(PA, T) IF (CVG.LT.ER20) GO TO 920 IF (CVG.LT.ER10) GO TO 910 VV = F51A04(PA, T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 RHOG = 1.0D00/VV DPT = G22A04(RHOG, T) IF (DPT.LT.ER20) GO TO 920 IF (DPT.LT.ER10) GO TO 910 F92A04 = (CPG-CVG)/(DPT*T*VV) RETURN 910 F92A04 = -1.0E10 RETURN 920 F92A04 = -1.0E20 RETURN END C ***************************************************** BVPT(PA,T) DOUBLE PRECISION FUNCTION F93A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ VV = F51A04(PA, T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 RHOG = 1.0D00/VV DPT = G22A04(RHOG, T) IF (DPT.LT.ER20) GO TO 920 IF (DPT.LT.ER10) GO TO 910 F93A04 = DPT/PA RETURN 910 F93A04 = -1.0E10 RETURN 920 F93A04 = -1.0E20 RETURN END C ***************************************************** AJTPT(PA,T) DOUBLE PRECISION FUNCTION F94A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA ER10/-1.0D08/, ER20/-1.0D18/ BP = F92A04(PA,T) IF (BP.LT.ER20) GO TO 920 IF (BP.LT.ER10) GO TO 910 CPG = F18A04(PA,T) IF (CPG.LT.ER20) GO TO 920 IF (CPG.LT.ER10) GO TO 910 VV = F51A04(PA,T) IF (VV.LT.ER20) GO TO 920 IF (VV.LT.ER10) GO TO 910 F94A04 = (BP*T-1.0D00)/(CPG*VV) RETURN 910 F94A04 = -1.0E10 RETURN 920 F94A04 = -1.0E20 RETURN END C****************************************************** GAMPDD(P) DOUBLE PRECISION FUNCTION F96A04(PA) IMPLICIT DOUBLE PRECISION (A-H,L-Z) DATA PMIN/1.47514D00/, ER10/-1.0D08/, ER20/-1.0D18/ DATA PMAX/2.164543D05/ IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 CPV=F17A04(PA) IF(CPV.LT.ER20) GO TO 920 IF(CPV.LT.ER10) GO TO 910 CVV=F76A04(PA) IF(CVV.LT.ER20) GO TO 920 IF(CVV.LT.ER10) GO TO 910 F96A04=CPV/CVV RETURN 910 F96A04=-1.0E+10 RETURN 920 F96A04=-1.0E+20 RETURN END C****************************************************** GAMTDD(T) DOUBLE PRECISION FUNCTION F97A04(TK) IMPLICIT DOUBLE PRECISION (A-H,L-Z) DATA TMIN/0.79999D00/, ER10/-1.0D08/, ER20/-1.0D18/ DATA TMAX/5.13001D00/ IF (TK.LT.TMIN.OR.TK.GT.TMAX) GO TO 920 CPV=F20A04(TK) IF(CPV.LT.ER20) GO TO 920 IF(CPV.LT.ER10) GO TO 910 CVV=F78A04(TK) IF(CVV.LT.ER20) GO TO 920 IF(CVV.LT.ER10) GO TO 910 F97A04=CPV/CVV RETURN 910 F97A04=-1.0E+10 RETURN 920 F97A04=-1.0E+20 RETURN END C*** F98A04 ******************************************* TPSEUP(P) DOUBLE PRECISION FUNCTION F98A04(PA) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION T(2),C(2),TL(3),TR(3),CL(3),CR(3) P1=F21A04('P') PP=DABS((PA-P1)/P1) T1=F21A04('T') IF (PA.LT.P1.OR.PA.GT.500.001D05) THEN F98A04=-1.0E20 RETURN ENDIF IF (PP.LT.1.0D-5) THEN F98A04=T1 RETURN ENDIF T2=0.85*T1 P2=F30A04(T2) TM0=T1+(T1-T2)*(PA-P1)/(P1-P2) C-----TM0 KINJICHI TC=T1 T(1)=TM0 DEL=T1*0.001D00 PBR=PA/1.0D05 IF (PA.GT.5.0D05.AND.PA.LE.50.0D05) THEN T(1)=0.362*PBR+4.5 DEL=T(1)*0.005D00 ELSE IF (PA.GT.50.0D05.AND.PA.LE.120.0D05) THEN T(1)=0.257*PBR+10.0 DEL=T(1)*0.01D00 ELSE T(1)=0.181*PBR+19.0 DEL=T(1)*0.01D00 ENDIF 150 EPS=1.0D-7 IREP=0 IREM=5000 KCONT=0 ICONT=0 C(1)=F18A04(PA,T(1)) T(2)=T(1)-DEL C(2)=F18A04(PA,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)=F18A04(PA,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=-F18A04(PA,TT) 3000 CONV=DABS((CC-C(2))/CC) IF(CONV.LT.EPS) THEN F98A04=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=-F18A04(PA,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=F18A04(PA,TA) CB=F18A04(PA,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=F18A04(PA,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)=F18A04(PA,TL(2)) CR(2)=F18A04(PA,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=F18A04(PA,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=F18A04(PA,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 F98A04=TC RETURN 8000 F98A04=-1.0E10 RETURN END C ********************************************* DD(PA,T) FOR RHO T<3.0 DOUBLE PRECISION FUNCTION G10A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.8D00, 3.0D00/ DATA PC, DLT, DELT7, MM / 2.2746D05, 1.0D-5, 1.0D-7, 5000 / IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 PM=F68A04(T) IF (PA.GT.PC) GO TO 100 PS=F30A04(T) DP = DABS((PA-PS)/PS) IF (DP.LE.DLT) GO TO 100 IF (PS.LT.-1.0E18) GO TO 920 IF (PS.LT.-1.0E08) GO TO 910 IF (PA.GT.PM.OR.PA.LT.PS) GO TO 920 100 IF (PA.LT.2.848D06) THEN D1 = 170. ELSE D1 = 210. ENDIF IC = 0 IK = 0 IN = 0 FM = 0. 8 F1 = G1A04(D1,T)-PA DF1 = G4A04(D1,T) D2 = D1-F1/DF1 DELT = DABS((D2-D1)/D2) IF (DELT.LT.DELT7.OR.DABS(F1).LT.DELT7) GO TO 200 IF ((F1*FM).LT.0.) THEN IN = IN +1 IF (IN.GT.10) GO TO 500 ELSE IN = 0 ENDIF FM = F1 D1 = D2 IC = IC+1 D1 = D2 IF (IC.LT.MM) THEN GO TO 8 ELSE GO TO 910 ENDIF 200 G10A04 = D2 RETURN 500 DM = (D1+D2)/2.D00 F2 = G1A04(DM,T)-PA IF (F1*F2.GT.0) THEN D1 = DM ELSE D2 = DM ENDIF IK = IK+1 IF (IK.GT.2000) GO TO 910 DLTT2 = DABS((D1-D2)/D2) IF (DLTT2.GT.DELT7) GO TO 500 GO TO 200 910 G10A04 = -1.E10 RETURN 920 G10A04 = -1.E20 RETURN END C ********************************************* DD(PA,T) FOR RHO T>2.3 DOUBLE PRECISION FUNCTION G11A04(PA, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / DATA PC /2.2746D05/ PI=0.33033259D-02 R=0.831431D-02 P0=PA*1.0D-6 RHO=200. IF (PA.LT.PC) THEN PS=F30A04(T) IF (PA.LT.PS) RHO=1.0 ENDIF IF (PA.GT.15.0D06) RHO=250. RHO=RHO/4.0026D00 IF (PA.LT.2.105D05.AND.T.GT.5.18.AND.T.LT.5.80) GO TO 500 IC=0 10 IC=IC+1 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH4*RH2 RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH9*RH2 RH12=RH11*RHO RH13=RH11*RH2 T12=T**0.5 T2=T*T T3=T2*T T4=T3*T DPP=DEXP(-PI*RH2) PIRH=2.D00*PI*RH2 C1=DN(1)*T+DN(2)*T12+DN(3)+DN(4)/T+DN(5)/T2 C6=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 C10=DN(10)*T+DN(11)+DN(12)/T C13=DN(13) C14=DN(14)/T+DN(15)/T2 C16=DN(16)/T C17=DN(17)/T+DN(18)/T2 C19=DN(19)/T2 C20=(DN(20)/T2+DN(21)/T3)*DPP C22=(DN(22)/T2+DN(23)/T4)*DPP C24=(DN(24)/T2+DN(25)/T3)*DPP C26=(DN(26)/T2+DN(27)/T4)*DPP C28=(DN(28)/T2+DN(29)/T3)*DPP C30=(DN(30)/T2+DN(31)/T3+DN(32)/T4)*DPP FN1=RH2*C1 FN2=RH3*C6 FN3=RH4*C10 FN4=RH5*C13 FN5=RH6*C14 FN6=RH7*C16 FN7=RH8*C17 FN8=RH9*C19 FN9=RH3*C20 FN10=RH5*C22 FN11=RH7*C24 FN12=RH9*C26 FN13=RH11*C28 FN14=RH13*C30 FUNC=RHO*R*T+FN1+FN2+FN3+FN4+FN5+FN6+FN7+FN8+FN9+FN10+FN11+FN12+ & FN13+FN14 FUN=FUNC-P0 DFN1=2.D00*RHO*C1 DFN2=3.D00*RH2*C6 DFN3=4.D00*RH3*C10 DFN4=5.D00*RH4*C13 DFN5=6.D00*RH5*C14 DFN6=7.D00*RH6*C16 DFN7=8.D00*RH7*C17 DFN8=9.D00*RH8*C19 DFN9=RH2*C20*(3.D00-PIRH) DFN10=RH4*C22*(5.D00-PIRH) DFN11=RH6*C24*(7.D00-PIRH) DFN12=RH8*C26*(9.D00-PIRH) DFN13=RH10*C28*(11.D00-PIRH) DFN14=RH12*C30*(13.D00-PIRH) DFN=R*T+DFN1+DFN2+DFN3+DFN4+DFN5+DFN6+DFN7+DFN8+DFN9 & +DFN10+DFN11+DFN12+DFN13+DFN14 RO2=RHO-FUN/DFN IF (IC.GT.5000) GO TO 910 IF (DABS((RHO-RO2)/RHO).GT.1.0D-8) THEN RHO=RO2 GO TO 10 ENDIF DD=4.0026D00*RHO G11A04=DD RETURN 500 G11A04 = G27A04(PA,T) RETURN 910 G11A04 = -1.0E10 RETURN END C ********************************************* UU(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G12A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B14=1.D00/4.D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C10=10.0D00 C12=12.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=C05*DN(2)*T05+DN(3)+C2*DN(4)/T+C3*DN(5)/T2 CC2=DN(7)+C2*DN(8)/T+C3*DN(9)/T2 CC3=DN(11)+C2*DN(12)/T CC4=DN(13) CC5=C2*DN(14)/T+C3*DN(15)/T2 CC6=C2*DN(16)/T CC7=C2*DN(17)/T+C3*DN(18)/T2 CC8=C3*DN(19)/T2 CC9=C3*DN(20)/T2+C4*DN(21)/T3 CC10=C3*DN(22)/T2+C5*DN(23)/T4 CC11=C3*DN(24)/T2+C4*DN(25)/T3 CC12=C3*DN(26)/T2+C5*DN(27)/T4 CC13=C3*DN(28)/T2+C4*DN(29)/T3 CC14=C3*DN(30)/T2+C4*DN(31)/T3+C5*DN(32)/T4 REX=DEXP(-PI*RH2) U1=CC1*RHO+C05*CC2*RH2+B13*CC3*RH3+B14*CC4*RH4 1 +C02*CC5*RH5+B16*CC6*RH6+B17*CC7*RH7+B18*CC8*RH8 2 -CC9*(REX-C1)/(PI*C2)-CC10*(REX*(PIRH2+C1)-C1)/(C2*PI2) 3 -CC11*(REX*(P2RH4+PIRH2+C1)-C1)/PI3 4 -C3*CC12*(REX*(P3RH6+P2RH4+PIRH2+C1)-C1)/PI4 5 -C12*CC13*(REX*(P4RH8+P3RH6+P2RH4+PIRH2+C1)-C1)/PI5 6 -C60*CC14*(REX*(P5RH10+P4RH8+P3RH6+P2RH4+PIRH2+C1) 7 -C1)/PI6 U0=61.132D0 U1=U1*1.0D03 U2=(C5/C2*R*T-R*T)*1.0D03 U=U0+U1+U2 UU=U/4.0026D00*1.0D03 G12A04=UU RETURN END C *********************************************** SS(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G13A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C10=10.0D00 C12=12.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 PREF=0.101325D00 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=-DN(1)-C05*DN(2)/T05+DN(4)/T2+C2*DN(5)/T3 CC2=-DN(6)+DN(8)/T2+C2*DN(9)/T3 CC3=-DN(10)+DN(12)/T2 CC4=DN(14)/T2+C2*DN(15)/T3 CC5=DN(16)/T2 CC6=DN(17)/T2+C2*DN(18)/T3 CC7=C2*DN(19)/T3 CC8=C2*DN(20)/T3+C3*DN(21)/T4 CC9=C2*DN(22)/T3+C4*DN(23)/T5 CC10=C2*DN(24)/T3+C3*DN(25)/T4 CC11=C2*DN(26)/T3+C4*DN(27)/T5 CC12=C2*DN(28)/T3+C3*DN(29)/T4 CC13=C2*DN(30)/T3+C3*DN(31)/T4+C4*DN(32)/T5 REX=DEXP(-PI*RH2) S1=CC1*RHO+C05*CC2*RH2+B13*CC3*RH3+C02*CC4*RH5 1 +B16*CC5*RH6+B17*CC6*RH7+B18*CC7*RH8 2 -CC8*(REX-C1)/(PI*C2)-CC9*(REX*(PIRH2+C1)-C1)/(C2*PI2) 3 -CC10*(REX*(P2RH4+PIRH2+C1)-C1)/PI3 4 -C3*CC11*(REX*(P3RH6+P2RH4+PIRH2+C1)-C1)/PI4 5 -C12*CC12*(REX*(P4RH8+P3RH6+P2RH4+PIRH2+C1)-C1)/PI5 6 -C60*CC13*(REX*(P5RH10+P4RH8+P3RH6+P2RH4+PIRH2+C1) 7 -C1)/PI6 S0=7.862D0+20.7857D0*DLOG(T) S1=S1*1.0D03 S2=R*DLOG(R*RHO*T/PREF)*1.0D03 S=S0-S2+S1 SS=S/4.0026D00*1.0D03 G13A04=SS RETURN END C *********************************************** CP(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G14A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B14=1.D00/4.D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C9=9.0D00 C10=10.0D00 C11=11.0D00 C12=12.0D00 C13=13.0D00 C20=20.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH10*RHO RH12=RH11*RHO RH13=RH12*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=-B14*DN(2)/T05+C2*DN(4)/T2+C6*DN(5)/T3 CC2=C2*DN(8)/T2+C6*DN(9)/T3 CC3=C2*DN(12)/T2 CC4=C2*DN(14)/T2+C6*DN(15)/T3 CC5=C2*DN(16)/T2 CC6=C2*DN(17)/T2+C6*DN(18)/T3 CC7=C6*DN(19)/T3 CC8=C6*DN(20)/T3+C12*DN(21)/T4 CC9=C6*DN(22)/T3+C20*DN(23)/T5 CC10=C6*DN(24)/T3+C12*DN(25)/T4 CC11=C6*DN(26)/T3+C20*DN(27)/T5 CC12=C6*DN(28)/T3+C12*DN(29)/T4 CC13=C6*DN(30)/T3+C12*DN(31)/T4+C20*DN(32)/T5 REX=DEXP(-PI*RH2) CV1=-CC1*RHO-C05*CC2*RH2-B13*CC3*RH3-C02*CC4*RH5 1 -B16*CC5*RH6-B17*CC6*RH7-B18*CC7*RH8 2 +CC8*(REX-C1)/(PI*C2)+CC9*(REX*(PIRH2+C1)-C1)/(C2*PI2) 3 +CC10*(REX*(P2RH4+PIRH2+C1)-C1)/PI3 4 +C3*CC11*(REX*(P3RH6+P2RH4+PIRH2+C1)-C1)/PI4 5 +C12*CC12*(REX*(P4RH8+P3RH6+P2RH4+PIRH2+C1)-C1)/PI5 6 +C60*CC13*(REX*(P5RH10+P4RH8+P3RH6+P2RH4+PIRH2+C1) 7 -C1)/PI6 CV1=CV1*1.0D03 CV0=C3/C2*R*1.0D03 CV=(CV0+CV1)/4.0026D00 CT1=DN(1)+C05*DN(2)/T05-DN(4)/T2-C2*DN(5)/T3 CT2=DN(6)-DN(8)/T2-C2*DN(9)/T3 CT3=DN(10)-DN(12)/T2 CT4=-DN(14)/T2-C2*DN(15)/T3 CT5=-DN(16)/T2 CT6=-DN(17)/T2-C2*DN(18)/T3 CT7=-C2*DN(19)/T3 CT8=-C2*DN(20)/T3-C3*DN(21)/T4 CT9=-C2*DN(22)/T3-C4*DN(23)/T5 CT10=-C2*DN(24)/T3-C3*DN(25)/T4 CT11=-C2*DN(26)/T3-C4*DN(27)/T5 CT12=-C2*DN(28)/T3-C3*DN(29)/T4 CT13=-C2*DN(30)/T3-C3*DN(31)/T4-C4*DN(32)/T5 CD1=DN(1)*T+DN(2)*T05+DN(3)+DN(4)/T+DN(5)/T2 CD2=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 CD3=DN(10)*T+DN(11)+DN(12)/T CD4=DN(13) CD5=DN(14)/T+DN(15)/T2 CD6=DN(16)/T CD7=DN(17)/T+DN(18)/T2 CD8=DN(19)/T2 CD9=DN(20)/T2+DN(21)/T3 CD10=DN(22)/T2+DN(23)/T4 CD11=DN(24)/T2+DN(25)/T3 CD12=DN(26)/T2+DN(27)/T4 CD13=DN(28)/T2+DN(29)/T3 CD14=DN(30)/T2+DN(31)/T3+DN(32)/T4 DPDT=RHO*R+CT1*RH2+CT2*RH3+CT3*RH4+CT4*RH6+CT5*RH7+CT6*RH8 1 +CT7*RH9+(CT8*RH3+CT9*RH5+CT10*RH7+CT11*RH9+CT12*RH11+CT13* 2 RH13)*REX DPDR=T*R+C2*CD1*RHO+C3*CD2*RH2+C4*CD3*RH3+C5*CD4*RH4+C6*CD5*RH5 1 +C7*CD6*RH6+C8*CD7*RH7+C9*CD8*RH8+(CD9*(C3-C2*PIRH2)*RH2+CD10* 2 (C5-C2*PIRH2)*RH4+CD11*(C7-C2*PIRH2)*RH6+CD12*(C9-C2*PIRH2)*RH8 3 +CD13*(C11-C2*PIRH2)*RH10+CD14*(C13-C2*PIRH2)*RH12)*REX CP1=T*DPDT*DPDT/DPDR/RH2/4.0026D00*1.0D03 CP=CV+CP1 G14A04=CP*1.0D03 RETURN END C *********************************************** CV(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G18A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B14=1.D00/4.D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C9=9.0D00 C10=10.0D00 C11=11.0D00 C12=12.0D00 C13=13.0D00 C20=20.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH10*RHO RH12=RH11*RHO RH13=RH12*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=-B14*DN(2)/T05+C2*DN(4)/T2+C6*DN(5)/T3 CC2=C2*DN(8)/T2+C6*DN(9)/T3 CC3=C2*DN(12)/T2 CC4=C2*DN(14)/T2+C6*DN(15)/T3 CC5=C2*DN(16)/T2 CC6=C2*DN(17)/T2+C6*DN(18)/T3 CC7=C6*DN(19)/T3 CC8=C6*DN(20)/T3+C12*DN(21)/T4 CC9=C6*DN(22)/T3+C20*DN(23)/T5 CC10=C6*DN(24)/T3+C12*DN(25)/T4 CC11=C6*DN(26)/T3+C20*DN(27)/T5 CC12=C6*DN(28)/T3+C12*DN(29)/T4 CC13=C6*DN(30)/T3+C12*DN(31)/T4+C20*DN(32)/T5 REX=DEXP(-PI*RH2) CV1=-CC1*RHO-C05*CC2*RH2-B13*CC3*RH3-C02*CC4*RH5 1 -B16*CC5*RH6-B17*CC6*RH7-B18*CC7*RH8 2 +CC8*(REX-C1)/(PI*C2)+CC9*(REX*(PIRH2+C1)-C1)/(C2*PI2) 3 +CC10*(REX*(P2RH4+PIRH2+C1)-C1)/PI3 4 +C3*CC11*(REX*(P3RH6+P2RH4+PIRH2+C1)-C1)/PI4 5 +C12*CC12*(REX*(P4RH8+P3RH6+P2RH4+PIRH2+C1)-C1)/PI5 6 +C60*CC13*(REX*(P5RH10+P4RH8+P3RH6+P2RH4+PIRH2+C1) 7 -C1)/PI6 CV1=CV1*1.0D03 CV0=C3/C2*R*1.0D03 CV=(CV0+CV1)/4.0026D00*1.0D03 G18A04=CV RETURN END C ********************************************* DP/DT(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G20A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B14=1.D00/4.D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C9=9.0D00 C10=10.0D00 C11=11.0D00 C12=12.0D00 C13=13.0D00 C20=20.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH10*RHO RH12=RH11*RHO RH13=RH12*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=-B14*DN(2)/T05+C2*DN(4)/T2+C6*DN(5)/T3 CC2=C2*DN(8)/T2+C6*DN(9)/T3 CC3=C2*DN(12)/T2 CC4=C2*DN(14)/T2+C6*DN(15)/T3 CC5=C2*DN(16)/T2 CC6=C2*DN(17)/T2+C6*DN(18)/T3 CC7=C6*DN(19)/T3 CC8=C6*DN(20)/T3+C12*DN(21)/T4 CC9=C6*DN(22)/T3+C20*DN(23)/T5 CC10=C6*DN(24)/T3+C12*DN(25)/T4 CC11=C6*DN(26)/T3+C20*DN(27)/T5 CC12=C6*DN(28)/T3+C12*DN(29)/T4 CC13=C6*DN(30)/T3+C12*DN(31)/T4+C20*DN(32)/T5 REX=DEXP(-PI*RH2) CT1=DN(1)+C05*DN(2)/T05-DN(4)/T2-C2*DN(5)/T3 CT2=DN(6)-DN(8)/T2-C2*DN(9)/T3 CT3=DN(10)-DN(12)/T2 CT4=-DN(14)/T2-C2*DN(15)/T3 CT5=-DN(16)/T2 CT6=-DN(17)/T2-C2*DN(18)/T3 CT7=-C2*DN(19)/T3 CT8=-C2*DN(20)/T3-C3*DN(21)/T4 CT9=-C2*DN(22)/T3-C4*DN(23)/T5 CT10=-C2*DN(24)/T3-C3*DN(25)/T4 CT11=-C2*DN(26)/T3-C4*DN(27)/T5 CT12=-C2*DN(28)/T3-C3*DN(29)/T4 CT13=-C2*DN(30)/T3-C3*DN(31)/T4-C4*DN(32)/T5 DPDT=RHO*R+CT1*RH2+CT2*RH3+CT3*RH4+CT4*RH6+CT5*RH7+CT6*RH8 1 +CT7*RH9+(CT8*RH3+CT9*RH5+CT10*RH7+CT11*RH9+CT12*RH11+CT13* 2 RH13)*REX G20A04 = DPDT*1.D06 RETURN END C ********************************************* DP/DR(DD,T) FOR T>2.3 DOUBLE PRECISION FUNCTION G21A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / C05=0.5D00 B13=1.D00/3.D00 C02=0.2D00 B14=1.D00/4.D00 B16=1.D00/6.D00 B17=1.D00/7.D00 B18=1.0D00/8.0D00 C1=1.0D00 C2=2.0D00 C3=3.0D00 C4=4.0D00 C5=5.0D00 C6=6.0D00 C7=7.0D00 C8=8.0D00 C9=9.0D00 C10=10.0D00 C11=11.0D00 C12=12.0D00 C13=13.0D00 C20=20.0D00 C30=30.0D00 C24=24.0D00 C60=60.0D00 C120=120.0D00 PI=0.33033259D-02 PIT2=C2*PI PI2=PI*PI PI3=PI2*PI PI4=PI3*PI PI5=PI4*PI PI6=PI5*PI R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH5*RHO RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH10*RHO RH12=RH11*RHO RH13=RH12*RHO PIRH2=PI*RH2 P2RH4=PIRH2*PIRH2/C2 P3RH6=P2RH4*PIRH2/C3 P4RH8=P3RH6*PIRH2/C4 P5RH10=P4RH8*PIRH2/C5 T05=T**0.5 T2=T*T T3=T2*T T4=T3*T T5=T4*T CC1=-B14*DN(2)/T05+C2*DN(4)/T2+C6*DN(5)/T3 CC2=C2*DN(8)/T2+C6*DN(9)/T3 CC3=C2*DN(12)/T2 CC4=C2*DN(14)/T2+C6*DN(15)/T3 CC5=C2*DN(16)/T2 CC6=C2*DN(17)/T2+C6*DN(18)/T3 CC7=C6*DN(19)/T3 CC8=C6*DN(20)/T3+C12*DN(21)/T4 CC9=C6*DN(22)/T3+C20*DN(23)/T5 CC10=C6*DN(24)/T3+C12*DN(25)/T4 CC11=C6*DN(26)/T3+C20*DN(27)/T5 CC12=C6*DN(28)/T3+C12*DN(29)/T4 CC13=C6*DN(30)/T3+C12*DN(31)/T4+C20*DN(32)/T5 REX=DEXP(-PI*RH2) CD1=DN(1)*T+DN(2)*T05+DN(3)+DN(4)/T+DN(5)/T2 CD2=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 CD3=DN(10)*T+DN(11)+DN(12)/T CD4=DN(13) CD5=DN(14)/T+DN(15)/T2 CD6=DN(16)/T CD7=DN(17)/T+DN(18)/T2 CD8=DN(19)/T2 CD9=DN(20)/T2+DN(21)/T3 CD10=DN(22)/T2+DN(23)/T4 CD11=DN(24)/T2+DN(25)/T3 CD12=DN(26)/T2+DN(27)/T4 CD13=DN(28)/T2+DN(29)/T3 CD14=DN(30)/T2+DN(31)/T3+DN(32)/T4 DPDR=T*R+C2*CD1*RHO+C3*CD2*RH2+C4*CD3*RH3+C5*CD4*RH4+C6*CD5*RH5 1 +C7*CD6*RH6+C8*CD7*RH7+C9*CD8*RH8+(CD9*(C3-C2*PIRH2)*RH2+CD10* 2 (C5-C2*PIRH2)*RH4+CD11*(C7-C2*PIRH2)*RH6+CD12*(C9-C2*PIRH2)*RH8 3 +CD13*(C11-C2*PIRH2)*RH10+CD14*(C13-C2*PIRH2)*RH12)*REX G21A04=DPDR/4.0026D00*1.0D06 RETURN END C ********************************************************* DP/DT(DD,T) DOUBLE PRECISION FUNCTION G22A04(DD, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.79999D00, 1500.0001D00/ DPDT2 = 5000. TA = 2.0 TM = 3.0 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (T.GT.2.3) DPDT2 = G20A04(DD, T) IF (T.LT.2.98.AND.DD.GT.140.AND.DD.LT.189.) THEN DPDT1 = G3A04(DD, T) TA=2.98-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) TM=2.53-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) IF (DD.GT.140.AND.DD.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN G22A04 = DPDT1*((TA-T)/(TA-TM))+DPDT2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN G22A04 = DPDT2 RETURN ELSEIF (T.LT.TM) THEN G22A04 = DPDT1 RETURN ENDIF ENDIF ELSE G22A04 = DPDT2 RETURN ENDIF 920 G22A04 = -1.E20 RETURN END C ******************************************************** DP/DR(DD,T) DOUBLE PRECISION FUNCTION G23A04(DD, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.79999D00, 1500.0001D00/ DPDR2 = 5000. TA = 2.0 TM = 3.0 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 IF (T.GT.2.3) DPDR2=G21A04(DD, T) IF (T.LT.2.98.AND.DD.GT.140.AND.DD.LT.189.) THEN DPDR1=G4A04(DD, T) TA=2.98-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) TM=2.53-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) IF (DD.GT.140.AND.DD.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN G23A04=DPDR1*((TA-T)/(TA-TM))+DPDR2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN G23A04 = DPDR2 RETURN ELSEIF (T.LT.TM) THEN G23A04 = DPDR1 RETURN ENDIF ENDIF ELSE G23A04=DPDR2 RETURN ENDIF 920 G23A04 = -1.E20 RETURN END C ******************************************* FUNC(DD,T) FOR TPV T>2.3 DOUBLE PRECISION FUNCTION G24A04(DD, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / PI=0.33033259D-02 R=0.831431D-02 RHO=DD/4.0026D00 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH4*RH2 RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH9*RH2 RH12=RH11*RHO RH13=RH11*RH2 T12=T**0.5 T2=T*T T3=T2*T T4=T3*T DPP=DEXP(-PI*RH2) PIRH=2.D00*PI*RH2 C1=DN(1)*T+DN(2)*T12+DN(3)+DN(4)/T+DN(5)/T2 C6=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 C10=DN(10)*T+DN(11)+DN(12)/T C13=DN(13) C14=DN(14)/T+DN(15)/T2 C16=DN(16)/T C17=DN(17)/T+DN(18)/T2 C19=DN(19)/T2 C20=(DN(20)/T2+DN(21)/T3)*DPP C22=(DN(22)/T2+DN(23)/T4)*DPP C24=(DN(24)/T2+DN(25)/T3)*DPP C26=(DN(26)/T2+DN(27)/T4)*DPP C28=(DN(28)/T2+DN(29)/T3)*DPP C30=(DN(30)/T2+DN(31)/T3+DN(32)/T4)*DPP FN1=RH2*C1 FN2=RH3*C6 FN3=RH4*C10 FN4=RH5*C13 FN5=RH6*C14 FN6=RH7*C16 FN7=RH8*C17 FN8=RH9*C19 FN9=RH3*C20 FN10=RH5*C22 FN11=RH7*C24 FN12=RH9*C26 FN13=RH11*C28 FN14=RH13*C30 FUNC=RHO*R*T+FN1+FN2+FN3+FN4+FN5+FN6+FN7+FN8+FN9+FN10+FN11+FN12+ & FN13+FN14 G24A04=FUNC RETURN END C ******************************************* DDL(PA,L) FOR VPD T>2.3 DOUBLE PRECISION FUNCTION G25A04(PA, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / PI=0.33033259D-02 R=0.831431D-02 P0=PA*1.0D-6 RHO=200. RHO=RHO/4.0026D00 IC=0 10 IC=IC+1 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH4*RH2 RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH9*RH2 RH12=RH11*RHO RH13=RH11*RH2 T12=T**0.5 T2=T*T T3=T2*T T4=T3*T DPP=DEXP(-PI*RH2) PIRH=2.D00*PI*RH2 C1=DN(1)*T+DN(2)*T12+DN(3)+DN(4)/T+DN(5)/T2 C6=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 C10=DN(10)*T+DN(11)+DN(12)/T C13=DN(13) C14=DN(14)/T+DN(15)/T2 C16=DN(16)/T C17=DN(17)/T+DN(18)/T2 C19=DN(19)/T2 C20=(DN(20)/T2+DN(21)/T3)*DPP C22=(DN(22)/T2+DN(23)/T4)*DPP C24=(DN(24)/T2+DN(25)/T3)*DPP C26=(DN(26)/T2+DN(27)/T4)*DPP C28=(DN(28)/T2+DN(29)/T3)*DPP C30=(DN(30)/T2+DN(31)/T3+DN(32)/T4)*DPP FN1=RH2*C1 FN2=RH3*C6 FN3=RH4*C10 FN4=RH5*C13 FN5=RH6*C14 FN6=RH7*C16 FN7=RH8*C17 FN8=RH9*C19 FN9=RH3*C20 FN10=RH5*C22 FN11=RH7*C24 FN12=RH9*C26 FN13=RH11*C28 FN14=RH13*C30 FUNC=RHO*R*T+FN1+FN2+FN3+FN4+FN5+FN6+FN7+FN8+FN9+FN10+FN11+FN12+ & FN13+FN14 FUN=FUNC-P0 DFN1=2.D00*RHO*C1 DFN2=3.D00*RH2*C6 DFN3=4.D00*RH3*C10 DFN4=5.D00*RH4*C13 DFN5=6.D00*RH5*C14 DFN6=7.D00*RH6*C16 DFN7=8.D00*RH7*C17 DFN8=9.D00*RH8*C19 DFN9=RH2*C20*(3.D00-PIRH) DFN10=RH4*C22*(5.D00-PIRH) DFN11=RH6*C24*(7.D00-PIRH) DFN12=RH8*C26*(9.D00-PIRH) DFN13=RH10*C28*(11.D00-PIRH) DFN14=RH12*C30*(13.D00-PIRH) DFN=R*T+DFN1+DFN2+DFN3+DFN4+DFN5+DFN6+DFN7+DFN8+DFN9 & +DFN10+DFN11+DFN12+DFN13+DFN14 RO2=RHO-FUN/DFN IF (IC.GT.5000) GO TO 910 IF (DABS((RHO-RO2)/RHO).GT.1.0D-8) THEN RHO=RO2 GO TO 10 ENDIF DD=4.0026D00*RHO G25A04=DD RETURN 910 G25A04= -1.E10 RETURN END C ********************************************* FUNC(PA, DD, T) FOR TPV DOUBLE PRECISION FUNCTION G26A04(PA, DD, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.79999D00, 1500.0001D00/ P0 = PA/1.0D06 TA = 2.0 TM = 3.0 IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 920 FUNC2=G24A04(DD, T)-P0 IF (T.LT.2.98.AND.DD.GT.140.AND.DD.LT.189.) THEN FUNC1=G1A04(DD, T)-PA TA=2.98-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) TM=2.53-5.6D-3*(DD-140.)-3.5D-2*(DMAX1(0.D00,DD-180.D00)) IF (DD.GT.140.AND.DD.LT.189.) THEN IF (T.LE.TA.AND.T.GE.TM) THEN G26A04=FUNC1*((TA-T)/(TA-TM))+FUNC2*((T-TM)/(TA-TM)) RETURN ELSEIF (T.GT.TA) THEN G26A04 = FUNC2 RETURN ELSEIF (T.LT.TM) THEN G26A04 = FUNC1 RETURN ENDIF ENDIF ELSE G26A04 = FUNC2 RETURN ENDIF 920 G26A04 = -1.E20 RETURN END C ************************************************** DD(PA,T) FOR G11A04 DOUBLE PRECISION FUNCTION G27A04(PA,T) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DATA DELT8 /1.0D-8/ IC = 0 TS = F40A04(PA) P0 = PA/1.0D06 T1 = TS+0.01 D1 = G28A04(PA,T1) IF (D1.LT.-1.0E08) GO TO 910 F1 = G24A04(D1,T)-P0 T2 = 6.0 D2 = G28A04(PA,T2) IF (D2.LT.-1.0E08) GO TO 910 50 DM = (D1+D2)/2.D00 F2 = G24A04(DM,T)-P0 IF (F1*F2.GT.0) THEN D1 = DM ELSE D2 = DM ENDIF IC = IC+1 IF (IC.GT.2000) GO TO 910 DLTT2 = DABS((D1-D2)/D2) IF (DLTT2.GT.DELT8) GO TO 50 200 G27A04 = DM RETURN 910 G27A04 = -1.0E+10 RETURN END C ************************************************ DD(PA,T) FOR G36A04 DOUBLE PRECISION FUNCTION G28A04(PA, T) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / DATA PC /2.2746D05/ PI=0.33033259D-02 R=0.831431D-02 P=PA*1.0D-6 PB=PA*1.0D-5 RHO=200. IF (PA.LT.PC) THEN PS=F30A04(T) IF (PA.LT.PS) RHO=1.0 ENDIF RHO=RHO/4.0026D00 IC=0 10 IC=IC+1 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH4*RH2 RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH9*RH2 RH12=RH11*RHO RH13=RH11*RH2 T12=T**0.5 T2=T*T T3=T2*T T4=T3*T DPP=DEXP(-PI*RH2) PIRH=2.D00*PI*RH2 C1=DN(1)*T+DN(2)*T12+DN(3)+DN(4)/T+DN(5)/T2 C6=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 C10=DN(10)*T+DN(11)+DN(12)/T C13=DN(13) C14=DN(14)/T+DN(15)/T2 C16=DN(16)/T C17=DN(17)/T+DN(18)/T2 C19=DN(19)/T2 C20=(DN(20)/T2+DN(21)/T3)*DPP C22=(DN(22)/T2+DN(23)/T4)*DPP C24=(DN(24)/T2+DN(25)/T3)*DPP C26=(DN(26)/T2+DN(27)/T4)*DPP C28=(DN(28)/T2+DN(29)/T3)*DPP C30=(DN(30)/T2+DN(31)/T3+DN(32)/T4)*DPP FN1=RH2*C1 FN2=RH3*C6 FN3=RH4*C10 FN4=RH5*C13 FN5=RH6*C14 FN6=RH7*C16 FN7=RH8*C17 FN8=RH9*C19 FN9=RH3*C20 FN10=RH5*C22 FN11=RH7*C24 FN12=RH9*C26 FN13=RH11*C28 FN14=RH13*C30 FUNC=RHO*R*T+FN1+FN2+FN3+FN4+FN5+FN6+FN7+FN8+FN9+FN10+FN11+FN12+ & FN13+FN14 FUN=FUNC-P DFN1=2.D00*RHO*C1 DFN2=3.D00*RH2*C6 DFN3=4.D00*RH3*C10 DFN4=5.D00*RH4*C13 DFN5=6.D00*RH5*C14 DFN6=7.D00*RH6*C16 DFN7=8.D00*RH7*C17 DFN8=9.D00*RH8*C19 DFN9=RH2*C20*(3.D00-PIRH) DFN10=RH4*C22*(5.D00-PIRH) DFN11=RH6*C24*(7.D00-PIRH) DFN12=RH8*C26*(9.D00-PIRH) DFN13=RH10*C28*(11.D00-PIRH) DFN14=RH12*C30*(13.D00-PIRH) DFN=R*T+DFN1+DFN2+DFN3+DFN4+DFN5+DFN6+DFN7+DFN8+DFN9 & +DFN10+DFN11+DFN12+DFN13+DFN14 RO2=RHO-FUN/DFN IF (IC.GT.5000) GO TO 910 IF (DABS((RHO-RO2)/RHO).GT.1.0D-8) THEN RHO=RO2 GO TO 10 ENDIF DD=4.0026D00*RHO G28A04=DD RETURN 910 G28A04 = -1.0E10 RETURN END C ******************************************************** DDV(PA) DOUBLE PRECISION FUNCTION G29A04(PA) IMPLICIT REAL*8 (A-H,O-Z) DIMENSION DN(32) DATA DN/ & 0.4558980227431D-04, 0.1260692007853D-02, -0.7139657549318D-02, & 0.9728903861441D-02, -0.1589302471562D-01, 0.1454229259623D-05, &-0.4708238429298D-04, 0.1132915232587D-02, 0.2410763742104D-02, &-0.5093547838381D-08, 0.2699726927900D-05, -0.3954146691114D-04, & 0.1551961438127D-08, 0.1050712335785D-07, -0.5501158366750D-07, &-0.1037673478521D-09, 0.6446881346448D-12, 0.3298960057071D-10, &-0.3555585738784D-12, -0.6885401367690D-02, 0.9166109232806D-02, &-0.6544314242937D-05, -0.3315398880031D-04, -0.2067693644676D-07, & 0.3850153114958D-07, -0.1399040626999D-10, -0.1888462892389D-11, &-0.4595138561035D-14, 0.6872567403738D-14, -0.6097223119177D-18, &-0.7636186157005D-17, 0.3848665703556D-17 / DATA PMIN, PMAX / 1.47514D00, 2.2746D05 / IF (PA.LT.PMIN.OR.PA.GT.PMAX) GO TO 920 TA = 2.0 TM = 3.0 T = F40A04(PA) PI=0.33033259D-02 R=0.831431D-02 P=PA*1.0D-6 RHO=0.5 RHO=RHO/4.0026D00 IC=0 10 IC=IC+1 RH2=RHO*RHO RH3=RH2*RHO RH4=RH3*RHO RH5=RH4*RHO RH6=RH4*RH2 RH7=RH6*RHO RH8=RH7*RHO RH9=RH8*RHO RH10=RH9*RHO RH11=RH9*RH2 RH12=RH11*RHO RH13=RH11*RH2 T12=T**0.5 T2=T*T T3=T2*T T4=T3*T DPP=DEXP(-PI*RH2) PIRH=2.D00*PI*RH2 C1=DN(1)*T+DN(2)*T12+DN(3)+DN(4)/T+DN(5)/T2 C6=DN(6)*T+DN(7)+DN(8)/T+DN(9)/T2 C10=DN(10)*T+DN(11)+DN(12)/T C13=DN(13) C14=DN(14)/T+DN(15)/T2 C16=DN(16)/T C17=DN(17)/T+DN(18)/T2 C19=DN(19)/T2 C20=(DN(20)/T2+DN(21)/T3)*DPP C22=(DN(22)/T2+DN(23)/T4)*DPP C24=(DN(24)/T2+DN(25)/T3)*DPP C26=(DN(26)/T2+DN(27)/T4)*DPP C28=(DN(28)/T2+DN(29)/T3)*DPP C30=(DN(30)/T2+DN(31)/T3+DN(32)/T4)*DPP FN1=RH2*C1 FN2=RH3*C6 FN3=RH4*C10 FN4=RH5*C13 FN5=RH6*C14 FN6=RH7*C16 FN7=RH8*C17 FN8=RH9*C19 FN9=RH3*C20 FN10=RH5*C22 FN11=RH7*C24 FN12=RH9*C26 FN13=RH11*C28 FN14=RH13*C30 FUNC=RHO*R*T+FN1+FN2+FN3+FN4+FN5+FN6+FN7+FN8+FN9+FN10+FN11+FN12+ & FN13+FN14 FUN=FUNC-P DFN1=2.D00*RHO*C1 DFN2=3.D00*RH2*C6 DFN3=4.D00*RH3*C10 DFN4=5.D00*RH4*C13 DFN5=6.D00*RH5*C14 DFN6=7.D00*RH6*C16 DFN7=8.D00*RH7*C17 DFN8=9.D00*RH8*C19 DFN9=RH2*C20*(3.D00-PIRH) DFN10=RH4*C22*(5.D00-PIRH) DFN11=RH6*C24*(7.D00-PIRH) DFN12=RH8*C26*(9.D00-PIRH) DFN13=RH10*C28*(11.D00-PIRH) DFN14=RH12*C30*(13.D00-PIRH) DFN=R*T+DFN1+DFN2+DFN3+DFN4+DFN5+DFN6+DFN7+DFN8+DFN9 & +DFN10+DFN11+DFN12+DFN13+DFN14 RO2=RHO-FUN/DFN IF (IC.GT.5000) GO TO 910 IF (DABS((RHO-RO2)/RHO).GT.1.0D-8) THEN RHO=RO2 GO TO 10 ENDIF RHO=4.0026D00*RHO G29A04 = RHO RETURN 910 G29A04 = -1.0E10 RETURN 920 G29A04 = -1.0E20 RETURN END C ******************************************** SUB. TPH(P,H,....) *** SUBROUTINE S4A04 (PA,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 (RDELT.LT.1.D-10.OR.DABS((T0-T1)/T1).LE.1.D-07) RETURN 300 IC = IC+1 IF (IC.GT.20000) GO TO 910 T0 = T1 H0 = H1 T1 = T1+DELT IF (T1.LT.TMIN) T1=TMIN H1 = F25A04(PA,T1) IF (T1.GT.0.OR.H1.GT.-1.D05) GO TO 100 910 T1 = -1.E+10 RETURN END C ******************************************** SUB. TPS(P,S,...) *** SUBROUTINE S5A04 (PA,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 (RDELT.LT.1.D-10.OR.DABS((T0-T1)/T1).LT.1.D-07) RETURN 300 IC = IC+1 IF (IC.GT.20000) GO TO 910 T0 = T1 S0 = S1 T1 = T1+DELT IF (T1.LT.TMIN) T1=TMIN S1 = F35A04(PA,T1) IF (T1.GT.0.OR.S1.GT.-1.D05) GO TO 100 910 T1 = -1.E+10 RETURN END C *************************************** THERMAL CONDUCTIVITY *** SUBROUTINE S6A04(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 S7A04(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 ************************************* COEFFICIENT OF VISCOSITY *** SUBROUTINE S7A04(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 ************************************************ CPPT DOUBLE PRECISION FUNCTION G17A04(DD, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) AA = G3A04 (DD, T) IF (AA.LT.-1.E08) GO TO 900 BB = G4A04 (DD, T) IF (BB.LT.-1.E08) GO TO 900 CCV = G5A04 (DD, T) IF (CCV.LT.-1.0E08) GO TO 900 F6 = 1.D00/(BB*DD) F4 = T*AA*F6 F5 = AA/(DD*CCV) F3 = 1.+F4*F5 CP = F3*CCV G17A04 = CP RETURN 900 G17A04 = -1.E20 RETURN END C*************F21A04 * CRP(A) QUANTITIES AT THE CRITICAL POINT DOUBLE PRECISION FUNCTION F21A04(A) CHARACTER*1 A,B(5) DATA B/'P','T','V','H','S'/ IF (A.EQ.B(1)) THEN F21A04=0.22746D06 ELSE IF (A.EQ.B(2)) THEN F21A04=5.1953D00 ELSE IF (A.EQ.B(3)) THEN F21A04=1.43596D-02 ELSE IF (A.EQ.B(4)) THEN F21A04=21.948D03 ELSE IF (A.EQ.B(5)) THEN F21A04=5.7685D03 ELSE F21A04=-1.E+20 ENDIF RETURN END C ************************************************ UPT DOUBLE PRECISION FUNCTION G15A04(DD, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX/ 0.8D00, 3.3D00/ IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 900 SS=G6A04(DD, T) IF (SS.LT.-1.0E08) GO TO 900 FF=G8A04(DD, T) IF (FF.LT.-1.0E08) GO TO 900 UU=FF+T*SS G15A04 = UU RETURN 900 G15A04 = -1.E20 RETURN END C ************************************************ PST DOUBLE PRECISION FUNCTION F30A04 (T) C SATURATION PRESSURE AND DP/DT AS A FUNCTION OF TEMPERATURE C ON THE T76 SCALE, FROM 0.5 TO 5.193 K; EQUATION IS GIVEN BY C DURIEUX AND RUSBY, METROLOGIA 19, 67 (1983). C V. ARP, NOV. 14, 1987 IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION C(11,2) DATA EPS5/1.0D-06/ DATA C/ -30.93285D0, 392.47361D0, -2328.04587D0, 1 8111.30347D0, -17809.80901D0, 25766.52747D0, -24601.4D0, 2 14944.65142D0, -5240.36518D0, 807.93168D0, 14.5333D0, 3 -7.41816D0, 5.42128D0, 9.903203D0, -9.617095D0, 4 6.804602D0, -3.0154606D0, 0.7461357D0, -0.0791791D0, 5 0.D0, 0.D0, 0.D0/ C* AT TL, P=5041.8 AND DPDT=12407.9; FIXED POINTS ON THE T76 SCALE DATA TL, TC, TMIN /2.1768D0, 5.1953D0, 0.8D0/ IF (DABS((T-TC)/TC).LE.EPS5) THEN F30A04 = 2.27460D05 RETURN ELSE IF (DABS((T-TMIN)/TMIN).LE.EPS5) THEN F30A04 = 1.47515D00 RETURN ENDIF IF ((T .GT. 5.1953001D00) .OR. (T .LT. 0.799999D00)) THEN F30A04 = -1.E20 RETURN ELSE IF (T .GT. TL) THEN X = T/TC IF (X .GT. 1.D0) X = 1.D0 M = 1 MX = 10 ELSE X = T M = 2 MX = 8 ENDIF Q0 = C(1,M)/X + C(2,M) TN = 1.D0 DO 10 J = 3, MX TN = TN*X 10 Q0 = Q0 + C(J,M)*TN X = 1.D0 - X IF ((M .EQ. 1) .AND. (X .GT. 0.D0)) THEN Y = X**0.9D0 Q0 = Q0 + X*C(11,M)*Y ENDIF F30A04 = DEXP(Q0) RETURN END C ************************************************ TSP DOUBLE PRECISION FUNCTION F40A04(P) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DATA PMAX, PMIN/ 2.2746D05, 1.47515D00/ DATA EPS, EPS5/ 1.0D-8, 1.0D-6/ IF (DABS((P-PMAX)/PMAX).LE.EPS5) THEN XM = 5.1953D00 GO TO 200 ELSE IF (DABS((P-PMIN)/PMIN).LE.EPS5) THEN XM = 0.8D00 GO TO 200 ELSE IF (P.LT.PMIN-1.0D-7.OR.P.GT.PMAX+1.0D-4) THEN F40A04 = -1.E20 RETURN ENDIF X1=0.8D00 X2=5.1953D00 Y=F30A04(X1)-P 10 XM=(X1+X2)/2.D00 IF (Y*(F30A04(XM)-P).GT.0) THEN X1=XM ELSE X2=XM ENDIF IF (DABS((X2-X1)/X1).GE.EPS) GO TO 10 200 F40A04 = XM RETURN END C******************F41A04 * TRPL(A) QUANTITIES AT THE TRIPLE POINT DOUBLE PRECISION FUNCTION F41A04(A) CHARACTER*1 A, B(2) DATA B/'P', 'T'/ IF (A.EQ.B(1)) THEN F41A04=5041.8D00 ELSEIF (A.EQ.B(2)) THEN F41A04=2.1768D00 ELSE F41A04=-1.E+20 ENDIF RETURN END C ************************************************ PLDT DOUBLE PRECISION FUNCTION F66A04 (T) C P = LAMBDA-LINE PRESSURE = 5041.8 PA AT T = 2.1768 K C = 30.134E+5 PA AT T = 1.7673 K C REFERENCE: C H.A. KIERSTEAD, PHYS REV 162, 153 (1967) (T58 SCALE) C REFITTED TO T76 SCALE, V. ARP, JAN. 22, 1988 IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION B(7) DATA B / .42774167D00, -94.820469D00, -85.817089D00, 1 -102.39597D00, -76.735240D00, -.37798315D00, 42.148155D00/ X = T - 2.1768D0 DATA TMIN, TMAX, EPS5/1.7673D00, 2.1768D00, 1.0D-6/ IF (DABS((T-TMAX)/TMAX).LE.EPS5) THEN F66A04 = 5014.8D00 RETURN ELSE IF (DABS((T-TMIN)/TMIN).LE.EPS5) THEN F66A04 = 30.134D05 RETURN ENDIF IF (T.LT.TMIN-1.0D-7.OR.T.GT.TMAX+1.0D-7) THEN F66A04=-1.E20 ELSE P=B(1)+(B(2)+(B(3)+(B(4)+B(5)*X)*X)*X)*X+B(6)*DEXP(B(7)*X) F66A04=P*101325.D00 ENDIF RETURN END C ************************************************ TLDP DOUBLE PRECISION FUNCTION F67A04 (P) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DATA PMAX, PMIN/ 30.134D05, 5041.8D00/ DATA EPS, EPS5/ 1.0D-8, 1.0D-6/ IF (DABS((P-PMAX)/PMAX).LE.EPS5) THEN XM = 1.7673D00 GO TO 200 ELSE IF (DABS((P-PMIN)/PMIN).LE.EPS5) THEN XM = 2.1768D00 GO TO 200 ENDIF IF (P.LT.PMIN-1.0D-3.OR.P.GT.PMAX+1.0D-2) THEN F67A04=-1.E20 RETURN ENDIF X1=1.7673D00 X2=2.1768D00 Y=F66A04(X1)-P 10 XM=(X1+X2)/2.D00 IF (Y*(F66A04(XM)-P).GT.0) THEN X1=XM ELSE X2=XM ENDIF IF (DABS((X2-X1)/X1).GE.EPS) GO TO 10 200 F67A04 = XM RETURN END C ************************************************ PMLT DOUBLE PRECISION FUNCTION F68A04 (TT) C MELTING PRESSURE AS A FUNCTION OF TEMPERATURE, T76 SCALE BELOW 5.1953 C F1 TO F3 HAVE BEEN REFITTED TO T58 EQUATIONS OF GRILLY AND MILLS. C F4 AND F5 HAVE BEEN FITTED TO GRILLY AND MILLS DATA. C V. ARP, FEB 22, 1988 IMPLICIT DOUBLE PRECISION (A-H, O-Z) DATA TMIN, TMAX, EPS5 /0.8D00, 13.8943D00, 1.0D-6/ DATA CONV /98066.5D00/, CON2 /101325.D00/, DELT /0.005D00/ F1 (X) = (-17.80D00+17.31457D00*X**1.555414D00)*CONV F2 (X) = (34.2097D00 + X*(-45.31231D00 + X*(32.26926D00 + X* & (-4.91282D00 + X*0.310795D00)))) * CONV F3 (X) = (17.8537D00 + X*(-15.49444D00 + X*12.57562D00)) * CON2 F4 (X) = (99.3328D00 + X*(-101.44970D00 + X*35.1175D00)) * CON2 F5 (X) = (24.997D00 + 3.36930D00 * (X-0.725D00)**4) * CON2 T = TT IF (DABS((T-TMAX)/TMAX).LE.EPS5) THEN F68A04 = 1000.0D05 RETURN ELSE IF (DABS((T-TMIN)/TMIN).LE.EPS5) THEN F68A04 = 25.328D05 RETURN ELSE IF (T.LT.0.799999D00.OR.T.GT.13.8943001D00) THEN F68A04=-1.E20 RETURN ENDIF C* CUT OFF AT THE UPPER LIMIT OF HEPROP IF (T .GT. 14.010) THEN PMFT = 1013.25E+05 ELSE IF (T .GE. 5.1953 + DELT) THEN PMFT = F1 (T) ELSE IF (T .GT. 5.1953 - DELT) THEN W = (T - 5.1953)/(2.*DELT) + 0.5 PMFT = W*F1(T) + (1.-W)*F2(T) ELSE IF (T .GE. 2.0044 + DELT) THEN PMFT = F2 (T) ELSE IF (T .GT. 2.0044 - DELT) THEN W = (T - 2.0044)/(2.*DELT) + 0.5 PMFT = W*F2(T) + (1.-W)*F3(T) ELSE IF (T .GE. 1.7660 + DELT) THEN PMFT = F3 (T) ELSE IF (T .GT. 1.7660 - DELT) THEN W = (T - 1.7660)/(2.*DELT) + 0.5 PMFT = W*F3(T) + (1.-W)*F4(T) ELSE IF (T .GE. 1.4676 + DELT) THEN PMFT = F4 (T) ELSE IF (T .GT. 1.4676 - DELT) THEN W = (T - 1.4676)/(2.*DELT) + 0.5 PMFT = W*F4(T) + (1.-W)*F5(T) ELSE IF (T .GE. 0.8D0) THEN PMFT = F5 (T) ELSE PMFT = -1.E20 ENDIF F68A04 = PMFT RETURN END C ************************************************ TMLP DOUBLE PRECISION FUNCTION F69A04(P) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DATA PMAX, PMIN/ 1000.0D05, 25.328D05/ DATA EPS, EPS5/ 1.0D-8, 1.0D-6/ IF (DABS((P-PMAX)/PMAX).LE.EPS5) THEN XM = 13.8943D00 GO TO 200 ELSE IF (DABS((P-PMIN)/PMIN).LE.EPS5) THEN XM = 0.8D00 GO TO 200 ELSE IF (P.LT.PMIN-1.0D-2.OR.P.GT.PMAX+1.0D-2) THEN F69A04=-1.E20 RETURN ENDIF X1=0.8D00 X2=13.8943D00 Y=F68A04(X1)-P 10 XM=(X1+X2)/2.D00 IF (Y*(F68A04(XM)-P).GT.0) THEN X1=XM ELSE X2=XM ENDIF IF (DABS((X2-X1)/X1).GE.EPS) GO TO 10 200 F69A04 = XM RETURN END C ************************************************ GAMPT DOUBLE PRECISION FUNCTION F95A04(PA, T) IMPLICIT DOUBLE PRECISION (A-H, O-Z) CCV = F77A04 (PA, T) IF (CCV.LT.-1.0E08) GO TO 920 CP = F18A04(PA, T) IF (CP.LT.-1.0E08) GO TO 920 F95A04 = CP/CCV RETURN 920 F95A04 = -1.E20 RETURN END C ************************************************ G1A04 FUNCTION G1A04 (D, TT) C PRESSURE [PA] AS A FUNCTION OF DENSITY (KG/M3) AND TEMPERATURE (K) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000.D0/D X = V - V0 CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF CALL S2A04 (DT) C I SPECIFIES THE ILOG INDEX IN CV (0 TO 2) C K IS THE COEFFICIENT NUMBER. (1 TO 5) C M SPECIFIES HEI OR HEII (1 OR 2) C N SPECIFIES THE EXPONENT ON T IN CV (3, 5, OR 7) K = 0 PRESS = 0. DO 30 N = 1, 5, 2 R0 = T**(N-1) DO 20 I = 0, 2 IF ((N .EQ. 5) .OR. (I .EQ. 0)) THEN K = K + 1 R = R0 A = 0. DO 10 J = 1, N A = A + DBLE(J) * CL(J,N) * FL(I+J) *R 10 R = R/T PRESS = PRESS + C(K,M) * A ENDIF 20 CONTINUE 30 CONTINUE PRESS = -DTLDV*PRESS + G2A04 (X, T, M) G1A04 = PRESS * 1.D06 END C ************************************************ G2A04 FUNCTION G2A04 (X, T, M) C "BACKGROUND" PRESSURE (MPA) AS A FUNCTION OF V-V0 (CM3/GM), AND T (K). C M SPECIFIES HEI OR HEII. V0 IS THE VOLUME AT THE LOWER LAMBDA POINT. IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(7), Q(5) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 T2 = T*T F(1) = 1. DO 40 J = 1, 5 40 Q(J) = 0. DO 60 K = 1, 6 DO 50 J = 1, 5 50 Q(J) = Q(J) + F(K)*C(6*J+K-1, M) 60 F(K+1) = F(K)*X G2A04 = Q(1)+T2*(Q(2)+T2*(Q(3)+T2*(Q(4)+T2*Q(5)))) END C ************************************************ G3A04 FUNCTION G3A04 (D, TT) C DP/DT AS A FUNCTION OF DENSITY AND TEMPERATURE (SI UNITS) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(7), Q(4) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000./D X = V - V0 IF (V .NE. VSAVE) CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF IF (DT .NE. DTSAVE) CALL S2A04 (DT) T2 = T*T K = 0 DPDT = 0. DO 30 N = 1, 5, 2 R0 = T**(N-1) DO 20 I = 0, 2 IF ((N .EQ. 5) .OR. (I .EQ. 0)) THEN K = K + 1 R = R0 A = 0. DO 10 J = 1, N A = A + CL(J,N) * FL(I+J-1) *R 10 R = R/T DPDT = DPDT + C(K,M) * A ENDIF 20 CONTINUE 30 CONTINUE DPDT = -DTLDV * DPDT F(1) = 1. DO 40 J = 1, 4 40 Q(J) = 0. DO 60 K = 1, 6 DO 50 J = 1, 4 50 Q(J) = Q(J) + F(K)*C(6*J+K+5, M) 60 F(K+1) = F(K)*X DPDT = DPDT+T*(2.*Q(1)+T2*(4.*Q(2)+T2*(6.*Q(3)+T2*8.*Q(4)))) G3A04 = DPDT * 1.D+06 END C ************************************************ G4A04 FUNCTION G4A04 (D, TT) C DP/DD AS A FUNCTION OF DENSITY AND TEMPERATURE (SI UNITS) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(6), Q(5) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000./D X = V - V0 IF (V .NE. VSAVE) CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF IF (DT .NE. DTSAVE) CALL S2A04 (DT) T2 = T*T ASUM = 0. BSUM = 0. K = 0 DO 30 N = 1, 5, 2 R0 = T**(N-1) DO 20 I = 0, 2 IF ((N .EQ. 5) .OR. (I .EQ. 0)) THEN K = K+1 R = R0 A = 0. B = 0. DO 10 J = 1, N E = DBLE(J)*CL(J,N)*R A = A + E*FL(I+J) B = B + E*FL(I+J-1) 10 R = R/T ASUM = ASUM + A*C(K,M) BSUM = BSUM + B*C(K,M) ENDIF 20 CONTINUE 30 CONTINUE DPDV = DTLDV*DTLDV*BSUM - D2TDV2*ASUM F(1) = 1. DO 40 J = 1, 5 40 Q(J) = 0. DO 60 K = 1, 5 FF = F(K)*DBLE(K) DO 50 J = 1, 5 50 Q(J) = Q(J) + FF*C(6*J+K, M) 60 F(K+1) = F(K)*X DPDV = DPDV+Q(1)+T2*(Q(2)+T2*(Q(3)+T2*(Q(4)+T2*Q(5)))) G4A04 = -1.E+09 * DPDV/(D*D) END C ************************************************ G5A04 FUNCTION G5A04 (D, TT) C CV AS A FUNCTION OF DENSITY AND TEMPERATURE (SI UNITS) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(7), Q(4) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000./D X = V - V0 IF (V .NE. VSAVE) CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF IF (DT .NE. DTSAVE) CALL S2A04 (DT) T2 = T*T CV = T * (FL(0)*C(1,M) + T2*(FL(0)*C(2,M) + T2* 1 (FL(0)*C(3,M) + FL(1)*C(4,M) + FL(2)*C(5,M)))) 2 + T*(C(36,M) + T2*(C(37,M) + T2*(C(38,M) + T2*C(39,M)))) F(1) = X DO 10 J = 1, 4 10 Q(J) = 0. DO 30 K = 1, 6 FF = F(K)/DBLE(K) DO 20 J = 1, 4 20 Q(J) = Q(J) + FF*C(6*J+K+5, M) 30 F(K+1) = F(K)*X CV = CV + T*(2.D0*Q(1) + T2*(12.D0*Q(2) + T2* A (30.D0*Q(3) + T2*56.D0*Q(4)))) G5A04 = CV * 1.D+03 END C ************************************************ G6A04 FUNCTION G6A04 (D, TT) C ENTROPY AS A FUNCTION OF DENSITY AND TEMPERATURE (SI UNITS) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000./D X = V - V0 IF (V .NE. VSAVE) CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF IF (DT .NE. DTSAVE) CALL S2A04 (DT) K = 0 ENTROP = 0. DO 30 N = 1, 5, 2 R0 = T**(N-1) DO 20 I = 0, 2 IF ((N .EQ. 5) .OR. (I .EQ. 0)) THEN K = K + 1 R = R0 A = 0. DO 10 J = 1, N A = A + CL(J,N) * FL(I+J) *R 10 R = R/T ENTROP = ENTROP + C(K,M) * A ENDIF 20 CONTINUE 30 CONTINUE ENTROP = ENTROP + G7A04 (X, T, M) G6A04 = ENTROP * 1000. END C ************************************************ G7A04 FUNCTION G7A04 (X, T, M) C"BACKGROUND" ENTROPY (J/GM-K) AS A FUNCTION OF V-V0 (CM3/GM), AND T (K) C M SPECIFIES HEI OR HEII. V0 IS THE VOLUME AT THE LOWER LAMBDA POINT. IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(7), Q(4) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 T2 = T*T T5 = T2*T2*T F(1) = X DO 40 J = 1, 4 40 Q(J) = 0. DO 60 K = 1, 6 FF = F(K)/REAL(K) DO 50 J = 1, 4 50 Q(J) = Q(J) + FF*C(6*J+K+5,M) 60 F(K+1) = F(K) * X G7A04 = T*(2.D0*Q(1) + T2*(4.D0*Q(2) + T2*(6.D0*Q(3) A + T2*8.D0*Q(4)))) + C(36,M)*T + C(37,M)*T2*T/3.D0 B + C(38,M)*T5/5.D0 + C(39,M)*T5*T2/7.D0 + C(40,M) END C ************************************************ G8A04 FUNCTION G8A04 (D, TT) C HELMHOLZ ENERGY AS A FUNCTION OF DENSITY AND TEMPERATURE (SI UNITS) C VALID IN COMPRESSED LIQUID BELOW ABOUT 3 K (THE ARP EQUATION). IMPLICIT DOUBLE PRECISION (A-H, O-Z) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 COMMON /SUBLAM/ TL, DTLDV, D2TDV2, VSAVE, FL(0:10), DTSAVE CALL S3A04 V = 1000./D X = V - V0 IF (V .NE. VSAVE) CALL S1A04 (V) IF (TT .GT. 0.7999D0) THEN T = TT DT = T - TL ELSE DT = TT T = TL + DT ENDIF IF (DT .GT. 0.) THEN M = 1 ELSE M = 2 ENDIF IF (DT .NE. DTSAVE) CALL S2A04 (DT) K = 0 K = 0 HELMH = 0. DO 30 N = 1, 5, 2 R0 = T**(N-1) DO 20 I = 0, 2 IF ((N .EQ. 5) .OR. (I .EQ. 0)) THEN K = K + 1 R = R0 A = 0. DO 10 J = 1, N A = A + DBLE(J) * CL(J,N) * FL(I+J+1) *R 10 R = R/T HELMH = HELMH + C(K,M) * A ENDIF 20 CONTINUE 30 CONTINUE HELMH = HELMH + G9A04 (X, T, M) G8A04 = -1000.*HELMH END C ************************************************ G9A04 FUNCTION G9A04 (X, T, M) C "BACKGROUND" HELMHOLZ ENERGY (J/GM) AS A FUNCTION OF V-V0 (CM3/GM),AND C M SPECIFIES HEI OR HEII. V0 IS THE VOLUME AT THE LOWER LAMBDA POINT. IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION F(7), Q(5) COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 T2 = T*T T4 = T2*T2 F(1) = X DO 40 J = 1, 5 40 Q(J) = 0. DO 60 K = 1, 6 FF = F(K)/REAL(K) DO 50 J = 1, 5 50 Q(J) = Q(J) + FF*C(6*J+K-1, M) 60 F(K+1) = F(K)*X G9A04 = Q(1)+T2*(Q(2)+T2*(Q(3)+T2*(Q(4)+T2*Q(5)))) & + C(36,M)*T2/2.D0 + C(37,M)*T4/12.D0 & + C(38,M)*T4*T2/30.D0+ C(39,M)*T4*T4/56.D0 & + C(40,M)*T + C(41,M) END C ************************************************ G19A04 DOUBLE PRECISION FUNCTION G19A04(PA) IMPLICIT DOUBLE PRECISION (A-H,O-Z) DIMENSION A(8) DATA A/ 8882.60065279D00, - 17126.2418267D00, 9753.63359086D00, & 717.791034485D00, -2429.6722001399D00, 605.812001399D00, & 49.64224399D00, -24.0102845047D00 / P0 = PA/1.D06 T = 10.0 5 F = A(1)+T*(A(2)+T*(A(3)+T*(A(4)+T*(A(5)+T*(A(6) & +T*(A(7)+A(8)*T))))))-P0 DF = A(2)+T*(2.D00*A(3)+T*(3.D00*A(4)+T*(4.D00*A(5) & +T*(5.D00*A(6)+T*(6.D00*A(7)+7.D00*A(8)*T))))) T1 = T-F/DF IF (DABS((T1-T)/T).LT.1.D-7) GO TO 100 T = T1 GO TO 5 100 G19A04 = T1 END C ************************************************ S1A04 SUBROUTINE S1A04 (V) C-----OUTPUT C T = LAMBDA TEMPERATURE (K) (T76 SCALE) C DTDV = DT/DV ( (K-GM)/CM3 ) C D2TDV2 = D2T/DV2 ( (K-GM2)/CM6 ) C-----INPUT C V = SPECIFIC VOLUME [CM3/GM] C-----VERSION 15 JAN 89, V. ARP IMPLICIT DOUBLE PRECISION (A-H, O-Z) COMMON /SUBLAM/ T, DTDV, D2TDV2, VSAVE, FL(11), DELT COMMON /SUBHEC/ C(41,2), CL(8,8), V0, T0 DIMENSION A(5) C FOLLOWING CONSTANTS FROM INDEPENDENT FIT TO KIERSTEAD DATA APRIL 88 C DATA A / 0.9163419802E-01, -0.7663982954E-01, 0.7930218537E-01, C & 0.4802206106E-01, 0.3327932733E-01/ C FOLLOWING CONSTANTS FROM AATZ, JAN 15, 1989 DATA A / 0.91672438D-01, -0.82840336D-01, & 0.71832749D-01, 0.48395170D-01, 0.39159012D-01/ X = V - V0 T = T0 + X*(A(1) + X*(A(2) + X*(A(3) + X*(A(4) + X*A(5))))) DTDV = A(1)+ X*(2.D0*A(2) + X*(3.D0*A(3) + X*(4.D0*A(4) A + X*5.D0*A(5)))) D2TDV2 = 2.D0*A(2)+X*(6.D0*A(3)+X*(12.D0*A(4)+X*20.D0*A(5))) END C ************************************************ S2A04 SUBROUTINE S2A04 (X) C OUTPUT: Y(1) = QUASI-LOGARITHMIC SINGULARITY C Y(I), I>1, = SUCCESSIVE INDEFINITE INTEGRALS C INPUT: X = DELTA-T TO THE LAMBDA LINE C (THE FIRST FOUR VARIABLES IN /SUBLAM/ ARE NOT USED BY THIS SUBROUTINE) IMPLICIT DOUBLE PRECISION (A-H, O-Z) COMMON /SUBLAM/ TL, DTDV, DT2DV2, V, Y(11), XSAVE DATA MAX /11/, ZERO /0.D0/ XSAVE = X IF (X .EQ. ZERO) THEN Y(1) = -100.D0 DO 10 I = 2, MAX 10 Y(I) = ZERO ELSE Z = DABS (X) Y(1) = DLOG (Z) FACT = 1.D0 DO 20 I = 2, MAX IF ((DABS(Y(I-1)) .LT. 1.D-25) .AND. (I .GT. 2)) THEN Y(I) = ZERO ELSE C Y(I) = (X*Y(I-1) - X**(I-1)/FACT) / (DBLE(I-1) - ALFA) Y(I) = (X*Y(I-1) - X**(I-1)/FACT) / DBLE(I-1) ENDIF 20 FACT = FACT*DBLE(I) ENDIF END C ************************************************ S3A04 SUBROUTINE S3A04 IMPLICIT DOUBLE PRECISION (A-H, O-Z) DIMENSION CD1(39),CD2(39),CDL(8,8) COMMON /SUBHEC/ C1(39),S01,A01, C2(39),S02,A02, CL(8,8), V0, T0 C* FOLLOWING FROM CYBER AAUI, JAN 15, 1989; ALFA=0.; BEST FIT TO 3.0 K DATA CD1 / A -0.55524231D+00, 0.18333872D+00,-0.42819388D-01,-0.49962336D-01, B -0.83818656D-01,-0.14081863D+00, 0.89901300D+00, 0.66841845D+01, C 0.98899347D+01, 0.73876336D+01, 0.20130513D+01,-0.25251119D-01, D -0.10649342D+01,-0.35520547D+01,-0.58465160D+01,-0.42352097D+01, E -0.12206228D+01, 0.16905918D-01, 0.21835833D+00, 0.74548499D+00, F 0.11932777D+01, 0.86121784D+00, 0.25519882D+00,-0.13224602D-02, G -0.20161157D-01,-0.65799115D-01,-0.10330654D+00,-0.75388298D-01, H -0.23788871D-01, 0.52810077D-04, 0.66823412D-03, 0.21297216D-02, I 0.32541554D-02, 0.24421110D-02, 0.82971274D-03,-0.19221178D+01, J 0.11203003D+01,-0.16430199D+00,-0.34805651D-02/ C* FOLLOWING CONSTANTS FROM CYBER AATZ, JAN. 15, 1989; ALFA = 0 DATA CD2 / -0.14629809D+01, 0.76365652D+00,-0.11943389D+00, A 0.30525707D-02,-0.10405828D+00,-0.14748127D+00,-0.96308013D+00, B 0.67943869D+00,-0.55276320D+00, 0.19535989D-01,-0.47928840D-01, C 0.75410775D-01,-0.80609070D-01,-0.76641800D+00,-0.35260824D+01, D -0.51006607D+01,-0.20845765D+01,-0.12991643D-01,-0.14392344D-01, E 0.72407645D+00, 0.35487925D+01, 0.49508203D+01, 0.20014366D+01, F 0.12086022D-02, 0.10894370D-01,-0.21841979D+00,-0.10697343D+01, G -0.14541070D+01,-0.56844527D+00,-0.82701229D-04,-0.10295381D-02, H 0.20777040D-01, 0.10085164D+00, 0.13267276D+00, 0.49701540D-01, I 0.95829687D+00,-0.14025436D+01, 0.63551829D+00,-0.55978553D-01/ DATA CDL/1., 0., 0., 0., 0., 0., 0., 0., A 1., -1., 0., 0., 0., 0., 0., 0., B 1., -2., 2., 0., 0., 0., 0., 0., C 1., -3., 6., -6., 0., 0., 0., 0., D 1., -4., 12., -24., 24., 0., 0., 0., E 1., -5., 20., -60., 120., -120., 0., 0., F 1., -6., 30.,-120., 360., -720., 720., 0., G 1., -7., 42.,-210., 840.,-2520., 5040.,-5040./ DATA VD0 ,TD0 /6.842285D0, 2.1768D0/ C* S(146.15,0.8) = 4.515 J/KG-K; A(146.15,0.8) = -0.606 J/KG C AETZ HEII CONSTANTS GIVE: C S(146.15,2.1768) = 1580.0 J/KG-K; A(146.15,2.1768) = -512.5 J/KG C* FOLLOWING ARE FITTED TO AAUI HEI AND AATZ HEII CONSTANTS. DATA SD1,SD2 / 3.63345D0, -.044009D0/ DATA AD1,AD2 /-4.32500D0, -.787770D0/ DO 100 I=1,39 C1(I) = CD1(I) 100 C2(I) = CD2(I) DO 150 I=1,8 DO 150 J=1,8 150 CL(I,J) = CDL(I,J) S01 = SD1 S02 = SD2 A01 = AD1 A02 = AD2 V0 = VD0 T0 = TD0 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