Line data Source code
1 : !> \file mo_mhm_interface_run.f90
2 : !> \brief \copybrief mo_mhm_interface_run
3 : !> \details \copydetails mo_mhm_interface_run
4 :
5 : !> \brief Module providing interfaces for running preconfigured mHM.
6 : !> \details Interfaces to control the mHM run from outside (prepare domain, do timestep, ...).
7 : !> \authors Sebastian Mueller, Matthias Kelbling
8 : !> \version 0.1
9 : !> \date Jan 2022
10 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
11 : !! mHM is released under the LGPLv3+ license \license_note
12 : !> \ingroup f_mhm
13 : module mo_mhm_interface_run
14 :
15 : ! forces
16 : use mo_kind, only: i4, dp
17 : use mo_message, only: message, error_message
18 : use mo_string_utils, only : num2str
19 : ! mhm
20 : use mo_common_run_variables, only : run_cfg
21 : use mo_optimization_types, only : optidata_sim
22 : use mo_common_datetime_type, only : datetimeinfo
23 : use mo_common_mHM_mRM_variables, only : &
24 : resolutionRouting, &
25 : LCyearId, &
26 : mhmFileRestartIn, &
27 : mrmFileRestartIn, &
28 : nTstepDay, &
29 : optimize, &
30 : read_restart, &
31 : simPer, &
32 : timeStep, &
33 : warmingDays, &
34 : c2TSTu, &
35 : restart_reset_fluxes_states
36 : use mo_common_variables, only : &
37 : global_parameters, &
38 : level0, &
39 : level1, &
40 : domainMeta, &
41 : processMatrix
42 : use mo_global_variables, only : &
43 : meteo_handler, &
44 : L1_Throughfall, &
45 : L1_aETCanopy, &
46 : L1_aETSealed, &
47 : L1_aETSoil, &
48 : L1_baseflow, &
49 : L1_fastRunoff, &
50 : L1_infilSoil, &
51 : L1_inter, &
52 : L1_melt, &
53 : L1_neutrons, &
54 : L1_percol, &
55 : L1_pet_calc, &
56 : L1_temp_calc, &
57 : L1_prec_calc, &
58 : L1_preEffect, &
59 : L1_rain, &
60 : L1_runoffSeal, &
61 : L1_satSTW, &
62 : L1_sealSTW, &
63 : L1_slowRunoff, &
64 : L1_snow, &
65 : L1_snowPack, &
66 : L1_soilMoist, &
67 : L1_total_runoff, &
68 : L1_unsatSTW, &
69 : L1_twsaObs, &
70 : L1_etObs, &
71 : L1_smObs, &
72 : L1_neutronsObs, &
73 : evap_coeff, &
74 : nSoilHorizons_sm_input, &
75 : neutron_integral_AFast, &
76 : outputFlxState, &
77 : timeStep_model_outputs, &
78 : BFI_qBF_sum, &
79 : BFI_qT_sum
80 : use mo_init_states, only : variables_default_init, fluxes_states_default_init
81 : use mo_julian, only : caldat, julday
82 : use mo_string_utils, only : num2str
83 : use mo_mhm, only : mhm
84 : use mo_restart, only : read_restart_states
85 : use mo_write_fluxes_states, only : mHM_updateDataset, mHM_OutputDataset
86 : use mo_mrm_write_fluxes_states, only : mRM_updateDataset, mRM_OutputDataset, GW_OutputDataset, GW_updateDataset
87 : use mo_constants, only : HourSecs
88 : use mo_common_variables, only : resolutionHydrology
89 : use mo_mrm_global_variables, only : &
90 : InflowGauge, &
91 : L11_C1, &
92 : L11_C2, &
93 : L11_L1_Id, &
94 : L11_TSrout, &
95 : L11_fromN, &
96 : L11_length, &
97 : L11_nLinkFracFPimp, &
98 : L11_nOutlets, &
99 : L11_netPerm, &
100 : L11_qMod, &
101 : L11_qOUT, &
102 : L11_qTIN, &
103 : L11_qTR, &
104 : L11_slope, &
105 : L11_toN, &
106 : L1_L11_Id, &
107 : domain_mrm, &
108 : level11, &
109 : mRM_runoff, &
110 : outputFlxState_mrm, &
111 : timeStep_model_outputs_mrm, &
112 : gw_coupling, &
113 : L0_river_head_mon_sum, &
114 : riv_temp_pcs
115 : use mo_mpr_global_variables, only : &
116 : L1_HarSamCoeff, &
117 : L1_PrieTayAlpha, &
118 : L1_aeroResist, &
119 : L1_alpha, &
120 : L1_degDay, &
121 : L1_degDayInc, &
122 : L1_degDayMax, &
123 : L1_degDayNoPre, &
124 : L1_fAsp, &
125 : L1_fRoots, &
126 : L1_fSealed, &
127 : L1_jarvis_thresh_c1, &
128 : L1_kBaseFlow, &
129 : L1_kPerco, &
130 : L1_kSlowFlow, &
131 : L1_karstLoss, &
132 : L1_kfastFlow, &
133 : L1_maxInter, &
134 : L1_petLAIcorFactor, &
135 : L1_sealedThresh, &
136 : L1_soilMoistExp, &
137 : L1_soilMoistFC, &
138 : L1_soilMoistSat, &
139 : L1_surfResist, &
140 : L1_tempThresh, &
141 : L1_unsatThresh, &
142 : L1_wiltingPoint, &
143 : L1_No_Count, &
144 : L1_bulkDens, &
145 : L1_latticeWater, &
146 : L1_COSMICL3, &
147 : HorizonDepth_mHM, &
148 : nSoilHorizons_mHM
149 : use mo_mrm_init, only : variables_default_init_routing, fluxes_states_default_init_routing
150 : use mo_mrm_mpr, only : mrm_update_param
151 : use mo_mrm_restart, only : mrm_read_restart_states
152 : use mo_mrm_routing, only : mrm_routing
153 : use mo_utils, only : ge
154 : use mo_mrm_river_head, only: calc_river_head
155 : use mo_mpr_eval, only : mpr_eval
156 :
157 : contains
158 :
159 : !> \brief prepare single run of mHM
160 122 : subroutine mhm_interface_run_prepare(parameterset, opti_domain_indices, runoff_present, BFI_present)
161 : implicit none
162 : !> a set of global parameter (gamma) to run mHM, DIMENSION [no. of global_Parameters]
163 : real(dp), dimension(:), optional, intent(in) :: parameterset
164 : !> selected domains for optimization
165 : integer(i4), dimension(:), optional, intent(in) :: opti_domain_indices
166 : !> whether runoff is present
167 : logical, optional, intent(in) :: runoff_present
168 : !> whether BFI is present
169 : logical, optional, intent(in) :: BFI_present
170 :
171 : integer(i4) :: i, iDomain, domainID
172 :
173 : ! run_cfg%output_runoff = .false. by default
174 61 : if (present(runoff_present)) run_cfg%output_runoff = runoff_present
175 : ! run_cfg%output_BFI = .false. by default
176 61 : if (present(BFI_present)) run_cfg%output_BFI = BFI_present
177 :
178 : ! store current parameter set
179 183 : allocate(run_cfg%parameterset(size(global_parameters, dim=1)))
180 61 : if (.not. present(parameterset) .and. optimize) then
181 0 : call error_message("mhm_interface_run_prepare: Can't optimize without parameter!")
182 61 : else if (.not. present(parameterset)) then
183 0 : run_cfg%parameterset = global_parameters(:, 3)
184 : else
185 3283 : run_cfg%parameterset = parameterset
186 : end if
187 :
188 : ! store domain indices for domain loop
189 61 : if (optimize .and. present(opti_domain_indices)) then
190 28 : run_cfg%nDomains = size(opti_domain_indices)
191 84 : allocate(run_cfg%domain_indices(run_cfg%nDomains))
192 98 : run_cfg%domain_indices = opti_domain_indices
193 : else
194 33 : run_cfg%nDomains = domainMeta%nDomains
195 99 : allocate(run_cfg%domain_indices(run_cfg%nDomains))
196 186 : run_cfg%domain_indices = [(i, i=1, run_cfg%nDomains)]
197 : end if
198 :
199 : !----------------------------------------------------------
200 : ! Check optionals and initialize
201 : !----------------------------------------------------------
202 61 : if (run_cfg%output_runoff) then
203 50 : do i = 1, run_cfg%nDomains
204 25 : iDomain = run_cfg%get_domain_index(i)
205 25 : domainID = domainMeta%indices(iDomain)
206 50 : if (.not. domainMeta%doRouting(iDomain)) then
207 : call error_message("***ERROR: runoff for domain", trim(num2str(domainID)),&
208 0 : "can not be produced, since routing process is off in Process Matrix")
209 : end if
210 : end do
211 : end if
212 :
213 : ! prepare BFI calculation
214 61 : if (run_cfg%output_BFI) then
215 0 : allocate(BFI_qBF_sum(run_cfg%nDomains))
216 0 : allocate(BFI_qT_sum(run_cfg%nDomains))
217 0 : BFI_qBF_sum = 0.0_dp
218 0 : BFI_qT_sum = 0.0_dp
219 : end if
220 :
221 61 : if (read_restart) then
222 2 : do i = 1, run_cfg%nDomains
223 1 : iDomain = run_cfg%get_domain_index(i)
224 1 : domainID = domainMeta%indices(iDomain)
225 : ! this reads the eff. parameters and optionally the states and fluxes
226 2 : call read_restart_states(iDomain, domainID, mhmFileRestartIn(iDomain))
227 : end do
228 1 : if (restart_reset_fluxes_states) then
229 0 : call message(' Resetting mHM states and fluxes from restart files ...')
230 0 : call fluxes_states_default_init()
231 : end if
232 : else
233 60 : call variables_default_init()
234 60 : call mpr_eval(run_cfg%parameterset)
235 60 : if (processMatrix(8, 1) > 0) call variables_default_init_routing()
236 : end if
237 :
238 305 : allocate(run_cfg%L1_fNotSealed(size(L1_fSealed, 1), size(L1_fSealed, 2), size(L1_fSealed, 3)))
239 15650 : run_cfg%L1_fNotSealed = 1.0_dp - L1_fSealed
240 :
241 61 : end subroutine mhm_interface_run_prepare
242 :
243 : !> \brief get number of domains for looping
244 61 : subroutine mhm_interface_run_get_ndomains(ndomains)
245 : implicit none
246 : integer(i4), intent(inout) :: ndomains !< number of domains
247 61 : ndomains = run_cfg%nDomains
248 61 : end subroutine mhm_interface_run_get_ndomains
249 :
250 : !> \brief prepare single domain to run mHM on
251 164 : subroutine mhm_interface_run_prepare_domain(domain, etOptiSim, twsOptiSim, neutronsOptiSim, smOptiSim)
252 : implicit none
253 : !> domain loop counter
254 : integer(i4), intent(in), optional :: domain
255 : !> returns soil moisture time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
256 : type(optidata_sim), dimension(:), optional, intent(inout) :: smOptiSim
257 : !> dim1=ncells, dim2=time
258 : type(optidata_sim), dimension(:), optional, intent(inout) :: neutronsOptiSim
259 : !> returns evapotranspiration time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
260 : type(optidata_sim), dimension(:), optional, intent(inout) :: etOptiSim
261 : !> returns tws time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
262 : type(optidata_sim), dimension(:), optional, intent(inout) :: twsOptiSim
263 :
264 : integer(i4) :: iDomain, domainID
265 :
266 : ! get domain index
267 82 : run_cfg%selected_domain = 1_i4
268 82 : if (present(domain)) run_cfg%selected_domain = domain
269 82 : iDomain = run_cfg%get_domain_index(run_cfg%selected_domain)
270 82 : domainID = domainMeta%indices(iDomain)
271 :
272 : ! evapotranspiration optimization
273 82 : if (present(etOptiSim)) call etOptiSim(iDomain)%init(L1_etObs(iDomain))
274 : ! total water storage optimization
275 82 : if (present(twsOptiSim)) call twsOptiSim(iDomain)%init(L1_twsaObs(iDomain))
276 : ! neutrons optimization
277 82 : if (present(neutronsOptiSim)) call neutronsOptiSim(iDomain)%init(L1_neutronsObs(iDomain))
278 : ! sm optimization
279 82 : if (present(smOptiSim)) call smOptiSim(iDomain)%init(L1_smObs(iDomain))
280 :
281 : ! get Domain information
282 82 : run_cfg%nCells = level1(iDomain)%nCells
283 82 : run_cfg%mask1 => level1(iDomain)%mask
284 82 : run_cfg%s1 = level1(iDomain)%iStart
285 82 : run_cfg%e1 = level1(iDomain)%iEnd
286 :
287 82 : if (domainMeta%doRouting(iDomain)) then
288 : ! ----------------------------------------
289 : ! initialize factor between routing resolution and hydrologic model resolution
290 : ! ----------------------------------------
291 46 : run_cfg%tsRoutFactor = 1_i4
292 138 : allocate(run_cfg%InflowDischarge(size(InflowGauge%Q, dim = 2)))
293 95 : run_cfg%InflowDischarge = 0._dp
294 :
295 : ! read states from restart
296 46 : if (read_restart) then
297 1 : call mrm_read_restart_states(iDomain, domainID, mrmFileRestartIn(iDomain))
298 1 : if (restart_reset_fluxes_states) then
299 0 : call message(' Resetting mRM states and fluxes from restart files for domain ', num2str(iDomain), ' ...')
300 0 : call fluxes_states_default_init_routing(iDomain)
301 : end if
302 : end if
303 : ! get Domain information at L11 if routing is activated
304 46 : run_cfg%s11 = level11(iDomain)%iStart
305 46 : run_cfg%e11 = level11(iDomain)%iEnd
306 46 : run_cfg%mask11 => level11(iDomain)%mask
307 :
308 : ! initialize routing parameters (has to be called for routing option 2)
309 46 : if ((processMatrix(8, 1) .eq. 2) .or. (processMatrix(8, 1) .eq. 3)) &
310 : call mrm_update_param( &
311 10 : iDomain, run_cfg%parameterset(processMatrix(8, 3) - processMatrix(8, 2) + 1 : processMatrix(8, 3)))
312 : ! initialize variable for runoff for routing
313 138 : allocate(run_cfg%RunToRout(run_cfg%e1 - run_cfg%s1 + 1))
314 2180 : run_cfg%RunToRout = 0._dp
315 :
316 46 : if ( riv_temp_pcs%active ) then
317 : ! set indices for current L11 domain
318 1 : riv_temp_pcs%s11 = run_cfg%s11
319 1 : riv_temp_pcs%e11 = run_cfg%e11
320 : ! allocate current L1 lateral components
321 1 : call riv_temp_pcs%alloc_lateral(run_cfg%nCells)
322 : end if
323 : end if
324 :
325 : ! init datetime variable
326 82 : call run_cfg%domainDateTime%init(iDomain)
327 :
328 : ! init time-loop
329 82 : run_cfg%time_step = 0_i4
330 :
331 143 : end subroutine mhm_interface_run_prepare_domain
332 :
333 : !> \brief check if current time loop is finished
334 2020272 : subroutine mhm_interface_run_finished(time_loop_finished)
335 : implicit none
336 : logical, intent(inout) :: time_loop_finished !< flag to indicate end of timeloop
337 2020272 : time_loop_finished = run_cfg%time_step == run_cfg%domainDateTime%nTimeSteps
338 82 : end subroutine mhm_interface_run_finished
339 :
340 : !> \brief do one time-step on current domain
341 2020272 : subroutine mhm_interface_run_do_time_step()
342 : implicit none
343 :
344 : integer(i4) :: iDomain, domainID, tt, jj, s1, e1
345 :
346 : ! increment time step count (first input is 0)
347 2020272 : run_cfg%time_step = run_cfg%time_step + 1_i4
348 : ! current time counter and domain indices
349 2020272 : tt = run_cfg%time_step
350 2020272 : s1 = run_cfg%s1
351 2020272 : e1 = run_cfg%e1
352 2020272 : iDomain = run_cfg%get_domain_index(run_cfg%selected_domain)
353 2020272 : domainID = domainMeta%indices(iDomain)
354 :
355 2020272 : call run_cfg%domainDateTime%update_LAI_timestep()
356 :
357 : ! update the meteo-handler
358 : call meteo_handler%update_timestep( &
359 : tt=tt, &
360 : time=run_cfg%domainDateTime%newTime - 0.5_dp, &
361 : iDomain=iDomain, &
362 : level1=level1, &
363 2020272 : simPer=simPer)
364 :
365 : ! get corrected pet
366 0 : call meteo_handler%get_corrected_pet(pet_calc=L1_pet_calc(s1 : e1), &
367 : ! pet calculation dependencies
368 0 : petLAIcorFactorL1 = L1_petLAIcorFactor(s1 : e1, run_cfg%domainDateTime%iLAI, run_cfg%domainDateTime%yId), &
369 0 : fAsp = L1_fAsp(s1 : e1, 1, 1), &
370 0 : HarSamCoeff = L1_HarSamCoeff(s1 : e1, 1, 1), &
371 0 : latitude = pack(level1(iDomain)%y, level1(iDomain)%mask), &
372 0 : PrieTayAlpha = L1_PrieTayAlpha(s1 : e1, run_cfg%domainDateTime%iLAI, 1), &
373 0 : aeroResist = L1_aeroResist(s1 : e1, run_cfg%domainDateTime%iLAI, run_cfg%domainDateTime%yId), &
374 2020272 : surfResist = L1_surfResist(s1 : e1, run_cfg%domainDateTime%iLAI, 1))
375 :
376 : ! get temperature and precipitation
377 2020272 : call meteo_handler%get_temp(temp_calc=L1_temp_calc(s1 : e1))
378 2020272 : call meteo_handler%get_prec(prec_calc=L1_prec_calc(s1 : e1))
379 :
380 : ! -------------------------------------------------------------------------
381 : ! ARGUMENT LIST KEY FOR mHM
382 : ! -------------------------------------------------------------------------
383 : ! C CONFIGURATION
384 : ! F FORCING DATA L2
385 : ! Q INFLOW FROM UPSTREAM AREAS
386 : ! L0 MORPHOLOGIC DATA L0
387 : ! L1 MORPHOLOGIC DATA L1
388 : ! L11 MORPHOLOGIC DATA L11
389 : ! P GLOBAL PARAMETERS
390 : ! E EFFECTIVE PARAMETER FIELDS (L1, L11 levels)
391 : ! S STATE VARIABLES L1
392 : ! X FLUXES (L1, L11 levels)
393 : ! --------------------------------------------------------------------------
394 : call mhm( &
395 : read_states = read_restart, & ! IN C
396 : tt = tt, &
397 : time = run_cfg%domainDateTime%newTime - 0.5_dp, &
398 : processMatrix = processMatrix, &
399 : horizon_depth = HorizonDepth_mHM, & ! IN C
400 : nCells1 = run_cfg%nCells, &
401 : nHorizons_mHM = nSoilHorizons_mHM, &
402 : c2TSTu = c2TSTu, & ! IN C
403 : neutron_integral_AFast = neutron_integral_AFast, & ! IN C
404 : evap_coeff = evap_coeff, &
405 0 : fSealed1 = L1_fSealed(s1 : e1, 1, run_cfg%domainDateTime%yId), & ! INOUT L1
406 0 : interc = L1_inter(s1 : e1), &
407 0 : snowpack = L1_snowPack(s1 : e1), &
408 0 : sealedStorage = L1_sealSTW(s1 : e1), & ! INOUT S
409 0 : soilMoisture = L1_soilMoist(s1 : e1, :), &
410 0 : unsatStorage = L1_unsatSTW(s1 : e1), &
411 0 : satStorage = L1_satSTW(s1 : e1), & ! INOUT S
412 0 : neutrons = L1_neutrons(s1 : e1), & ! INOUT S
413 0 : pet_calc = L1_pet_calc(s1 : e1), & ! INOUT X
414 0 : temp_calc = L1_temp_calc(s1 : e1), & ! INOUT X
415 0 : prec_calc = L1_prec_calc(s1 : e1), & ! INOUT X
416 0 : aet_soil = L1_aETSoil(s1 : e1, :), &
417 0 : aet_canopy = L1_aETCanopy(s1 : e1), &
418 0 : aet_sealed = L1_aETSealed(s1 : e1), & ! INOUT X
419 0 : baseflow = L1_baseflow(s1 : e1), &
420 0 : infiltration = L1_infilSoil(s1 : e1, :), &
421 0 : fast_interflow = L1_fastRunoff(s1 : e1), & ! INOUT X
422 0 : melt = L1_melt(s1 : e1), &
423 0 : perc = L1_percol(s1 : e1), &
424 0 : prec_effect = L1_preEffect(s1 : e1), &
425 0 : rain = L1_rain(s1 : e1), & ! INOUT X
426 0 : runoff_sealed = L1_runoffSeal(s1 : e1), &
427 0 : slow_interflow = L1_slowRunoff(s1 : e1), &
428 0 : snow = L1_snow(s1 : e1), & ! INOUT X
429 0 : throughfall = L1_Throughfall(s1 : e1), &
430 0 : total_runoff = L1_total_runoff(s1 : e1), & ! INOUT X
431 : ! MPR
432 0 : alpha = L1_alpha(s1 : e1, 1, run_cfg%domainDateTime%yId), &
433 0 : deg_day_incr = L1_degDayInc(s1 : e1, 1, run_cfg%domainDateTime%yId), &
434 0 : deg_day_max = L1_degDayMax(s1 : e1, 1, run_cfg%domainDateTime%yId), & ! INOUT E1
435 0 : deg_day_noprec = L1_degDayNoPre(s1 : e1, 1, run_cfg%domainDateTime%yId), &
436 0 : deg_day = L1_degDay(s1 : e1, 1, 1), & ! INOUT E1
437 0 : frac_roots = L1_fRoots(s1 : e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
438 0 : interc_max = L1_maxInter(s1 : e1, run_cfg%domainDateTime%iLAI, 1), &
439 0 : karst_loss = L1_karstLoss(s1 : e1, 1, 1), & ! INOUT E1
440 0 : k0 = L1_kFastFlow(s1 : e1, 1, run_cfg%domainDateTime%yId), &
441 0 : k1 = L1_kSlowFlow(s1 : e1, 1, run_cfg%domainDateTime%yId), & ! INOUT E1
442 0 : k2 = L1_kBaseFlow(s1 : e1, 1, run_cfg%domainDateTime%yId), &
443 0 : kp = L1_kPerco(s1 : e1, 1, run_cfg%domainDateTime%yId), & ! INOUT E1
444 0 : soil_moist_FC = L1_soilMoistFC(s1 : e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
445 0 : soil_moist_sat = L1_soilMoistSat(s1 : e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
446 0 : soil_moist_exponen = L1_soilMoistExp(s1 : e1, :, run_cfg%domainDateTime%yId), &
447 0 : jarvis_thresh_c1 = L1_jarvis_thresh_c1(s1 : e1, 1, 1), & ! INOUT E1
448 0 : temp_thresh = L1_tempThresh(s1 : e1, 1, run_cfg%domainDateTime%yId), &
449 0 : unsat_thresh = L1_unsatThresh(s1 : e1, 1, 1), & ! INOUT E1
450 0 : water_thresh_sealed = L1_sealedThresh(s1 : e1, 1, 1), & ! INOUT E1
451 0 : wilting_point = L1_wiltingPoint(s1 : e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
452 : ! >> neutron count
453 0 : No_count = L1_No_Count(s1:e1, 1, 1), & ! INOUT E1
454 0 : bulkDens = L1_bulkDens(s1:e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
455 0 : latticeWater = L1_latticeWater(s1:e1, :, run_cfg%domainDateTime%yId), & ! INOUT E1
456 0 : COSMICL3 = L1_COSMICL3(s1:e1, :, run_cfg%domainDateTime%yId) & ! INOUT E1
457 2020272 : )
458 :
459 : ! call mRM routing
460 2020272 : run_cfg%doRoute = .false.
461 2020272 : if (domainMeta%doRouting(iDomain)) then
462 : ! set discharge timestep
463 779928 : run_cfg%iDischargeTS = ceiling(real(tt, dp) / real(nTstepDay, dp))
464 : ! set input variables for routing
465 779928 : if (processMatrix(8, 1) .eq. 1) then
466 : ! >>>
467 : ! >>> original Muskingum routing, executed every time
468 : ! >>>
469 495096 : run_cfg%doRoute = .True.
470 495096 : run_cfg%tsRoutFactorIn = 1._dp
471 495096 : run_cfg%timestep_rout = timestep
472 25781256 : run_cfg%RunToRout = L1_total_runoff(s1 : e1) ! runoff [mm TS-1] mm per timestep
473 1524600 : run_cfg%InflowDischarge = InflowGauge%Q(run_cfg%iDischargeTS, :) ! inflow discharge in [m3 s-1]
474 : !
475 284832 : else if ((processMatrix(8, 1) .eq. 2) .or. &
476 : (processMatrix(8, 1) .eq. 3)) then
477 : ! >>>
478 : ! >>> adaptive timestep
479 : ! >>>
480 : run_cfg%doRoute = .False.
481 : ! calculate factor
482 284832 : run_cfg%tsRoutFactor = L11_tsRout(iDomain) / (timestep * HourSecs)
483 : ! print *, 'routing factor: ', tsRoutFactor
484 : ! prepare routing call
485 284832 : if (run_cfg%tsRoutFactor .lt. 1._dp) then
486 : ! ----------------------------------------------------------------
487 : ! routing timesteps are shorter than hydrologic time steps
488 : ! ----------------------------------------------------------------
489 : ! set all input variables
490 0 : run_cfg%tsRoutFactorIn = run_cfg%tsRoutFactor
491 0 : run_cfg%RunToRout = L1_total_runoff(s1 : e1) ! runoff [mm TS-1] mm per timestep
492 0 : run_cfg%InflowDischarge = InflowGauge%Q(run_cfg%iDischargeTS, :) ! inflow discharge in [m3 s-1]
493 0 : run_cfg%timestep_rout = timestep
494 0 : run_cfg%doRoute = .True.
495 : else
496 : ! ----------------------------------------------------------------
497 : ! routing timesteps are longer than hydrologic time steps
498 : ! ----------------------------------------------------------------
499 : ! set all input variables
500 284832 : run_cfg%tsRoutFactorIn = run_cfg%tsRoutFactor
501 : ! Runoff is accumulated in [mm]
502 10253952 : run_cfg%RunToRout = run_cfg%RunToRout + L1_total_runoff(s1 : e1)
503 854496 : run_cfg%InflowDischarge = run_cfg%InflowDischarge + InflowGauge%Q(run_cfg%iDischargeTS, :)
504 : ! reset tsRoutFactorIn if last period did not cover full period
505 284832 : if ((tt == run_cfg%domainDateTime%nTimeSteps) .and. (mod(tt, nint(run_cfg%tsRoutFactorIn)) /= 0_i4)) &
506 0 : run_cfg%tsRoutFactorIn = mod(tt, nint(run_cfg%tsRoutFactorIn))
507 284832 : if ((mod(tt, nint(run_cfg%tsRoutFactorIn)) .eq. 0_i4) .or. (tt .eq. run_cfg%domainDateTime%nTimeSteps)) then
508 : ! Inflow discharge is given as flow-rate and has to be converted to [m3]
509 278280 : run_cfg%InflowDischarge = run_cfg%InflowDischarge / run_cfg%tsRoutFactorIn
510 139140 : run_cfg%timestep_rout = timestep * nint(run_cfg%tsRoutFactorIn, i4)
511 139140 : run_cfg%doRoute = .True.
512 : end if
513 : end if
514 : end if
515 : ! prepare temperature routing
516 779928 : if ( riv_temp_pcs%active ) then
517 : ! init riv-temp from current air temp
518 13104 : if ( tt .eq. 1_i4 ) call riv_temp_pcs%init_riv_temp( &
519 0 : temp_air = L1_temp_calc(s1 : e1), &
520 0 : efecarea = level1(iDomain)%CellArea * 1.E-6_dp, &
521 0 : L1_L11_Id = L1_L11_Id(s1 : e1), &
522 0 : L11_areacell = level11(iDomain)%CellArea * 1.E-6_dp, &
523 0 : L11_L1_Id = L11_L1_Id(run_cfg%s11 : run_cfg%e11), &
524 0 : map_flag = ge(resolutionRouting(iDomain), resolutionHydrology(iDomain)) &
525 69 : )
526 : ! get riv-temp specific meteo arrays
527 13104 : call meteo_handler%get_ssrd(riv_temp_pcs%L1_ssrd_calc)
528 13104 : call meteo_handler%get_strd(riv_temp_pcs%L1_strd_calc)
529 13104 : call meteo_handler%get_tann(riv_temp_pcs%L1_tann_calc)
530 : ! accumulate source Energy at L1 level
531 : call riv_temp_pcs%acc_source_E( &
532 0 : fSealed_area_fraction = L1_fSealed(s1 : e1, 1, run_cfg%domainDateTime%yId), &
533 0 : fast_interflow = L1_fastRunoff(s1 : e1), &
534 0 : slow_interflow = L1_slowRunoff(s1 : e1), &
535 0 : baseflow = L1_baseflow(s1 : e1), &
536 0 : direct_runoff = L1_runoffSeal(s1 : e1), &
537 0 : temp_air = L1_temp_calc(s1 : e1) &
538 13104 : )
539 : ! if routing should be performed, scale source energy to L11 level
540 13104 : if ( run_cfg%doRoute ) call riv_temp_pcs%finalize_source_E( &
541 0 : efecarea = level1(iDomain)%CellArea * 1.E-6_dp, &
542 0 : L1_L11_Id = L1_L11_Id(s1 : e1), &
543 0 : L11_areacell = level11(iDomain)%CellArea * 1.E-6_dp, &
544 0 : L11_L1_Id = L11_L1_Id(run_cfg%s11 : run_cfg%e11), &
545 : timestep = run_cfg%timestep_rout, &
546 0 : map_flag = ge(resolutionRouting(iDomain), resolutionHydrology(iDomain)) &
547 452088 : )
548 : end if
549 : ! -------------------------------------------------------------------
550 : ! execute routing
551 : ! -------------------------------------------------------------------
552 779928 : if (run_cfg%doRoute) call mRM_routing(&
553 : ! general INPUT variables
554 : read_restart, &
555 : processMatrix(8, 1), & ! parse process Case to be used
556 0 : run_cfg%parameterset(processMatrix(8, 3) - processMatrix(8, 2) + 1 : processMatrix(8, 3)), & ! routing par.
557 : run_cfg%RunToRout, & ! runoff [mm TS-1] mm per timestep old: L1_total_runoff_in(run_cfg%s1:run_cfg%e1, tt), &
558 0 : level1(iDomain)%CellArea * 1.E-6_dp, &
559 0 : L1_L11_Id(s1 : e1), &
560 0 : level11(iDomain)%CellArea * 1.E-6_dp, &
561 0 : L11_L1_Id(run_cfg%s11 : run_cfg%e11), &
562 0 : L11_netPerm(run_cfg%s11 : run_cfg%e11), & ! routing order at L11
563 0 : L11_fromN(run_cfg%s11 : run_cfg%e11), & ! link source at L11
564 0 : L11_toN(run_cfg%s11 : run_cfg%e11), & ! link target at L11
565 0 : L11_nOutlets(iDomain), & ! number of outlets
566 : run_cfg%timestep_rout, & ! timestep of runoff to rout [h]
567 : run_cfg%tsRoutFactorIn, & ! simulate timestep in [h]
568 0 : level11(iDomain)%nCells, & ! number of Nodes
569 0 : domain_mrm(iDomain)%nInflowGauges, &
570 : domain_mrm(iDomain)%InflowGaugeIndexList(:), &
571 : domain_mrm(iDomain)%InflowGaugeHeadwater(:), &
572 : domain_mrm(iDomain)%InflowGaugeNodeList(:), &
573 : run_cfg%InflowDischarge, &
574 : domain_mrm(iDomain)%nGauges, &
575 : domain_mrm(iDomain)%gaugeIndexList(:), &
576 : domain_mrm(iDomain)%gaugeNodeList(:), &
577 0 : ge(resolutionRouting(iDomain), resolutionHydrology(iDomain)), &
578 : ! original routing specific input variables
579 0 : L11_length(run_cfg%s11 : run_cfg%e11 - 1), & ! link length
580 0 : L11_slope(run_cfg%s11 : run_cfg%e11 - 1), &
581 0 : L11_nLinkFracFPimp(run_cfg%s11 : run_cfg%e11, run_cfg%domainDateTime%yId), & ! fraction of impervious layer at L11 scale
582 : ! general INPUT/OUTPUT variables
583 0 : L11_C1(run_cfg%s11 : run_cfg%e11), & ! first muskingum parameter
584 0 : L11_C2(run_cfg%s11 : run_cfg%e11), & ! second muskigum parameter
585 0 : L11_qOUT(run_cfg%s11 : run_cfg%e11), & ! routed runoff flowing out of L11 cell
586 0 : L11_qTIN(run_cfg%s11 : run_cfg%e11, :), & ! inflow water into the reach at L11
587 0 : L11_qTR(run_cfg%s11 : run_cfg%e11, :), & !
588 0 : L11_qMod(run_cfg%s11 : run_cfg%e11), &
589 : mRM_runoff(tt, :) &
590 58038084 : )
591 : ! -------------------------------------------------------------------
592 : ! groundwater coupling
593 : ! -------------------------------------------------------------------
594 779928 : if (gw_coupling) then
595 : call calc_river_head(iDomain, L11_Qmod, L0_river_head_mon_sum)
596 : end if
597 : ! -------------------------------------------------------------------
598 : ! reset variables
599 : ! -------------------------------------------------------------------
600 779928 : if (processMatrix(8, 1) .eq. 1) then
601 : ! reset Input variables
602 1029504 : run_cfg%InflowDischarge = 0._dp
603 25286160 : run_cfg%RunToRout = 0._dp
604 284832 : else if ((processMatrix(8, 1) .eq. 2) .or. (processMatrix(8, 1) .eq. 3)) then
605 284832 : if ((.not. (run_cfg%tsRoutFactorIn .lt. 1._dp)) .and. run_cfg%doRoute) then
606 423972 : do jj = 1, nint(run_cfg%tsRoutFactorIn) ! BUG: this should start at 2
607 708804 : mRM_runoff(tt - jj + 1, :) = mRM_runoff(tt, :)
608 : end do
609 : ! reset Input variables
610 278280 : run_cfg%InflowDischarge = 0._dp
611 4869900 : run_cfg%RunToRout = 0._dp
612 : ! reset lateral fluxes and time-step counter if routing was done
613 139140 : if ( riv_temp_pcs%active ) call riv_temp_pcs%reset_timestep()
614 : end if
615 : ! if routing is done every time-step, reset river-temp time step
616 284832 : if (run_cfg%tsRoutFactor .lt. 1._dp .and. riv_temp_pcs%active ) call riv_temp_pcs%reset_timestep()
617 : end if
618 : end if
619 :
620 : ! output only for evaluation period
621 2020272 : run_cfg%domainDateTime%tIndex_out = (tt - warmingDays(iDomain) * nTstepDay) ! tt if write out of warming period
622 :
623 2020272 : call run_cfg%domainDateTime%increment()
624 :
625 : ! update the year-dependent domainDateTime%yId (land cover id)
626 2020272 : if (run_cfg%domainDateTime%is_new_year .and. tt < run_cfg%domainDateTime%nTimeSteps) then
627 165 : run_cfg%domainDateTime%yId = LCyearId(run_cfg%domainDateTime%year, iDomain)
628 : end if
629 :
630 : ! calculate BFI releated after warming days if wanted
631 2020272 : if ( run_cfg%output_BFI .and. (run_cfg%domainDateTime%tIndex_out > 0_i4) ) then
632 0 : BFI_qBF_sum(iDomain) = BFI_qBF_sum(iDomain) &
633 0 : + sum(L1_baseflow(s1 : e1) * level1(iDomain)%CellArea) / level1(iDomain)%nCells
634 0 : BFI_qT_sum(iDomain) = BFI_qT_sum(iDomain) &
635 0 : + sum(L1_total_runoff(s1 : e1) * level1(iDomain)%CellArea) / level1(iDomain)%nCells
636 : end if
637 :
638 4040544 : end subroutine mhm_interface_run_do_time_step
639 :
640 : !> \brief write output after current time-step
641 2020272 : subroutine mhm_interface_run_write_output()
642 : implicit none
643 : integer(i4) :: iDomain, tt, s0, e0
644 :
645 : ! get time step
646 2020272 : tt = run_cfg%time_step
647 :
648 : ! get domain index
649 2020272 : iDomain = run_cfg%get_domain_index(run_cfg%selected_domain)
650 :
651 2020272 : if ( (.not. optimize) .and. (run_cfg%domainDateTime%tIndex_out > 0_i4)) then
652 350760 : if (any(outputFlxState_mrm) .AND. (domainMeta%doRouting(iDomain))) then
653 :
654 35040 : if (run_cfg%domainDateTime%tIndex_out == 1) then
655 4 : run_cfg%nc_mrm = mRM_OutputDataset(iDomain, run_cfg%mask11)
656 : end if
657 :
658 : ! update Dataset (riv-temp as optional input)
659 35040 : if ( riv_temp_pcs%active ) then
660 : call mRM_updateDataset(run_cfg%nc_mrm, &
661 8760 : L11_qmod(run_cfg%s11 : run_cfg%e11), riv_temp_pcs%river_temp(riv_temp_pcs%s11 : riv_temp_pcs%e11))
662 : else
663 26280 : call mRM_updateDataset(run_cfg%nc_mrm, L11_qmod(run_cfg%s11 : run_cfg%e11))
664 : end if
665 :
666 : ! write data
667 35040 : if (run_cfg%domainDateTime%writeout(timeStep_model_outputs_mrm, tt)) then
668 732 : call run_cfg%nc_mrm%writeTimestep(run_cfg%domainDateTime%tIndex_out * timestep)
669 : end if
670 :
671 35040 : if(tt == run_cfg%domainDateTime%nTimeSteps) then
672 4 : call run_cfg%nc_mrm%close()
673 : end if
674 :
675 35040 : if ( gw_coupling ) then
676 : ! create
677 0 : if (run_cfg%domainDateTime%tIndex_out == 1) then
678 0 : run_cfg%nc_gw = GW_OutputDataset(iDomain, level0(iDomain)%mask)
679 : end if
680 : ! add data
681 0 : s0 = level0(iDomain)%iStart
682 0 : e0 = level0(iDomain)%iEnd
683 0 : call GW_updateDataset(run_cfg%nc_gw, L0_river_head_mon_sum(s0 : e0))
684 : ! write
685 0 : if (run_cfg%domainDateTime%writeout(-2, tt)) then ! -2 for monthly
686 0 : call run_cfg%nc_gw%writeTimestep(run_cfg%domainDateTime%tIndex_out * timestep)
687 : end if
688 : ! close
689 0 : if(tt == run_cfg%domainDateTime%nTimeSteps) then
690 0 : call run_cfg%nc_gw%close()
691 : end if
692 : end if
693 :
694 : end if
695 :
696 587280 : if (any(outputFlxState)) then
697 :
698 131520 : if (run_cfg%domainDateTime%tIndex_out == 1) then
699 38 : run_cfg%nc_mhm = mHM_OutputDataset(iDomain, run_cfg%mask1)
700 : end if
701 :
702 : call mHM_updateDataset(&
703 : run_cfg%nc_mhm, &
704 : run_cfg%s1, run_cfg%e1, &
705 : L1_fSealed(:, 1, run_cfg%domainDateTime%yId), &
706 : run_cfg%L1_fNotSealed(:, 1, run_cfg%domainDateTime%yId), &
707 : L1_inter, &
708 : L1_snowPack, &
709 : L1_soilMoist, &
710 : L1_soilMoistSat(:, :, run_cfg%domainDateTime%yId), &
711 : L1_sealSTW, &
712 : L1_unsatSTW, &
713 : L1_satSTW, &
714 : L1_neutrons, &
715 : L1_pet_calc, &
716 : L1_aETSoil, &
717 : L1_aETCanopy, &
718 : L1_aETSealed, &
719 : L1_total_runoff, &
720 : L1_runoffSeal, &
721 : L1_fastRunoff, &
722 : L1_slowRunoff, &
723 : L1_baseflow, &
724 : L1_percol, &
725 : L1_infilSoil, &
726 : L1_preEffect, &
727 : L1_melt &
728 131520 : )
729 :
730 : ! write data
731 131520 : if (run_cfg%domainDateTime%writeout(timeStep_model_outputs, tt)) then
732 93 : call run_cfg%nc_mhm%writeTimestep(run_cfg%domainDateTime%tIndex_out * timestep)
733 : end if
734 :
735 131520 : if(tt == run_cfg%domainDateTime%nTimeSteps) then
736 15 : call run_cfg%nc_mhm%close()
737 : end if
738 :
739 : end if
740 : end if ! <-- if (.not. optimize)
741 :
742 2020272 : end subroutine mhm_interface_run_write_output
743 :
744 : !> \brief add simulation data to optimization data types
745 4040544 : subroutine mhm_interface_run_update_optisim(etOptiSim, twsOptiSim, neutronsOptiSim, smOptiSim)
746 : implicit none
747 : !> returns soil moisture time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
748 : type(optidata_sim), dimension(:), optional, intent(inout) :: smOptiSim
749 : !> dim1=ncells, dim2=time
750 : type(optidata_sim), dimension(:), optional, intent(inout) :: neutronsOptiSim
751 : !> returns evapotranspiration time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
752 : type(optidata_sim), dimension(:), optional, intent(inout) :: etOptiSim
753 : !> returns tws time series for all grid cells (of multiple Domains concatenated),DIMENSION [nCells, nTimeSteps]
754 : type(optidata_sim), dimension(:), optional, intent(inout) :: twsOptiSim
755 :
756 : integer(i4) :: gg, s1, e1, lcId
757 :
758 : integer(i4) :: iDomain, tt
759 :
760 : ! get time step
761 2020272 : tt = run_cfg%time_step
762 :
763 : ! get domain index
764 2020272 : iDomain = run_cfg%get_domain_index(run_cfg%selected_domain)
765 :
766 : ! slice settings
767 2020272 : s1 = run_cfg%s1
768 2020272 : e1 = run_cfg%e1
769 2020272 : lcId = run_cfg%domainDateTime%yId
770 :
771 : !----------------------------------------------------------------------
772 : ! FOR SOIL MOISTURE
773 : ! NOTE:: modeled soil moisture is averaged according to input time step
774 : ! soil moisture (timeStep_sm_input)
775 : !----------------------------------------------------------------------
776 2020272 : if (present(smOptiSim)) then
777 : ! only for evaluation period - ignore warming days
778 78624 : if ((tt - warmingDays(iDomain) * nTstepDay) .GT. 0) then
779 : ! decide for daily, monthly or yearly aggregation
780 52560 : call smOptiSim(iDomain)%average_per_timestep(L1_smObs(iDomain)%timeStepInput, &
781 105120 : run_cfg%domainDateTime%is_new_day, run_cfg%domainDateTime%is_new_month, run_cfg%domainDateTime%is_new_year)
782 : ! last timestep is already done - write_counter exceeds size(smOptiSim(iDomain)%dataSim, dim=2)
783 52560 : if (tt /= run_cfg%domainDateTime%nTimeSteps) then
784 : ! aggregate soil moisture to needed time step for optimization
785 52554 : call smOptiSim(iDomain)%average_add(sum(L1_soilMoist(:, 1 : nSoilHorizons_sm_input), dim = 2) / &
786 5465616 : sum(L1_soilMoistSat(:, 1 : nSoilHorizons_sm_input, lcId), dim = 2))
787 : end if
788 : end if
789 : end if
790 :
791 : !----------------------------------------------------------------------
792 : ! FOR NEUTRONS
793 : ! NOTE:: modeled neutrons are averaged daily
794 : !----------------------------------------------------------------------
795 2020272 : if (present(neutronsOptiSim)) then
796 : ! only for evaluation period - ignore warming days
797 78624 : if ((tt - warmingDays(iDomain) * nTstepDay) .GT. 0) then
798 : ! decide for daily, monthly or yearly aggregation
799 : ! daily
800 52560 : if (run_cfg%domainDateTime%is_new_day) then
801 2190 : call neutronsOptiSim(iDomain)%average()
802 : end if
803 :
804 : ! last timestep is already done - write_counter exceeds size(sm_opti, dim=2)
805 52560 : if (tt /= run_cfg%domainDateTime%nTimeSteps) then
806 : ! aggregate neutrons to needed time step for optimization
807 52554 : call neutronsOptiSim(iDomain)%average_add(L1_neutrons(s1 : e1))
808 : end if
809 : end if
810 : end if
811 :
812 : !----------------------------------------------------------------------
813 : ! FOR EVAPOTRANSPIRATION
814 : ! NOTE:: modeled evapotranspiration is averaged according to input time step
815 : ! evapotranspiration (timeStep_et_input)
816 : !----------------------------------------------------------------------
817 2020272 : if (present(etOptiSim)) then
818 : ! only for evaluation period - ignore warming days
819 736344 : if ((tt - warmingDays(iDomain) * nTstepDay) .GT. 0) then
820 : ! decide for daily, monthly or yearly aggregation
821 736344 : call etOptiSim(iDomain)%increment_counter(L1_etObs(iDomain)%timeStepInput, &
822 1472688 : run_cfg%domainDateTime%is_new_day, run_cfg%domainDateTime%is_new_month, run_cfg%domainDateTime%is_new_year)
823 :
824 : ! last timestep is already done - write_counter exceeds size(etOptiSim(iDomain)%dataSim, dim=2)
825 736344 : if (tt /= run_cfg%domainDateTime%nTimeSteps) then
826 : ! aggregate evapotranspiration to needed time step for optimization
827 0 : call etOptiSim(iDomain)%add(sum(L1_aETSoil(s1 : e1, :), dim = 2) * &
828 0 : run_cfg%L1_fNotSealed(s1 : e1, 1, lcId) + &
829 0 : L1_aETCanopy(s1 : e1) + &
830 75841269 : L1_aETSealed(s1 : e1) * L1_fSealed(s1 : e1, 1, lcId))
831 : end if
832 : end if
833 : end if
834 :
835 : !----------------------------------------------------------------------
836 : ! FOR TWS
837 : ! NOTE:: modeled tws is averaged according to input time step
838 : ! (timeStepInput)
839 : !----------------------------------------------------------------------
840 2020272 : if (present(twsOptiSim)) then
841 : ! only for evaluation period - ignore warming days
842 841464 : if ((tt - warmingDays(iDomain) * nTstepDay) > 0) then
843 : ! decide for daily, monthly or yearly aggregation
844 841464 : call twsOptiSim(iDomain)%average_per_timestep(L1_twsaObs(iDomain)%timeStepInput, &
845 1682928 : run_cfg%domainDateTime%is_new_day, run_cfg%domainDateTime%is_new_month, run_cfg%domainDateTime%is_new_year)
846 :
847 : ! last timestep is already done - write_counter exceeds size(twsOptiSim(iDomain)%dataSim, dim=2)
848 841464 : if (tt /= run_cfg%domainDateTime%nTimeSteps) then
849 : ! aggregate evapotranspiration to needed time step for optimization
850 0 : call twsOptiSim(iDomain)%average_add(L1_inter(s1 : e1) + L1_snowPack(s1 : e1) + L1_sealSTW(s1 : e1) + &
851 29450295 : L1_unsatSTW(s1 : e1) + L1_satSTW(s1 : e1))
852 2524311 : do gg = 1, nSoilHorizons_mHM
853 2524311 : call twsOptiSim(iDomain)%add(L1_soilMoist (s1 : e1, gg))
854 : end do
855 : end if
856 : end if
857 : end if
858 :
859 : ! TODO-RIV-TEMP: add OptiSim for river temperature
860 :
861 2020272 : end subroutine mhm_interface_run_update_optisim
862 :
863 : !> \brief finalize current domain after running
864 82 : subroutine mhm_interface_run_finalize_domain()
865 : implicit none
866 :
867 : integer(i4) :: iDomain
868 :
869 : ! get domain index
870 82 : iDomain = run_cfg%get_domain_index(run_cfg%selected_domain)
871 :
872 46 : if (allocated(run_cfg%InflowDischarge)) deallocate(run_cfg%InflowDischarge)
873 82 : if (domainMeta%doRouting(iDomain)) then
874 : ! clean runoff variable
875 46 : deallocate(run_cfg%RunToRout)
876 46 : if ( riv_temp_pcs%active ) call riv_temp_pcs%dealloc_lateral()
877 : end if
878 :
879 2020272 : end subroutine mhm_interface_run_finalize_domain
880 :
881 : !> \brief finalize run
882 61 : subroutine mhm_interface_run_finalize(runoff, BFI)
883 : implicit none
884 : !> returns runoff time series, DIMENSION [nTimeSteps, nGaugesTotal]
885 : real(dp), dimension(:, :), allocatable, optional, intent(out) :: runoff
886 : !> baseflow index, dim1=domainID
887 : real(dp), dimension(:), allocatable, optional, intent(out) :: BFI
888 :
889 508027 : if (present(runoff) .and. (processMatrix(8, 1) > 0)) runoff = mRM_runoff
890 61 : if (present(BFI)) then
891 0 : BFI = BFI_qBF_sum / BFI_qT_sum
892 0 : deallocate(BFI_qBF_sum)
893 0 : deallocate(BFI_qT_sum)
894 : end if
895 :
896 61 : if (allocated(run_cfg%parameterset)) deallocate(run_cfg%parameterset)
897 61 : if (allocated(run_cfg%L1_fNotSealed)) deallocate(run_cfg%L1_fNotSealed)
898 61 : if (allocated(run_cfg%domain_indices)) deallocate(run_cfg%domain_indices)
899 :
900 82 : end subroutine mhm_interface_run_finalize
901 :
902 : end module mo_mhm_interface_run
|
