--- trunk/libf/phylmd/suphec.f90 2008/08/05 13:31:32 17 +++ trunk/libf/phylmd/suphec.f90 2011/08/24 11:43:14 49 @@ -2,167 +2,159 @@ implicit none + ! A1.0 Fundamental constants + REAL RPI + real, parameter:: RCLUM = 299792458. + real, parameter:: RHPLA = 6.6260755E-34 + real, parameter:: RKBOL = 1.380658E-23 ! Boltzmann constant, in J K-1 + real, parameter:: RNAVO = 6.0221367E+23 ! Avogadro number, in mol-1 + + ! A1.1 Astronomical constants + REAL RSIYEA, RSIDAY, ROMEGA + real, parameter:: RDAY = 86400. + real, parameter:: REA = 149597870000. + real, parameter:: REPSM = 0.409093 + + ! A1.2 Geoide + real, parameter:: RG = 9.80665 ! acceleration of gravity, in m s-2 + real, parameter:: RA = 6371229. + + ! A1.3 Radiation + REAL RSIGMA + + ! A1.4 Thermodynamic gas phase + REAL, parameter:: R = RNAVO * RKBOL ! ideal gas constant, in J K-1 mol-1 + real RV, RCPD, RCPV, RCVD, RCVV + real, parameter:: RMD = 28.9644 ! molar mass of dry air, in g mol-1 + + real, parameter:: RD = 1000. * R / RMD + ! specific ideal gas constant for dry air, in J K-1 kg-1 + ! (factor 1000: conversion from g to kg) + + real, parameter:: RMO3 = 47.9942 + real, parameter:: RMV = 18.0153 + REAL RKAPPA, RETV + + ! A1.5, 6 Thermodynamic liquid, solid phases + REAL RCW, RCS + + ! A1.7 Thermodynamic transition of phase + REAL RLMLT + real, parameter:: RTT = 273.16 + real, parameter:: RLVTT = 2.5008E+6 + real, parameter:: RLSTT = 2.8345E+6 + real, parameter:: RATM = 100000. + + ! A1.8 Curve of saturation + REAL RALPW, RBETW, RGAMW, RALPS, RBETS, RGAMS + real, parameter:: RESTT = 611.14 + REAL RALPD, RBETD, RGAMD + + save + contains SUBROUTINE suphec - ! From phylmd/suphec.F,v 1.2 2005/06/06 13:16:33 - + ! From phylmd/suphec.F, v 1.2 2005/06/06 13:16:33 ! Initialise certaines constantes et parametres physiques. - use YOMCST, only: rpi, rclum, rhpla, rkbol, rnavo, rday, rea, repsm, & - rsiyea, rsiday,romega, rg, ra, r1sa, rsigma, r, rmd, rmo3, rmv, rd, & - rv, rcpd, rcvd, rcpv, rcvv, rkappa, retv, rcw, rcs, rtt, rlvtt, & - rlstt, rlmlt, ratm, restt, rgamw, rbetw, ralpw, rgams, rbets, ralps, & - rgamd, rbetd, ralpd - use yoethf + !------------------------------------------ - LOGICAL:: firstcall = .TRUE. + PRINT *, 'Call sequence information: suphec' - !------------------------------------------ + ! 1. DEFINE FUNDAMENTAL CONSTANTS - IF (firstcall) THEN - PRINT *, 'suphec initialise les constantes du GCM' - firstcall = .FALSE. - - !* 1. DEFINE FUNDAMENTAL CONSTANTS. - - WRITE(UNIT=6,FMT='(''0*** Constants of the ICM ***'')') - RPI=2.*ASIN(1.) - RCLUM=299792458. - RHPLA=6.6260755E-34 - RKBOL=1.380658E-23 - RNAVO=6.0221367E+23 - WRITE(UNIT=6,FMT='('' *** Fundamental constants ***'')') - WRITE(UNIT=6,FMT='('' PI = '',E13.7,'' -'')')RPI - WRITE(UNIT=6,FMT='('' c = '',E13.7,''m s-1'')') & - RCLUM - WRITE(UNIT=6,FMT='('' h = '',E13.7,''J s'')') & - RHPLA - WRITE(UNIT=6,FMT='('' K = '',E13.7,''J K-1'')') & - RKBOL - WRITE(UNIT=6,FMT='('' N = '',E13.7,''mol-1'')') & - RNAVO - - !* 2. DEFINE ASTRONOMICAL CONSTANTS. - - RDAY=86400. - REA=149597870000. - REPSM=0.409093 - - RSIYEA=365.25*RDAY*2.*RPI/6.283076 - RSIDAY=RDAY/(1.+RDAY/RSIYEA) - ROMEGA=2.*RPI/RSIDAY - - WRITE(UNIT=6,FMT='('' *** Astronomical constants ***'')') - WRITE(UNIT=6,FMT='('' day = '',E13.7,'' s'')')RDAY - WRITE(UNIT=6,FMT='('' half g. axis = '',E13.7,'' m'')')REA - WRITE(UNIT=6,FMT='('' mean anomaly = '',E13.7,'' -'')')REPSM - WRITE(UNIT=6,FMT='('' sideral year = '',E13.7,'' s'')')RSIYEA - WRITE(UNIT=6,FMT='('' sideral day = '',E13.7,'' s'')')RSIDAY - WRITE(UNIT=6,FMT='('' omega = '',E13.7,'' s-1'')') & - ROMEGA - - !* 3. DEFINE GEOIDE. - - RG=9.80665 - RA=6371229. - R1SA=SNGL(1.D0/DBLE(RA)) - WRITE(UNIT=6,FMT='('' *** Geoide ***'')') - WRITE(UNIT=6,FMT='('' Gravity = '',E13.7,'' m s-2'')') & - RG - WRITE(UNIT=6,FMT='('' Earth radius = '',E13.7,'' m'')')RA - WRITE(UNIT=6,FMT='('' Inverse E.R. = '',E13.7,'' m'')')R1SA - - !* 4. DEFINE RADIATION CONSTANTS. - - rsigma = 2.*rpi**5 * (rkbol/rhpla)**3 * rkbol/rclum/rclum/15. - !IM init. dans conf_phys.F90 RI0=1365. - WRITE(UNIT=6,FMT='('' *** Radiation ***'')') - WRITE(UNIT=6,FMT='('' Stefan-Bol. = '',E13.7,'' W m-2 K-4'')') RSIGMA - - !* 5. DEFINE THERMODYNAMIC CONSTANTS, GAS PHASE. - - R=RNAVO*RKBOL - RMD=28.9644 - RMO3=47.9942 - RMV=18.0153 - RD=1000.*R/RMD - RV=1000.*R/RMV - RCPD=3.5*RD - RCVD=RCPD-RD - RCPV=4. *RV - RCVV=RCPV-RV - RKAPPA=RD/RCPD - RETV=RV/RD-1. - WRITE(UNIT=6,FMT='('' *** Thermodynamic, gas ***'')') - WRITE(UNIT=6,FMT='('' Perfect gas = '',e13.7)') R - WRITE(UNIT=6,FMT='('' Dry air mass = '',e13.7)') RMD - WRITE(UNIT=6,FMT='('' Ozone mass = '',e13.7)') RMO3 - WRITE(UNIT=6,FMT='('' Vapour mass = '',e13.7)') RMV - WRITE(UNIT=6,FMT='('' Dry air cst. = '',e13.7)') RD - WRITE(UNIT=6,FMT='('' Vapour cst. = '',e13.7)') RV - WRITE(UNIT=6,FMT='('' Cpd = '',e13.7)') RCPD - WRITE(UNIT=6,FMT='('' Cvd = '',e13.7)') RCVD - WRITE(UNIT=6,FMT='('' Cpv = '',e13.7)') RCPV - WRITE(UNIT=6,FMT='('' Cvv = '',e13.7)') RCVV - WRITE(UNIT=6,FMT='('' Rd/Cpd = '',e13.7)') RKAPPA - WRITE(UNIT=6,FMT='('' Rv/Rd-1 = '',e13.7)') RETV - - !* 6. DEFINE THERMODYNAMIC CONSTANTS, LIQUID PHASE. - - RCW=RCPV - WRITE(UNIT=6,FMT='('' *** Thermodynamic, liquid ***'')') - WRITE(UNIT=6,FMT='('' Cw = '',E13.7)') RCW - - !* 7. DEFINE THERMODYNAMIC CONSTANTS, SOLID PHASE. - - RCS=RCPV - WRITE(UNIT=6,FMT='('' *** thermodynamic, solid ***'')') - WRITE(UNIT=6,FMT='('' Cs = '',E13.7)') RCS - - !* 8. DEFINE THERMODYNAMIC CONSTANTS, TRANSITION OF PHASE. - - RTT=273.16 - RLVTT=2.5008E+6 - RLSTT=2.8345E+6 - RLMLT=RLSTT-RLVTT - RATM=100000. - WRITE(UNIT=6,FMT='('' *** Thermodynamic, trans. ***'')') - WRITE(UNIT=6,FMT='('' Fusion point = '',E13.7)') RTT - WRITE(UNIT=6,FMT='('' RLvTt = '',E13.7)') RLVTT - WRITE(UNIT=6,FMT='('' RLsTt = '',E13.7)') RLSTT - WRITE(UNIT=6,FMT='('' RLMlt = '',E13.7)') RLMLT - WRITE(UNIT=6,FMT='('' Normal press. = '',E13.7)') RATM - WRITE(UNIT=6,FMT='('' Latent heat : '')') - - !* 9. SATURATED VAPOUR PRESSURE. - - RESTT=611.14 - RGAMW=(RCW-RCPV)/RV - RBETW=RLVTT/RV+RGAMW*RTT - RALPW=LOG(RESTT)+RBETW/RTT+RGAMW*LOG(RTT) - RGAMS=(RCS-RCPV)/RV - RBETS=RLSTT/RV+RGAMS*RTT - RALPS=LOG(RESTT)+RBETS/RTT+RGAMS*LOG(RTT) - RGAMD=RGAMS-RGAMW - RBETD=RBETS-RBETW - RALPD=RALPS-RALPW - - ! calculer les constantes pour les fonctions thermodynamiques - - RVTMP2=RCPV/RCPD-1. - RHOH2O=RATM/100. - R2ES=RESTT*RD/RV - R3LES=17.269 - R3IES=21.875 - R4LES=35.86 - R4IES=7.66 - R5LES=R3LES*(RTT-R4LES) - R5IES=R3IES*(RTT-R4IES) - ELSE - PRINT *, 'suphec DEJA APPELE ' - ENDIF + print *, 'Constants of the ICM' + RPI = 2.*ASIN(1.) + print *, 'Fundamental constants ' + print '('' PI = '', E13.7, '' -'')', RPI + print '('' c = '', E13.7, ''m s-1'')', RCLUM + print '('' h = '', E13.7, ''J s'')', RHPLA + print '('' K = '', E13.7, ''J K-1'')', RKBOL + print '('' N = '', E13.7, ''mol-1'')', RNAVO + + ! 2. DEFINE ASTRONOMICAL CONSTANTS + + RSIYEA = 365.25*RDAY*2.*RPI/6.283076 + RSIDAY = RDAY/(1.+RDAY/RSIYEA) + ROMEGA = 2.*RPI/RSIDAY + + print *, 'Astronomical constants ' + print '('' day = '', E13.7, '' s'')', RDAY + print '('' half g. axis = '', E13.7, '' m'')', REA + print '('' mean anomaly = '', E13.7, '' -'')', REPSM + print '('' sideral year = '', E13.7, '' s'')', RSIYEA + print '('' sideral day = '', E13.7, '' s'')', RSIDAY + print '('' omega = '', E13.7, '' s-1'')', ROMEGA + + ! 3. DEFINE GEOIDE. + + print *, ' Geoide ' + print '('' Gravity = '', E13.7, '' m s-2'')', RG + print '('' Earth radius = '', E13.7, '' m'')', RA + + ! 4. DEFINE RADIATION CONSTANTS. + + rsigma = 2.*rpi**5 * (rkbol/rhpla)**3 * rkbol/rclum/rclum/15. + print *, ' Radiation ' + print '('' Stefan-Bol. = '', E13.7, '' W m-2 K-4'')', RSIGMA + + ! 5. DEFINE THERMODYNAMIC CONSTANTS, GAS PHASE. + + RV = 1000.*R/RMV + RCPD = 3.5*RD + RCVD = RCPD-RD + RCPV = 4. *RV + RCVV = RCPV-RV + RKAPPA = RD/RCPD + RETV = RV/RD-1. + print *, 'Thermodynamic, gas ' + print '('' Perfect gas = '', e13.7)', R + print '('' Ozone mass = '', e13.7)', RMO3 + print '('' Vapour mass = '', e13.7)', RMV + print '('' Dry air constant = '', e13.7)', RD + print '('' Vapour constant = '', e13.7)', RV + print '('' Cpd = '', e13.7)', RCPD + print '('' Cvd = '', e13.7)', RCVD + print '('' Cpv = '', e13.7)', RCPV + print '('' Cvv = '', e13.7)', RCVV + print '('' Rd/Cpd = '', e13.7)', RKAPPA + print '('' Rv/Rd-1 = '', e13.7)', RETV + + ! 6. DEFINE THERMODYNAMIC CONSTANTS, LIQUID PHASE. + + RCW = RCPV + print *, 'Thermodynamic, liquid ' + print '('' Cw = '', E13.7)', RCW + + ! 7. DEFINE THERMODYNAMIC CONSTANTS, SOLID PHASE. + + RCS = RCPV + print *, 'thermodynamic, solid' + print '('' Cs = '', E13.7)', RCS + + ! 8. DEFINE THERMODYNAMIC CONSTANTS, TRANSITION OF PHASE. + + RLMLT = RLSTT-RLVTT + print *, 'Thermodynamic, trans. ' + print '('' Fusion point = '', E13.7)', RTT + print '('' RLvTt = '', E13.7)', RLVTT + print '('' RLsTt = '', E13.7)', RLSTT + print '('' RLMlt = '', E13.7)', RLMLT + print '('' Normal press. = '', E13.7)', RATM + + ! 9. SATURATED VAPOUR PRESSURE. + + RGAMW = (RCW-RCPV)/RV + RBETW = RLVTT/RV+RGAMW*RTT + RALPW = LOG(RESTT)+RBETW/RTT+RGAMW*LOG(RTT) + RGAMS = (RCS-RCPV)/RV + RBETS = RLSTT/RV+RGAMS*RTT + RALPS = LOG(RESTT)+RBETS/RTT+RGAMS*LOG(RTT) + RGAMD = RGAMS-RGAMW + RBETD = RBETS-RBETW + RALPD = RALPS-RALPW END SUBROUTINE suphec