ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE_old/src/mdtools/libmdCode/simulation_module.F90
Revision: 330
Committed: Wed Mar 12 23:15:46 2003 UTC (21 years, 5 months ago) by gezelter
File size: 13406 byte(s)
Log Message:
forceGlobals is gone (part3)

File Contents

# User Rev Content
1 chuckv 290 !! Fortran interface to C entry plug.
2    
3 mmeineke 270 module simulation
4 gezelter 325 use definitions
5     use neighborLists
6 gezelter 328 use atype_module
7 mmeineke 270 #ifdef IS_MPI
8     use mpiSimulation
9     #endif
10    
11     implicit none
12     PRIVATE
13    
14 mmeineke 285 #define __FORTRAN90
15 gezelter 309 #include "fSimulation.h"
16 mmeineke 270
17     type (simtype), public :: thisSim
18    
19 gezelter 325 logical, save :: simulation_setup_complete = .false.
20 mmeineke 270
21 gezelter 325 integer :: natoms
22 gezelter 330 integer, public, save :: nExcludes_Global = 0
23     integer, public, save :: nExcludes_Local = 0
24 gezelter 325 real(kind=dp), save :: rcut2 = 0.0_DP
25     real(kind=dp), save :: rcut6 = 0.0_DP
26     real(kind=dp), save :: rlist2 = 0.0_DP
27 gezelter 330 real(kind=dp), public, dimension(3), save :: box
28 mmeineke 270
29 gezelter 325 #ifdef IS_MPI
30     real( kind = dp ), allocatable, dimension(:,:), public :: q_Row
31     real( kind = dp ), allocatable, dimension(:,:), public :: q_Col
32     real( kind = dp ), allocatable, dimension(:,:), public :: u_l_Row
33     real( kind = dp ), allocatable, dimension(:,:), public :: u_l_Col
34     real( kind = dp ), allocatable, dimension(:,:), public :: A_Row
35     real( kind = dp ), allocatable, dimension(:,:), public :: A_Col
36    
37     real( kind = dp ), allocatable, dimension(:), public :: pot_Row
38     real( kind = dp ), allocatable, dimension(:), public :: pot_Col
39     real( kind = dp ), allocatable, dimension(:), public :: pot_Temp
40     real( kind = dp ), allocatable, dimension(:,:), public :: f_Row
41     real( kind = dp ), allocatable, dimension(:,:), public :: f_Col
42     real( kind = dp ), allocatable, dimension(:,:), public :: f_Temp
43     real( kind = dp ), allocatable, dimension(:,:), public :: t_Row
44     real( kind = dp ), allocatable, dimension(:,:), public :: t_Col
45     real( kind = dp ), allocatable, dimension(:,:), public :: t_Temp
46 gezelter 329 real( kind = dp ), allocatable, dimension(:,:), public :: rf_Row
47     real( kind = dp ), allocatable, dimension(:,:), public :: rf_Col
48 gezelter 330 real( kind = dp ), allocatable, dimension(:,:), public :: rf_Temp
49 gezelter 329
50 gezelter 325 integer, allocatable, dimension(:), public :: atid_Row
51     integer, allocatable, dimension(:), public :: atid_Col
52     #else
53     integer, allocatable, dimension(:), public :: atid
54     #endif
55 gezelter 330 real( kind = dp ), allocatable, dimension(:,:), public :: rf
56 gezelter 329 real(kind = dp), dimension(9), public :: tau_Temp = 0.0_dp
57     real(kind = dp), public :: virial_Temp = 0.0_dp
58 gezelter 325
59     public :: SimulationSetup
60 mmeineke 270 public :: getBox
61     public :: getRcut
62     public :: getRlist
63 gezelter 312 public :: getRrf
64     public :: getRt
65 gezelter 329 public :: getDielect
66 mmeineke 270 public :: getNlocal
67 gezelter 325 public :: SimUsesPBC
68     public :: SimUsesLJ
69     public :: SimUsesDipoles
70     public :: SimUsesSticky
71     public :: SimUsesRF
72     public :: SimUsesGB
73     public :: SimUsesEAM
74     public :: SimRequiresPrepairCalc
75     public :: SimRequiresPostpairCalc
76 mmeineke 285
77 mmeineke 270 interface getBox
78     module procedure getBox_3d
79     module procedure getBox_dim
80     end interface
81 gezelter 325
82 mmeineke 270 contains
83 gezelter 325
84     subroutine SimulationSetup(setThisSim, nComponents, c_idents, &
85     nExcludes_local, excludesLocal, nExcludes_global, excludesGlobal, &
86     status)
87    
88 chuckv 290 type (simtype) :: setThisSim
89 gezelter 325 integer, intent(inout) :: nComponents
90     integer, dimension(nComponents),intent(inout) :: c_idents
91     integer :: nExcludes_local
92     integer, dimension(nExcludes_local),intent(inout) :: excludesLocal
93     integer :: nExcludes_global
94     integer, dimension(nExcludes_global),intent(inout) :: excludesGlobal
95     !! Result status, success = 0, status = -1
96     integer, intent(out) :: status
97     integer :: i, me, thisStat, alloc_stat, myNode
98 chuckv 290
99 gezelter 325 #ifdef IS_MPI
100     integer, allocatable, dimension(:) :: c_idents_Row
101     integer, allocatable, dimension(:) :: c_idents_Col
102     integer :: nrow
103     integer :: ncol
104     #endif
105    
106     simulation_setup_complete = .false.
107     status = 0
108    
109 gezelter 312 ! copy C struct into fortran type
110 chuckv 290 thisSim = setThisSim
111 gezelter 325 natoms = nComponents
112     rcut2 = thisSim%rcut * thisSim%rcut
113     rcut6 = rcut2 * rcut2 * rcut2
114     rlist2 = thisSim%rlist * thisSim%rlist
115 gezelter 330 box = thisSim%box
116 mmeineke 270
117 gezelter 325 #ifdef IS_MPI
118     ! We can only set up forces if mpiSimulation has been setup.
119     if (.not. isMPISimSet()) then
120     write(default_error,*) "MPI is not set"
121     status = -1
122     return
123     endif
124     nrow = getNrow(plan_row)
125     ncol = getNcol(plan_col)
126     mynode = getMyNode()
127    
128     allocate(c_idents_Row(nrow),stat=alloc_stat)
129     if (alloc_stat /= 0 ) then
130     status = -1
131     return
132     endif
133 mmeineke 270
134 gezelter 325 allocate(c_idents_Col(ncol),stat=alloc_stat)
135     if (alloc_stat /= 0 ) then
136     status = -1
137     return
138     endif
139    
140     call gather(c_idents, c_idents_Row, plan_row)
141     call gather(c_idents, c_idents_Col, plan_col)
142    
143     allocate(atid_Row(nrow),stat=alloc_stat)
144     if (alloc_stat /= 0 ) then
145     status = -1
146     return
147     endif
148    
149     allocate(atid_Col(ncol),stat=alloc_stat)
150     if (alloc_stat /= 0 ) then
151     status = -1
152     return
153     endif
154    
155     do i = 1, nrow
156     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Row(i))
157     atid_Row(i) = me
158     enddo
159    
160     do i = 1, ncol
161     me = getFirstMatchingElement(atypes, "c_ident", c_idents_Col(i))
162     atid_Col(i) = me
163     enddo
164    
165     !! free temporary ident arrays
166     if (allocated(c_idents_Col)) then
167     deallocate(c_idents_Col)
168     end if
169     if (allocated(c_idents_Row)) then
170     deallocate(c_idents_Row)
171     endif
172    
173     #else
174     do i = 1, nComponents
175    
176     me = getFirstMatchingElement(atypes, "c_ident", c_idents(i))
177     atid(i) = me
178    
179     enddo
180     #endif
181    
182     call setupGlobals(thisStat)
183     if (thisStat /= 0) then
184     status = -1
185     return
186     endif
187    
188     !! Create neighbor lists
189     call expandNeighborList(nComponents, thisStat)
190     if (thisStat /= 0) then
191     status = -1
192     return
193     endif
194    
195     if (status == 0) simulation_setup_complete = .true.
196    
197     end subroutine SimulationSetup
198    
199 mmeineke 270 subroutine change_box_size(new_box_size)
200     real(kind=dp), dimension(3) :: new_box_size
201     thisSim%box = new_box_size
202 gezelter 330 box = thisSim%box
203 mmeineke 270 end subroutine change_box_size
204    
205     function getBox_3d() result(thisBox)
206     real( kind = dp ), dimension(3) :: thisBox
207     thisBox = thisSim%box
208     end function getBox_3d
209    
210     function getBox_dim(dim) result(thisBox)
211     integer, intent(in) :: dim
212     real( kind = dp ) :: thisBox
213 gezelter 312
214 mmeineke 270 thisBox = thisSim%box(dim)
215     end function getBox_dim
216 gezelter 312
217 gezelter 325 subroutine getRcut(thisrcut,rc2,rc6,status)
218 mmeineke 270 real( kind = dp ), intent(out) :: thisrcut
219 gezelter 325 real( kind = dp ), intent(out), optional :: rc2
220     real( kind = dp ), intent(out), optional :: rc6
221 mmeineke 270 integer, optional :: status
222    
223     if (present(status)) status = 0
224 gezelter 325
225     if (.not.simulation_setup_complete ) then
226 mmeineke 270 if (present(status)) status = -1
227     return
228     end if
229    
230     thisrcut = thisSim%rcut
231 gezelter 325 if(present(rc2)) rc2 = rcut2
232     if(present(rc6)) rc6 = rcut6
233 mmeineke 270 end subroutine getRcut
234    
235 gezelter 325 subroutine getRlist(thisrlist,rl2,status)
236 mmeineke 270 real( kind = dp ), intent(out) :: thisrlist
237 gezelter 325 real( kind = dp ), intent(out), optional :: rl2
238 mmeineke 270
239     integer, optional :: status
240    
241     if (present(status)) status = 0
242    
243 gezelter 325 if (.not.simulation_setup_complete ) then
244 mmeineke 270 if (present(status)) status = -1
245     return
246     end if
247    
248     thisrlist = thisSim%rlist
249 gezelter 325 if(present(rl2)) rl2 = rlist2
250 gezelter 312 end subroutine getRlist
251 mmeineke 270
252 gezelter 312 function getRrf() result(rrf)
253     real( kind = dp ) :: rrf
254     rrf = thisSim%rrf
255     end function getRrf
256 mmeineke 270
257 gezelter 312 function getRt() result(rt)
258     real( kind = dp ) :: rt
259     rt = thisSim%rt
260     end function getRt
261 gezelter 329
262     function getDielect() result(dielect)
263     real( kind = dp ) :: dielect
264     dielect = thisSim%dielect
265     end function getDielect
266 gezelter 312
267     pure function getNlocal() result(nlocal)
268 mmeineke 270 integer :: nlocal
269 gezelter 325 nlocal = natoms
270 mmeineke 270 end function getNlocal
271 gezelter 325
272     function SimUsesPBC() result(doesit)
273     logical :: doesit
274     doesit = thisSim%SIM_uses_PBC
275     end function SimUsesPBC
276    
277     function SimUsesLJ() result(doesit)
278     logical :: doesit
279     doesit = thisSim%SIM_uses_LJ
280     end function SimUsesLJ
281    
282     function SimUsesSticky() result(doesit)
283     logical :: doesit
284     doesit = thisSim%SIM_uses_sticky
285     end function SimUsesSticky
286    
287     function SimUsesDipoles() result(doesit)
288     logical :: doesit
289     doesit = thisSim%SIM_uses_dipoles
290     end function SimUsesDipoles
291    
292     function SimUsesRF() result(doesit)
293     logical :: doesit
294     doesit = thisSim%SIM_uses_RF
295     end function SimUsesRF
296    
297     function SimUsesGB() result(doesit)
298     logical :: doesit
299     doesit = thisSim%SIM_uses_GB
300     end function SimUsesGB
301    
302     function SimUsesEAM() result(doesit)
303     logical :: doesit
304     doesit = thisSim%SIM_uses_EAM
305     end function SimUsesEAM
306    
307 gezelter 329 function SimUsesDirectionalAtoms() result(doesit)
308     logical :: doesit
309     doesit = thisSim%SIM_uses_dipoles .or. thisSim%SIM_uses_sticky .or. &
310     thisSim%SIM_uses_GB .or. thisSim%SIM_uses_RF
311 gezelter 330 end function SimUsesDirectionalAtoms
312 gezelter 329
313 gezelter 325 function SimRequiresPrepairCalc() result(doesit)
314     logical :: doesit
315     doesit = thisSim%SIM_uses_EAM
316     end function SimRequiresPrepairCalc
317    
318     function SimRequiresPostpairCalc() result(doesit)
319     logical :: doesit
320     doesit = thisSim%SIM_uses_RF
321     end function SimRequiresPostpairCalc
322 gezelter 312
323 gezelter 325 subroutine setupGlobals(thisStat)
324     integer, intent(out) :: thisStat
325     integer :: nrow
326     integer :: ncol
327     integer :: nlocal
328     integer :: ndim = 3
329     integer :: alloc_stat
330    
331     thisStat = 0
332    
333     #ifdef IS_MPI
334     nrow = getNrow(plan_row)
335     ncol = getNcol(plan_col)
336     #endif
337     nlocal = getNlocal()
338    
339     call freeGlobals()
340    
341     #ifdef IS_MPI
342    
343     allocate(q_Row(ndim,nrow),stat=alloc_stat)
344     if (alloc_stat /= 0 ) then
345     thisStat = 0
346     return
347     endif
348    
349     allocate(q_Col(ndim,ncol),stat=alloc_stat)
350     if (alloc_stat /= 0 ) then
351     thisStat = 0
352     return
353     endif
354    
355     allocate(u_l_Row(ndim,nrow),stat=alloc_stat)
356     if (alloc_stat /= 0 ) then
357     thisStat = 0
358     return
359     endif
360    
361     allocate(u_l_Col(ndim,ncol),stat=alloc_stat)
362     if (alloc_stat /= 0 ) then
363     thisStat = 0
364     return
365     endif
366    
367     allocate(A_row(9,nrow),stat=alloc_stat)
368     if (alloc_stat /= 0 ) then
369     thisStat = 0
370     return
371     endif
372    
373     allocate(A_Col(9,ncol),stat=alloc_stat)
374     if (alloc_stat /= 0 ) then
375     thisStat = 0
376     return
377     endif
378    
379     allocate(pot_row(nrow),stat=alloc_stat)
380     if (alloc_stat /= 0 ) then
381     thisStat = 0
382     return
383     endif
384    
385     allocate(pot_Col(ncol),stat=alloc_stat)
386     if (alloc_stat /= 0 ) then
387     thisStat = 0
388     return
389     endif
390    
391     allocate(pot_Temp(nlocal),stat=alloc_stat)
392     if (alloc_stat /= 0 ) then
393     thisStat = 0
394     return
395     endif
396    
397     allocate(f_Row(ndim,nrow),stat=alloc_stat)
398     if (alloc_stat /= 0 ) then
399     thisStat = 0
400     return
401     endif
402    
403     allocate(f_Col(ndim,ncol),stat=alloc_stat)
404     if (alloc_stat /= 0 ) then
405     thisStat = 0
406     return
407     endif
408    
409     allocate(f_Temp(ndim,nlocal),stat=alloc_stat)
410     if (alloc_stat /= 0 ) then
411     thisStat = 0
412     return
413     endif
414    
415     allocate(t_Row(ndim,nrow),stat=alloc_stat)
416     if (alloc_stat /= 0 ) then
417     thisStat = 0
418     return
419     endif
420    
421     allocate(t_Col(ndim,ncol),stat=alloc_stat)
422     if (alloc_stat /= 0 ) then
423     thisStat = 0
424     return
425     endif
426    
427     allocate(t_temp(ndim,nlocal),stat=alloc_stat)
428     if (alloc_stat /= 0 ) then
429     thisStat = 0
430     return
431     endif
432    
433     allocate(atid_Row(nrow),stat=alloc_stat)
434     if (alloc_stat /= 0 ) then
435     thisStat = 0
436     return
437     endif
438    
439     allocate(atid_Col(ncol),stat=alloc_stat)
440     if (alloc_stat /= 0 ) then
441     thisStat = 0
442     return
443     endif
444    
445 gezelter 329 allocate(rf_Row(ndim,nrow),stat=alloc_stat)
446     if (alloc_stat /= 0 ) then
447     thisStat = 0
448     return
449     endif
450    
451     allocate(rf_Col(ndim,ncol),stat=alloc_stat)
452     if (alloc_stat /= 0 ) then
453     thisStat = 0
454     return
455     endif
456    
457 gezelter 330 allocate(rf_Temp(ndim,nlocal),stat=alloc_stat)
458     if (alloc_stat /= 0 ) then
459     thisStat = 0
460     return
461     endif
462 gezelter 329
463 gezelter 330
464 gezelter 325 #else
465    
466     allocate(atid(nlocal),stat=alloc_stat)
467     if (alloc_stat /= 0 ) then
468     thisStat = 0
469     return
470     end if
471    
472 gezelter 330 #endif
473    
474 gezelter 329 allocate(rf(ndim,nlocal),stat=alloc_stat)
475     if (alloc_stat /= 0 ) then
476     thisStat = 0
477     return
478     endif
479    
480 gezelter 325 end subroutine setupGlobals
481 gezelter 312
482 gezelter 325 subroutine freeGlobals()
483    
484     !We free in the opposite order in which we allocate in.
485 gezelter 330
486     if (allocated(rf)) deallocate(rf)
487 gezelter 325 #ifdef IS_MPI
488 gezelter 330 if (allocated(rf_Temp)) deallocate(rf_Temp)
489 gezelter 329 if (allocated(rf_Col)) deallocate(rf_Col)
490     if (allocated(rf_Row)) deallocate(rf_Row)
491 gezelter 325 if (allocated(atid_Col)) deallocate(atid_Col)
492     if (allocated(atid_Row)) deallocate(atid_Row)
493     if (allocated(t_Temp)) deallocate(t_Temp)
494     if (allocated(t_Col)) deallocate(t_Col)
495     if (allocated(t_Row)) deallocate(t_Row)
496     if (allocated(f_Temp)) deallocate(f_Temp)
497     if (allocated(f_Col)) deallocate(f_Col)
498     if (allocated(f_Row)) deallocate(f_Row)
499     if (allocated(pot_Temp)) deallocate(pot_Temp)
500     if (allocated(pot_Col)) deallocate(pot_Col)
501     if (allocated(pot_Row)) deallocate(pot_Row)
502     if (allocated(A_Col)) deallocate(A_Col)
503     if (allocated(A_Row)) deallocate(A_Row)
504     if (allocated(u_l_Col)) deallocate(u_l_Col)
505     if (allocated(u_l_Row)) deallocate(u_l_Row)
506     if (allocated(q_Col)) deallocate(q_Col)
507 gezelter 330 if (allocated(q_Row)) deallocate(q_Row)
508     #else
509     if (allocated(atid)) deallocate(atid)
510 gezelter 325 #endif
511    
512     end subroutine freeGlobals
513 gezelter 309
514 mmeineke 270 end module simulation