CMake netCDF Compatibility with WPS (#2121)
[WRF.git] / phys / module_cam_physconst.F
blobcdb33ec7cd36bd8cf4058ad75751580c6a3808a5
1 #define WRF_PORT
2 #define MODAL_AERO
3 !------------------------------------------------------------------------
4 ! Based on physconst.F90 from CAM
5 ! Ported to WRF by William.Gustafson@pnl.gov, Nov. 2009
6 ! Updated to version from CESM 1.0.1 Nov. 2010
7 ! Updated to CESM1.0.3 (CAM5.1.01) by Balwinder.Singh@pnnl.gov
8 !------------------------------------------------------------------------
9 module physconst
11    ! Physical constants.  Use CCSM shared values whenever available.
13    use shr_kind_mod, only: r8 => shr_kind_r8
14    use shr_const_mod, only: shr_const_g,      shr_const_stebol, shr_const_tkfrz,  &
15                             shr_const_mwdair, shr_const_rdair,  shr_const_mwwv,   &
16                             shr_const_latice, shr_const_latvap, shr_const_cpdair, &
17                             shr_const_rhofw,  shr_const_cpwv,   shr_const_rgas,   &
18                             shr_const_karman, shr_const_pstd,   shr_const_rhodair,&
19                             shr_const_avogad, shr_const_boltz,  shr_const_cpfw,   &
20                             shr_const_rwv,    shr_const_zvir,   shr_const_pi,     &
21                             shr_const_rearth, shr_const_sday,   shr_const_cday,   &
22                             shr_const_spval         
23    implicit none
24    private
25    public  :: physconst_readnl
26    save
27    ! Constants based off share code or defined in physconst
28    real(r8), public, parameter :: avogad      = shr_const_avogad     ! Avogadro's number (molecules/kmole)
29    real(r8), public, parameter :: boltz       = shr_const_boltz      ! Boltzman's constant (J/K/molecule)
30    real(r8), public, parameter :: cday        = shr_const_cday       ! sec in calendar day ~ sec
31    real(r8), public, parameter :: cpair       = shr_const_cpdair     ! specific heat of dry air (J/K/kg)
32    real(r8), public, parameter :: cpliq       = shr_const_cpfw       ! specific heat of fresh h2o (J/K/kg)
33    real(r8), public, parameter :: karman      = shr_const_karman     ! Von Karman constant
34    real(r8), public, parameter :: latice      = shr_const_latice     ! Latent heat of fusion (J/kg)
35    real(r8), public, parameter :: latvap      = shr_const_latvap     ! Latent heat of vaporization (J/kg)
36    real(r8), public, parameter :: pi          = shr_const_pi         ! 3.14...
37    real(r8), public, parameter :: pstd        = shr_const_pstd       ! Standard pressure (Pascals)
38    real(r8), public, parameter :: r_universal = shr_const_rgas       ! Universal gas constant (J/K/kmol)
39    real(r8), public, parameter :: rhoh2o      = shr_const_rhofw      ! Density of liquid water (STP)
40    real(r8), public, parameter :: spval       = shr_const_spval      !special value 
41    real(r8), public, parameter :: stebol      = shr_const_stebol     ! Stefan-Boltzmann's constant (W/m^2/K^4)
43    real(r8), public, parameter :: c0          = 2.99792458e8_r8      ! Speed of light in a vacuum (m/s)
44    real(r8), public, parameter :: planck      = 6.6260755e-34_r8     ! Planck's constant (J.s)
46    ! Molecular weights
47    real(r8), public, parameter :: mwco2       =  44._r8             ! molecular weight co2
48    real(r8), public, parameter :: mwn2o       =  44._r8             ! molecular weight n2o
49    real(r8), public, parameter :: mwch4       =  16._r8             ! molecular weight ch4
50    real(r8), public, parameter :: mwf11       = 136._r8             ! molecular weight cfc11
51    real(r8), public, parameter :: mwf12       = 120._r8             ! molecular weight cfc12
52    real(r8), public, parameter :: mwo3        =  48._r8             ! molecular weight O3
53    real(r8), public, parameter :: mwso2       =  64._r8
54    real(r8), public, parameter :: mwso4       =  96._r8
55    real(r8), public, parameter :: mwh2o2      =  34._r8
56    real(r8), public, parameter :: mwdms       =  62._r8
59    ! modifiable physical constants for aquaplanet
61    real(r8), public           :: gravit       = shr_const_g     ! gravitational acceleration (m/s**2)
62    real(r8), public           :: sday         = shr_const_sday  ! sec in siderial day ~ sec
63    real(r8), public           :: mwh2o        = shr_const_mwwv  ! molecular weight h2o
64    real(r8), public           :: cpwv         = shr_const_cpwv  ! specific heat of water vapor (J/K/kg)
65    real(r8), public           :: mwdry        = shr_const_mwdair! molecular weight dry air
66    real(r8), public           :: rearth       = shr_const_rearth! radius of earth (m)
67    real(r8), public           :: tmelt        = shr_const_tkfrz ! Freezing point of water (K)
69 !---------------  Variables below here are derived from those above -----------------------
71    real(r8), public           :: rga          = 1._r8/shr_const_g                 ! reciprocal of gravit
72    real(r8), public           :: ra           = 1._r8/shr_const_rearth            ! reciprocal of earth radius
73    real(r8), public           :: omega        = 2.0_R8*shr_const_pi/shr_const_sday! earth rot ~ rad/sec
74    real(r8), public           :: rh2o         = shr_const_rgas/shr_const_mwwv     ! Water vapor gas constant ~ J/K/kg
75    real(r8), public           :: rair         = shr_const_rdair   ! Dry air gas constant     ~ J/K/kg
76    real(r8), public           :: epsilo       = shr_const_mwwv/shr_const_mwdair   ! ratio of h2o to dry air molecular weights 
77    real(r8), public           :: zvir         = (shr_const_rwv/shr_const_rdair)-1.0_R8 ! (rh2o/rair) - 1
78    real(r8), public           :: cpvir        = (shr_const_cpwv/shr_const_cpdair)-1.0_R8 ! CPWV/CPDAIR - 1.0
79    real(r8), public           :: rhodair      = shr_const_pstd/(shr_const_rdair*shr_const_tkfrz)
80    real(r8), public           :: cappa        = (shr_const_rgas/shr_const_mwdair)/shr_const_cpdair  ! R/Cp
81    real(r8), public           :: ez           ! Coriolis expansion coeff -> omega/sqrt(0.375)   
82    real(r8), public           :: Cpd_on_Cpv   = shr_const_cpdair/shr_const_cpwv
83                          
84 !================================================================================================
85 contains
86 !================================================================================================
88    ! Read namelist variables.
89    subroutine physconst_readnl(nlfile)
90 #ifndef WRF_PORT
91       use namelist_utils,  only: find_group_name
92       use units,           only: getunit, freeunit
93       use mpishorthand
94       use spmd_utils,      only: masterproc
95       use abortutils,      only: endrun
96       use cam_logfile,     only: iulog
97 #endif   
98       character(len=*), intent(in) :: nlfile  ! filepath for file containing namelist input
99 #ifndef WRF_PORT
100       ! Local variables
101       integer :: unitn, ierr
102       character(len=*), parameter :: subname = 'physconst_readnl'
103       logical       newg, newsday, newmwh2o, newcpwv, newmwdry, newrearth, newtmelt
105       ! Physical constants needing to be reset (ie. for aqua planet experiments)
106       namelist /physconst_nl/  cpwv, gravit, mwdry, mwh2o, rearth, sday, tmelt
108       !-----------------------------------------------------------------------------
110       if (masterproc) then
111          unitn = getunit()
112          open( unitn, file=trim(nlfile), status='old' )
113          call find_group_name(unitn, 'physconst_nl', status=ierr)
114          if (ierr == 0) then
115             read(unitn, physconst_nl, iostat=ierr)
116             if (ierr /= 0) then
117                call endrun(subname // ':: ERROR reading namelist')
118             end if
119          end if
120          close(unitn)
121          call freeunit(unitn)
122       end if
124 #ifdef SPMD
125       ! Broadcast namelist variables
126       call mpibcast(cpwv,      1,                   mpir8,   0, mpicom)
127       call mpibcast(gravit,    1,                   mpir8,   0, mpicom)
128       call mpibcast(mwdry,     1,                   mpir8,   0, mpicom)
129       call mpibcast(mwh2o,     1,                   mpir8,   0, mpicom)
130       call mpibcast(rearth,    1,                   mpir8,   0, mpicom)
131       call mpibcast(sday,      1,                   mpir8,   0, mpicom)
132       call mpibcast(tmelt,     1,                   mpir8,   0, mpicom)
133 #endif
136       
137       newg     =  gravit .ne. shr_const_g 
138       newsday  =  sday   .ne. shr_const_sday
139       newmwh2o =  mwh2o  .ne. shr_const_mwwv
140       newcpwv  =  cpwv   .ne. shr_const_cpwv
141       newmwdry =  mwdry  .ne. shr_const_mwdair
142       newrearth=  rearth .ne. shr_const_rearth
143       newtmelt =  tmelt  .ne. shr_const_tkfrz
144       
145       
146       
147       if (newg .or. newsday .or. newmwh2o .or. newcpwv .or. newmwdry .or. newrearth .or. newtmelt) then
148          if (masterproc) then
149             write(iulog,*)'****************************************************************************'
150 #ifdef WRF_PORT
151             call wrf_message(iulog)
152 #endif 
153             write(iulog,*)'***    New Physical Constant Values set via namelist                     ***'
154 #ifdef WRF_PORT
155             call wrf_message(iulog)
156 #endif 
157             write(iulog,*)'***                                                                      ***'
158 #ifdef WRF_PORT
159             call wrf_message(iulog)
160 #endif 
161             write(iulog,*)'***    Physical Constant    Old Value                  New Value         ***'
162 #ifdef WRF_PORT
163             call wrf_message(iulog)
164 #endif 
165             if (newg)       write(iulog,*)'***       GRAVITY   ',shr_const_g,gravit,'***'
166             if (newsday)    write(iulog,*)'***       SDAY      ',shr_const_sday,sday,'***'
167             if (newmwh2o)   write(iulog,*)'***       MWH20     ',shr_const_mwwv,mwh2o,'***'
168             if (newcpwv)    write(iulog,*)'***       CPWV      ',shr_const_cpwv,cpwv,'***'
169             if (newmwdry)   write(iulog,*)'***       MWDRY     ',shr_const_mwdair,mwdry,'***'
170             if (newrearth)  write(iulog,*)'***       REARTH    ',shr_const_rearth,rearth,'***'
171             if (newtmelt)   write(iulog,*)'***       TMELT     ',shr_const_tkfrz,tmelt,'***'
172             write(iulog,*)'****************************************************************************'
173 #ifdef WRF_PORT
174             call wrf_message(iulog)
175 #endif 
176          end if
177          rga         = 1._r8/gravit 
178          ra          = 1._r8/rearth
179          omega       = 2.0_R8*pi/sday
180          cpvir       = cpwv/cpair - 1._r8
181          epsilo      = mwh2o/mwdry      
182          
183          !  rair and rh2o have to be defined before any of the variables that use them
184          
185          rair        = r_universal/mwdry
186          rh2o        = r_universal/mwh2o  
187          
188          cappa       = rair/cpair       
189          rhodair     = pstd/(rair*tmelt)
190          zvir        =  (rh2o/rair)-1.0_R8
191          ez          = omega / sqrt(0.375_r8)
192          Cpd_on_Cpv  = cpair/cpwv
193          
194       else      
195          ez          = omega / sqrt(0.375_r8)
196       end if
197 #endif     
198     end subroutine physconst_readnl
199   end module physconst