ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/neighborLists.F90
Revision: 872
Committed: Fri Nov 21 19:31:05 2003 UTC (20 years, 7 months ago) by chrisfen
File size: 6926 byte(s)
Log Message:
Fixed a bug in SimInfo ordering of radii

File Contents

# Content
1 !! Module neighborLists
2 !! Impliments verlet neighbor lists for force modules.
3 !! Automagically expands neighbor list if size too small
4 !! up to maxAllocations times. If after maxAllocations we try to
5 !! expand the neighbor list, we get an error message and quit.
6 !! @author Charles F. Vardeman II
7 !! @author Matthew Meineke
8 !! @author J. Daniel Gezelter
9 !! @version $Id: neighborLists.F90,v 1.9 2003-11-21 19:31:05 chrisfen Exp $,
10
11 module neighborLists
12
13 use definitions
14 #ifdef IS_MPI
15 use mpiSimulation
16 #endif
17
18 implicit none
19 PRIVATE
20
21 !--------------MODULE VARIABLES---------------------->
22 !! Parameter for size > # of long range particles neighbor list
23 !! should be.
24 integer, parameter :: listMultiplier = 80
25 !! Maximum number of times we should reallocate neighbor list.
26 integer, parameter :: maxAllocations = 5
27 !! Number of times we have allocated the neighbor list.
28 integer, save :: nAllocations = 0
29 !! Pointer array to location in list for atom i.
30 integer, dimension(:),public, pointer :: point => null()
31 !! Neighbor list for atom i.
32 integer, dimension(:),public, pointer :: list => null()
33 !! Position array of previous positions for check. Allocated first time
34 !! into saveNeighborList.
35 real( kind = dp ), dimension(:,:), allocatable, save :: q0
36 !! Current list size
37 integer, save :: listSize
38 !--------------MODULE ACCESS-------------------------->
39 public :: expandNeighborList
40 public :: checkNeighborList
41 public :: saveNeighborList
42 public :: getNeighborListSize
43
44 contains
45
46
47 subroutine expandNeighborList(natoms, error)
48 integer, intent(out) :: error
49 integer :: natoms
50 integer :: alloc_error
51 integer :: oldSize = 0
52 integer :: newSize = 0
53 integer :: i
54 integer, dimension(:), pointer :: newList => null()
55 error = 0
56
57 !! First time through we should allocate point and list.
58 !! If one is associated and one is not, something is wrong
59 !! and return a error.
60
61 #ifdef IS_MPI !! // MPI
62 if (.not. associated(point) .and. &
63 .not. associated(list) ) then
64 allocate(point(getNRow(plan_row)),stat=alloc_error)
65 if (alloc_error /= 0) then
66 write(default_error,*) &
67 "ExpandNeighborLists: Error in allocating MPI point"
68 error = -1
69 return
70 end if
71 allocate(list(listMultiplier * getNCol(plan_col)),stat=alloc_error)
72 if (alloc_error /= 0) then
73 write(default_error,*) &
74 "ExpandNeighborLists: Error in allocating MPI list"
75 error = -1
76 return
77 end if
78 listSize = size(list)
79 nAllocations = nAllocations + 1
80 return
81 end if
82 #else !! // NONMPI
83 if (.not. associated(point) .and. &
84 .not. associated(list) ) then
85 allocate(point(natoms),stat=alloc_error)
86 if (alloc_error /= 0) then
87 write(default_error,*) &
88 "ExpandNeighborLists: Error in allocating point"
89 error = -1
90 return
91 end if
92 allocate(list(listMultiplier * natoms),stat=alloc_error)
93 if (alloc_error /= 0) then
94 write(default_error,*) &
95 "ExpandNeighborLists: Error in allocating list"
96 error = -1
97 return
98 end if
99 listSize = size(list)
100 nAllocations = nAllocations + 1
101 return
102 end if
103 #endif !! //MPI
104
105 ! Expand the neighbor list
106
107 ! Check to see if we have exceeded the maximum number of allocations.
108 if (nAllocations > maxAllocations) then
109 write(default_error,*) &
110 "ExpandNeighborList: exceeded maximum number of re-allocations"
111 error = -1
112 return
113 else !! Expand the list.
114 oldSize = size(list)
115
116
117 #ifdef IS_MPI !! MPI
118 newSize = listMultiplier * getNCol(plan_col) + oldSize
119 #else
120 newSize = listMultiplier * natoms + oldSize
121 #endif !! MPI
122
123
124
125 allocate(newList(newSize), stat=alloc_error)
126 if (alloc_error /= 0) then
127 write(*,*) "Error allocating new neighborlist"
128 error = -1
129 return
130 end if
131
132 !! Copy old list to new list
133 do i = 1, oldSize
134 newList(i) = list(i)
135 end do
136 !! Free old list
137 if(associated(list)) deallocate(list,stat=alloc_error)
138 if (alloc_error /= 0) then
139 error = -1
140 return
141 end if
142
143 !! Point list at new list
144
145 list => newList
146 end if
147
148 nAllocations = nAllocations + 1
149 listSize = size(list)
150 end subroutine expandNeighborList
151
152 !! checks to see if any long range particle has moved
153 !! through the neighbor list skin thickness.
154 subroutine checkNeighborList(natoms, q, listSkin, update_nlist)
155 integer :: natoms
156 real(kind = dp), intent(in) :: listSkin
157 real( kind = dp ), dimension(:,:) :: q
158 integer :: i
159 real( kind = DP ) :: dispmx
160 logical, intent(out) :: update_nlist
161 real( kind = DP ) :: dispmx_tmp
162 integer :: nlocal
163
164 nlocal = natoms
165 dispmx = 0.0E0_DP
166 !! calculate the largest displacement of any atom in any direction
167
168 !! If we have changed the particle idents, then we need to update
169 if (.not. allocated(q0)) then
170 update_nlist = .true.
171 return
172 end if
173
174 if (size(q0,2) /= nlocal) then
175 update_nlist = .true.
176 return
177 end if
178
179
180 #ifdef MPI
181
182 dispmx_tmp = 0.0E0_DP
183 do i = 1, nlocal
184 dispmx_tmp = max( abs ( q(1,i) - q0(1,i) ), dispmx_tmp )
185 dispmx_tmp = max( abs ( q(2,i) - q0(2,i) ), dispmx_tmp )
186 dispmx_tmp = max( abs ( q(3,i) - q0(3,i) ), dispmx_tmp )
187 end do
188 call mpi_allreduce(dispmx_tmp,dispmx,1,mpi_double_precision, &
189 mpi_max,mpi_comm_world,mpi_err)
190
191 #else
192
193 dispmx = 0.0_DP
194 do i = 1, nlocal
195 dispmx = max( abs ( q(1,i) - q0(1,i) ), dispmx )
196 dispmx = max( abs ( q(2,i) - q0(2,i) ), dispmx )
197 dispmx = max( abs ( q(3,i) - q0(3,i) ), dispmx )
198 end do
199
200 #endif
201
202 !! a conservative test of list skin crossings
203 dispmx = 2.0E0_DP * sqrt (3.0E0_DP * dispmx * dispmx)
204
205 update_nlist = (dispmx.gt.listSkin)
206
207 end subroutine checkNeighborList
208
209
210 !! Saves neighbor list for comparison in check.
211 !! Save_neighborList will work even if the number of
212 !! local atoms has changed.
213 subroutine saveNeighborList(natoms, q)
214
215 integer :: natoms
216 real(kind = dp ), dimension(3,natoms), intent(in) :: q
217 integer :: list_size
218
219 !! get size of list
220 list_size = natoms
221
222 if (.not. allocated(q0)) then
223 allocate(q0(3,list_size))
224 else if( list_size > size(q0,2)) then
225 deallocate(q0)
226 allocate(q0(3,list_size))
227 endif
228 q0 = q
229 end subroutine saveNeighborList
230
231
232 function getNeighborListSize() result(returnListSize)
233 integer :: returnListSize
234 returnListSize = listSize
235 end function getNeighborListSize
236
237 end module neighborLists