********>Bugfix 35: Author: Ross Walker Date: April 29, 2010 Programs: sander Description: QMMM simulations with small boxes can cause certain atoms to be shifted to outside the central box when building the coordinates for the QM region. This can lead to triggering of the cut off vs box size check leading to the calculation stopping. The solution is to shift the QM coordinates to have the center of their coordinates to be at the origin instead of just shifting the first QM atom to be at the origin. Use this patch in amber10/src/sander/qm_mm.f ------------------------------------------------------------------------------ --- qm_mm.f 2010-04-29 01:31:38.000000000 -0700 +++ qm_mm.f 2010-04-29 01:35:18.000000000 -0700 @@ -623,21 +623,40 @@ _REAL_, intent(in), dimension(natom) :: scaled_mm_chrgs _REAL_, intent(out), dimension(3,natom) :: real_scratch - integer j,iqm_one, jqmatom, ier + integer j, jqmatom, ier _REAL_ :: xbnd0,xbnd1,ybnd0,ybnd1,zbnd0,zbnd1 - _REAL_ :: offset(3), frac(3), xx,yy,zz + _REAL_ :: offset(3), frac(3), xx,yy,zz,xtmp,ytmp,ztmp,one_nquant _REAL_ , dimension(6) :: bxbnd integer :: m,n, n1 _REAL_ :: fbndx0, fbndx1, fbndy0, fbndy1, fbndz0, fbndz1, dx2 logical :: include_atom - real_scratch(1:3,1) = zero - !Move the first QM atom to be at the origin. - iqm_one = iqmatoms(1) - frac(1:3) = x(1,iqm_one)*recip(1,1:3)+x(2,iqm_one)*recip(2,1:3)+ & - x(3,iqm_one)*recip(3,1:3) +! iqm_one = iqmatoms(1) +! frac(1:3) = x(1,iqm_one)*recip(1,1:3)+x(2,iqm_one)*recip(2,1:3)+ & +! x(3,iqm_one)*recip(3,1:3) + +!Moving the first QM atom can be an issue if the box size is small +!since it could end up with atoms that can stick outside the box. +!Thus we can change this to calculate the center of the coordinates of +!the QM region and move this instead. This could cause problems if the QM +!region diffuses significantly but if that happens we are in a whole +!world of problems anyway and it should trigger the cutoff too large error. + xtmp=zero; ytmp=zero; ztmp=zero + do j = 1, qmmm_struct%nquant + jqmatom = iqmatoms(j) + xtmp = xtmp + x(1,jqmatom) + ytmp = ytmp + x(2,jqmatom) + ztmp = ztmp + x(3,jqmatom) + end do + one_nquant = 1.0d0/real(qmmm_struct%nquant) + xtmp = xtmp * one_nquant + ytmp = ytmp * one_nquant + ztmp = ztmp * one_nquant + + frac(1:3) = xtmp*recip(1,1:3)+ytmp*recip(2,1:3)+ & + ztmp*recip(3,1:3) frac(1:3) = frac(1:3) - anint(frac(1:3)) offset(1) = frac(1)*ucell(1,1) + frac(2)*ucell(1,2) + frac(3)*ucell(1,3) offset(2) = frac(1)*ucell(2,1) + frac(2)*ucell(2,2) + frac(3)*ucell(2,3) ------------------------------------------------------------------------------ Temporary workarounds: Use a bigger solvent box.