Line data Source code
1 : !> \file mo_mhm_interface.f90
2 : !> \brief \copybrief mo_mhm_interface
3 : !> \details \copydetails mo_mhm_interface
4 :
5 : !> \brief Module providing interfaces for mHM.
6 : !> \details Interfaces to control the mHM workflow from outside (init, run, get infos, etc.).
7 : !> \authors Sebastian Mueller
8 : !> \version 0.1
9 : !> \date Oct 2021
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
14 :
15 : use mo_kind, only: i4, dp
16 : use mo_message, only: message, error_message
17 : use mo_string_utils, only: num2str
18 :
19 : #ifdef MPI
20 : use mpi_f08
21 : #endif
22 :
23 : implicit none
24 :
25 : private
26 :
27 : public :: mhm_interface_init
28 : public :: mhm_interface_get_parameter
29 : public :: mhm_interface_get_parameter_number
30 : public :: mhm_interface_run
31 : public :: mhm_interface_run_optimization
32 : public :: mhm_interface_finalize
33 :
34 : contains
35 :
36 : !> \brief initialize mHM from given namelist paths.
37 28 : subroutine mhm_interface_init(namelist_mhm, namelist_mhm_param, namelist_mhm_output, namelist_mrm_output, cwd)
38 : use mo_file, only: &
39 : file_namelist_mhm, &
40 : file_namelist_mhm_param, &
41 : file_defOutput
42 : use mo_mrm_file, only: mrm_file_defOutput => file_defOutput
43 : use mo_common_read_config, only: &
44 : common_read_config
45 : use mo_common_mHM_mRM_read_config, only : &
46 : common_mHM_mRM_read_config, &
47 : check_optimization_settings
48 : use mo_mpr_read_config, only: mpr_read_config
49 : use mo_mhm_read_config, only: mhm_read_config
50 : use mo_read_wrapper, only : read_data
51 : use mo_mrm_init, only: &
52 : mrm_init, &
53 : mrm_configuration
54 : use mo_common_variables, only: &
55 : level0, &
56 : level1, &
57 : itimer, &
58 : domainMeta, &
59 : processMatrix
60 : use mo_common_mHM_mRM_variables, only : &
61 : timeStep, &
62 : simPer, &
63 : optimize, &
64 : opti_function, &
65 : mrm_coupling_mode, &
66 : read_restart
67 : use mo_mhm_messages, only: &
68 : startup_message, &
69 : domain_dir_check_message
70 : use mo_timer, only: &
71 : timers_init, &
72 : timer_start, &
73 : timer_stop, &
74 : timer_get
75 : use mo_startup, only: mhm_initialize
76 : use mo_global_variables, only: &
77 : couple_cfg, &
78 : meteo_handler, &
79 : L1_twsaObs, &
80 : L1_etObs, &
81 : L1_neutronsObs, &
82 : L1_smObs, &
83 : BFI_calc
84 : use mo_read_optional_data, only: readOptidataObs
85 : use mo_write_ascii, only: write_configfile
86 : use mo_mhm_bfi, only: calculate_BFI
87 : use mo_os, only: change_dir
88 :
89 : implicit none
90 :
91 : character(*), optional, intent(in) :: namelist_mhm !< path to mHM configuration namelist
92 : character(*), optional, intent(in) :: namelist_mhm_param !< path to mHM parameter namelist
93 : character(*), optional, intent(in) :: namelist_mhm_output !< path to mHM output namelist
94 : character(*), optional, intent(in) :: namelist_mrm_output !< path to mRM output namelist
95 : character(*), optional, intent(in) :: cwd !< desired working directory
96 :
97 : integer(i4) :: domainID, iDomain
98 :
99 : #ifdef MPI
100 : integer :: ierror
101 : integer(i4) :: nproc, rank
102 : #endif
103 :
104 : ! reset nml paths if wanted
105 0 : if (present(namelist_mhm)) file_namelist_mhm = namelist_mhm
106 14 : if (present(namelist_mhm_param)) file_namelist_mhm_param = namelist_mhm_param
107 14 : if (present(namelist_mhm_output)) file_defOutput = namelist_mhm_output
108 14 : if (present(namelist_mrm_output)) mrm_file_defOutput = namelist_mrm_output
109 : ! change working directory
110 14 : if (present(cwd)) call change_dir(cwd)
111 :
112 : ! startup message
113 14 : call startup_message()
114 :
115 : ! coupling configuration
116 14 : call couple_cfg%read_config(file_namelist_mhm)
117 : ! read configs
118 14 : call common_read_config(file_namelist_mhm)
119 14 : call mpr_read_config(file_namelist_mhm, file_namelist_mhm_param)
120 14 : call common_mHM_mRM_read_config(file_namelist_mhm)
121 14 : call mhm_read_config(file_namelist_mhm)
122 14 : call couple_cfg%check(domainMeta, optimize)
123 14 : call meteo_handler%config(file_namelist_mhm, optimize, domainMeta, processMatrix, timestep, couple_cfg)
124 14 : mrm_coupling_mode = 2_i4 ! TODO: this shouldn't be needed
125 14 : call mrm_configuration(file_namelist_mhm, file_namelist_mhm_param)
126 14 : if (optimize) call check_optimization_settings()
127 :
128 : ! Message about input directories
129 14 : call domain_dir_check_message()
130 :
131 : ! Start timings
132 14 : call timers_init
133 :
134 : ! --------------------------------------------------------------------------
135 : ! READ AND INITIALIZE
136 : ! --------------------------------------------------------------------------
137 14 : itimer = 1
138 : #ifdef MPI
139 : call MPI_Comm_size(domainMeta%comMaster, nproc, ierror)
140 : ! find the number the process is referred to, called rank
141 : call MPI_Comm_rank(domainMeta%comMaster, rank, ierror)
142 : ! ComLocal is a communicator, i.e. a group of processes assigned to the same
143 : ! domain, with a master and subprocesses. Only the master processes of these
144 : ! groups need to read the data. The master process with rank 0 only
145 : ! coordinates the other processes and does not need to read the data.
146 : if (rank > 0 .and. domainMeta%isMasterInComLocal) then
147 : #endif
148 14 : call message()
149 :
150 14 : if (.not. read_restart) then
151 13 : call message(' Read data ...')
152 13 : call timer_start(itimer)
153 : ! for DEM, slope, ... define nGvar local
154 : ! read_data has a domain loop inside
155 13 : call read_data(simPer)
156 13 : call timer_stop(itimer)
157 13 : call message(' in ', trim(num2str(timer_get(itimer), '(F9.3)')), ' seconds.')
158 : end if
159 :
160 : ! read data for every domain
161 14 : itimer = itimer + 1
162 14 : call message(' Initialize domains ...')
163 14 : call timer_start(itimer)
164 14 : call mhm_initialize()
165 14 : call meteo_handler%init_level2(level0, level1)
166 14 : call timer_stop(itimer)
167 14 : call message(' in ', trim(num2str(timer_get(itimer), '(F9.3)')), ' seconds.')
168 27 : if (processMatrix(8, 1) > 0) &
169 13 : call mrm_init(file_namelist_mhm, file_namelist_mhm_param)
170 :
171 14 : itimer = itimer + 1
172 14 : call message(' Read forcing and optional data ...')
173 14 : call timer_start(itimer)
174 :
175 40 : do iDomain = 1, domainMeta%nDomains
176 26 : domainID = domainMeta%indices(iDomain)
177 : ! read meteorology now, if it should be loaded in one go
178 26 : if (meteo_handler%single_read(iDomain)) call meteo_handler%prepare_data(1, iDomain, level1, simPer)
179 :
180 : ! read optional optional data if necessary
181 40 : if (optimize) then
182 1 : select case (opti_function)
183 : case(10 : 13, 28)
184 : ! read optional spatio-temporal soil mositure data
185 1 : call readOptidataObs(iDomain, domainID, L1_smObs(iDomain))
186 : case(17)
187 : ! read optional spatio-temporal neutrons data
188 1 : call readOptidataObs(iDomain, domainID, L1_neutronsObs(iDomain))
189 : case(27, 29, 30)
190 : ! read optional spatio-temporal evapotranspiration data
191 0 : call readOptidataObs(iDomain, domainID, L1_etObs(iDomain))
192 : case(15)
193 : ! read optional spatio-temporal tws data
194 1 : call readOptidataObs(iDomain, domainID, L1_twsaObs(iDomain))
195 : case(33)
196 : ! read optional spatio-temporal evapotranspiration data
197 6 : if (domainMeta%optidata(iDomain) == 0 .or. domainMeta%optidata(iDomain) == 5 .or. &
198 : domainMeta%optidata(iDomain) == 6 ) then
199 3 : call readOptidataObs(iDomain, domainID, L1_etObs(iDomain))
200 : end if
201 : ! read optional spatio-temporal tws data
202 6 : if (domainMeta%optidata(iDomain) == 0 .or. domainMeta%optidata(iDomain) == 3 .or. &
203 13 : domainMeta%optidata(iDomain) == 6 ) then
204 3 : call readOptidataObs(iDomain, domainID, L1_twsaObs(iDomain))
205 : end if
206 : end select
207 : end if
208 : end do
209 :
210 : ! calculate observed BFI if wanted
211 14 : if ( optimize .and. opti_function==34 .and. BFI_calc ) call calculate_BFI()
212 :
213 14 : call timer_stop(itimer)
214 14 : call message(' in ', trim(num2str(timer_get(itimer), '(F9.3)')), ' seconds.')
215 :
216 : !this call may be moved to another position as it writes the master config out file for all domains
217 : call write_configfile(meteo_handler%dirPrecipitation, meteo_handler%dirReferenceET, meteo_handler%dirTemperature)
218 :
219 : #ifdef MPI
220 : end if
221 : #endif
222 :
223 28 : end subroutine mhm_interface_init
224 :
225 : !> \brief Get current global parameter value of mHM.
226 0 : subroutine mhm_interface_get_parameter(para)
227 14 : use mo_common_variables, only: global_parameters
228 :
229 : implicit none
230 :
231 : real(dp), dimension(:), allocatable, intent(out) :: para !< global parameter values of mHM
232 :
233 0 : allocate(para(size(global_parameters, dim=1)))
234 :
235 0 : para = global_parameters(:, 3)
236 :
237 0 : end subroutine mhm_interface_get_parameter
238 :
239 : !> \brief Get number of current global parameter value of mHM.
240 0 : subroutine mhm_interface_get_parameter_number(n)
241 0 : use mo_common_variables, only: global_parameters
242 :
243 : implicit none
244 :
245 : integer(i4), intent(out) :: n !< number of global parameter values of mHM
246 :
247 0 : n = size(global_parameters, dim=1)
248 :
249 0 : end subroutine mhm_interface_get_parameter_number
250 :
251 : !> \brief Run mHM with current settings.
252 9 : subroutine mhm_interface_run()
253 : #ifdef MPI
254 : use mo_common_variables, only: domainMeta
255 : #endif
256 : use mo_common_variables, only: &
257 0 : itimer, &
258 : global_parameters
259 : use mo_timer, only: &
260 : timer_start, &
261 : timer_stop, &
262 : timer_get
263 : use mo_mhm_eval, only: mhm_eval
264 :
265 : implicit none
266 :
267 : #ifdef MPI
268 : integer :: ierror
269 : integer(i4) :: nproc, rank
270 :
271 : call MPI_Comm_size(domainMeta%comMaster, nproc, ierror)
272 : ! find the number the process is referred to, called rank
273 : call MPI_Comm_rank(domainMeta%comMaster, rank, ierror)
274 : #endif
275 :
276 9 : itimer = itimer + 1
277 :
278 : ! --------------------------------------------------------------------------
279 : ! call mHM
280 : ! get runoff timeseries if possible (i.e. when domainMeta%doRouting,
281 : ! processMatrix(8,1) > 0)
282 : ! get other model outputs (i.e. gridded fields of model output)
283 : ! --------------------------------------------------------------------------
284 :
285 : #ifdef MPI
286 : if (rank > 0 .and. domainMeta%isMasterInComLocal) then
287 : #endif
288 :
289 9 : call message(' Run mHM')
290 9 : call timer_start(itimer)
291 9 : call mhm_eval(global_parameters(:, 3))
292 9 : call timer_stop(itimer)
293 9 : call message(' in ', trim(num2str(timer_get(itimer), '(F12.3)')), ' seconds.')
294 :
295 : #ifdef MPI
296 : endif
297 : #endif
298 :
299 9 : end subroutine mhm_interface_run
300 :
301 : !> \brief Run mHM optimization with current settings.
302 5 : subroutine mhm_interface_run_optimization()
303 : #ifdef MPI
304 : use mo_common_variables, only: domainMeta
305 : use mo_objective_function, only: &
306 : objective_subprocess, &
307 : objective_master
308 : use mo_mrm_objective_function_runoff, only: &
309 : single_objective_runoff_master, &
310 : single_objective_runoff_subprocess
311 : #endif
312 : use mo_common_variables, only: &
313 : itimer, &
314 : dirConfigOut, &
315 : global_parameters,&
316 9 : global_parameters_name, &
317 : processMatrix
318 : use mo_common_mHM_mRM_variables, only : &
319 : opti_function
320 : use mo_timer, only: &
321 : timer_start, &
322 : timer_stop, &
323 : timer_get
324 : use mo_mhm_eval, only: mhm_eval
325 : use mo_objective_function, only: objective
326 : use mo_optimization, only: optimization
327 : use mo_mrm_objective_function_runoff, only: single_objective_runoff
328 : use mo_write_ascii, only: &
329 : write_optifile, & ! Writing optimized parameter set and objective
330 : write_optinamelist ! Writing optimized parameter set to a namelist
331 :
332 : implicit none
333 :
334 : procedure(mhm_eval), pointer :: eval
335 : procedure(objective), pointer :: obj_func
336 :
337 5 : real(dp) :: funcbest ! best objective function achivied during optimization
338 5 : logical, dimension(:), allocatable :: maskpara ! true = parameter will be optimized, = parameter(i,4) = 1
339 : ! ! false = parameter will not be optimized = parameter(i,4) = 0
340 :
341 : #ifdef MPI
342 : integer :: ierror
343 : integer(i4) :: nproc, rank
344 :
345 : call MPI_Comm_size(domainMeta%comMaster, nproc, ierror)
346 : ! find the number the process is referred to, called rank
347 : call MPI_Comm_rank(domainMeta%comMaster, rank, ierror)
348 : #endif
349 :
350 5 : itimer = itimer + 1
351 5 : call message(' Run mHM optimization')
352 5 : call timer_start(itimer)
353 :
354 5 : eval => mhm_eval
355 :
356 1 : select case(opti_function)
357 : case(1 : 9, 14, 31 : 32)
358 : ! call optimization against only runoff (no other variables)
359 1 : obj_func => single_objective_runoff
360 : #ifdef MPI
361 : if (rank == 0 .and. domainMeta%isMasterInComLocal) then
362 : obj_func => single_objective_runoff_master
363 : call optimization(eval, obj_func, dirConfigOut, funcBest, maskpara)
364 : else if (domainMeta%isMasterInComLocal) then
365 : ! In case of a master process from ComLocal, i.e. a master of a group of
366 : ! processes that are assigned to a single domain, this process calls the
367 : ! objective subroutine directly. The master over all processes collects
368 : ! the data and runs the dds/sce/other opti method.
369 : call single_objective_runoff_subprocess(eval)
370 : end if
371 : #else
372 1 : call optimization(eval, obj_func, dirConfigOut, funcBest, maskpara)
373 : #endif
374 :
375 : case(10 : 13, 15, 17, 27, 28, 29, 30, 33, 34)
376 : ! call optimization for other variables
377 4 : obj_func => objective
378 : #ifdef MPI
379 : if (rank == 0 .and. domainMeta%isMasterInComLocal) then
380 : obj_func => objective_master
381 : call optimization(eval, obj_func, dirConfigOut, funcBest, maskpara)
382 : else if (domainMeta%isMasterInComLocal) then
383 : ! In case of a master process from ComLocal, i.e. a master of a group of
384 : ! processes that are assigned to a single domain, this process calls the
385 : ! objective subroutine directly. The master over all processes collects
386 : ! the data and runs the dds/sce/other opti method.
387 : call objective_subprocess(eval)
388 : end if
389 : #else
390 4 : call optimization(eval, obj_func, dirConfigOut, funcBest, maskpara)
391 : #endif
392 :
393 : case default
394 : call error_message('***ERROR: mhm_driver: The given objective function number ', &
395 10 : trim(adjustl(num2str(opti_function))), ' in mhm.nml is not valid!')
396 : end select
397 :
398 : #ifdef MPI
399 : if (rank == 0 .and. domainMeta%isMasterInComLocal) then
400 : #endif
401 :
402 : ! write a file with final objective function and the best parameter set
403 5 : call write_optifile(funcbest, global_parameters(:, 3), global_parameters_name(:))
404 : ! write a file with final best parameter set in a namlist format
405 5 : call write_optinamelist(processMatrix, global_parameters, maskpara, global_parameters_name(:))
406 5 : deallocate(maskpara)
407 :
408 : #ifdef MPI
409 : end if
410 : #endif
411 :
412 5 : call timer_stop(itimer)
413 5 : call message(' in ', trim(num2str(timer_get(itimer), '(F12.3)')), ' seconds.')
414 :
415 5 : end subroutine mhm_interface_run_optimization
416 :
417 : !> \brief Write mHM restart.
418 14 : subroutine mhm_interface_finalize()
419 : #ifdef MPI
420 : use mo_common_variables, only: domainMeta
421 : #endif
422 : use mo_common_variables, only: &
423 : itimer, &
424 : mhmFileRestartOut, &
425 5 : write_restart, &
426 : processMatrix
427 : use mo_common_mHM_mRM_variables, only : &
428 : optimize
429 : use mo_timer, only: &
430 : timer_start, &
431 : timer_stop, &
432 : timer_get
433 : use mo_restart, only: write_restart_files
434 : use mo_mrm_write, only : mrm_write
435 : use mo_mhm_messages, only: finish_message
436 : use mo_clean_up, only: deallocate_global_variables
437 :
438 : implicit none
439 :
440 : #ifdef MPI
441 : integer :: ierror
442 : integer(i4) :: nproc, rank
443 :
444 : call MPI_Comm_size(domainMeta%comMaster, nproc, ierror)
445 : ! find the number the process is referred to, called rank
446 : call MPI_Comm_rank(domainMeta%comMaster, rank, ierror)
447 : if (rank > 0 .and. domainMeta%isMasterInComLocal) then
448 : #endif
449 :
450 : ! --------------------------------------------------------------------------
451 : ! WRITE RESTART files
452 : ! --------------------------------------------------------------------------
453 23 : if (write_restart .AND. (.NOT. optimize)) then
454 9 : itimer = itimer + 1
455 9 : call message()
456 9 : call message(' Write restart file')
457 9 : call timer_start(itimer)
458 : call write_restart_files(mhmFileRestartOut)
459 9 : call timer_stop(itimer)
460 9 : call message(' in ', trim(num2str(timer_get(itimer), '(F9.3)')), ' seconds.')
461 : end if
462 :
463 : ! --------------------------------------------------------------------------
464 : ! WRITE RUNOFF (INCLUDING RESTART FILES, has to be called after mHM restart
465 : ! files are written)
466 : ! --------------------------------------------------------------------------
467 14 : if (processMatrix(8, 1) > 0) call mrm_write()
468 :
469 : #ifdef MPI
470 : end if
471 : #endif
472 :
473 14 : call finish_message()
474 :
475 : ! clean up all allocated variables
476 14 : call deallocate_global_variables()
477 :
478 14 : end subroutine mhm_interface_finalize
479 :
480 : end module mo_mhm_interface
|
