Thanks Ross!
It would appear that "noshakemask" DOES work when applied to QM atoms. For
some reason I did not expect this behavior. Is this a bug that has become a
feature? :P
My reason for using SHAKE is not to increase the time step. Of course that
could be a good idea, but I'm not really sure what the constrained C-H bond
lengths should be (although I'd imagine most reasonable guesses would give
rather similar results). In any event I am trying to keep apples to apples
comparison with my previous simulations and I would think that new
constraints or changes in the time step should give one pause in that
regard.
Brian
On Tue, Jul 26, 2011 at 12:40 PM, Ross Walker <ross.rosswalker.co.uk> wrote:
> Hi Brian,
>
> I have not looked at the code in a while but there is the option for
> noshakemask in the &ctrl namelist which allows you to specify atoms to NOT
> be shaken. This applies to ALL atoms in a MM simulation. I am not sure if
> this will work for QM or not although if you have things setup already it
> would be quick to test. Just add the QM atoms you do not want shaken to
> this
> mask and then run a short MD simulation and take a look at the bond
> lengths,
> are they changing or not?
>
> Essentially what happens in the code is that all the bonds in the QM
> region,
> between QM-QM atoms are deleted during two calls to setbon during setup.
> One
> for bonds without hydrogen and one with. You can find this code in set.f
> inside sander. If the above does not work you could very easily go in and
> hack this code for your specific system, i.e. inside the loop just have an
> 'if atom number = blah then ...' That way you can have it delete the bond,
> instead of preserving the bond if shake is in use, and I think you will
> then
> not be shaking that atom.
>
> Of all the suggestions tweaking the code is probably going to be the
> easiest
> if noshakemask does not work.
>
> Note, you don't gain from having atoms unshaken. E.g. if you have just a
> single hydrogen that is not shaken you will need to run with 1fs timestep
> or
> less. The main reason for wanting to keep shake on some of the atoms is
> because things like TIP3P were parameterized to be used ONLY with shake.
>
> Good luck,
>
> All the best
> Ross
>
> > -----Original Message-----
> > From: Brian Radak [mailto:radak004.umn.edu]
> > Sent: Tuesday, July 26, 2011 8:57 AM
> > To: AMBER Mailing List
> > Subject: [AMBER] selective qmshake?
> >
> > I would like to run QM/MM MD where only *part* of the QM region (waters
> > in
> > TIP3P geometry) is constrained by SHAKE. That is, there are several C-
> > H
> > bonds which I would like to leave unconstrained if I can help it. The
> > reason being, and perhaps people have opinions on this as well, that I
> > would
> > like to compare to other simulations with no QM waters in which I did
> > not
> > use SHAKE.
> >
> > My understanding from the manual is that qmshake is all or none wrt
> > bonds
> > with hydrogen. I can think of two doable, but arguably outlandish,
> > ways to
> > accomplish this. Will either of these work? Is there a better way to
> > do
> > this?
> >
> > 1.) Manually go into the parm7 and change the BONDS_INC_HYDROGEN and
> > BONDS_WITHOUT_HYDROGEN tables... yuck! and will probably fail at first
> >
> > 2.) Create a new atom type for hydrogens bound to carbons (maybe I can
> > steal
> > one of the halide symbols?) so that leap does not recognize them as
> > bonds
> > with hydrogen. Of course this new atom type would have the same mass
> > as
> > hydrogen. Is that enough to trick leap?
> >
> > Thanks,
> > Brian
> >
> >
> > --
> > ================================ 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.
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Tue Jul 26 2011 - 12:00:02 PDT
