ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/UseTheForce/doForces.F90
(Generate patch)

Comparing trunk/OOPSE-4/src/UseTheForce/doForces.F90 (file contents):
Revision 2398 by chrisfen, Wed Oct 26 23:31:18 2005 UTC vs.
Revision 2432 by chuckv, Tue Nov 15 16:01:06 2005 UTC

# Line 45 | Line 45
45  
46   !! @author Charles F. Vardeman II
47   !! @author Matthew Meineke
48 < !! @version $Id: doForces.F90,v 1.63 2005-10-26 23:31:18 chrisfen Exp $, $Date: 2005-10-26 23:31:18 $, $Name: not supported by cvs2svn $, $Revision: 1.63 $
48 > !! @version $Id: doForces.F90,v 1.68 2005-11-15 16:01:06 chuckv Exp $, $Date: 2005-11-15 16:01:06 $, $Name: not supported by cvs2svn $, $Revision: 1.68 $
49  
50  
51   module doForces
# Line 62 | Line 62 | module doForces
62    use shapes
63    use vector_class
64    use eam
65 +  use suttonchen
66    use status
67   #ifdef IS_MPI
68    use mpiSimulation
# Line 92 | Line 93 | module doForces
93    logical, save :: FF_uses_Dipoles
94    logical, save :: FF_uses_GayBerne
95    logical, save :: FF_uses_EAM
96 +  logical, save :: FF_uses_SC
97 +  logical, save :: FF_uses_MEAM
98 +
99  
100    logical, save :: SIM_uses_DirectionalAtoms
101    logical, save :: SIM_uses_EAM
102 +  logical, save :: SIM_uses_SC
103 +  logical, save :: SIM_uses_MEAM
104    logical, save :: SIM_requires_postpair_calc
105    logical, save :: SIM_requires_prepair_calc
106    logical, save :: SIM_uses_PBC
# Line 158 | Line 164 | contains
164      logical :: i_is_GB
165      logical :: i_is_EAM
166      logical :: i_is_Shape
167 +    logical :: i_is_SC
168 +    logical :: i_is_MEAM
169      logical :: j_is_LJ
170      logical :: j_is_Elect
171      logical :: j_is_Sticky
# Line 165 | Line 173 | contains
173      logical :: j_is_GB
174      logical :: j_is_EAM
175      logical :: j_is_Shape
176 +    logical :: j_is_SC
177 +    logical :: j_is_MEAM
178      real(kind=dp) :: myRcut
179  
180 +
181      status = 0  
182  
183      if (.not. associated(atypes)) then
# Line 204 | Line 215 | contains
215         call getElementProperty(atypes, i, "is_GayBerne", i_is_GB)
216         call getElementProperty(atypes, i, "is_EAM", i_is_EAM)
217         call getElementProperty(atypes, i, "is_Shape", i_is_Shape)
218 +       call getElementProperty(atypes, i, "is_SC", i_is_SC)
219 +       call getElementProperty(atypes, i, "is_MEAM", i_is_MEAM)
220  
221         do j = i, nAtypes
222  
# Line 217 | Line 230 | contains
230            call getElementProperty(atypes, j, "is_GayBerne", j_is_GB)
231            call getElementProperty(atypes, j, "is_EAM", j_is_EAM)
232            call getElementProperty(atypes, j, "is_Shape", j_is_Shape)
233 +          call getElementProperty(atypes, j, "is_SC", j_is_SC)
234 +          call getElementProperty(atypes, j, "is_MEAM", j_is_MEAM)
235  
236            if (i_is_LJ .and. j_is_LJ) then
237               iHash = ior(iHash, LJ_PAIR)            
# Line 238 | Line 253 | contains
253               iHash = ior(iHash, EAM_PAIR)
254            endif
255  
256 +          if (i_is_SC .and. j_is_SC) then
257 +             iHash = ior(iHash, SC_PAIR)
258 +          endif
259 +
260            if (i_is_GB .and. j_is_GB) iHash = ior(iHash, GAYBERNE_PAIR)
261            if (i_is_GB .and. j_is_LJ) iHash = ior(iHash, GAYBERNE_LJ)
262            if (i_is_LJ .and. j_is_GB) iHash = ior(iHash, GAYBERNE_LJ)
# Line 569 | Line 588 | contains
588    subroutine setSimVariables()
589      SIM_uses_DirectionalAtoms = SimUsesDirectionalAtoms()
590      SIM_uses_EAM = SimUsesEAM()
591 +    SIM_uses_SC  = SimUsesSC()
592      SIM_requires_postpair_calc = SimRequiresPostpairCalc()
593      SIM_requires_prepair_calc = SimRequiresPrepairCalc()
594      SIM_uses_PBC = SimUsesPBC()
# Line 644 | Line 664 | contains
664      integer, intent(out) :: thisStat  
665      integer :: my_status, nMatches
666      integer, pointer :: MatchList(:) => null()
647 –    real(kind=dp) :: rcut, rrf, rt, dielect
667  
668      !! assume things are copacetic, unless they aren't
669      thisStat = 0
# Line 756 | Line 775 | contains
775      integer :: loopStart, loopEnd, loop
776      integer :: iHash
777      integer :: i1
759 –    logical :: indirect_only
778    
779  
780      !! initialize local variables  
# Line 904 | Line 922 | contains
922  
923                     list(nlist) = j
924                  endif
925 <
926 <                if (loop .eq. PAIR_LOOP) then
909 <                   vij = 0.0d0
910 <                   fij(1:3) = 0.0d0
911 <                endif
925 >                
926 >                if (rgrpsq < gtypeCutoffMap(groupToGtypeRow(i),groupToGtypeCol(j))%rCutsq) then
927  
928 <                call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
929 <                     in_switching_region)
930 <
931 <                n_in_j = groupStartCol(j+1) - groupStartCol(j)
932 <
933 <                do ia = groupStartRow(i), groupStartRow(i+1)-1
934 <
935 <                   atom1 = groupListRow(ia)
936 <
937 <                   inner: do jb = groupStartCol(j), groupStartCol(j+1)-1
938 <
939 <                      atom2 = groupListCol(jb)
940 <
941 <                      indirect_only = .false.
942 <    
943 <                      if (skipThisPair(atom1, atom2)) then
944 <                         if (electrostaticSummationMethod.ne.REACTION_FIELD) then
945 <                            cycle inner
946 <                         else
947 <                            indirect_only = .true.
948 <                         endif
949 <                      endif
950 <    
951 <
937 <                      if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
938 <                         d_atm(1:3) = d_grp(1:3)
939 <                         ratmsq = rgrpsq
940 <                      else
928 >                   if (loop .eq. PAIR_LOOP) then
929 >                      vij = 0.0d0
930 >                      fij(1:3) = 0.0d0
931 >                   endif
932 >                  
933 >                   call get_switch(rgrpsq, sw, dswdr, rgrp, group_switch, &
934 >                        in_switching_region)
935 >                  
936 >                   n_in_j = groupStartCol(j+1) - groupStartCol(j)
937 >                  
938 >                   do ia = groupStartRow(i), groupStartRow(i+1)-1
939 >                      
940 >                      atom1 = groupListRow(ia)
941 >                      
942 >                      inner: do jb = groupStartCol(j), groupStartCol(j+1)-1
943 >                        
944 >                         atom2 = groupListCol(jb)
945 >                        
946 >                         if (skipThisPair(atom1, atom2))  cycle inner
947 >                        
948 >                         if ((n_in_i .eq. 1).and.(n_in_j .eq. 1)) then
949 >                            d_atm(1:3) = d_grp(1:3)
950 >                            ratmsq = rgrpsq
951 >                         else
952   #ifdef IS_MPI
953 <                         call get_interatomic_vector(q_Row(:,atom1), &
954 <                              q_Col(:,atom2), d_atm, ratmsq)
953 >                            call get_interatomic_vector(q_Row(:,atom1), &
954 >                                 q_Col(:,atom2), d_atm, ratmsq)
955   #else
956 <                         call get_interatomic_vector(q(:,atom1), &
957 <                              q(:,atom2), d_atm, ratmsq)
956 >                            call get_interatomic_vector(q(:,atom1), &
957 >                                 q(:,atom2), d_atm, ratmsq)
958   #endif
959 <                      endif
960 <
961 <                      if (loop .eq. PREPAIR_LOOP) then
959 >                         endif
960 >                        
961 >                         if (loop .eq. PREPAIR_LOOP) then
962   #ifdef IS_MPI                      
963 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
964 <                              rgrpsq, d_grp, do_pot, do_stress, &
965 <                              eFrame, A, f, t, pot_local)
966 < #else
967 <                         call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
968 <                              rgrpsq, d_grp, do_pot, do_stress, &
969 <                              eFrame, A, f, t, pot)
963 >                            call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
964 >                                 rgrpsq, d_grp, do_pot, do_stress, &
965 >                                 eFrame, A, f, t, pot_local)
966 > #else
967 >                            call do_prepair(atom1, atom2, ratmsq, d_atm, sw, &
968 >                                 rgrpsq, d_grp, do_pot, do_stress, &
969 >                                 eFrame, A, f, t, pot)
970   #endif                                              
971 <                      else
971 >                         else
972   #ifdef IS_MPI                      
973 <                         call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
974 <                              do_pot, eFrame, A, f, t, pot_local, vpair, &
975 <                              fpair, d_grp, rgrp, indirect_only)
973 >                            call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
974 >                                 do_pot, eFrame, A, f, t, pot_local, vpair, &
975 >                                 fpair, d_grp, rgrp)
976   #else
977 <                         call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
978 <                              do_pot, eFrame, A, f, t, pot, vpair, fpair, &
979 <                              d_grp, rgrp, indirect_only)
977 >                            call do_pair(atom1, atom2, ratmsq, d_atm, sw, &
978 >                                 do_pot, eFrame, A, f, t, pot, vpair, fpair, &
979 >                                 d_grp, rgrp)
980   #endif
981 +                            vij = vij + vpair
982 +                            fij(1:3) = fij(1:3) + fpair(1:3)
983 +                         endif
984 +                      enddo inner
985 +                   enddo
986  
987 <                         vij = vij + vpair
988 <                         fij(1:3) = fij(1:3) + fpair(1:3)
989 <                      endif
990 <                   enddo inner
991 <                enddo
992 <
993 <                if (loop .eq. PAIR_LOOP) then
994 <                   if (in_switching_region) then
995 <                      swderiv = vij*dswdr/rgrp
996 <                      fij(1) = fij(1) + swderiv*d_grp(1)
981 <                      fij(2) = fij(2) + swderiv*d_grp(2)
982 <                      fij(3) = fij(3) + swderiv*d_grp(3)
983 <
984 <                      do ia=groupStartRow(i), groupStartRow(i+1)-1
985 <                         atom1=groupListRow(ia)
986 <                         mf = mfactRow(atom1)
987 >                   if (loop .eq. PAIR_LOOP) then
988 >                      if (in_switching_region) then
989 >                         swderiv = vij*dswdr/rgrp
990 >                         fij(1) = fij(1) + swderiv*d_grp(1)
991 >                         fij(2) = fij(2) + swderiv*d_grp(2)
992 >                         fij(3) = fij(3) + swderiv*d_grp(3)
993 >                        
994 >                         do ia=groupStartRow(i), groupStartRow(i+1)-1
995 >                            atom1=groupListRow(ia)
996 >                            mf = mfactRow(atom1)
997   #ifdef IS_MPI
998 <                         f_Row(1,atom1) = f_Row(1,atom1) + swderiv*d_grp(1)*mf
999 <                         f_Row(2,atom1) = f_Row(2,atom1) + swderiv*d_grp(2)*mf
1000 <                         f_Row(3,atom1) = f_Row(3,atom1) + swderiv*d_grp(3)*mf
998 >                            f_Row(1,atom1) = f_Row(1,atom1) + swderiv*d_grp(1)*mf
999 >                            f_Row(2,atom1) = f_Row(2,atom1) + swderiv*d_grp(2)*mf
1000 >                            f_Row(3,atom1) = f_Row(3,atom1) + swderiv*d_grp(3)*mf
1001   #else
1002 <                         f(1,atom1) = f(1,atom1) + swderiv*d_grp(1)*mf
1003 <                         f(2,atom1) = f(2,atom1) + swderiv*d_grp(2)*mf
1004 <                         f(3,atom1) = f(3,atom1) + swderiv*d_grp(3)*mf
1002 >                            f(1,atom1) = f(1,atom1) + swderiv*d_grp(1)*mf
1003 >                            f(2,atom1) = f(2,atom1) + swderiv*d_grp(2)*mf
1004 >                            f(3,atom1) = f(3,atom1) + swderiv*d_grp(3)*mf
1005   #endif
1006 <                      enddo
1007 <
1008 <                      do jb=groupStartCol(j), groupStartCol(j+1)-1
1009 <                         atom2=groupListCol(jb)
1010 <                         mf = mfactCol(atom2)
1006 >                         enddo
1007 >                        
1008 >                         do jb=groupStartCol(j), groupStartCol(j+1)-1
1009 >                            atom2=groupListCol(jb)
1010 >                            mf = mfactCol(atom2)
1011   #ifdef IS_MPI
1012 <                         f_Col(1,atom2) = f_Col(1,atom2) - swderiv*d_grp(1)*mf
1013 <                         f_Col(2,atom2) = f_Col(2,atom2) - swderiv*d_grp(2)*mf
1014 <                         f_Col(3,atom2) = f_Col(3,atom2) - swderiv*d_grp(3)*mf
1012 >                            f_Col(1,atom2) = f_Col(1,atom2) - swderiv*d_grp(1)*mf
1013 >                            f_Col(2,atom2) = f_Col(2,atom2) - swderiv*d_grp(2)*mf
1014 >                            f_Col(3,atom2) = f_Col(3,atom2) - swderiv*d_grp(3)*mf
1015   #else
1016 <                         f(1,atom2) = f(1,atom2) - swderiv*d_grp(1)*mf
1017 <                         f(2,atom2) = f(2,atom2) - swderiv*d_grp(2)*mf
1018 <                         f(3,atom2) = f(3,atom2) - swderiv*d_grp(3)*mf
1016 >                            f(1,atom2) = f(1,atom2) - swderiv*d_grp(1)*mf
1017 >                            f(2,atom2) = f(2,atom2) - swderiv*d_grp(2)*mf
1018 >                            f(3,atom2) = f(3,atom2) - swderiv*d_grp(3)*mf
1019   #endif
1020 <                      enddo
1021 <                   endif
1020 >                         enddo
1021 >                      endif
1022  
1023 <                   if (do_stress) call add_stress_tensor(d_grp, fij)
1023 >                      if (do_stress) call add_stress_tensor(d_grp, fij)
1024 >                   endif
1025                  endif
1026 <             end if
1026 >             endif
1027            enddo
1028 <
1028 >          
1029         enddo outer
1030  
1031         if (update_nlist) then
# Line 1110 | Line 1121 | contains
1121  
1122            if ( iand(iHash, ELECTROSTATIC_PAIR).ne.0 ) then
1123   #ifdef IS_MPI
1124 <             call rf_self_self(i, eFrame, pot_local(ELECTROSTATIC_POT), &
1124 >             call self_self(i, eFrame, pot_local(ELECTROSTATIC_POT), &
1125                    t, do_pot)
1126   #else
1127 <             call rf_self_self(i, eFrame, pot(ELECTROSTATIC_POT), &
1127 >             call self_self(i, eFrame, pot(ELECTROSTATIC_POT), &
1128                    t, do_pot)
1129   #endif
1130            endif
1131    
1132 <          ! loop over the excludes to accumulate any additional RF components
1133 <
1123 <          do i1 = 1, nSkipsForAtom(i)
1124 <             j = skipsForAtom(i, i1)
1132 >          
1133 >          if (electrostaticSummationMethod.eq.REACTION_FIELD) then
1134              
1135 <             ! prevent overcounting of the skips
1136 <             if (i.lt.j) then
1137 <                call get_interatomic_vector(q(:,i), &
1138 <                     q(:,j), d_atm, ratmsq)
1139 <                rVal = dsqrt(ratmsq)
1140 <                call get_switch(ratmsq, sw, dswdr, rVal, group_switch, &
1141 <                     in_switching_region)
1135 >             ! loop over the excludes to accumulate RF stuff we've
1136 >             ! left out of the normal pair loop
1137 >            
1138 >             do i1 = 1, nSkipsForAtom(i)
1139 >                j = skipsForAtom(i, i1)
1140 >                
1141 >                ! prevent overcounting of the skips
1142 >                if (i.lt.j) then
1143 >                   call get_interatomic_vector(q(:,i), &
1144 >                        q(:,j), d_atm, ratmsq)
1145 >                   rVal = dsqrt(ratmsq)
1146 >                   call get_switch(ratmsq, sw, dswdr, rVal, group_switch, &
1147 >                        in_switching_region)
1148   #ifdef IS_MPI
1149 <                call rf_self_excludes(i, j, sw, eFrame, d_atm, rVal, vpair, &
1150 <                     pot_local(ELECTROSTATIC_POT), f, t, do_pot)
1149 >                   call rf_self_excludes(i, j, sw, eFrame, d_atm, rVal, &
1150 >                        vpair, pot_local(ELECTROSTATIC_POT), f, t, do_pot)
1151   #else
1152 <                call rf_self_excludes(i, j, sw, eFrame, d_atm, rVal, vpair, &
1153 <                     pot(ELECTROSTATIC_POT), f, t, do_pot)
1152 >                   call rf_self_excludes(i, j, sw, eFrame, d_atm, rVal, &
1153 >                        vpair, pot(ELECTROSTATIC_POT), f, t, do_pot)
1154   #endif
1155 <             endif
1156 <          enddo
1157 <       enddo      
1155 >                endif
1156 >             enddo
1157 >          endif
1158 >       enddo
1159      endif
1160      
1161   #ifdef IS_MPI
# Line 1168 | Line 1184 | contains
1184    end subroutine do_force_loop
1185  
1186    subroutine do_pair(i, j, rijsq, d, sw, do_pot, &
1187 <       eFrame, A, f, t, pot, vpair, fpair, d_grp, r_grp, indirect_only)
1187 >       eFrame, A, f, t, pot, vpair, fpair, d_grp, r_grp)
1188  
1189      real( kind = dp ) :: vpair, sw
1190      real( kind = dp ), dimension(LR_POT_TYPES) :: pot
# Line 1180 | Line 1196 | contains
1196      real( kind = dp ), dimension(3,nLocal) :: t
1197  
1198      logical, intent(inout) :: do_pot
1183 –    logical, intent(inout) :: indirect_only
1199      integer, intent(in) :: i, j
1200      real ( kind = dp ), intent(inout) :: rijsq
1201      real ( kind = dp ), intent(inout) :: r_grp
# Line 1204 | Line 1219 | contains
1219   #endif
1220  
1221      iHash = InteractionHash(me_i, me_j)
1222 <
1223 <    if (indirect_only) then
1224 <       if ( iand(iHash, ELECTROSTATIC_PAIR).ne.0 ) then
1225 <          call doElectrostaticPair(i, j, d, r, rijsq, sw, vpair, fpair, &
1226 <               pot(ELECTROSTATIC_POT), eFrame, f, t, do_pot, indirect_only)
1227 <       endif
1228 <    else
1229 <          
1230 <       if ( iand(iHash, LJ_PAIR).ne.0 ) then
1231 <          call do_lj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1232 <               pot(VDW_POT), f, do_pot)
1233 <       endif
1222 >    
1223 >    if ( iand(iHash, LJ_PAIR).ne.0 ) then
1224 >       call do_lj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1225 >            pot(VDW_POT), f, do_pot)
1226 >    endif
1227 >    
1228 >    if ( iand(iHash, ELECTROSTATIC_PAIR).ne.0 ) then
1229 >       call doElectrostaticPair(i, j, d, r, rijsq, sw, vpair, fpair, &
1230 >            pot(ELECTROSTATIC_POT), eFrame, f, t, do_pot)
1231 >    endif
1232 >    
1233 >    if ( iand(iHash, STICKY_PAIR).ne.0 ) then
1234 >       call do_sticky_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1235 >            pot(HB_POT), A, f, t, do_pot)
1236 >    endif
1237 >    
1238 >    if ( iand(iHash, STICKYPOWER_PAIR).ne.0 ) then
1239 >       call do_sticky_power_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1240 >            pot(HB_POT), A, f, t, do_pot)
1241 >    endif
1242 >    
1243 >    if ( iand(iHash, GAYBERNE_PAIR).ne.0 ) then
1244 >       call do_gb_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1245 >            pot(VDW_POT), A, f, t, do_pot)
1246 >    endif
1247 >    
1248 >    if ( iand(iHash, GAYBERNE_LJ).ne.0 ) then
1249 >       call do_gb_lj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1250 >            pot(VDW_POT), A, f, t, do_pot)
1251 >    endif
1252 >    
1253 >    if ( iand(iHash, EAM_PAIR).ne.0 ) then      
1254 >       call do_eam_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1255 >            pot(METALLIC_POT), f, do_pot)
1256 >    endif
1257 >    
1258 >    if ( iand(iHash, SHAPE_PAIR).ne.0 ) then      
1259 >       call do_shape_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1260 >            pot(VDW_POT), A, f, t, do_pot)
1261 >    endif
1262 >    
1263 >    if ( iand(iHash, SHAPE_LJ).ne.0 ) then      
1264 >       call do_shape_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1265 >            pot(VDW_POT), A, f, t, do_pot)
1266 >    endif
1267  
1268 <       if ( iand(iHash, ELECTROSTATIC_PAIR).ne.0 ) then
1269 <          call doElectrostaticPair(i, j, d, r, rijsq, sw, vpair, fpair, &
1270 <               pot(ELECTROSTATIC_POT), eFrame, f, t, do_pot, indirect_only)
1223 <       endif
1224 <      
1225 <       if ( iand(iHash, STICKY_PAIR).ne.0 ) then
1226 <          call do_sticky_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1227 <               pot(HB_POT), A, f, t, do_pot)
1228 <       endif
1229 <      
1230 <       if ( iand(iHash, STICKYPOWER_PAIR).ne.0 ) then
1231 <          call do_sticky_power_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1232 <               pot(HB_POT), A, f, t, do_pot)
1233 <       endif
1234 <      
1235 <       if ( iand(iHash, GAYBERNE_PAIR).ne.0 ) then
1236 <          call do_gb_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1237 <               pot(VDW_POT), A, f, t, do_pot)
1238 <       endif
1239 <      
1240 <       if ( iand(iHash, GAYBERNE_LJ).ne.0 ) then
1241 <          call do_gb_lj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1242 <               pot(VDW_POT), A, f, t, do_pot)
1243 <       endif
1244 <      
1245 <       if ( iand(iHash, EAM_PAIR).ne.0 ) then      
1246 <          call do_eam_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1247 <               pot(METALLIC_POT), f, do_pot)
1248 <       endif
1249 <      
1250 <       if ( iand(iHash, SHAPE_PAIR).ne.0 ) then      
1251 <          call do_shape_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1252 <               pot(VDW_POT), A, f, t, do_pot)
1253 <       endif
1254 <      
1255 <       if ( iand(iHash, SHAPE_LJ).ne.0 ) then      
1256 <          call do_shape_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1257 <               pot(VDW_POT), A, f, t, do_pot)
1258 <       endif
1268 >    if ( iand(iHash, SC_PAIR).ne.0 ) then      
1269 >       call do_SC_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1270 >            pot(METALLIC_POT), f, do_pot)
1271      endif
1272 +
1273      
1274 +    
1275    end subroutine do_pair
1276  
1277    subroutine do_prepair(i, j, rijsq, d, sw, rcijsq, dc, &
# Line 1293 | Line 1307 | contains
1307      if ( iand(iHash, EAM_PAIR).ne.0 ) then      
1308              call calc_EAM_prepair_rho(i, j, d, r, rijsq )
1309      endif
1310 +
1311 +    if ( iand(iHash, SC_PAIR).ne.0 ) then      
1312 +            call calc_SC_prepair_rho(i, j, d, r, rijsq )
1313 +    endif
1314      
1315    end subroutine do_prepair
1316  
# Line 1303 | Line 1321 | contains
1321  
1322      if (FF_uses_EAM .and. SIM_uses_EAM) then
1323         call calc_EAM_preforce_Frho(nlocal,pot(METALLIC_POT))
1324 +    endif
1325 +    if (FF_uses_SC .and. SIM_uses_SC) then
1326 +       call calc_SC_preforce_Frho(nlocal,pot(METALLIC_POT))
1327      endif
1328  
1329  
# Line 1388 | Line 1409 | contains
1409      pot_Col = 0.0_dp
1410      pot_Temp = 0.0_dp
1411  
1391 –    rf_Row = 0.0_dp
1392 –    rf_Col = 0.0_dp
1393 –    rf_Temp = 0.0_dp
1394 –
1412   #endif
1413  
1414      if (FF_uses_EAM .and. SIM_uses_EAM) then
1415         call clean_EAM()
1416      endif
1417  
1401 –    rf = 0.0_dp
1418      tau_Temp = 0.0_dp
1419      virial_Temp = 0.0_dp
1420    end subroutine zero_work_arrays
# Line 1492 | Line 1508 | contains
1508  
1509    function FF_RequiresPrepairCalc() result(doesit)
1510      logical :: doesit
1511 <    doesit = FF_uses_EAM
1511 >    doesit = FF_uses_EAM .or. FF_uses_SC &
1512 >         .or. FF_uses_MEAM
1513    end function FF_RequiresPrepairCalc
1514  
1515   #ifdef PROFILE

Diff Legend

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