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

Comparing trunk/mdtools/md_code/f_longRange_module.F90 (file contents):
Revision 185 by chuckv, Fri Nov 22 20:39:33 2002 UTC vs.
Revision 194 by chuckv, Wed Dec 4 21:19:38 2002 UTC

# Line 1 | Line 1
1   module long_range_forces
2    use definitions, ONLY : dp, ndim
3    use simulation
4 + #ifdef MPI
5 +  use mpi_force_module
6 + #endif
7    implicit none
8 +  PRIVATE
9  
10 +  type, public :: ForceSelecter
11 +     sequence
12 +    
13 +     logical :: is_dipole
14 +     logical :: is_vdw
15 +     logical :: is_lj
16 +     logical :: is_ssd
17  
18 +  end type ForceSelecter
19  
20  
9 –  real(kind = dp ) :: rlistsq
10 –  real(kind = dp ) :: rcutsq
21  
22  
23  
24 + contains
25  
15 –  integer, allocatable, dimension(:)              :: point
16 –  integer, allocatable, dimension(:)              :: list
26  
27 +  subroutine long_range_force(nComponents,nDim,q,rField,t,rF,A,TotPE)
28  
29 + !-------------------- arguments-------------------------->
30 + !   integer components
31 +    integer, intent(in) :: nComponents  ! number of Particles
32 +    integer, intent(in) :: nDim         ! number of dimensions
33  
34 + !   real scalars
35 +    real( kind = dp ), intent(out), optional :: TotPE   ! potential energy
36  
37 < contains
37 > !   real vectors
38 >    real( kind = dp ),dimension(nDim,nComponents), &
39 >         intent(inout) :: q                   ! position vector
40 >    real( kind = dp ),dimension(nDim,nComponents), &
41 >         intent(inout) :: f                   ! force vector
42 >    real( kind = dp ),dimension(nDim,nComponents), &
43 >         intent(inout) :: t                   ! torque vector
44 >    real( kind = dp ),dimension(nDim,nComponents), &
45 >         intent(inout) :: rField              ! reaction field vector
46 >    real( kind = dp ),dimension(nDim,nDim), &
47 >         intent(inout) :: A                   ! rotation matrix
48 > !------------------- end arguments------------------------->
49  
50  
51 + !------------------- local variables----------------------->
52 +    real( kind = dp ), dimension(nDim)  :: q_i
53 +    real( kind = dp ), dimension(nDim)  :: q_ij
54 +    
55  
56 +    logical :: do_pot = .false.
57  
26 –  subroutine long_range_force ( check_exclude, &
27 –       pair_i, pair_j, n_pairs, natoms, sigma, epslon, &
28 –       box, &
29 –       rx, ry, rz, fx, fy, fz, &
30 –       v, &
31 –       update, point, list, maxnab, &
32 –       rx0, ry0, rz0, &
33 –       ux, uy, uz, &
34 –       tx, ty, tz, &
35 –       ex, ey, ez, &
36 –       mu, rrf, rtaper, &
37 –       dielectric, &
38 –       is_dipole, is_vdw, is_lj, is_ssd, &
39 –       Axx, Axy, Axz, &
40 –       Ayx, Ayy, Ayz, &
41 –       Azx, Azy, Azz )
58  
59 <    implicit none
59 >
60 > !!$
61 > !!$  subroutine long_range_force ( check_exclude, &
62 > !!$       pair_i, pair_j, n_pairs, natoms, sigma, epslon, &
63 > !!$       box, &
64 > !!$       rx, ry, rz, fx, fy, fz, &
65 > !!$       v, &
66 > !!$       update, point, list, maxnab, &
67 > !!$       rx0, ry0, rz0, &
68 > !!$       ux, uy, uz, &
69 > !!$       tx, ty, tz, &
70 > !!$       ex, ey, ez, &
71 > !!$       mu, rrf, rtaper, &
72 > !!$       dielectric, &
73 > !!$       is_dipole, is_vdw, is_lj, is_ssd, &
74 > !!$       Axx, Axy, Axz, &
75 > !!$       Ayx, Ayy, Ayz, &
76 > !!$       Azx, Azy, Azz )
77 >
78 >
79      
80      ! This routine implements the Verlet neighbor list for all of
81      ! the long range force calculations. After a pair has been
82      ! determined to lie within the neighbor list, the pair is passed
83      ! to the appropriate force functions.
84      
85 <    
86 <    
87 <    ! Passed Arguments
88 <    
89 <    integer :: n_pairs  ! the number of excluded pairs
90 <    integer :: natoms   ! the number of atoms in the system
91 <    integer :: maxnab   ! the max number of neighbors for the verlet list
92 <    
93 <    logical(kind=2) :: check_exclude ! boolean to check for exclusion pairs
94 <    logical(kind=2) :: update        ! boolean to update the neighbor list
95 <    
96 <    double precision  rcut                ! the VDW/LJ cutoff radius
97 <    double precision  rlist               ! Verlet list cutoff
98 <    double precision  box_x, box_y, box_z ! periodic boundry conditions
99 <    double precision  rrf                 ! dipole reaction field cutoff
100 <    double precision  rtaper              ! the taper radius of the rxn field
101 <    double precision  dielectric          ! the dielectric cnst. of the rxn field
102 <    double precision  v                   ! the potential energy
103 <    
104 <    
105 <    
106 <    ! Passed Arrays
107 <    
108 <    integer, dimension(n_pairs)::  pair_i, pair_j  ! the excluded pairs i and j
109 <    integer, dimension(natoms) ::  point           ! the Verlet list ptr array
110 <    integer, dimension(maxnab) ::  list            ! the verlet list
111 <
112 <    logical(kind=2), dimension(natoms) :: is_dipole  ! dipole boolean array
113 <    logical(kind=2), dimension(natoms) :: is_vdw     ! VDW boolean array
114 <    logical(kind=2), dimension(natoms) :: is_lj      ! LJ boolean array
115 <    logical(kind=2), dimension(natoms) :: is_ssd     ! ssd boolean array
116 <
117 <    double precision, dimension(natoms) :: sigma         ! VDW/LJ distance prmtr.
118 <    double precision, dimension(natoms) :: epslon        ! VDW/LJ well depth prmtr.
119 <    double precision, dimension(natoms) :: rx, ry, rz    ! positions
120 <    double precision, dimension(natoms) :: fx, fy, fz    ! forces
121 <    double precision, dimension(natoms) :: rx0, ry0, rz0 ! intial verlet positions
122 <    double precision, dimension(natoms) :: ux, uy, uz    ! dipole unit vectors
123 <    double precision, dimension(natoms) :: tx, ty, tz    ! torsions
124 <    double precision, dimension(natoms) :: ex, ey, ez    ! reacion field
125 <    double precision, dimension(natoms) :: mu            ! dipole moments
126 <    double precision, dimension(natoms) :: Axx, Axy, Axz ! rotational array
127 <    double precision, dimension(natoms) :: Ayx, Ayy, Ayz !
128 <    double precision, dimension(natoms) :: Azx, Azy, Azz !
129 <
130 <    ! local variables
131 <
132 <    double precision rxi, ryi, rzi
133 <    double precision rcutsq, rlstsq, rrfsq, rijsq, rij
134 <    double precision rxij, ryij, rzij
135 <    double precision prerf ! a reaction field preterm
136 <
137 <    double precision, dimension(9,natoms) :: A
138 <
139 <    integer :: i, j, k, ix, jx, nlist
140 <    integer :: jbeg, jend, jnab
141 <
142 <    logical :: exclude_temp1, exclude_temp2, exclude
85 > !!$    
86 > !!$    
87 > !!$    ! Passed Arguments
88 > !!$    
89 > !!$    integer :: n_pairs  ! the number of excluded pairs
90 > !!$    integer :: natoms   ! the number of atoms in the system
91 > !!$    integer :: maxnab   ! the max number of neighbors for the verlet list
92 > !!$    
93 > !!$    logical(kind=2) :: check_exclude ! boolean to check for exclusion pairs
94 > !!$    logical(kind=2) :: update        ! boolean to update the neighbor list
95 > !!$    
96 > !!$    double precision  rcut                ! the VDW/LJ cutoff radius
97 > !!$    double precision  rlist               ! Verlet list cutoff
98 > !!$    double precision  box_x, box_y, box_z ! periodic boundry conditions
99 > !!$    double precision  rrf                 ! dipole reaction field cutoff
100 > !!$    double precision  rtaper              ! the taper radius of the rxn field
101 > !!$    double precision  dielectric          ! the dielectric cnst. of the rxn field
102 > !!$    double precision  v                   ! the potential energy
103 > !!$    
104 > !!$    
105 > !!$    
106 > !!$    ! Passed Arrays
107 > !!$    
108 > !!$    integer, dimension(n_pairs)::  pair_i, pair_j  ! the excluded pairs i and j
109 > !!$    integer, dimension(natoms) ::  point           ! the Verlet list ptr array
110 > !!$    integer, dimension(maxnab) ::  list            ! the verlet list
111 > !!$
112 > !!$    logical(kind=2), dimension(natoms) :: is_dipole  ! dipole boolean array
113 > !!$    logical(kind=2), dimension(natoms) :: is_vdw     ! VDW boolean array
114 > !!$    logical(kind=2), dimension(natoms) :: is_lj      ! LJ boolean array
115 > !!$    logical(kind=2), dimension(natoms) :: is_ssd     ! ssd boolean array
116 > !!$
117 > !!$    double precision, dimension(natoms) :: sigma         ! VDW/LJ distance prmtr.
118 > !!$    double precision, dimension(natoms) :: epslon        ! VDW/LJ well depth prmtr.
119 > !!$    double precision, dimension(natoms) :: rx, ry, rz    ! positions
120 > !!$    double precision, dimension(natoms) :: fx, fy, fz    ! forces
121 > !!$    double precision, dimension(natoms) :: rx0, ry0, rz0 ! intial verlet positions
122 > !!$    double precision, dimension(natoms) :: ux, uy, uz    ! dipole unit vectors
123 > !!$    double precision, dimension(natoms) :: tx, ty, tz    ! torsions
124 > !!$    double precision, dimension(natoms) :: ex, ey, ez    ! reacion field
125 > !!$    double precision, dimension(natoms) :: mu            ! dipole moments
126 > !!$    double precision, dimension(natoms) :: Axx, Axy, Axz ! rotational array
127 > !!$    double precision, dimension(natoms) :: Ayx, Ayy, Ayz !
128 > !!$    double precision, dimension(natoms) :: Azx, Azy, Azz !
129 > !!$
130 > !!$    ! local variables
131 > !!$
132 > !!$    double precision rxi, ryi, rzi
133 > !!$    double precision rcutsq, rlstsq, rrfsq, rijsq, rij
134 > !!$    double precision rxij, ryij, rzij
135 > !!$    double precision prerf ! a reaction field preterm
136 > !!$
137 > !!$    double precision, dimension(9,natoms) :: A
138 > !!$
139 > !!$    integer :: i, j, k, ix, jx, nlist
140 > !!$    integer :: jbeg, jend, jnab
141 > !!$
142 > !!$    logical :: exclude_temp1, exclude_temp2, exclude
143      
144  
145      !*******************************************************************
146 <    
147 <    do i=1, natoms
148 <       A(1,i) = Axx(i)
149 <       A(2,i) = Axy(i)
150 <       A(3,i) = Axz(i)
151 <      
152 <       A(4,i) = Ayx(i)
153 <       A(5,i) = Ayy(i)
154 <       A(6,i) = Ayz(i)
155 <      
156 <       A(7,i) = Azx(i)
157 <       A(8,i) = Azy(i)
158 <       A(9,i) = Azz(i)
159 <    end do
146 > !!$    
147 > !!$    do i=1, natoms
148 > !!$       A(1,i) = Axx(i)
149 > !!$       A(2,i) = Axy(i)
150 > !!$       A(3,i) = Axz(i)
151 > !!$      
152 > !!$       A(4,i) = Ayx(i)
153 > !!$       A(5,i) = Ayy(i)
154 > !!$       A(6,i) = Ayz(i)
155 > !!$      
156 > !!$       A(7,i) = Azx(i)
157 > !!$       A(8,i) = Azy(i)
158 > !!$       A(9,i) = Azz(i)
159 > !!$    end do
160  
161  
162 +
163 +    if (present(TotPE)) do_pot = .true.
164 +
165      rlstsq = rlist * rlist
166      rcutsq = rcut * rcut
167      rrfsq  = rrf * rrf
# Line 131 | Line 169 | contains
169      prerf = 14.39257d0 * 2.0d0 * ( dielectric - 1.0d0 ) / &
170           ( ( 2.0d0 * dielectric + 1.0d0 ) * rrfsq * rrf )
171  
134 –    ! ** zero forces **
172  
173 <    do I = 1, natoms
174 <
175 <       fx(i) = 0.0d0
176 <       fy(i) = 0.0d0
140 <       fz(i) = 0.0d0
173 >    ! ** zero forces **
174 >    f  = 0.0_dp
175 >    t  = 0.0_dp
176 >    rF = 0.0_dp
177  
178 <       tx(i) = 0.0d0
143 <       ty(i) = 0.0d0
144 <       tz(i) = 0.0d0
145 <
146 <       ex(i) = 0.0d0
147 <       ey(i) = 0.0d0
148 <       ez(i) = 0.0d0
178 >    TotPE = 0.0_dp
179  
150 –    end do
180  
152 –    V  = 0.0d0
181  
182      if ( update ) then
183  
# Line 165 | Line 193 | contains
193  
194            point(i) = nlist + 1
195  
196 +          q_i(1:nDim) = q(1:nDim,i)
197 +
198            rxi      = rx(i)
199            ryi      = ry(i)
200            rzi      = rz(i)
201  
202            do j = i + 1, natoms
203  
204 +             q_ij(1:nDim) = q(1:nDim,j) - q_i(1:nDim)
205               rxij  = rx(j) - rxi
206               ryij  = ry(j) - ryi
207               rzij  = rz(j) - rzi

Diff Legend

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