ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/mdtools/md_code/simulation_module.F90
(Generate patch)

Comparing trunk/mdtools/md_code/simulation_module.F90 (file contents):
Revision 222 by chuckv, Thu Jan 2 21:45:45 2003 UTC vs.
Revision 284 by chuckv, Tue Feb 25 21:30:09 2003 UTC

# Line 1 | Line 1 | module simulation
1   module simulation
2    use definitions, ONLY :dp
3 <  use force_wrappers, ONLY : alloc_force_wrappers
3 > #ifdef IS_MPI
4 >  use mpiSimulation
5 > #endif
6  
7    implicit none
8    PRIVATE
9  
10 + #define __FORTRAN90
11 + #include "../headers/fsimulation.h"
12  
13 <  integer :: nLR
13 >  type (simtype), public :: thisSim
14 > !! Tag for MPI calculations  
15 >  integer, allocatable, dimension(:) :: tag
16  
17 + #ifdef IS_MPI
18 +  integer, allocatable, dimension(:) :: tag_row
19 +  integer, allocatable, dimension(:) :: tag_column
20 + #endif
21  
22 <  real ( kind = dp ), dimension(3) :: box
22 > !! WARNING: use_pbc hardcoded, fixme
23 >   logical :: setSim = .false.
24  
25 + !! array for saving previous positions for neighbor lists.  
26 +  real( kind = dp ), allocatable,dimension(:,:),save :: q0
27  
15 –  real ( kind = dp ), public :: rlist
16 –  real ( kind = dp ), public :: rcut
17 –  real ( kind = dp ), public :: rlistsq
18 –  real ( kind = dp ), public :: rcutsq
28  
20 –  integer,public, allocatable, dimension(:)   :: point
21 –  integer,public, allocatable, dimension(:)   :: list
22 –
23 –
24 –
25 –
26 –
27 –  public :: check
28 –  public :: save_nlist
29    public :: wrap
30    public :: getBox
31 +  public :: getRcut
32 +  public :: getRlist
33 +  public :: getNlocal
34 +  public :: setSimulation
35 +  public :: isEnsemble
36 +  public :: isPBC
37 +  public :: getStringLen
38 +  public :: returnMixingRules
39  
40 + !  public :: setRcut
41 +
42    interface wrap
43       module procedure wrap_1d
44       module procedure wrap_3d
# Line 40 | Line 50 | module simulation
50    end interface
51  
52  
43 – !MPI dependent routines
53  
45 – #ifdef IS_MPI
46 – ! Universal routines: All types of force calculations will need these arrays
47 – ! Arrays specific to a type of force calculation should be declared in that module.
48 –  real( kind = dp ), allocatable, dimension(:,:) :: qRow
49 –  real( kind = dp ), allocatable, dimension(:,:) :: qColumn
54  
51 –  real( kind = dp ), allocatable, dimension(:,:) :: fRow
52 –  real( kind = dp ), allocatable, dimension(:,:) :: fColumn
55  
54 –  real( kind = dp ), allocatable, dimension(:,:) :: tRow
55 –  real( kind = dp ), allocatable, dimension(:,:) :: tColumn
56  
57 – #endif
57  
58  
59  
60  
62 –
63 –
61   contains
62  
63 < #ifdef MPI
64 < ! Allocated work arrays for MPI
65 <  subroutine allocate_mpi_arrays(nDimensions,numComponents)
66 <    integer, intent(in) :: nDimensions
67 <    integer, intent(in) :: numComponents
63 >  subroutine setSimulation(nLRParticles,box,rlist,rcut,ensemble,mixingRule,use_pbc)
64 >    integer, intent(in) :: nLRParticles
65 >    real(kind = dp ), intent(in), dimension(3) :: box
66 >    real(kind = dp ), intent(in) :: rlist
67 >    real(kind = dp ), intent(in) :: rcut
68 >    character( len = stringLen), intent(in)  :: ensemble
69 >    character( len = stringLen), intent(in)  :: mixingRule
70 >    logical, intent(in) :: use_pbc
71 >    integer :: alloc_stat
72 >    if( setsim ) return  ! simulation is already initialized
73 >    setSim = .true.
74  
75 <
76 <
77 <
78 <
79 <  end subroutine allocate_mpi_arrays
80 < #endif
81 <
79 <  subroutine set_simulation(box,rlist,rcut)
75 >    thisSim%nLRParticles = nLRParticles
76 >    thisSim%box          = box
77 >    thisSim%rlist        = rlist
78 >    thisSIm%rlistsq      = rlist * rlist
79 >    thisSim%rcut         = rcut
80 >    thisSim%rcutsq       = rcut * rcut
81 >    thisSim%rcut6        = thisSim%rcutsq * thisSim%rcutsq * thisSim%rcutsq
82      
83 +    thisSim%ensemble = ensemble
84 +    thisSim%mixingRule = mixingRule
85 +    thisSim%use_pbc = use_pbc
86  
87 +    if (.not. allocated(q0)) then
88 +       allocate(q0(3,nLRParticles),stat=alloc_stat)
89 +    endif
90 +  end subroutine setSimulation
91  
92 <  end subroutine set_simulation
92 >  function getNparticles() result(nparticles)
93 >    integer :: nparticles
94 >    nparticles = thisSim%nLRparticles
95 >  end function getNparticles
96  
97  
86 –
98    subroutine change_box_size(new_box_size)
99      real(kind=dp), dimension(3) :: new_box_size
100  
101 <    box = new_box_size
101 >    thisSim%box = new_box_size
102  
103    end subroutine change_box_size
104  
105  
106 <  elemental function getBox_3d() result(thisBox)
106 >  function getBox_3d() result(thisBox)
107      real( kind = dp ), dimension(3) :: thisBox
108 <    thisBox = box
108 >    thisBox = thisSim%box
109    end function getBox_3d
110  
111    function getBox_dim(dim) result(thisBox)
112      integer, intent(in) :: dim
113      real( kind = dp ) :: thisBox
114  
115 <    thisBox = box(dim)
115 >    thisBox = thisSim%box(dim)
116    end function getBox_dim
106 –    
107 –  subroutine check(update_nlist)
117    
109 –    integer :: i
110 –    real( kind = DP ) :: dispmx
111 –    logical, intent(out) :: update_nlist
112 –    real( kind = DP ) :: dispmx_tmp
113 –    
114 –    dispmx = 0.0E0_DP
115 –    
116 –    !! calculate the largest displacement of any atom in any direction
117 –    
118 – #ifdef MPI
119 –    dispmx_tmp = 0.0E0_DP
120 –    do i = 1, nlocal
121 –       dispmx_tmp = max( abs ( q(1,i) - q0(1,i) ), dispmx )
122 –       dispmx_tmp = max( abs ( q(2,i) - q0(2,i) ), dispmx )
123 –       dispmx_tmp = max( abs ( q(3,i) - q0(3,i) ), dispmx )
124 –    end do
125 –    call mpi_allreduce(dispmx_tmp,dispmx,1,mpi_double_precision, &
126 –       mpi_max,mpi_comm_world,mpi_err)
127 – #else
128 –    do i = 1, natoms
129 –       dispmx = max( abs ( q(1,i) - q0(1,i) ), dispmx )
130 –       dispmx = max( abs ( q(2,i) - q0(2,i) ), dispmx )
131 –       dispmx = max( abs ( q(3,i) - q0(3,i) ), dispmx )
132 –    end do
133 – #endif
134 –  
135 –    !! a conservative test of list skin crossings
136 –    dispmx = 2.0E0_DP * sqrt (3.0E0_DP * dispmx * dispmx)
137 –    
138 –    update_nlist = (dispmx.gt.(skin_thickness))
139 –    
140 –  end subroutine check
141 –  
142 –  subroutine save_nlist()
143 –    integer :: i
144 – #ifdef MPI
145 –    do i = 1, nlocal
146 – #else
147 –       do i = 1, natoms
148 – #endif
149 –          q0(1,i) = q(1,i)
150 –          q0(2,i) = q(2,i)
151 –          q0(3,i) = q(3,i)
152 –       end do
118  
154 –  end subroutine save_nlist
155 –  
156 –
119    function wrap_1d(r,dim) result(this_wrap)
120      
121      
# Line 163 | Line 125 | contains
125      
126      if (use_pbc) then
127         !     this_wrap = r - box(dim)*dsign(1.0E0_DP,r)*int(abs(r/box(dim)) + 0.5E0_DP)
128 <       this_wrap = r - box(dim)*nint(r/box(dim))
128 >       this_wrap = r - thisSim%box(dim)*nint(r/thisSim%box(dim))
129      else
130         this_wrap = r
131      endif
# Line 171 | Line 133 | contains
133      return
134    end function wrap_1d
135  
136 <  elemental function wrap_3d(r) result(this_wrap)
136 >  function wrap_3d(r) result(this_wrap)
137      real( kind = dp ), dimension(3), intent(in) :: r
138      real( kind = dp ), dimension(3) :: this_wrap
139  
140      
141 <    if (use_pbc) then
141 >    if (this_sim%use_pbc) then
142         !     this_wrap = r - box(dim)*dsign(1.0E0_DP,r)*int(abs(r/box(dim)) + 0.5E0_DP)
143 <       this_wrap(1:3) = r(1:3) - box(1:3)*nint(r(1:3)/box(1:3))
143 >       this_wrap = r - thisSim%box*nint(r/thisSim%box)
144      else
145         this_wrap = r
146      endif
147    end function wrap_3d
148  
149    
150 +
151 +  subroutine getRcut(thisrcut,rcut2,rcut6,status)
152 +    real( kind = dp ), intent(out) :: thisrcut
153 +    real( kind = dp ), intent(out), optional :: rcut2
154 +    real( kind = dp ), intent(out), optional :: rcut6
155 +    integer, optional :: status
156 +
157 +    if (present(status)) status = 0
158 +
159 +    if (.not.setSim ) then
160 +       if (present(status)) status = -1
161 +       return
162 +    end if
163 +    
164 +    thisrcut = thisSim%rcut
165 +    if(present(rcut2)) rcut2 = thisSim%rcutsq
166 +    if(present(rcut6)) rcut6 = thisSim%rcut6
167 +
168 +  end subroutine getRcut
169 +  
170 +  
171 +  
172 +
173 +  subroutine getRlist(thisrlist,rlist2,status)
174 +    real( kind = dp ), intent(out) :: thisrlist
175 +    real( kind = dp ), intent(out), optional :: rlist2
176 +
177 +    integer, optional :: status
178 +
179 +    if (present(status)) status = 0
180 +
181 +    if (.not.setSim ) then
182 +       if (present(status)) status = -1
183 +       return
184 +    end if
185 +    
186 +    thisrlist = thisSim%rlist
187 +    if(present(rlist2)) rlist2 = thisSim%rlistsq
188 +
189 +
190 +  end subroutine getRlist
191 +  
192 +  
193 +
194 + pure function getNlocal() result(nlocal)
195 +    integer :: nlocal
196 +    nlocal = thisSim%nLRparticles
197 +  end function getNlocal
198 +
199 +
200 +  function isEnsemble(this_ensemble) result(is_this_ensemble)
201 +    character(len = *) :: this_ensemble
202 +    logical :: is_this_enemble
203 +    is_this_ensemble = .false.
204 +    if (this_ensemble == thisSim%ensemble) is_this_ensemble = .true.
205 +  end function isEnsemble
206 +
207 +  function returnEnsemble() result(thisEnsemble)
208 +    character (len = len(thisSim%ensemble)) :: thisEnsemble
209 +    thisEnsemble = thisSim%ensemble
210 +  end function returnEnsemble
211 +
212 +  function returnMixingRules() result(thisMixingRule)
213 +    character (len = len(thisSim%ensemble)) :: thisMixingRule
214 +    thisMixingRule = thisSim%MixingRule
215 +  end function returnMixingRules
216 +
217 +  function isPBC() result(PBCset)
218 +    logical :: PBCset
219 +    PBCset = .false.
220 +    if (thisSim%use_pbc) PBCset = .true.
221 +  end function isPBC
222 +
223 +  pure function getStringLen() result (thislen)
224 +    integer :: thislen    
225 +    thislen = string_len
226 +  end function setStringLen
227 +
228   end module simulation

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines