LCOV - code coverage report
Current view: top level - src/common - cp_min_heap.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:d1f8d1b) Lines: 110 121 90.9 %
Date: 2024-11-29 06:42:44 Functions: 14 19 73.7 %

          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             : MODULE cp_min_heap
       9             :    USE kinds,                           ONLY: int_4,&
      10             :                                               int_8
      11             : #include "../base/base_uses.f90"
      12             : 
      13             :    IMPLICIT NONE
      14             :    PRIVATE
      15             : 
      16             :    PUBLIC :: cp_heap_type, keyt, valt
      17             :    PUBLIC :: cp_heap_pop, cp_heap_reset_node, cp_heap_fill
      18             :    PUBLIC :: cp_heap_new, cp_heap_release
      19             :    PUBLIC :: cp_heap_get_first, cp_heap_reset_first
      20             : 
      21             :    ! Sets the types
      22             :    INTEGER, PARAMETER         :: keyt = int_4
      23             :    INTEGER, PARAMETER         :: valt = int_8
      24             : 
      25             :    TYPE cp_heap_node
      26             :       INTEGER(KIND=keyt) :: key = -1_keyt
      27             :       INTEGER(KIND=valt) :: value = -1_valt
      28             :    END TYPE cp_heap_node
      29             : 
      30             :    TYPE cp_heap_node_e
      31             :       TYPE(cp_heap_node) :: node = cp_heap_node()
      32             :    END TYPE cp_heap_node_e
      33             : 
      34             :    TYPE cp_heap_type
      35             :       INTEGER :: n = -1
      36             :       INTEGER, DIMENSION(:), POINTER           :: index => NULL()
      37             :       TYPE(cp_heap_node_e), DIMENSION(:), POINTER :: nodes => NULL()
      38             :    END TYPE cp_heap_type
      39             : 
      40             : CONTAINS
      41             : 
      42             :    ! Lookup functions
      43             : 
      44             : ! **************************************************************************************************
      45             : !> \brief ...
      46             : !> \param n ...
      47             : !> \return ...
      48             : ! **************************************************************************************************
      49     1295697 :    ELEMENTAL FUNCTION get_parent(n) RESULT(parent)
      50             :       INTEGER, INTENT(IN)                                :: n
      51             :       INTEGER                                            :: parent
      52             : 
      53     1295697 :       parent = INT(n/2)
      54     1295697 :    END FUNCTION get_parent
      55             : 
      56             : ! **************************************************************************************************
      57             : !> \brief ...
      58             : !> \param n ...
      59             : !> \return ...
      60             : ! **************************************************************************************************
      61      602309 :    ELEMENTAL FUNCTION get_left_child(n) RESULT(child)
      62             :       INTEGER, INTENT(IN)                                :: n
      63             :       INTEGER                                            :: child
      64             : 
      65      602309 :       child = 2*n
      66      602309 :    END FUNCTION get_left_child
      67             : 
      68             : ! **************************************************************************************************
      69             : !> \brief ...
      70             : !> \param n ...
      71             : !> \return ...
      72             : ! **************************************************************************************************
      73           0 :    ELEMENTAL FUNCTION get_right_child(n) RESULT(child)
      74             :       INTEGER, INTENT(IN)                                :: n
      75             :       INTEGER                                            :: child
      76             : 
      77           0 :       child = 2*n + 1
      78           0 :    END FUNCTION get_right_child
      79             : 
      80             : ! **************************************************************************************************
      81             : !> \brief ...
      82             : !> \param heap ...
      83             : !> \param n ...
      84             : !> \return ...
      85             : ! **************************************************************************************************
      86     1204657 :    ELEMENTAL FUNCTION get_value(heap, n) RESULT(value)
      87             :       TYPE(cp_heap_type), INTENT(IN)                     :: heap
      88             :       INTEGER, INTENT(IN)                                :: n
      89             :       INTEGER(KIND=valt)                                 :: value
      90             : 
      91     1204657 :       value = heap%nodes(n)%node%value
      92     1204657 :    END FUNCTION get_value
      93             : 
      94             : ! **************************************************************************************************
      95             : !> \brief ...
      96             : !> \param heap ...
      97             : !> \param key ...
      98             : !> \return ...
      99             : ! **************************************************************************************************
     100           0 :    ELEMENTAL FUNCTION get_value_by_key(heap, key) RESULT(value)
     101             :       TYPE(cp_heap_type), INTENT(IN)                     :: heap
     102             :       INTEGER(KIND=keyt), INTENT(IN)                     :: key
     103             :       INTEGER(KIND=valt)                                 :: value
     104             : 
     105             :       INTEGER                                            :: n
     106             : 
     107           0 :       n = heap%index(key)
     108           0 :       value = get_value(heap, n)
     109           0 :    END FUNCTION get_value_by_key
     110             : 
     111             :    ! Initialization functions
     112             : 
     113             : ! **************************************************************************************************
     114             : !> \brief ...
     115             : !> \param heap ...
     116             : !> \param n ...
     117             : ! **************************************************************************************************
     118       28339 :    SUBROUTINE cp_heap_new(heap, n)
     119             :       TYPE(cp_heap_type), INTENT(OUT)                    :: heap
     120             :       INTEGER, INTENT(IN)                                :: n
     121             : 
     122       28339 :       heap%n = n
     123       85017 :       ALLOCATE (heap%index(n))
     124      138928 :       ALLOCATE (heap%nodes(n))
     125       28339 :    END SUBROUTINE cp_heap_new
     126             : 
     127             : ! **************************************************************************************************
     128             : !> \brief ...
     129             : !> \param heap ...
     130             : ! **************************************************************************************************
     131       28339 :    SUBROUTINE cp_heap_release(heap)
     132             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     133             : 
     134       28339 :       DEALLOCATE (heap%index)
     135       28339 :       DEALLOCATE (heap%nodes)
     136       28339 :       heap%n = 0
     137       28339 :    END SUBROUTINE cp_heap_release
     138             : 
     139             : ! **************************************************************************************************
     140             : !> \brief Fill heap with given values
     141             : !> \param heap ...
     142             : !> \param values ...
     143             : ! **************************************************************************************************
     144       28339 :    SUBROUTINE cp_heap_fill(heap, values)
     145             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     146             :       INTEGER(KIND=valt), DIMENSION(:), INTENT(IN)       :: values
     147             : 
     148             :       INTEGER                                            :: first, i, n
     149             : 
     150       28339 :       n = SIZE(values)
     151       28339 :       CPASSERT(heap%n >= n)
     152             : 
     153       82250 :       DO i = 1, n
     154       53911 :          heap%index(i) = i
     155       53911 :          heap%nodes(i)%node%key = i
     156       82250 :          heap%nodes(i)%node%value = values(i)
     157             :       END DO
     158             :       ! Sort from the last full subtree
     159       28339 :       first = get_parent(n)
     160       53900 :       DO i = first, 1, -1
     161       53900 :          CALL bubble_down(heap, i)
     162             :       END DO
     163             : 
     164       28339 :    END SUBROUTINE cp_heap_fill
     165             : 
     166             : ! **************************************************************************************************
     167             : !> \brief Returns the first heap element without removing it.
     168             : !> \param heap ...
     169             : !> \param key ...
     170             : !> \param value ...
     171             : !> \param found ...
     172             : ! **************************************************************************************************
     173      639460 :    SUBROUTINE cp_heap_get_first(heap, key, value, found)
     174             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     175             :       INTEGER(KIND=keyt), INTENT(OUT)                    :: key
     176             :       INTEGER(KIND=valt), INTENT(OUT)                    :: value
     177             :       LOGICAL, INTENT(OUT)                               :: found
     178             : 
     179      639460 :       IF (heap%n .LT. 1) THEN
     180           0 :          found = .FALSE.
     181             :       ELSE
     182      639460 :          found = .TRUE.
     183      639460 :          key = heap%nodes(1)%node%key
     184      639460 :          value = heap%nodes(1)%node%value
     185             :       END IF
     186      639460 :    END SUBROUTINE cp_heap_get_first
     187             : 
     188             : ! **************************************************************************************************
     189             : !> \brief Returns and removes the first heap element and rebalances
     190             : !>        the heap.
     191             : !> \param heap ...
     192             : !> \param key ...
     193             : !> \param value ...
     194             : !> \param found ...
     195             : ! **************************************************************************************************
     196          85 :    SUBROUTINE cp_heap_pop(heap, key, value, found)
     197             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     198             :       INTEGER(KIND=keyt), INTENT(OUT)                    :: key
     199             :       INTEGER(KIND=valt), INTENT(OUT)                    :: value
     200             :       LOGICAL, INTENT(OUT)                               :: found
     201             : 
     202             : !
     203             : 
     204          85 :       CALL cp_heap_get_first(heap, key, value, found)
     205          85 :       IF (found) THEN
     206          85 :          IF (heap%n .GT. 1) THEN
     207          48 :             CALL cp_heap_copy_node(heap, 1, heap%n)
     208          48 :             heap%n = heap%n - 1
     209          48 :             CALL bubble_down(heap, 1)
     210             :          ELSE
     211          37 :             heap%n = heap%n - 1
     212             :          END IF
     213             :       END IF
     214          85 :    END SUBROUTINE cp_heap_pop
     215             : 
     216             : ! **************************************************************************************************
     217             : !> \brief Changes the value of the heap element with given key and
     218             : !>        rebalances the heap.
     219             : !> \param heap ...
     220             : !> \param key ...
     221             : !> \param value ...
     222             : ! **************************************************************************************************
     223         102 :    SUBROUTINE cp_heap_reset_node(heap, key, value)
     224             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     225             :       INTEGER(KIND=keyt), INTENT(IN)                     :: key
     226             :       INTEGER(KIND=valt), INTENT(IN)                     :: value
     227             : 
     228             :       INTEGER                                            :: n, new_pos
     229             : 
     230          51 :       CPASSERT(heap%n > 0)
     231             : 
     232          51 :       n = heap%index(key)
     233          51 :       CPASSERT(heap%nodes(n)%node%key == key)
     234          51 :       heap%nodes(n)%node%value = value
     235          51 :       CALL bubble_up(heap, n, new_pos)
     236          51 :       CALL bubble_down(heap, new_pos)
     237          51 :    END SUBROUTINE cp_heap_reset_node
     238             : 
     239             : ! **************************************************************************************************
     240             : !> \brief Changes the value of the minimum heap element and rebalances the heap.
     241             : !> \param heap ...
     242             : !> \param value ...
     243             : ! **************************************************************************************************
     244      639375 :    SUBROUTINE cp_heap_reset_first(heap, value)
     245             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     246             :       INTEGER(KIND=valt), INTENT(IN)                     :: value
     247             : 
     248      639375 :       CPASSERT(heap%n > 0)
     249      639375 :       heap%nodes(1)%node%value = value
     250      639375 :       CALL bubble_down(heap, 1)
     251      639375 :    END SUBROUTINE cp_heap_reset_first
     252             : 
     253             : ! **************************************************************************************************
     254             : !> \brief ...
     255             : !> \param heap ...
     256             : !> \param e1 ...
     257             : !> \param e2 ...
     258             : ! **************************************************************************************************
     259      602316 :    PURE SUBROUTINE cp_heap_swap(heap, e1, e2)
     260             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     261             :       INTEGER, INTENT(IN)                                :: e1, e2
     262             : 
     263             :       INTEGER(KIND=keyt)                                 :: key1, key2
     264             :       TYPE(cp_heap_node)                                 :: tmp_node
     265             : 
     266      602316 :       key1 = heap%nodes(e1)%node%key
     267      602316 :       key2 = heap%nodes(e2)%node%key
     268             : 
     269      602316 :       tmp_node = heap%nodes(e1)%node
     270      602316 :       heap%nodes(e1)%node = heap%nodes(e2)%node
     271      602316 :       heap%nodes(e2)%node = tmp_node
     272             : 
     273      602316 :       heap%index(key1) = e2
     274      602316 :       heap%index(key2) = e1
     275      602316 :    END SUBROUTINE cp_heap_swap
     276             : 
     277             : ! **************************************************************************************************
     278             : !> \brief Sets node e1 to e2
     279             : !> \param heap ...
     280             : !> \param e1 ...
     281             : !> \param e2 ...
     282             : ! **************************************************************************************************
     283          48 :    PURE SUBROUTINE cp_heap_copy_node(heap, e1, e2)
     284             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     285             :       INTEGER, INTENT(IN)                                :: e1, e2
     286             : 
     287             :       INTEGER(KIND=keyt)                                 :: key1, key2
     288             : 
     289          48 :       key1 = heap%nodes(e1)%node%key
     290          48 :       key2 = heap%nodes(e2)%node%key
     291             : 
     292          48 :       heap%nodes(e1)%node = heap%nodes(e2)%node
     293             : 
     294          48 :       heap%index(key1) = 0
     295          48 :       heap%index(key2) = e1
     296          48 :    END SUBROUTINE cp_heap_copy_node
     297             : 
     298             : ! **************************************************************************************************
     299             : !> \brief Balances a heap by bubbling down from the given element.
     300             : !> \param heap ...
     301             : !> \param first ...
     302             : ! **************************************************************************************************
     303      665035 :    SUBROUTINE bubble_down(heap, first)
     304             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     305             :       INTEGER, INTENT(IN)                                :: first
     306             : 
     307             :       INTEGER                                            :: e, left_child, right_child, smallest
     308             :       INTEGER(kind=valt)                                 :: left_child_value, min_value, &
     309             :                                                             right_child_value
     310             :       LOGICAL                                            :: all_done
     311             : 
     312             : !
     313      665035 :       CPASSERT(0 < first .AND. first <= heap%n)
     314             : 
     315      665035 :       e = first
     316      665035 :       all_done = .FALSE.
     317             :       ! Check whether we are finished, i.e,. whether the element to
     318             :       ! bubble down is childless.
     319     1267344 :       DO WHILE (e .LE. get_parent(heap%n) .AND. .NOT. all_done)
     320             :          ! Determines which node (current, left, or right child) has the
     321             :          ! smallest value.
     322      602309 :          smallest = e
     323      602309 :          min_value = get_value(heap, e)
     324      602309 :          left_child = get_left_child(e)
     325      602309 :          IF (left_child .LE. heap%n) THEN
     326      602309 :             left_child_value = get_value(heap, left_child)
     327      602309 :             IF (left_child_value .LT. min_value) THEN
     328      359346 :                min_value = left_child_value
     329      359346 :                smallest = left_child
     330             :             END IF
     331             :          END IF
     332      602309 :          right_child = left_child + 1
     333      602309 :          IF (right_child .LE. heap%n) THEN
     334          11 :             right_child_value = get_value(heap, right_child)
     335          11 :             IF (right_child_value .LT. min_value) THEN
     336           0 :                min_value = right_child_value
     337           0 :                smallest = right_child
     338             :             END IF
     339             :          END IF
     340             :          !
     341      602309 :          CALL cp_heap_swap(heap, e, smallest)
     342     1267344 :          IF (smallest .EQ. e) THEN
     343             :             all_done = .TRUE.
     344             :          ELSE
     345      359346 :             e = smallest
     346             :          END IF
     347             :       END DO
     348      665035 :    END SUBROUTINE bubble_down
     349             : 
     350             : ! **************************************************************************************************
     351             : !> \brief Balances a heap by bubbling up from the given element.
     352             : !> \param heap ...
     353             : !> \param first ...
     354             : !> \param new_pos ...
     355             : ! **************************************************************************************************
     356          51 :    SUBROUTINE bubble_up(heap, first, new_pos)
     357             :       TYPE(cp_heap_type), INTENT(INOUT)                  :: heap
     358             :       INTEGER, INTENT(IN)                                :: first
     359             :       INTEGER, INTENT(OUT)                               :: new_pos
     360             : 
     361             :       INTEGER                                            :: e, parent
     362             :       INTEGER(kind=valt)                                 :: my_value, parent_value
     363             :       LOGICAL                                            :: all_done
     364             : 
     365          51 :       CPASSERT(0 < first .AND. first <= heap%n)
     366             : 
     367          51 :       e = first
     368          51 :       all_done = .FALSE.
     369          51 :       IF (e .GT. 1) THEN
     370          14 :          my_value = get_value(heap, e)
     371             :       END IF
     372             :       ! Check whether we are finished, i.e,. whether the element to
     373             :       ! bubble up is an orphan.
     374          51 :       new_pos = e
     375          65 :       DO WHILE (e .GT. 1 .AND. .NOT. all_done)
     376             :          ! Switches the parent and the current element if the current
     377             :          ! element's value is greater than the parent's value.
     378          14 :          parent = get_parent(e)
     379          14 :          parent_value = get_value(heap, parent)
     380          65 :          IF (my_value .LT. parent_value) THEN
     381           7 :             CALL cp_heap_swap(heap, e, parent)
     382           7 :             e = parent
     383             :          ELSE
     384             :             all_done = .TRUE.
     385             :          END IF
     386             :       END DO
     387          51 :       new_pos = e
     388          51 :    END SUBROUTINE bubble_up
     389             : 
     390           0 : END MODULE cp_min_heap

Generated by: LCOV version 1.15