source: LMDZ4/trunk/libf/phylmd/suphec.F @ 679

Last change on this file since 679 was 652, checked in by Laurent Fairhead, 20 years ago

Correction du bug sur l'ozone (unite, calcul). On peut retrouver le
bug, qui est présent dans les simulations IPCC 2005, en positionnant
le flag bug_ozone à true dans physiq.def.
MPL/LF
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.0 KB
RevLine 
[524]1!
2! $Header$
3!
4      SUBROUTINE suphec
5C
6#include "YOMCST.h"
7#include "YOETHF.h"
8cIM cf. JLD
9       LOGICAL firstcall
10       SAVE firstcall
11       DATA firstcall /.TRUE./
12       IF (firstcall) THEN
13         PRINT*, 'suphec initialise les constantes du GCM'
14         firstcall = .FALSE.
15       ELSE
16         PRINT*, 'suphec DEJA APPELE '
17         RETURN
18       ENDIF
19C      -----------------------------------------------------------------
20C
21C*       1.    DEFINE FUNDAMENTAL CONSTANTS.
22C              -----------------------------
23C
24      WRITE(UNIT=6,FMT='(''0*** Constants of the ICM   ***'')')
25      RPI=2.*ASIN(1.)
26      RCLUM=299792458.
27      RHPLA=6.6260755E-34
28      RKBOL=1.380658E-23
29      RNAVO=6.0221367E+23
30      WRITE(UNIT=6,FMT='('' *** Fundamental constants ***'')')
31      WRITE(UNIT=6,FMT='(''           PI = '',E13.7,'' -'')')RPI
32      WRITE(UNIT=6,FMT='(''            c = '',E13.7,''m s-1'')')
33     S RCLUM
34      WRITE(UNIT=6,FMT='(''            h = '',E13.7,''J s'')')
35     S RHPLA
36      WRITE(UNIT=6,FMT='(''            K = '',E13.7,''J K-1'')')
37     S RKBOL
38      WRITE(UNIT=6,FMT='(''            N = '',E13.7,''mol-1'')')
39     S RNAVO
40C
41C     ----------------------------------------------------------------
42C
43C*       2.    DEFINE ASTRONOMICAL CONSTANTS.
44C              ------------------------------
45C
46      RDAY=86400.
47      REA=149597870000.
48      REPSM=0.409093
49C
50      RSIYEA=365.25*RDAY*2.*RPI/6.283076
51      RSIDAY=RDAY/(1.+RDAY/RSIYEA)
52      ROMEGA=2.*RPI/RSIDAY
53c
54c exp1      R_ecc = 0.05
55c exp1      R_peri = 102.04
56c exp1      R_incl = 22.5
57c exp1      print*, 'Parametres orbitaux modifies'
58c ref      R_ecc = 0.016724
59c ref      R_peri = 102.04
60c ref      R_incl = 23.5
61c
62cIM 161002 : pour avoir les ctes AMIP II
63cIM 161002   R_ecc = 0.016724
64cIM 161002   R_peri = 102.04
65cIM 161002   R_incl = 23.5
66cIM on mets R_ecc, R_peri, R_incl dans conf_phys.F90
67c     R_ecc = 0.016715
68c     R_peri = 102.7
69c     R_incl = 23.441
70c
71      WRITE(UNIT=6,FMT='('' *** Astronomical constants ***'')')
72      WRITE(UNIT=6,FMT='(''          day = '',E13.7,'' s'')')RDAY
73      WRITE(UNIT=6,FMT='('' half g. axis = '',E13.7,'' m'')')REA
74      WRITE(UNIT=6,FMT='('' mean anomaly = '',E13.7,'' -'')')REPSM
75      WRITE(UNIT=6,FMT='('' sideral year = '',E13.7,'' s'')')RSIYEA
76      WRITE(UNIT=6,FMT='(''  sideral day = '',E13.7,'' s'')')RSIDAY
77      WRITE(UNIT=6,FMT='(''        omega = '',E13.7,'' s-1'')')
78     S                  ROMEGA
79c     write(unit=6,fmt='('' excentricite = '',e13.7,''-'')')R_ecc
80c     write(unit=6,fmt='(''     equinoxe = '',e13.7,''-'')')R_peri
81c     write(unit=6,fmt='(''  inclinaison = '',e13.7,''-'')')R_incl
82C
83C     ------------------------------------------------------------------
84C
85C*       3.    DEFINE GEOIDE.
86C              --------------
87C
88      RG=9.80665
89      RA=6371229.
90      R1SA=SNGL(1.D0/DBLE(RA))
91      WRITE(UNIT=6,FMT='('' ***         Geoide         ***'')')
92      WRITE(UNIT=6,FMT='(''      Gravity = '',E13.7,'' m s-2'')')
93     S      RG
94      WRITE(UNIT=6,FMT='('' Earth radius = '',E13.7,'' m'')')RA
95      WRITE(UNIT=6,FMT='('' Inverse E.R. = '',E13.7,'' m'')')R1SA
96C
97C     -----------------------------------------------------------------
98C
99C*       4.    DEFINE RADIATION CONSTANTS.
100C              ---------------------------
101C
102c z.x.li      RSIGMA=2. * RPI**5 * RKBOL**4 /(15.* RCLUM**2 * RHPLA**3)
103      rsigma = 2.*rpi**5 * (rkbol/rhpla)**3 * rkbol/rclum/rclum/15.
104cIM init. dans conf_phys.F90   RI0=1365.
105      WRITE(UNIT=6,FMT='('' ***        Radiation       ***'')')
106      WRITE(UNIT=6,FMT='('' Stefan-Bol.  = '',E13.7,'' W m-2 K-4''
107     S )')  RSIGMA
108cIM init. dans conf_phys.F90   WRITE(UNIT=6,FMT='('' Solar const. = '',E13.7,'' W m-2'')')
109cIM init. dans conf_phys.F90  S      RI0
110C
111C     -----------------------------------------------------------------
112C
113C*       5.    DEFINE THERMODYNAMIC CONSTANTS, GAS PHASE.
114C              ------------------------------------------
115C
116      R=RNAVO*RKBOL
117      RMD=28.9644
[652]118      RMO3=47.9942
[524]119      RMV=18.0153
120      RD=1000.*R/RMD
121      RV=1000.*R/RMV
122      RCPD=3.5*RD
123      RCVD=RCPD-RD
124      RCPV=4. *RV
125      RCVV=RCPV-RV
126      RKAPPA=RD/RCPD
127      RETV=RV/RD-1.
128      WRITE(UNIT=6,FMT='('' *** Thermodynamic, gas     ***'')')
129      WRITE(UNIT=6,FMT='('' Perfect gas  = '',e13.7)') R
130      WRITE(UNIT=6,FMT='('' Dry air mass = '',e13.7)') RMD
[652]131      WRITE(UNIT=6,FMT='('' Ozone   mass = '',e13.7)') RMO3
[524]132      WRITE(UNIT=6,FMT='('' Vapour  mass = '',e13.7)') RMV
133      WRITE(UNIT=6,FMT='('' Dry air cst. = '',e13.7)') RD
134      WRITE(UNIT=6,FMT='('' Vapour  cst. = '',e13.7)') RV
135      WRITE(UNIT=6,FMT='(''         Cpd  = '',e13.7)') RCPD
136      WRITE(UNIT=6,FMT='(''         Cvd  = '',e13.7)') RCVD
137      WRITE(UNIT=6,FMT='(''         Cpv  = '',e13.7)') RCPV
138      WRITE(UNIT=6,FMT='(''         Cvv  = '',e13.7)') RCVV
139      WRITE(UNIT=6,FMT='(''      Rd/Cpd  = '',e13.7)') RKAPPA
140      WRITE(UNIT=6,FMT='(''     Rv/Rd-1  = '',e13.7)') RETV
141C
142C     ----------------------------------------------------------------
143C
144C*       6.    DEFINE THERMODYNAMIC CONSTANTS, LIQUID PHASE.
145C              ---------------------------------------------
146C
147      RCW=RCPV
148      WRITE(UNIT=6,FMT='('' *** Thermodynamic, liquid  ***'')')
149      WRITE(UNIT=6,FMT='(''         Cw   = '',E13.7)') RCW
150C
151C     ----------------------------------------------------------------
152C
153C*       7.    DEFINE THERMODYNAMIC CONSTANTS, SOLID PHASE.
154C              --------------------------------------------
155C
156      RCS=RCPV
157      WRITE(UNIT=6,FMT='('' *** thermodynamic, solid   ***'')')
158      WRITE(UNIT=6,FMT='(''         Cs   = '',E13.7)') RCS
159C
160C     ----------------------------------------------------------------
161C
162C*       8.    DEFINE THERMODYNAMIC CONSTANTS, TRANSITION OF PHASE.
163C              ----------------------------------------------------
164C
165      RTT=273.16
166      RLVTT=2.5008E+6
167      RLSTT=2.8345E+6
168      RLMLT=RLSTT-RLVTT
169      RATM=100000.
170      WRITE(UNIT=6,FMT='('' *** Thermodynamic, trans.  ***'')')
171      WRITE(UNIT=6,FMT='('' Fusion point  = '',E13.7)') RTT
172      WRITE(UNIT=6,FMT='(''        RLvTt  = '',E13.7)') RLVTT
173      WRITE(UNIT=6,FMT='(''        RLsTt  = '',E13.7)') RLSTT
174      WRITE(UNIT=6,FMT='(''        RLMlt  = '',E13.7)') RLMLT
175      WRITE(UNIT=6,FMT='('' Normal press. = '',E13.7)') RATM
176      WRITE(UNIT=6,FMT='('' Latent heat :  '')')
177C
178C     ----------------------------------------------------------------
179C
180C*       9.    SATURATED VAPOUR PRESSURE.
181C              --------------------------
182C
183      RESTT=611.14
184      RGAMW=(RCW-RCPV)/RV
185      RBETW=RLVTT/RV+RGAMW*RTT
186      RALPW=LOG(RESTT)+RBETW/RTT+RGAMW*LOG(RTT)
187      RGAMS=(RCS-RCPV)/RV
188      RBETS=RLSTT/RV+RGAMS*RTT
189      RALPS=LOG(RESTT)+RBETS/RTT+RGAMS*LOG(RTT)
190      RGAMD=RGAMS-RGAMW
191      RBETD=RBETS-RBETW
192      RALPD=RALPS-RALPW
193C
194C     ------------------------------------------------------------------
195c
196c calculer les constantes pour les fonctions thermodynamiques
197c
198      RVTMP2=RCPV/RCPD-1.
199      RHOH2O=RATM/100.
200      R2ES=RESTT*RD/RV
201      R3LES=17.269
202      R3IES=21.875
203      R4LES=35.86
204      R4IES=7.66
205      R5LES=R3LES*(RTT-R4LES)
206      R5IES=R3IES*(RTT-R4IES)
207C
208      RETURN
209      END
Note: See TracBrowser for help on using the repository browser.