source: trunk/LMDZ.GENERIC/libf/phystd/watercommon_h.F90 @ 728

Last change on this file since 728 was 728, checked in by jleconte, 13 years ago

18/07/2012 == JL

  • New water cycle scheme:
    • largescale now in F90. Robustness increased by i) including evap inside largescale ii) computing the

condensed water amount iteratively

  • same improvements in moistadj.
  • Water thermodynamical data and saturation curves centralized in module watercommn_h
    • The saturation curves used are now Tetens formula as they are analyticaly inversible (Ts(P)-> Ps(T)).

New saturation curve yields very good agreement with the former one.

  • Saturation curves are now generalized for arbitrary water amount (not just q<<1)
  • The old watersat should be removed soon.
  • The effect of water vapor on total (surface) pressure can be taken into account by setting

mass_redistrib=.true. in callphys.def (routine mass_redistribution inspired from co2_condense in martian
model but with a different scheme as many routines evaporate/condense water vapor).

  • New cloud and precipitation scheme (JL + BC):
    • The default recovery assumption for computing the total cloud fraction has been changed (total random gave too

large cloud fractions). See totalcloudfrac.F90 for details and to change this.

  • Totalcloudfraction now set the total cloud fraction to the fraction of the

optically thickest cloud and totalcloudfrac is thus called in aeropacity.

  • Only the total cloud fraction is used to compute optical depth in aeropacity (no more effective

optical depth with exponential formula).

  • 4 precipitation schemes are now available (see rain.F90 for details). The choice can be made using precip_scheme

in callphys.def. Usage of the more physically based model of Boucher et al 95 (precip_scheme=4) is recommended.
default behavior is set to the former "simple scheme" (precip_scheme=1).

  • See rain.f90 to determine the parameter to be defined in callphys.def as a function of the precipitation scheme used.
  • Physiq.F90 now written in a matricial (more F90) way.
  • Radii (H2O and CO2 cloud particles, aerosols, duts, ...) calculations now centralized in module radii_mod.F90

and work with the new aerosol scheme implemented by Laura K. Some inconsistency may remain in callsedim.


Implementation compiled with ifort and pgf90.
gcm.e runs in Earth and Early Mars case with CO2 and H2O cycle + dust.

File size: 5.0 KB
Line 
1module watercommon_h
2
3      implicit none
4
5      real, parameter :: T_coup = 234.0
6      real, parameter :: T_h2O_ice_liq = 273.16
7      real, parameter :: T_h2O_ice_clouds = T_h2O_ice_liq-15.
8      real, parameter :: mH2O = 18.01528   
9
10      ! benjamin additions
11      real, parameter :: RLVTT = 2.257E+6 ! Latent heat of vaporization (J kg-1)
12      real, parameter :: RLSTT = 2.257E+6 ! 2.591E+6 in reality ! Latent heat of sublimation (J kg-1)
13
14      real, parameter :: RLFTT = 3.334E+5 ! Latent heat of fusion (J kg-1) ! entails an energy sink but better description of albedo
15      real, parameter :: rhowater = 1.0E+3 ! mass of water (kg/m^3)
16      real, parameter :: rhowaterice = 9.2E+2 ! mass of water (kg/m^3)
17      real, parameter :: capcal_h2o_liq = 4181.3 ! specific heat capacity of liquid water J/kg/K
18      real, parameter :: mx_eau_sol = 150 ! mass of water (kg/m^2)
19
20      real, save :: epsi, RCPD, RCPV, RV, RVTMP2
21     
22      contains
23
24     
25!==================================================================
26      subroutine Psat_water(T,p,psat,qsat)
27
28         implicit none
29
30!==================================================================
31!     Purpose
32!     -------
33!     Compute the saturation vapor pressure and mass mixing ratio at saturation (kg/kg)
34!     for a given pressure (Pa) and temperature (K)
35!     Based on the Tetens formula from L.Li physical parametrization manual
36!
37!     Authors
38!     -------
39!     Jeremy Leconte (2012)
40!
41!==================================================================
42
43!        input
44         real, intent(in) :: T, p
45 
46!        output
47         real psat,qsat
48
49! JL12 variables for tetens formula
50         real,parameter :: Pref_solid_liquid=611.14
51         real,parameter :: Trefvaporization=35.86
52         real,parameter :: Trefsublimation=7.66
53         real,parameter :: r3vaporization=17.269
54         real,parameter :: r3sublimation=21.875
55
56! checked vs. old watersat data 14/05/2012 by JL.
57
58         if (T.lt.T_h2O_ice_liq) then ! solid / vapour
59            psat = Pref_solid_liquid*Exp(r3sublimation*(T-T_h2O_ice_liq)/(T-Trefsublimation))
60         else                 ! liquid / vapour
61            psat = Pref_solid_liquid*Exp(r3vaporization*(T-T_h2O_ice_liq)/(T-Trefvaporization))
62         endif
63         if(psat.gt.p) then
64            qsat=1.
65         else
66            qsat=epsi*psat/(p-(1.-epsi)*psat)
67         endif
68         return
69      end subroutine Psat_water
70
71
72
73
74!==================================================================
75      subroutine Lcpdqsat_water(T,p,psat,qsat,dqsat)
76
77         implicit none
78
79!==================================================================
80!     Purpose
81!     -------
82!     Compute L/cp*d (q_sat)/d T
83!     for a given temperature (K)!
84!     Based on the Tetens formula from L.Li physical parametrization manual
85!
86!     Authors
87!     -------
88!     Jeremy Leconte (2012)
89!
90!==================================================================
91
92!        input
93         real T, p, psat, qsat
94 
95!        output
96         real dqsat
97
98! JL12 variables for tetens formula
99         real,parameter :: Pref_solid_liquid=611.14
100         real,parameter :: Trefvaporization=35.86
101         real,parameter :: Trefsublimation=7.66
102         real,parameter :: r3vaporization=17.269
103         real,parameter :: r3sublimation=21.875
104
105         real :: dummy
106
107         if (psat.gt.p) then
108            dqsat=0.
109            return
110         endif
111
112         if (T.lt.T_h2O_ice_liq) then ! solid / vapour
113            dummy = r3sublimation*(T_h2O_ice_liq-Trefsublimation)/(T-Trefsublimation)**2
114         else                 ! liquid / vapour
115            dummy = r3vaporization*(T_h2O_ice_liq-Trefvaporization)/(T-Trefvaporization)**2
116         endif
117
118         dqsat=RLVTT/RCPD*qsat**2*(p/(epsi*psat))*dummy
119         return
120      end subroutine Lcpdqsat_water
121
122
123
124
125!==================================================================
126      subroutine Tsat_water(p,Tsat)
127
128         implicit none
129
130!==================================================================
131!     Purpose
132!     -------
133!     Compute the saturation temperature
134!     for a given pressure (Pa)
135!     Based on the Tetens formula from L.Li physical parametrization manual
136!
137!     Authors
138!     -------
139!     Jeremy Leconte (2012)
140!
141!==================================================================
142
143!        input
144         real p
145 
146!        output
147         real Tsat
148
149! JL12 variables for tetens formula
150         real,parameter :: Pref_solid_liquid=611.14
151         real,parameter :: Trefvaporization=35.86
152         real,parameter :: Trefsublimation=7.66
153         real,parameter :: r3vaporization=17.269
154         real,parameter :: r3sublimation=21.875
155
156         if (p.lt.Pref_solid_liquid) then ! solid / vapour
157            Tsat =(T_h2O_ice_liq*r3sublimation- Trefsublimation*Log(p/Pref_solid_liquid))/(r3sublimation-Log(p/Pref_solid_liquid))
158         else                 ! liquid / vapour
159            Tsat =(T_h2O_ice_liq*r3vaporization- Trefvaporization*Log(p/Pref_solid_liquid))/(r3vaporization-Log(p/Pref_solid_liquid))
160         endif
161
162         return
163      end subroutine Tsat_water
164
165
166end module watercommon_h
Note: See TracBrowser for help on using the repository browser.