********>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.