source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/config/ppsrc/phys/calc_cpp_mugaz.f90 @ 224

Last change on this file since 224 was 224, checked in by ymipsl, 10 years ago
File size: 10.1 KB
Line 
1      subroutine calc_cpp_mugaz
2
3!==================================================================
4!     Purpose
5!     -------
6!     Check to see if the atmospheric specific heat capacity and
7!     mean molar mass for the gas mixture defined in gases.def
8!     corresponds to what we're using. If it doesn't, abort run
9!     unless option 'check_cpp_match' is set to false in
10!     callphys.def.
11!
12!     Authors
13!     -------
14!     Robin Wordsworth (2009)
15!     A. Spiga: make the routine OK with latest changes in rcm1d
16!
17!==================================================================
18
19      use gases_h
20      implicit none
21
22!-----------------------------------------------------------------------
23!   INCLUDE 'dimensions.h'
24!
25!   dimensions.h contient les dimensions du modele
26!   ndm est tel que iim=2**ndm
27!-----------------------------------------------------------------------
28
29      INTEGER iim,jjm,llm,ndm
30
31      PARAMETER (iim= 128,jjm=96,llm=64,ndm=1)
32
33!-----------------------------------------------------------------------
34!-----------------------------------------------------------------------
35!   INCLUDE 'dimphys.h'
36
37! ngridmx : number of horizontal grid points
38! note: the -1/jjm term will be 0; unless jj=1
39      integer, parameter :: ngridmx = (2+(jjm-1)*iim - 1/jjm)   
40! nlayermx : number of atmospheric layers
41      integer, parameter :: nlayermx = llm 
42! nsoilmx : number of subterranean layers ! nsoilmx is now in comsoil_h
43      !integer, parameter :: nsoilmx = 4 ! for a test
44      !integer, parameter :: nsoilmx = 18 ! for z1=0.0002 m, depth = 18 m => mars case
45      !integer, parameter :: nsoilmx = 13 ! for z1=0.03 m, depth = 104.8 m => earth case
46!-----------------------------------------------------------------------
47
48!-----------------------------------------------------------------------
49! INCLUDE "comcstfi.h"
50
51      common/comcstfi/pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg
52      common/comcstfi/avocado!,molrad,visc
53     
54      real pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg
55      real avocado!,molrad,visc
56
57!
58! For Fortran 77/Fortran 90 compliance always use line continuation
59! symbols '&' in columns 73 and 6
60!
61! Group commons according to their type for minimal performance impact
62
63      COMMON/callkeys_l/callrad,corrk,calldifv,UseTurbDiff,calladj      &
64     &   , co2cond,callsoil                                             &
65     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
66     &   , callstats,calleofdump                                        &
67     &   , enertest                                                     &
68     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
69     &   , radfixed                                                     &
70     &   , meanOLR, specOLR                                             &
71     &   , kastprof                                                     &
72     &   , nosurf, oblate                                               &     
73     &   , newtonian, testradtimes                                      &
74     &   , check_cpp_match, force_cpp                                   &
75     &   , rayleigh                                                     &
76     &   , stelbbody                                                    &
77     &   , nearco2cond                                                  &
78     &   , tracer, mass_redistrib, varactive, varfixed                  &
79     &   , sedimentation,water,watercond,waterrain                      &
80     &   , aeroco2,aeroh2o,aeroh2so4,aeroback2lay                       &
81     &   , aerofixco2,aerofixh2o                                        &
82     &   , hydrology, sourceevol                                        &
83     &   , CLFvarying                                                   &
84     &   , strictboundcorrk                                             &                                       
85     &   , ok_slab_ocean                                                &
86     &   , ok_slab_sic                                                  &
87     &   , ok_slab_heat_transp                                         
88
89
90      COMMON/callkeys_i/iaervar,iddist,iradia,startype
91     
92      COMMON/callkeys_r/topdustref,Nmix_co2,dusttau,Fat1AU,stelTbb,     &
93     &                  Tstrat,tplanet,obs_tau_col_tropo,               &
94     &                  obs_tau_col_strato,pres_bottom_tropo,           &
95     &                  pres_top_tropo,pres_bottom_strato,              &
96     &                  pres_top_strato,size_tropo,size_strato,satval,  &
97     &                  CLFfixval,n2mixratio,co2supsat,pceil,albedosnow,&
98     &                  maxicethick,Tsaldiff,tau_relax,cloudlvl,        &
99     &                  icetstep,intheat,flatten,Rmean,J2,MassPlanet
100     
101      logical callrad,corrk,calldifv,UseTurbDiff                        &
102     &   , calladj,co2cond,callsoil                                     &
103     &   , season,diurnal,tlocked,rings_shadow,lwrite                   &
104     &   , callstats,calleofdump                                        &
105     &   , callgasvis,continuum,H2Ocont_simple,graybody                 &
106     &   , strictboundcorrk                                             
107
108      logical enertest
109      logical nonideal
110      logical meanOLR
111      logical specOLR
112      logical kastprof
113      logical newtonian
114      logical check_cpp_match
115      logical force_cpp
116      logical testradtimes
117      logical rayleigh
118      logical stelbbody
119      logical ozone
120      logical nearco2cond
121      logical tracer
122      logical mass_redistrib
123      logical varactive
124      logical varfixed
125      logical radfixed
126      logical sedimentation
127      logical water,watercond,waterrain
128      logical aeroco2,aeroh2o,aeroh2so4,aeroback2lay
129      logical aerofixco2,aerofixh2o
130      logical hydrology
131      logical sourceevol
132      logical CLFvarying
133      logical nosurf
134      logical oblate
135      logical ok_slab_ocean
136      logical ok_slab_sic
137      logical ok_slab_heat_transp
138
139      integer iddist
140      integer iaervar
141      integer iradia
142      integer startype
143
144      real topdustref
145      real Nmix_co2
146      real dusttau
147      real Fat1AU
148      real stelTbb
149      real Tstrat
150      real tplanet
151      real obs_tau_col_tropo
152      real obs_tau_col_strato
153      real pres_bottom_tropo
154      real pres_top_tropo
155      real pres_bottom_strato
156      real pres_top_strato
157      real size_tropo
158      real size_strato
159      real satval
160      real CLFfixval
161      real n2mixratio
162      real co2supsat
163      real pceil
164      real albedosnow
165      real maxicethick
166      real Tsaldiff
167      real tau_relax
168      real cloudlvl
169      real icetstep
170      real intheat
171      real flatten
172      real Rmean
173      real J2
174      real MassPlanet
175
176      real cpp_c   
177      real mugaz_c
178
179      integer igas
180
181      cpp_c   = 0.0
182      mugaz_c = 0.0
183
184      do igas=1,ngasmx
185
186         if(igas.eq.vgas)then
187            ! ignore variable gas in cpp calculation
188         else
189            ! all values at 300 K from Engineering Toolbox
190            if(igas.eq.igas_CO2)then
191               !cpp_c   = cpp_c   + 0.744*gfrac(igas) ! @ ~210 K (better for Mars conditions)
192               cpp_c   = cpp_c   + 0.846*gfrac(igas)
193               mugaz_c = mugaz_c + 44.01*gfrac(igas)
194            elseif(igas.eq.igas_N2)then
195               cpp_c   = cpp_c   + 1.040*gfrac(igas)
196               mugaz_c = mugaz_c + 28.01*gfrac(igas)
197            elseif(igas.eq.igas_H2)then
198               cpp_c   = cpp_c   + 14.31*gfrac(igas)
199               mugaz_c = mugaz_c + 2.01*gfrac(igas)
200            elseif(igas.eq.igas_He)then
201               cpp_c   = cpp_c   + 5.19*gfrac(igas)
202               mugaz_c = mugaz_c + 4.003*gfrac(igas)
203            elseif(igas.eq.igas_H2O)then
204               cpp_c   = cpp_c   + 1.864*gfrac(igas)
205               mugaz_c = mugaz_c + 18.02*gfrac(igas)
206            elseif(igas.eq.igas_SO2)then
207               cpp_c   = cpp_c   + 0.64*gfrac(igas)
208               mugaz_c = mugaz_c + 64.066*gfrac(igas)
209            elseif(igas.eq.igas_H2S)then
210               cpp_c   = cpp_c   + 1.003*gfrac(igas) ! from wikipedia...
211               mugaz_c = mugaz_c + 34.08*gfrac(igas)
212            elseif(igas.eq.igas_CH4)then
213               cpp_c   = cpp_c   + 2.226*gfrac(igas)
214               mugaz_c = mugaz_c + 16.04*gfrac(igas)
215            elseif(igas.eq.igas_NH3)then
216               cpp_c   = cpp_c   + 2.175*gfrac(igas)
217               mugaz_c = mugaz_c + 17.03*gfrac(igas)
218               print*,'WARNING, cpp for NH3 may be for liquid'
219            elseif(igas.eq.igas_C2H6)then 
220               ! C2H6 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=28
221               cpp_c   = cpp_c   + 1.763*gfrac(igas)
222               mugaz_c = mugaz_c + 30.07*gfrac(igas)
223            elseif(igas.eq.igas_C2H2)then
224               ! C2H2 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=1
225               cpp_c   = cpp_c   + 1.575*gfrac(igas)
226               mugaz_c = mugaz_c + 26.04*gfrac(igas)
227            else
228               print*,'Error in calc_cpp_mugaz: Gas species not recognised!'
229               call abort
230            endif
231         endif
232
233      enddo
234
235      cpp_c = 1000.0*cpp_c
236
237      print*,'Cp in calc_cpp_mugaz is ',cpp_c,'J kg^-1 K^-1'
238      print*,'Mg in calc_cpp_mugaz is ',mugaz_c,'amu'
239      print*,'Predefined Cp in physics is ',cpp,'J kg^-1 K^-1'
240      print*,'Predefined Mg in physics is ',mugaz,'amu'
241
242      if (check_cpp_match) then
243         print*,'REQUEST TO CHECK cpp_match :'
244         if((abs(1.-cpp/cpp_c).gt.1.e-6) .or.  &
245              (abs(1.-mugaz/mugaz_c).gt.1.e-6)) then
246            ! Ehouarn: tolerate a small mismatch between computed/stored values
247            print*,'--> Values do not match!'
248            print*,'    Either adjust cpp / mugaz via newstart to calculated values,'
249            print*,'    or set check_cpp_match to .false. in callphys.def.'
250            stop
251         else
252            print*,'--> OK. Settings match composition.'
253         endif
254      endif
255
256      if (.not.force_cpp) then
257          print*,'*** Setting cpp & mugaz to computations in calc_cpp_mugaz.'
258          mugaz = mugaz_c
259          cpp = cpp_c
260      else
261          print*,'*** Setting cpp & mugaz to predefined values.'
262      endif
263
264
265      return
266    end subroutine calc_cpp_mugaz
Note: See TracBrowser for help on using the repository browser.