1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2019 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5! ************************************************************************************************** 6MODULE integration_grid_types 7 8 USE kinds, ONLY: dp 9#include "./base/base_uses.f90" 10 11 IMPLICIT NONE 12 13 PRIVATE 14 15 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'integration_grid_types' 16 17 TYPE grid_batch_val_1d_type 18 INTEGER :: np1 19 REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: val1d 20 END TYPE grid_batch_val_1d_type 21 22 TYPE grid_batch_val_2d_type 23 INTEGER :: np1, np2 24 REAL(KIND=dp), DIMENSION(:, :), ALLOCATABLE :: val2d 25 END TYPE grid_batch_val_2d_type 26 27 TYPE gnlist_type 28 INTEGER, DIMENSION(:), ALLOCATABLE :: atom_list 29 REAL(KIND=dp), DIMENSION(:, :), ALLOCATABLE :: atom_pos 30 END TYPE gnlist_type 31 32 TYPE grid_batch_info_type 33 INTEGER :: np 34 INTEGER :: ref_atom 35 INTEGER :: ibatch 36 TYPE(gnlist_type) :: gnlist 37 REAL(KIND=dp), DIMENSION(3) :: rcenter 38 REAL(KIND=dp) :: radius 39 REAL(dp), DIMENSION(:, :), ALLOCATABLE :: rco 40 REAL(dp), DIMENSION(:), ALLOCATABLE :: weight 41 REAL(dp), DIMENSION(:), ALLOCATABLE :: wref 42 REAL(dp), DIMENSION(:), ALLOCATABLE :: wsum 43 END TYPE grid_batch_info_type 44 45 TYPE integration_grid_type 46 INTEGER :: nbatch 47 TYPE(grid_batch_info_type), DIMENSION(:), ALLOCATABLE :: grid_batch 48 END TYPE integration_grid_type 49 50 TYPE integration_grid_value_type 51 INTEGER :: nbatch 52 TYPE(grid_batch_val_1d_type), DIMENSION(:), ALLOCATABLE :: grid_val_1d 53 TYPE(grid_batch_val_2d_type), DIMENSION(:), ALLOCATABLE :: grid_val_2d 54 END TYPE integration_grid_value_type 55 56 PUBLIC :: integration_grid_type, allocate_intgrid, deallocate_intgrid 57 PUBLIC :: integration_grid_value_type, allocate_intgrid_val, deallocate_intgrid_val 58 59! ************************************************************************************************** 60 61CONTAINS 62 63! ************************************************************************************************** 64!> \brief Initialize integration_grid_type 65!> \param int_grid ... 66!> \date 02.2018 67!> \param 68!> \author JGH 69!> \version 1.0 70! ************************************************************************************************** 71 SUBROUTINE allocate_intgrid(int_grid) 72 73 TYPE(integration_grid_type), POINTER :: int_grid 74 75 CHARACTER(len=*), PARAMETER :: routineN = 'allocate_intgrid', & 76 routineP = moduleN//':'//routineN 77 78 IF (ASSOCIATED(int_grid)) CALL deallocate_intgrid(int_grid) 79 ALLOCATE (int_grid) 80 int_grid%nbatch = 0 81 82 END SUBROUTINE allocate_intgrid 83 84! ************************************************************************************************** 85!> \brief Deallocate integration_grid_type 86!> \param int_grid ... 87!> \date 02.2018 88!> \param 89!> \author JGH 90!> \version 1.0 91! ************************************************************************************************** 92 SUBROUTINE deallocate_intgrid(int_grid) 93 TYPE(integration_grid_type), POINTER :: int_grid 94 95 CHARACTER(len=*), PARAMETER :: routineN = 'deallocate_intgrid', & 96 routineP = moduleN//':'//routineN 97 98 INTEGER :: i 99 100 IF (ASSOCIATED(int_grid)) THEN 101 IF (ALLOCATED(int_grid%grid_batch)) THEN 102 DO i = 1, int_grid%nbatch 103 IF (ALLOCATED(int_grid%grid_batch(i)%rco)) DEALLOCATE (int_grid%grid_batch(i)%rco) 104 IF (ALLOCATED(int_grid%grid_batch(i)%weight)) DEALLOCATE (int_grid%grid_batch(i)%weight) 105 IF (ALLOCATED(int_grid%grid_batch(i)%wref)) DEALLOCATE (int_grid%grid_batch(i)%wref) 106 IF (ALLOCATED(int_grid%grid_batch(i)%wsum)) DEALLOCATE (int_grid%grid_batch(i)%wsum) 107 ! 108 IF (ALLOCATED(int_grid%grid_batch(i)%gnlist%atom_list)) DEALLOCATE (int_grid%grid_batch(i)%gnlist%atom_list) 109 IF (ALLOCATED(int_grid%grid_batch(i)%gnlist%atom_pos)) DEALLOCATE (int_grid%grid_batch(i)%gnlist%atom_pos) 110 END DO 111 DEALLOCATE (int_grid%grid_batch) 112 END IF 113 DEALLOCATE (int_grid) 114 ELSE 115 CALL cp_abort(__LOCATION__, & 116 "The pointer int_grid is not associated and "// & 117 "cannot be deallocated") 118 END IF 119 END SUBROUTINE deallocate_intgrid 120 121! ************************************************************************************************** 122!> \brief Initialize integration_grid_value_type 123!> \param int_grid ... 124!> \date 02.2018 125!> \param 126!> \author JGH 127!> \version 1.0 128! ************************************************************************************************** 129 SUBROUTINE allocate_intgrid_val(int_grid) 130 131 TYPE(integration_grid_value_type), POINTER :: int_grid 132 133 CHARACTER(len=*), PARAMETER :: routineN = 'allocate_intgrid_val', & 134 routineP = moduleN//':'//routineN 135 136 IF (ASSOCIATED(int_grid)) CALL deallocate_intgrid_val(int_grid) 137 ALLOCATE (int_grid) 138 int_grid%nbatch = 0 139 140 END SUBROUTINE allocate_intgrid_val 141 142! ************************************************************************************************** 143!> \brief Deallocate integration_grid_value_type 144!> \param int_grid ... 145!> \date 02.2018 146!> \param 147!> \author JGH 148!> \version 1.0 149! ************************************************************************************************** 150 SUBROUTINE deallocate_intgrid_val(int_grid) 151 TYPE(integration_grid_value_type), POINTER :: int_grid 152 153 CHARACTER(len=*), PARAMETER :: routineN = 'deallocate_intgrid_val', & 154 routineP = moduleN//':'//routineN 155 156 INTEGER :: i 157 158 IF (ASSOCIATED(int_grid)) THEN 159 IF (ALLOCATED(int_grid%grid_val_1d)) THEN 160 DO i = 1, int_grid%nbatch 161 IF (ALLOCATED(int_grid%grid_val_1d(i)%val1d)) DEALLOCATE (int_grid%grid_val_1d(i)%val1d) 162 END DO 163 DEALLOCATE (int_grid%grid_val_1d) 164 END IF 165 IF (ALLOCATED(int_grid%grid_val_2d)) THEN 166 DO i = 1, int_grid%nbatch 167 IF (ALLOCATED(int_grid%grid_val_2d(i)%val2d)) DEALLOCATE (int_grid%grid_val_2d(i)%val2d) 168 END DO 169 DEALLOCATE (int_grid%grid_val_2d) 170 END IF 171 DEALLOCATE (int_grid) 172 ELSE 173 CALL cp_abort(__LOCATION__, & 174 "The pointer int_grid is not associated and "// & 175 "cannot be deallocated") 176 END IF 177 END SUBROUTINE deallocate_intgrid_val 178 179END MODULE integration_grid_types 180