I am running many QM/MM NVT simulations of small organic molecules (QM) in
water (MM) using rhombic dodecahedral PBCs. The starting structures come
from 1 ns NpT MM equilibrations with the volume reset to a rounded value
(these are for TI calculations and I thought it sensible and not
unreasonable to set the volume constant amongst all lambda values, the
change is usually less than 0.5%). Sometimes, but not always, sander bombs
out with the following error, either after many steps of QM/MM MD or
sometimes after no steps at all:
****************************************************
ERROR: QM region + cutoff larger than box dimension:
QM-MM Cutoff = 10.0000
Coord Lower Upper Size Radius of largest sphere inside unit
cell
X -11.486 11.276 22.762 12.092
Y -12.143 12.058 24.200 12.092
Z -11.619 11.293 22.911 12.092
****************************************************
SANDER BOMB in subroutine QM_CHECK_PERIODIC<qm_mm.f>
QM region + cutoff larger than box
cannot continue, need larger box.
As a specific example this happens to isopropanol (generously ~6 A wide) in
a 34 A RHDO. I don't think I am wrong in concluding that there is easily 12
A (qm_cut + skinnb) to the nearest RHDO face. In any event I got no such
error during the MM simulation.
I looked in the code and I don't fully understand how the table of values is
arrived at. It would seem that the system is shifted by the geometric
center of mass of the QM region (which is why "Lower" is negative and
"Upper" is positive?), but I'm not exactly sure what the reason for this is
(I'm sure it's a very good one, I just don't know what it is). Anyway,
somewhere in that shifting, the molecule seems to get split across a
periodic wall and thus the "size" of the molecule is calculated
erroneously. There is also a comment in the code that something bad (this?)
would happen if the QM region diffuses extensively.
Can anyone explain what's going on here? Has anyone seen this problem
before? I'm willing to give fixing this a try, but I think I need a little
bit broader vision of the code as this extends across sqm, the image code,
and the ewald code.
Thanks!
Brian
P.S. On a related note, 12.092 A seems too small for the largest sphere in
the "familiar" form of a RHDO. Is this a consequence of the triclinic form
used in the simulation?
--
================================ Current Address =======================
Brian Radak : BioMaPS
Institute for Quantitative Biology
PhD candidate - York Research Group : Rutgers, The State
University of New Jersey
University of Minnesota - Twin Cities : Wright-Rieman Hall 101
Graduate Program in Chemical Physics : 610 Taylor Road,
Department of Chemistry : Piscataway, NJ
08854-8066
radak004.umn.edu :
radakb.biomaps.rutgers.edu
====================================================================
Sorry for the multiple e-mail addresses, just use the institute appropriate
address.
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Sep 07 2011 - 15:00:03 PDT
