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

Last change on this file since 1157 was 869, checked in by aslmd, 12 years ago

24/01/2013 == AS + JL

A more robust way to refer to gas type.

  • Gas names with an arbitrary number of characters (<20) can be used This is good for C2H2, C2H6, H2SO4, C17H21NO4, etc... !!! Remember this must be compliant with Q.dat in corrk_data !!!
  • igas_... labels are assigned once for all in su_gases Then using igas_... everywhere instead of gnom (except for kcm stuff)
  • Users can still use e.g. H2_ but H2 also works
  • Simplified condense_cloud so that igas_CO2 is used directly
File size: 4.2 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#include "dimensions.h"
23#include "dimphys.h"
24#include "comcstfi.h"
25#include "callkeys.h"
26
27      real cpp_c   
28      real mugaz_c
29
30      integer igas
31
32      cpp_c   = 0.0
33      mugaz_c = 0.0
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(igas.eq.igas_CO2)then
42               !cpp_c   = cpp_c   + 0.744*gfrac(igas) ! @ ~210 K (better for Mars conditions)
43               cpp_c   = cpp_c   + 0.846*gfrac(igas)
44               mugaz_c = mugaz_c + 44.01*gfrac(igas)
45            elseif(igas.eq.igas_N2)then
46               cpp_c   = cpp_c   + 1.040*gfrac(igas)
47               mugaz_c = mugaz_c + 28.01*gfrac(igas)
48            elseif(igas.eq.igas_H2)then
49               cpp_c   = cpp_c   + 14.31*gfrac(igas)
50               mugaz_c = mugaz_c + 2.01*gfrac(igas)
51            elseif(igas.eq.igas_He)then
52               cpp_c   = cpp_c   + 5.19*gfrac(igas)
53               mugaz_c = mugaz_c + 4.003*gfrac(igas)
54            elseif(igas.eq.igas_H2O)then
55               cpp_c   = cpp_c   + 1.864*gfrac(igas)
56               mugaz_c = mugaz_c + 18.02*gfrac(igas)
57            elseif(igas.eq.igas_SO2)then
58               cpp_c   = cpp_c   + 0.64*gfrac(igas)
59               mugaz_c = mugaz_c + 64.066*gfrac(igas)
60            elseif(igas.eq.igas_H2S)then
61               cpp_c   = cpp_c   + 1.003*gfrac(igas) ! from wikipedia...
62               mugaz_c = mugaz_c + 34.08*gfrac(igas)
63            elseif(igas.eq.igas_CH4)then
64               cpp_c   = cpp_c   + 2.226*gfrac(igas)
65               mugaz_c = mugaz_c + 16.04*gfrac(igas)
66            elseif(igas.eq.igas_NH3)then
67               cpp_c   = cpp_c   + 2.175*gfrac(igas)
68               mugaz_c = mugaz_c + 17.03*gfrac(igas)
69               print*,'WARNING, cpp for NH3 may be for liquid'
70            elseif(igas.eq.igas_C2H6)then
71               ! C2H6 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=28
72               cpp_c   = cpp_c   + 1.763*gfrac(igas)
73               mugaz_c = mugaz_c + 30.07*gfrac(igas)
74            elseif(igas.eq.igas_C2H2)then
75               ! C2H2 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=1
76               cpp_c   = cpp_c   + 1.575*gfrac(igas)
77               mugaz_c = mugaz_c + 26.04*gfrac(igas)
78            else
79               print*,'Error in calc_cpp_mugaz: Gas species not recognised!'
80               call abort
81            endif
82         endif
83
84      enddo
85
86      cpp_c = 1000.0*cpp_c
87
88      print*,'Cp in calc_cpp_mugaz is ',cpp_c,'J kg^-1 K^-1'
89      print*,'Mg in calc_cpp_mugaz is ',mugaz_c,'amu'
90      print*,'Predefined Cp in physics is ',cpp,'J kg^-1 K^-1'
91      print*,'Predefined Mg in physics is ',mugaz,'amu'
92
93      if (check_cpp_match) then
94         print*,'REQUEST TO CHECK cpp_match :'
95         if((abs(1.-cpp/cpp_c).gt.1.e-6) .or.  &
96              (abs(1.-mugaz/mugaz_c).gt.1.e-6)) then
97            ! Ehouarn: tolerate a small mismatch between computed/stored values
98            print*,'--> Values do not match!'
99            print*,'    Either adjust cpp / mugaz via newstart to calculated values,'
100            print*,'    or set check_cpp_match to .false. in callphys.def.'
101            stop
102         else
103            print*,'--> OK. Settings match composition.'
104         endif
105      endif
106
107      if (.not.force_cpp) then
108          print*,'*** Setting cpp & mugaz to computations in calc_cpp_mugaz.'
109          mugaz = mugaz_c
110          cpp = cpp_c
111      else
112          print*,'*** Setting cpp & mugaz to predefined values.'
113      endif
114
115
116      return
117    end subroutine calc_cpp_mugaz
Note: See TracBrowser for help on using the repository browser.