1!--------------------------------------------------------------------------------------------------!
2!   CP2K: A general program to perform molecular dynamics simulations                              !
3!   Copyright (C) 2000 - 2020  CP2K developers group                                               !
4!--------------------------------------------------------------------------------------------------!
5
6! **************************************************************************************************
7!> \brief module contains the worker routine handling the communication and
8!>        the calculation / creation of the configurations
9!>        - WORKER these are all TMC cores, instead of master core
10!>          and maybe some idle cores
11!>        - divided in groups, in every group exists group master
12!>          - there can be two kind of groups, one for exact energy calculation
13!>            and one calculating configurational change using an approximate
14!>            potential
15!>        - Algorithm:
16!>          - group master receive messages and decide what to do,
17!>          - (if nessesary) broadcast of working task
18!>            to all other group members (needed for parallel CP2K)
19!>          - process task, calculations of energy or configurational change
20!>          - result, exist on group master, sent to master core
21!>        Communication structure (master->worker, worker->master):
22!>        - message structure is defined in TMC message module
23!> \par History
24!>      11.2012 created [Mandes Schoenherr]
25!> \author Mandes
26! **************************************************************************************************
27
28MODULE tmc_worker
29   USE cell_types,                      ONLY: cell_copy,&
30                                              cell_type,&
31                                              init_cell
32   USE cp_external_control,             ONLY: set_external_comm
33   USE cp_log_handling,                 ONLY: cp_to_string
34   USE cp_para_types,                   ONLY: cp_para_env_type
35   USE cp_result_methods,               ONLY: cp_results_erase,&
36                                              put_results
37   USE cp_result_types,                 ONLY: cp_result_type
38   USE cp_subsys_types,                 ONLY: cp_subsys_get,&
39                                              cp_subsys_type
40   USE f77_interface,                   ONLY: f_env_get_from_id,&
41                                              f_env_type,&
42                                              get_natom,&
43                                              get_pos,&
44                                              get_result_r1
45   USE force_env_types,                 ONLY: force_env_get,&
46                                              force_env_get_natom
47   USE kinds,                           ONLY: default_string_length,&
48                                              dp
49   USE molecule_list_types,             ONLY: molecule_list_type
50   USE particle_list_types,             ONLY: particle_list_type
51   USE tmc_analysis,                    ONLY: analysis_init,&
52                                              analysis_restart_print,&
53                                              analysis_restart_read,&
54                                              analyze_file_configurations,&
55                                              do_tmc_analysis,&
56                                              finalize_tmc_analysis
57   USE tmc_analysis_types,              ONLY: tmc_ana_list_type
58   USE tmc_calculations,                ONLY: calc_potential_energy
59   USE tmc_messages,                    ONLY: bcast_group,&
60                                              check_if_group_master,&
61                                              communicate_atom_types,&
62                                              master_comm_id,&
63                                              recv_msg,&
64                                              send_msg,&
65                                              stop_whole_group,&
66                                              tmc_message
67   USE tmc_move_handle,                 ONLY: clear_move_probs,&
68                                              prob_update,&
69                                              select_random_move_type
70   USE tmc_move_types,                  ONLY: mv_type_MD,&
71                                              mv_type_NMC_moves
72   USE tmc_moves,                       ONLY: change_pos
73   USE tmc_stati,                       ONLY: &
74        TMC_CANCELING_MESSAGE, TMC_CANCELING_RECEIPT, TMC_STATUS_CALCULATING, TMC_STATUS_FAILED, &
75        TMC_STATUS_STOP_RECEIPT, TMC_STATUS_WAIT_FOR_NEW_TASK, TMC_STATUS_WORKER_INIT, &
76        TMC_STAT_ANALYSIS_REQUEST, TMC_STAT_ANALYSIS_RESULT, TMC_STAT_APPROX_ENERGY_REQUEST, &
77        TMC_STAT_APPROX_ENERGY_RESULT, TMC_STAT_ENERGY_REQUEST, TMC_STAT_ENERGY_RESULT, &
78        TMC_STAT_INIT_ANALYSIS, TMC_STAT_MD_REQUEST, TMC_STAT_MD_RESULT, TMC_STAT_NMC_REQUEST, &
79        TMC_STAT_NMC_RESULT, TMC_STAT_SCF_STEP_ENER_RECEIVE, TMC_STAT_START_CONF_REQUEST, &
80        TMC_STAT_START_CONF_RESULT, task_type_MC, task_type_ideal_gas
81   USE tmc_tree_acceptance,             ONLY: acceptance_check
82   USE tmc_tree_build,                  ONLY: allocate_new_sub_tree_node,&
83                                              deallocate_sub_tree_node
84   USE tmc_tree_types,                  ONLY: tree_type
85   USE tmc_types,                       ONLY: allocate_tmc_atom_type,&
86                                              tmc_atom_type,&
87                                              tmc_env_type,&
88                                              tmc_param_type
89#include "../base/base_uses.f90"
90
91   IMPLICIT NONE
92
93   PRIVATE
94
95   CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'tmc_worker'
96
97   PUBLIC :: do_tmc_worker
98   PUBLIC :: get_initial_conf, get_atom_kinds_and_cell
99
100   INTEGER, PARAMETER :: DEBUG = 0
101
102CONTAINS
103
104! **************************************************************************************************
105!> \brief worker get tasks form master and fulfill them
106!> \param tmc_env structure for storing all the tmc parameters
107!> \param ana_list ...
108!> \author Mandes 11.2012
109! **************************************************************************************************
110   SUBROUTINE do_tmc_worker(tmc_env, ana_list)
111      TYPE(tmc_env_type), POINTER                        :: tmc_env
112      TYPE(tmc_ana_list_type), DIMENSION(:), OPTIONAL, &
113         POINTER                                         :: ana_list
114
115      CHARACTER(LEN=*), PARAMETER :: routineN = 'do_tmc_worker', routineP = moduleN//':'//routineN
116
117      CHARACTER(LEN=default_string_length)               :: c_tmp
118      INTEGER                                            :: calc_stat, handle, i1, i2, ierr, itmp, &
119                                                            num_dim, work_stat
120      INTEGER, DIMENSION(:), POINTER                     :: ana_restart_conf
121      LOGICAL                                            :: flag, master
122      TYPE(cp_para_env_type), POINTER                    :: para_env_m_w
123      TYPE(tree_type), POINTER                           :: conf
124
125      master = .FALSE.
126      i1 = -1
127      i2 = -1
128      NULLIFY (conf, para_env_m_w, ana_restart_conf)
129
130      CPASSERT(ASSOCIATED(tmc_env))
131
132      ! start the timing
133      CALL timeset(routineN, handle)
134
135      ! initialize
136      IF (tmc_env%tmc_comp_set%group_nr .GT. 0) THEN
137         CPASSERT(ASSOCIATED(tmc_env%tmc_comp_set%para_env_sub_group))
138         IF (tmc_env%w_env%env_id_ener .GT. 0) THEN
139            itmp = tmc_env%w_env%env_id_ener
140         ELSE
141            itmp = tmc_env%w_env%env_id_approx
142         END IF
143
144         CALL get_atom_kinds_and_cell(env_id=itmp, &
145                                      atoms=tmc_env%params%atoms, cell=tmc_env%params%cell)
146         para_env_m_w => tmc_env%tmc_comp_set%para_env_m_w
147         master = check_if_group_master(tmc_env%tmc_comp_set%para_env_sub_group)
148      ELSE
149         ! analysis group
150         CPASSERT(ASSOCIATED(tmc_env%tmc_comp_set%para_env_m_ana))
151         para_env_m_w => tmc_env%tmc_comp_set%para_env_m_ana
152         master = .TRUE.
153      END IF
154
155      !-- GROUP MASTER only --------------
156      ! get messages from master and handle them
157      IF (master) THEN
158         ! NOT the analysis group
159         IF (tmc_env%tmc_comp_set%group_nr .GT. 0) THEN
160            IF (tmc_env%w_env%env_id_ener .GT. 0) THEN
161               itmp = tmc_env%w_env%env_id_ener
162            ELSE
163               itmp = tmc_env%w_env%env_id_approx
164            END IF
165            ! set the communicator in the external control for receiving exit tags
166            !  and sending additional information (e.g. the intermediate scf energies)
167            IF (tmc_env%params%use_scf_energy_info) &
168               CALL set_intermediate_info_comm(env_id=itmp, &
169                                               comm=tmc_env%tmc_comp_set%para_env_m_w%group)
170            IF (tmc_env%params%SPECULATIVE_CANCELING) &
171               CALL set_external_comm(comm=tmc_env%tmc_comp_set%para_env_m_w%group, &
172                                      in_external_master_id=MASTER_COMM_ID, &
173                                      in_exit_tag=TMC_CANCELING_MESSAGE)
174         END IF
175         !-- WORKING LOOP --!
176         master_work_time: DO
177            work_stat = TMC_STATUS_WAIT_FOR_NEW_TASK
178            ! -- receive message from master
179            ! check for new task (wait for it)
180            itmp = MASTER_COMM_ID
181            CALL tmc_message(msg_type=work_stat, send_recv=recv_msg, &
182                             dest=itmp, &
183                             para_env=para_env_m_w, &
184                             result_count=ana_restart_conf, &
185                             tmc_params=tmc_env%params, elem=conf)
186
187            IF (DEBUG .GE. 1 .AND. work_stat .NE. TMC_STATUS_WAIT_FOR_NEW_TASK) &
188               WRITE (tmc_env%w_env%io_unit, *) "worker: group master of group ", &
189               tmc_env%tmc_comp_set%group_nr, "got task ", work_stat
190            calc_stat = TMC_STATUS_CALCULATING
191            SELECT CASE (work_stat)
192            CASE (TMC_STATUS_WAIT_FOR_NEW_TASK)
193            CASE (TMC_STATUS_WORKER_INIT)
194               CALL init_cell(cell=tmc_env%params%cell)
195               itmp = bcast_group
196               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
197                                dest=itmp, &
198                                para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
199                                tmc_params=tmc_env%params)
200            CASE (TMC_CANCELING_MESSAGE)
201               work_stat = TMC_CANCELING_RECEIPT
202               itmp = MASTER_COMM_ID
203               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
204                                dest=itmp, &
205                                para_env=para_env_m_w, &
206                                tmc_params=tmc_env%params)
207            CASE (TMC_STATUS_FAILED)
208               IF (DEBUG .GE. 1) &
209                  WRITE (tmc_env%w_env%io_unit, *) "master worker of group", &
210                  tmc_env%tmc_comp_set%group_nr, " exit work time."
211               EXIT master_work_time
212               !-- group master read the CP2K input file, and write data to master
213            CASE (TMC_STAT_START_CONF_REQUEST)
214               IF (tmc_env%w_env%env_id_ener .GT. 0) THEN
215                  itmp = tmc_env%w_env%env_id_ener
216               ELSE
217                  itmp = tmc_env%w_env%env_id_approx
218               END IF
219               CALL get_initial_conf(tmc_params=tmc_env%params, init_conf=conf, &
220                                     env_id=itmp)
221               ! send start configuration back to master
222               work_stat = TMC_STAT_START_CONF_RESULT
223               itmp = MASTER_COMM_ID
224               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
225                                dest=itmp, &
226                                para_env=para_env_m_w, &
227                                tmc_params=tmc_env%params, elem=conf, &
228                                wait_for_message=.TRUE.)
229
230               IF (ASSOCIATED(tmc_env%tmc_comp_set%para_env_m_first_w)) &
231                  CALL communicate_atom_types(atoms=tmc_env%params%atoms, &
232                                              source=1, &
233                                              para_env=tmc_env%tmc_comp_set%para_env_m_first_w)
234               !-- calculate the approximate energy
235            CASE (TMC_STAT_APPROX_ENERGY_REQUEST)
236               CPASSERT(tmc_env%w_env%env_id_approx .GT. 0)
237               itmp = bcast_group
238               !-- DISTRIBUTING WORK (group master) to all other group members
239               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
240                                dest=itmp, &
241                                para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
242                                tmc_params=tmc_env%params, elem=conf)
243               CALL calc_potential_energy(conf=conf, &
244                                          env_id=tmc_env%w_env%env_id_approx, &
245                                          exact_approx_pot=.FALSE., &
246                                          tmc_env=tmc_env)
247               work_stat = TMC_STAT_APPROX_ENERGY_RESULT
248               itmp = MASTER_COMM_ID
249               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
250                                dest=itmp, &
251                                para_env=para_env_m_w, &
252                                tmc_params=tmc_env%params, elem=conf)
253               ! -- Nested Monte Carlo routines
254            CASE (TMC_STAT_MD_REQUEST, TMC_STAT_NMC_REQUEST)
255               CALL clear_move_probs(tmc_env%params%nmc_move_types)
256               itmp = bcast_group
257               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
258                                dest=itmp, &
259                                para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
260                                tmc_params=tmc_env%params, elem=conf)
261               !-- collective calculation for MD/NMC steps
262               IF (work_stat .EQ. TMC_STAT_NMC_REQUEST) THEN
263                  !-- calculate MD steps, in case of 2 different potentials do nested Monte Carlo
264                  CALL nested_markov_chain_MC(conf=conf, &
265                                              env_id=tmc_env%w_env%env_id_approx, &
266                                              tmc_env=tmc_env, calc_status=calc_stat)
267               ELSEIF (work_stat .EQ. TMC_STAT_MD_REQUEST) THEN
268                  !TODO Hybrid MC routine
269                  CPABORT("there is no Hybrid MC implemented yet.")
270
271               ELSE
272                  CPABORT("unknown task type for workers.")
273               END IF
274               !-- in case of cancelation send receipt
275               itmp = MASTER_COMM_ID
276               CALL tmc_message(msg_type=calc_stat, send_recv=recv_msg, &
277                                dest=itmp, &
278                                para_env=para_env_m_w, &
279                                tmc_params=tmc_env%params, &
280                                success=flag)
281               SELECT CASE (calc_stat)
282               CASE (TMC_STATUS_CALCULATING)
283                  SELECT CASE (work_stat)
284                  CASE (TMC_STAT_MD_REQUEST)
285                     work_stat = TMC_STAT_MD_RESULT
286                  CASE (TMC_STAT_NMC_REQUEST)
287                     work_stat = TMC_STAT_NMC_RESULT
288                  CASE DEFAULT
289                     CALL cp_abort(__LOCATION__, &
290                                   "unknown work status after possible NMC subgroup "// &
291                                   "cancelation, work_stat="//cp_to_string(work_stat))
292                  END SELECT
293               CASE (TMC_CANCELING_MESSAGE)
294                  work_stat = TMC_CANCELING_RECEIPT
295               CASE DEFAULT
296                  CALL cp_abort(__LOCATION__, &
297                                "unknown calc status before sending NMC result "// &
298                                cp_to_string(calc_stat))
299               END SELECT
300               ! send message back to master
301               itmp = MASTER_COMM_ID
302               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
303                                dest=itmp, &
304                                para_env=para_env_m_w, &
305                                tmc_params=tmc_env%params, elem=conf)
306            CASE (TMC_STAT_ENERGY_REQUEST)
307               CPASSERT(tmc_env%w_env%env_id_ener .GT. 0)
308               !-- DISTRIBUTING WORK (group master) to all other group members
309               itmp = bcast_group
310               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
311                                dest=itmp, &
312                                para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
313                                tmc_params=tmc_env%params, elem=conf)
314
315               CALL calc_potential_energy(conf=conf, &
316                                          env_id=tmc_env%w_env%env_id_ener, &
317                                          exact_approx_pot=.TRUE., &
318                                          tmc_env=tmc_env)
319               !-- in case of cancelation send receipt
320               flag = .FALSE.
321               itmp = MASTER_COMM_ID
322               CALL tmc_message(msg_type=calc_stat, send_recv=recv_msg, &
323                                dest=itmp, &
324                                para_env=para_env_m_w, &
325                                tmc_params=tmc_env%params, success=flag)
326               SELECT CASE (calc_stat)
327               CASE (TMC_STATUS_CALCULATING)
328                  SELECT CASE (work_stat)
329                  CASE (TMC_STAT_ENERGY_REQUEST)
330                     work_stat = TMC_STAT_ENERGY_RESULT
331                     !-- if nessesary get the exact dipoles (for e.g. quantum potential)
332                     IF (tmc_env%params%print_dipole) THEN
333                        c_tmp = "[DIPOLE]"
334                        CALL get_result_r1(env_id=tmc_env%w_env%env_id_ener, &
335                                           description=c_tmp, N=3, RESULT=conf%dipole, &
336                                           res_exist=flag, ierr=ierr)
337                        IF (.NOT. flag) tmc_env%params%print_dipole = .FALSE.
338                        ! TODO maybe let run with the changed option, but inform user properly
339                        IF (.NOT. flag) &
340                           CALL cp_abort(__LOCATION__, &
341                                         "TMC: The requested dipoles are not porvided by the "// &
342                                         "force environment.")
343                     END IF
344                  CASE DEFAULT
345                     CALL cp_abort(__LOCATION__, &
346                                   "energy worker should handle unknown stat "// &
347                                   cp_to_string(work_stat))
348                  END SELECT
349               CASE (TMC_CANCELING_MESSAGE)
350                  work_stat = TMC_CANCELING_RECEIPT
351               CASE DEFAULT
352                  CALL cp_abort(__LOCATION__, &
353                                "worker while energy calc is in unknown state "// &
354                                cp_to_string(work_stat))
355               END SELECT
356
357               !-- send information back to master
358               IF (DEBUG .GE. 1) &
359                  WRITE (tmc_env%w_env%io_unit, *) "worker group ", &
360                  tmc_env%tmc_comp_set%group_nr, &
361                  "calculations done, send result energy", conf%potential
362               itmp = MASTER_COMM_ID
363               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
364                                dest=itmp, &
365                                para_env=para_env_m_w, &
366                                tmc_params=tmc_env%params, elem=conf)
367            CASE (TMC_STAT_INIT_ANALYSIS)
368               CPASSERT(ASSOCIATED(ana_restart_conf))
369               CPASSERT(SIZE(ana_restart_conf) .EQ. tmc_env%params%nr_temp)
370               CPASSERT(PRESENT(ana_list))
371               CPASSERT(ASSOCIATED(ana_list))
372               itmp = MASTER_COMM_ID
373               CALL communicate_atom_types(atoms=tmc_env%params%atoms, &
374                                           source=itmp, para_env=tmc_env%tmc_comp_set%para_env_m_ana)
375
376               num_dim = SIZE(conf%pos)
377               DO itmp = 1, tmc_env%params%nr_temp
378                  ! do not forget to nullify the pointer at the end, deallcoated at tmc_env%params
379                  ana_list(itmp)%temp%temperature = tmc_env%params%Temp(itmp)
380                  ana_list(itmp)%temp%atoms => tmc_env%params%atoms
381                  ana_list(itmp)%temp%cell => tmc_env%params%cell
382!              ana_list(itmp)%temp%io_unit     = tmc_env%w_env%io_unit
383
384                  CALL analysis_init(ana_env=ana_list(itmp)%temp, nr_dim=num_dim)
385                  ana_list(itmp)%temp%print_test_output = tmc_env%params%print_test_output
386                  IF (.NOT. ASSOCIATED(conf)) &
387                     CALL allocate_new_sub_tree_node(tmc_params=tmc_env%params, &
388                                                     next_el=conf, nr_dim=num_dim)
389                  CALL analysis_restart_read(ana_env=ana_list(itmp)%temp, &
390                                             elem=conf)
391                  !check if we have the read the file
392                  flag = .FALSE.
393                  IF ((.NOT. ASSOCIATED(ana_list(itmp)%temp%last_elem)) .AND. &
394                      ana_restart_conf(itmp) .GT. 0) THEN
395                     flag = .TRUE.
396                     i1 = 0
397                     i2 = ana_restart_conf(itmp)
398                     CALL cp_warn(__LOCATION__, &
399                                  "analysis old trajectory up to "// &
400                                  "elem "//cp_to_string(ana_restart_conf(itmp))// &
401                                  ". Read trajectory file.")
402                  ELSE IF (ASSOCIATED(ana_list(itmp)%temp%last_elem)) THEN
403                     IF (.NOT. (ana_list(itmp)%temp%last_elem%nr .EQ. ana_restart_conf(itmp))) THEN
404                        flag = .TRUE.
405                        i1 = ana_list(itmp)%temp%last_elem%nr
406                        i2 = ana_restart_conf(itmp)
407                        CALL cp_warn(__LOCATION__, &
408                                     "analysis restart with the incorrect configuration "// &
409                                     "TMC "//cp_to_string(ana_restart_conf(itmp))// &
410                                     " ana "//cp_to_string(ana_list(itmp)%temp%last_elem%nr)// &
411                                     ". REread trajectory file.")
412                     END IF
413                  END IF
414                  IF (flag) THEN
415                     CALL analyze_file_configurations(start_id=i1, &
416                                                      end_id=i2, &
417                                                      ana_env=ana_list(itmp)%temp, &
418                                                      tmc_params=tmc_env%params)
419                  END IF
420               END DO
421            CASE (TMC_STAT_ANALYSIS_REQUEST)
422               CPASSERT(PRESENT(ana_list))
423               CPASSERT(ASSOCIATED(ana_list(conf%sub_tree_nr)%temp))
424               CALL do_tmc_analysis(elem=conf, &
425                                    ana_env=ana_list(conf%sub_tree_nr)%temp)
426               work_stat = TMC_STAT_ANALYSIS_RESULT
427               itmp = MASTER_COMM_ID
428               CALL tmc_message(msg_type=work_stat, send_recv=send_msg, &
429                                dest=itmp, &
430                                para_env=para_env_m_w, &
431                                tmc_params=tmc_env%params, elem=conf)
432            CASE DEFAULT
433               CALL cp_abort(__LOCATION__, &
434                             "worker received unknown message task type "// &
435                             cp_to_string(work_stat))
436            END SELECT
437
438            IF (DEBUG .GE. 1 .AND. work_stat .NE. TMC_STATUS_WAIT_FOR_NEW_TASK) &
439               WRITE (tmc_env%w_env%io_unit, *) "worker: group ", &
440               tmc_env%tmc_comp_set%group_nr, &
441               "send back status:", work_stat
442            IF (ASSOCIATED(conf)) &
443               CALL deallocate_sub_tree_node(tree_elem=conf)
444         END DO master_work_time
445         !-- every other group paricipants----------------------------------------
446      ELSE
447         worker_work_time: DO
448            work_stat = TMC_STATUS_WAIT_FOR_NEW_TASK
449            flag = .FALSE.
450            itmp = bcast_group
451            CALL tmc_message(msg_type=work_stat, send_recv=recv_msg, &
452                             dest=itmp, &
453                             para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
454                             tmc_params=tmc_env%params, elem=conf)
455            calc_stat = TMC_STATUS_CALCULATING
456            SELECT CASE (work_stat)
457            CASE (TMC_STATUS_WORKER_INIT)
458               CALL init_cell(cell=tmc_env%params%cell)
459            CASE (TMC_CANCELING_MESSAGE)
460               ! error message
461            CASE (TMC_STATUS_FAILED)
462               EXIT worker_work_time
463               ! all group members have to calculate the (MD potential) energy together
464            CASE (TMC_STAT_START_CONF_RESULT)
465               CPASSERT(tmc_env%w_env%env_id_approx .GT. 0)
466               !-- collective calculation of the potential energy of MD potential
467               SELECT CASE (tmc_env%params%task_type)
468               CASE (task_type_MC, task_type_ideal_gas)
469                  IF (tmc_env%params%NMC_inp_file .NE. "") THEN
470                     conf%box_scale(:) = 1.0_dp
471                     CALL calc_potential_energy(conf=conf, &
472                                                env_id=tmc_env%w_env%env_id_approx, &
473                                                exact_approx_pot=.FALSE., &
474                                                tmc_env=tmc_env)
475                  END IF
476               CASE DEFAULT
477                  CALL cp_abort(__LOCATION__, &
478                                "unknown task_type for participants in "// &
479                                "START_CONF_RESULT request ")
480               END SELECT
481               !-- HMC - calculating MD steps
482            CASE (TMC_STAT_NMC_REQUEST, TMC_STAT_MD_REQUEST)
483               !-- collective calculation for MD/NMC steps
484               IF (work_stat .EQ. TMC_STAT_NMC_REQUEST) THEN
485                  !-- calculate MD steps, in case of 2 different potentials do nested Monte Carlo
486                  CALL nested_markov_chain_MC(conf=conf, &
487                                              env_id=tmc_env%w_env%env_id_approx, &
488                                              tmc_env=tmc_env, calc_status=calc_stat)
489               ELSEIF (work_stat .EQ. TMC_STAT_MD_REQUEST) THEN
490                  !TODO Hybrid MC routine
491                  CPABORT("there is no Hybrid MC implemented yet.")
492
493               ELSE
494                  CPABORT("unknown task type for workers.")
495               END IF
496               !-- energy calculations
497            CASE (TMC_STAT_APPROX_ENERGY_REQUEST)
498               !--- do calculate energy
499               CPASSERT(tmc_env%w_env%env_id_approx .GT. 0)
500               CALL calc_potential_energy(conf=conf, &
501                                          env_id=tmc_env%w_env%env_id_approx, &
502                                          exact_approx_pot=.FALSE., &
503                                          tmc_env=tmc_env)
504            CASE (TMC_STAT_ENERGY_REQUEST)
505               !--- do calculate energy
506               CPASSERT(tmc_env%w_env%env_id_ener .GT. 0)
507               CALL calc_potential_energy(conf=conf, &
508                                          env_id=tmc_env%w_env%env_id_ener, &
509                                          exact_approx_pot=.TRUE., &
510                                          tmc_env=tmc_env)
511            CASE DEFAULT
512               CALL cp_abort(__LOCATION__, &
513                             "group participant got unknown working type "// &
514                             cp_to_string(work_stat))
515            END SELECT
516            IF (ASSOCIATED(conf)) &
517               CALL deallocate_sub_tree_node(tree_elem=conf)
518         END DO worker_work_time
519      END IF
520      ! --------------------------------------------------------------------
521      ! finalizing analysis, writing files etc.
522      IF (ASSOCIATED(tmc_env%tmc_comp_set%para_env_m_ana)) THEN
523         DO itmp = 1, tmc_env%params%nr_temp
524            CALL analysis_restart_print(ana_env=ana_list(itmp)%temp)
525            IF (ASSOCIATED(conf)) &
526               CALL deallocate_sub_tree_node(tree_elem=ana_list(itmp)%temp%last_elem)
527            CALL finalize_tmc_analysis(ana_list(itmp)%temp)
528         END DO
529      END IF
530      !-- stopping and finalizing
531      ! sending back receipt for stopping
532      IF (master) THEN
533         ! NOT the analysis group
534         IF (tmc_env%tmc_comp_set%group_nr .GT. 0) THEN
535            ! remove the communicator in the external control for receiving exit tags
536            !  and sending additional information (e.g. the intermediate scf energies)
537            IF (tmc_env%params%use_scf_energy_info) THEN
538               IF (tmc_env%w_env%env_id_ener .GT. 0) THEN
539                  itmp = tmc_env%w_env%env_id_ener
540               ELSE
541                  itmp = tmc_env%w_env%env_id_approx
542               END IF
543               CALL remove_intermediate_info_comm(env_id=itmp)
544            END IF
545         END IF
546         IF (ASSOCIATED(tmc_env%tmc_comp_set%para_env_sub_group)) &
547            CALL stop_whole_group(para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
548                                  tmc_params=tmc_env%params)
549
550         work_stat = TMC_STATUS_STOP_RECEIPT
551         itmp = MASTER_COMM_ID
552         CALL tmc_message(msg_type=work_stat, send_recv=send_msg, dest=itmp, &
553                          para_env=para_env_m_w, &
554                          tmc_params=tmc_env%params)
555      ELSE IF (ASSOCIATED(tmc_env%tmc_comp_set%para_env_sub_group)) THEN
556         work_stat = TMC_STATUS_STOP_RECEIPT
557         itmp = MASTER_COMM_ID
558         CALL tmc_message(msg_type=work_stat, send_recv=send_msg, dest=itmp, &
559                          para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
560                          tmc_params=tmc_env%params)
561      END IF
562
563      IF (DEBUG .GE. 5) &
564         WRITE (tmc_env%w_env%io_unit, *) "worker ", &
565         tmc_env%tmc_comp_set%para_env_sub_group%mepos, "of group ", &
566         tmc_env%tmc_comp_set%group_nr, "stops working!!!!!!!!!!!!!!!!!!"
567
568      IF (PRESENT(ana_list)) THEN
569         DO itmp = 1, tmc_env%params%nr_temp
570            ana_list(itmp)%temp%atoms => NULL()
571            ana_list(itmp)%temp%cell => NULL()
572         END DO
573      END IF
574      IF (ASSOCIATED(conf)) &
575         CALL deallocate_sub_tree_node(tree_elem=conf)
576      IF (ASSOCIATED(ana_restart_conf)) DEALLOCATE (ana_restart_conf)
577
578      ! end the timing
579      CALL timestop(handle)
580   END SUBROUTINE do_tmc_worker
581
582! **************************************************************************************************
583!> \brief Nested Monte Carlo (NMC), do several Markov Chain Monte Carlo steps
584!>        usually using the approximate potential, could be also Hybrid MC.
585!>        The amount of steps are predefined by the user, but should be huge
586!>        enough to reach the equilibrium state for this potential
587!> \param conf ...
588!> \param env_id ...
589!> \param tmc_env ...
590!> \param calc_status ...
591!> \param
592!> \author Mandes 11.2012
593! **************************************************************************************************
594   SUBROUTINE nested_markov_chain_MC(conf, env_id, tmc_env, calc_status)
595      TYPE(tree_type), POINTER                           :: conf
596      INTEGER, INTENT(IN)                                :: env_id
597      TYPE(tmc_env_type), POINTER                        :: tmc_env
598      INTEGER, INTENT(OUT)                               :: calc_status
599
600      CHARACTER(LEN=*), PARAMETER :: routineN = 'nested_markov_chain_MC', &
601         routineP = moduleN//':'//routineN
602
603      INTEGER                                            :: comm_dest, handle, substeps
604      LOGICAL                                            :: accept, change_rejected, flag
605      REAL(KIND=dp)                                      :: rnd_nr
606      TYPE(tree_type), POINTER                           :: last_acc_conf
607
608      NULLIFY (last_acc_conf)
609
610      CPASSERT(ASSOCIATED(tmc_env))
611      CPASSERT(ASSOCIATED(tmc_env%params))
612      CPASSERT(ASSOCIATED(tmc_env%tmc_comp_set))
613      CPASSERT(ALLOCATED(tmc_env%rng_stream))
614      CPASSERT(ASSOCIATED(conf))
615      CPASSERT(conf%temp_created .GT. 0)
616      CPASSERT(conf%temp_created .LE. tmc_env%params%nr_temp)
617      CPASSERT(env_id .GT. 0)
618
619      ! start the timing
620      CALL timeset(routineN, handle)
621
622      CALL allocate_new_sub_tree_node(tmc_params=tmc_env%params, &
623                                      next_el=last_acc_conf, nr_dim=SIZE(conf%pos))
624
625      last_acc_conf%pos = conf%pos
626      last_acc_conf%box_scale = conf%box_scale
627
628      ! energy of the last accepted configuration
629      CALL calc_potential_energy(conf=last_acc_conf, &
630                                 env_id=tmc_env%w_env%env_id_approx, exact_approx_pot=.FALSE., &
631                                 tmc_env=tmc_env)
632
633      NMC_steps: DO substeps = 1, INT(tmc_env%params%move_types%mv_size(mv_type_NMC_moves, 1))
634         ! check for canceling message
635         IF (ASSOCIATED(tmc_env%tmc_comp_set%para_env_m_w)) THEN
636            flag = .FALSE.
637            comm_dest = MASTER_COMM_ID
638            ! check for new canceling message
639            CALL tmc_message(msg_type=calc_status, send_recv=recv_msg, &
640                             dest=comm_dest, &
641                             para_env=tmc_env%tmc_comp_set%para_env_m_w, &
642                             tmc_params=tmc_env%params, success=flag)
643         END IF
644         comm_dest = bcast_group
645         CALL tmc_message(msg_type=calc_status, send_recv=send_msg, &
646                          dest=comm_dest, &
647                          para_env=tmc_env%tmc_comp_set%para_env_sub_group, &
648                          tmc_params=tmc_env%params)
649         SELECT CASE (calc_status)
650         CASE (TMC_STATUS_CALCULATING)
651            ! keep on working
652         CASE (TMC_CANCELING_MESSAGE)
653            ! nothing to do, because calculation CANCELING, exit with cancel status
654            EXIT NMC_steps
655         CASE DEFAULT
656            CALL cp_abort(__LOCATION__, &
657                          "unknown status "//cp_to_string(calc_status)// &
658                          "in the NMC routine, expect only caneling status. ")
659         END SELECT
660
661         ! set move type
662         CALL tmc_env%rng_stream%set( &
663            bg=conf%rng_seed(:, :, 1), cg=conf%rng_seed(:, :, 2), &
664            ig=conf%rng_seed(:, :, 3))
665         conf%move_type = select_random_move_type( &
666                          move_types=tmc_env%params%nmc_move_types, &
667                          rnd=tmc_env%rng_stream%next())
668         CALL tmc_env%rng_stream%get( &
669            bg=conf%rng_seed(:, :, 1), cg=conf%rng_seed(:, :, 2), &
670            ig=conf%rng_seed(:, :, 3))
671
672         ! do move
673         CALL change_pos(tmc_params=tmc_env%params, &
674                         move_types=tmc_env%params%nmc_move_types, &
675                         rng_stream=tmc_env%rng_stream, &
676                         elem=conf, mv_conf=1, new_subbox=.FALSE., &
677                         move_rejected=change_rejected)
678         ! for Hybrid MC the change_pos is only velocity change,
679         !   the actual MD step hast to be done in this module for communication reason
680         IF (conf%move_type .EQ. mv_type_MD) THEN
681            !TODO implement the MD part
682            !CALL calc_MD_step(...)
683            !CALL calc_calc_e_kin(...)
684            CALL cp_abort(__LOCATION__, &
685                          "Hybrid MC is not implemented yet, "// &
686                          "(no MD section in TMC yet). ")
687         END IF
688
689         ! update the subbox acceptance probabilities
690         CALL prob_update(move_types=tmc_env%params%nmc_move_types, elem=conf, &
691                          acc=.NOT. change_rejected, subbox=.TRUE., &
692                          prob_opt=tmc_env%params%esimate_acc_prob)
693
694         ! calculate potential energy if necessary
695         IF (.NOT. change_rejected) THEN
696            CALL calc_potential_energy(conf=conf, &
697                                       env_id=tmc_env%w_env%env_id_approx, exact_approx_pot=.FALSE., &
698                                       tmc_env=tmc_env)
699         ELSE
700            conf%e_pot_approx = HUGE(conf%e_pot_approx)
701         END IF
702
703         !check NMC step
704         CALL tmc_env%rng_stream%set( &
705            bg=conf%rng_seed(:, :, 1), cg=conf%rng_seed(:, :, 2), &
706            ig=conf%rng_seed(:, :, 3))
707         rnd_nr = tmc_env%rng_stream%next()
708         CALL tmc_env%rng_stream%get( &
709            bg=conf%rng_seed(:, :, 1), cg=conf%rng_seed(:, :, 2), &
710            ig=conf%rng_seed(:, :, 3))
711
712         IF (.NOT. change_rejected) THEN
713            CALL acceptance_check(tree_element=conf, parent_element=last_acc_conf, &
714                                  tmc_params=tmc_env%params, &
715                                  temperature=tmc_env%params%Temp(conf%temp_created), &
716                                  diff_pot_check=.FALSE., &
717                                  accept=accept, approx_ener=.TRUE., rnd_nr=rnd_nr)
718         ELSE
719            accept = .FALSE.
720         END IF
721         ! update the NMC accpetance per move
722         CALL prob_update(move_types=tmc_env%params%nmc_move_types, elem=conf, &
723                          acc=accept, prob_opt=tmc_env%params%esimate_acc_prob)
724
725         ! update last accepted configuration or actual configuration
726         IF (accept .AND. (.NOT. change_rejected)) THEN
727            last_acc_conf%pos = conf%pos
728            last_acc_conf%vel = conf%vel
729            last_acc_conf%e_pot_approx = conf%e_pot_approx
730            last_acc_conf%ekin = conf%ekin
731            last_acc_conf%ekin_before_md = conf%ekin_before_md
732            last_acc_conf%box_scale = conf%box_scale
733         ELSE
734            conf%pos = last_acc_conf%pos
735            conf%vel = last_acc_conf%vel
736            conf%box_scale = last_acc_conf%box_scale
737         END IF
738      END DO NMC_steps
739
740      ! result values of Nested Monte Carlo (NMC) steps
741      !   regard that the calculated potential energy is the one of the approximated potential
742      conf%pos = last_acc_conf%pos
743      conf%vel = last_acc_conf%vel
744      conf%e_pot_approx = last_acc_conf%e_pot_approx
745      conf%potential = 0.0_dp
746      conf%ekin = last_acc_conf%ekin
747      conf%ekin_before_md = last_acc_conf%ekin_before_md
748
749      CALL deallocate_sub_tree_node(tree_elem=last_acc_conf)
750
751      ! end the timing
752      CALL timestop(handle)
753   END SUBROUTINE nested_markov_chain_MC
754
755! **************************************************************************************************
756!> \brief get the initial confuguration (pos,...)
757!> \param tmc_params ...
758!> \param init_conf the structure the data should be stored
759!> force_env
760!> \param env_id ...
761!> \author Mandes 11.2012
762! **************************************************************************************************
763   SUBROUTINE get_initial_conf(tmc_params, init_conf, env_id)
764      TYPE(tmc_param_type), POINTER                      :: tmc_params
765      TYPE(tree_type), POINTER                           :: init_conf
766      INTEGER                                            :: env_id
767
768      CHARACTER(LEN=*), PARAMETER :: routineN = 'get_initial_conf', &
769         routineP = moduleN//':'//routineN
770
771      INTEGER                                            :: handle, ierr, mol, ndim, nr_atoms
772      TYPE(cp_subsys_type), POINTER                      :: subsys
773      TYPE(f_env_type), POINTER                          :: f_env
774      TYPE(molecule_list_type), POINTER                  :: molecule_new
775
776      CPASSERT(.NOT. ASSOCIATED(init_conf))
777
778      ! start the timing
779      CALL timeset(routineN, handle)
780
781      ! get positions
782      CALL get_natom(env_id=env_id, n_atom=nr_atoms, ierr=ierr)
783      CPASSERT(ierr .EQ. 0)
784      ndim = 3*nr_atoms
785      CALL allocate_new_sub_tree_node(tmc_params=tmc_params, &
786                                      next_el=init_conf, nr_dim=ndim)
787      CALL get_pos(env_id=env_id, pos=init_conf%pos, n_el=SIZE(init_conf%pos), &
788                   ierr=ierr)
789
790      ! get the molecule info
791      CALL f_env_get_from_id(env_id, f_env)
792      CALL force_env_get(f_env%force_env, subsys=subsys)
793
794      CALL cp_subsys_get(subsys=subsys, molecules=molecule_new)
795      loop_mol: DO mol = 1, SIZE(molecule_new%els(:))
796         init_conf%mol(molecule_new%els(mol)%first_atom: &
797                       molecule_new%els(mol)%last_atom) = mol
798      END DO loop_mol
799
800      ! end the timing
801      CALL timestop(handle)
802
803   END SUBROUTINE get_initial_conf
804
805! **************************************************************************************************
806!> \brief get the pointer to the atoms, for easy handling
807!> \param env_id ...
808!> \param atoms pointer to atomic_kind
809!> \param cell ...
810!> \author Mandes 01.2013
811! **************************************************************************************************
812   SUBROUTINE get_atom_kinds_and_cell(env_id, atoms, cell)
813      INTEGER                                            :: env_id
814      TYPE(tmc_atom_type), DIMENSION(:), POINTER         :: atoms
815      TYPE(cell_type), POINTER                           :: cell
816
817      CHARACTER(LEN=*), PARAMETER :: routineN = 'get_atom_kinds_and_cell', &
818         routineP = moduleN//':'//routineN
819
820      INTEGER                                            :: iparticle, nr_atoms, nunits_tot
821      TYPE(cell_type), POINTER                           :: cell_tmp
822      TYPE(cp_subsys_type), POINTER                      :: subsys
823      TYPE(f_env_type), POINTER                          :: f_env
824      TYPE(particle_list_type), POINTER                  :: particles
825
826      NULLIFY (f_env, subsys, particles)
827      nr_atoms = 0
828
829      CPASSERT(env_id .GT. 0)
830      CPASSERT(.NOT. ASSOCIATED(atoms))
831      CPASSERT(.NOT. ASSOCIATED(cell))
832
833      CALL f_env_get_from_id(env_id, f_env)
834      nr_atoms = force_env_get_natom(f_env%force_env)
835      CALL force_env_get(f_env%force_env, subsys=subsys, cell=cell_tmp)
836      ALLOCATE (cell)
837      CALL cell_copy(cell_in=cell_tmp, cell_out=cell)
838
839      !get atom kinds
840      CALL allocate_tmc_atom_type(atoms, nr_atoms)
841      CALL cp_subsys_get(subsys, particles=particles)
842      nunits_tot = SIZE(particles%els(:))
843      IF (nunits_tot .GT. 0) THEN
844         DO iparticle = 1, nunits_tot
845            atoms(iparticle)%name = particles%els(iparticle)%atomic_kind%name
846            atoms(iparticle)%mass = particles%els(iparticle)%atomic_kind%mass
847         END DO
848         CPASSERT(iparticle - 1 .EQ. nr_atoms)
849      ENDIF
850   END SUBROUTINE get_atom_kinds_and_cell
851
852! **************************************************************************************************
853!> \brief set the communicator in the SCF environment
854!>        to receive the intermediate energies on the (global) master side
855!> \param comm the master-worker communicator
856!> \param env_id the ID of the related force environment
857!> \author Mandes 10.2013
858! **************************************************************************************************
859   SUBROUTINE set_intermediate_info_comm(comm, env_id)
860      INTEGER, INTENT(IN)                                :: comm
861      INTEGER                                            :: env_id
862
863      CHARACTER(LEN=*), PARAMETER :: routineN = 'set_intermediate_info_comm', &
864         routineP = moduleN//':'//routineN
865
866      CHARACTER(LEN=default_string_length)               :: description
867      REAL(KIND=dp), DIMENSION(3)                        :: values
868      TYPE(cp_result_type), POINTER                      :: results
869      TYPE(cp_subsys_type), POINTER                      :: subsys
870      TYPE(f_env_type), POINTER                          :: f_env
871
872      NULLIFY (results, subsys)
873      CPASSERT(env_id .GT. 0)
874
875      CALL f_env_get_from_id(env_id, f_env)
876
877      CPASSERT(ASSOCIATED(f_env))
878      CPASSERT(ASSOCIATED(f_env%force_env))
879      IF (.NOT. ASSOCIATED(f_env%force_env%qs_env)) &
880         CALL cp_abort(__LOCATION__, &
881                       "the intermediate SCF energy request can not be set "// &
882                       "employing this force environment! ")
883
884      ! set the information
885      values(1) = REAL(comm, KIND=dp)
886      values(2) = REAL(MASTER_COMM_ID, KIND=dp)
887      values(3) = REAL(TMC_STAT_SCF_STEP_ENER_RECEIVE, KIND=dp)
888      description = "[EXT_SCF_ENER_COMM]"
889
890      ! set the communicator information in the qs_env result container
891      CALL force_env_get(f_env%force_env, subsys=subsys)
892      CALL cp_subsys_get(subsys, results=results)
893      CALL put_results(results, description=description, values=values)
894   END SUBROUTINE set_intermediate_info_comm
895
896! **************************************************************************************************
897!> \brief set the communicator in the SCF environment
898!>        to receive the intermediate energies on the (global) master side
899!> \param env_id the ID of the related force environment
900!> \author Mandes 10.2013
901! **************************************************************************************************
902   SUBROUTINE remove_intermediate_info_comm(env_id)
903      INTEGER                                            :: env_id
904
905      CHARACTER(LEN=*), PARAMETER :: routineN = 'remove_intermediate_info_comm', &
906         routineP = moduleN//':'//routineN
907
908      CHARACTER(LEN=default_string_length)               :: description
909      TYPE(cp_result_type), POINTER                      :: results
910      TYPE(cp_subsys_type), POINTER                      :: subsys
911      TYPE(f_env_type), POINTER                          :: f_env
912
913      NULLIFY (subsys, results)
914      CPASSERT(env_id .GT. 0)
915
916      CALL f_env_get_from_id(env_id, f_env)
917
918      CPASSERT(ASSOCIATED(f_env))
919      CPASSERT(ASSOCIATED(f_env%force_env))
920      IF (.NOT. ASSOCIATED(f_env%force_env%qs_env)) &
921         CALL cp_abort(__LOCATION__, &
922                       "the SCF intermediate energy communicator can not be "// &
923                       "removed! ")
924
925      description = "[EXT_SCF_ENER_COMM]"
926
927      ! set the communicator information in the qs_env result container
928      CALL force_env_get(f_env%force_env, subsys=subsys)
929      CALL cp_subsys_get(subsys, results=results)
930      CALL cp_results_erase(results, description=description)
931   END SUBROUTINE remove_intermediate_info_comm
932
933END MODULE tmc_worker
934