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 2727 by chrisfen, Fri Apr 21 19:32:07 2006 UTC vs.
Revision 2917 by chrisfen, Mon Jul 3 13:18:43 2006 UTC

# Line 45 | Line 45
45  
46   !! @author Charles F. Vardeman II
47   !! @author Matthew Meineke
48 < !! @version $Id: doForces.F90,v 1.80 2006-04-21 19:32:07 chrisfen Exp $, $Date: 2006-04-21 19:32:07 $, $Name: not supported by cvs2svn $, $Revision: 1.80 $
48 > !! @version $Id: doForces.F90,v 1.84 2006-07-03 13:18:43 chrisfen Exp $, $Date: 2006-07-03 13:18:43 $, $Name: not supported by cvs2svn $, $Revision: 1.84 $
49  
50  
51   module doForces
# Line 89 | Line 89 | module doForces
89    logical, save :: haveElectrostaticSummationMethod = .false.
90    logical, save :: haveCutoffPolicy = .false.
91    logical, save :: VisitCutoffsAfterComputing = .false.
92 +  logical, save :: do_box_dipole = .false.
93  
94    logical, save :: FF_uses_DirectionalAtoms
95    logical, save :: FF_uses_Dipoles
# Line 117 | Line 118 | module doForces
118    public :: setCutoffs
119    public :: cWasLame
120    public :: setElectrostaticMethod
121 +  public :: setBoxDipole
122 +  public :: getBoxDipole
123    public :: setCutoffPolicy
124    public :: setSkinThickness
125    public :: do_force_loop
# Line 147 | Line 150 | module doForces
150    end type gtypeCutoffs
151    type(gtypeCutoffs), dimension(:,:), allocatable :: gtypeCutoffMap
152  
153 +  real(kind=dp), dimension(3) :: boxDipole
154  
155   contains
156  
# Line 423 | Line 427 | contains
427      !! largest cutoff for any atypes present in this group.  We also
428      !! create gtypes at this point.
429      
430 <    tol = 1.0d-6
430 >    tol = 1.0e-6_dp
431      nGroupTypesRow = 0
432      nGroupTypesCol = 0
433      do i = istart, iend      
# Line 610 | Line 614 | contains
614      
615     end subroutine setCutoffPolicy
616    
617 +   subroutine setBoxDipole()
618 +
619 +     do_box_dipole = .true.
620 +    
621 +   end subroutine setBoxDipole
622 +
623 +   subroutine getBoxDipole( box_dipole )
624 +
625 +     real(kind=dp), intent(inout), dimension(3) :: box_dipole
626 +
627 +     box_dipole = boxDipole
628 +
629 +   end subroutine getBoxDipole
630 +
631     subroutine setElectrostaticMethod( thisESM )
632  
633       integer, intent(in) :: thisESM
# Line 636 | Line 654 | contains
654       SIM_requires_prepair_calc = SimRequiresPrepairCalc()
655       SIM_uses_PBC = SimUsesPBC()
656       SIM_uses_SC = SimUsesSC()
657 <    
657 >
658       haveSIMvariables = .true.
659      
660       return
# Line 810 | Line 828 | contains
828      integer :: loopStart, loopEnd, loop
829      integer :: iHash
830      integer :: i1
831 <  
831 >
832 >    !! the variables for the box dipole moment
833 > #ifdef IS_MPI
834 >    integer :: pChgCount_local
835 >    integer :: nChgCount_local
836 >    real(kind=dp) :: pChg_local
837 >    real(kind=dp) :: nChg_local
838 >    real(kind=dp), dimension(3) :: pChgPos_local
839 >    real(kind=dp), dimension(3) :: nChgPos_local
840 >    real(kind=dp), dimension(3) :: dipVec_local
841 > #endif
842 >    integer :: pChgCount
843 >    integer :: nChgCount
844 >    real(kind=dp) :: pChg
845 >    real(kind=dp) :: nChg
846 >    real(kind=dp) :: chg_value
847 >    real(kind=dp), dimension(3) :: pChgPos
848 >    real(kind=dp), dimension(3) :: nChgPos
849 >    real(kind=dp), dimension(3) :: dipVec
850 >    real(kind=dp), dimension(3) :: chgVec
851 >
852 >    !! initialize box dipole variables
853 >    if (do_box_dipole) then
854 > #ifdef IS_MPI
855 >       pChg_local = 0.0_dp
856 >       nChg_local = 0.0_dp
857 >       pChgCount_local = 0
858 >       nChgCount_local = 0
859 >       do i=1, 3
860 >          pChgPos_local = 0.0_dp
861 >          nChgPos_local = 0.0_dp
862 >          dipVec_local = 0.0_dp
863 >       enddo
864 > #endif
865 >       pChg = 0.0_dp
866 >       nChg = 0.0_dp
867 >       pChgCount = 0
868 >       nChgCount = 0
869 >       chg_value = 0.0_dp
870 >      
871 >       do i=1, 3
872 >          pChgPos(i) = 0.0_dp
873 >          nChgPos(i) = 0.0_dp
874 >          dipVec(i) = 0.0_dp
875 >          chgVec(i) = 0.0_dp
876 >          boxDipole(i) = 0.0_dp
877 >       enddo
878 >    endif
879  
880      !! initialize local variables  
881  
# Line 962 | Line 1027 | contains
1027  
1028                     rCut = gtypeCutoffMap(groupToGtypeRow(i),groupToGtypeCol(j))%rCut
1029                     if (loop .eq. PAIR_LOOP) then
1030 <                      vij = 0.0d0
1030 >                      vij = 0.0_dp
1031                        fij(1) = 0.0_dp
1032                        fij(2) = 0.0_dp
1033                        fij(3) = 0.0_dp
# Line 1182 | Line 1247 | contains
1247                  ! prevent overcounting of the skips
1248                  if (i.lt.j) then
1249                     call get_interatomic_vector(q(:,i), q(:,j), d_atm, ratmsq)
1250 <                   rVal = dsqrt(ratmsq)
1250 >                   rVal = sqrt(ratmsq)
1251                     call get_switch(ratmsq, sw, dswdr, rVal,in_switching_region)
1252   #ifdef IS_MPI
1253                     call rf_self_excludes(i, j, sw, eFrame, d_atm, rVal, &
# Line 1194 | Line 1259 | contains
1259                  endif
1260               enddo
1261            endif
1262 +
1263 +          if (do_box_dipole) then
1264 + #ifdef IS_MPI
1265 +             call accumulate_box_dipole(i, eFrame, q(:,i), pChg_local, &
1266 +                  nChg_local, pChgPos_local, nChgPos_local, dipVec_local, &
1267 +                  pChgCount_local, nChgCount_local)
1268 + #else
1269 +             call accumulate_box_dipole(i, eFrame, q(:,i), pChg, nChg, &
1270 +                  pChgPos, nChgPos, dipVec, pChgCount, nChgCount)
1271 + #endif
1272 +          endif
1273         enddo
1274      endif
1275 <    
1275 >
1276   #ifdef IS_MPI
1201 –    
1277      if (do_pot) then
1278 <       call mpi_allreduce(pot_local, pot, LR_POT_TYPES,mpi_double_precision,mpi_sum, &
1278 > #ifdef SINGLE_PRECISION
1279 >       call mpi_allreduce(pot_local, pot, LR_POT_TYPES,mpi_real,mpi_sum, &
1280              mpi_comm_world,mpi_err)            
1281 + #else
1282 +       call mpi_allreduce(pot_local, pot, LR_POT_TYPES,mpi_double_precision, &
1283 +            mpi_sum, mpi_comm_world,mpi_err)            
1284 + #endif
1285      endif
1286      
1287      if (do_stress) then
1288 + #ifdef SINGLE_PRECISION
1289 +       call mpi_allreduce(tau_Temp, tau, 9,mpi_real,mpi_sum, &
1290 +            mpi_comm_world,mpi_err)
1291 +       call mpi_allreduce(virial_Temp, virial,1,mpi_real,mpi_sum, &
1292 +            mpi_comm_world,mpi_err)
1293 + #else
1294         call mpi_allreduce(tau_Temp, tau, 9,mpi_double_precision,mpi_sum, &
1295              mpi_comm_world,mpi_err)
1296         call mpi_allreduce(virial_Temp, virial,1,mpi_double_precision,mpi_sum, &
1297              mpi_comm_world,mpi_err)
1298 + #endif
1299      endif
1300      
1301 +    if (do_box_dipole) then
1302 +
1303 + #ifdef SINGLE_PRECISION
1304 +       call mpi_allreduce(pChg_local, pChg, 1, mpi_real, mpi_sum, &
1305 +            mpi_comm_world, mpi_err)
1306 +       call mpi_allreduce(nChg_local, nChg, 1, mpi_real, mpi_sum, &
1307 +            mpi_comm_world, mpi_err)
1308 +       call mpi_allreduce(pChgCount_local, pChgCount, 1, mpi_integer, mpi_sum,&
1309 +            mpi_comm_world, mpi_err)
1310 +       call mpi_allreduce(nChgCount_local, nChgCount, 1, mpi_integer, mpi_sum,&
1311 +            mpi_comm_world, mpi_err)
1312 +       call mpi_allreduce(pChgPos_local, pChgPos, 3, mpi_real, mpi_sum, &
1313 +            mpi_comm_world, mpi_err)
1314 +       call mpi_allreduce(nChgPos_local, nChgPos, 3, mpi_real, mpi_sum, &
1315 +            mpi_comm_world, mpi_err)
1316 +       call mpi_allreduce(dipVec_local, dipVec, 3, mpi_real, mpi_sum, &
1317 +            mpi_comm_world, mpi_err)
1318   #else
1319 +       call mpi_allreduce(pChg_local, pChg, 1, mpi_double_precision, mpi_sum, &
1320 +            mpi_comm_world, mpi_err)
1321 +       call mpi_allreduce(nChg_local, nChg, 1, mpi_double_precision, mpi_sum, &
1322 +            mpi_comm_world, mpi_err)
1323 +       call mpi_allreduce(pChgCount_local, pChgCount, 1, mpi_integer,&
1324 +            mpi_sum, mpi_comm_world, mpi_err)
1325 +       call mpi_allreduce(nChgCount_local, nChgCount, 1, mpi_integer,&
1326 +            mpi_sum, mpi_comm_world, mpi_err)
1327 +       call mpi_allreduce(pChgPos_local, pChgPos, 3, mpi_double_precision, &
1328 +            mpi_sum, mpi_comm_world, mpi_err)
1329 +       call mpi_allreduce(nChgPos_local, nChgPos, 3, mpi_double_precision, &
1330 +            mpi_sum, mpi_comm_world, mpi_err)
1331 +       call mpi_allreduce(dipVec_local, dipVec, 3, mpi_double_precision, &
1332 +            mpi_sum, mpi_comm_world, mpi_err)
1333 + #endif
1334 +
1335 +    endif
1336 +
1337 + #else
1338      
1339      if (do_stress) then
1340         tau = tau_Temp
# Line 1219 | Line 1342 | contains
1342      endif
1343      
1344   #endif
1345 <    
1345 >
1346 >    if (do_box_dipole) then
1347 >       ! first load the accumulated dipole moment (if dipoles were present)
1348 >       boxDipole(1) = dipVec(1)
1349 >       boxDipole(2) = dipVec(2)
1350 >       boxDipole(3) = dipVec(3)
1351 >
1352 >       ! now include the dipole moment due to charges
1353 >       ! use the lesser of the positive and negative charge totals
1354 >       if (nChg .le. pChg) then
1355 >          chg_value = nChg
1356 >       else
1357 >          chg_value = pChg
1358 >       endif
1359 >      
1360 >       ! find the average positions
1361 >       if (pChgCount .gt. 0 .and. nChgCount .gt. 0) then
1362 >          pChgPos = pChgPos / pChgCount
1363 >          nChgPos = nChgPos / nChgCount
1364 >       endif
1365 >
1366 >       ! dipole is from the negative to the positive (physics notation)
1367 >       chgVec(1) = pChgPos(1) - nChgPos(1)
1368 >       chgVec(2) = pChgPos(2) - nChgPos(2)
1369 >       chgVec(3) = pChgPos(3) - nChgPos(3)
1370 >
1371 >       boxDipole(1) = boxDipole(1) + chgVec(1) * chg_value
1372 >       boxDipole(2) = boxDipole(2) + chgVec(2) * chg_value
1373 >       boxDipole(3) = boxDipole(3) + chgVec(3) * chg_value
1374 >
1375 >    endif
1376 >
1377    end subroutine do_force_loop
1378  
1379    subroutine do_pair(i, j, rijsq, d, sw, do_pot, &
# Line 1250 | Line 1404 | contains
1404  
1405      r = sqrt(rijsq)
1406      
1407 <    vpair = 0.0d0
1408 <    fpair(1:3) = 0.0d0
1407 >    vpair = 0.0_dp
1408 >    fpair(1:3) = 0.0_dp
1409  
1410   #ifdef IS_MPI
1411      me_i = atid_row(i)
# Line 1289 | Line 1443 | contains
1443      endif
1444      
1445      if ( iand(iHash, GAYBERNE_LJ).ne.0 ) then
1446 <       call do_gb_lj_pair(i, j, d, r, rijsq, rcut, sw, vpair, fpair, &
1446 >       call do_gb_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1447              pot(VDW_POT), A, f, t, do_pot)
1448      endif
1449      
# Line 1394 | Line 1548 | contains
1548            
1549            ! wrap the scaled coordinates
1550  
1551 <          scaled(1) = scaled(1) - dnint(scaled(1))
1552 <          scaled(2) = scaled(2) - dnint(scaled(2))
1553 <          scaled(3) = scaled(3) - dnint(scaled(3))
1551 >          scaled(1) = scaled(1) - anint(scaled(1), kind=dp)
1552 >          scaled(2) = scaled(2) - anint(scaled(2), kind=dp)
1553 >          scaled(3) = scaled(3) - anint(scaled(3), kind=dp)
1554  
1555            ! calc the wrapped real coordinates from the wrapped scaled
1556            ! coordinates
# Line 1414 | Line 1568 | contains
1568            
1569            ! wrap the scaled coordinates
1570            
1571 <          scaled(1) = scaled(1) - dnint(scaled(1))
1572 <          scaled(2) = scaled(2) - dnint(scaled(2))
1573 <          scaled(3) = scaled(3) - dnint(scaled(3))
1571 >          scaled(1) = scaled(1) - anint(scaled(1), kind=dp)
1572 >          scaled(2) = scaled(2) - anint(scaled(2), kind=dp)
1573 >          scaled(3) = scaled(3) - anint(scaled(3), kind=dp)
1574  
1575            ! calc the wrapped real coordinates from the wrapped scaled
1576            ! coordinates

Diff Legend

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