ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/simulation_module.F90
Revision: 1198
Committed: Thu May 27 00:48:12 2004 UTC (20 years, 1 month ago) by tim
File size: 14374 byte(s)
Log Message:
in the progress of fixing MPI version of cutoff group

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
228 groupStartCol(nGroupsInCol+1) = nAtomsInCol
229 call gather(groupListLocal, groupListRow, plan_atom_row)
230 call gather(groupListLocal, groupListCol, plan_atom_col)
231 call gather(mfactLocal, mfactRow, plan_atom_row)
232 call gather(mfactLocal, mfactCol, plan_atom_col)
233
234 if (allocated(mfactLocal)) then
235 deallocate(mfactLocal)
236 end if
237 if (allocated(groupListLocal)) then
238 deallocate(groupListLocal)
239 endif
240 if (allocated(groupStartLocal)) then
241 deallocate(groupStartLocal)
242 endif
243 #else
244 allocate(groupStartRow(nGroups+1),stat=alloc_stat)
245 if (alloc_stat /= 0 ) then
246 status = -1
247 return
248 endif
249 allocate(groupStartCol(nGroups+1),stat=alloc_stat)
250 if (alloc_stat /= 0 ) then
251 status = -1
252 return
253 endif
254 allocate(groupListRow(nLocal),stat=alloc_stat)
255 if (alloc_stat /= 0 ) then
256 status = -1
257 return
258 endif
259 allocate(groupListCol(nLocal),stat=alloc_stat)
260 if (alloc_stat /= 0 ) then
261 status = -1
262 return
263 endif
264 allocate(mfactRow(nLocal),stat=alloc_stat)
265 if (alloc_stat /= 0 ) then
266 status = -1
267 return
268 endif
269 allocate(mfactCol(nLocal),stat=alloc_stat)
270 if (alloc_stat /= 0 ) then
271 status = -1
272 return
273 endif
274 do i = 1, nGroups
275 groupStartRow(i) = CgroupStart(i)
276 groupStartCol(i) = CgroupStart(i)
277 end do
278 groupStartRow(nGroups+1) = nLocal
279 groupStartCol(nGroups+1) = nLocal
280 do i = 1, nLocal
281 groupListRow(i) = CgroupList(i)
282 groupListCol(i) = CgroupList(i)
283 mfactRow(i) = Cmfact(i)
284 mfactCol(i) = Cmfact(i)
285 end do
286
287 #endif
288
289
290 ! We build the local atid's for both mpi and nonmpi
291 do i = 1, nLocal
292
293 me = getFirstMatchingElement(atypes, "c_ident", c_idents(i))
294 atid(i) = me
295
296 enddo
297
298 do i = 1, nExcludes_Local
299 excludesLocal(1,i) = CexcludesLocal(1,i)
300 excludesLocal(2,i) = CexcludesLocal(2,i)
301 enddo
302
303 maxSkipsForAtom = 0
304 do j = 1, nLocal
305 nSkipsForAtom(j) = 0
306 #ifdef IS_MPI
307 id1 = AtomRowToGlobal(j)
308 #else
309 id1 = j
310 #endif
311 do i = 1, nExcludes_Local
312 if (excludesLocal(1,i) .eq. id1 ) then
313 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
314
315 if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
316 maxSkipsForAtom = nSkipsForAtom(j)
317 endif
318 endif
319 if (excludesLocal(2,i) .eq. id1 ) then
320 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
321
322 if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
323 maxSkipsForAtom = nSkipsForAtom(j)
324 endif
325 endif
326 end do
327 enddo
328
329 allocate(skipsForAtom(nLocal, maxSkipsForAtom), stat=alloc_stat)
330 if (alloc_stat /= 0 ) then
331 write(*,*) 'Could not allocate skipsForAtom array'
332 return
333 endif
334
335 do j = 1, nLocal
336 nSkipsForAtom(j) = 0
337 #ifdef IS_MPI
338 id1 = AtomRowToGlobal(j)
339 #else
340 id1 = j
341 #endif
342 do i = 1, nExcludes_Local
343 if (excludesLocal(1,i) .eq. id1 ) then
344 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
345 #ifdef IS_MPI
346 id2 = AtomColToGlobal(excludesLocal(2,i))
347 #else
348 id2 = excludesLocal(2,i)
349 #endif
350 skipsForAtom(j, nSkipsForAtom(j)) = id2
351 endif
352 if (excludesLocal(2, i) .eq. id2 ) then
353 nSkipsForAtom(j) = nSkipsForAtom(j) + 1
354 #ifdef IS_MPI
355 id2 = AtomColToGlobal(excludesLocal(1,i))
356 #else
357 id2 = excludesLocal(1,i)
358 #endif
359 skipsForAtom(j, nSkipsForAtom(j)) = id2
360 endif
361 end do
362 enddo
363
364 do i = 1, nExcludes_Global
365 excludesGlobal(i) = CexcludesGlobal(i)
366 enddo
367
368 do i = 1, nGlobal
369 molMemberShipList(i) = CmolMembership(i)
370 enddo
371
372 if (status == 0) simulation_setup_complete = .true.
373
374 end subroutine SimulationSetup
375
376 subroutine setBox(cHmat, cHmatInv, cBoxIsOrthorhombic)
377 real(kind=dp), dimension(3,3) :: cHmat, cHmatInv
378 integer :: cBoxIsOrthorhombic
379 integer :: smallest, status, i
380
381 Hmat = cHmat
382 HmatInv = cHmatInv
383 if (cBoxIsOrthorhombic .eq. 0 ) then
384 boxIsOrthorhombic = .false.
385 else
386 boxIsOrthorhombic = .true.
387 endif
388
389 return
390 end subroutine setBox
391
392 function getDielect() result(dielect)
393 real( kind = dp ) :: dielect
394 dielect = thisSim%dielect
395 end function getDielect
396
397 function SimUsesPBC() result(doesit)
398 logical :: doesit
399 doesit = thisSim%SIM_uses_PBC
400 end function SimUsesPBC
401
402 function SimUsesLJ() result(doesit)
403 logical :: doesit
404 doesit = thisSim%SIM_uses_LJ
405 end function SimUsesLJ
406
407 function SimUsesSticky() result(doesit)
408 logical :: doesit
409 doesit = thisSim%SIM_uses_sticky
410 end function SimUsesSticky
411
412 function SimUsesCharges() result(doesit)
413 logical :: doesit
414 doesit = thisSim%SIM_uses_charges
415 end function SimUsesCharges
416
417 function SimUsesDipoles() result(doesit)
418 logical :: doesit
419 doesit = thisSim%SIM_uses_dipoles
420 end function SimUsesDipoles
421
422 function SimUsesRF() result(doesit)
423 logical :: doesit
424 doesit = thisSim%SIM_uses_RF
425 end function SimUsesRF
426
427 function SimUsesGB() result(doesit)
428 logical :: doesit
429 doesit = thisSim%SIM_uses_GB
430 end function SimUsesGB
431
432 function SimUsesEAM() result(doesit)
433 logical :: doesit
434 doesit = thisSim%SIM_uses_EAM
435 end function SimUsesEAM
436
437 function SimUsesDirectionalAtoms() result(doesit)
438 logical :: doesit
439 doesit = thisSim%SIM_uses_dipoles .or. thisSim%SIM_uses_sticky .or. &
440 thisSim%SIM_uses_GB .or. thisSim%SIM_uses_RF
441 end function SimUsesDirectionalAtoms
442
443 function SimRequiresPrepairCalc() result(doesit)
444 logical :: doesit
445 doesit = thisSim%SIM_uses_EAM
446 end function SimRequiresPrepairCalc
447
448 function SimRequiresPostpairCalc() result(doesit)
449 logical :: doesit
450 doesit = thisSim%SIM_uses_RF
451 end function SimRequiresPostpairCalc
452
453 subroutine InitializeSimGlobals(thisStat)
454 integer, intent(out) :: thisStat
455 integer :: alloc_stat
456
457 thisStat = 0
458
459 call FreeSimGlobals()
460
461 allocate(nSkipsForAtom(nLocal), stat=alloc_stat)
462 if (alloc_stat /= 0 ) then
463 thisStat = -1
464 return
465 endif
466
467 allocate(excludesLocal(2,nExcludes_Local), stat=alloc_stat)
468 if (alloc_stat /= 0 ) then
469 thisStat = -1
470 return
471 endif
472
473 allocate(excludesGlobal(nExcludes_Global), stat=alloc_stat)
474 if (alloc_stat /= 0 ) then
475 thisStat = -1
476 return
477 endif
478
479 allocate(molMembershipList(nGlobal), stat=alloc_stat)
480 if (alloc_stat /= 0 ) then
481 thisStat = -1
482 return
483 endif
484
485 end subroutine InitializeSimGlobals
486
487 subroutine FreeSimGlobals()
488
489 !We free in the opposite order in which we allocate in.
490
491 if (allocated(skipsForAtom)) deallocate(skipsForAtom)
492 if (allocated(mfactCol)) deallocate(mfactCol)
493 if (allocated(mfactRow)) deallocate(mfactRow)
494 if (allocated(groupListCol)) deallocate(groupListCol)
495 if (allocated(groupListRow)) deallocate(groupListRow)
496 if (allocated(groupStartCol)) deallocate(groupStartCol)
497 if (allocated(groupStartRow)) deallocate(groupStartRow)
498 if (allocated(molMembershipList)) deallocate(molMembershipList)
499 if (allocated(excludesGlobal)) deallocate(excludesGlobal)
500 if (allocated(excludesLocal)) deallocate(excludesLocal)
501
502 end subroutine FreeSimGlobals
503
504 pure function getNlocal() result(n)
505 integer :: n
506 n = nLocal
507 end function getNlocal
508
509
510 end module simulation