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 Contains utilities for a Dimer Method calculations
      10              : !> \par History
      11              : !>      none
      12              : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
      13              : ! **************************************************************************************************
      14              : MODULE dimer_utils
      15              :    USE cp_log_handling,                 ONLY: cp_logger_get_default_io_unit
      16              :    USE dimer_types,                     ONLY: dimer_env_type
      17              :    USE input_section_types,             ONLY: section_vals_get_subs_vals,&
      18              :                                               section_vals_remove_values,&
      19              :                                               section_vals_type,&
      20              :                                               section_vals_val_set
      21              :    USE kinds,                           ONLY: dp
      22              :    USE memory_utilities,                ONLY: reallocate
      23              : #include "../base/base_uses.f90"
      24              : 
      25              :    IMPLICIT NONE
      26              :    PRIVATE
      27              : 
      28              :    LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
      29              :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dimer_utils'
      30              : 
      31              :    PUBLIC :: rotate_dimer, update_dimer_vec, get_theta
      32              :    REAL(KIND=dp), PARAMETER, PUBLIC     :: dimer_thrs = EPSILON(0.0_dp)*1.0E4_dp
      33              : 
      34              : CONTAINS
      35              : 
      36              : ! **************************************************************************************************
      37              : !> \brief Performs a rotation of the unit dimer vector
      38              : !> \param nvec ...
      39              : !> \param theta ...
      40              : !> \param dt ...
      41              : !> \par History
      42              : !>      none
      43              : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
      44              : ! **************************************************************************************************
      45         1998 :    SUBROUTINE rotate_dimer(nvec, theta, dt)
      46              :       REAL(KIND=dp), DIMENSION(:), POINTER               :: nvec, theta
      47              :       REAL(KIND=dp)                                      :: dt
      48              : 
      49              :       INTEGER                                            :: output_unit
      50              :       LOGICAL                                            :: check
      51              : 
      52         1998 :       output_unit = cp_logger_get_default_io_unit()
      53              : 
      54              :       ! Orthogonality check for the rotation..
      55        36348 :       check = ABS(DOT_PRODUCT(nvec, theta)) < MAX(1.0E-9_dp, dimer_thrs)
      56         1998 :       IF (.NOT. check .AND. (output_unit > 0)) THEN
      57            0 :          WRITE (output_unit, *) "NVEC and THETA should be orthogonal! Residue: ", &
      58            0 :             ABS(DOT_PRODUCT(nvec, theta))
      59              :       END IF
      60         1998 :       CPASSERT(check)
      61        70698 :       nvec = nvec*COS(dt) + theta*SIN(dt)
      62              : 
      63         1998 :    END SUBROUTINE rotate_dimer
      64              : 
      65              : ! **************************************************************************************************
      66              : !> \brief Updates the orientation of the dimer vector in the input file
      67              : !> \param dimer_env ...
      68              : !> \param motion_section ...
      69              : !> \par History
      70              : !>      none
      71              : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
      72              : ! **************************************************************************************************
      73          858 :    SUBROUTINE update_dimer_vec(dimer_env, motion_section)
      74              :       TYPE(dimer_env_type), POINTER                      :: dimer_env
      75              :       TYPE(section_vals_type), POINTER                   :: motion_section
      76              : 
      77              :       INTEGER                                            :: i, i_rep_val, isize, j, size_array
      78          858 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: array
      79              :       TYPE(section_vals_type), POINTER                   :: nvec_section
      80              : 
      81              :       nvec_section => section_vals_get_subs_vals(motion_section, &
      82         1716 :                                                  "GEO_OPT%TRANSITION_STATE%DIMER%DIMER_VECTOR")
      83              :       ! Clean the content of the section first..
      84          858 :       CALL section_vals_remove_values(nvec_section)
      85              :       ! Fill in the section with the present values..
      86          858 :       size_array = 6
      87          858 :       isize = 0
      88          858 :       i_rep_val = 0
      89         2354 :       Main_Loop: DO i = 1, SIZE(dimer_env%nvec), size_array
      90         2354 :          ALLOCATE (array(size_array))
      91         2354 :          i_rep_val = i_rep_val + 1
      92        15386 :          DO j = 1, size_array
      93        13890 :             isize = isize + 1
      94        13890 :             array(j) = dimer_env%nvec(isize)
      95        15386 :             IF (isize == SIZE(dimer_env%nvec)) THEN
      96          858 :                CALL reallocate(array, 1, j)
      97              :                CALL section_vals_val_set(nvec_section, "_DEFAULT_KEYWORD_", r_vals_ptr=array, &
      98          858 :                                          i_rep_val=i_rep_val)
      99          858 :                EXIT Main_Loop
     100              :             END IF
     101              :          END DO
     102              :          CALL section_vals_val_set(nvec_section, "_DEFAULT_KEYWORD_", r_vals_ptr=array, &
     103         1496 :                                    i_rep_val=i_rep_val)
     104              :       END DO Main_Loop
     105          858 :       CPASSERT(isize == SIZE(dimer_env%nvec))
     106          858 :    END SUBROUTINE update_dimer_vec
     107              : 
     108              : ! **************************************************************************************************
     109              : !> \brief This function orthonormalize the vector for the rotational search
     110              : !> \param gradient ...
     111              : !> \param dimer_env ...
     112              : !> \param norm ...
     113              : !> \par History
     114              : !>      none
     115              : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
     116              : ! **************************************************************************************************
     117         1428 :    SUBROUTINE get_theta(gradient, dimer_env, norm)
     118              :       REAL(KIND=dp), DIMENSION(:)                        :: gradient
     119              :       TYPE(dimer_env_type), POINTER                      :: dimer_env
     120              :       REAL(KIND=dp), INTENT(OUT)                         :: norm
     121              : 
     122        49668 :       gradient = gradient - DOT_PRODUCT(gradient, dimer_env%nvec)*dimer_env%nvec
     123        25548 :       norm = SQRT(DOT_PRODUCT(gradient, gradient))
     124         1428 :       IF (norm < EPSILON(0.0_dp)) THEN
     125              :          ! This means that NVEC is totally aligned with minimum curvature mode
     126            0 :          gradient = 0.0_dp
     127              :       ELSE
     128        25548 :          gradient = gradient/norm
     129              :       END IF
     130         1428 :    END SUBROUTINE get_theta
     131              : 
     132              : END MODULE dimer_utils