ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/UseTheForce/DarkSide/simulation.F90
Revision: 1633
Committed: Fri Oct 22 20:22:48 2004 UTC (19 years, 8 months ago) by gezelter
File size: 17294 byte(s)
Log Message:
Added un-busticated fortran files and c/Fortran interfaces

File Contents

# User Rev Content
1 gezelter 1608 !! 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 "brains/fSimulation.h"
19     #include "UseTheForce/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 gezelter 1633
53     public :: SimUsesDirectionalAtoms
54     public :: SimUsesLennardJones
55     public :: SimUsesElectrostatics
56 gezelter 1608 public :: SimUsesCharges
57     public :: SimUsesDipoles
58     public :: SimUsesSticky
59 gezelter 1633 public :: SimUsesGayBerne
60     public :: SimUsesEAM
61     public :: SimUsesShapes
62     public :: SimUsesFLARB
63 gezelter 1608 public :: SimUsesRF
64     public :: SimRequiresPrepairCalc
65     public :: SimRequiresPostpairCalc
66 gezelter 1633
67 gezelter 1608
68     contains
69    
70     subroutine SimulationSetup(setThisSim, CnGlobal, CnLocal, c_idents, &
71     CnLocalExcludes, CexcludesLocal, CnGlobalExcludes, CexcludesGlobal, &
72     CmolMembership, Cmfact, CnGroups, CglobalGroupMembership, &
73     status)
74    
75     type (simtype) :: setThisSim
76     integer, intent(inout) :: CnGlobal, CnLocal
77     integer, dimension(CnLocal),intent(inout) :: c_idents
78    
79     integer :: CnLocalExcludes
80     integer, dimension(2,CnLocalExcludes), intent(in) :: CexcludesLocal
81     integer :: CnGlobalExcludes
82     integer, dimension(CnGlobalExcludes), intent(in) :: CexcludesGlobal
83     integer, dimension(CnGlobal),intent(in) :: CmolMembership
84     !! Result status, success = 0, status = -1
85     integer, intent(out) :: status
86     integer :: i, j, me, thisStat, alloc_stat, myNode, id1, id2
87     integer :: ia
88    
89     !! mass factors used for molecular cutoffs
90     real ( kind = dp ), dimension(CnLocal) :: Cmfact
91     integer, intent(in):: CnGroups
92     integer, dimension(CnGlobal), intent(in):: CglobalGroupMembership
93     integer :: maxSkipsForAtom, glPointer
94    
95     #ifdef IS_MPI
96     integer, allocatable, dimension(:) :: c_idents_Row
97     integer, allocatable, dimension(:) :: c_idents_Col
98     integer :: nAtomsInRow, nGroupsInRow, aid
99     integer :: nAtomsInCol, nGroupsInCol, gid
100     #endif
101    
102     simulation_setup_complete = .false.
103     status = 0
104    
105     ! copy C struct into fortran type
106    
107     nLocal = CnLocal
108     nGlobal = CnGlobal
109     nGroups = CnGroups
110    
111     thisSim = setThisSim
112    
113     nExcludes_Global = CnGlobalExcludes
114     nExcludes_Local = CnLocalExcludes
115    
116     call InitializeForceGlobals(nLocal, thisStat)
117     if (thisStat /= 0) then
118     write(default_error,*) "SimSetup: InitializeForceGlobals error"
119     status = -1
120     return
121     endif
122    
123     call InitializeSimGlobals(thisStat)
124     if (thisStat /= 0) then
125     write(default_error,*) "SimSetup: InitializeSimGlobals error"
126     status = -1
127     return
128     endif
129    
130     #ifdef IS_MPI
131     ! We can only set up forces if mpiSimulation has been setup.
132     if (.not. isMPISimSet()) then
133     write(default_error,*) "MPI is not set"
134     status = -1
135     return
136     endif
137     nAtomsInRow = getNatomsInRow(plan_atom_row)
138     nAtomsInCol = getNatomsInCol(plan_atom_col)
139     nGroupsInRow = getNgroupsInRow(plan_group_row)
140     nGroupsInCol = getNgroupsInCol(plan_group_col)
141     mynode = getMyNode()
142    
143     allocate(c_idents_Row(nAtomsInRow),stat=alloc_stat)
144     if (alloc_stat /= 0 ) then
145     status = -1
146     return
147     endif
148    
149     allocate(c_idents_Col(nAtomsInCol),stat=alloc_stat)
150     if (alloc_stat /= 0 ) then
151     status = -1
152     return
153     endif
154    
155     call gather(c_idents, c_idents_Row, plan_atom_row)
156     call gather(c_idents, c_idents_Col, plan_atom_col)
157    
158     do i = 1, nAtomsInRow
159     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Row(i))
160     atid_Row(i) = me
161     enddo
162    
163     do i = 1, nAtomsInCol
164     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Col(i))
165     atid_Col(i) = me
166     enddo
167    
168     !! free temporary ident arrays
169     if (allocated(c_idents_Col)) then
170     deallocate(c_idents_Col)
171     end if
172     if (allocated(c_idents_Row)) then
173     deallocate(c_idents_Row)
174     endif
175    
176     #endif
177    
178     #ifdef IS_MPI
179     allocate(groupStartRow(nGroupsInRow+1),stat=alloc_stat)
180     if (alloc_stat /= 0 ) then
181     status = -1
182     return
183     endif
184     allocate(groupStartCol(nGroupsInCol+1),stat=alloc_stat)
185     if (alloc_stat /= 0 ) then
186     status = -1
187     return
188     endif
189     allocate(groupListRow(nAtomsInRow),stat=alloc_stat)
190     if (alloc_stat /= 0 ) then
191     status = -1
192     return
193     endif
194     allocate(groupListCol(nAtomsInCol),stat=alloc_stat)
195     if (alloc_stat /= 0 ) then
196     status = -1
197     return
198     endif
199     allocate(mfactRow(nAtomsInRow),stat=alloc_stat)
200     if (alloc_stat /= 0 ) then
201     status = -1
202     return
203     endif
204     allocate(mfactCol(nAtomsInCol),stat=alloc_stat)
205     if (alloc_stat /= 0 ) then
206     status = -1
207     return
208     endif
209     allocate(mfactLocal(nLocal),stat=alloc_stat)
210     if (alloc_stat /= 0 ) then
211     status = -1
212     return
213     endif
214    
215     glPointer = 1
216    
217     do i = 1, nGroupsInRow
218    
219     gid = GroupRowToGlobal(i)
220     groupStartRow(i) = glPointer
221    
222     do j = 1, nAtomsInRow
223     aid = AtomRowToGlobal(j)
224     if (CglobalGroupMembership(aid) .eq. gid) then
225     groupListRow(glPointer) = j
226     glPointer = glPointer + 1
227     endif
228     enddo
229     enddo
230     groupStartRow(nGroupsInRow+1) = nAtomsInRow + 1
231    
232     glPointer = 1
233    
234     do i = 1, nGroupsInCol
235    
236     gid = GroupColToGlobal(i)
237     groupStartCol(i) = glPointer
238    
239     do j = 1, nAtomsInCol
240     aid = AtomColToGlobal(j)
241     if (CglobalGroupMembership(aid) .eq. gid) then
242     groupListCol(glPointer) = j
243     glPointer = glPointer + 1
244     endif
245     enddo
246     enddo
247     groupStartCol(nGroupsInCol+1) = nAtomsInCol + 1
248    
249     mfactLocal = Cmfact
250    
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     #else
258     allocate(groupStartRow(nGroups+1),stat=alloc_stat)
259     if (alloc_stat /= 0 ) then
260     status = -1
261     return
262     endif
263     allocate(groupStartCol(nGroups+1),stat=alloc_stat)
264     if (alloc_stat /= 0 ) then
265     status = -1
266     return
267     endif
268     allocate(groupListRow(nLocal),stat=alloc_stat)
269     if (alloc_stat /= 0 ) then
270     status = -1
271     return
272     endif
273     allocate(groupListCol(nLocal),stat=alloc_stat)
274     if (alloc_stat /= 0 ) then
275     status = -1
276     return
277     endif
278     allocate(mfactRow(nLocal),stat=alloc_stat)
279     if (alloc_stat /= 0 ) then
280     status = -1
281     return
282     endif
283     allocate(mfactCol(nLocal),stat=alloc_stat)
284     if (alloc_stat /= 0 ) then
285     status = -1
286     return
287     endif
288     allocate(mfactLocal(nLocal),stat=alloc_stat)
289     if (alloc_stat /= 0 ) then
290     status = -1
291     return
292     endif
293    
294     glPointer = 1
295     do i = 1, nGroups
296     groupStartRow(i) = glPointer
297     groupStartCol(i) = glPointer
298     do j = 1, nLocal
299     if (CglobalGroupMembership(j) .eq. i) then
300     groupListRow(glPointer) = j
301     groupListCol(glPointer) = j
302     glPointer = glPointer + 1
303     endif
304     enddo
305     enddo
306     groupStartRow(nGroups+1) = nLocal + 1
307     groupStartCol(nGroups+1) = nLocal + 1
308    
309     do i = 1, nLocal
310     mfactRow(i) = Cmfact(i)
311     mfactCol(i) = Cmfact(i)
312     end do
313    
314     #endif
315    
316    
317     ! We build the local atid's for both mpi and nonmpi
318     do i = 1, nLocal
319    
320     me = getFirstMatchingElement(atypes, "c_ident", c_idents(i))
321     atid(i) = me
322    
323     enddo
324    
325     do i = 1, nExcludes_Local
326     excludesLocal(1,i) = CexcludesLocal(1,i)
327     excludesLocal(2,i) = CexcludesLocal(2,i)
328     enddo
329    
330     #ifdef IS_MPI
331     allocate(nSkipsForAtom(nAtomsInRow), stat=alloc_stat)
332     #else
333     allocate(nSkipsForAtom(nLocal), stat=alloc_stat)
334     #endif
335     if (alloc_stat /= 0 ) then
336     thisStat = -1
337     write(*,*) 'Could not allocate nSkipsForAtom array'
338     return
339     endif
340    
341     maxSkipsForAtom = 0
342     #ifdef IS_MPI
343     do j = 1, nAtomsInRow
344     #else
345     do j = 1, nLocal
346     #endif
347     nSkipsForAtom(j) = 0
348     #ifdef IS_MPI
349     id1 = AtomRowToGlobal(j)
350     #else
351     id1 = j
352     #endif
353     do i = 1, nExcludes_Local
354     if (excludesLocal(1,i) .eq. id1 ) then
355     nSkipsForAtom(j) = nSkipsForAtom(j) + 1
356    
357     if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
358     maxSkipsForAtom = nSkipsForAtom(j)
359     endif
360     endif
361     if (excludesLocal(2,i) .eq. id1 ) then
362     nSkipsForAtom(j) = nSkipsForAtom(j) + 1
363    
364     if (nSkipsForAtom(j) .gt. maxSkipsForAtom) then
365     maxSkipsForAtom = nSkipsForAtom(j)
366     endif
367     endif
368     end do
369     enddo
370    
371     #ifdef IS_MPI
372     allocate(skipsForAtom(nAtomsInRow, maxSkipsForAtom), stat=alloc_stat)
373     #else
374     allocate(skipsForAtom(nLocal, maxSkipsForAtom), stat=alloc_stat)
375     #endif
376     if (alloc_stat /= 0 ) then
377     write(*,*) 'Could not allocate skipsForAtom array'
378     return
379     endif
380    
381     #ifdef IS_MPI
382     do j = 1, nAtomsInRow
383     #else
384     do j = 1, nLocal
385     #endif
386     nSkipsForAtom(j) = 0
387     #ifdef IS_MPI
388     id1 = AtomRowToGlobal(j)
389     #else
390     id1 = j
391     #endif
392     do i = 1, nExcludes_Local
393     if (excludesLocal(1,i) .eq. id1 ) then
394     nSkipsForAtom(j) = nSkipsForAtom(j) + 1
395     ! exclude lists have global ID's so this line is
396     ! the same in MPI and non-MPI
397     id2 = excludesLocal(2,i)
398     skipsForAtom(j, nSkipsForAtom(j)) = id2
399     endif
400     if (excludesLocal(2, i) .eq. id1 ) then
401     nSkipsForAtom(j) = nSkipsForAtom(j) + 1
402     ! exclude lists have global ID's so this line is
403     ! the same in MPI and non-MPI
404     id2 = excludesLocal(1,i)
405     skipsForAtom(j, nSkipsForAtom(j)) = id2
406     endif
407     end do
408     enddo
409    
410     do i = 1, nExcludes_Global
411     excludesGlobal(i) = CexcludesGlobal(i)
412     enddo
413    
414     do i = 1, nGlobal
415     molMemberShipList(i) = CmolMembership(i)
416     enddo
417    
418     if (status == 0) simulation_setup_complete = .true.
419    
420     end subroutine SimulationSetup
421    
422     subroutine setBox(cHmat, cHmatInv, cBoxIsOrthorhombic)
423     real(kind=dp), dimension(3,3) :: cHmat, cHmatInv
424     integer :: cBoxIsOrthorhombic
425     integer :: smallest, status, i
426    
427     Hmat = cHmat
428     HmatInv = cHmatInv
429     if (cBoxIsOrthorhombic .eq. 0 ) then
430     boxIsOrthorhombic = .false.
431     else
432     boxIsOrthorhombic = .true.
433     endif
434    
435     return
436     end subroutine setBox
437    
438     function getDielect() result(dielect)
439     real( kind = dp ) :: dielect
440     dielect = thisSim%dielect
441     end function getDielect
442    
443     function SimUsesPBC() result(doesit)
444     logical :: doesit
445     doesit = thisSim%SIM_uses_PBC
446     end function SimUsesPBC
447    
448 gezelter 1633 function SimUsesDirectionalAtoms() result(doesit)
449 gezelter 1608 logical :: doesit
450 gezelter 1633 doesit = thisSim%SIM_uses_dipoles .or. thisSim%SIM_uses_sticky .or. &
451     thisSim%SIM_uses_GayBerne .or. thisSim%SIM_uses_Shapes
452     end function SimUsesDirectionalAtoms
453 gezelter 1608
454 gezelter 1633 function SimUsesLennardJones() result(doesit)
455 gezelter 1608 logical :: doesit
456 gezelter 1633 doesit = thisSim%SIM_uses_LennardJones
457     end function SimUsesLennardJones
458 gezelter 1608
459 gezelter 1633 function SimUsesElectrostatics() result(doesit)
460     logical :: doesit
461     doesit = thisSim%SIM_uses_Electrostatics
462     end function SimUsesElectrostatics
463    
464 gezelter 1608 function SimUsesCharges() result(doesit)
465     logical :: doesit
466 gezelter 1633 doesit = thisSim%SIM_uses_Charges
467 gezelter 1608 end function SimUsesCharges
468    
469     function SimUsesDipoles() result(doesit)
470     logical :: doesit
471 gezelter 1633 doesit = thisSim%SIM_uses_Dipoles
472 gezelter 1608 end function SimUsesDipoles
473    
474 gezelter 1633 function SimUsesSticky() result(doesit)
475 gezelter 1608 logical :: doesit
476 gezelter 1633 doesit = thisSim%SIM_uses_Sticky
477     end function SimUsesSticky
478 gezelter 1608
479 gezelter 1633 function SimUsesGayBerne() result(doesit)
480 gezelter 1608 logical :: doesit
481 gezelter 1633 doesit = thisSim%SIM_uses_GayBerne
482     end function SimUsesGayBerne
483    
484 gezelter 1608 function SimUsesEAM() result(doesit)
485     logical :: doesit
486     doesit = thisSim%SIM_uses_EAM
487     end function SimUsesEAM
488    
489 gezelter 1633 function SimUsesShapes() result(doesit)
490 gezelter 1608 logical :: doesit
491 gezelter 1633 doesit = thisSim%SIM_uses_Shapes
492     end function SimUsesShapes
493 gezelter 1608
494 gezelter 1633 function SimUsesFLARB() result(doesit)
495     logical :: doesit
496     doesit = thisSim%SIM_uses_FLARB
497     end function SimUsesFLARB
498    
499     function SimUsesRF() result(doesit)
500     logical :: doesit
501     doesit = thisSim%SIM_uses_RF
502     end function SimUsesRF
503    
504 gezelter 1608 function SimRequiresPrepairCalc() result(doesit)
505     logical :: doesit
506     doesit = thisSim%SIM_uses_EAM
507     end function SimRequiresPrepairCalc
508    
509     function SimRequiresPostpairCalc() result(doesit)
510     logical :: doesit
511     doesit = thisSim%SIM_uses_RF
512     end function SimRequiresPostpairCalc
513    
514     subroutine InitializeSimGlobals(thisStat)
515     integer, intent(out) :: thisStat
516     integer :: alloc_stat
517    
518     thisStat = 0
519    
520     call FreeSimGlobals()
521    
522     allocate(excludesLocal(2,nExcludes_Local), stat=alloc_stat)
523     if (alloc_stat /= 0 ) then
524     thisStat = -1
525     return
526     endif
527    
528     allocate(excludesGlobal(nExcludes_Global), stat=alloc_stat)
529     if (alloc_stat /= 0 ) then
530     thisStat = -1
531     return
532     endif
533    
534     allocate(molMembershipList(nGlobal), stat=alloc_stat)
535     if (alloc_stat /= 0 ) then
536     thisStat = -1
537     return
538     endif
539    
540     end subroutine InitializeSimGlobals
541    
542     subroutine FreeSimGlobals()
543    
544     !We free in the opposite order in which we allocate in.
545    
546     if (allocated(skipsForAtom)) deallocate(skipsForAtom)
547     if (allocated(nSkipsForAtom)) deallocate(nSkipsForAtom)
548     if (allocated(mfactLocal)) deallocate(mfactLocal)
549     if (allocated(mfactCol)) deallocate(mfactCol)
550     if (allocated(mfactRow)) deallocate(mfactRow)
551     if (allocated(groupListCol)) deallocate(groupListCol)
552     if (allocated(groupListRow)) deallocate(groupListRow)
553     if (allocated(groupStartCol)) deallocate(groupStartCol)
554     if (allocated(groupStartRow)) deallocate(groupStartRow)
555     if (allocated(molMembershipList)) deallocate(molMembershipList)
556     if (allocated(excludesGlobal)) deallocate(excludesGlobal)
557     if (allocated(excludesLocal)) deallocate(excludesLocal)
558    
559     end subroutine FreeSimGlobals
560    
561     pure function getNlocal() result(n)
562     integer :: n
563     n = nLocal
564     end function getNlocal
565    
566    
567     end module simulation
568    
569    
570     subroutine setFortranSim(setThisSim, CnGlobal, CnLocal, c_idents, &
571     CnLocalExcludes, CexcludesLocal, CnGlobalExcludes, CexcludesGlobal, &
572     CmolMembership, Cmfact, CnGroups, CglobalGroupMembership, &
573     status)
574     use definitions, ONLY : dp
575     use simulation
576    
577     type (simtype) :: setThisSim
578     integer, intent(inout) :: CnGlobal, CnLocal
579     integer, dimension(CnLocal),intent(inout) :: c_idents
580    
581     integer :: CnLocalExcludes
582     integer, dimension(2,CnLocalExcludes), intent(inout) :: CexcludesLocal
583     integer :: CnGlobalExcludes
584     integer, dimension(CnGlobalExcludes), intent(inout) :: CexcludesGlobal
585     integer, dimension(CnGlobal),intent(inout) :: CmolMembership
586     !! Result status, success = 0, status = -1
587     integer, intent(inout) :: status
588    
589     !! mass factors used for molecular cutoffs
590     real ( kind = dp ), dimension(CnLocal) :: Cmfact
591     integer, intent(in):: CnGroups
592     integer, dimension(CnGlobal), intent(inout):: CglobalGroupMembership
593    
594     call SimulationSetup(setThisSim, CnGlobal, CnLocal, c_idents, &
595     CnLocalExcludes, CexcludesLocal, CnGlobalExcludes, CexcludesGlobal, &
596     CmolMembership, Cmfact, CnGroups, CglobalGroupMembership, &
597     status)
598     end subroutine setFortranSim
599    
600     subroutine setFortranBox(cHmat, cHmatInv, cBoxIsOrthorhombic)
601     use simulation, only : setBox
602     use definitions, ONLY : dp
603     real(kind=dp), dimension(3,3) :: cHmat, cHmatInv
604     integer :: cBoxIsOrthorhombic
605    
606     call setBox(cHmat, cHmatInv, cBoxIsOrthorhombic)
607    
608 gezelter 1609 end subroutine setFortranBox