-
Notifications
You must be signed in to change notification settings - Fork 10
/
Copy pathsupersat.F90
107 lines (83 loc) · 3.47 KB
/
supersat.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
! Include shortname defintions, so that the F77 code does not have to be modified to
! reference the CARMA structure.
#include "carma_globaer.h"
!! This routine evaluates supersaturations <supsatl> and <supsati> for all gases.
!!
!! @author Andy Ackerman, Chuck Bardeen
!! @version Dec-1995, Aug-2010
subroutine supersat(carma, cstate, iz, igas, rc)
! types
use carma_precision_mod
use carma_enums_mod
use carma_constants_mod
use carma_types_mod
use carmastate_mod
use carma_mod
implicit none
type(carma_type), intent(in) :: carma !! the carma object
type(carmastate_type), intent(inout) :: cstate !! the carma state object
integer, intent(in) :: iz !! z index
integer, intent(in) :: igas !! gas index
integer, intent(inout) :: rc !! return code, negative indicates failure
! Local declarations
real(kind=f) :: rvap
real(kind=f) :: gc_cgs
real(kind=f) :: alpha
! Calculate vapor pressures.
call vaporp(carma, cstate, iz, igas, rc)
! Define gas constant for this gas
rvap = RGAS / gwtmol(igas)
gc_cgs = gc(iz,igas) / zmet(iz)
supsatl(iz,igas) = (gc_cgs * rvap * t(iz) - pvapl(iz,igas)) / pvapl(iz,igas)
supsati(iz,igas) = (gc_cgs * rvap * t(iz) - pvapi(iz,igas)) / pvapi(iz,igas)
! For subgrid scale clouds, the supersaturation needs to be increased be scaled
! based upon cloud fraction. This approach is similar to Wilson and Ballard (1999),
! except that only the water vapor (no liquid water) is used to determine the available
! water.
!
! NOTE: This assumes that the cloud is an ice cloud.
if (do_incloud) then
alpha = rhcrit(iz) * (1._f - cldfrc(iz)) + cldfrc(iz)
supsatl(iz,igas) = (gc_cgs * rvap * t(iz) - alpha * pvapl(iz,igas)) / pvapl(iz,igas)
supsati(iz,igas) = (gc_cgs * rvap * t(iz) - alpha * pvapi(iz,igas)) / pvapi(iz,igas)
! Limit supersaturation to liquid saturation.
supsatl(iz,igas) = min(supsatl(iz,igas), 0._f)
supsati(iz,igas) = min(supsati(iz,igas), (pvapl(iz,igas) - alpha * pvapi(iz,igas)) / pvapi(iz,igas))
end if
return
end
!! This routine evaluates supersaturations <supsatl> and <supsati> for all gases, but
!! thus version of the routine does not scale the supersaturation based on the cloud
!! fraction. It also assumes that vaporp has already been called.
!!
!! @author Andy Ackerman, Chuck Bardeen
!! @version Dec-1995, Aug-2010
subroutine supersat_nocldf(carma, cstate, iz, igas, ssi, ssl, rc)
! types
use carma_precision_mod
use carma_enums_mod
use carma_constants_mod
use carma_types_mod
use carmastate_mod
use carma_mod
implicit none
type(carma_type), intent(in) :: carma !! the carma object
type(carmastate_type), intent(inout) :: cstate !! the carma state object
integer, intent(in) :: iz !! z index
integer, intent(in) :: igas !! gas index
real(kind=f), intent(out) :: ssl
real(kind=f), intent(out) :: ssi
integer, intent(inout) :: rc !! return code, negative indicates failure
! Local declarations
real(kind=f) :: rvap
real(kind=f) :: gc_cgs
real(kind=f) :: alpha
! Calculate vapor pressures.
call vaporp(carma, cstate, iz, igas, rc)
! Define gas constant for this gas
rvap = RGAS / gwtmol(igas)
gc_cgs = gc(iz,igas) / zmet(iz)
ssl = (gc_cgs * rvap * t(iz) - pvapl(iz,igas)) / pvapl(iz,igas)
ssi = (gc_cgs * rvap * t(iz) - pvapi(iz,igas)) / pvapi(iz,igas)
return
end