On 23 May 2014 12:18, Massimiliano Porrini <m.porrini.iecb.u-bordeaux.fr>wrote:
> Hi Jason,
>
> Very useful appreciated hints!
>
> Below my few comments.
>
>
> On 22 May 2014 23:42, Jason Swails <jason.swails.gmail.com> wrote:
>
>>
>> On May 22, 2014, at 3:52 PM, Massimiliano Porrini <
>> m.porrini.iecb.u-bordeaux.fr> wrote:
>>
>> > Hi all,
>> >
>> > I am attempting to run in vacuo simulations and I would like to know
>> which
>> > would be
>> > the best setting with sander module.
>> >
>> > As stated at p.60 of Amber14 user guide, igb=6 would correspond to an in
>> > vacuo
>> > simulation and I quote here what it is written "This option is logically
>> > equivalent to
>> > setting igb=0 and eedmeth=4, although the implementation (and
>> computational
>> > efficiency)
>> > is not the same".
>> >
>> > It was not clear to me in which of the 2 cases (igb=6 or ntb=0, igb=0,
>> > eedmeth=4)
>> > the computational efficiency should improve, therefore I made some tests
>> > with the
>> > 3 different temperature couplings available in Amber12, here are the
>> > results:
>> >
>> > *igb=6*
>> > Berendsen: 113.7 ns/day
>> > Andersen(vrand=1000): 70.12 ns/day
>> > Langevin(gamma_ln=2): 61.79 ns/day
>> >
>> > *ntb=0, igb=0, eedmeth=4*
>> > Berendsen: 66.48 ns/day
>> > Andersen(vrand=1000): 68.94 ns/day
>> > Langevin(gamma_ln=2): 72.70 ns/day
>> >
>> > I ran a simulation of a 758 atoms molecule, on a blade with 8 CPUs dual
>> core
>> > (hence 16 cores) and hyper-threading (hence 32 threads), however I
>> > utilised only 16 cores of such a blade and with the following
>> > sander input file (example from the simulation with Berendsen and
>> igb=6):
>>
>> First a few words of caution about running vacuum dynamics -- the
>> parameters (including the standard charge sets and am1-bcc charge
>> derivation schemes) are designed to reproduce solvent-induced polarization
>> (particularly polar solvents). So you should be careful about how you
>> approach vacuum dynamics in general (for instance, dominant protonation
>> states are completely different -- zwitterions are unstable in vacuum but
>> very stable in solution).
>>
>
> Yes, I am aware of this risk, as removing the solvent (the very polar
> water) the potential energy surface is totally
> different and in principle one should recalculate the whole set of
> parameters for the solute.
> However, I just want the solute to explore the accessible microstates in
> absence of water, and compare
> the associated collision cross sections with ion mobility mass
> spectrometry experimental values.
>
>
>
>> >
>> > dyn in vacuo (igb=6)
>> > &cntrl
>> > imin = 0, irest = 1, ntx = 5,
>> > nstlim = 100000, dt = 0.001,
>> > ntc = 2, ntf = 2,
>> > cut = 1000.0,
>> > igb = 6,
>> > ntpr = 5000, ntwr = 10000, ntwx = 2000,
>> > ntt = 1,
>> > temp0 = 300.0,
>> > ioutfm = 1
>> > /
>> >
>> > Therefore I would obviously be tempted to utilise Berendsen thermostat
>> with
>> > igb=6,
>> > as it is the fastest one (113.7 ns/day)
>>
>> I question the use of any thermostat at all. The idea behind a
>> thermostat is that there is some "unlimited heat bath at a constant
>> temperature" surrounding the system of interest that is able to provide or
>> absorb heat energy as needed. This makes far more sense in a condensed
>> phase simulation than it does in vacuum or gas. If I was running a gas
>> phase simulation, for whatever reason, I would probably use a short time
>> step, stringent SHAKE tolerance (if I was using SHAKE at all), and if I
>> didn't observe acceptable energy conservation, a Berendsen thermostat with
>> a very weak coupling constant (on the order of 10 to 20 /ps) to smooth over
>> integration errors.
>>
>
> I totally forgot to mention 2 things:
> 1) the simulation in microcanonical ensemble was in my "to-do" list, I
> have just sent this email prematurely and I agree that in vacuo
> (no external bath) one should implement pure Newton dynamics, however (and
> here it comes the 2nd thing)
> 2) the simulation in not actually in the "high" vacuum, but there is a
> buffer gas at (relatively) low pressure ranges
> (with respect to the atmospheric one) and temperature T (which, in my
> case, is set equal to 300K), hence I would
> pursue the Berendsen coupling choice, with the settings suggested by you
> (short time-step, very low SHAKE tolerance
> and weak coupling).
>
Sorry! One errata:
In the 2nd point I meant to say that ion mobility experiments are conducted
with a buffer gas at low P and
temperature T.
Whereas I will carry out my simulations in vacuo with Berendsen weak T
coupling.
>
>
>> > But, provided that Berendsen coupling is kind of deprecated for several
>> > reasons (Amber14
>> > user guide, p.290), especially for simulations in implicit solvent, in
>> > which a uniform
>> > temperature is not guaranteed because collisions with water molecules
>> are
>> > missing, do you think that it would be a sensible choice to use
>> Berendsen
>> > coupling
>> > for in vacuo simulations?
>> >
>> > My answer would be yes for two reasons:
>> > 1) performance efficiency
>> > 2) with in vacuo simulations one does not actually want
>> > any "random kick" (like it occurs with Andersen or Langevin
>> > coupling) that would emulate the collisions with the solvent.
>>
>> The weak coupling algorithm (for either temperature or pressure) has
>> well-documented shortcomings and is fundamentally _known_ not to sample
>> from the desired distribution. For one thing, fluctuations are too small,
>> and the thermostat doesn't even keep correct partitioning between
>> translational, rotational, and vibrational motion (flying block of ice
>> effect in explicit solvent). I would personally not use the Berendsen
>> thermostat here except as I described above.
>>
>
> I was aware of these Berendsen thermostat drawbacks too, this is one main
> reason why I wanted any suggestion about
> running dynamics in gas phase (not actually vacuum).
>
>
>
>>
>> >
>> > Am I right/wrong?
>> > Any comments/hints on this issue would be very appreciated.
>> >
>> > Many thanks!
>> >
>> > PS Just to be meticulous (running 1 step of minimisation) I have also
>> > checked if the potential functions utilised were
>> > actually the same between igb=6 and igb=0, ntb=0, eedmeth=4 settings and
>> > they
>> > perfectly match :) as you can see below (therefore igb=6 with Berendsen
>> is
>> > actually
>> > the most efficient implementation):
>>
>> The difference is, as the manual states, in the implementation. I would
>> use igb=6 when doing gas-phase energy calculations since the electrostatic
>> code used to compute the energy in that case respects the exclusion list
>> stored in the topology file whereas the PME code (eedmeth=4, igb=0)
>> computes the exclusions from the bond list. This only really matters if
>> you add custom exclusions with ParmEd.
>>
>
> I am slightly puzzled here, in that I can not use PME where there are no
> periodic boundaries
> (like for molecules in gas phase).
> In fact, if igb>0 then ntb=0 by default (thus sander module does not
> implement PME algorithm).
> Therefore, if one wants to utilise exactly the same setting, one needs to
> impose ntb=0 when [eedmeth=4, igb=0],
> otherwise ntb would be equal to 1 by default (igb=0).
> As a matter of fact, my simulation with [eedmeth=4, igb=0] terminated
> because pmemd complained about the missing
> box.
>
>
>
>> Anyway, that's my take on this.
>>
>> HTH,
>> Jason
>>
>
> Thanks a lot again, TDH (this definitely helps)!
>
>
>>
>> --
>> Jason M. Swails
>> BioMaPS,
>> Rutgers University
>> Postdoctoral Researcher
>>
>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>
>
>
> --
> Dr Massimiliano Porrini
> Valérie Gabelica Team
> U869 ARNA - Inserm / Bordeaux University
> Institut Européen de Chimie et Biologie (IECB)
> 2, rue Robert Escarpit
> 33607 Pessac Cedex
> FRANCE
>
> Tel : 33 (0)5 40 00 63 31
>
> http://www.iecb.u-bordeaux.fr/teams/GABELICA
> Emails: m.porrini.iecb.u-bordeaux.fr
> mozz76.gmail.com
>
--
Dr Massimiliano Porrini
Valérie Gabelica Team
U869 ARNA - Inserm / Bordeaux University
Institut Européen de Chimie et Biologie (IECB)
2, rue Robert Escarpit
33607 Pessac Cedex
FRANCE
Tel : 33 (0)5 40 00 63 31
http://www.iecb.u-bordeaux.fr/teams/GABELICA
Emails: m.porrini.iecb.u-bordeaux.fr
mozz76.gmail.com
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Fri May 23 2014 - 04:00:02 PDT
