--- trunk/OOPSE/libmdtools/notifyCutoffs.F90 2003/08/07 00:47:33 669 +++ trunk/OOPSE/libmdtools/notifyCutoffs.F90 2004/06/04 02:38:23 1233 @@ -2,37 +2,87 @@ module notifyCutoffs use definitions use do_Forces, only: setRlistDF - use dipole_dipole, only: setCutoffsDipole use reaction_field, only: setCutoffsRF use lj, only: setCutoffLJ use eam, only: setCutoffEAM + use switcheroo, only: set_switch implicit none PRIVATE + +#define __FORTRAN90 +#include "fSwitchingFunction.h" +#include "simError.h" public::cutoffNotify contains - subroutine cutoffNotify( this_rcut, this_rlist, this_ecr, this_est ) + subroutine cutoffNotify( this_rcut, this_rsw, this_rlist ) - real(kind=dp), intent(in) :: this_rcut, this_rlist, this_ecr, this_est + real(kind=dp), intent(in) :: this_rcut, this_rsw, this_rlist - real(kind=dp) :: rtaper, rcut, rlist, ecr + real(kind=dp) :: rsw, rcut, rlist integer :: localError - - + logical :: do_shift + rcut = this_rcut + rsw = this_rsw rlist = this_rlist - ecr = this_ecr - rtaper = this_ecr - this_est + if (rcut .lt. rsw) then + write(painCave%errMsg, *) 'cutoffRadius is ', rcut, & + achar(10) // achar(9), & + 'but switchingRadius is set larger at ', rsw , & + achar(10) // achar(9) , & + 'That is probably not what you wanted to do!', & + achar(10) // achar(0) + + painCave%severity = OOPSE_WARNING + painCave%isFatal = .false. + call c_simError(painCave) + + endif + + if (rlist .lt. rcut) then + + write(painCave%errMsg, *) 'neighborListRadius is ', rlist, & + achar(10) // achar(9), & + 'but cutoffRadius is set larger at ', rcut , & + achar(10) // achar(9) , & + 'That is probably a programming error!', & + achar(10) // achar(0) + painCave%severity = OOPSE_WARNING + painCave%isFatal = .false. + call c_simError(painCave) + + endif + + do_shift = .false. + if (abs(rcut-rsw) .lt. 0.0001) then + + write(painCave%errMsg,*) 'cutoffRadius and switchingRadius', & + achar(10) // achar(9), & + 'are set to the same value. OOPSE will use', & + achar(10) // achar(9), & + 'shifted Lennard-Jones potentials instead of', & + achar(10) // achar(9), & + 'switching functions.', & + achar(10) // achar(0) + painCave%severity = OOPSE_INFO + painCave%isFatal = .false. + call c_simError(painCave) + + do_shift = .true. + + endif + call setRlistDF( rlist ) - call setCutoffsDipole( ecr, rtaper ) - call setCutoffsRF( ecr, rtaper ) - call setCutoffLJ( rcut, localError ) - call setCutoffEAM(rcut,localError) + call setCutoffsRF( rcut, rsw ) + call setCutoffLJ( rcut, do_shift, localError ) + call setCutoffEAM( rcut, localError) + call set_switch(GROUP_SWITCH, rsw, rcut) end subroutine cutoffNotify