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 2787 by gezelter, Mon Jun 5 18:24:45 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.83 2006-06-05 18:24:45 gezelter Exp $, $Date: 2006-06-05 18:24:45 $, $Name: not supported by cvs2svn $, $Revision: 1.83 $
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 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
813 –  
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  
882   #ifdef IS_MPI
# Line 1193 | Line 1258 | contains
1258   #endif
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   #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,mpi_sum, &
1283 <            mpi_comm_world,mpi_err)            
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      
# Line 1223 | Line 1298 | contains
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
# Line 1231 | 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, &

Diff Legend

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