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 1192 by gezelter, Mon May 24 21:03:30 2004 UTC vs.
Revision 1197 by gezelter, Wed May 26 16:41:23 2004 UTC

# Line 4 | Line 4
4  
5   !! @author Charles F. Vardeman II
6   !! @author Matthew Meineke
7 < !! @version $Id: do_Forces.F90,v 1.61 2004-05-24 21:03:25 gezelter Exp $, $Date: 2004-05-24 21:03:25 $, $Name: not supported by cvs2svn $, $Revision: 1.61 $
7 > !! @version $Id: do_Forces.F90,v 1.62 2004-05-26 16:41:23 gezelter Exp $, $Date: 2004-05-26 16:41:23 $, $Name: not supported by cvs2svn $, $Revision: 1.62 $
8  
9   module do_Forces
10    use force_globals
# Line 32 | Line 32 | module do_Forces
32   #define __FORTRAN90
33   #include "fForceField.h"
34   #include "fSwitchingFunction.h"
35 +
36 +  INTEGER, PARAMETER:: PREPAIR_LOOP = 1
37 +  INTEGER, PARAMETER:: PAIR_LOOP    = 2
38  
39    logical, save :: haveRlist = .false.
40    logical, save :: haveNeighborList = .false.
# Line 399 | Line 402 | contains
402   #endif
403      integer :: natoms    
404      logical :: update_nlist  
405 <    integer :: i, j, jbeg, jend, jnab
406 <    integer :: istart, iend, jstart
405 >    integer :: i, j, jstart, jend, jnab
406 >    integer :: istart, iend
407      integer :: ia, jb, atom1, atom2
408      integer :: nlist
409      real( kind = DP ) :: ratmsq, rgrpsq, rgrp, vpair, vij
# Line 413 | Line 416 | contains
416      integer :: listerror, error
417      integer :: localError
418      integer :: propPack_i, propPack_j
419 +    integer :: loopStart, loopEnd, loop
420  
421      real(kind=dp) :: listSkin = 1.0  
422      
# Line 465 | Line 469 | contains
469      nloops = nloops + 1
470   #endif
471      
472 +    loopEnd = PAIR_LOOP
473      if (FF_RequiresPrepairCalc() .and. SIM_requires_prepair_calc) then
474 <       !! See if we need to update neighbor lists
474 >       loopStart = PREPAIR_LOOP
475 >    else
476 >       loopStart = PAIR_LOOP
477 >    endif
478  
479 <       call checkNeighborList(nGroup, q_group, listSkin, update_nlist)
479 >    do loop = loopStart, loopEnd
480  
481 <       !! if_mpi_gather_stuff_for_prepair
482 <       !! do_prepair_loop_if_needed
483 <       !! if_mpi_scatter_stuff_from_prepair
484 <       !! if_mpi_gather_stuff_from_prepair_to_main_loop
481 >       ! See if we need to update neighbor lists
482 >       ! (but only on the first time through):
483 >       if (loop .eq. loopStart) then
484 >          call checkNeighborList(nGroup, q_group, listSkin, update_nlist)
485 >       endif
486        
478       !--------------------PREFORCE LOOP----------->>>>>>>>>>>>>>>>>>>>>>>>>>>
479
487         if (update_nlist) then
488 <          
489 <          !! save current configuration, construct neighbor list,
490 <          !! and calculate forces
491 <          
492 <          call saveNeighborList(nGroup, q_group)
493 <          
494 <          neighborListSize = size(list)
488 <          nlist = 0      
489 <          
490 <          istart = 1
488 >          !! save current configuration and construct neighbor list
489 >          call saveNeighborList(nGroup, q_group)          
490 >          neighborListSize = size(list)
491 >          nlist = 0
492 >       endif
493 >      
494 >       istart = 1
495   #ifdef IS_MPI
496 <          iend = nrow_group
496 >       iend = nrow_group
497   #else
498 <          iend = nGroup - 1
498 >       iend = nGroup - 1
499   #endif
500 <          do i = istart, iend
501 <            
502 <             point(i) = nlist + 1
503 <            
504 <             n_in_i = groupStart(i+1) - groupStart(i)
505 <            
500 >       outer: do i = istart, iend
501 >          
502 >          if (update_nlist) point(i) = nlist + 1
503 >          
504 >          n_in_i = groupStart(i+1) - groupStart(i)
505 >          
506 >          if (update_nlist) then
507   #ifdef IS_MPI
508               jstart = 1
509               jend = ncol_group
# Line 506 | Line 511 | contains
511               jstart = i+1
512               jend = nGroup
513   #endif
514 <             do j = jstart, jend
515 <                
516 < #ifdef IS_MPI
517 <                call get_interatomic_vector(q_group_Row(:,i), &
518 <                     q_group_Col(:,j), d_grp, rgrpsq)
514 >          else            
515 >             jstart = point(i)
516 >             jend = point(i+1) - 1
517 >             ! make sure group i has neighbors
518 >             if (jstart .gt. jend) cycle outer
519 >          endif
520 >          
521 >          do jnab = jstart, jend
522 >             if (update_nlist) then
523 >                j = jnab
524 >             else
525 >                j = list(jnab)
526 >             endif
527 > #ifdef IS_MPI
528 >             call get_interatomic_vector(q_group_Row(:,i), &
529 >                  q_group_Col(:,j), d_grp, rgrpsq)
530   #else
531 <                call get_interatomic_vector(q_group(:,i), &
532 <                     q_group(:,j), d_grp, rgrpsq)
533 < #endif            
534 <                if (rgrpsq < rlistsq) then
531 >             call get_interatomic_vector(q_group(:,i), &
532 >                  q_group(:,j), d_grp, rgrpsq)
533 > #endif
534 >             if (rgrpsq < rlistsq) then
535 >                if (update_nlist) then
536                     nlist = nlist + 1
537                    
538                     if (nlist > neighborListSize) then
# Line 528 | Line 545 | contains
545                        neighborListSize = size(list)
546                     endif
547                    
548 <                   list(nlist) = j                  
532 <                  
533 <                   call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
534 <                        in_switching_region)
535 <                  
536 <                   n_in_j = groupStart(j+1) - groupStart(j)
537 <                  
538 <                   do ia = groupStart(i), groupStart(i+1)-1
539 <                      atom1 = groupList(ia)
540 <                      
541 <                      prepair_inner1: do jb = groupStart(j), groupStart(j+1)-1
542 <                         atom2 = groupList(jb)                    
543 <                        
544 <                         if (skipThisPair(atom1, atom2)) cycle prepair_inner1
545 <                        
546 <                         if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
547 <                            d_atm(1:3) = d_grp(1:3)
548 <                            ratmsq = rgrpsq
549 <                         else
550 < #ifdef IS_MPI
551 <                            call get_interatomic_vector(q_Row(:,atom1), &
552 <                                 q_Col(:,atom2), d_atm, ratmsq)
553 < #else
554 <                            call get_interatomic_vector(q(:,atom1), &
555 <                                 q(:,atom2), d_atm, ratmsq)
556 < #endif
557 <                         endif
558 < #ifdef IS_MPI                      
559 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
560 <                              rgrpsq, d_grp, do_pot, do_stress, &
561 <                              u_l, A, f, t, pot_local)
562 < #else
563 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
564 <                              rgrpsq, d_grp, do_pot, do_stress, &
565 <                              u_l, A, f, t, pot)
566 < #endif                                              
567 <                      enddo prepair_inner1
568 <                   enddo
569 <                  
570 <                end if
571 <             enddo
572 <          enddo
573 <          
574 < #ifdef IS_MPI
575 <          point(nrow_group + 1) = nlist + 1
576 < #else
577 <          point(nGroup) = nlist + 1
578 < #endif
579 <          
580 <       else  !! (of update_check)
581 <          
582 <          ! use the list to find the neighbors
583 <          
584 <          istart = 1
585 < #ifdef IS_MPI
586 <          iend = nrow_group
587 < #else
588 <          iend = nGroup - 1
589 < #endif
590 <          
591 <          do i = istart, iend
592 <            
593 <             n_in_i = groupStart(i+1) - groupStart(i)
594 <            
595 <             JBEG = POINT(i)
596 <             JEND = POINT(i+1) - 1
597 <             ! check that group i has neighbors
598 <             if (jbeg .le. jend) then
599 <                
600 <                do jnab = jbeg, jend
601 <                   j = list(jnab)
602 <                  
603 < #ifdef IS_MPI
604 <                   call get_interatomic_vector(q_group_Row(:,i), &
605 <                        q_group_Col(:,j), d_grp, rgrpsq)
606 < #else
607 <                   call get_interatomic_vector(q_group(:,i), &
608 <                        q_group(:,j), d_grp, rgrpsq)
609 < #endif                                  
610 <                   call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
611 <                        in_switching_region)
612 <                  
613 <                   n_in_j = groupStart(j+1) - groupStart(j)
614 <                  
615 <                   do ia = groupStart(i), groupStart(i+1)-1
616 <                      atom1 = groupList(ia)
617 <                      
618 <                      prepair_inner2: do jb = groupStart(j), groupStart(j+1)-1
619 <                        
620 <                         atom2 = groupList(jb)                        
621 <                        
622 <                         if (skipThisPair(atom1, atom2)) cycle prepair_inner2
623 <                        
624 <                         if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
625 <                            d_atm(1:3) = d_grp(1:3)
626 <                            ratmsq = rgrpsq
627 <                         else
628 < #ifdef IS_MPI
629 <                            call get_interatomic_vector(q_Row(:,atom1), &
630 <                                 q_Col(:,atom2), d_atm, ratmsq)
631 < #else
632 <                            call get_interatomic_vector(q(:,atom1), &
633 <                                 q(:,atom2), d_atm, ratmsq)
634 < #endif
635 <                         endif
636 <
637 < #ifdef IS_MPI                      
638 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
639 <                              rgrpsq, d_grp, do_pot, do_stress, &
640 <                              u_l, A, f, t, pot_local)
641 < #else
642 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
643 <                              rgrpsq, d_grp, do_pot, do_stress, &
644 <                              u_l, A, f, t, pot)
645 < #endif                                              
646 <                      enddo prepair_inner2
647 <                   enddo
648 <                enddo
649 <             endif
650 <          enddo
651 <       endif
652 <       !! Do rest of preforce calculations
653 <       !! do necessary preforce calculations  
654 <       call do_preforce(nlocal,pot)
655 <       ! we have already updated the neighbor list set it to false...
656 <       update_nlist = .false.
657 <    else
658 <       !! See if we need to update neighbor lists for non pre-pair
659 <       call checkNeighborList(nGroup, q_group, listSkin, update_nlist)  
660 <    end if
661 <    
662 <    !---------------------------------MAIN Pair LOOP->>>>>>>>>>>>>
663 <
664 <    if (update_nlist) then
665 <      
666 <       !! save current configuration, construct neighbor list,
667 <       !! and calculate forces
668 <      
669 <       call saveNeighborList(nGroup, q_group)
670 <      
671 <       neighborListSize = size(list)
672 <       nlist = 0      
673 <      
674 <       istart = 1
675 < #ifdef IS_MPI
676 <       iend = nrow_group
677 < #else
678 <       iend = nGroup - 1
679 < #endif
680 <       do i = istart, iend
681 <
682 <          point(i) = nlist + 1
683 <
684 <          n_in_i = groupStart(i+1) - groupStart(i)
685 <          
686 < #ifdef IS_MPI
687 <          jstart = 1
688 <          jend = ncol_group
689 < #else
690 <          jstart = i+1
691 <          jend = nGroup
692 < #endif
693 <          do j = jstart, jend
694 <            
695 < #ifdef IS_MPI
696 <             call get_interatomic_vector(q_group_Row(:,i), &
697 <                  q_group_Col(:,j), d_grp, rgrpsq)
698 < #else
699 <             call get_interatomic_vector(q_group(:,i), &
700 <                  q_group(:,j), d_grp, rgrpsq)
701 < #endif            
702 <             if (rgrpsq < rlistsq) then
703 <                nlist = nlist + 1
704 <                
705 <                if (nlist > neighborListSize) then
706 <                   call expandNeighborList(nGroup, listerror)
707 <                   if (listerror /= 0) then
708 <                      error = -1
709 <                      write(DEFAULT_ERROR,*) "ERROR: nlist > list size and max allocations exceeded."
710 <                      return
711 <                   end if
712 <                   neighborListSize = size(list)
548 >                   list(nlist) = j
549                  endif
550                  
551 <                list(nlist) = j
552 <                
553 <                vij = 0.0d0
554 <                fij(1:3) = 0.0d0
551 >                if (loop .eq. PAIR_LOOP) then
552 >                   vij = 0.0d0
553 >                   fij(1:3) = 0.0d0
554 >                endif
555                  
556                  call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
557                       in_switching_region)
558 <
558 >                
559                  n_in_j = groupStart(j+1) - groupStart(j)
560                  
561                  do ia = groupStart(i), groupStart(i+1)-1
562 +                  
563                     atom1 = groupList(ia)
564                    
565 <                   inner1: do jb = groupStart(j), groupStart(j+1)-1
729 <                      atom2 = groupList(jb)                    
730 <
731 <                      if (skipThisPair(atom1, atom2)) cycle inner1
565 >                   inner: do jb = groupStart(j), groupStart(j+1)-1
566                        
567 +                      atom2 = groupList(jb)
568 +                      
569 +                      if (skipThisPair(atom1, atom2)) cycle inner
570 +                      
571                        if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
572                           d_atm(1:3) = d_grp(1:3)
573                           ratmsq = rgrpsq
# Line 742 | Line 580 | contains
580                                q(:,atom2), d_atm, ratmsq)
581   #endif
582                        endif
583 +                      if (loop .eq. PREPAIR_LOOP) then
584   #ifdef IS_MPI                      
585 <                      call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
586 <                           do_pot, &
587 <                           u_l, A, f, t, pot_local, vpair, fpair)
585 >                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
586 >                              rgrpsq, d_grp, do_pot, do_stress, &
587 >                              u_l, A, f, t, pot_local)
588   #else
589 <                      call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
590 <                           do_pot,  &
591 <                           u_l, A, f, t, pot, vpair, fpair)
592 < #endif                      
593 <                      vij = vij + vpair
594 <                      fij(1:3) = fij(1:3) + fpair(1:3)
595 <
596 <                   enddo inner1
589 >                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
590 >                              rgrpsq, d_grp, do_pot, do_stress, &
591 >                              u_l, A, f, t, pot)
592 > #endif                                              
593 >                      else
594 > #ifdef IS_MPI                      
595 >                         call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
596 >                              do_pot, &
597 >                              u_l, A, f, t, pot_local, vpair, fpair)
598 > #else
599 >                         call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
600 >                              do_pot,  &
601 >                              u_l, A, f, t, pot, vpair, fpair)
602 > #endif
603 >                         vij = vij + vpair
604 >                         fij(1:3) = fij(1:3) + fpair(1:3)
605 >                      endif
606 >                   enddo inner
607                  enddo
608 <
609 <                if (in_switching_region) then
610 <                   swderiv = vij*dswdr/rgrp
611 <                   fij(1) = fij(1) + swderiv*d_grp(1)
612 <                   fij(2) = fij(2) + swderiv*d_grp(2)
613 <                   fij(3) = fij(3) + swderiv*d_grp(3)
614 <                  
615 <                   do ia=groupStart(i), groupStart(i+1)-1
616 <                      atom1=groupList(ia)
617 <                      mf = mfact(atom1)                  
608 >                
609 >                if (loop .eq. PAIR_LOOP) then
610 >                   if (in_switching_region) then
611 >                      swderiv = vij*dswdr/rgrp
612 >                      fij(1) = fij(1) + swderiv*d_grp(1)
613 >                      fij(2) = fij(2) + swderiv*d_grp(2)
614 >                      fij(3) = fij(3) + swderiv*d_grp(3)
615 >                      
616 >                      do ia=groupStart(i), groupStart(i+1)-1
617 >                         atom1=groupList(ia)
618 >                         mf = mfact(atom1)
619   #ifdef IS_MPI
620 <                      f_Row(1,atom1) = f_Row(1,atom1) + swderiv*d_grp(1)*mf
621 <                      f_Row(2,atom1) = f_Row(2,atom1) + swderiv*d_grp(2)*mf
622 <                      f_Row(3,atom1) = f_Row(3,atom1) + swderiv*d_grp(3)*mf
620 >                         f_Row(1,atom1) = f_Row(1,atom1) + swderiv*d_grp(1)*mf
621 >                         f_Row(2,atom1) = f_Row(2,atom1) + swderiv*d_grp(2)*mf
622 >                         f_Row(3,atom1) = f_Row(3,atom1) + swderiv*d_grp(3)*mf
623   #else
624 <                      f(1,atom1) = f(1,atom1) + swderiv*d_grp(1)*mf
625 <                      f(2,atom1) = f(2,atom1) + swderiv*d_grp(2)*mf
626 <                      f(3,atom1) = f(3,atom1) + swderiv*d_grp(3)*mf
627 < #endif                      
628 <                   enddo
629 <                  
630 <                   do jb=groupStart(j), groupStart(j+1)-1
631 <                      atom2=groupList(jb)
632 <                      mf = mfact(atom2)
624 >                         f(1,atom1) = f(1,atom1) + swderiv*d_grp(1)*mf
625 >                         f(2,atom1) = f(2,atom1) + swderiv*d_grp(2)*mf
626 >                         f(3,atom1) = f(3,atom1) + swderiv*d_grp(3)*mf
627 > #endif
628 >                      enddo
629 >                      
630 >                      do jb=groupStart(j), groupStart(j+1)-1
631 >                         atom2=groupList(jb)
632 >                         mf = mfact(atom2)
633   #ifdef IS_MPI
634 <                      f_Col(1,atom2) = f_Col(1,atom2) - swderiv*d_grp(1)*mf
635 <                      f_Col(2,atom2) = f_Col(2,atom2) - swderiv*d_grp(2)*mf
636 <                      f_Col(3,atom2) = f_Col(3,atom2) - swderiv*d_grp(3)*mf
634 >                         f_Col(1,atom2) = f_Col(1,atom2) - swderiv*d_grp(1)*mf
635 >                         f_Col(2,atom2) = f_Col(2,atom2) - swderiv*d_grp(2)*mf
636 >                         f_Col(3,atom2) = f_Col(3,atom2) - swderiv*d_grp(3)*mf
637   #else
638 <                      f(1,atom2) = f(1,atom2) - swderiv*d_grp(1)*mf
639 <                      f(2,atom2) = f(2,atom2) - swderiv*d_grp(2)*mf
640 <                      f(3,atom2) = f(3,atom2) - swderiv*d_grp(3)*mf
638 >                         f(1,atom2) = f(1,atom2) - swderiv*d_grp(1)*mf
639 >                         f(2,atom2) = f(2,atom2) - swderiv*d_grp(2)*mf
640 >                         f(3,atom2) = f(3,atom2) - swderiv*d_grp(3)*mf
641   #endif
642 <                   enddo                        
642 >                      enddo
643 >                   endif
644 >                  
645 >                   if (do_stress) call add_stress_tensor(d_grp, fij)
646                  endif
794
795                if (do_stress) call add_stress_tensor(d_grp, fij)
796
647               end if
648            enddo
649 <       enddo
650 <
649 >       enddo outer
650 >      
651 >       if (update_nlist) then
652   #ifdef IS_MPI
653 <       point(nrow_group + 1) = nlist + 1
653 >          point(nrow_group + 1) = nlist + 1
654   #else
655 <       point(nGroup) = nlist + 1
655 >          point(nGroup) = nlist + 1
656   #endif
657 <      
658 <    else  !! (of update_check)
659 <      
660 <       ! use the list to find the neighbors
661 <
662 <       istart = 1
663 < #ifdef IS_MPI
813 <       iend = nrow_group
814 < #else
815 <       iend = nGroup - 1
816 < #endif
817 <
818 <       do i = istart, iend
819 <
820 <          n_in_i = groupStart(i+1) - groupStart(i)
821 <          
822 <          JBEG = POINT(i)
823 <          JEND = POINT(i+1) - 1
824 <          ! check that group i has neighbors
825 <          if (jbeg .le. jend) then
657 >          if (loop .eq. PREPAIR_LOOP) then
658 >             ! we just did the neighbor list update on the first
659 >             ! pass, so we don't need to do it
660 >             ! again on the second pass
661 >             update_nlist = .false.                              
662 >          endif
663 >       endif
664              
665 <             do jnab = jbeg, jend
666 <                j = list(jnab)
667 <                              
668 < #ifdef IS_MPI
669 <                call get_interatomic_vector(q_group_Row(:,i), &
832 <                     q_group_Col(:,j), d_grp, rgrpsq)
833 < #else
834 <                call get_interatomic_vector(q_group(:,i), &
835 <                     q_group(:,j), d_grp, rgrpsq)
836 < #endif                
837 <                vij = 0.0d0
838 <                fij(1:3) = 0.0d0
839 <
840 <                call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
841 <                     in_switching_region)
842 <                
843 <                n_in_j = groupStart(j+1) - groupStart(j)
844 <                
845 <                do ia = groupStart(i), groupStart(i+1)-1
846 <                   atom1 = groupList(ia)
847 <                  
848 <                   inner2: do jb = groupStart(j), groupStart(j+1)-1
849 <
850 <                      atom2 = groupList(jb)                        
851 <                      
852 <                      if (skipThisPair(atom1, atom2)) cycle inner2
853 <
854 <                      if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
855 <                         d_atm(1:3) = d_grp(1:3)
856 <                         ratmsq = rgrpsq
857 <                      else
858 < #ifdef IS_MPI
859 <                         call get_interatomic_vector(q_Row(:,atom1), &
860 <                              q_Col(:,atom2), d_atm, ratmsq)
861 < #else
862 <                         call get_interatomic_vector(q(:,atom1), &
863 <                              q(:,atom2), d_atm, ratmsq)
864 < #endif
865 <                      endif
866 < #ifdef IS_MPI                      
867 <                      call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
868 <                           do_pot,  &
869 <                           u_l, A, f, t, pot_local, vpair, fpair)
870 < #else
871 <                      call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
872 <                           do_pot,  &
873 <                           u_l, A, f, t, pot, vpair, fpair)
874 < #endif                      
875 <                      vij = vij + vpair
876 <                      fij(1:3) = fij(1:3) + fpair(1:3)
877 <                                            
878 <                   enddo inner2
879 <                enddo
880 <                
881 <                if (in_switching_region) then
882 <                   swderiv = vij*dswdr/rgrp
883 <                   fij(1) = fij(1) + swderiv*d_grp(1)
884 <                   fij(2) = fij(2) + swderiv*d_grp(2)
885 <                   fij(3) = fij(3) + swderiv*d_grp(3)
886 <  
887 <                   do ia=groupStart(i), groupStart(i+1)-1
888 <                      atom1=groupList(ia)
889 <                      mf = mfact(atom1)
890 < #ifdef IS_MPI
891 <                      f_Row(1,atom1) = f_Row(1,atom1) + swderiv*d_grp(1)*mf
892 <                      f_Row(2,atom1) = f_Row(2,atom1) + swderiv*d_grp(2)*mf
893 <                      f_Row(3,atom1) = f_Row(3,atom1) + swderiv*d_grp(3)*mf
894 < #else
895 <                      f(1,atom1) = f(1,atom1) + swderiv*d_grp(1)*mf
896 <                      f(2,atom1) = f(2,atom1) + swderiv*d_grp(2)*mf
897 <                      f(3,atom1) = f(3,atom1) + swderiv*d_grp(3)*mf
898 < #endif
899 <                   enddo
900 <                  
901 <                   do jb=groupStart(j), groupStart(j+1)-1
902 <                      atom2=groupList(jb)
903 <                      mf = mfact(atom2)
904 < #ifdef IS_MPI
905 <                      f_Col(1,atom2) = f_Col(1,atom2) - swderiv*d_grp(1)*mf
906 <                      f_Col(2,atom2) = f_Col(2,atom2) - swderiv*d_grp(2)*mf
907 <                      f_Col(3,atom2) = f_Col(3,atom2) - swderiv*d_grp(3)*mf
908 < #else
909 <                      f(1,atom2) = f(1,atom2) - swderiv*d_grp(1)*mf
910 <                      f(2,atom2) = f(2,atom2) - swderiv*d_grp(2)*mf
911 <                      f(3,atom2) = f(3,atom2) - swderiv*d_grp(3)*mf
912 < #endif
913 <                   enddo
914 <                endif
915 <
916 <                if (do_stress) call add_stress_tensor(d_grp, fij)
917 <
918 <             enddo
919 <          endif
920 <       enddo
921 <    endif
922 <    
923 <    ! phew, done with main loop.
665 >       if (loop .eq. PREPAIR_LOOP) then
666 >          call do_preforce(nlocal, pot)
667 >       endif
668 >      
669 >    enddo
670      
671      !! Do timing
672   #ifdef PROFILE
# Line 1043 | Line 789 | contains
789      endif
790      
791   #endif
792 <    
1047 <    
792 >      
793    end subroutine do_force_loop
1049
794    
795    subroutine do_pair(i, j, rijsq, d, sw, do_pot, &
796         u_l, A, f, t, pot, vpair, fpair)

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines