AMBER Archive (2009)

Subject: Re: [AMBER] Maxwell-Boltzmann distribution and ig=-1

From: Carlos Simmerling (carlos.simmerling_at_gmail.com)
Date: Wed Apr 01 2009 - 06:01:51 CDT


if you assign temperatures at the start and do not use a thermostat you will
almost certainly not get the temperature that you want, mostly due to shake.
you should look at the runs more closely. did those articles advocate
constant energy runs? even if they did, short equilibration with a
thermostat is probably a good idea before your production run. otherwise it
is difficult to know how to compare the different runs in terms of their
effective temperature.

On Tue, Mar 31, 2009 at 10:43 PM, Naser Alijabbari <na3m_at_virginia.edu>wrote:

> Last week I came across the following papers:
> M. A. Balsera, W. Wriggers, Y. Oono, and K. Schulten, “Principal Component
> Analysis and Long Time Protein Dynamics”, J. Phys. Chem. 1996, 100,
> 2567-2572.
> L. S. Caves, J. D. Evanseck, and M. Karplus, “Locally accessible
> conformations of proteins: Multiple molecular dynamics simulations of
> crambin”,* *Protein Science* *(1998), 7649-666.
>
> In the last paper (which is old) they advocate using multiple short
> simulations with new velocity assigments. I think the only way to do this
> in
> AMBER is by setting ntx=0, irest=1, and ig=-1 otherwise the same velocities
> will be assigned for all simulations. Some results have interesting
> features
> and I was wondering what seeds were used, and hence what initial velocities
> were being assigned.
>
> On Tue, Mar 31, 2009 at 10:24 PM, Carlos Simmerling <
> carlos.simmerling_at_gmail.com> wrote:
>
> > I don't think any of this is relevant since you have ntt=0 and therefore
> no
> > thermostat. is that what you really want?
> >
> > On Tue, Mar 31, 2009 at 10:15 PM, Naser Alijabbari <na3m_at_virginia.edu
> > >wrote:
> >
> > > I am using Amber 10 and these are the parameters for a test run in
> > sander.
> > >
> > > &cntrl
> > > imin = 0, irest = 0, ntx = 1,
> > > ntb = 1,
> > > cut = 10, ntr = 0,
> > > ntc = 2, ntf = 2, tol=.000001,
> > > tempi = 300, temp0 = 300.0, ig=-1,
> > > ntt = 0, nscm=0,
> > > nstlim = 95000, dt = 0.001,
> > > ntpr = 100, ntwx = 50, ntwr = 1000,
> > > &end
> > > &ewald
> > >
> > > ew_type=0,
> > > dsum_tol=0.00001,
> > > &end
> > >
> > > END
> > >
> > >
> > > I am not sure where ig comes into setvel. In mdread there is only 'if(
> > > ig==-1 ) call microsec(ig)'
> > >
> > > ! Assign velocities from a Maxwellian distribution.
> > >
> > > subroutine setvel(nr,v,winv,tempi,init,iscale,scalm)
> > >
> > > use amoeba_mdin, only : iamoeba,beeman_integrator
> > >
> > > implicit none
> > >
> > >
> > >
> > > double precision v(*),winv(*),tempi,scalm
> > >
> > > integer nr,init,iscale
> > >
> > >
> > >
> > > integer i,j,nr3, icopy
> > >
> > > double precision boltz,sd
> > >
> > >
> > >
> > > nr3 = 3*nr
> > >
> > >
> > >
> > > ! ----- Assign velocities from a Maxwellian distribution
> > >
> > >
> > >
> > > if (tempi < 1.d-6) then
> > >
> > > v(1:nr3+iscale) = 0.d0
> > >
> > > return
> > >
> > > end if
> > >
> > >
> > >
> > > boltz = 8.31441d-3*tempi/4.184d0
> > >
> > > i = 0
> > >
> > > do j=1,nr
> > >
> > > sd = sqrt(boltz*winv(j))
> > >
> > > call gauss(0.d0,sd,v(i+1))
> > >
> > > call gauss(0.d0,sd,v(i+2))
> > >
> > > call gauss(0.d0,sd,v(i+3))
> > >
> > > i = i+3
> > >
> > > end do
> > >
> > > if (iscale > 0) then
> > >
> > > sd = sqrt(boltz/scalm)
> > >
> > > do j=1,iscale
> > >
> > > call gauss(0.d0,sd,v(i+j))
> > >
> > > end do
> > >
> > > end if
> > >
> > > if ( iamoeba == 1 .and. beeman_integrator == 1 )then
> > >
> > > do j=1,nr3
> > >
> > > v(j) = 20.455d0*v(j) ! beeman uses time in ps units
> > >
> > > enddo
> > >
> > > endif
> > >
> > > return
> > >
> > > end subroutine setvel
> > >
> > > On Tue, Mar 31, 2009 at 9:41 PM, Carlos Simmerling <
> > > carlos.simmerling_at_gmail.com> wrote:
> > >
> > > > in dynlib.f:
> > > >
> > > > ! Assign velocities from a Maxwellian distribution.
> > > > subroutine setvel(nr,v,winv,tempi,init,iscale,scalm)
> > > >
> > > > look in mdread.f- ig should be written to the mdout file for ntt=2 or
> > 3.
> > > >
> > > > this assumes amber10, so you'll need to be more specific if this
> > doesn't
> > > > solve your problem.
> > > >
> > > >
> > > > On Tue, Mar 31, 2009 at 9:35 PM, Naser Alijabbari <na3m_at_virginia.edu
> >
> > > > wrote:
> > > >
> > > > > The default value for ig=71277, but is there a way to find what was
> > > used
> > > > as
> > > > > a seed when ig=-1?
> > > > > Also where is the routine for assigning velocities
> > > > > according Maxwell-Boltzmann distribution? fgrep of the sander
> source
> > > code
> > > > > only returns poisson_boltzmann routine.
> > > > > Naser
> > > > > _______________________________________________
> > > > > AMBER mailing list
> > > > > AMBER_at_ambermd.org
> > > > > http://lists.ambermd.org/mailman/listinfo/amber
> > > > >
> > > > _______________________________________________
> > > > AMBER mailing list
> > > > AMBER_at_ambermd.org
> > > > http://lists.ambermd.org/mailman/listinfo/amber
> > > >
> > > _______________________________________________
> > > AMBER mailing list
> > > AMBER_at_ambermd.org
> > > http://lists.ambermd.org/mailman/listinfo/amber
> > >
> > _______________________________________________
> > AMBER mailing list
> > AMBER_at_ambermd.org
> > http://lists.ambermd.org/mailman/listinfo/amber
> >
> _______________________________________________
> AMBER mailing list
> AMBER_at_ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
_______________________________________________
AMBER mailing list
AMBER_at_ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber