ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE_old/src/mdtools/libmdCode/simulation_module.F90
Revision: 332
Committed: Thu Mar 13 15:28:43 2003 UTC (21 years, 4 months ago) by gezelter
File size: 14067 byte(s)
Log Message:
initialization routines, bug fixes, exclude lists

File Contents

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