Re: [AMBER] Temperature and TI with qm/mm

From: Fabrício Bracht <>
Date: Sat, 9 Nov 2013 16:01:18 -0200

Hello. Indeed, qmshake was off, but turning it on didn't do a thing. But I
have noticed something else, that might be more related to my system than
anything else. The simulation started ok (except for the temperature
thing), but after a few hundreds of steps, the system in V1 presented a
very abnormal structural behaviour. What I mean by abnormal is, when
comparing to pure qm/mm simulations (no TI), the hydroxyl group bound to
the zinc center went completely crazy. This is actually something that I
have noticed over and over again. If I happen to exchange V0 and V1, the
other structure goes crazy. So, no matter what choice I make (appear or
disappearing of atoms) V1 always goes crazy. V0 stays ok though.
I took the liberty of testing this on behaviour on the malonaldehyde
structure from the EVB test. I chose these structures because it seemed
like a good starting point for a small test on something that you can also
easily check (if you wish, of course).
If I run a qm/mm simulation (solo, no TI), with either or,
both simulations run fine. Temperature stays around the 300K that was set
and the structure holds like I would expect it to hold (no bonds going
crazy). If I turn icfe and ifcs on, than the same thing happens here. It
does not matter which structure I choose to be V1, it always goes
bananas.As for the temperature, it is hard to say if the same thing happens
here, since the temperature oscillates wildly for such a small system, but
the average seems to be a lot higher (maybe due to the wild fluctuations in
V1?), around 500K or so.
In case you want to test things for yourselves there, here is the mdin file
I used

  imin = 0,
  irest = 1,
  ntx = 7,
  ntb = 0,
  cut = 9999.0,
  tempi = 300.0,
  temp0 = 300.0,
  ntt = 3,
  gamma_ln = 1.0,
  nstlim = 20000, dt = 0.0005, nscm = 1000,
  ntpr = 1, ntwx = 100, ntwr = 100
  icfe = 1, ifsc = 1, clambda = 0.01592,
  noshakemask = '', scmask = ':1', idecomp = 0,
  ifqnt=1, nmropt = 0,
 dftb_3rd_order = 'PA',
 qmcut=12.0, vsolv=0,adjust_q=2,
The input files I used are the, mr.crd and,mp.crd in the
$AMBERHOME/test/evb/malon_bond_umb_GAFF/ folder.

Thank you

2013/11/9 case <>

> On Sat, Nov 09, 2013, Fabrício Bracht wrote:
> > The problem is
> > that the temperature of the system falls from the previous 298K to around
> > 213 K and it stays like that for as long as I have simulated.
> This looks like a problem with SHAKE inputs: you are shaking the MM bonds
> involving hydrogen, but have set qmshake=0, which turns off shake in the
> quantum region. I'm guessing that the code is computing the number of
> degrees of freedom assuming shake is turned on, and is not correcting that
> number when it processes the qmshake=0 input.
> Try adjusting either qmshake or the noshakemask variables so that you are
> doing the same amount of SHAKEing in both qm and mm calculations.
> (Simplest
> is to use the default value for qmshake, which should be OK unless bonds to
> hydrogen are being made or broken.)
> ...good luck...dac
> _______________________________________________
> AMBER mailing list
AMBER mailing list
Received on Sat Nov 09 2013 - 10:30:03 PST
Custom Search