ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/simulation_module.F90
Revision: 1203
Committed: Thu May 27 18:59:17 2004 UTC (20 years, 1 month ago) by gezelter
File size: 14978 byte(s)
Log Message:
Cutoff group changes under MPI

File Contents

# Content
1 !! Fortran interface to C entry plug.
2
3 module simulation
4 use definitions
5 use neighborLists
6 use force_globals
7 use vector_class
8 use atype_module
9 use switcheroo
10 #ifdef IS_MPI
11 use mpiSimulation
12 #endif
13
14 implicit none
15 PRIVATE
16
17 #define __FORTRAN90
18 #include "fSimulation.h"
19 #include "fSwitchingFunction.h"
20
21 type (simtype), public, save :: thisSim
22
23 logical, save :: simulation_setup_complete = .false.
24
25 integer, public, save :: nLocal, nGlobal
26 integer, public, save :: nGroups, nGroupGlobal
27 integer, public, save :: nExcludes_Global = 0
28 integer, public, save :: nExcludes_Local = 0
29 integer, allocatable, dimension(:,:), public :: excludesLocal
30 integer, allocatable, dimension(:), public :: excludesGlobal
31 integer, allocatable, dimension(:), public :: molMembershipList
32 integer, allocatable, dimension(:), public :: groupListRow
33 integer, allocatable, dimension(:), public :: groupStartRow
34 integer, allocatable, dimension(:), public :: groupListCol
35 integer, allocatable, dimension(:), public :: groupStartCol
36 integer, allocatable, dimension(:), public :: groupListLocal
37 integer, allocatable, dimension(:), public :: groupStartLocal
38 integer, allocatable, dimension(:), public :: nSkipsForAtom
39 integer, allocatable, dimension(:,:), public :: skipsForAtom
40 real(kind=dp), allocatable, dimension(:), public :: mfactRow
41 real(kind=dp), allocatable, dimension(:), public :: mfactCol
42 real(kind=dp), allocatable, dimension(:), public :: mfactLocal
43
44 real(kind=dp), public, dimension(3,3), save :: Hmat, HmatInv
45 logical, public, save :: boxIsOrthorhombic
46
47 public :: SimulationSetup
48 public :: getNlocal
49 public :: setBox
50 public :: getDielect
51 public :: SimUsesPBC
52 public :: SimUsesLJ
53 public :: SimUsesCharges
54 public :: SimUsesDipoles
55 public :: SimUsesSticky
56 public :: SimUsesRF
57 public :: SimUsesGB
58 public :: SimUsesEAM
59 public :: SimRequiresPrepairCalc
60 public :: SimRequiresPostpairCalc
61 public :: SimUsesDirectionalAtoms
62
63 contains
64
65 subroutine SimulationSetup(setThisSim, CnGlobal, CnLocal, c_idents, &
66 CnLocalExcludes, CexcludesLocal, CnGlobalExcludes, CexcludesGlobal, &
67 CmolMembership, Cmfact, CnGroups, CgroupList, CgroupStart, &
68 status)
69
70 type (simtype) :: setThisSim
71 integer, intent(inout) :: CnGlobal, CnLocal
72 integer, dimension(CnLocal),intent(inout) :: c_idents
73
74 integer :: CnLocalExcludes
75 integer, dimension(2,CnLocalExcludes), intent(in) :: CexcludesLocal
76 integer :: CnGlobalExcludes
77 integer, dimension(CnGlobalExcludes), intent(in) :: CexcludesGlobal
78 integer, dimension(CnGlobal),intent(in) :: CmolMembership
79 !! Result status, success = 0, status = -1
80 integer, intent(out) :: status
81 integer :: i, j, me, thisStat, alloc_stat, myNode, id1, id2
82 integer :: ia
83
84 !! mass factors used for molecular cutoffs
85 real ( kind = dp ), dimension(CnLocal) :: Cmfact
86 integer, intent(in):: CnGroups
87 integer, dimension(CnLocal),intent(in) :: CgroupList
88 integer, dimension(CnGroups),intent(in) :: CgroupStart
89 integer :: maxSkipsForAtom
90
91 #ifdef IS_MPI
92 integer, allocatable, dimension(:) :: c_idents_Row
93 integer, allocatable, dimension(:) :: c_idents_Col
94 integer :: nAtomsInRow, nGroupsInRow
95 integer :: nAtomsInCol, nGroupsInCol
96 #endif
97
98 simulation_setup_complete = .false.
99 status = 0
100
101 ! copy C struct into fortran type
102
103 nLocal = CnLocal
104 nGlobal = CnGlobal
105 nGroups = CnGroups
106
107 thisSim = setThisSim
108
109 nExcludes_Global = CnGlobalExcludes
110 nExcludes_Local = CnLocalExcludes
111
112 call InitializeForceGlobals(nLocal, thisStat)
113 if (thisStat /= 0) then
114 write(default_error,*) "SimSetup: InitializeForceGlobals error"
115 status = -1
116 return
117 endif
118
119 call InitializeSimGlobals(thisStat)
120 if (thisStat /= 0) then
121 write(default_error,*) "SimSetup: InitializeSimGlobals error"
122 status = -1
123 return
124 endif
125
126 #ifdef IS_MPI
127 ! We can only set up forces if mpiSimulation has been setup.
128 if (.not. isMPISimSet()) then
129 write(default_error,*) "MPI is not set"
130 status = -1
131 return
132 endif
133 nAtomsInRow = getNatomsInRow(plan_atom_row)
134 nAtomsInCol = getNatomsInCol(plan_atom_col)
135 nGroupsInRow = getNgroupsInRow(plan_group_row)
136 nGroupsInCol = getNgroupsInCol(plan_group_col)
137 mynode = getMyNode()
138
139 allocate(c_idents_Row(nAtomsInRow),stat=alloc_stat)
140 if (alloc_stat /= 0 ) then
141 status = -1
142 return
143 endif
144
145 allocate(c_idents_Col(nAtomsInCol),stat=alloc_stat)
146 if (alloc_stat /= 0 ) then
147 status = -1
148 return
149 endif
150
151 call gather(c_idents, c_idents_Row, plan_atom_row)
152 call gather(c_idents, c_idents_Col, plan_atom_col)
153
154 do i = 1, nAtomsInRow
155 me = getFirstMatchingElement(atypes, "c_ident", c_idents_Row(i))
156 atid_Row(i) = me
157 enddo
158
159 do i = 1, nAtomsInCol
160 me = getFirstMatchingElement(atypes, "c_ident", c_idents_Col(i))
161 atid_Col(i) = me
162 enddo
163
164 !! free temporary ident arrays
165 if (allocated(c_idents_Col)) then
166 deallocate(c_idents_Col)
167 end if
168 if (allocated(c_idents_Row)) then
169 deallocate(c_idents_Row)
170 endif
171
172 #endif
173
174 #ifdef IS_MPI
175 allocate(groupStartRow(nGroupsInRow+1),stat=alloc_stat)
176 if (alloc_stat /= 0 ) then
177 status = -1
178 return
179 endif
180 allocate(groupStartCol(nGroupsInCol+1),stat=alloc_stat)
181 if (alloc_stat /= 0 ) then
182 status = -1
183 return
184 endif
185 allocate(groupListRow(nAtomsInRow),stat=alloc_stat)
186 if (alloc_stat /= 0 ) then
187 status = -1
188 return
189 endif
190 allocate(groupListCol(nAtomsInCol),stat=alloc_stat)
191 if (alloc_stat /= 0 ) then
192 status = -1
193 return
194 endif
195 allocate(mfactRow(nAtomsInRow),stat=alloc_stat)
196 if (alloc_stat /= 0 ) then
197 status = -1
198 return
199 endif
200 allocate(mfactCol(nAtomsInCol),stat=alloc_stat)
201 if (alloc_stat /= 0 ) then
202 status = -1
203 return
204 endif
205 allocate(groupStartLocal(nGroups), stat=alloc_stat)
206 if (alloc_stat /= 0 ) then
207 status = -1
208 return
209 endif
210 allocate(groupListLocal(nLocal), stat=alloc_stat)
211 if (alloc_stat /= 0 ) then
212 status = -1
213 return
214 endif
215 allocate(mfactLocal(nLocal), stat=alloc_stat)
216 if (alloc_stat /= 0 ) then
217 status = -1
218 return
219 endif
220
221 groupStartLocal = CgroupStart
222 groupListLocal = CgroupList
223 mfactLocal = Cmfact
224
225 call gather(groupStartLocal, groupStartRow, plan_group_row)
226 call gather(groupStartLocal, groupStartCol, plan_group_col)
227 groupStartRow(nGroupsInRow+1) = nAtomsInRow + 1
228 groupStartCol(nGroupsInCol+1) = nAtomsInCol + 1
229 call gather(groupListLocal, groupListRow, plan_atom_row)
230 call gather(groupListLocal, groupListCol, plan_atom_col)
231
232 ! C passes us groupList as globalID numbers for the atoms
233 ! we need to remap these onto the row and column ids on this
234 ! processor. This is a linear search, but is only done once
235 ! (we hope)
236
237 do i = 1, nAtomsInRow
238 do j = 1, nAtomsInRow
239 if (AtomRowToGlobal(j) .eq. groupListRow(i)) then
240 groupListRow(i) = j
241 endif
242 enddo
243 enddo
244 do i = 1, nAtomsInCol
245 do j = 1, nAtomsInCol
246 if (AtomColToGlobal(j) .eq. groupListCol(i)) then
247 groupListCol(i) = j
248 endif
249 enddo
250 enddo
251 call gather(mfactLocal, mfactRow, plan_atom_row)
252 call gather(mfactLocal, mfactCol, plan_atom_col)
253
254 if (allocated(mfactLocal)) then
255 deallocate(mfactLocal)
256 end if
257 if (allocated(groupListLocal)) then
258 deallocate(groupListLocal)
259 endif
260 if (allocated(groupStartLocal)) then
261 deallocate(groupStartLocal)
262 endif
263 #else
264 allocate(groupStartRow(nGroups+1),stat=alloc_stat)
265 if (alloc_stat /= 0 ) then
266 status = -1
267 return
268 endif
269 allocate(groupStartCol(nGroups+1),stat=alloc_stat)
270 if (alloc_stat /= 0 ) then
271 status = -1
272 return
273 endif
274 allocate(groupListRow(nLocal),stat=alloc_stat)
275 if (alloc_stat /= 0 ) then
276 status = -1
277 return
278 endif
279 allocate(groupListCol(nLocal),stat=alloc_stat)
280 if (alloc_stat /= 0 ) then
281 status = -1
282 return
283 endif
284 allocate(mfactRow(nLocal),stat=alloc_stat)
285 if (alloc_stat /= 0 ) then
286 status = -1
287 return
288 endif
289 allocate(mfactCol(nLocal),stat=alloc_stat)
290 if (alloc_stat /= 0 ) then
291 status = -1
292 return
293 endif
294 do i = 1, nGroups
295 groupStartRow(i) = CgroupStart(i)
296 groupStartCol(i) = CgroupStart(i)
297 end do
298 groupStartRow(nGroups+1) = nLocal + 1
299 groupStartCol(nGroups+1) = nLocal + 1
300 do i = 1, nLocal
301 groupListRow(i) = CgroupList(i)
302 groupListCol(i) = CgroupList(i)
303 mfactRow(i) = Cmfact(i)
304 mfactCol(i) = Cmfact(i)
305 end do
306
307 #endif
308
309
310 ! We build the local atid's for both mpi and nonmpi
311 do i = 1, nLocal
312
313 me = getFirstMatchingElement(atypes, "c_ident", c_idents(i))
314 atid(i) = me
315
316 enddo
317
318 do i = 1, nExcludes_Local
319 excludesLocal(1,i) = CexcludesLocal(1,i)
320 excludesLocal(2,i) = CexcludesLocal(2,i)
321 enddo
322
323 maxSkipsForAtom = 0
324 do j = 1, nLocal
325 nSkipsForAtom(j) = 0
326 #ifdef IS_MPI
327 id1 = AtomRowToGlobal(j)
328 #else
329 id1 = j
330 #endif
331 do i = 1, nExcludes_Local
332 if (excludesLocal(1,i) .eq. id1 ) then
333 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
334
335 if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
336 maxSkipsForAtom = nSkipsForAtom(j)
337 endif
338 endif
339 if (excludesLocal(2,i) .eq. id1 ) then
340 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
341
342 if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
343 maxSkipsForAtom = nSkipsForAtom(j)
344 endif
345 endif
346 end do
347 enddo
348
349 allocate(skipsForAtom(nLocal, maxSkipsForAtom), stat=alloc_stat)
350 if (alloc_stat /= 0 ) then
351 write(*,*) 'Could not allocate skipsForAtom array'
352 return
353 endif
354
355 do j = 1, nLocal
356 nSkipsForAtom(j) = 0
357 #ifdef IS_MPI
358 id1 = AtomRowToGlobal(j)
359 #else
360 id1 = j
361 #endif
362 do i = 1, nExcludes_Local
363 if (excludesLocal(1,i) .eq. id1 ) then
364 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
365 #ifdef IS_MPI
366 id2 = AtomColToGlobal(excludesLocal(2,i))
367 #else
368 id2 = excludesLocal(2,i)
369 #endif
370 skipsForAtom(j, nSkipsForAtom(j)) = id2
371 endif
372 if (excludesLocal(2, i) .eq. id2 ) then
373 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
374 #ifdef IS_MPI
375 id2 = AtomColToGlobal(excludesLocal(1,i))
376 #else
377 id2 = excludesLocal(1,i)
378 #endif
379 skipsForAtom(j, nSkipsForAtom(j)) = id2
380 endif
381 end do
382 enddo
383
384 do i = 1, nExcludes_Global
385 excludesGlobal(i) = CexcludesGlobal(i)
386 enddo
387
388 do i = 1, nGlobal
389 molMemberShipList(i) = CmolMembership(i)
390 enddo
391
392 if (status == 0) simulation_setup_complete = .true.
393
394 end subroutine SimulationSetup
395
396 subroutine setBox(cHmat, cHmatInv, cBoxIsOrthorhombic)
397 real(kind=dp), dimension(3,3) :: cHmat, cHmatInv
398 integer :: cBoxIsOrthorhombic
399 integer :: smallest, status, i
400
401 Hmat = cHmat
402 HmatInv = cHmatInv
403 if (cBoxIsOrthorhombic .eq. 0 ) then
404 boxIsOrthorhombic = .false.
405 else
406 boxIsOrthorhombic = .true.
407 endif
408
409 return
410 end subroutine setBox
411
412 function getDielect() result(dielect)
413 real( kind = dp ) :: dielect
414 dielect = thisSim%dielect
415 end function getDielect
416
417 function SimUsesPBC() result(doesit)
418 logical :: doesit
419 doesit = thisSim%SIM_uses_PBC
420 end function SimUsesPBC
421
422 function SimUsesLJ() result(doesit)
423 logical :: doesit
424 doesit = thisSim%SIM_uses_LJ
425 end function SimUsesLJ
426
427 function SimUsesSticky() result(doesit)
428 logical :: doesit
429 doesit = thisSim%SIM_uses_sticky
430 end function SimUsesSticky
431
432 function SimUsesCharges() result(doesit)
433 logical :: doesit
434 doesit = thisSim%SIM_uses_charges
435 end function SimUsesCharges
436
437 function SimUsesDipoles() result(doesit)
438 logical :: doesit
439 doesit = thisSim%SIM_uses_dipoles
440 end function SimUsesDipoles
441
442 function SimUsesRF() result(doesit)
443 logical :: doesit
444 doesit = thisSim%SIM_uses_RF
445 end function SimUsesRF
446
447 function SimUsesGB() result(doesit)
448 logical :: doesit
449 doesit = thisSim%SIM_uses_GB
450 end function SimUsesGB
451
452 function SimUsesEAM() result(doesit)
453 logical :: doesit
454 doesit = thisSim%SIM_uses_EAM
455 end function SimUsesEAM
456
457 function SimUsesDirectionalAtoms() result(doesit)
458 logical :: doesit
459 doesit = thisSim%SIM_uses_dipoles .or. thisSim%SIM_uses_sticky .or. &
460 thisSim%SIM_uses_GB .or. thisSim%SIM_uses_RF
461 end function SimUsesDirectionalAtoms
462
463 function SimRequiresPrepairCalc() result(doesit)
464 logical :: doesit
465 doesit = thisSim%SIM_uses_EAM
466 end function SimRequiresPrepairCalc
467
468 function SimRequiresPostpairCalc() result(doesit)
469 logical :: doesit
470 doesit = thisSim%SIM_uses_RF
471 end function SimRequiresPostpairCalc
472
473 subroutine InitializeSimGlobals(thisStat)
474 integer, intent(out) :: thisStat
475 integer :: alloc_stat
476
477 thisStat = 0
478
479 call FreeSimGlobals()
480
481 allocate(nSkipsForAtom(nLocal), stat=alloc_stat)
482 if (alloc_stat /= 0 ) then
483 thisStat = -1
484 return
485 endif
486
487 allocate(excludesLocal(2,nExcludes_Local), stat=alloc_stat)
488 if (alloc_stat /= 0 ) then
489 thisStat = -1
490 return
491 endif
492
493 allocate(excludesGlobal(nExcludes_Global), stat=alloc_stat)
494 if (alloc_stat /= 0 ) then
495 thisStat = -1
496 return
497 endif
498
499 allocate(molMembershipList(nGlobal), stat=alloc_stat)
500 if (alloc_stat /= 0 ) then
501 thisStat = -1
502 return
503 endif
504
505 end subroutine InitializeSimGlobals
506
507 subroutine FreeSimGlobals()
508
509 !We free in the opposite order in which we allocate in.
510
511 if (allocated(skipsForAtom)) deallocate(skipsForAtom)
512 if (allocated(mfactCol)) deallocate(mfactCol)
513 if (allocated(mfactRow)) deallocate(mfactRow)
514 if (allocated(groupListCol)) deallocate(groupListCol)
515 if (allocated(groupListRow)) deallocate(groupListRow)
516 if (allocated(groupStartCol)) deallocate(groupStartCol)
517 if (allocated(groupStartRow)) deallocate(groupStartRow)
518 if (allocated(molMembershipList)) deallocate(molMembershipList)
519 if (allocated(excludesGlobal)) deallocate(excludesGlobal)
520 if (allocated(excludesLocal)) deallocate(excludesLocal)
521
522 end subroutine FreeSimGlobals
523
524 pure function getNlocal() result(n)
525 integer :: n
526 n = nLocal
527 end function getNlocal
528
529
530 end module simulation