Line data Source code
1 : !> \file mo_mpr_runoff.f90
2 : !> \brief \copybrief mo_mpr_runoff
3 : !> \details \copydetails mo_mpr_runoff
4 :
5 : !> \brief multiscale parameter regionalization for runoff generation
6 : !> \details This contains the routine for multiscale parameter regionalization of the runoff parametrization.
7 : !> \authors Stephan Thober, Rohini Kumar
8 : !> \date Dec 2012
9 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
10 : !! mHM is released under the LGPLv3+ license \license_note
11 : !> \ingroup f_mpr
12 : module mo_mpr_runoff
13 :
14 : use mo_kind, only : i4, dp
15 :
16 : implicit none
17 :
18 : private
19 :
20 : public :: mpr_runoff
21 :
22 : ! ----------------------------------------------------------------------------
23 :
24 : contains
25 :
26 : ! ----------------------------------------------------------------------------
27 :
28 : ! NAME
29 : ! mpr_runoff
30 :
31 : ! PURPOSE
32 : !> \brief multiscale parameter regionalization for runoff parameters
33 :
34 : !> \details Perform the multiscale parameter regionalization for runoff
35 : !> global parameters (see mhm_parameter.nml). These are the following five
36 : !> parameters:
37 : !> - param(1) = interflowStorageCapacityFactor
38 : !> - param(2) = interflowRecession_slope
39 : !> - param(3) = fastInterflowRecession_forest
40 : !> - param(4) = slowInterflowRecession_Ks
41 : !> - param(5) = exponentSlowInterflow
42 :
43 : ! INTENT(IN)
44 : !> \param[in] "integer(i4), dimension(:) :: LCOVER0" land cover at level 0
45 : !> \param[in] "logical, dimension(:, :) :: mask0" mask at Level 0
46 : !> \param[in] "real(dp), dimension(:) :: SMs_FC0" [-] soil mositure deficit from field
47 : !> \param[in] "real(dp), dimension(:) :: slope_emp0" empirical quantile values F(slope)
48 : !> \param[in] "real(dp), dimension(:) :: KsVar_H0" [-] relative variability of saturated
49 : !> \param[in] "real(dp), dimension(5) :: param" global parameters
50 : !> \param[in] "integer(i4), dimension(:) :: cell_id0" Cell ids of hi res field
51 : !> \param[in] "integer(i4), dimension(:) :: upp_row_L1" Upper row of hi res block
52 : !> \param[in] "integer(i4), dimension(:) :: low_row_L1" Lower row of hi res block
53 : !> \param[in] "integer(i4), dimension(:) :: lef_col_L1" Left column of hi res block
54 : !> \param[in] "integer(i4), dimension(:) :: rig_col_L1" Right column of hi res block
55 : !> \param[in] "integer(i4), dimension(:) :: nL0_in_L1" Number of L0 cells within a L1 cell
56 : !> \param[in] "real(dp) :: c2TSTu" unit transformations
57 :
58 : ! INTENT(OUT)
59 : !> \param[out] "real(dp), dimension(:) :: L1_HL1" [10^-3 m] Threshhold water depth
60 : !> \param[out] "real(dp), dimension(:) :: L1_K0" [10^-3 m] Recession coefficient
61 : !> \param[out] "real(dp), dimension(:) :: L1_K1" [10^-3 m] Recession coefficient
62 : !> \param[out] "real(dp), dimension(:) :: L1_alpha" [1] Exponent for the upper reservoir
63 :
64 : ! HISTORY
65 : !> \authors Stephan Thober, Rohini Kumar
66 :
67 : !> \date Dec 2012
68 :
69 : ! Modifications:
70 : ! Stephan Thober Jan 2013 - updated calling sequence for upscaling operators
71 : ! Stephan Thober Dec 2013 - made header conform with mo_template
72 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
73 :
74 1242 : subroutine mpr_runoff(LCOVER0, mask0, SMs_FC0, slope_emp0, KsVar_H0, param, cell_id0, upp_row_L1, low_row_L1, &
75 1242 : lef_col_L1, rig_col_L1, nL0_in_L1, L1_HL1, L1_K0, L1_K1, L1_alpha)
76 :
77 : use mo_common_constants, only : nodata_dp, nodata_i4
78 : use mo_upscaling_operators, only : upscale_arithmetic_mean
79 :
80 : implicit none
81 :
82 : ! global parameters
83 : real(dp), dimension(5), intent(in) :: param
84 :
85 : ! [-] soil mositure deficit from field
86 : real(dp), dimension(:), intent(in) :: SMs_FC0
87 :
88 : ! empirical quantile values F(slope)
89 : real(dp), dimension(:), intent(in) :: slope_emp0
90 :
91 : ! [-] relative variability of saturated
92 : real(dp), dimension(:), intent(in) :: KsVar_H0
93 :
94 : ! land cover at level 0
95 : integer(i4), dimension(:), intent(in) :: LCOVER0
96 :
97 : ! mask at Level 0
98 : logical, dimension(:, :), intent(in) :: mask0
99 :
100 : ! Cell ids of hi res field
101 : integer(i4), dimension(:), intent(in) :: cell_id0
102 :
103 : ! Upper row of hi res block
104 : integer(i4), dimension(:), intent(in) :: upp_row_L1
105 :
106 : ! Lower row of hi res block
107 : integer(i4), dimension(:), intent(in) :: low_row_L1
108 :
109 : ! Left column of hi res block
110 : integer(i4), dimension(:), intent(in) :: lef_col_L1
111 :
112 : ! Right column of hi res block
113 : integer(i4), dimension(:), intent(in) :: rig_col_L1
114 :
115 : ! Number of L0 cells within a L1 cell
116 : integer(i4), dimension(:), intent(in) :: nL0_in_L1
117 :
118 : ! [10^-3 m] Threshhold water depth
119 : real(dp), dimension(:), intent(out) :: L1_HL1
120 :
121 : ! [10^-3 m] Recession coefficient
122 : real(dp), dimension(:), intent(out) :: L1_K0
123 :
124 : ! [10^-3 m] Recession coefficient
125 : real(dp), dimension(:), intent(out) :: L1_K1
126 :
127 : ! [1] Exponent for the upper reservoir
128 : real(dp), dimension(:), intent(out) :: L1_alpha
129 :
130 : ! temporal variable
131 19168378 : real(dp), dimension(size(SMs_FC0, 1)) :: tmp
132 :
133 :
134 : !-----------------------------
135 : ! FAST INTERFLOW
136 : !-----------------------------
137 : ! HL1 = f(soil properties; No reference found)
138 : ! Based on the saturation deficit from the field capacity status
139 : ! seems more reasonable and intutative.
140 : ! NOTE: This value for the sandy soils will have higher value of HL1, as compared to
141 : ! to clayey soil and so these soils can hold larger amount of amount.
142 19167964 : tmp = merge(param(1) * SMs_FC0, nodata_dp, cell_id0 .ne. nodata_i4)
143 : L1_HL1 = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, &
144 414 : Lef_col_L1, Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
145 :
146 : ! 1/K0 = f(terrain slope) [Booij, et. al.(2005), JoH]
147 : ! Steeper slopes resists (1/K0) fast water flows lesser as
148 : ! compared to that on the flater slope areas.
149 : ! Assuming that above relationship holds for all kind of land cover classes
150 :
151 : ! In the forested area surface resistance to fast interflow is higher as compared
152 : ! to the permeable land surface
153 :
154 19167964 : tmp = merge(param(2) * (2.0_dp - slope_emp0), nodata_dp, cell_id0 .ne. nodata_i4)
155 19167964 : tmp = merge(tmp * param(3), tmp, LCOVER0 .eq. 1)
156 : L1_K0 = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, &
157 414 : Lef_col_L1, Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
158 :
159 : ! To avoid numerical error in fully impervious areas (K0 == 0)
160 : ! minimum value of K0 is 1-day
161 15630 : L1_K0 = merge(1.0_dp, L1_K0, L1_K0 .lt. 1.0_dp)
162 :
163 : ! ------------------------------------------------------------------
164 : ! SLOW INTERFLOW
165 : ! ------------------------------------------------------------------
166 : ! K1 = f(terrian slope, Booij, et. al.(2005), JoH)
167 : ! = f(soil properties, LC & New modification)
168 : ! K1 = K0 + K1(soil-Ks)
169 :
170 0 : tmp = merge(param(2) * (2.0_dp - slope_emp0) + param(4) * (1.0_dp + KsVar_H0), &
171 19167964 : nodata_dp, cell_id0 .ne. nodata_i4)
172 : L1_K1 = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, &
173 414 : Lef_col_L1, Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
174 :
175 : ! minimum value of K1 is 1-day
176 15630 : L1_K1 = merge(2.0_dp, L1_K1, L1_K1 .lt. 2.0_dp)
177 :
178 :
179 : ! alpha = f(soil type; variabitity of Ks)
180 : ! Lower the alpha (exponent of slow interflow) means lower amount of
181 : ! water released from the storage to contribute for slow interflow.
182 : ! For instance sandy soils will have lower value of alpha as comapred to
183 : ! the clayey soils.
184 : ! This assumption is quite realistic in physical sense...
185 828 : tmp = merge(param(5) * (1.0_dp / KsVar_H0) * (1.0_dp / (1.0_dp + SMs_FC0)), &
186 19168378 : nodata_dp, cell_id0 .ne. nodata_i4)
187 : L1_alpha = upscale_arithmetic_mean(nL0_in_L1, Upp_row_L1, Low_row_L1, Lef_col_L1, &
188 414 : Rig_col_L1, cell_id0, mask0, nodata_dp, tmp)
189 :
190 : ! constraints and unit transformation
191 16044 : L1_K0 = merge(L1_K1, L1_K0, L1_K0 .gt. L1_K1)
192 :
193 :
194 414 : end subroutine mpr_runoff
195 :
196 : end module mo_mpr_runoff
|