Line data    Source code

       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 Callback used by global geometry optimization schemes
      10             : !> \author Ole Schuett
      11             : ! **************************************************************************************************
      12             : MODULE glbopt_callback
      13             :    USE cp_subsys_types,                 ONLY: cp_subsys_get,&
      14             :                                               cp_subsys_type,&
      15             :                                               pack_subsys_particles
      16             :    USE force_env_types,                 ONLY: force_env_get,&
      17             :                                               force_env_type
      18             :    USE kinds,                           ONLY: dp
      19             :    USE md_ener_types,                   ONLY: md_ener_type
      20             :    USE md_environment_types,            ONLY: get_md_env,&
      21             :                                               md_environment_type
      22             :    USE mdctrl_types,                    ONLY: glbopt_mdctrl_data_type
      23             : #include "../base/base_uses.f90"
      24             : 
      25             :    IMPLICIT NONE
      26             :    PRIVATE
      27             : 
      28             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'glbopt_callback'
      29             : 
      30             :    PUBLIC :: glbopt_md_callback
      31             : 
      32             : CONTAINS
      33             : 
      34             : ! **************************************************************************************************
      35             : !> \brief Callback used to hook into the main MD-loop.
      36             : !>        It recognizes and counts bumps in the potential energy.
      37             : !>        When MD_BUMPS_MAX is reached, the MD simulations is stoped.
      38             : !> \param mdctrl_data ...
      39             : !> \param md_env ...
      40             : !> \param should_stop ...
      41             : !> \author Ole Schuett
      42             : ! **************************************************************************************************
      43        4075 :    SUBROUTINE glbopt_md_callback(mdctrl_data, md_env, should_stop)
      44             :       TYPE(glbopt_mdctrl_data_type), POINTER             :: mdctrl_data
      45             :       TYPE(md_environment_type), POINTER                 :: md_env
      46             :       LOGICAL, INTENT(inout)                             :: should_stop
      47             : 
      48             :       INTEGER                                            :: i, iw, n_atoms
      49             :       INTEGER, POINTER                                   :: itimes
      50             :       LOGICAL                                            :: passed_minimum
      51        4075 :       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:)           :: positions
      52             :       TYPE(cp_subsys_type), POINTER                      :: subsys
      53             :       TYPE(force_env_type), POINTER                      :: force_env
      54             :       TYPE(md_ener_type), POINTER                        :: md_ener
      55             : 
      56        4075 :       CPASSERT(ASSOCIATED(mdctrl_data))
      57        4075 :       CPASSERT(ASSOCIATED(md_env))
      58             : 
      59        4075 :       iw = mdctrl_data%output_unit
      60             : 
      61             :       ! add new potential energy value to history
      62        4075 :       NULLIFY (md_ener, itimes)
      63        4075 :       CALL get_md_env(md_env=md_env, md_ener=md_ener, itimes=itimes, force_env=force_env)
      64        4075 :       mdctrl_data%itimes = itimes
      65             : 
      66       24450 :       mdctrl_data%epot_history(:) = EOSHIFT(mdctrl_data%epot_history, shift=-1)
      67        4075 :       mdctrl_data%epot_history(1) = md_ener%epot
      68             : 
      69             :       ! check if we passed a minimum
      70        4075 :       passed_minimum = .TRUE.
      71       12225 :       DO i = 1, mdctrl_data%bump_steps_upwards
      72        8150 :          IF (mdctrl_data%epot_history(i) <= mdctrl_data%epot_history(i + 1)) &
      73        7547 :             passed_minimum = .FALSE.
      74             :       END DO
      75             : 
      76       12225 :       DO i = mdctrl_data%bump_steps_upwards + 1, mdctrl_data%bump_steps_upwards + mdctrl_data%bump_steps_downwards
      77        8150 :          IF (mdctrl_data%epot_history(i) >= mdctrl_data%epot_history(i + 1)) &
      78        8849 :             passed_minimum = .FALSE.
      79             :       END DO
      80             : 
      81             :       ! count the passed bumps and stop md_run when md_bumps_max is reached.
      82        4075 :       IF (passed_minimum) &
      83         144 :          mdctrl_data%md_bump_counter = mdctrl_data%md_bump_counter + 1
      84             : 
      85        4075 :       IF (mdctrl_data%md_bump_counter >= mdctrl_data%md_bumps_max) THEN
      86          48 :          should_stop = .TRUE.
      87          48 :          IF (iw > 0) WRITE (iw, "(A)") " GLBOPT| Stopping MD because of MD_BUMPS_MAX."
      88             :       END IF
      89             : 
      90        4075 :       CALL force_env_get(force_env, subsys=subsys)
      91        4075 :       CALL cp_subsys_get(subsys, natom=n_atoms)
      92       12225 :       ALLOCATE (positions(3*n_atoms))
      93        4075 :       CALL pack_subsys_particles(subsys, r=positions)
      94             : 
      95        8150 :    END SUBROUTINE glbopt_md_callback
      96             : 
      97             : END MODULE glbopt_callback
      98             :