LCOV - code coverage report
Current view: top level - src/swarm - glbopt_history.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:42dac4a) Lines: 83.3 % 120 100
Test Date: 2025-07-25 12:55:17 Functions: 64.3 % 14 9

            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 History of minima, calculates, stores and compares fingerprints of minima.
      10              : !>        Used by Minima Hopping and Minima Crawling.
      11              : !> \author Ole Schuett
      12              : ! **************************************************************************************************
      13              : MODULE glbopt_history
      14              :    USE input_section_types,             ONLY: section_vals_type,&
      15              :                                               section_vals_val_get
      16              :    USE kinds,                           ONLY: dp
      17              : #include "../base/base_uses.f90"
      18              : 
      19              :    IMPLICIT NONE
      20              :    PRIVATE
      21              : 
      22              :    TYPE history_fingerprint_type
      23              :       PRIVATE
      24              :       REAL(KIND=dp)                            :: Epot = 0.0
      25              :       REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: goedecker
      26              :    END TYPE history_fingerprint_type
      27              : 
      28              :    TYPE history_entry_type
      29              :       TYPE(history_fingerprint_type), POINTER :: p => Null()
      30              :       INTEGER                                 :: id = -1
      31              :    END TYPE history_entry_type
      32              : 
      33              :    TYPE history_type
      34              :       PRIVATE
      35              :       TYPE(history_entry_type), DIMENSION(:), POINTER :: entries => Null()
      36              :       INTEGER                              :: length = 0
      37              :       INTEGER                              :: iw = -1
      38              :       REAL(KIND=dp)                        :: E_precision = 0.0
      39              :       REAL(KIND=dp)                        :: FP_precision = 0.0
      40              :    END TYPE history_type
      41              : 
      42              :    PUBLIC :: history_type, history_fingerprint_type
      43              :    PUBLIC :: history_init, history_finalize
      44              :    PUBLIC :: history_add, history_lookup
      45              :    PUBLIC :: history_fingerprint
      46              :    PUBLIC :: history_fingerprint_match
      47              : 
      48              :    LOGICAL, PARAMETER                     :: debug = .FALSE.
      49              :    INTEGER, PARAMETER                     :: history_grow_unit = 1000
      50              : CONTAINS
      51              : 
      52              : ! **************************************************************************************************
      53              : !> \brief Initializes a history.
      54              : !> \param history ...
      55              : !> \param history_section ...
      56              : !> \param iw ...
      57              : !> \author Ole Schuett
      58              : ! **************************************************************************************************
      59            3 :    SUBROUTINE history_init(history, history_section, iw)
      60              :       TYPE(history_type), INTENT(INOUT)                  :: history
      61              :       TYPE(section_vals_type), POINTER                   :: history_section
      62              :       INTEGER                                            :: iw
      63              : 
      64         3003 :       ALLOCATE (history%entries(history_grow_unit))
      65            3 :       history%iw = iw
      66              :       CALL section_vals_val_get(history_section, "ENERGY_PRECISION", &
      67            3 :                                 r_val=history%E_precision)
      68              :       CALL section_vals_val_get(history_section, "FINGERPRINT_PRECISION", &
      69            3 :                                 r_val=history%FP_precision)
      70              : 
      71            3 :       IF (iw > 0) THEN
      72              :          WRITE (iw, '(A,T66,E15.3)') &
      73            3 :             " GLBOPT| History energy precision", history%E_precision
      74              :          WRITE (iw, '(A,T66,E15.3)') &
      75            3 :             " GLBOPT| History fingerprint precision", history%FP_precision
      76              :       END IF
      77            3 :    END SUBROUTINE history_init
      78              : 
      79              : ! **************************************************************************************************
      80              : !> \brief Calculates a fingerprint for a given configuration.
      81              : !> \param Epot ...
      82              : !> \param pos ...
      83              : !> \return ...
      84              : !> \author Ole Schuett
      85              : ! **************************************************************************************************
      86           56 :    FUNCTION history_fingerprint(Epot, pos) RESULT(fp)
      87              :       REAL(KIND=dp), INTENT(IN)                          :: Epot
      88              :       REAL(KIND=dp), DIMENSION(:), INTENT(IN)            :: pos
      89              :       TYPE(history_fingerprint_type)                     :: fp
      90              : 
      91              :       INTEGER                                            :: handle
      92           28 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: tmp
      93              : 
      94           28 :       CALL timeset("glbopt_history_fingerprint", handle)
      95              : 
      96           28 :       NULLIFY (tmp)
      97           28 :       fp%Epot = Epot
      98           28 :       CALL goedecker_fingerprint(pos, tmp)
      99              : 
     100              :       !copy pointer to allocatable
     101           84 :       ALLOCATE (fp%goedecker(SIZE(tmp)))
     102          308 :       fp%goedecker(:) = tmp
     103           28 :       DEALLOCATE (tmp)
     104              : 
     105           28 :       CALL timestop(handle)
     106           28 :    END FUNCTION history_fingerprint
     107              : 
     108              : ! **************************************************************************************************
     109              : !> \brief Helper routine for history_fingerprint.
     110              : !>        Calculates a fingerprint based on inter-atomic distances.
     111              : !> \param pos ...
     112              : !> \param res ...
     113              : !> \author Stefan Goedecker
     114              : ! **************************************************************************************************
     115           28 :    SUBROUTINE goedecker_fingerprint(pos, res)
     116              :       REAL(KIND=dp), DIMENSION(:), INTENT(IN)            :: pos
     117              :       REAL(KIND=dp), DIMENSION(:), POINTER               :: res
     118              : 
     119              :       INTEGER                                            :: i, info, j, N
     120              :       REAL(KIND=dp)                                      :: d2, t
     121           28 :       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :)        :: matrix, work
     122              :       REAL(KIND=dp), DIMENSION(3)                        :: d
     123              : 
     124           28 :       IF (ASSOCIATED(res)) CPABORT("goedecker_fingerprint: res already allocated")
     125           28 :       N = SIZE(pos)/3 ! number of atoms
     126              : 
     127          168 :       ALLOCATE (matrix(N, N), work(N, N))
     128          308 :       DO i = 1, N
     129          280 :          matrix(i, i) = 1.0
     130         1568 :          DO j = i + 1, N
     131         5040 :             d = pos(3*i - 2:3*i) - pos(3*j - 2:3*j)
     132         5040 :             d2 = SUM(d**2)
     133         1260 :             t = EXP(-0.5*d2)
     134         1260 :             matrix(i, j) = t
     135         1540 :             matrix(j, i) = t
     136              :          END DO
     137              :       END DO
     138           84 :       ALLOCATE (res(N))
     139              :       ! matrix values are garbage on exit because of jobz='N'
     140           28 :       CALL dsyev('N', 'U', N, matrix, N, res, work, N**2, info)
     141           28 :       IF (info /= 0) CPABORT("goedecker_fingerprint: DSYEV failed")
     142           28 :    END SUBROUTINE goedecker_fingerprint
     143              : 
     144              : ! **************************************************************************************************
     145              : !> \brief Checks if two given fingerprints match.
     146              : !> \param history ...
     147              : !> \param fp1 ...
     148              : !> \param fp2 ...
     149              : !> \return ...
     150              : !> \author Ole Schuett
     151              : ! **************************************************************************************************
     152           11 :    FUNCTION history_fingerprint_match(history, fp1, fp2) RESULT(res)
     153              :       TYPE(history_type), INTENT(IN)                     :: history
     154              :       TYPE(history_fingerprint_type), INTENT(IN)         :: fp1, fp2
     155              :       LOGICAL                                            :: res
     156              : 
     157              :       res = (ABS(fp1%Epot - fp2%Epot) < history%E_precision) .AND. &
     158           11 :             (fingerprint_distance(fp1, fp2) < history%fp_precision)
     159              : 
     160           11 :    END FUNCTION history_fingerprint_match
     161              : 
     162              : ! **************************************************************************************************
     163              : !> \brief Helper routine for history_fingerprint_match
     164              : !>        Calculates the distance between two given fingerprints.
     165              : !> \param fp1 ...
     166              : !> \param fp2 ...
     167              : !> \return ...
     168              : !> \author Stefan Goedecker
     169              : ! **************************************************************************************************
     170           12 :    PURE FUNCTION fingerprint_distance(fp1, fp2) RESULT(res)
     171              :       TYPE(history_fingerprint_type), INTENT(IN)         :: fp1, fp2
     172              :       REAL(KIND=dp)                                      :: res
     173              : 
     174          132 :       res = SQRT(SUM((fp1%goedecker - fp2%goedecker)**2)/SIZE(fp1%goedecker))
     175           12 :    END FUNCTION fingerprint_distance
     176              : 
     177              : ! **************************************************************************************************
     178              : !> \brief Addes a new fingerprints to the history.
     179              : !>        Optionally, an abitrary id can be stored alongside the fingerprint.
     180              : !> \param history ...
     181              : !> \param fingerprint ...
     182              : !> \param id ...
     183              : !> \author Ole Schuett
     184              : ! **************************************************************************************************
     185           17 :    SUBROUTINE history_add(history, fingerprint, id)
     186              :       TYPE(history_type), INTENT(INOUT)                  :: history
     187              :       TYPE(history_fingerprint_type), INTENT(IN)         :: fingerprint
     188              :       INTEGER, INTENT(IN), OPTIONAL                      :: id
     189              : 
     190              :       INTEGER                                            :: handle, i, k, n
     191           17 :       TYPE(history_entry_type), DIMENSION(:), POINTER    :: tmp
     192              : 
     193           17 :       CALL timeset("glbopt_history_add", handle)
     194              : 
     195           17 :       n = SIZE(history%entries)
     196           17 :       IF (n == history%length) THEN
     197              :          ! grow history%entries array
     198            0 :          tmp => history%entries
     199            0 :          ALLOCATE (history%entries(n + history_grow_unit))
     200            0 :          history%entries(1:n) = tmp(:)
     201            0 :          DEALLOCATE (tmp)
     202            0 :          n = n + history_grow_unit
     203              :       END IF
     204              : 
     205           17 :       k = interpolation_search(history, fingerprint%Epot)
     206              : 
     207              :       !history%entries(k+1:) = history%entries(k:n-1)
     208              :       !Workaround for an XLF bug - pointer array copy does
     209              :       !not work correctly
     210        16981 :       DO i = n, k + 1, -1
     211        16981 :          history%entries(i) = history%entries(i - 1)
     212              :       END DO
     213              : 
     214           17 :       ALLOCATE (history%entries(k)%p)
     215           17 :       history%entries(k)%p = fingerprint
     216           17 :       IF (PRESENT(id)) &
     217           10 :          history%entries(k)%id = id
     218           17 :       history%length = history%length + 1
     219              : 
     220              :       IF (debug) THEN
     221              :          ! check history for correct order
     222              :          DO k = 1, history%length
     223              :             !WRITE(*,*) "history: ", k, "Epot",history%entries(k)%p%Epot
     224              :             IF (k > 1) THEN
     225              :                IF (history%entries(k - 1)%p%Epot > history%entries(k)%p%Epot) &
     226              :                   CPABORT("history_add: history in wrong order")
     227              :             END IF
     228              :          END DO
     229              :       END IF
     230              : 
     231           17 :       CALL timestop(handle)
     232           17 :    END SUBROUTINE history_add
     233              : 
     234              : ! **************************************************************************************************
     235              : !> \brief Checks if a given fingerprints is contained in the history.
     236              : !> \param history ...
     237              : !> \param fingerprint ...
     238              : !> \param found ...
     239              : !> \param id ...
     240              : !> \author Ole Schuett
     241              : ! **************************************************************************************************
     242           25 :    SUBROUTINE history_lookup(history, fingerprint, found, id)
     243              :       TYPE(history_type), INTENT(IN)                     :: history
     244              :       TYPE(history_fingerprint_type), INTENT(IN)         :: fingerprint
     245              :       LOGICAL, INTENT(OUT)                               :: found
     246              :       INTEGER, INTENT(OUT), OPTIONAL                     :: id
     247              : 
     248              :       INTEGER                                            :: found_i, handle, i, k, k_max, k_min
     249              :       REAL(KIND=dp)                                      :: best_match, dist, Epot
     250              : 
     251           25 :       CALL timeset("glbopt_history_lookup", handle)
     252              : 
     253           25 :       found = .FALSE.
     254           25 :       IF (PRESENT(id)) id = -1
     255           25 :       best_match = HUGE(1.0_dp)
     256              : 
     257           25 :       IF (history%length > 0) THEN
     258           22 :          Epot = fingerprint%Epot
     259           22 :          k = interpolation_search(history, fingerprint%Epot)
     260              : 
     261           28 :          DO k_min = k - 1, 1, -1
     262           28 :             IF (history%entries(k_min)%p%Epot < Epot - history%E_precision) EXIT
     263              :          END DO
     264              : 
     265           25 :          DO k_max = k, history%length
     266           25 :             IF (history%entries(k_max)%p%Epot > Epot + history%E_precision) EXIT
     267              :          END DO
     268              : 
     269           22 :          k_min = MAX(k_min + 1, 1)
     270           22 :          k_max = MIN(k_max - 1, history%length)
     271              : 
     272              :          IF (debug) found_i = -1
     273              : 
     274           31 :          DO i = k_min, k_max
     275            9 :             dist = fingerprint_distance(fingerprint, history%entries(i)%p)
     276              :             !WRITE(*,*) "entry ", i, " dist: ",dist
     277           31 :             IF (dist < history%fp_precision .AND. dist < best_match) THEN
     278            8 :                best_match = dist
     279            8 :                found = .TRUE.
     280            8 :                IF (PRESENT(id)) id = history%entries(i)%id
     281              :                IF (debug) found_i = i
     282              :             END IF
     283              :          END DO
     284              : 
     285              :          IF (debug) CALL verify_history_lookup(history, fingerprint, found_i)
     286              :       END IF
     287              : 
     288           25 :       CALL timestop(handle)
     289              : 
     290           25 :    END SUBROUTINE history_lookup
     291              : 
     292              : ! **************************************************************************************************
     293              : !> \brief Helper routine for history_lookup
     294              : !> \param history ...
     295              : !> \param Efind ...
     296              : !> \return ...
     297              : !> \author Ole Schuett
     298              : ! **************************************************************************************************
     299           39 :    FUNCTION interpolation_search(history, Efind) RESULT(res)
     300              :       TYPE(history_type), INTENT(IN)                     :: history
     301              :       REAL(KIND=dp), INTENT(IN)                          :: Efind
     302              :       INTEGER                                            :: res
     303              : 
     304              :       INTEGER                                            :: high, low, mid
     305              :       REAL(KIND=dp)                                      :: slope
     306              : 
     307           39 :       low = 1
     308           39 :       high = history%length
     309              : 
     310           84 :       DO WHILE (low < high)
     311              :          !linear interpolation
     312           45 :          slope = REAL(high - low, KIND=dp)/(history%entries(high)%p%Epot - history%entries(low)%p%Epot)
     313           45 :          mid = low + INT(slope*(Efind - history%entries(low)%p%Epot))
     314           45 :          mid = MIN(MAX(mid, low), high)
     315              : 
     316           84 :          IF (history%entries(mid)%p%Epot < Efind) THEN
     317           21 :             low = mid + 1
     318              :          ELSE
     319           24 :             high = mid - 1
     320              :          END IF
     321              :       END DO
     322              : 
     323           39 :       IF (0 < low .AND. low <= history%length) THEN
     324           36 :          IF (Efind > history%entries(low)%p%Epot) low = low + 1
     325              :       END IF
     326              : 
     327           39 :       res = low
     328           39 :    END FUNCTION interpolation_search
     329              : 
     330              : ! **************************************************************************************************
     331              : !> \brief Debugging routine, performs a slow (but robust) linear search.
     332              : !> \param history ...
     333              : !> \param fingerprint ...
     334              : !> \param found_i_ref ...
     335              : !> \author Ole Schuett
     336              : ! **************************************************************************************************
     337            0 :    SUBROUTINE verify_history_lookup(history, fingerprint, found_i_ref)
     338              :       TYPE(history_type), INTENT(IN)                     :: history
     339              :       TYPE(history_fingerprint_type), INTENT(IN)         :: fingerprint
     340              :       INTEGER, INTENT(IN)                                :: found_i_ref
     341              : 
     342              :       INTEGER                                            :: found_i, i
     343              :       REAL(KIND=dp)                                      :: best_fp_match, Epot_dist, fp_dist
     344              : 
     345            0 :       found_i = -1
     346            0 :       best_fp_match = HUGE(1.0_dp)
     347              : 
     348            0 :       DO i = 1, history%length
     349            0 :          Epot_dist = ABS(fingerprint%Epot - history%entries(i)%p%Epot)
     350            0 :          IF (Epot_dist > history%E_precision) CYCLE
     351            0 :          fp_dist = fingerprint_distance(fingerprint, history%entries(i)%p)
     352              :          !WRITE(*,*) "entry ", i, " dist: ",dist
     353            0 :          IF (fp_dist < history%fp_precision .AND. fp_dist < best_fp_match) THEN
     354            0 :             best_fp_match = fp_dist
     355            0 :             found_i = i
     356              :          END IF
     357              :       END DO
     358              : 
     359            0 :       IF (found_i /= found_i_ref) THEN
     360            0 :          WRITE (*, *) found_i, found_i_ref
     361            0 :          CPABORT("verify_history_lookup failed")
     362              :       END IF
     363              : 
     364            0 :    END SUBROUTINE verify_history_lookup
     365              : 
     366              : ! **************************************************************************************************
     367              : !> \brief Finalizes a history.
     368              : !> \param history ...
     369              : !> \author Ole Schuett
     370              : ! **************************************************************************************************
     371            3 :    SUBROUTINE history_finalize(history)
     372              :       TYPE(history_type)                                 :: history
     373              : 
     374              :       INTEGER                                            :: i
     375              : 
     376           20 :       DO i = 1, history%length
     377           17 :          IF (ASSOCIATED(history%entries(i)%p)) &
     378           20 :             DEALLOCATE (history%entries(i)%p)
     379              :       END DO
     380              : 
     381            3 :       DEALLOCATE (history%entries)
     382              : 
     383            3 :    END SUBROUTINE history_finalize
     384              : 
     385            0 : END MODULE glbopt_history
        

Generated by: LCOV version 2.0-1