ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE_old/src/mdtools/libmdCode/do_Forces.F90
(Generate patch)

Comparing trunk/OOPSE_old/src/mdtools/libmdCode/do_Forces.F90 (file contents):
Revision 295 by chuckv, Thu Mar 6 19:57:03 2003 UTC vs.
Revision 298 by chuckv, Fri Mar 7 18:26:30 2003 UTC

# Line 1 | Line 1
1   !! Calculates Long Range forces.
2   !! @author Charles F. Vardeman II
3   !! @author Matthew Meineke
4 < !! @version $Id: do_Forces.F90,v 1.3 2003-03-06 19:57:03 chuckv Exp $, $Date: 2003-03-06 19:57:03 $, $Name: not supported by cvs2svn $, $Revision: 1.3 $
4 > !! @version $Id: do_Forces.F90,v 1.5 2003-03-07 18:26:30 chuckv Exp $, $Date: 2003-03-07 18:26:30 $, $Name: not supported by cvs2svn $, $Revision: 1.5 $
5  
6  
7  
8   module do_Forces
9    use simulation
10    use definitions
11 <  use generic_atypes
11 >  use forceGlobals
12 >  use atype_typedefs
13    use neighborLists
14 +
15    
16    use lj_FF
17    use sticky_FF
# Line 22 | Line 24 | module do_Forces
24    implicit none
25    PRIVATE
26  
27 < !! Number of lj_atypes in lj_atype_list
26 <  integer, save :: n_atypes = 0
27 > public :: do_force_loop
28  
28 – !! Global list of lj atypes in simulation
29 –  type (atype), pointer :: ListHead => null()
30 –  type (atype), pointer :: ListTail => null()
31 –
32 –
33 –
34 –
35 –  logical, save :: firstTime = .True.
36 –
37 – !! Atype identity pointer lists
38 – #ifdef IS_MPI
39 – !! Row lj_atype pointer list
40 –  type (identPtrList), dimension(:), pointer :: identPtrListRow => null()
41 – !! Column lj_atype pointer list
42 –  type (identPtrList), dimension(:), pointer :: identPtrListColumn => null()
43 – #else
44 –  type(identPtrList ), dimension(:), pointer :: identPtrList => null()
45 – #endif
46 –
47 –
48 – !! Logical has lj force field module been initialized?
49 –  logical, save :: isFFinit = .false.
50 –
51 – !! Use periodic boundry conditions
52 –  logical :: wrap = .false.
53 –
54 – !! Potential energy global module variables
55 – #ifdef IS_MPI
56 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: qRow = 0.0_dp
57 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: qCol = 0.0_dp
58 –
59 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: muRow = 0.0_dp
60 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: muCol = 0.0_dp
61 –
62 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: u_lRow = 0.0_dp
63 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: u_lCol = 0.0_dp
64 –
65 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: ARow = 0.0_dp
66 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: ACol = 0.0_dp
67 –
68 –  
69 –
70 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: fRow = 0.0_dp
71 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: fCol = 0.0_dp
72 –
73 –
74 –  real(kind = dp), dimension(3,getNrow(plan_row)) :: tRow = 0.0_dp
75 –  real(kind = dp), dimension(3,getNcol(plan_col)) :: tCol = 0.0_dp
76 –
77 –
78 –
79 –  real(kind = dp), dimension(getNrow(plan_row)) :: eRow = 0.0_dp
80 –  real(kind = dp), dimension(getNcol(plan_col)) :: eCol = 0.0_dp
81 –
82 –  real(kind = dp), dimension(getNlocal()) :: eTemp = 0.0_dp
83 – #endif
84 –  real(kind = dp) :: pe = 0.0_dp
85 –  real(kind = dp), dimension(3,getNlocal()) :: fTemp = 0.0_dp
86 –  real(kind = dp), dimension(3,getNlocal()) :: tTemp = 0.0_dp
87 –  real(kind = dp), dimension(9) :: tauTemp = 0.0_dp
88 –
89 –  logical :: do_preForce  = .false.
90 –  logical :: do_postForce = .false.
91 –
92 –
93 –
94 – !! Public methods and data
95 –  public :: new_atype
96 –  public :: do_forceLoop
97 –  public :: init_FF
98 –
99 –  
100 –
101 –
29   contains
103 –
104 – !! Adds a new lj_atype to the list.
105 –  subroutine new_atype(ident,mass,epsilon,sigma, &
106 –       is_LJ,is_Sticky,is_DP,is_GB,w0,v0,dipoleMoment,status)
107 –    real( kind = dp ), intent(in) :: mass
108 –    real( kind = dp ), intent(in) :: epsilon
109 –    real( kind = dp ), intent(in) :: sigma
110 –    real( kind = dp ), intent(in) :: w0
111 –    real( kind = dp ), intent(in) :: v0
112 –    real( kind = dp ), intent(in) :: dipoleMoment
113 –
114 –    integer, intent(in) :: ident
115 –    integer, intent(out) :: status
116 –    integer, intent(in) :: is_Sticky
117 –    integer, intent(in) :: is_DP
118 –    integer, intent(in) :: is_GB
119 –    integer, intent(in) :: is_LJ
120 –
121 –
122 –    type (atype), pointer :: the_new_atype
123 –    integer :: alloc_error
124 –    integer :: atype_counter = 0
125 –    integer :: alloc_size
126 –    integer :: err_stat
127 –    status = 0
128 –
129 –
130 –
131 – ! allocate a new atype    
132 –    allocate(the_new_atype,stat=alloc_error)
133 –    if (alloc_error /= 0 ) then
134 –       status = -1
135 –       return
136 –    end if
137 –
138 – ! assign our new atype information
139 –    the_new_atype%mass        = mass
140 –    the_new_atype%epsilon     = epsilon
141 –    the_new_atype%sigma       = sigma
142 –    the_new_atype%sigma2      = sigma * sigma
143 –    the_new_atype%sigma6      = the_new_atype%sigma2 * the_new_atype%sigma2 &
144 –         * the_new_atype%sigma2
145 –    the_new_atype%w0       = w0
146 –    the_new_atype%v0       = v0
147 –    the_new_atype%dipoleMoment       = dipoleMoment
148 –
149 –    
150 – ! assume that this atype will be successfully added
151 –    the_new_atype%atype_ident = ident
152 –    the_new_atype%atype_number = n_lj_atypes + 1
153 –    
154 –    if ( is_Sticky /= 0 )    the_new_atype%is_Sticky   = .true.
155 –    if ( is_GB /= 0 )        the_new_atype%is_GB       = .true.
156 –    if ( is_LJ /= 0 )        the_new_atype%is_LJ       = .true.
157 –    if ( is_DP /= 0 )        the_new_atype%is_DP       = .true.
158 –
159 –    call add_atype(the_new_atype,ListHead,ListTail,err_stat)
160 –    if (err_stat /= 0 ) then
161 –       status = -1
162 –       return
163 –    endif
164 –
165 –    n_atypes = n_atypes + 1
166 –
167 –
168 –  end subroutine new_atype
169 –
30  
171 –  subroutine init_FF(nComponents,ident, status)
172 – !! Number of components in ident array
173 –    integer, intent(inout) :: nComponents
174 – !! Array of identities nComponents long corresponding to
175 – !! ljatype ident.
176 –    integer, dimension(nComponents),intent(inout) :: ident
177 – !!  Result status, success = 0, error = -1
178 –    integer, intent(out) :: Status
179 –
180 –    integer :: alloc_stat
181 –
182 –    integer :: thisStat
183 –    integer :: i
184 –
185 –    integer :: myNode
186 – #ifdef IS_MPI
187 –    integer, allocatable, dimension(:) :: identRow
188 –    integer, allocatable, dimension(:) :: identCol
189 –    integer :: nrow
190 –    integer :: ncol
191 – #endif
192 –    status = 0
193 –  
194 –
195 –    
196 –
197 – !! if were're not in MPI, we just update ljatypePtrList
198 – #ifndef IS_MPI
199 –    call create_IdentPtrlst(ident,ListHead,identPtrList,thisStat)
200 –    if ( thisStat /= 0 ) then
201 –       status = -1
202 –       return
203 –    endif
204 –
205 –    
206 – ! if were're in MPI, we also have to worry about row and col lists    
207 – #else
208 –  
209 – ! We can only set up forces if mpiSimulation has been setup.
210 –    if (.not. isMPISimSet()) then
211 –       write(default_error,*) "MPI is not set"
212 –       status = -1
213 –       return
214 –    endif
215 –    nrow = getNrow(plan_row)
216 –    ncol = getNcol(plan_col)
217 –    mynode = getMyNode()
218 – !! Allocate temperary arrays to hold gather information
219 –    allocate(identRow(nrow),stat=alloc_stat)
220 –    if (alloc_stat /= 0 ) then
221 –       status = -1
222 –       return
223 –    endif
224 –
225 –    allocate(identCol(ncol),stat=alloc_stat)
226 –    if (alloc_stat /= 0 ) then
227 –       status = -1
228 –       return
229 –    endif
230 –
231 – !! Gather idents into row and column idents
232 –
233 –    call gather(ident,identRow,plan_row)
234 –    call gather(ident,identCol,plan_col)
235 –    
236 –  
237 – !! Create row and col pointer lists
238 –  
239 –    call create_IdentPtrlst(identRow,ListHead,identPtrListRow,thisStat)
240 –    if (thisStat /= 0 ) then
241 –       status = -1
242 –       return
243 –    endif
244 –  
245 –    call create_IdentPtrlst(identCol,ListHead,identPtrListColumn,thisStat)
246 –    if (thisStat /= 0 ) then
247 –       status = -1
248 –       return
249 –    endif
250 –
251 – !! free temporary ident arrays
252 –    if (allocated(identCol)) then
253 –       deallocate(identCol)
254 –    end if
255 –    if (allocated(identCol)) then
256 –       deallocate(identRow)
257 –    endif
258 –
259 – #endif
260 –    
261 –    call initForce_Modules(thisStat)
262 –    if (thisStat /= 0) then
263 –       status = -1
264 –       return
265 –    endif
266 –
267 – !! Create neighbor lists
268 –    call expandList(thisStat)
269 –    if (thisStat /= 0) then
270 –       status = -1
271 –       return
272 –    endif
273 –
274 –    isFFinit = .true.
275 –
276 –
277 –  end subroutine init_FF
278 –
279 –
280 –
281 –
282 –  subroutine initForce_Modules(thisStat)
283 –    integer, intent(out) :: thisStat
284 –    integer :: my_status
285 –    
286 –    thisStat = 0
287 –    call init_lj_FF(ListHead,my_status)
288 –    if (my_status /= 0) then
289 –       thisStat = -1
290 –       return
291 –    end if
292 –
293 –  end subroutine initForce_Modules
294 –
295 –
296 –
297 –
31   !! FORCE routine Calculates Lennard Jones forces.
32   !------------------------------------------------------------->
33 <  subroutine do__force_loop(q,A,mu,u_l,f,t,tau,potE,do_pot,FFerror)
33 >  subroutine do_force_loop(q,A,mu,u_l,f,t,tau,potE,do_pot,FFerror)
34   !! Position array provided by C, dimensioned by getNlocal
35      real ( kind = dp ), dimension(3,getNlocal()) :: q
36    !! Rotation Matrix for each long range particle in simulation.
# Line 393 | Line 126 | contains
126    if (isPBC()) wrap = .true.
127  
128  
396 – #ifndef IS_MPI
397 –  nrow = natoms - 1
398 –  ncol = natoms
399 – #else
400 –  nrow = getNrow(plan_row)
401 –  ncol = getNcol(plan_col)
402 –  nlocal = natoms
403 –  j_start = 1
404 – #endif
129  
130    
131   !! See if we need to update neighbor lists
# Line 431 | Line 155 | contains
155   #endif
156  
157  
158 <  if (update_nlist) then
158 > #ifdef IS_MPI
159 >    
160 >    if (update_nlist) then
161 >      
162 >       ! save current configuration, contruct neighbor list,
163 >       ! and calculate forces
164 >       call save_neighborList(q)
165 >      
166 >       neighborListSize = getNeighborListSize()
167 >       nlist = 0
168 >      
169 >       nrow = getNrow(plan_row)
170 >       ncol = getNcol(plan_col)
171 >       nlocal = getNlocal()
172 >      
173 >       do i = 1, nrow
174 >          point(i) = nlist + 1
175 >          Atype_i => identPtrListRow(i)%this
176 >          
177 >          inner: do j = 1, ncol
178 >             Atype_j => identPtrListColumn(j)%this
179 >            
180 >             call get_interatomic_vector(i,j,qRow(:,i),qCol(:,j),&
181 >                  rxij,ryij,rzij,rijsq,r)
182 >            
183 >             ! skip the loop if the atoms are identical
184 >             if (mpi_cycle_jLoop(i,j)) cycle inner:
185 >            
186 >             if (rijsq <  rlistsq) then            
187 >                
188 >                nlist = nlist + 1
189 >                
190 >                if (nlist > neighborListSize) then
191 >                   call expandList(listerror)
192 >                   if (listerror /= 0) then
193 >                      FFerror = -1
194 >                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
195 >                      return
196 >                   end if
197 >                endif
198 >                
199 >                list(nlist) = j
200 >                
201 >                
202 >                if (rijsq <  rcutsq) then
203 >                   call do_pair(Atype_i,Atype_j,i,j,r,rxij,ryij,rzij)
204 >                endif
205 >             endif
206 >          enddo inner
207 >       enddo
208  
209 <     ! save current configuration, contruct neighbor list,
210 <     ! and calculate forces
211 <     call save_neighborList(q)
439 <    
440 <     neighborListSize = getNeighborListSize()
441 <     nlist = 0
442 <    
443 <    
209 >       point(nrow + 1) = nlist + 1
210 >      
211 >    else !! (update)
212  
213 <     do i = 1, nrow
214 <        point(i) = nlist + 1
215 < #ifdef IS_MPI
216 <        Atype_i => identPtrListRow(i)%this
217 <        tag_i = tagRow(i)
218 < #else
451 <        Atype_i   => identPtrList(i)%this
452 <        j_start = i + 1
453 < #endif
213 >       ! use the list to find the neighbors
214 >       do i = 1, nrow
215 >          JBEG = POINT(i)
216 >          JEND = POINT(i+1) - 1
217 >          ! check thiat molecule i has neighbors
218 >          if (jbeg .le. jend) then
219  
220 <        inner: do j = j_start, ncol
221 < ! Assign identity pointers and tags
222 < #ifdef IS_MPI
223 <           Atype_j => identPtrListColumn(j)%this
224 <          
225 <           call get_interatomic_vector(i,j,qRow(:,i),qCol(:,j),&
226 <                rxij,ryij,rzij,rijsq,r)
227 < !! For mpi, use newtons 3rd law when building neigbor list
228 < !! Also check to see the particle i != j.
229 <           if (mpi_cycle_jLoop(i,j)) cycle inner:
220 >             Atype_i => identPtrListRow(i)%this
221 >             do jnab = jbeg, jend
222 >                j = list(jnab)
223 >                Atype_j = identPtrListColumn(j)%this
224 >                call get_interatomic_vector(i,j,qRow(:,i),qCol(:,j),&
225 >                     rxij,ryij,rzij,rijsq,r)
226 >                
227 >                call do_pair(Atype_i,Atype_j,i,j,r,rxij,ryij,rzij)
228 >             enddo
229 >          endif
230 >       enddo
231 >    endif
232  
233 < #else          
234 <           Atype_j   => identPtrList(j)%this
235 <           call get_interatomic_vector(i,j,q(:,i),q(:,j),&
469 <                rxij,ryij,rzij,rijsq,r)
233 > #else
234 >    
235 >    if (update_nlist) then
236        
237 < #endif          
237 >       ! save current configuration, contruct neighbor list,
238 >       ! and calculate forces
239 >       call save_neighborList(q)
240 >      
241 >       neighborListSize = getNeighborListSize()
242 >       nlist = 0
243 >      
244 >    
245 >       do i = 1, natoms-1
246 >          point(i) = nlist + 1
247 >          Atype_i   => identPtrList(i)%this
248 >
249 >          inner: do j = i+1, natoms
250 >             Atype_j   => identPtrList(j)%this
251 >             call get_interatomic_vector(i,j,q(:,i),q(:,j),&
252 >                  rxij,ryij,rzij,rijsq,r)
253            
254 <           if (rijsq <  rlistsq) then
254 >             if (rijsq <  rlistsq) then
255  
256 <              nlist = nlist + 1
256 >                nlist = nlist + 1
257 >                
258 >                if (nlist > neighborListSize) then
259 >                   call expandList(listerror)
260 >                   if (listerror /= 0) then
261 >                      FFerror = -1
262 >                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
263 >                      return
264 >                   end if
265 >                endif
266 >                
267 >                list(nlist) = j
268  
477 –              if (nlist > neighborListSize) then
478 –                 call expandList(listerror)
479 –                 if (listerror /= 0) then
480 –                    FFerror = -1
481 –                    write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
482 –                    return
483 –                 end if
484 –              endif
485 –
486 –              list(nlist) = j
487 –
269      
270 <              if (rijsq <  rcutsq) then
271 <                 call do_pair(Atype_i,Atype_j,i,j,r,rxij,ryij,rzij)
272 <              endif
270 >                if (rijsq <  rcutsq) then
271 >                   call do_pair(Atype_i,Atype_j,i,j,r,rxij,ryij,rzij)
272 >                endif
273 >             endif
274            enddo inner
275 <     enddo
275 >       enddo
276 >      
277 >       point(natoms) = nlist + 1
278 >      
279 >    else !! (update)
280 >      
281 >       ! use the list to find the neighbors
282 >       do i = 1, nrow
283 >          JBEG = POINT(i)
284 >          JEND = POINT(i+1) - 1
285 >          ! check thiat molecule i has neighbors
286 >          if (jbeg .le. jend) then
287 >            
288 >             Atype_i => identPtrList(i)%this
289 >             do jnab = jbeg, jend
290 >                j = list(jnab)
291 >                Atype_j = identPtrList(j)%this
292 >                call get_interatomic_vector(i,j,q(:,i),q(:,j),&
293 >                     rxij,ryij,rzij,rijsq,r)
294 >                call do_pair(Atype_i,Atype_j,i,j,r,rxij,ryij,rzij)
295 >             enddo
296 >          endif
297 >       enddo
298 >    endif
299  
300 + #endif
301 +
302 +
303   #ifdef IS_MPI
304 <     point(nrow + 1) = nlist + 1
305 < #else
306 <     point(natoms) = nlist + 1
304 >    !! distribute all reaction field stuff (or anything for post-pair):
305 >    call scatter(rflRow,rflTemp1,plan_row3d)
306 >    call scatter(rflCol,rflTemp2,plan_col3d)
307 >    do i = 1,nlocal
308 >       rflTemp(1:3,i) = rflTemp1(1:3,i) + rflTemp2(1:3,i)
309 >    end do
310   #endif
311  
312 <  else !! (update)
312 > ! This is the post-pair loop:
313 > #ifdef IS_MPI
314  
315 <     ! use the list to find the neighbors
316 <     do i = 1, nrow
317 <        JBEG = POINT(i)
318 <        JEND = POINT(i+1) - 1
319 <        ! check thiat molecule i has neighbors
320 <        if (jbeg .le. jend) then
315 >    if (system_has_postpair_atoms) then
316 >       do i = 1, nlocal
317 >          Atype_i => identPtrListRow(i)%this
318 >          call do_postpair(i, Atype_i)
319 >       enddo
320 >    endif
321  
510 – #ifdef IS_MPI
511 –           ljAtype_i => identPtrListRow(i)%this
322   #else
323 <           ljAtype_i => identPtrList(i)%this
323 >
324 >    if (system_has_postpair_atoms) then
325 >       do i = 1, natoms
326 >          Atype_i => identPtr(i)%this
327 >          call do_postpair(i, Atype_i)
328 >       enddo
329 >    endif
330 >
331   #endif
515 –           do jnab = jbeg, jend
516 –              j = list(jnab)
517 – #ifdef IS_MPI
518 –              ljAtype_j = identPtrListColumn(j)%this
519 –              call get_interatomic_vector(i,j,qRow(:,i),qCol(:,j),&
520 –                   rxij,ryij,rzij,rijsq,r)
521 –              
522 – #else
523 –              ljAtype_j = identPtrList(j)%this
524 –              call get_interatomic_vector(i,j,q(:,i),q(:,j),&
525 –                   rxij,ryij,rzij,rijsq,r)
526 – #endif
527 –              call do_pair(i,j,r,rxij,ryij,rzij)
528 –          enddo
529 –       endif
530 –    enddo
531 – endif
532 –
332  
333  
334 +
335 +
336   #ifdef IS_MPI
337      !!distribute forces
338  
339 <    call scatter(fRow,f,plan_row3d)
340 <    call scatter(fCol,fTemp,plan_col3d)
339 >    call scatter(fRow,fTemp1,plan_row3d)
340 >    call scatter(fCol,fTemp2,plan_col3d)
341  
342 +
343      do i = 1,nlocal
344 <       f(1:3,i) = f(1:3,i) + fTemp(1:3,i)
344 >       fTemp(1:3,i) = fTemp1(1:3,i) + fTemp2(1:3,i)
345      end do
346  
347      if (do_torque) then
348 <       call scatter(tRow,t,plan_row3d)
349 <       call scatter(tCol,tTemp,plan_col3d)
348 >       call scatter(tRow,tTemp1,plan_row3d)
349 >       call scatter(tCol,tTemp2,plan_col3d)
350      
351         do i = 1,nlocal
352 <          t(1:3,i) = t(1:3,i) + tTemp(1:3,i)
352 >          tTemp(1:3,i) = tTemp1(1:3,i) + tTemp2(1:3,i)
353         end do
354      endif
355  
# Line 571 | Line 373 | contains
373      endif
374   #else
375   ! Copy local array into return array for c
376 <    f = fTemp
377 <    t = tTemp
376 >    f = f+fTemp
377 >    t = t+tTemp
378   #endif
379  
380      potE = pe
# Line 654 | Line 456 | contains
456      real( kind = dp ) :: fx = 0.0_dp
457      real( kind = dp ) :: fy = 0.0_dp
458      real( kind = dp ) :: fz = 0.0_dp  
459 <
459 >  
460      real( kind = dp ) ::  drdx = 0.0_dp
461      real( kind = dp ) ::  drdy = 0.0_dp
462      real( kind = dp ) ::  drdz = 0.0_dp
463      
464 +
465 +    if (Atype_i%is_LJ .and. Atype_j%is_LJ) then
466 +       call getLJForce(r,pot,dudr,ljAtype_i,ljAtype_j,fx,fy,fz)
467 +    endif
468 +
469 +    if (Atype_i%is_dp .and. Atype_j%is_dp) then
470 +
471 + #ifdef IS_MPI
472 +       call dipole_dipole(i, j, atype_i, atype_j, rx_ij, ry_ij, rz_ij, r_ij, &
473 +            ulRow(:,i), ulCol(:,j), rt, rrf, pot)
474 + #else
475 +       call dipole_dipole(i, j, atype_i, atype_j, rx_ij, ry_ij, rz_ij, r_ij, &
476 +            ul(:,i), ul(:,j), rt, rrf, pot)
477 + #endif
478  
479 +       if (do_reaction_field) then
480 + #ifdef IS_MPI
481 +          call accumulate_rf(i, j, r_ij, rflRow(:,i), rflCol(:j), &
482 +               ulRow(:i), ulCol(:,j), rt, rrf)
483 + #else
484 +          call accumulate_rf(i, j, r_ij, rfl(:,i), rfl(:j), &
485 +               ul(:,i), ul(:,j), rt, rrf)
486 + #endif
487 +       endif
488  
489  
490 +    endif
491  
492 +    if (Atype_i%is_sticky .and. Atype_j%is_sticky) then
493 +       call getstickyforce(r,pot,dudr,ljAtype_i,ljAtype_j)
494 +    endif
495  
667 –    call getLJForce(r,pot,dudr,ljAtype_i,ljAtype_j)
496        
497   #ifdef IS_MPI
498                  eRow(i) = eRow(i) + pot*0.5

Diff Legend

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