ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/do_Forces.F90
(Generate patch)

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents):
Revision 631 by chuckv, Thu Jul 17 19:25:51 2003 UTC vs.
Revision 650 by chuckv, Thu Jul 24 19:57:35 2003 UTC

# Line 4 | Line 4
4  
5   !! @author Charles F. Vardeman II
6   !! @author Matthew Meineke
7 < !! @version $Id: do_Forces.F90,v 1.22 2003-07-17 19:25:51 chuckv Exp $, $Date: 2003-07-17 19:25:51 $, $Name: not supported by cvs2svn $, $Revision: 1.22 $
7 > !! @version $Id: do_Forces.F90,v 1.25 2003-07-24 19:57:35 chuckv Exp $, $Date: 2003-07-24 19:57:35 $, $Name: not supported by cvs2svn $, $Revision: 1.25 $
8  
9   module do_Forces
10    use force_globals
# Line 18 | Line 18 | module do_Forces
18    use reaction_field
19    use gb_pair
20    use vector_class
21 +  use eam
22   #ifdef IS_MPI
23    use mpiSimulation
24   #endif
# Line 256 | Line 257 | contains
257         !! do_prepair_loop_if_needed
258         !! if_mpi_scatter_stuff_from_prepair
259         !! if_mpi_gather_stuff_from_prepair_to_main_loop
260 +
261 + !--------------------PREFORCE LOOP----------->>>>>>>>>>>>>>>>>>>>>>>>>>>
262 + #ifdef IS_MPI
263 +    
264 +    if (update_nlist) then
265 +      
266 +       !! save current configuration, construct neighbor list,
267 +       !! and calculate forces
268 +       call saveNeighborList(nlocal, q)
269 +      
270 +       neighborListSize = size(list)
271 +       nlist = 0      
272 +      
273 +       do i = 1, nrow
274 +          point(i) = nlist + 1
275 +          
276 +          prepair_inner: do j = 1, ncol
277 +            
278 +             if (skipThisPair(i,j)) cycle prepair_inner
279 +            
280 +             call get_interatomic_vector(q_Row(:,i), q_Col(:,j), d, rijsq)
281 +            
282 +             if (rijsq < rlistsq) then            
283 +                
284 +                nlist = nlist + 1
285 +                
286 +                if (nlist > neighborListSize) then
287 +                   call expandNeighborList(nlocal, listerror)
288 +                   if (listerror /= 0) then
289 +                      error = -1
290 +                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
291 +                      return
292 +                   end if
293 +                   neighborListSize = size(list)
294 +                endif
295 +                
296 +                list(nlist) = j
297 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot_local)                      
298 +             endif
299 +          enddo prepair_inner
300 +       enddo
301 +
302 +       point(nrow + 1) = nlist + 1
303 +      
304 +    else  !! (of update_check)
305 +
306 +       ! use the list to find the neighbors
307 +       do i = 1, nrow
308 +          JBEG = POINT(i)
309 +          JEND = POINT(i+1) - 1
310 +          ! check thiat molecule i has neighbors
311 +          if (jbeg .le. jend) then
312 +            
313 +             do jnab = jbeg, jend
314 +                j = list(jnab)
315 +
316 +                call get_interatomic_vector(q_Row(:,i), q_Col(:,j), d, rijsq)
317 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
318 +                     u_l, A, f, t, pot_local)
319 +
320 +             enddo
321 +          endif
322 +       enddo
323 +    endif
324 +    
325 + #else
326 +    
327 +    if (update_nlist) then
328 +      
329 +       ! save current configuration, contruct neighbor list,
330 +       ! and calculate forces
331 +       call saveNeighborList(natoms, q)
332 +      
333 +       neighborListSize = size(list)
334 +  
335 +       nlist = 0
336 +      
337 +       do i = 1, natoms-1
338 +          point(i) = nlist + 1
339 +          
340 +          prepair_inner: do j = i+1, natoms
341 +            
342 +             if (skipThisPair(i,j))  cycle prepair_inner
343 +                          
344 +             call get_interatomic_vector(q(:,i), q(:,j), d, rijsq)
345 +          
346 +
347 +             if (rijsq < rlistsq) then
348 +                
349 +                nlist = nlist + 1
350 +              
351 +                if (nlist > neighborListSize) then
352 +                   call expandNeighborList(natoms, listerror)
353 +                   if (listerror /= 0) then
354 +                      error = -1
355 +                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
356 +                      return
357 +                   end if
358 +                   neighborListSize = size(list)
359 +                endif
360 +                
361 +                list(nlist) = j
362 +                
363 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
364 +                        u_l, A, f, t, pot)
365 +                
366 +             endif
367 +          enddo prepair_inner
368 +       enddo
369 +      
370 +       point(natoms) = nlist + 1
371 +      
372 +    else !! (update)
373 +      
374 +       ! use the list to find the neighbors
375 +       do i = 1, natoms-1
376 +          JBEG = POINT(i)
377 +          JEND = POINT(i+1) - 1
378 +          ! check thiat molecule i has neighbors
379 +          if (jbeg .le. jend) then
380 +            
381 +             do jnab = jbeg, jend
382 +                j = list(jnab)
383 +
384 +                call get_interatomic_vector(q(:,i), q(:,j), d, rijsq)
385 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
386 +                     u_l, A, f, t, pot)
387 +
388 +             enddo
389 +          endif
390 +       enddo
391 +    endif    
392 + #endif
393 +    !! Do rest of preforce calculations
394 +   call do_preforce(nlocal,pot)
395      else
396 <       !! See if we need to update neighbor lists
396 >       !! See if we need to update neighbor lists for non pre-pair
397         call checkNeighborList(nlocal, q, listSkin, update_nlist)  
398      endif
399 <    
399 >
400 >
401 >
402 >
403 >
404 > !---------------------------------MAIN Pair LOOP->>>>>>>>>>>>>>>>>>>>>>>>>>>>
405 >
406 >
407 >
408 >
409 >  
410   #ifdef IS_MPI
411      
412      if (update_nlist) then
# Line 531 | Line 677 | contains
677      logical :: is_LJ_i, is_LJ_j
678      logical :: is_DP_i, is_DP_j
679      logical :: is_GB_i, is_GB_j
680 +    logical :: is_EAM_i,is_EAM_j
681      logical :: is_Sticky_i, is_Sticky_j
682      integer :: me_i, me_j
683  
# Line 599 | Line 746 | contains
746      endif
747      
748  
749 +  
750 +   if (FF_uses_EAM .and. SimUsesEAM()) then
751 +      call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
752 +      call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
753 +      
754 +      if ( is_EAM_i .and. is_EAM_j ) &
755 +           call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
756 +   endif
757 +
758  
759 +
760 +
761    end subroutine do_pair
762  
763  
764  
765 <  subroutine do_preforce(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
765 >  subroutine do_prepair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
766     real( kind = dp ) :: pot
767     real( kind = dp ), dimension(3,getNlocal()) :: u_l
768     real (kind=dp), dimension(9,getNlocal()) :: A
# Line 643 | Line 801 | contains
801        call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
802        
803        if ( is_EAM_i .and. is_EAM_j ) &
804 <           call calc_EAM_prepair(i, j, d, r, rijsq )
804 >           call calc_EAM_prepair_rho(i, j, d, r, rijsq )
805 >   endif
806 >  end subroutine do_prepair
807 >
808 >
809 >
810 >
811 >  subroutine do_preforce(nlocal,pot)
812 >    integer :: nlocal
813 >    real( kind = dp ) :: pot
814 >
815 >   if (FF_uses_EAM .and. SimUsesEAM()) then
816 >      call calc_EAM_preforce_Frho(nlocal,pot)
817     endif
818 +
819 +
820    end subroutine do_preforce
821    
822    

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines