ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/simulation_module.F90
Revision: 490
Committed: Fri Apr 11 15:16:59 2003 UTC (21 years, 3 months ago) by gezelter
File size: 9749 byte(s)
Log Message:
Bug fix in progress for NPT

File Contents

# User Rev Content
1 mmeineke 377 !! 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     #ifdef IS_MPI
10     use mpiSimulation
11     #endif
12    
13     implicit none
14     PRIVATE
15    
16     #define __FORTRAN90
17     #include "fSimulation.h"
18    
19     type (simtype), public :: thisSim
20    
21     logical, save :: simulation_setup_complete = .false.
22    
23     integer, public, save :: natoms
24     integer, public, save :: nExcludes_Global = 0
25     integer, public, save :: nExcludes_Local = 0
26     integer, allocatable, dimension(:,:), public :: excludesLocal
27     integer, allocatable, dimension(:), public :: excludesGlobal
28 chuckv 482 integer, allocatable, dimension(:), public :: molMembershipList
29 mmeineke 377
30     real(kind=dp), save :: rcut2 = 0.0_DP
31     real(kind=dp), save :: rcut6 = 0.0_DP
32     real(kind=dp), save :: rlist2 = 0.0_DP
33     real(kind=dp), public, dimension(3), save :: box
34    
35    
36     public :: SimulationSetup
37     public :: getNlocal
38     public :: setBox
39     public :: setBox_3d
40     public :: getBox
41     public :: setRcut
42     public :: getRcut
43     public :: getRlist
44     public :: getRrf
45     public :: getRt
46     public :: getDielect
47     public :: SimUsesPBC
48     public :: SimUsesLJ
49     public :: SimUsesDipoles
50     public :: SimUsesSticky
51     public :: SimUsesRF
52     public :: SimUsesGB
53     public :: SimUsesEAM
54     public :: SimRequiresPrepairCalc
55     public :: SimRequiresPostpairCalc
56     public :: SimUsesDirectionalAtoms
57    
58     interface getBox
59     module procedure getBox_3d
60     module procedure getBox_1d
61     end interface
62    
63     interface setBox
64     module procedure setBox_3d
65     module procedure setBox_1d
66     end interface
67    
68     contains
69    
70 gezelter 490 subroutine SimulationSetup(setThisSim, nGlobal, nLocal, c_idents, &
71 gezelter 483 CnLocalExcludes, CexcludesLocal, CnGlobalExcludes, CexcludesGlobal, &
72     CmolMembership, &
73 mmeineke 377 status)
74    
75     type (simtype) :: setThisSim
76 gezelter 490 integer, intent(inout) :: nGlobal, nLocal
77     integer, dimension(nLocal),intent(inout) :: c_idents
78 mmeineke 377
79     integer :: CnLocalExcludes
80     integer, dimension(2,CnLocalExcludes), intent(in) :: CexcludesLocal
81     integer :: CnGlobalExcludes
82     integer, dimension(CnGlobalExcludes), intent(in) :: CexcludesGlobal
83 gezelter 490 integer, dimension(nGlobal),intent(in) :: CmolMembership
84 mmeineke 377 !! Result status, success = 0, status = -1
85     integer, intent(out) :: status
86     integer :: i, me, thisStat, alloc_stat, myNode
87     #ifdef IS_MPI
88     integer, allocatable, dimension(:) :: c_idents_Row
89     integer, allocatable, dimension(:) :: c_idents_Col
90     integer :: nrow
91     integer :: ncol
92     #endif
93    
94     simulation_setup_complete = .false.
95     status = 0
96    
97     ! copy C struct into fortran type
98     thisSim = setThisSim
99 gezelter 490 natoms = nLocal
100 mmeineke 377 rcut2 = thisSim%rcut * thisSim%rcut
101     rcut6 = rcut2 * rcut2 * rcut2
102     rlist2 = thisSim%rlist * thisSim%rlist
103     box = thisSim%box
104    
105     nExcludes_Global = CnGlobalExcludes
106     nExcludes_Local = CnLocalExcludes
107    
108     call InitializeForceGlobals(natoms, thisStat)
109     if (thisStat /= 0) then
110 chuckv 480 write(default_error,*) "SimSetup: InitializeForceGlobals error"
111 mmeineke 377 status = -1
112     return
113     endif
114    
115     call InitializeSimGlobals(thisStat)
116     if (thisStat /= 0) then
117 chuckv 480 write(default_error,*) "SimSetup: InitializeSimGlobals error"
118 mmeineke 377 status = -1
119     return
120     endif
121    
122     #ifdef IS_MPI
123     ! We can only set up forces if mpiSimulation has been setup.
124     if (.not. isMPISimSet()) then
125     write(default_error,*) "MPI is not set"
126     status = -1
127     return
128     endif
129     nrow = getNrow(plan_row)
130     ncol = getNcol(plan_col)
131     mynode = getMyNode()
132    
133     allocate(c_idents_Row(nrow),stat=alloc_stat)
134     if (alloc_stat /= 0 ) then
135     status = -1
136     return
137     endif
138    
139     allocate(c_idents_Col(ncol),stat=alloc_stat)
140     if (alloc_stat /= 0 ) then
141     status = -1
142     return
143     endif
144    
145     call gather(c_idents, c_idents_Row, plan_row)
146     call gather(c_idents, c_idents_Col, plan_col)
147    
148     do i = 1, nrow
149     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Row(i))
150     atid_Row(i) = me
151     enddo
152    
153     do i = 1, ncol
154     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Col(i))
155     atid_Col(i) = me
156     enddo
157    
158     !! free temporary ident arrays
159     if (allocated(c_idents_Col)) then
160     deallocate(c_idents_Col)
161     end if
162     if (allocated(c_idents_Row)) then
163     deallocate(c_idents_Row)
164     endif
165    
166     #else
167 gezelter 490 do i = 1, nLocal
168 mmeineke 377
169     me = getFirstMatchingElement(atypes, "c_ident", c_idents(i))
170     atid(i) = me
171    
172     enddo
173     #endif
174    
175    
176 chuckv 388
177 mmeineke 377 do i = 1, nExcludes_Local
178     excludesLocal(1,i) = CexcludesLocal(1,i)
179     excludesLocal(2,i) = CexcludesLocal(2,i)
180     enddo
181    
182     do i = 1, nExcludes_Global
183     excludesGlobal(i) = CexcludesGlobal(i)
184     enddo
185 mmeineke 435
186 gezelter 490 do i = 1, nGlobal
187 gezelter 483 molMemberShipList(i) = CmolMembership(i)
188     ! write(0,*) 'molMembershipList(',i,') = ', molMemberShipList(i)
189     enddo
190 chuckv 482
191 mmeineke 377 if (status == 0) simulation_setup_complete = .true.
192    
193     end subroutine SimulationSetup
194    
195     subroutine setBox_3d(new_box_size)
196     real(kind=dp), dimension(3) :: new_box_size
197     integer :: smallest, status, i
198    
199     thisSim%box = new_box_size
200     box = thisSim%box
201    
202     return
203     end subroutine setBox_3d
204    
205     subroutine setBox_1d(dim, new_box_size)
206     integer :: dim, status
207     real(kind=dp) :: new_box_size
208     thisSim%box(dim) = new_box_size
209     box(dim) = thisSim%box(dim)
210     end subroutine setBox_1d
211    
212     subroutine setRcut(new_rcut, status)
213     real(kind = dp) :: new_rcut
214     integer :: myStatus, status
215     thisSim%rcut = new_rcut
216     rcut2 = thisSim%rcut * thisSim%rcut
217     rcut6 = rcut2 * rcut2 * rcut2
218     status = 0
219     return
220     end subroutine setRcut
221    
222     function getBox_3d() result(thisBox)
223     real( kind = dp ), dimension(3) :: thisBox
224     thisBox = thisSim%box
225     end function getBox_3d
226    
227     function getBox_1d(dim) result(thisBox)
228     integer, intent(in) :: dim
229     real( kind = dp ) :: thisBox
230    
231     thisBox = thisSim%box(dim)
232     end function getBox_1d
233    
234     subroutine getRcut(thisrcut,rc2,rc6,status)
235     real( kind = dp ), intent(out) :: thisrcut
236     real( kind = dp ), intent(out), optional :: rc2
237     real( kind = dp ), intent(out), optional :: rc6
238     integer, optional :: status
239    
240     if (present(status)) status = 0
241    
242     if (.not.simulation_setup_complete ) then
243     if (present(status)) status = -1
244     return
245     end if
246    
247     thisrcut = thisSim%rcut
248     if(present(rc2)) rc2 = rcut2
249     if(present(rc6)) rc6 = rcut6
250     end subroutine getRcut
251    
252     subroutine getRlist(thisrlist,rl2,status)
253     real( kind = dp ), intent(out) :: thisrlist
254     real( kind = dp ), intent(out), optional :: rl2
255    
256     integer, optional :: status
257    
258     if (present(status)) status = 0
259    
260     if (.not.simulation_setup_complete ) then
261     if (present(status)) status = -1
262     return
263     end if
264    
265     thisrlist = thisSim%rlist
266     if(present(rl2)) rl2 = rlist2
267     end subroutine getRlist
268    
269     function getRrf() result(rrf)
270     real( kind = dp ) :: rrf
271     rrf = thisSim%rrf
272 gezelter 394 write(*,*) 'getRrf = ', rrf, thisSim%rrf
273 mmeineke 377 end function getRrf
274    
275     function getRt() result(rt)
276     real( kind = dp ) :: rt
277     rt = thisSim%rt
278     end function getRt
279    
280     function getDielect() result(dielect)
281     real( kind = dp ) :: dielect
282     dielect = thisSim%dielect
283     end function getDielect
284    
285     function SimUsesPBC() result(doesit)
286     logical :: doesit
287     doesit = thisSim%SIM_uses_PBC
288     end function SimUsesPBC
289    
290     function SimUsesLJ() result(doesit)
291     logical :: doesit
292     doesit = thisSim%SIM_uses_LJ
293     end function SimUsesLJ
294    
295     function SimUsesSticky() result(doesit)
296     logical :: doesit
297     doesit = thisSim%SIM_uses_sticky
298     end function SimUsesSticky
299    
300     function SimUsesDipoles() result(doesit)
301     logical :: doesit
302     doesit = thisSim%SIM_uses_dipoles
303     end function SimUsesDipoles
304    
305     function SimUsesRF() result(doesit)
306     logical :: doesit
307     doesit = thisSim%SIM_uses_RF
308     end function SimUsesRF
309    
310     function SimUsesGB() result(doesit)
311     logical :: doesit
312     doesit = thisSim%SIM_uses_GB
313     end function SimUsesGB
314    
315     function SimUsesEAM() result(doesit)
316     logical :: doesit
317     doesit = thisSim%SIM_uses_EAM
318     end function SimUsesEAM
319    
320     function SimUsesDirectionalAtoms() result(doesit)
321     logical :: doesit
322     doesit = thisSim%SIM_uses_dipoles .or. thisSim%SIM_uses_sticky .or. &
323     thisSim%SIM_uses_GB .or. thisSim%SIM_uses_RF
324     end function SimUsesDirectionalAtoms
325    
326     function SimRequiresPrepairCalc() result(doesit)
327     logical :: doesit
328     doesit = thisSim%SIM_uses_EAM
329     end function SimRequiresPrepairCalc
330    
331     function SimRequiresPostpairCalc() result(doesit)
332     logical :: doesit
333     doesit = thisSim%SIM_uses_RF
334     end function SimRequiresPostpairCalc
335    
336     subroutine InitializeSimGlobals(thisStat)
337     integer, intent(out) :: thisStat
338     integer :: alloc_stat
339    
340     thisStat = 0
341    
342     call FreeSimGlobals()
343    
344     allocate(excludesLocal(2,nExcludes_Local), stat=alloc_stat)
345     if (alloc_stat /= 0 ) then
346     thisStat = -1
347     return
348     endif
349    
350     allocate(excludesGlobal(nExcludes_Global), stat=alloc_stat)
351     if (alloc_stat /= 0 ) then
352     thisStat = -1
353     return
354     endif
355 chuckv 482
356     allocate(molMembershipList(getNlocal()), stat=alloc_stat)
357     if (alloc_stat /= 0 ) then
358     thisStat = -1
359     return
360     endif
361 mmeineke 377
362     end subroutine InitializeSimGlobals
363    
364     subroutine FreeSimGlobals()
365    
366     !We free in the opposite order in which we allocate in.
367 gezelter 483
368     if (allocated(molMembershipList)) deallocate(molMembershipList)
369 mmeineke 377 if (allocated(excludesGlobal)) deallocate(excludesGlobal)
370     if (allocated(excludesLocal)) deallocate(excludesLocal)
371 gezelter 483
372 mmeineke 377 end subroutine FreeSimGlobals
373    
374     pure function getNlocal() result(nlocal)
375     integer :: nlocal
376     nlocal = natoms
377     end function getNlocal
378    
379    
380     end module simulation