-
Notifications
You must be signed in to change notification settings - Fork 1
/
pme_tools.F
171 lines (148 loc) · 6.74 KB
/
pme_tools.F
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
!--------------------------------------------------------------------------------------------------!
! CP2K: A general program to perform molecular dynamics simulations !
! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
! !
! SPDX-License-Identifier: GPL-2.0-or-later !
!--------------------------------------------------------------------------------------------------!
! **************************************************************************************************
!> \brief Tools common both to PME and SPME
!> \par History
!> JGH (03-May-2001) : first correctly working version
!> teo (Feb-2007) : Merging common routines to spme and pme
!> \author JGH (21-Mar-2001)
! **************************************************************************************************
MODULE pme_tools
USE atomic_kind_types, ONLY: atomic_kind_type,&
get_atomic_kind
USE cell_types, ONLY: cell_type,&
real_to_scaled
USE kinds, ONLY: dp
USE particle_types, ONLY: particle_type
USE realspace_grid_types, ONLY: realspace_grid_type
#include "./base/base_uses.f90"
IMPLICIT NONE
PRIVATE
PUBLIC :: get_center, set_list
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pme_tools'
CONTAINS
! **************************************************************************************************
!> \brief ...
!> \param part ...
!> \param npart ...
!> \param center ...
!> \param p1 ...
!> \param rs ...
!> \param ipart ...
!> \param core_center ...
! **************************************************************************************************
SUBROUTINE set_list(part, npart, center, p1, rs, ipart, core_center)
TYPE(particle_type), DIMENSION(:), INTENT(IN) :: part
INTEGER, INTENT(IN) :: npart
INTEGER, DIMENSION(:, :), INTENT(IN) :: center
INTEGER, INTENT(OUT) :: p1
TYPE(realspace_grid_type), INTENT(IN) :: rs
INTEGER, INTENT(INOUT) :: ipart
INTEGER, DIMENSION(:, :), OPTIONAL, POINTER :: core_center
INTEGER :: ndim, npos
INTEGER, DIMENSION(3) :: lb, ub
REAL(KIND=dp) :: charge
TYPE(atomic_kind_type), POINTER :: atomic_kind
p1 = 0
lb = rs%lb_real
ub = rs%ub_real
DO
ipart = ipart + 1
IF (ipart > npart) EXIT
atomic_kind => part(ipart)%atomic_kind
CALL get_atomic_kind(atomic_kind=atomic_kind, qeff=charge)
IF (charge == 0.0_dp .AND. part(ipart)%shell_index == 0) CYCLE
IF (rs%desc%parallel) THEN
! check if the rs grid is distributed or not
IF (ALL(rs%desc%group_dim == 1)) THEN
ndim = rs%desc%group_size
npos = rs%desc%my_pos
! All processors work on the same grid
IF (MOD(ipart, ndim) == npos) THEN
p1 = ipart
EXIT
END IF
ELSE
! First check if this atom is on my grid
IF (part(ipart)%shell_index /= 0 .AND. PRESENT(core_center)) THEN
IF (in_slice(core_center(:, part(ipart)%shell_index), lb, ub)) THEN
p1 = ipart
END IF
ELSE
IF (in_slice(center(:, ipart), lb, ub)) THEN
p1 = ipart
EXIT
END IF
END IF
END IF
ELSE
p1 = ipart
EXIT
END IF
END DO
END SUBROUTINE set_list
! **************************************************************************************************
!> \brief ...
!> \param pos ...
!> \param lb ...
!> \param ub ...
!> \return ...
! **************************************************************************************************
FUNCTION in_slice(pos, lb, ub) RESULT(internal)
INTEGER, DIMENSION(3), INTENT(IN) :: pos, lb, ub
LOGICAL :: internal
IF (ALL(pos >= lb) .AND. ALL(pos <= ub)) THEN
internal = .TRUE.
ELSE
internal = .FALSE.
END IF
END FUNCTION in_slice
! **************************************************************************************************
!> \brief ...
!> \param part ...
!> \param box ...
!> \param center ...
!> \param delta ...
!> \param npts ...
!> \param n ...
! **************************************************************************************************
SUBROUTINE get_center(part, box, center, delta, npts, n)
TYPE(particle_type), DIMENSION(:), INTENT(IN) :: part
TYPE(cell_type), POINTER :: box
INTEGER, DIMENSION(:, :), INTENT(OUT) :: center
REAL(KIND=dp), DIMENSION(:, :), INTENT(OUT) :: delta
INTEGER, DIMENSION(:), INTENT(IN) :: npts
INTEGER, INTENT(IN) :: n
INTEGER :: ipart, mp
REAL(KIND=dp) :: rmp
REAL(KIND=dp), DIMENSION(3) :: ca, gp, s
! The pbc algorithm is sensitive to numeric noise and compiler optimization because of ANINT.
! Therefore center and delta have to be computed simultaneously to ensure they are consistent.
mp = MAXVAL(npts(:))
rmp = REAL(mp, KIND=dp)
DO ipart = 1, SIZE(part)
! compute the scaled coordinate of atom ipart
CALL real_to_scaled(s, part(ipart)%r, box)
s = s - ANINT(s)
! find the continuous ``grid'' point
gp = REAL(npts, KIND=dp)*s
! find the closest grid point (on big grid)
IF (MOD(n, 2) == 0) THEN
center(:, ipart) = INT(gp + rmp) - mp
ca(:) = REAL(center(:, ipart), KIND=dp) + 0.5_dp
ELSE
center(:, ipart) = NINT(gp)
ca(:) = REAL(center(:, ipart), KIND=dp)
END IF
center(:, ipart) = center(:, ipart) + npts(:)/2
center(:, ipart) = MODULO(center(:, ipart), npts(:))
center(:, ipart) = center(:, ipart) - npts(:)/2
! find the distance vector
delta(:, ipart) = gp - ca(:)
END DO
END SUBROUTINE get_center
END MODULE pme_tools