Sorry, apparently the trajectory was too large. You probably wanted the
output instead anyway.
---------- Forwarded message ----------
From: Brian Radak <radak004.umn.edu>
Date: Thu, Sep 8, 2011 at 10:32 AM
Subject: Re: [AMBER] false positive on too large QM region?
To: AMBER Mailing List <amber.ambermd.org>
Hi Gustavo,
That was sensible advice; apparently I was missing 10 integration steps
before the crash due to a large print step. Attached is a parm7 and a
NetCDF file with the 11 steps before the crash. Below is the input file I
am using. Some probably very important clarifications:
1.) These are actually QM/MM-MM TI simulations; the error occurs in ALL of
them in differing number of steps (they use the same starting structure but
different random seeds)
2.) The error occurs at lambda=0 (all QM) with and without coupling to the
MM potential (i.e. I tried again with icfe = 0 in both inputs)
3.) I have used this exact same protocol with other similar molecules with
success, but have also had other systems where the error presents itself
later in the simulation (say ~1ps in)
4.) Suspecting an imaging problem, I tried iwrap = 0 and iwrap = 1, the
problem seems to be identical in both
I inspected the QM region pdb of one of the simulations and the molecule is
intact.
Brian
QM/MM constant volume equilibration/production
&cntrl
nstlim = 500000, dt = 0.001,
ntpr = 1, ntave = 500000, ! frequency of
printing energies and averages
ntwr = 500000, ntwx = 1, ioutfm = 1, ! frequency of
restart and coord print in NetCDF format
ntb = 1, iwrap = 0, ! const vol PBCs
with wrapping
ntc = 2, ntf = 2, tol = 1e-8, ! Shake all
hydrogens, ignore all bonded terms
cut = 10, ! direct space and
LJ MM cutoff
ntt = 2, temp0 = 310.0, tempi = 310.0, vrand = 2000, ! Andersen
thermostat
ig = 1315422706,
icfe = 1, klambda = 1,
clambda = 0.00,
ifqnt = 1,
noshakemask = ':1' ! Needs to match
MM input for FEP-type calcs!
/
&ewald
ew_type = 0, order = 6,
nfft1 = 36, nfft2 = 36, nfft3 = 36
/
&qmmm
qmmask = ':1', qmcharge = 0,
qm_theory = 'AM1', parameter_file = 'AM1D_HCNOPMGD.PAR',
qm_ewald = 1, qmmm_switch = 1, qmcut = 10,
qmshake = 0,
writepdb = 1
/
On Wed, Sep 7, 2011 at 8:19 PM, Gustavo Seabra <gustavo.seabra.gmail.com>wrote:
> Hi Brian,
>
> As you mentioned, yes, the coordinates are shifted for the QM calculation.
> In fact, the QM zone is put in the center of the box and the whole system is
> wrapped around it, but that only happens for the QM calculation, without
> changing the actual coordinates.
>
> As to what is happening, it would be helpful to see your input and output
> files. Also, it may be a good idea to run the calculation printing the
> coordinates at every step to follow your system. There is also an option
> (can't recall now, but it is in the manual) that tells sander to print the
> QM zone as a PDB file at every step. My guess is that somehow sander is
> including in the QM zone some atom(s) that should not be there.
>
> Good luck,
> Gustavo.
> --
> Sent from my iPad.
>
> On 07/09/2011, at 18:37, Brian Radak <radak004.umn.edu> wrote:
>
> > 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
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
--
================================ 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.
--
================================ 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 Thu Sep 08 2011 - 08:00:03 PDT
