(git:0de0cc2)
metadynamics_utils.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Performs the metadynamics calculation
10 !> \par History
11 !> 01.2005 created [fawzi and ale]
12 !> 11.2007 Teodoro Laino [tlaino] - University of Zurich
13 ! **************************************************************************************************
15  USE cp_files, ONLY: close_file,&
16  open_file
18  cp_logger_type,&
19  cp_to_string
22  USE cp_subsys_types, ONLY: cp_subsys_type
23  USE force_env_types, ONLY: force_env_get,&
24  force_env_type
25  USE input_constants, ONLY: do_fe_meta,&
27  do_wall_m,&
28  do_wall_none,&
29  do_wall_p,&
35  enumeration_type
36  USE input_keyword_types, ONLY: keyword_get,&
37  keyword_type
41  section_type,&
44  section_vals_type,&
46  USE kinds, ONLY: default_path_length,&
47  dp
48  USE machine, ONLY: m_mov
49  USE message_passing, ONLY: mp_para_env_type
50  USE metadynamics_types, ONLY: hills_env_type,&
51  meta_env_type,&
53  metavar_type,&
54  multiple_walkers_type
55  USE physcon, ONLY: kelvin
56 #include "./base/base_uses.f90"
57 
58  IMPLICIT NONE
59  PRIVATE
60 
61  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'metadynamics_utils'
62 
63  PUBLIC :: metadyn_read, &
66  restart_hills, &
69 
70 CONTAINS
71 
72 ! **************************************************************************************************
73 !> \brief reads metadynamics section
74 !> \param meta_env ...
75 !> \param force_env ...
76 !> \param root_section ...
77 !> \param para_env ...
78 !> \param fe_section ...
79 !> \par History
80 !> 04.2004 created
81 !> \author Teodoro Laino [tlaino] - University of Zurich. 11.2007
82 ! **************************************************************************************************
83  SUBROUTINE metadyn_read(meta_env, force_env, root_section, para_env, fe_section)
84  TYPE(meta_env_type), POINTER :: meta_env
85  TYPE(force_env_type), POINTER :: force_env
86  TYPE(section_vals_type), POINTER :: root_section
87  TYPE(mp_para_env_type), POINTER :: para_env
88  TYPE(section_vals_type), OPTIONAL, POINTER :: fe_section
89 
90  CHARACTER(len=*), PARAMETER :: routinen = 'metadyn_read'
91 
92  CHARACTER(LEN=default_path_length) :: walkers_file_name
93  INTEGER :: handle, i, id_method, n_colvar, n_rep, &
94  number_allocated_colvars
95  INTEGER, DIMENSION(:), POINTER :: walkers_status
96  LOGICAL :: check, explicit
97  REAL(kind=dp) :: dt
98  TYPE(cp_subsys_type), POINTER :: subsys
99  TYPE(section_vals_type), POINTER :: md_section, metadyn_section, &
100  metavar_section, walkers_section
101 
102  NULLIFY (subsys)
103  CALL timeset(routinen, handle)
104 
105  CALL section_vals_get(fe_section, explicit=explicit)
106  IF (explicit) THEN
107  number_allocated_colvars = 0
108  CALL force_env_get(force_env, subsys=subsys)
109  IF (ASSOCIATED(subsys%colvar_p)) THEN
110  number_allocated_colvars = SIZE(subsys%colvar_p)
111  END IF
112  CALL section_vals_val_get(fe_section, "METHOD", i_val=id_method)
113  IF (id_method /= do_fe_meta) THEN
114  CALL timestop(handle)
115  RETURN
116  END IF
117  metadyn_section => section_vals_get_subs_vals(fe_section, "METADYN")
118  cpassert(.NOT. ASSOCIATED(meta_env))
119 
120  md_section => section_vals_get_subs_vals(root_section, "MOTION%MD")
121  CALL section_vals_val_get(md_section, "TIMESTEP", r_val=dt)
122 
123  metavar_section => section_vals_get_subs_vals(metadyn_section, "METAVAR")
124  CALL section_vals_get(metavar_section, n_repetition=n_colvar)
125  ALLOCATE (meta_env)
126  CALL metadyn_create(meta_env, n_colvar=n_colvar, &
127  dt=dt, para_env=para_env, metadyn_section=metadyn_section)
128 
129  !Check if using plumed. If so, only get the file name and read nothing else
130  CALL section_vals_val_get(metadyn_section, "USE_PLUMED", l_val=meta_env%use_plumed)
131  IF (meta_env%use_plumed .EQV. .true.) THEN
132  CALL section_vals_val_get(metadyn_section, "PLUMED_INPUT_FILE", c_val=meta_env%plumed_input_file)
133  meta_env%plumed_input_file = trim(meta_env%plumed_input_file)//char(0)
134  meta_env%langevin = .false.
135  CALL timestop(handle)
136  RETURN
137  END IF
138 
139  CALL section_vals_val_get(metadyn_section, "DO_HILLS", l_val=meta_env%do_hills)
140  CALL section_vals_val_get(metadyn_section, "LAGRANGE", l_val=meta_env%extended_lagrange)
141  CALL section_vals_val_get(metadyn_section, "TAMCSteps", i_val=meta_env%TAMCSteps)
142  IF (meta_env%TAMCSteps < 0) THEN
143  cpabort("TAMCSteps must be positive!")
144  END IF
145  CALL section_vals_val_get(metadyn_section, "Timestep", r_val=meta_env%zdt)
146  IF (meta_env%zdt <= 0.0_dp) THEN
147  cpabort("Timestep must be positive!")
148  END IF
149  CALL section_vals_val_get(metadyn_section, "WW", r_val=meta_env%hills_env%ww)
150  CALL section_vals_val_get(metadyn_section, "NT_HILLS", i_val=meta_env%hills_env%nt_hills)
151  CALL section_vals_val_get(metadyn_section, "MIN_NT_HILLS", i_val=meta_env%hills_env%min_nt_hills)
152  IF (meta_env%hills_env%nt_hills <= 0) THEN
153  meta_env%hills_env%min_nt_hills = meta_env%hills_env%nt_hills
154  CALL cp_warn(__location__, &
155  "NT_HILLS has a value <= 0; "// &
156  "Setting MIN_NT_HILLS to the same value! "// &
157  "Overriding input specification!")
158  END IF
159  check = meta_env%hills_env%nt_hills >= meta_env%hills_env%min_nt_hills
160  IF (.NOT. check) &
161  CALL cp_abort(__location__, "MIN_NT_HILLS must have a value smaller or equal to NT_HILLS! "// &
162  "Cross check with the input reference!")
163  !RG Adaptive hills
164  CALL section_vals_val_get(metadyn_section, "MIN_DISP", r_val=meta_env%hills_env%min_disp)
165  CALL section_vals_val_get(metadyn_section, "OLD_HILL_NUMBER", i_val=meta_env%hills_env%old_hill_number)
166  CALL section_vals_val_get(metadyn_section, "OLD_HILL_STEP", i_val=meta_env%hills_env%old_hill_step)
167 
168  !Hills tail damping
169  CALL section_vals_val_get(metadyn_section, "HILL_TAIL_CUTOFF", r_val=meta_env%hills_env%tail_cutoff)
170  CALL section_vals_val_get(metadyn_section, "P_EXPONENT", i_val=meta_env%hills_env%p_exp)
171  CALL section_vals_val_get(metadyn_section, "Q_EXPONENT", i_val=meta_env%hills_env%q_exp)
172 
173  CALL section_vals_val_get(metadyn_section, "SLOW_GROWTH", l_val=meta_env%hills_env%slow_growth)
174 
175  !RG Adaptive hills
176  CALL section_vals_val_get(metadyn_section, "STEP_START_VAL", i_val=meta_env%n_steps)
177  cpassert(meta_env%n_steps >= 0)
178  CALL section_vals_val_get(metadyn_section, "NHILLS_START_VAL", &
179  i_val=meta_env%hills_env%n_hills)
180  CALL section_vals_val_get(metadyn_section, "TEMPERATURE", r_val=meta_env%temp_wanted)
181  CALL section_vals_val_get(metadyn_section, "LANGEVIN", l_val=meta_env%langevin)
182  CALL section_vals_val_get(metadyn_section, "TEMP_TOL", explicit=meta_env%tempcontrol, &
183  r_val=meta_env%toll_temp)
184  CALL section_vals_val_get(metadyn_section, "WELL_TEMPERED", l_val=meta_env%well_tempered)
185  CALL section_vals_val_get(metadyn_section, "DELTA_T", explicit=meta_env%hills_env%wtcontrol, &
186  r_val=meta_env%delta_t)
187  CALL section_vals_val_get(metadyn_section, "WTGAMMA", explicit=check, &
188  r_val=meta_env%wtgamma)
189  IF (meta_env%well_tempered) THEN
190  meta_env%hills_env%wtcontrol = meta_env%hills_env%wtcontrol .OR. check
191  check = meta_env%hills_env%wtcontrol
192  IF (.NOT. check) &
193  CALL cp_abort(__location__, "When using Well-Tempered metadynamics, "// &
194  "DELTA_T (or WTGAMMA) should be explicitly specified.")
195  IF (meta_env%extended_lagrange) &
196  CALL cp_abort(__location__, &
197  "Well-Tempered metadynamics not possible with extended-lagrangian formulation.")
198  IF (meta_env%hills_env%min_disp > 0.0_dp) &
199  CALL cp_abort(__location__, &
200  "Well-Tempered metadynamics not possible with Adaptive hills.")
201  END IF
202 
203  CALL section_vals_val_get(metadyn_section, "COLVAR_AVG_TEMPERATURE_RESTART", &
204  r_val=meta_env%avg_temp)
205  ! Parsing Metavar Section
206  DO i = 1, n_colvar
207  CALL metavar_read(meta_env%metavar(i), meta_env%extended_lagrange, &
208  meta_env%langevin, i, metavar_section)
209  check = (meta_env%metavar(i)%icolvar <= number_allocated_colvars)
210  IF (.NOT. check) &
211  CALL cp_abort(__location__, &
212  "An error occurred in the specification of COLVAR for METAVAR. "// &
213  "Specified COLVAR #("//trim(adjustl(cp_to_string(meta_env%metavar(i)%icolvar)))//") "// &
214  "is larger than the maximum number of COLVARS defined in the SUBSYS ("// &
215  trim(adjustl(cp_to_string(number_allocated_colvars)))//") !")
216  END DO
217 
218  ! Parsing the Multiple Walkers Info
219  IF (meta_env%do_multiple_walkers) THEN
220  NULLIFY (walkers_status)
221  walkers_section => section_vals_get_subs_vals(metadyn_section, "MULTIPLE_WALKERS")
222 
223  ! General setup for walkers
224  CALL section_vals_val_get(walkers_section, "WALKER_ID", &
225  i_val=meta_env%multiple_walkers%walker_id)
226  CALL section_vals_val_get(walkers_section, "NUMBER_OF_WALKERS", &
227  i_val=meta_env%multiple_walkers%walkers_tot_nr)
228  CALL section_vals_val_get(walkers_section, "WALKER_COMM_FREQUENCY", &
229  i_val=meta_env%multiple_walkers%walkers_freq_comm)
230 
231  ! Handle status and file names
232  ALLOCATE (meta_env%multiple_walkers%walkers_status(meta_env%multiple_walkers%walkers_tot_nr))
233  ALLOCATE (meta_env%multiple_walkers%walkers_file_name(meta_env%multiple_walkers%walkers_tot_nr))
234  CALL section_vals_val_get(walkers_section, "WALKERS_STATUS", explicit=explicit)
235  IF (explicit) THEN
236  CALL section_vals_val_get(walkers_section, "WALKERS_STATUS", i_vals=walkers_status)
237  check = (SIZE(walkers_status) == meta_env%multiple_walkers%walkers_tot_nr)
238  IF (.NOT. check) &
239  CALL cp_abort(__location__, &
240  "Number of Walkers specified in the input does not match with the "// &
241  "size of the WALKERS_STATUS. Please check your input and in case "// &
242  "this is a restart run consider the possibility to switch off the "// &
243  "RESTART_WALKERS in the EXT_RESTART section! ")
244  meta_env%multiple_walkers%walkers_status = walkers_status
245  ELSE
246  meta_env%multiple_walkers%walkers_status = 0
247  END IF
248  meta_env%multiple_walkers%n_hills_local = &
249  meta_env%multiple_walkers%walkers_status(meta_env%multiple_walkers%walker_id)
250 
251  CALL section_vals_val_get(walkers_section, "WALKERS_FILE_NAME%_DEFAULT_KEYWORD_", &
252  n_rep_val=n_rep)
253  check = (n_rep == meta_env%multiple_walkers%walkers_tot_nr)
254  IF (.NOT. check) &
255  CALL cp_abort(__location__, &
256  "Number of Walkers specified in the input does not match with the "// &
257  "number of Walkers File names provided. Please check your input and in case "// &
258  "this is a restart run consider the possibility to switch off the "// &
259  "RESTART_WALKERS in the EXT_RESTART section! ")
260  DO i = 1, n_rep
261  CALL section_vals_val_get(walkers_section, "WALKERS_FILE_NAME%_DEFAULT_KEYWORD_", &
262  i_rep_val=i, c_val=walkers_file_name)
263  meta_env%multiple_walkers%walkers_file_name(i) = walkers_file_name
264  END DO
265  END IF
266 
267  ! Print Metadynamics Info
268  CALL print_metadyn_info(meta_env, n_colvar, metadyn_section)
269  END IF
270 
271  CALL timestop(handle)
272 
273  END SUBROUTINE metadyn_read
274 
275 ! **************************************************************************************************
276 !> \brief prints information on the metadynamics run
277 !> \param meta_env ...
278 !> \param n_colvar ...
279 !> \param metadyn_section ...
280 !> \author Teodoro Laino [tlaino] - University of Zurich. 10.2008
281 ! **************************************************************************************************
282  SUBROUTINE print_metadyn_info(meta_env, n_colvar, metadyn_section)
283  TYPE(meta_env_type), POINTER :: meta_env
284  INTEGER, INTENT(IN) :: n_colvar
285  TYPE(section_vals_type), POINTER :: metadyn_section
286 
287  CHARACTER(len=*), PARAMETER :: routinen = 'print_metadyn_info'
288 
289  CHARACTER(LEN=10) :: my_id, my_tag
290  INTEGER :: handle, i, iw, j
291  TYPE(cp_logger_type), POINTER :: logger
292  TYPE(enumeration_type), POINTER :: enum
293  TYPE(keyword_type), POINTER :: keyword
294  TYPE(section_type), POINTER :: section, wall_section, work_section
295 
296  CALL timeset(routinen, handle)
297 
298  logger => cp_get_default_logger()
299  iw = cp_print_key_unit_nr(logger, metadyn_section, &
300  "PRINT%PROGRAM_RUN_INFO", extension=".metadynLog")
301  NULLIFY (section, enum, keyword)
302  CALL create_metavar_section(section)
303  wall_section => section_get_subsection(section, "WALL")
304  IF (iw > 0) THEN
305  WRITE (iw, '( /A )') ' METADYN| Meta Dynamics Protocol '
306  WRITE (iw, '( A,T71,I10)') ' METADYN| Number of interval time steps to spawn hills', &
307  meta_env%hills_env%nt_hills
308  WRITE (iw, '( A,T71,I10)') ' METADYN| Number of previously spawned hills', &
309  meta_env%hills_env%n_hills
310  IF (meta_env%extended_lagrange) THEN
311  WRITE (iw, '( A )') ' METADYN| Extended Lagrangian Scheme '
312  IF (meta_env%tempcontrol) WRITE (iw, '( A,T71,F10.2)') &
313  ' METADYN| Collective Variables Temperature control', meta_env%toll_temp
314  IF (meta_env%langevin) THEN
315  WRITE (iw, '(A,T71)') ' METADYN| Langevin Thermostat in use for COLVAR '
316  WRITE (iw, '(A,T71,F10.4)') ' METADYN| Langevin Thermostat. Target Temperature = ', &
317  meta_env%temp_wanted*kelvin
318  END IF
319  WRITE (iw, '(A,T71,F10.4)') ' METADYN| COLVARS restarted average temperature ', &
320  meta_env%avg_temp
321  END IF
322  IF (meta_env%do_hills) THEN
323  WRITE (iw, '( A )') ' METADYN| Spawning the Hills '
324  WRITE (iw, '( A,T71,F10.3)') ' METADYN| Height of the Spawned Gaussian', meta_env%hills_env%ww
325  !RG Adaptive hills
326  IF (meta_env%hills_env%min_disp .GT. 0.0_dp) THEN
327  WRITE (iw, '(A)') ' METADYN| Adapative meta time step is activated'
328  WRITE (iw, '(A,T71,F10.4)') ' METADYN| Minimum displacement for next hill', &
329  meta_env%hills_env%min_disp
330  END IF
331  !RG Adaptive hills
332  END IF
333 
334  IF (meta_env%well_tempered) THEN
335  WRITE (iw, '( A )') ' METADYN| Well-Tempered metadynamics '
336  IF (meta_env%delta_t > epsilon(1._dp)) THEN
337  WRITE (iw, '( A,T71,F10.3)') ' METADYN| Temperature parameter (Delta T) [K]', meta_env%delta_t*kelvin
338  ELSE
339  WRITE (iw, '( A,T71,F10.3)') ' METADYN| Temperature parameter (gamma)', meta_env%wtgamma
340  END IF
341  END IF
342 
343  IF (meta_env%do_multiple_walkers) THEN
344  WRITE (iw, '( A,T71,A10)') ' METADYN| Multiple Walkers', ' ENABLED'
345  WRITE (iw, '( A,T71,I10)') ' METADYN| Number of Multiple Walkers', &
346  meta_env%multiple_walkers%walkers_tot_nr
347  WRITE (iw, '( A,T71,I10)') ' METADYN| Local Walker ID', &
348  meta_env%multiple_walkers%walker_id
349  WRITE (iw, '( A,T71,I10)') ' METADYN| Walker Communication Frequency', &
350  meta_env%multiple_walkers%walkers_freq_comm
351  DO i = 1, meta_env%multiple_walkers%walkers_tot_nr
352  my_tag = ""
353  IF (i == meta_env%multiple_walkers%walker_id) my_tag = " ( Local )"
354  my_id = '( '//trim(adjustl(cp_to_string(i)))//' )'
355  WRITE (iw, '(/,A,T71,A10)') ' WALKERS| Walker ID'//trim(my_tag), adjustr(my_id)
356  WRITE (iw, '( A,T71,I10)') ' WALKERS| Number of Hills communicated', &
357  meta_env%multiple_walkers%walkers_status(i)
358  WRITE (iw, '( A,T24,A57)') ' WALKERS| Base Filename', &
359  adjustr(meta_env%multiple_walkers%walkers_file_name(i) (1:57))
360  END DO
361  WRITE (iw, '(/)')
362  END IF
363 
364  WRITE (iw, '( A,T71,I10)') ' METADYN| Number of collective variables', meta_env%n_colvar
365  DO i = 1, n_colvar
366  WRITE (iw, '( A )') ' '//'----------------------------------------------------------------------'
367  WRITE (iw, '( A,T71,I10)') ' METAVARS| Collective Variable Number', meta_env%metavar(i)%icolvar
368  IF (meta_env%extended_lagrange) THEN
369  WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Lambda Parameter', meta_env%metavar(i)%lambda
370  WRITE (iw, '( A,T66,F15.6)') ' METAVARS| Collective Variable Mass', meta_env%metavar(i)%mass
371  END IF
372  WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Scaling factor', meta_env%metavar(i)%delta_s
373  IF (meta_env%langevin) WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Friction for Langevin Thermostat', &
374  meta_env%metavar(i)%gamma
375  IF (meta_env%metavar(i)%do_wall) THEN
376  WRITE (iw, '( A,T71,I10)') ' METAVARS| Number of Walls present', SIZE(meta_env%metavar(i)%walls)
377  DO j = 1, SIZE(meta_env%metavar(i)%walls)
378  keyword => section_get_keyword(wall_section, "TYPE")
379  CALL keyword_get(keyword, enum=enum)
380  WRITE (iw, '(/,A,5X,I10,T50,A,T70,A11)') ' METAVARS| Wall Number:', j, 'Type of Wall:', &
381  adjustr(trim(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_type)))
382  ! Type of wall IO
383  SELECT CASE (meta_env%metavar(i)%walls(j)%id_type)
384  CASE (do_wall_none)
385  ! Do Nothing
386  cycle
387  CASE (do_wall_reflective)
388  work_section => section_get_subsection(wall_section, "REFLECTIVE")
389  keyword => section_get_keyword(work_section, "DIRECTION")
390  CALL keyword_get(keyword, enum=enum)
391  WRITE (iw, '(A,T70,A11)') ' METAVARS| Wall direction', &
392  adjustr(trim(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_direction)))
393  CASE (do_wall_quadratic)
394  work_section => section_get_subsection(wall_section, "QUADRATIC")
395  keyword => section_get_keyword(work_section, "DIRECTION")
396  CALL keyword_get(keyword, enum=enum)
397  WRITE (iw, '(A,T70,A11)') ' METAVARS| Wall direction', &
398  adjustr(trim(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_direction)))
399  WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Constant K of the quadratic potential', &
400  meta_env%metavar(i)%walls(j)%k_quadratic
401  CASE (do_wall_gaussian)
402  WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Height of the Wall Gaussian', &
403  meta_env%metavar(i)%walls(j)%ww_gauss
404  WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Scale of the Wall Gaussian', &
405  meta_env%metavar(i)%walls(j)%sigma_gauss
406  END SELECT
407  WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Wall location', &
408  meta_env%metavar(i)%walls(j)%pos
409  END DO
410  END IF
411  WRITE (iw, '( A )') ' '//'----------------------------------------------------------------------'
412  END DO
413  END IF
414  CALL section_release(section)
415  CALL cp_print_key_finished_output(iw, logger, metadyn_section, "PRINT%PROGRAM_RUN_INFO")
416 
417  CALL timestop(handle)
418 
419  END SUBROUTINE print_metadyn_info
420 
421 ! **************************************************************************************************
422 !> \brief reads metavar section
423 !> \param metavar ...
424 !> \param extended_lagrange ...
425 !> \param langevin ...
426 !> \param icol ...
427 !> \param metavar_section ...
428 !> \par History
429 !> 04.2004 created
430 !> \author alessandro laio and fawzi mohamed
431 !> Teodoro Laino [tlaino] - University of Zurich. 11.2007
432 ! **************************************************************************************************
433  SUBROUTINE metavar_read(metavar, extended_lagrange, langevin, icol, metavar_section)
434  TYPE(metavar_type), INTENT(INOUT) :: metavar
435  LOGICAL, INTENT(IN) :: extended_lagrange, langevin
436  INTEGER, INTENT(IN) :: icol
437  TYPE(section_vals_type), OPTIONAL, POINTER :: metavar_section
438 
439  CHARACTER(len=*), PARAMETER :: routinen = 'metavar_read'
440 
441  INTEGER :: handle, i, n_walls
442  TYPE(section_vals_type), POINTER :: wall_section, work_section
443 
444  CALL timeset(routinen, handle)
445 
446  CALL section_vals_val_get(metavar_section, "COLVAR", i_rep_section=icol, i_val=metavar%icolvar)
447  CALL section_vals_val_get(metavar_section, "SCALE", i_rep_section=icol, r_val=metavar%delta_s)
448  ! Walls
449  wall_section => section_vals_get_subs_vals(metavar_section, "WALL", i_rep_section=icol)
450  CALL section_vals_get(wall_section, n_repetition=n_walls)
451  IF (n_walls /= 0) THEN
452  metavar%do_wall = .true.
453  ALLOCATE (metavar%walls(n_walls))
454  DO i = 1, n_walls
455  CALL section_vals_val_get(wall_section, "TYPE", i_rep_section=i, i_val=metavar%walls(i)%id_type)
456  CALL section_vals_val_get(wall_section, "POSITION", i_rep_section=i, r_val=metavar%walls(i)%pos)
457  SELECT CASE (metavar%walls(i)%id_type)
458  CASE (do_wall_none)
459  ! Just cycle..
460  cycle
461  CASE (do_wall_reflective)
462  work_section => section_vals_get_subs_vals(wall_section, "REFLECTIVE", i_rep_section=i)
463  CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
464  CASE (do_wall_quadratic)
465  work_section => section_vals_get_subs_vals(wall_section, "QUADRATIC", i_rep_section=i)
466  CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
467  CALL section_vals_val_get(work_section, "K", r_val=metavar%walls(i)%k_quadratic)
468  CASE (do_wall_quartic)
469  work_section => section_vals_get_subs_vals(wall_section, "QUARTIC", i_rep_section=i)
470  CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
471  CALL section_vals_val_get(work_section, "K", r_val=metavar%walls(i)%k_quartic)
472  SELECT CASE (metavar%walls(i)%id_direction)
473  CASE (do_wall_m)
474  metavar%walls(i)%pos0 = metavar%walls(i)%pos + (0.05_dp/metavar%walls(i)%k_quartic**(0.25_dp))
475  CASE (do_wall_p)
476  metavar%walls(i)%pos0 = metavar%walls(i)%pos - (0.05_dp/metavar%walls(i)%k_quartic**(0.25_dp))
477  END SELECT
478  CASE (do_wall_gaussian)
479  work_section => section_vals_get_subs_vals(wall_section, "GAUSSIAN", i_rep_section=i)
480  CALL section_vals_val_get(work_section, "WW", r_val=metavar%walls(i)%ww_gauss)
481  CALL section_vals_val_get(work_section, "SIGMA", r_val=metavar%walls(i)%sigma_gauss)
482  END SELECT
483  END DO
484  END IF
485  ! Setup few more parameters for extended lagrangian
486  IF (extended_lagrange) THEN
487  CALL section_vals_val_get(metavar_section, "MASS", i_rep_section=icol, r_val=metavar%mass)
488  CALL section_vals_val_get(metavar_section, "LAMBDA", i_rep_section=icol, r_val=metavar%lambda)
489  IF (langevin) THEN
490  CALL section_vals_val_get(metavar_section, "GAMMA", i_rep_section=icol, r_val=metavar%gamma)
491  END IF
492  END IF
493 
494  CALL timestop(handle)
495 
496  END SUBROUTINE metavar_read
497 
498 ! **************************************************************************************************
499 !> \brief Synchronize with the rest of the walkers
500 !> \param multiple_walkers ...
501 !> \param hills_env ...
502 !> \param colvars ...
503 !> \param n_colvar ...
504 !> \param metadyn_section ...
505 !> \author Teodoro Laino [tlaino] - University of Zurich - 10.2008
506 ! **************************************************************************************************
507  SUBROUTINE synchronize_multiple_walkers(multiple_walkers, hills_env, colvars, &
508  n_colvar, metadyn_section)
509  TYPE(multiple_walkers_type), POINTER :: multiple_walkers
510  TYPE(hills_env_type), POINTER :: hills_env
511  TYPE(metavar_type), DIMENSION(:), POINTER :: colvars
512  INTEGER, INTENT(IN) :: n_colvar
513  TYPE(section_vals_type), POINTER :: metadyn_section
514 
515  CHARACTER(len=*), PARAMETER :: routinen = 'synchronize_multiple_walkers'
516 
517  CHARACTER(LEN=default_path_length) :: filename, tmpname
518  INTEGER :: delta_hills, handle, i, i_hills, ih, iw, &
519  unit_nr
520  LOGICAL :: exist
521  REAL(kind=dp) :: invdt, ww
522  REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: delta_s_save, ss0_save
523  REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: delta_s_ss0_buf
524  TYPE(cp_logger_type), POINTER :: logger
525  TYPE(mp_para_env_type), POINTER :: para_env
526 
527  CALL timeset(routinen, handle)
528 
529  logger => cp_get_default_logger()
530  para_env => logger%para_env
531 
532  ! Locally dump information on file..
533  IF (para_env%is_source()) THEN
534  ! Generate file name for the specific Hill
535  i = multiple_walkers%walker_id
536  filename = trim(multiple_walkers%walkers_file_name(i))//"_"// &
537  trim(adjustl(cp_to_string(multiple_walkers%n_hills_local)))
538  tmpname = trim(filename)//".tmp"
539  CALL open_file(file_name=tmpname, file_status="UNKNOWN", &
540  file_form="FORMATTED", file_action="WRITE", &
541  file_position="APPEND", unit_number=unit_nr)
542  WRITE (unit_nr, *) hills_env%ww_history(hills_env%n_hills)
543  DO ih = 1, n_colvar
544  WRITE (unit_nr, *) hills_env%ss_history(ih, hills_env%n_hills)
545  WRITE (unit_nr, *) hills_env%delta_s_history(ih, hills_env%n_hills)
546  END DO
547  IF (hills_env%wtcontrol) WRITE (unit_nr, *) hills_env%invdt_history(hills_env%n_hills)
548  CALL close_file(unit_nr)
549  CALL m_mov(tmpname, filename)
550  END IF
551 
552  IF (modulo(multiple_walkers%n_hills_local, multiple_walkers%walkers_freq_comm) == 0) THEN
553  ! Store colvars information
554  ALLOCATE (ss0_save(n_colvar))
555  ALLOCATE (delta_s_save(n_colvar))
556  ALLOCATE (delta_s_ss0_buf(2, 0:n_colvar))
557  delta_s_ss0_buf = 0
558  DO i = 1, n_colvar
559  ss0_save(i) = colvars(i)%ss0
560  delta_s_save(i) = colvars(i)%delta_s
561  END DO
562 
563  ! Watch for other walkers's file and update
564  DO i = 1, multiple_walkers%walkers_tot_nr
565  IF (i == multiple_walkers%walker_id) THEN
566  ! Update local counter
567  multiple_walkers%walkers_status(i) = multiple_walkers%n_hills_local
568  cycle
569  END IF
570 
571  i_hills = multiple_walkers%walkers_status(i) + 1
572  filename = trim(multiple_walkers%walkers_file_name(i))//"_"// &
573  trim(adjustl(cp_to_string(i_hills)))
574 
575  IF (para_env%is_source()) THEN
576  INQUIRE (file=trim(filename), exist=exist)
577  END IF
578  CALL para_env%bcast(exist)
579  DO WHILE (exist)
580  ! Read information from the walker's file
581  ! We shouldn't care too much about the concurrency of these I/O instructions..
582  ! In case, they can be fixed in the future..
583  IF (para_env%is_source()) THEN
584  CALL open_file(file_name=filename, file_status="OLD", &
585  file_form="FORMATTED", file_action="READ", &
586  file_position="REWIND", unit_number=unit_nr)
587  READ (unit_nr, *) delta_s_ss0_buf(1, 0)
588  DO ih = 1, n_colvar
589  READ (unit_nr, *) delta_s_ss0_buf(1, ih)
590  READ (unit_nr, *) delta_s_ss0_buf(2, ih)
591  END DO
592  IF (hills_env%wtcontrol) READ (unit_nr, *) delta_s_ss0_buf(2, 0)
593  CALL close_file(unit_nr)
594  END IF
595  CALL para_env%bcast(delta_s_ss0_buf)
596  ww = delta_s_ss0_buf(1, 0)
597  IF (hills_env%wtcontrol) invdt = delta_s_ss0_buf(2, 0)
598  DO ih = 1, n_colvar
599  colvars(ih)%ss0 = delta_s_ss0_buf(1, ih)
600  colvars(ih)%delta_s = delta_s_ss0_buf(2, ih)
601  END DO
602 
603  ! Add this hill to the history dependent terms
604  IF (hills_env%wtcontrol) THEN
605  CALL add_hill_single(hills_env, colvars, ww, hills_env%n_hills, n_colvar, invdt=invdt)
606  ELSE
607  CALL add_hill_single(hills_env, colvars, ww, hills_env%n_hills, n_colvar)
608  END IF
609 
610  i_hills = i_hills + 1
611  filename = trim(multiple_walkers%walkers_file_name(i))//"_"// &
612  trim(adjustl(cp_to_string(i_hills)))
613  IF (para_env%is_source()) THEN
614  INQUIRE (file=trim(filename), exist=exist)
615  END IF
616  CALL para_env%bcast(exist)
617  END DO
618 
619  delta_hills = i_hills - 1 - multiple_walkers%walkers_status(i)
620  multiple_walkers%walkers_status(i) = i_hills - 1
621  iw = cp_print_key_unit_nr(logger, metadyn_section, "PRINT%PROGRAM_RUN_INFO", &
622  extension=".metadynLog")
623  IF (iw > 0) THEN
624  WRITE (iw, '(T2,A,I0,A,I0,A,I0,A)') 'WALKERS| Walker #', i, '. Reading [', delta_hills, &
625  '] Hills. Total number of Hills acquired [', multiple_walkers%walkers_status(i), ']'
626  END IF
627  CALL cp_print_key_finished_output(iw, logger, metadyn_section, &
628  "PRINT%PROGRAM_RUN_INFO")
629  END DO
630 
631  ! Restore colvars information
632  DO i = 1, n_colvar
633  colvars(i)%ss0 = ss0_save(i)
634  colvars(i)%delta_s = delta_s_save(i)
635  END DO
636  DEALLOCATE (ss0_save)
637  DEALLOCATE (delta_s_save)
638  END IF
639 
640  CALL timestop(handle)
641 
642  END SUBROUTINE synchronize_multiple_walkers
643 
644 ! **************************************************************************************************
645 !> \brief Add a single Hill
646 !> \param hills_env ...
647 !> \param colvars ...
648 !> \param ww ...
649 !> \param n_hills ...
650 !> \param n_colvar ...
651 !> \param invdt ...
652 !> \author Teodoro Laino [tlaino] - University of Zurich - 10.2008
653 ! **************************************************************************************************
654  SUBROUTINE add_hill_single(hills_env, colvars, ww, n_hills, n_colvar, invdt)
655  TYPE(hills_env_type), POINTER :: hills_env
656  TYPE(metavar_type), DIMENSION(:), POINTER :: colvars
657  REAL(kind=dp), INTENT(IN) :: ww
658  INTEGER, INTENT(INOUT) :: n_hills
659  INTEGER, INTENT(IN) :: n_colvar
660  REAL(kind=dp), INTENT(IN), OPTIONAL :: invdt
661 
662  CHARACTER(len=*), PARAMETER :: routinen = 'add_hill_single'
663 
664  INTEGER :: handle, i
665  LOGICAL :: wtcontrol
666  REAL(kind=dp), DIMENSION(:), POINTER :: tnp
667  REAL(kind=dp), DIMENSION(:, :), POINTER :: tmp
668 
669  CALL timeset(routinen, handle)
670 
671  wtcontrol = PRESENT(invdt)
672  NULLIFY (tmp, tnp)
673  IF (SIZE(hills_env%ss_history, 2) < n_hills + 1) THEN
674  ALLOCATE (tmp(n_colvar, n_hills + 100))
675  tmp(:, :n_hills) = hills_env%ss_history
676  tmp(:, n_hills + 1:) = 0.0_dp
677  DEALLOCATE (hills_env%ss_history)
678  hills_env%ss_history => tmp
679  NULLIFY (tmp)
680  END IF
681  IF (SIZE(hills_env%delta_s_history, 2) < n_hills + 1) THEN
682  ALLOCATE (tmp(n_colvar, n_hills + 100))
683  tmp(:, :n_hills) = hills_env%delta_s_history
684  tmp(:, n_hills + 1:) = 0.0_dp
685  DEALLOCATE (hills_env%delta_s_history)
686  hills_env%delta_s_history => tmp
687  NULLIFY (tmp)
688  END IF
689  IF (SIZE(hills_env%ww_history) < n_hills + 1) THEN
690  ALLOCATE (tnp(n_hills + 100))
691  tnp(1:n_hills) = hills_env%ww_history
692  tnp(n_hills + 1:) = 0.0_dp
693  DEALLOCATE (hills_env%ww_history)
694  hills_env%ww_history => tnp
695  NULLIFY (tnp)
696  END IF
697  IF (wtcontrol) THEN
698  IF (SIZE(hills_env%invdt_history) < n_hills + 1) THEN
699  ALLOCATE (tnp(n_hills + 100))
700  tnp(1:n_hills) = hills_env%invdt_history
701  tnp(n_hills + 1:) = 0.0_dp
702  DEALLOCATE (hills_env%invdt_history)
703  hills_env%invdt_history => tnp
704  NULLIFY (tnp)
705  END IF
706  END IF
707  n_hills = n_hills + 1
708  ! Now add the hill
709  DO i = 1, n_colvar
710  hills_env%ss_history(i, n_hills) = colvars(i)%ss0
711  hills_env%delta_s_history(i, n_hills) = colvars(i)%delta_s
712  END DO
713  hills_env%ww_history(n_hills) = ww
714  IF (wtcontrol) hills_env%invdt_history(n_hills) = invdt
715 
716  CALL timestop(handle)
717 
718  END SUBROUTINE add_hill_single
719 
720 ! **************************************************************************************************
721 !> \brief Restart Hills Information
722 !> \param ss_history ...
723 !> \param delta_s_history ...
724 !> \param ww_history ...
725 !> \param ww ...
726 !> \param n_hills ...
727 !> \param n_colvar ...
728 !> \param colvars ...
729 !> \param metadyn_section ...
730 !> \param invdt_history ...
731 !> \author Teodoro Laino [tlaino] - University of Zurich - 10.2008
732 ! **************************************************************************************************
733  SUBROUTINE restart_hills(ss_history, delta_s_history, ww_history, ww, &
734  n_hills, n_colvar, colvars, metadyn_section, invdt_history)
735  REAL(kind=dp), DIMENSION(:, :), POINTER :: ss_history, delta_s_history
736  REAL(kind=dp), DIMENSION(:), POINTER :: ww_history
737  REAL(kind=dp) :: ww
738  INTEGER, INTENT(IN) :: n_hills, n_colvar
739  TYPE(metavar_type), DIMENSION(:), POINTER :: colvars
740  TYPE(section_vals_type), POINTER :: metadyn_section
741  REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: invdt_history
742 
743  CHARACTER(len=*), PARAMETER :: routinen = 'restart_hills'
744 
745  INTEGER :: handle, i, j, ndum
746  LOGICAL :: explicit, wtcontrol
747  REAL(kind=dp) :: rval
748  REAL(kind=dp), DIMENSION(:), POINTER :: rvals
749  TYPE(section_vals_type), POINTER :: hills_history
750 
751  CALL timeset(routinen, handle)
752 
753  wtcontrol = PRESENT(invdt_history)
754  NULLIFY (rvals)
755  hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_POS")
756  CALL section_vals_get(hills_history, explicit=explicit)
757  IF (explicit) THEN
758  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", n_rep_val=ndum)
759  ! ss_history, delta_s_history, ww_history, invdt_history : deallocate and reallocate with the proper size
760  DEALLOCATE (ss_history)
761  DEALLOCATE (delta_s_history)
762  DEALLOCATE (ww_history)
763  IF (wtcontrol) THEN
764  DEALLOCATE (invdt_history)
765  END IF
766  !
767  cpassert(n_hills == ndum)
768  ALLOCATE (ss_history(n_colvar, n_hills))
769  ALLOCATE (delta_s_history(n_colvar, n_hills))
770  ALLOCATE (ww_history(n_hills))
771  IF (wtcontrol) THEN
772  ALLOCATE (invdt_history(n_hills))
773  END IF
774  !
775  DO i = 1, n_hills
776  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
777  i_rep_val=i, r_vals=rvals)
778  cpassert(SIZE(rvals) == n_colvar)
779  ss_history(1:n_colvar, i) = rvals
780  END DO
781  !
782  hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_SCALE")
783  CALL section_vals_get(hills_history, explicit=explicit)
784  IF (explicit) THEN
785  ! delta_s_history
786  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", n_rep_val=ndum)
787  cpassert(n_hills == ndum)
788  DO i = 1, n_hills
789  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
790  i_rep_val=i, r_vals=rvals)
791  cpassert(SIZE(rvals) == n_colvar)
792  delta_s_history(1:n_colvar, i) = rvals
793  END DO
794  ELSE
795  CALL cp_warn(__location__, &
796  "Section SPAWNED_HILLS_SCALE is not present! Setting the scales of the "// &
797  "restarted hills according the parameters specified in the input file.")
798  DO i = 1, n_hills
799  DO j = 1, n_colvar
800  delta_s_history(j, i) = colvars(i)%delta_s
801  END DO
802  END DO
803  END IF
804  !
805  hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_HEIGHT")
806  CALL section_vals_get(hills_history, explicit=explicit)
807  IF (explicit) THEN
808  ! ww_history
809  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
810  n_rep_val=ndum)
811  cpassert(n_hills == ndum)
812  DO i = 1, n_hills
813  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
814  i_rep_val=i, r_val=rval)
815  cpassert(SIZE(rvals) == n_colvar)
816  ww_history(i) = rval
817  END DO
818  ELSE
819  CALL cp_warn(__location__, &
820  "Section SPAWNED_HILLS_HEIGHT is not present! Setting the height of the"// &
821  " restarted hills according the parameters specified in the input file. ")
822  ww_history = ww
823  END IF
824  !
825  hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_INVDT")
826  CALL section_vals_get(hills_history, explicit=explicit)
827  IF (wtcontrol) THEN
828  IF (explicit) THEN
829  ! invdt_history
830  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
831  n_rep_val=ndum)
832  cpassert(n_hills == ndum)
833  DO i = 1, n_hills
834  CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
835  i_rep_val=i, r_val=rval)
836  cpassert(SIZE(rvals) == n_colvar)
837  invdt_history(i) = rval
838  END DO
839  ELSE
840  CALL cp_warn(__location__, &
841  "Section SPAWNED_HILLS_INVDT is not present! Restarting from standard"// &
842  " metadynamics run i.e. setting 1/(Delta T) equal to zero. ")
843  invdt_history = 0._dp
844  END IF
845  ELSE
846  IF (explicit) THEN
847  CALL cp_abort(__location__, &
848  "Found section SPAWNED_HILLS_INVDT while restarting a standard metadynamics run..."// &
849  " Cannot restart metadynamics from well-tempered MetaD runs. ")
850  END IF
851  END IF
852  END IF
853 
854  CALL timestop(handle)
855 
856  END SUBROUTINE restart_hills
857 
858 ! **************************************************************************************************
859 !> \brief Retrieves the iteration level for the metadynamics loop
860 !> \param meta_env ...
861 !> \param iter_nr ...
862 !> \author Teodoro Laino [tlaino] - University of Zurich - 10.2008
863 ! **************************************************************************************************
864  SUBROUTINE get_meta_iter_level(meta_env, iter_nr)
865  TYPE(meta_env_type), POINTER :: meta_env
866  INTEGER, INTENT(OUT) :: iter_nr
867 
868  IF (meta_env%do_multiple_walkers) THEN
869  iter_nr = meta_env%multiple_walkers%n_hills_local
870  ELSE
871  iter_nr = meta_env%hills_env%n_hills
872  END IF
873 
874  END SUBROUTINE get_meta_iter_level
875 
876 ! **************************************************************************************************
877 !> \brief ...
878 !> \param meta_env ...
879 !> \par History
880 !> 11.2007 [created] [tlaino]
881 !> \author Teodoro Laino - University of Zurich - 11.2007
882 ! **************************************************************************************************
883  SUBROUTINE meta_walls(meta_env)
884  TYPE(meta_env_type), POINTER :: meta_env
885 
886  INTEGER :: ih, iwall
887  REAL(dp) :: ddp, delta_s, dfunc, diff_ss, dp2, &
888  efunc, ww
889  TYPE(metavar_type), DIMENSION(:), POINTER :: colvars
890 
891  colvars => meta_env%metavar
892  ! Forces from the Walls
893  DO ih = 1, SIZE(colvars)
894  IF (colvars(ih)%do_wall) THEN
895  colvars(ih)%epot_walls = 0.0_dp
896  colvars(ih)%ff_walls = 0.0_dp
897  DO iwall = 1, SIZE(colvars(ih)%walls)
898  SELECT CASE (colvars(ih)%walls(iwall)%id_type)
900  ! Do Nothing.. treated in the main metadyn function
901  cycle
902  CASE (do_wall_quadratic)
903  diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos
904  IF (colvars(ih)%periodic) THEN
905  ! The difference of a periodic COLVAR is always within [-pi,pi]
906  diff_ss = sign(1.0_dp, asin(sin(diff_ss)))*acos(cos(diff_ss))
907  END IF
908  efunc = colvars(ih)%walls(iwall)%k_quadratic*diff_ss**2
909  dfunc = 2.0_dp*colvars(ih)%walls(iwall)%k_quadratic*diff_ss
910  SELECT CASE (colvars(ih)%walls(iwall)%id_direction)
911  CASE (do_wall_p)
912  IF (diff_ss > 0.0_dp) THEN
913  colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
914  colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
915  END IF
916  CASE (do_wall_m)
917  IF (diff_ss < 0.0_dp) THEN
918  colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
919  colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
920  END IF
921  END SELECT
922  CASE (do_wall_quartic)
923  diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos0
924  IF (colvars(ih)%periodic) THEN
925  ! The difference of a periodic COLVAR is always within [-pi,pi]
926  diff_ss = sign(1.0_dp, asin(sin(diff_ss)))*acos(cos(diff_ss))
927  END IF
928  efunc = colvars(ih)%walls(iwall)%k_quartic*diff_ss*diff_ss**4
929  dfunc = 4.0_dp*colvars(ih)%walls(iwall)%k_quartic*diff_ss**3
930  SELECT CASE (colvars(ih)%walls(iwall)%id_direction)
931  CASE (do_wall_p)
932  IF (diff_ss > 0.0_dp) THEN
933  colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
934  colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
935  END IF
936  CASE (do_wall_m)
937  IF (diff_ss < 0.0_dp) THEN
938  colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
939  colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
940  END IF
941  END SELECT
942  CASE (do_wall_gaussian)
943  diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos
944  IF (colvars(ih)%periodic) THEN
945  ! The difference of a periodic COLVAR is always within [-pi,pi]
946  diff_ss = sign(1.0_dp, asin(sin(diff_ss)))*acos(cos(diff_ss))
947  END IF
948  ww = colvars(ih)%walls(iwall)%ww_gauss
949  delta_s = colvars(ih)%walls(iwall)%sigma_gauss
950  ddp = (diff_ss)/delta_s
951  dp2 = ddp**2
952  efunc = ww*exp(-0.5_dp*dp2)
953  dfunc = -efunc*ddp/delta_s
954  colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
955  colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
956  END SELECT
957  END DO
958  END IF
959  END DO
960  END SUBROUTINE meta_walls
961 
962 END MODULE metadynamics_utils
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Definition: grid_common.h:117
Utility routines to open and close files. Tracking of preconnections.
Definition: cp_files.F:16
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
Definition: cp_files.F:308
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Definition: cp_files.F:119
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
types that represent a subsys, i.e. a part of the system
Interface for the force calculations.
recursive subroutine, public force_env_get(force_env, in_use, fist_env, qs_env, meta_env, fp_env, subsys, para_env, potential_energy, additional_potential, kinetic_energy, harmonic_shell, kinetic_shell, cell, sub_force_env, qmmm_env, qmmmx_env, eip_env, pwdft_env, globenv, input, force_env_section, method_name_id, root_section, mixed_env, nnp_env, embed_env)
returns various attributes about the force environment
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_wall_p
integer, parameter, public do_wall_m
integer, parameter, public do_wall_none
integer, parameter, public do_wall_gaussian
integer, parameter, public do_wall_quartic
integer, parameter, public do_wall_reflective
integer, parameter, public do_wall_quadratic
integer, parameter, public do_fe_meta
subroutine, public create_metavar_section(section)
creates the metavar section
represents an enumeration, i.e. a mapping between integers and strings
character(len=default_string_length) function, public enum_i2c(enum, i)
maps an integer to a string
represents keywords in an input
subroutine, public keyword_get(keyword, names, usage, description, type_of_var, n_var, default_value, lone_keyword_value, repeats, enum, citations)
...
objects that represent the structure of input sections and the data contained in an input section
type(section_type) function, pointer, public section_get_subsection(section, subsection_name)
returns the requested subsection
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
recursive subroutine, public section_release(section)
releases the given keyword list (see doc/ReferenceCounting.html)
recursive type(keyword_type) function, pointer, public section_get_keyword(section, keyword_name)
returns the requested keyword
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
integer, parameter, public default_path_length
Definition: kinds.F:58
Machine interface based on Fortran 2003 and POSIX.
Definition: machine.F:17
subroutine, public m_mov(source, TARGET)
...
Definition: machine.F:595
Interface to the message passing library MPI.
defines types for metadynamics calculation
subroutine, public metadyn_create(meta_env, n_colvar, dt, para_env, metadyn_section)
allocates a metadynamic environment (performs only minimal initialization)
Performs the metadynamics calculation.
subroutine, public synchronize_multiple_walkers(multiple_walkers, hills_env, colvars, n_colvar, metadyn_section)
Synchronize with the rest of the walkers.
subroutine, public metadyn_read(meta_env, force_env, root_section, para_env, fe_section)
reads metadynamics section
subroutine, public get_meta_iter_level(meta_env, iter_nr)
Retrieves the iteration level for the metadynamics loop.
subroutine, public meta_walls(meta_env)
...
subroutine, public restart_hills(ss_history, delta_s_history, ww_history, ww, n_hills, n_colvar, colvars, metadyn_section, invdt_history)
Restart Hills Information.
subroutine, public add_hill_single(hills_env, colvars, ww, n_hills, n_colvar, invdt)
Add a single Hill.
Definition of physical constants:
Definition: physcon.F:68
real(kind=dp), parameter, public kelvin
Definition: physcon.F:165