-
Notifications
You must be signed in to change notification settings - Fork 3
/
HMC_Module_Phys_Retention.f90
162 lines (131 loc) · 7.8 KB
/
HMC_Module_Phys_Retention.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
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
!------------------------------------------------------------------------------------------
! File: HMC_Module_Phys_Retention.f90
! Author: fabio
!
! Created on April 2, 2014, 5:19 PM
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Module Header
module HMC_Module_Phys_Retention
!------------------------------------------------------------------------------------------
! External module(s) for all subroutine in this module
use HMC_Module_Namelist, only: oHMC_Namelist
use HMC_Module_Vars_Loader, only: oHMC_Vars
use HMC_Module_Tools_Debug
! Implicit none for all subroutines in this module
implicit none
!------------------------------------------------------------------------------------------
contains
!------------------------------------------------------------------------------------------
! Subroutine for calculating volume retention
subroutine HMC_Phys_Retention_Cpl(iID, iRows, iCols)
!------------------------------------------------------------------------------------------
! Variable(s) declaration
integer(kind = 4) :: iID
integer(kind = 4) :: iRows, iCols
integer(kind = 4) :: iFlagLAI
real(kind = 4), dimension (iRows, iCols) :: a2dVarRain
real(kind = 4), dimension (iRows, iCols) :: a2dVarLAI
integer(kind = 4), dimension (iRows, iCols) :: a2iVarMask
real(kind = 4), dimension (iRows, iCols) :: a2dVarDEM, a2dVarS
real(kind = 4), dimension (iRows, iCols) :: a2dVarVRet, a2dVarVRetMax, a2dVarVRetStep
real(kind = 4), dimension (iRows, iCols) :: a2dVarRainStep
real(kind = 4), dimension (iRows, iCols) :: a2dVarRetErr
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Initialize variable(s)
iFlagLAI = 0
a2dVarDEM = 0.0; a2dVarS = 0.0;
a2dVarRain = 0.0; a2dVarRainStep = 0.0; a2dVarLAI = 0.0;
a2dVarVRet = 0.0; a2dVarVRetStep = 0.0; a2dVarRetErr = 0.0; a2dVarVRetMax = 0.0
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Loading flag(s)
iFlagLAI = oHMC_Namelist(iID)%iFlagLAI
! Loading static data
a2iVarMask = oHMC_Vars(iID)%a2iMask
a2dVarDEM = oHMC_Vars(iID)%a2dDem
a2dVarS = oHMC_Vars(iID)%a2dS
! Loading dynamic data
a2dVarVRet = oHMC_Vars(iID)%a2dVRet
a2dVarRain = oHMC_Vars(iID)%a2dRain
a2dVarLAI = oHMC_Vars(iID)%a2dLAI
! Temporary variable(s) for retention subroutine
a2dVarRainStep = a2dVarRain
a2dVarVRetStep = a2dVarVRet
! Info start
call mprintf(.true., iINFO_Verbose, ' Phys :: Retention ... ' )
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Debug
if (iDEBUG.gt.0) then
call mprintf(.true., iINFO_Extra, ' ========= RETENTION START =========== ')
call mprintf(.true., iINFO_Extra, checkvar(a2dVarRainStep, a2iVarMask, 'RAIN START ') )
call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRetStep, a2iVarMask, 'VRET START ') )
call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRetMax, a2iVarMask, 'VRETMAX START ') )
call mprintf(.true., iINFO_Extra, checkvar(a2dVarS, a2iVarMask, 'S ') )
call mprintf(.true., iINFO_Extra, '')
endif
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Calculating volume retention max and retrieving Vret and Rain
if (iFlagLAI .eq. 0) then
! Without LAI data (all LAI values are undefined)
where (a2dVarDEM.gt.0.0)
a2dVarVRetMax = (0.038*a2dVarS + 0.4909)*1 ! magic numbers :)
elsewhere
a2dVarVRetMax = 0.0
endwhere
elseif (iFlagLAI .eq. 1) then
! With LAI data
where (a2dVarDEM.gt.0.0.and.a2dVarLAI.gt.0.0)
a2dVarVRetMax = 0.95 + 0.5*a2dVarLAI-0.006*a2dVarLAI**2
elsewhere
a2dVarVRetMax = 0.0
endwhere
where (a2dVarVRetMax.gt.50.0) a2dVarVRetMax = 0.0
else
! Exit error
call mprintf(.true., iERROR, ' Phys :: Retention :: Incorrect LAI type definition!' )
endif
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Calculating rain and retention volume
where( (a2dVarVRet + a2dVarRain) .le. a2dVarVRetMax .and. a2dVarDEM.gt.0.0 )
a2dVarVRet = a2dVarVRet + a2dVarRain
a2dVarRain = 0.0
endwhere
where( (a2dVarVRet + a2dVarRain) .gt. a2dVarVRetMax .and. a2dVarDEM.gt.0.0 )
a2dVarRain = (a2dVarVRet + a2dVarRain) - a2dVarVRetMax
a2dVarVRet = a2dVarVRetMax
endwhere
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Calculating retention subroutine error
a2dVarRetErr = a2dVarRainStep - (a2dVarVRet - a2dVarVRetStep) - a2dVarRain
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Checking for volume retention less then 0.0
where (a2dVarVRet.lt.0.0) a2dVarVRet = 0.0
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Debug
if (iDEBUG.gt.0) then
call mprintf(.true., iINFO_Extra, '')
call mprintf(.true., iINFO_Extra, checkvar(a2dVarRain, a2iVarMask, 'RAIN END ') )
call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRet, a2iVarMask, 'VRET END ') )
call mprintf(.true., iINFO_Extra, checkvar(a2dVarVRetMax, a2iVarMask, 'VRETMAX END ') )
call mprintf(.true., iINFO_Extra, ' ========= RETENTION END =========== ')
endif
!------------------------------------------------------------------------------------------
!------------------------------------------------------------------------------------------
! Passing variable(s) to global declaration
oHMC_Vars(iID)%a2dVRet = a2dVarVRet
oHMC_Vars(iID)%a2dRain = a2dVarRain
! Info end
call mprintf(.true., iINFO_Verbose, ' Phys :: Retention ... OK' )
!------------------------------------------------------------------------------------------
end subroutine HMC_Phys_Retention_Cpl
!------------------------------------------------------------------------------------------
end module HMC_Module_Phys_Retention
!------------------------------------------------------------------------------------------