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 634 by mmeineke, Thu Jul 17 19:38:23 2003 UTC vs.
Revision 669 by chuckv, Thu Aug 7 00:47:33 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.23 2003-07-17 19:38:23 mmeineke Exp $, $Date: 2003-07-17 19:38:23 $, $Name: not supported by cvs2svn $, $Revision: 1.23 $
7 > !! @version $Id: do_Forces.F90,v 1.27 2003-08-07 00:47:33 chuckv Exp $, $Date: 2003-08-07 00:47:33 $, $Name: not supported by cvs2svn $, $Revision: 1.27 $
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 +  use status
23   #ifdef IS_MPI
24    use mpiSimulation
25   #endif
# Line 136 | Line 138 | contains
138  
139      if (FF_uses_sticky) then
140         call check_sticky_FF(my_status)
141 +       if (my_status /= 0) then
142 +          thisStat = -1
143 +          return
144 +       end if
145 +    endif
146 +
147 +
148 +    if (FF_uses_EAM) then
149 +       call init_EAM_FF(my_status)
150         if (my_status /= 0) then
151            thisStat = -1
152            return
153         end if
154      endif
155 +
156 +
157      
158      if (FF_uses_GB) then
159         call check_gb_pair_FF(my_status)
# Line 161 | Line 174 | contains
174            return
175         endif
176      endif
177 +    
178  
179      havePolicies = .true.
180      if( haveRlist ) do_forces_initialized = .true.
181 <    
181 >
182    end subroutine init_FF
183    
184  
# Line 221 | Line 235 | contains
235      nlocal = getNlocal()
236      natoms = nlocal
237   #endif
238 <  
238 >
239      call check_initialization(localError)
240      if ( localError .ne. 0 ) then
241 +       call handleError("do_force_loop","Not Initialized")
242         error = -1
243         return
244      end if
# Line 232 | Line 247 | contains
247      do_pot = do_pot_c
248      do_stress = do_stress_c
249  
250 +
251      ! Gather all information needed by all force loops:
252      
253   #ifdef IS_MPI    
# Line 248 | Line 264 | contains
264      endif
265      
266   #endif
267 <    
267 >  
268      if (FF_RequiresPrepairCalc() .and. SimRequiresPrepairCalc()) then
269         !! See if we need to update neighbor lists
270         call checkNeighborList(nlocal, q, listSkin, update_nlist)  
# Line 256 | Line 272 | contains
272         !! do_prepair_loop_if_needed
273         !! if_mpi_scatter_stuff_from_prepair
274         !! if_mpi_gather_stuff_from_prepair_to_main_loop
275 +
276 + !--------------------PREFORCE LOOP----------->>>>>>>>>>>>>>>>>>>>>>>>>>>
277 + #ifdef IS_MPI
278 +    
279 +    if (update_nlist) then
280 +      
281 +       !! save current configuration, construct neighbor list,
282 +       !! and calculate forces
283 +       call saveNeighborList(nlocal, q)
284 +      
285 +       neighborListSize = size(list)
286 +       nlist = 0      
287 +      
288 +       do i = 1, nrow
289 +          point(i) = nlist + 1
290 +          
291 +          prepair_inner: do j = 1, ncol
292 +            
293 +             if (skipThisPair(i,j)) cycle prepair_inner
294 +            
295 +             call get_interatomic_vector(q_Row(:,i), q_Col(:,j), d, rijsq)
296 +            
297 +             if (rijsq < rlistsq) then            
298 +                
299 +                nlist = nlist + 1
300 +                
301 +                if (nlist > neighborListSize) then
302 +                   call expandNeighborList(nlocal, listerror)
303 +                   if (listerror /= 0) then
304 +                      error = -1
305 +                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
306 +                      return
307 +                   end if
308 +                   neighborListSize = size(list)
309 +                endif
310 +                
311 +                list(nlist) = j
312 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot_local)                      
313 +             endif
314 +          enddo prepair_inner
315 +       enddo
316 +
317 +       point(nrow + 1) = nlist + 1
318 +      
319 +    else  !! (of update_check)
320 +
321 +       ! use the list to find the neighbors
322 +       do i = 1, nrow
323 +          JBEG = POINT(i)
324 +          JEND = POINT(i+1) - 1
325 +          ! check thiat molecule i has neighbors
326 +          if (jbeg .le. jend) then
327 +            
328 +             do jnab = jbeg, jend
329 +                j = list(jnab)
330 +
331 +                call get_interatomic_vector(q_Row(:,i), q_Col(:,j), d, rijsq)
332 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
333 +                     u_l, A, f, t, pot_local)
334 +
335 +             enddo
336 +          endif
337 +       enddo
338 +    endif
339 +    
340 + #else
341 +    
342 +    if (update_nlist) then
343 +      
344 +       ! save current configuration, contruct neighbor list,
345 +       ! and calculate forces
346 +       call saveNeighborList(natoms, q)
347 +      
348 +       neighborListSize = size(list)
349 +  
350 +       nlist = 0
351 +      
352 +       do i = 1, natoms-1
353 +          point(i) = nlist + 1
354 +          
355 +          prepair_inner: do j = i+1, natoms
356 +            
357 +             if (skipThisPair(i,j))  cycle prepair_inner
358 +                          
359 +             call get_interatomic_vector(q(:,i), q(:,j), d, rijsq)
360 +          
361 +
362 +             if (rijsq < rlistsq) then
363 +                
364 +                nlist = nlist + 1
365 +              
366 +                if (nlist > neighborListSize) then
367 +                   call expandNeighborList(natoms, listerror)
368 +                   if (listerror /= 0) then
369 +                      error = -1
370 +                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
371 +                      return
372 +                   end if
373 +                   neighborListSize = size(list)
374 +                endif
375 +                
376 +                list(nlist) = j
377 +                
378 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
379 +                        u_l, A, f, t, pot)
380 +                
381 +             endif
382 +          enddo prepair_inner
383 +       enddo
384 +      
385 +       point(natoms) = nlist + 1
386 +      
387 +    else !! (update)
388 +      
389 +       ! use the list to find the neighbors
390 +       do i = 1, natoms-1
391 +          JBEG = POINT(i)
392 +          JEND = POINT(i+1) - 1
393 +          ! check thiat molecule i has neighbors
394 +          if (jbeg .le. jend) then
395 +            
396 +             do jnab = jbeg, jend
397 +                j = list(jnab)
398 +
399 +                call get_interatomic_vector(q(:,i), q(:,j), d, rijsq)
400 +                call do_prepair(i, j, rijsq, d, do_pot, do_stress, &
401 +                     u_l, A, f, t, pot)
402 +
403 +             enddo
404 +          endif
405 +       enddo
406 +    endif    
407 + #endif
408 +    !! Do rest of preforce calculations
409 +   call do_preforce(nlocal,pot)
410      else
411 <       !! See if we need to update neighbor lists
411 >       !! See if we need to update neighbor lists for non pre-pair
412         call checkNeighborList(nlocal, q, listSkin, update_nlist)  
413      endif
414 <    
414 >
415 >
416 >
417 >
418 >
419 > !---------------------------------MAIN Pair LOOP->>>>>>>>>>>>>>>>>>>>>>>>>>>>
420 >
421 >
422 >
423 >
424 >  
425   #ifdef IS_MPI
426      
427      if (update_nlist) then
# Line 531 | Line 692 | contains
692      logical :: is_LJ_i, is_LJ_j
693      logical :: is_DP_i, is_DP_j
694      logical :: is_GB_i, is_GB_j
695 +    logical :: is_EAM_i,is_EAM_j
696      logical :: is_Sticky_i, is_Sticky_j
697      integer :: me_i, me_j
698  
# Line 599 | Line 761 | contains
761      endif
762      
763  
764 +  
765 +   if (FF_uses_EAM .and. SimUsesEAM()) then
766 +      call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
767 +      call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
768 +      
769 +      if ( is_EAM_i .and. is_EAM_j ) &
770 +           call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
771 +   endif
772 +
773  
774 +
775 +
776    end subroutine do_pair
777  
778  
779  
780 <  subroutine do_preforce(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
780 >  subroutine do_prepair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
781     real( kind = dp ) :: pot
782     real( kind = dp ), dimension(3,getNlocal()) :: u_l
783     real (kind=dp), dimension(9,getNlocal()) :: A
# Line 623 | Line 796 | contains
796    
797     r = sqrt(rijsq)
798    
799 +
800   #ifdef IS_MPI
801     if (tagRow(i) .eq. tagColumn(j)) then
802        write(0,*) 'do_pair is doing', i , j, tagRow(i), tagColumn(j)
# Line 642 | Line 816 | contains
816        call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
817        call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
818        
819 < !!$      if ( is_EAM_i .and. is_EAM_j ) &
820 < !!$           call calc_EAM_prepair(i, j, d, r, rijsq )
819 >      if ( is_EAM_i .and. is_EAM_j ) &
820 >           call calc_EAM_prepair_rho(i, j, d, r, rijsq )
821     endif
822 +  end subroutine do_prepair
823 +
824 +
825 +
826 +
827 +  subroutine do_preforce(nlocal,pot)
828 +    integer :: nlocal
829 +    real( kind = dp ) :: pot
830 +
831 +    if (FF_uses_EAM .and. SimUsesEAM()) then
832 +       call calc_EAM_preforce_Frho(nlocal,pot)
833 +    endif
834 +
835 +
836    end subroutine do_preforce
837    
838    
# Line 705 | Line 893 | contains
893      error = 0
894      ! Make sure we are properly initialized.
895      if (.not. do_forces_initialized) then
896 +       write(*,*) "Forces not initialized"
897         error = -1
898         return
899      endif

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines