source: trunk/LMDZ.GENERIC/libf/phystd/calc_cpp_mugaz.F90 @ 724

Last change on this file since 724 was 716, checked in by rwordsworth, 13 years ago

Mainly updates to radiative transfer and gas management scheme.
Most CIA data now read from standard HITRAN datafiles. For the H2O
continuum, two options have been added: the standard CKD continuum,
and the empirical formula in PPC (Pierrehumbert 2010). Use the toggle
'H2Ocont_simple' in callphys.def to choose.

Note to Martians: I've changed the default values of 'sedimentation' and
'co2cond' in inifis to false. Both these are defined in the standard deftank
callphys.def file, so there shouldn't be any compatibility problems.

File size: 3.7 KB
RevLine 
[253]1      subroutine calc_cpp_mugaz
2
3!==================================================================
4!     Purpose
5!     -------
[538]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.
[253]11!
12!     Authors
13!     -------
14!     Robin Wordsworth (2009)
[589]15!     A. Spiga: make the routine OK with latest changes in rcm1d
[253]16!
17!==================================================================
18
[471]19      use gases_h
[253]20      implicit none
21
[538]22#include "dimensions.h"
23#include "dimphys.h"
[253]24#include "comcstfi.h"
[538]25#include "callkeys.h"
[253]26
[589]27      real cpp_c   
[538]28      real mugaz_c
29
[253]30      integer igas
31
[538]32      cpp_c   = 0.0
33      mugaz_c = 0.0
[253]34
35      do igas=1,ngasmx
36
37         if(igas.eq.vgas)then
38            ! ignore variable gas in cpp calculation
39         else
40            ! all values at 300 K from Engineering Toolbox
41            if(gnom(igas).eq.'CO2')then
[716]42               !cpp_c   = cpp_c   + 0.744*gfrac(igas) ! @ ~210 K (better for Mars conditions)
[538]43               cpp_c   = cpp_c   + 0.846*gfrac(igas)
44               mugaz_c = mugaz_c + 44.01*gfrac(igas)
[253]45            elseif(gnom(igas).eq.'N2_')then
[538]46               cpp_c   = cpp_c   + 1.040*gfrac(igas)
47               mugaz_c = mugaz_c + 28.01*gfrac(igas)
[253]48            elseif(gnom(igas).eq.'H2_')then
[538]49               cpp_c   = cpp_c   + 14.31*gfrac(igas)
50               mugaz_c = mugaz_c + 2.01*gfrac(igas)
[716]51            elseif(gnom(igas).eq.'He_')then
52               cpp_c   = cpp_c   + 5.19*gfrac(igas)
53               mugaz_c = mugaz_c + 4.003*gfrac(igas)
[253]54            elseif(gnom(igas).eq.'H2O')then
[538]55               cpp_c   = cpp_c   + 1.864*gfrac(igas)
56               mugaz_c = mugaz_c + 18.02*gfrac(igas)
[716]57            elseif(gnom(igas).eq.'SO2')then
58               cpp_c   = cpp_c   + 0.64*gfrac(igas)
59               mugaz_c = mugaz_c + 64.066*gfrac(igas)
60            elseif(gnom(igas).eq.'H2S')then
61               cpp_c   = cpp_c   + 1.003*gfrac(igas) ! from wikipedia...
62               mugaz_c = mugaz_c + 34.08*gfrac(igas)
[253]63            elseif(gnom(igas).eq.'CH4')then
[538]64               cpp_c   = cpp_c   + 2.226*gfrac(igas)
65               mugaz_c = mugaz_c + 16.04*gfrac(igas)
[253]66            elseif(gnom(igas).eq.'NH3')then
[538]67               cpp_c   = cpp_c   + 2.175*gfrac(igas)
68               mugaz_c = mugaz_c + 17.03*gfrac(igas)
[253]69               print*,'WARNING, cpp for NH3 may be for liquid'
70            else
71               print*,'Error in calc_cpp_mugaz: Gas species not recognised!'
72               call abort
73            endif
74         endif
75
76      enddo
77
[538]78      cpp_c = 1000.0*cpp_c
[253]79
[538]80      print*,'Cp in calc_cpp_mugaz is ',cpp_c,'J kg^-1 K^-1'
81      print*,'Mg in calc_cpp_mugaz is ',mugaz_c,'amu'
82      print*,'Predefined Cp in physics is ',cpp,'J kg^-1 K^-1'
83      print*,'Predefined Mg in physics is ',mugaz,'amu'
[253]84
[589]85      if (check_cpp_match) then
86         print*,'REQUEST TO CHECK cpp_match :'
87         if((abs(1.-cpp/cpp_c).gt.1.e-6) .or.  &
[588]88              (abs(1.-mugaz/mugaz_c).gt.1.e-6)) then
[589]89            ! Ehouarn: tolerate a small mismatch between computed/stored values
90            print*,'--> Values do not match!'
91            print*,'    Either adjust cpp / mugaz via newstart to calculated values,'
92            print*,'    or set check_cpp_match to .false. in callphys.def.'
[538]93            stop
[589]94         else
95            print*,'--> OK. Settings match composition.'
[538]96         endif
97      endif
[253]98
[589]99      if (.not.force_cpp) then
100          print*,'*** Setting cpp & mugaz to computations in calc_cpp_mugaz.'
101          mugaz = mugaz_c
102          cpp = cpp_c
103      else
104          print*,'*** Setting cpp & mugaz to predefined values.'
105      endif
106
107
[253]108      return
109    end subroutine calc_cpp_mugaz
Note: See TracBrowser for help on using the repository browser.