LCOV - code coverage report
Current view: top level - src - metadynamics_utils.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:936074a) Lines: 93.3 % 464 433
Test Date: 2025-12-04 06:27:48 Functions: 100.0 % 8 8

            Line data    Source code
       1              : !--------------------------------------------------------------------------------------------------!
       2              : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3              : !   Copyright 2000-2025 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              : ! **************************************************************************************************
      14              : MODULE metadynamics_utils
      15              :    USE cp_files,                        ONLY: close_file,&
      16              :                                               open_file
      17              :    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
      18              :                                               cp_logger_type,&
      19              :                                               cp_to_string
      20              :    USE cp_output_handling,              ONLY: cp_print_key_finished_output,&
      21              :                                               cp_print_key_unit_nr
      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,&
      26              :                                               do_wall_gaussian,&
      27              :                                               do_wall_m,&
      28              :                                               do_wall_none,&
      29              :                                               do_wall_p,&
      30              :                                               do_wall_quadratic,&
      31              :                                               do_wall_quartic,&
      32              :                                               do_wall_reflective
      33              :    USE input_cp2k_free_energy,          ONLY: create_metavar_section
      34              :    USE input_enumeration_types,         ONLY: enum_i2c,&
      35              :                                               enumeration_type
      36              :    USE input_keyword_types,             ONLY: keyword_get,&
      37              :                                               keyword_type
      38              :    USE input_section_types,             ONLY: section_get_keyword,&
      39              :                                               section_get_subsection,&
      40              :                                               section_release,&
      41              :                                               section_type,&
      42              :                                               section_vals_get,&
      43              :                                               section_vals_get_subs_vals,&
      44              :                                               section_vals_type,&
      45              :                                               section_vals_val_get
      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,&
      52              :                                               metadyn_create,&
      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, &
      64              :              synchronize_multiple_walkers, &
      65              :              add_hill_single, &
      66              :              restart_hills, &
      67              :              get_meta_iter_level, &
      68              :              meta_walls
      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        18718 :    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         9359 :       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         9359 :       NULLIFY (subsys)
     103         9359 :       CALL timeset(routineN, handle)
     104              : 
     105         9359 :       CALL section_vals_get(fe_section, explicit=explicit)
     106         9359 :       IF (explicit) THEN
     107          166 :          number_allocated_colvars = 0
     108          166 :          CALL force_env_get(force_env, subsys=subsys)
     109          166 :          IF (ASSOCIATED(subsys%colvar_p)) THEN
     110          166 :             number_allocated_colvars = SIZE(subsys%colvar_p)
     111              :          END IF
     112          166 :          CALL section_vals_val_get(fe_section, "METHOD", i_val=id_method)
     113          166 :          IF (id_method /= do_fe_meta) THEN
     114           14 :             CALL timestop(handle)
     115           16 :             RETURN
     116              :          END IF
     117          152 :          metadyn_section => section_vals_get_subs_vals(fe_section, "METADYN")
     118          152 :          CPASSERT(.NOT. ASSOCIATED(meta_env))
     119              : 
     120          152 :          md_section => section_vals_get_subs_vals(root_section, "MOTION%MD")
     121          152 :          CALL section_vals_val_get(md_section, "TIMESTEP", r_val=dt)
     122              : 
     123          152 :          metavar_section => section_vals_get_subs_vals(metadyn_section, "METAVAR")
     124          152 :          CALL section_vals_get(metavar_section, n_repetition=n_colvar)
     125          152 :          ALLOCATE (meta_env)
     126              :          CALL metadyn_create(meta_env, n_colvar=n_colvar, &
     127          152 :                              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          152 :          CALL section_vals_val_get(metadyn_section, "USE_PLUMED", l_val=meta_env%use_plumed)
     131          152 :          IF (meta_env%use_plumed .EQV. .TRUE.) THEN
     132            2 :             CALL section_vals_val_get(metadyn_section, "PLUMED_INPUT_FILE", c_val=meta_env%plumed_input_file)
     133            2 :             meta_env%plumed_input_file = TRIM(meta_env%plumed_input_file)//CHAR(0)
     134            2 :             meta_env%langevin = .FALSE.
     135            2 :             CALL timestop(handle)
     136            2 :             RETURN
     137              :          END IF
     138              : 
     139          150 :          CALL section_vals_val_get(metadyn_section, "DO_HILLS", l_val=meta_env%do_hills)
     140          150 :          CALL section_vals_val_get(metadyn_section, "LAGRANGE", l_val=meta_env%extended_lagrange)
     141          150 :          CALL section_vals_val_get(metadyn_section, "TAMCSteps", i_val=meta_env%TAMCSteps)
     142          150 :          IF (meta_env%TAMCSteps < 0) THEN
     143            0 :             CPABORT("TAMCSteps must be positive!")
     144              :          END IF
     145          150 :          CALL section_vals_val_get(metadyn_section, "Timestep", r_val=meta_env%zdt)
     146          150 :          IF (meta_env%zdt <= 0.0_dp) THEN
     147            0 :             CPABORT("Timestep must be positive!")
     148              :          END IF
     149          150 :          CALL section_vals_val_get(metadyn_section, "WW", r_val=meta_env%hills_env%ww)
     150          150 :          CALL section_vals_val_get(metadyn_section, "NT_HILLS", i_val=meta_env%hills_env%nt_hills)
     151          150 :          CALL section_vals_val_get(metadyn_section, "MIN_NT_HILLS", i_val=meta_env%hills_env%min_nt_hills)
     152          150 :          IF (meta_env%hills_env%nt_hills <= 0) THEN
     153            4 :             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            4 :                          "Overriding input specification!")
     158              :          END IF
     159          150 :          check = meta_env%hills_env%nt_hills >= meta_env%hills_env%min_nt_hills
     160          150 :          IF (.NOT. check) &
     161              :             CALL cp_abort(__LOCATION__, "MIN_NT_HILLS must have a value smaller or equal to NT_HILLS! "// &
     162            0 :                           "Cross check with the input reference!")
     163              :          !RG Adaptive hills
     164          150 :          CALL section_vals_val_get(metadyn_section, "MIN_DISP", r_val=meta_env%hills_env%min_disp)
     165          150 :          CALL section_vals_val_get(metadyn_section, "OLD_HILL_NUMBER", i_val=meta_env%hills_env%old_hill_number)
     166          150 :          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          150 :          CALL section_vals_val_get(metadyn_section, "HILL_TAIL_CUTOFF", r_val=meta_env%hills_env%tail_cutoff)
     170          150 :          CALL section_vals_val_get(metadyn_section, "P_EXPONENT", i_val=meta_env%hills_env%p_exp)
     171          150 :          CALL section_vals_val_get(metadyn_section, "Q_EXPONENT", i_val=meta_env%hills_env%q_exp)
     172              : 
     173          150 :          CALL section_vals_val_get(metadyn_section, "SLOW_GROWTH", l_val=meta_env%hills_env%slow_growth)
     174              : 
     175              :          !RG Adaptive hills
     176          150 :          CALL section_vals_val_get(metadyn_section, "STEP_START_VAL", i_val=meta_env%n_steps)
     177          150 :          CPASSERT(meta_env%n_steps >= 0)
     178              :          CALL section_vals_val_get(metadyn_section, "NHILLS_START_VAL", &
     179          150 :                                    i_val=meta_env%hills_env%n_hills)
     180          150 :          CALL section_vals_val_get(metadyn_section, "TEMPERATURE", r_val=meta_env%temp_wanted)
     181          150 :          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          150 :                                    r_val=meta_env%toll_temp)
     184          150 :          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          150 :                                    r_val=meta_env%delta_t)
     187              :          CALL section_vals_val_get(metadyn_section, "WTGAMMA", explicit=check, &
     188          150 :                                    r_val=meta_env%wtgamma)
     189          150 :          IF (meta_env%well_tempered) THEN
     190            2 :             meta_env%hills_env%wtcontrol = meta_env%hills_env%wtcontrol .OR. check
     191            2 :             check = meta_env%hills_env%wtcontrol
     192            2 :             IF (.NOT. check) &
     193              :                CALL cp_abort(__LOCATION__, "When using Well-Tempered metadynamics, "// &
     194            0 :                              "DELTA_T (or WTGAMMA) should be explicitly specified.")
     195            2 :             IF (meta_env%extended_lagrange) &
     196              :                CALL cp_abort(__LOCATION__, &
     197            0 :                              "Well-Tempered metadynamics not possible with extended-lagrangian formulation.")
     198            2 :             IF (meta_env%hills_env%min_disp > 0.0_dp) &
     199              :                CALL cp_abort(__LOCATION__, &
     200            0 :                              "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          150 :                                    r_val=meta_env%avg_temp)
     205              :          ! Parsing Metavar Section
     206          342 :          DO i = 1, n_colvar
     207              :             CALL metavar_read(meta_env%metavar(i), meta_env%extended_lagrange, &
     208          192 :                               meta_env%langevin, i, metavar_section)
     209          192 :             check = (meta_env%metavar(i)%icolvar <= number_allocated_colvars)
     210          192 :             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          150 :                              TRIM(ADJUSTL(cp_to_string(number_allocated_colvars)))//") !")
     216              :          END DO
     217              : 
     218              :          ! Parsing the Multiple Walkers Info
     219          150 :          IF (meta_env%do_multiple_walkers) THEN
     220            8 :             NULLIFY (walkers_status)
     221            8 :             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            8 :                                       i_val=meta_env%multiple_walkers%walker_id)
     226              :             CALL section_vals_val_get(walkers_section, "NUMBER_OF_WALKERS", &
     227            8 :                                       i_val=meta_env%multiple_walkers%walkers_tot_nr)
     228              :             CALL section_vals_val_get(walkers_section, "WALKER_COMM_FREQUENCY", &
     229            8 :                                       i_val=meta_env%multiple_walkers%walkers_freq_comm)
     230              : 
     231              :             ! Handle status and file names
     232           24 :             ALLOCATE (meta_env%multiple_walkers%walkers_status(meta_env%multiple_walkers%walkers_tot_nr))
     233           24 :             ALLOCATE (meta_env%multiple_walkers%walkers_file_name(meta_env%multiple_walkers%walkers_tot_nr))
     234            8 :             CALL section_vals_val_get(walkers_section, "WALKERS_STATUS", explicit=explicit)
     235            8 :             IF (explicit) THEN
     236            4 :                CALL section_vals_val_get(walkers_section, "WALKERS_STATUS", i_vals=walkers_status)
     237            4 :                check = (SIZE(walkers_status) == meta_env%multiple_walkers%walkers_tot_nr)
     238            4 :                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            0 :                                 "RESTART_WALKERS in the EXT_RESTART section! ")
     244           20 :                meta_env%multiple_walkers%walkers_status = walkers_status
     245              :             ELSE
     246           12 :                meta_env%multiple_walkers%walkers_status = 0
     247              :             END IF
     248              :             meta_env%multiple_walkers%n_hills_local = &
     249            8 :                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            8 :                                       n_rep_val=n_rep)
     253            8 :             check = (n_rep == meta_env%multiple_walkers%walkers_tot_nr)
     254            8 :             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            0 :                              "RESTART_WALKERS in the EXT_RESTART section! ")
     260           40 :             DO i = 1, n_rep
     261              :                CALL section_vals_val_get(walkers_section, "WALKERS_FILE_NAME%_DEFAULT_KEYWORD_", &
     262           16 :                                          i_rep_val=i, c_val=walkers_file_name)
     263           24 :                meta_env%multiple_walkers%walkers_file_name(i) = walkers_file_name
     264              :             END DO
     265              :          END IF
     266              : 
     267              :          ! Print Metadynamics Info
     268          150 :          CALL print_metadyn_info(meta_env, n_colvar, metadyn_section)
     269              :       END IF
     270              : 
     271         9343 :       CALL timestop(handle)
     272              : 
     273         9359 :    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          150 :    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          150 :       CALL timeset(routineN, handle)
     297              : 
     298          150 :       logger => cp_get_default_logger()
     299              :       iw = cp_print_key_unit_nr(logger, metadyn_section, &
     300          150 :                                 "PRINT%PROGRAM_RUN_INFO", extension=".metadynLog")
     301          150 :       NULLIFY (section, enum, keyword)
     302          150 :       CALL create_metavar_section(section)
     303          150 :       wall_section => section_get_subsection(section, "WALL")
     304          150 :       IF (iw > 0) THEN
     305           75 :          WRITE (iw, '( /A )') ' METADYN| Meta Dynamics Protocol '
     306           75 :          WRITE (iw, '( A,T71,I10)') ' METADYN| Number of interval time steps to spawn hills', &
     307          150 :             meta_env%hills_env%nt_hills
     308           75 :          WRITE (iw, '( A,T71,I10)') ' METADYN| Number of previously spawned hills', &
     309          150 :             meta_env%hills_env%n_hills
     310           75 :          IF (meta_env%extended_lagrange) THEN
     311           28 :             WRITE (iw, '( A )') ' METADYN| Extended Lagrangian Scheme '
     312           28 :             IF (meta_env%tempcontrol) WRITE (iw, '( A,T71,F10.2)') &
     313           10 :                ' METADYN| Collective Variables Temperature control', meta_env%toll_temp
     314           28 :             IF (meta_env%langevin) THEN
     315            3 :                WRITE (iw, '(A,T71)') ' METADYN| Langevin Thermostat in use for COLVAR '
     316            3 :                WRITE (iw, '(A,T71,F10.4)') ' METADYN| Langevin Thermostat. Target Temperature = ', &
     317            6 :                   meta_env%temp_wanted*kelvin
     318              :             END IF
     319           28 :             WRITE (iw, '(A,T71,F10.4)') ' METADYN| COLVARS restarted average temperature ', &
     320           56 :                meta_env%avg_temp
     321              :          END IF
     322           75 :          IF (meta_env%do_hills) THEN
     323           59 :             WRITE (iw, '( A )') ' METADYN| Spawning the Hills '
     324           59 :             WRITE (iw, '( A,T71,F10.3)') ' METADYN| Height of the Spawned Gaussian', meta_env%hills_env%ww
     325              :             !RG Adaptive hills
     326           59 :             IF (meta_env%hills_env%min_disp > 0.0_dp) THEN
     327            2 :                WRITE (iw, '(A)') ' METADYN| Adapative meta time step is activated'
     328            2 :                WRITE (iw, '(A,T71,F10.4)') ' METADYN| Minimum displacement for next hill', &
     329            4 :                   meta_env%hills_env%min_disp
     330              :             END IF
     331              :             !RG Adaptive hills
     332              :          END IF
     333              : 
     334           75 :          IF (meta_env%well_tempered) THEN
     335            1 :             WRITE (iw, '( A )') ' METADYN| Well-Tempered metadynamics '
     336            1 :             IF (meta_env%delta_t > EPSILON(1._dp)) THEN
     337            1 :                WRITE (iw, '( A,T71,F10.3)') ' METADYN| Temperature parameter (Delta T) [K]', meta_env%delta_t*kelvin
     338              :             ELSE
     339            0 :                WRITE (iw, '( A,T71,F10.3)') ' METADYN| Temperature parameter (gamma)', meta_env%wtgamma
     340              :             END IF
     341              :          END IF
     342              : 
     343           75 :          IF (meta_env%do_multiple_walkers) THEN
     344            4 :             WRITE (iw, '( A,T71,A10)') ' METADYN| Multiple Walkers', '   ENABLED'
     345            4 :             WRITE (iw, '( A,T71,I10)') ' METADYN| Number of Multiple Walkers', &
     346            8 :                meta_env%multiple_walkers%walkers_tot_nr
     347            4 :             WRITE (iw, '( A,T71,I10)') ' METADYN| Local Walker ID', &
     348            8 :                meta_env%multiple_walkers%walker_id
     349            4 :             WRITE (iw, '( A,T71,I10)') ' METADYN| Walker Communication Frequency', &
     350            8 :                meta_env%multiple_walkers%walkers_freq_comm
     351           12 :             DO i = 1, meta_env%multiple_walkers%walkers_tot_nr
     352            8 :                my_tag = ""
     353            8 :                IF (i == meta_env%multiple_walkers%walker_id) my_tag = " ( Local )"
     354            8 :                my_id = '( '//TRIM(ADJUSTL(cp_to_string(i)))//' )'
     355            8 :                WRITE (iw, '(/,A,T71,A10)') ' WALKERS| Walker ID'//TRIM(my_tag), ADJUSTR(my_id)
     356            8 :                WRITE (iw, '(  A,T71,I10)') ' WALKERS| Number of Hills communicated', &
     357           16 :                   meta_env%multiple_walkers%walkers_status(i)
     358            8 :                WRITE (iw, '(  A,T24,A57)') ' WALKERS| Base Filename', &
     359           20 :                   ADJUSTR(meta_env%multiple_walkers%walkers_file_name(i) (1:57))
     360              :             END DO
     361            4 :             WRITE (iw, '(/)')
     362              :          END IF
     363              : 
     364           75 :          WRITE (iw, '( A,T71,I10)') ' METADYN| Number of collective variables', meta_env%n_colvar
     365          171 :          DO i = 1, n_colvar
     366           96 :             WRITE (iw, '( A )') '          '//'----------------------------------------------------------------------'
     367           96 :             WRITE (iw, '( A,T71,I10)') ' METAVARS| Collective Variable Number', meta_env%metavar(i)%icolvar
     368           96 :             IF (meta_env%extended_lagrange) THEN
     369           39 :                WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Lambda Parameter', meta_env%metavar(i)%lambda
     370           39 :                WRITE (iw, '( A,T66,F15.6)') ' METAVARS| Collective Variable Mass', meta_env%metavar(i)%mass
     371              :             END IF
     372           96 :             WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Scaling factor', meta_env%metavar(i)%delta_s
     373          101 :             IF (meta_env%langevin) WRITE (iw, '( A,T71,F10.6)') ' METAVARS| Friction for Langevin Thermostat', &
     374           10 :                meta_env%metavar(i)%gamma
     375           96 :             IF (meta_env%metavar(i)%do_wall) THEN
     376           18 :                WRITE (iw, '( A,T71,I10)') ' METAVARS| Number of Walls present', SIZE(meta_env%metavar(i)%walls)
     377           41 :                DO j = 1, SIZE(meta_env%metavar(i)%walls)
     378           23 :                   keyword => section_get_keyword(wall_section, "TYPE")
     379           23 :                   CALL keyword_get(keyword, enum=enum)
     380           23 :                   WRITE (iw, '(/,A,5X,I10,T50,A,T70,A11)') ' METAVARS| Wall Number:', j, 'Type of Wall:', &
     381           46 :                      ADJUSTR(TRIM(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_type)))
     382              :                   ! Type of wall IO
     383           23 :                   SELECT CASE (meta_env%metavar(i)%walls(j)%id_type)
     384              :                   CASE (do_wall_none)
     385              :                      ! Do Nothing
     386            4 :                      CYCLE
     387              :                   CASE (do_wall_reflective)
     388            4 :                      work_section => section_get_subsection(wall_section, "REFLECTIVE")
     389            4 :                      keyword => section_get_keyword(work_section, "DIRECTION")
     390            4 :                      CALL keyword_get(keyword, enum=enum)
     391            4 :                      WRITE (iw, '(A,T70,A11)') ' METAVARS| Wall direction', &
     392            8 :                         ADJUSTR(TRIM(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_direction)))
     393              :                   CASE (do_wall_quadratic)
     394           12 :                      work_section => section_get_subsection(wall_section, "QUADRATIC")
     395           12 :                      keyword => section_get_keyword(work_section, "DIRECTION")
     396           12 :                      CALL keyword_get(keyword, enum=enum)
     397           12 :                      WRITE (iw, '(A,T70,A11)') ' METAVARS| Wall direction', &
     398           24 :                         ADJUSTR(TRIM(enum_i2c(enum, meta_env%metavar(i)%walls(j)%id_direction)))
     399           12 :                      WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Constant K of the quadratic potential', &
     400           24 :                         meta_env%metavar(i)%walls(j)%k_quadratic
     401              :                   CASE (do_wall_gaussian)
     402            3 :                      WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Height of the Wall Gaussian', &
     403            6 :                         meta_env%metavar(i)%walls(j)%ww_gauss
     404            3 :                      WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Scale of the Wall Gaussian', &
     405           29 :                         meta_env%metavar(i)%walls(j)%sigma_gauss
     406              :                   END SELECT
     407           21 :                   WRITE (iw, '(A,T70,F11.6)') ' METAVARS| Wall location', &
     408           60 :                      meta_env%metavar(i)%walls(j)%pos
     409              :                END DO
     410              :             END IF
     411          171 :             WRITE (iw, '( A )') '          '//'----------------------------------------------------------------------'
     412              :          END DO
     413              :       END IF
     414          150 :       CALL section_release(section)
     415          150 :       CALL cp_print_key_finished_output(iw, logger, metadyn_section, "PRINT%PROGRAM_RUN_INFO")
     416              : 
     417          150 :       CALL timestop(handle)
     418              : 
     419          150 :    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          384 :    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          192 :       CALL timeset(routineN, handle)
     445              : 
     446          192 :       CALL section_vals_val_get(metavar_section, "COLVAR", i_rep_section=icol, i_val=metavar%icolvar)
     447          192 :       CALL section_vals_val_get(metavar_section, "SCALE", i_rep_section=icol, r_val=metavar%delta_s)
     448              :       ! Walls
     449          192 :       wall_section => section_vals_get_subs_vals(metavar_section, "WALL", i_rep_section=icol)
     450          192 :       CALL section_vals_get(wall_section, n_repetition=n_walls)
     451          192 :       IF (n_walls /= 0) THEN
     452           36 :          metavar%do_wall = .TRUE.
     453          154 :          ALLOCATE (metavar%walls(n_walls))
     454           82 :          DO i = 1, n_walls
     455           46 :             CALL section_vals_val_get(wall_section, "TYPE", i_rep_section=i, i_val=metavar%walls(i)%id_type)
     456           46 :             CALL section_vals_val_get(wall_section, "POSITION", i_rep_section=i, r_val=metavar%walls(i)%pos)
     457           36 :             SELECT CASE (metavar%walls(i)%id_type)
     458              :             CASE (do_wall_none)
     459              :                ! Just cycle..
     460            8 :                CYCLE
     461              :             CASE (do_wall_reflective)
     462            8 :                work_section => section_vals_get_subs_vals(wall_section, "REFLECTIVE", i_rep_section=i)
     463            8 :                CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
     464              :             CASE (do_wall_quadratic)
     465           24 :                work_section => section_vals_get_subs_vals(wall_section, "QUADRATIC", i_rep_section=i)
     466           24 :                CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
     467           24 :                CALL section_vals_val_get(work_section, "K", r_val=metavar%walls(i)%k_quadratic)
     468              :             CASE (do_wall_quartic)
     469            4 :                work_section => section_vals_get_subs_vals(wall_section, "QUARTIC", i_rep_section=i)
     470            4 :                CALL section_vals_val_get(work_section, "DIRECTION", i_val=metavar%walls(i)%id_direction)
     471            4 :                CALL section_vals_val_get(work_section, "K", r_val=metavar%walls(i)%k_quartic)
     472           12 :                SELECT CASE (metavar%walls(i)%id_direction)
     473              :                CASE (do_wall_m)
     474            2 :                   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            4 :                   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            6 :                work_section => section_vals_get_subs_vals(wall_section, "GAUSSIAN", i_rep_section=i)
     480            6 :                CALL section_vals_val_get(work_section, "WW", r_val=metavar%walls(i)%ww_gauss)
     481           52 :                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          192 :       IF (extended_lagrange) THEN
     487           78 :          CALL section_vals_val_get(metavar_section, "MASS", i_rep_section=icol, r_val=metavar%mass)
     488           78 :          CALL section_vals_val_get(metavar_section, "LAMBDA", i_rep_section=icol, r_val=metavar%lambda)
     489           78 :          IF (langevin) THEN
     490           10 :             CALL section_vals_val_get(metavar_section, "GAMMA", i_rep_section=icol, r_val=metavar%gamma)
     491              :          END IF
     492              :       END IF
     493              : 
     494          192 :       CALL timestop(handle)
     495              : 
     496          192 :    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          132 :    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          132 :       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:)           :: delta_s_save, ss0_save
     523          132 :       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          132 :       CALL timeset(routineN, handle)
     528              : 
     529          132 :       logger => cp_get_default_logger()
     530          132 :       para_env => logger%para_env
     531              : 
     532              :       ! Locally dump information on file..
     533          132 :       IF (para_env%is_source()) THEN
     534              :          ! Generate file name for the specific Hill
     535           66 :          i = multiple_walkers%walker_id
     536              :          filename = TRIM(multiple_walkers%walkers_file_name(i))//"_"// &
     537           66 :                     TRIM(ADJUSTL(cp_to_string(multiple_walkers%n_hills_local)))
     538           66 :          tmpname = TRIM(filename)//".tmp"
     539              :          CALL open_file(file_name=tmpname, file_status="UNKNOWN", &
     540              :                         file_form="FORMATTED", file_action="WRITE", &
     541           66 :                         file_position="APPEND", unit_number=unit_nr)
     542           66 :          WRITE (unit_nr, *) hills_env%ww_history(hills_env%n_hills)
     543          132 :          DO ih = 1, n_colvar
     544           66 :             WRITE (unit_nr, *) hills_env%ss_history(ih, hills_env%n_hills)
     545          132 :             WRITE (unit_nr, *) hills_env%delta_s_history(ih, hills_env%n_hills)
     546              :          END DO
     547           66 :          IF (hills_env%wtcontrol) WRITE (unit_nr, *) hills_env%invdt_history(hills_env%n_hills)
     548           66 :          CALL close_file(unit_nr)
     549           66 :          CALL m_mov(tmpname, filename)
     550              :       END IF
     551              : 
     552          132 :       IF (MODULO(multiple_walkers%n_hills_local, multiple_walkers%walkers_freq_comm) == 0) THEN
     553              :          ! Store colvars information
     554          396 :          ALLOCATE (ss0_save(n_colvar))
     555          264 :          ALLOCATE (delta_s_save(n_colvar))
     556          396 :          ALLOCATE (delta_s_ss0_buf(2, 0:n_colvar))
     557          924 :          delta_s_ss0_buf = 0
     558          264 :          DO i = 1, n_colvar
     559          132 :             ss0_save(i) = colvars(i)%ss0
     560          264 :             delta_s_save(i) = colvars(i)%delta_s
     561              :          END DO
     562              : 
     563              :          ! Watch for other walkers's file and update
     564          396 :          DO i = 1, multiple_walkers%walkers_tot_nr
     565          264 :             IF (i == multiple_walkers%walker_id) THEN
     566              :                ! Update local counter
     567          132 :                multiple_walkers%walkers_status(i) = multiple_walkers%n_hills_local
     568          132 :                CYCLE
     569              :             END IF
     570              : 
     571          132 :             i_hills = multiple_walkers%walkers_status(i) + 1
     572              :             filename = TRIM(multiple_walkers%walkers_file_name(i))//"_"// &
     573          132 :                        TRIM(ADJUSTL(cp_to_string(i_hills)))
     574              : 
     575          132 :             IF (para_env%is_source()) THEN
     576           66 :                INQUIRE (FILE=TRIM(filename), EXIST=exist)
     577              :             END IF
     578          132 :             CALL para_env%bcast(exist)
     579          230 :             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           98 :                IF (para_env%is_source()) THEN
     584              :                   CALL open_file(file_name=filename, file_status="OLD", &
     585              :                                  file_form="FORMATTED", file_action="READ", &
     586           49 :                                  file_position="REWIND", unit_number=unit_nr)
     587           49 :                   READ (unit_nr, *) delta_s_ss0_buf(1, 0)
     588           98 :                   DO ih = 1, n_colvar
     589           49 :                      READ (unit_nr, *) delta_s_ss0_buf(1, ih)
     590           98 :                      READ (unit_nr, *) delta_s_ss0_buf(2, ih)
     591              :                   END DO
     592           49 :                   IF (hills_env%wtcontrol) READ (unit_nr, *) delta_s_ss0_buf(2, 0)
     593           49 :                   CALL close_file(unit_nr)
     594              :                END IF
     595           98 :                CALL para_env%bcast(delta_s_ss0_buf)
     596           98 :                ww = delta_s_ss0_buf(1, 0)
     597           98 :                IF (hills_env%wtcontrol) invdt = delta_s_ss0_buf(2, 0)
     598          196 :                DO ih = 1, n_colvar
     599           98 :                   colvars(ih)%ss0 = delta_s_ss0_buf(1, ih)
     600          196 :                   colvars(ih)%delta_s = delta_s_ss0_buf(2, ih)
     601              :                END DO
     602              : 
     603              :                ! Add this hill to the history dependent terms
     604           98 :                IF (hills_env%wtcontrol) THEN
     605            0 :                   CALL add_hill_single(hills_env, colvars, ww, hills_env%n_hills, n_colvar, invdt=invdt)
     606              :                ELSE
     607           98 :                   CALL add_hill_single(hills_env, colvars, ww, hills_env%n_hills, n_colvar)
     608              :                END IF
     609              : 
     610           98 :                i_hills = i_hills + 1
     611              :                filename = TRIM(multiple_walkers%walkers_file_name(i))//"_"// &
     612           98 :                           TRIM(ADJUSTL(cp_to_string(i_hills)))
     613           98 :                IF (para_env%is_source()) THEN
     614           49 :                   INQUIRE (FILE=TRIM(filename), EXIST=exist)
     615              :                END IF
     616          230 :                CALL para_env%bcast(exist)
     617              :             END DO
     618              : 
     619          132 :             delta_hills = i_hills - 1 - multiple_walkers%walkers_status(i)
     620          132 :             multiple_walkers%walkers_status(i) = i_hills - 1
     621              :             iw = cp_print_key_unit_nr(logger, metadyn_section, "PRINT%PROGRAM_RUN_INFO", &
     622          132 :                                       extension=".metadynLog")
     623          132 :             IF (iw > 0) THEN
     624           66 :                WRITE (iw, '(T2,A,I0,A,I0,A,I0,A)') 'WALKERS| Walker #', i, '. Reading [', delta_hills, &
     625          132 :                   '] 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          264 :                                               "PRINT%PROGRAM_RUN_INFO")
     629              :          END DO
     630              : 
     631              :          ! Restore colvars information
     632          264 :          DO i = 1, n_colvar
     633          132 :             colvars(i)%ss0 = ss0_save(i)
     634          264 :             colvars(i)%delta_s = delta_s_save(i)
     635              :          END DO
     636          132 :          DEALLOCATE (ss0_save)
     637          132 :          DEALLOCATE (delta_s_save)
     638              :       END IF
     639              : 
     640          132 :       CALL timestop(handle)
     641              : 
     642          132 :    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         1214 :    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         1214 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: tnp
     667         1214 :       REAL(KIND=dp), DIMENSION(:, :), POINTER            :: tmp
     668              : 
     669         1214 :       CALL timeset(routineN, handle)
     670              : 
     671         1214 :       wtcontrol = PRESENT(invdt)
     672         1214 :       NULLIFY (tmp, tnp)
     673         1214 :       IF (SIZE(hills_env%ss_history, 2) < n_hills + 1) THEN
     674          456 :          ALLOCATE (tmp(n_colvar, n_hills + 100))
     675          836 :          tmp(:, :n_hills) = hills_env%ss_history
     676        26514 :          tmp(:, n_hills + 1:) = 0.0_dp
     677          114 :          DEALLOCATE (hills_env%ss_history)
     678          114 :          hills_env%ss_history => tmp
     679          114 :          NULLIFY (tmp)
     680              :       END IF
     681         1214 :       IF (SIZE(hills_env%delta_s_history, 2) < n_hills + 1) THEN
     682          456 :          ALLOCATE (tmp(n_colvar, n_hills + 100))
     683          836 :          tmp(:, :n_hills) = hills_env%delta_s_history
     684        26514 :          tmp(:, n_hills + 1:) = 0.0_dp
     685          114 :          DEALLOCATE (hills_env%delta_s_history)
     686          114 :          hills_env%delta_s_history => tmp
     687          114 :          NULLIFY (tmp)
     688              :       END IF
     689         1214 :       IF (SIZE(hills_env%ww_history) < n_hills + 1) THEN
     690          342 :          ALLOCATE (tnp(n_hills + 100))
     691          500 :          tnp(1:n_hills) = hills_env%ww_history
     692        11514 :          tnp(n_hills + 1:) = 0.0_dp
     693          114 :          DEALLOCATE (hills_env%ww_history)
     694          114 :          hills_env%ww_history => tnp
     695          114 :          NULLIFY (tnp)
     696              :       END IF
     697         1214 :       IF (wtcontrol) THEN
     698            4 :          IF (SIZE(hills_env%invdt_history) < n_hills + 1) THEN
     699            6 :             ALLOCATE (tnp(n_hills + 100))
     700            4 :             tnp(1:n_hills) = hills_env%invdt_history
     701          202 :             tnp(n_hills + 1:) = 0.0_dp
     702            2 :             DEALLOCATE (hills_env%invdt_history)
     703            2 :             hills_env%invdt_history => tnp
     704            2 :             NULLIFY (tnp)
     705              :          END IF
     706              :       END IF
     707         1214 :       n_hills = n_hills + 1
     708              :       ! Now add the hill
     709         2708 :       DO i = 1, n_colvar
     710         1494 :          hills_env%ss_history(i, n_hills) = colvars(i)%ss0
     711         2708 :          hills_env%delta_s_history(i, n_hills) = colvars(i)%delta_s
     712              :       END DO
     713         1214 :       hills_env%ww_history(n_hills) = ww
     714         1214 :       IF (wtcontrol) hills_env%invdt_history(n_hills) = invdt
     715              : 
     716         1214 :       CALL timestop(handle)
     717              : 
     718         1214 :    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          236 :    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          118 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: rvals
     749              :       TYPE(section_vals_type), POINTER                   :: hills_history
     750              : 
     751          118 :       CALL timeset(routineN, handle)
     752              : 
     753          118 :       wtcontrol = PRESENT(invdt_history)
     754          118 :       NULLIFY (rvals)
     755          118 :       hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_POS")
     756          118 :       CALL section_vals_get(hills_history, explicit=explicit)
     757          118 :       IF (explicit) THEN
     758           18 :          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           18 :          DEALLOCATE (ss_history)
     761           18 :          DEALLOCATE (delta_s_history)
     762           18 :          DEALLOCATE (ww_history)
     763           18 :          IF (wtcontrol) THEN
     764            0 :             DEALLOCATE (invdt_history)
     765              :          END IF
     766              :          !
     767           18 :          CPASSERT(n_hills == ndum)
     768           72 :          ALLOCATE (ss_history(n_colvar, n_hills))
     769           72 :          ALLOCATE (delta_s_history(n_colvar, n_hills))
     770           54 :          ALLOCATE (ww_history(n_hills))
     771           18 :          IF (wtcontrol) THEN
     772            0 :             ALLOCATE (invdt_history(n_hills))
     773              :          END IF
     774              :          !
     775          162 :          DO i = 1, n_hills
     776              :             CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     777          144 :                                       i_rep_val=i, r_vals=rvals)
     778          144 :             CPASSERT(SIZE(rvals) == n_colvar)
     779          658 :             ss_history(1:n_colvar, i) = rvals
     780              :          END DO
     781              :          !
     782           18 :          hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_SCALE")
     783           18 :          CALL section_vals_get(hills_history, explicit=explicit)
     784           18 :          IF (explicit) THEN
     785              :             ! delta_s_history
     786           18 :             CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", n_rep_val=ndum)
     787           18 :             CPASSERT(n_hills == ndum)
     788          162 :             DO i = 1, n_hills
     789              :                CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     790          144 :                                          i_rep_val=i, r_vals=rvals)
     791          144 :                CPASSERT(SIZE(rvals) == n_colvar)
     792          658 :                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            0 :                          "restarted hills according the parameters specified in the input file.")
     798            0 :             DO i = 1, n_hills
     799            0 :                DO j = 1, n_colvar
     800            0 :                   delta_s_history(j, i) = colvars(i)%delta_s
     801              :                END DO
     802              :             END DO
     803              :          END IF
     804              :          !
     805           18 :          hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_HEIGHT")
     806           18 :          CALL section_vals_get(hills_history, explicit=explicit)
     807           18 :          IF (explicit) THEN
     808              :             ! ww_history
     809              :             CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     810           18 :                                       n_rep_val=ndum)
     811           18 :             CPASSERT(n_hills == ndum)
     812          162 :             DO i = 1, n_hills
     813              :                CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     814          144 :                                          i_rep_val=i, r_val=rval)
     815          144 :                CPASSERT(SIZE(rvals) == n_colvar)
     816          306 :                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            0 :                          " restarted hills according the parameters specified in the input file. ")
     822            0 :             ww_history = ww
     823              :          END IF
     824              :          !
     825           18 :          hills_history => section_vals_get_subs_vals(metadyn_section, "SPAWNED_HILLS_INVDT")
     826           18 :          CALL section_vals_get(hills_history, explicit=explicit)
     827           72 :          IF (wtcontrol) THEN
     828            0 :             IF (explicit) THEN
     829              :                ! invdt_history
     830              :                CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     831            0 :                                          n_rep_val=ndum)
     832            0 :                CPASSERT(n_hills == ndum)
     833            0 :                DO i = 1, n_hills
     834              :                   CALL section_vals_val_get(hills_history, "_DEFAULT_KEYWORD_", &
     835            0 :                                             i_rep_val=i, r_val=rval)
     836            0 :                   CPASSERT(SIZE(rvals) == n_colvar)
     837            0 :                   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            0 :                             " metadynamics run i.e. setting 1/(Delta T) equal to zero. ")
     843            0 :                invdt_history = 0._dp
     844              :             END IF
     845              :          ELSE
     846           18 :             IF (explicit) THEN
     847              :                CALL cp_abort(__LOCATION__, &
     848              :                              "Found section SPAWNED_HILLS_INVDT while restarting a standard metadynamics run..."// &
     849            0 :                              " Cannot restart metadynamics from well-tempered MetaD runs. ")
     850              :             END IF
     851              :          END IF
     852              :       END IF
     853              : 
     854          118 :       CALL timestop(handle)
     855              : 
     856          118 :    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        12394 :    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        12394 :       IF (meta_env%do_multiple_walkers) THEN
     869          540 :          iter_nr = meta_env%multiple_walkers%n_hills_local
     870              :       ELSE
     871        11854 :          iter_nr = meta_env%hills_env%n_hills
     872              :       END IF
     873              : 
     874        12394 :    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        13928 :    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        13928 :       TYPE(metavar_type), DIMENSION(:), POINTER          :: colvars
     890              : 
     891        13928 :       colvars => meta_env%metavar
     892              :       ! Forces from the Walls
     893        28600 :       DO ih = 1, SIZE(colvars)
     894        28600 :          IF (colvars(ih)%do_wall) THEN
     895        11276 :             colvars(ih)%epot_walls = 0.0_dp
     896        11276 :             colvars(ih)%ff_walls = 0.0_dp
     897        23162 :             DO iwall = 1, SIZE(colvars(ih)%walls)
     898        11276 :                SELECT CASE (colvars(ih)%walls(iwall)%id_type)
     899              :                CASE (do_wall_reflective, do_wall_none)
     900              :                   ! Do Nothing.. treated in the main metadyn function
     901        10564 :                   CYCLE
     902              :                CASE (do_wall_quadratic)
     903        10564 :                   diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos
     904        10564 :                   IF (colvars(ih)%periodic) THEN
     905              :                      ! The difference of a periodic COLVAR is always within [-pi,pi]
     906            0 :                      diff_ss = SIGN(1.0_dp, ASIN(SIN(diff_ss)))*ACOS(COS(diff_ss))
     907              :                   END IF
     908        10564 :                   efunc = colvars(ih)%walls(iwall)%k_quadratic*diff_ss**2
     909        10564 :                   dfunc = 2.0_dp*colvars(ih)%walls(iwall)%k_quadratic*diff_ss
     910        21084 :                   SELECT CASE (colvars(ih)%walls(iwall)%id_direction)
     911              :                   CASE (do_wall_p)
     912        10316 :                      IF (diff_ss > 0.0_dp) THEN
     913          170 :                         colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
     914          170 :                         colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
     915              :                      END IF
     916              :                   CASE (do_wall_m)
     917        10564 :                      IF (diff_ss < 0.0_dp) THEN
     918           30 :                         colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
     919           30 :                         colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
     920              :                      END IF
     921              :                   END SELECT
     922              :                CASE (do_wall_quartic)
     923          204 :                   diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos0
     924          204 :                   IF (colvars(ih)%periodic) THEN
     925              :                      ! The difference of a periodic COLVAR is always within [-pi,pi]
     926            0 :                      diff_ss = SIGN(1.0_dp, ASIN(SIN(diff_ss)))*ACOS(COS(diff_ss))
     927              :                   END IF
     928          204 :                   efunc = colvars(ih)%walls(iwall)%k_quartic*diff_ss*diff_ss**4
     929          204 :                   dfunc = 4.0_dp*colvars(ih)%walls(iwall)%k_quartic*diff_ss**3
     930          812 :                   SELECT CASE (colvars(ih)%walls(iwall)%id_direction)
     931              :                   CASE (do_wall_p)
     932          102 :                      IF (diff_ss > 0.0_dp) THEN
     933           46 :                         colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
     934           46 :                         colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
     935              :                      END IF
     936              :                   CASE (do_wall_m)
     937          204 :                      IF (diff_ss < 0.0_dp) THEN
     938           54 :                         colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
     939           54 :                         colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
     940              :                      END IF
     941              :                   END SELECT
     942              :                CASE (do_wall_gaussian)
     943          506 :                   diff_ss = colvars(ih)%ss0 - colvars(ih)%walls(iwall)%pos
     944          506 :                   IF (colvars(ih)%periodic) THEN
     945              :                      ! The difference of a periodic COLVAR is always within [-pi,pi]
     946            0 :                      diff_ss = SIGN(1.0_dp, ASIN(SIN(diff_ss)))*ACOS(COS(diff_ss))
     947              :                   END IF
     948          506 :                   ww = colvars(ih)%walls(iwall)%ww_gauss
     949          506 :                   delta_s = colvars(ih)%walls(iwall)%sigma_gauss
     950          506 :                   ddp = (diff_ss)/delta_s
     951          506 :                   dp2 = ddp**2
     952          506 :                   efunc = ww*EXP(-0.5_dp*dp2)
     953          506 :                   dfunc = -efunc*ddp/delta_s
     954          506 :                   colvars(ih)%ff_walls = colvars(ih)%ff_walls - dfunc
     955        12392 :                   colvars(ih)%epot_walls = colvars(ih)%epot_walls + efunc
     956              :                END SELECT
     957              :             END DO
     958              :          END IF
     959              :       END DO
     960        13928 :    END SUBROUTINE meta_walls
     961              : 
     962              : END MODULE metadynamics_utils
        

Generated by: LCOV version 2.0-1