**From:** Giacomo Fiorin (*giacomo.fiorin_at_gmail.com*)

**Date:** Fri Feb 17 2017 - 12:22:05 CST

**Next message:**Life Sciences Inc: "Re: Npt lipid simulation with electrostatics 2"**Previous message:**Life Sciences Inc: "Re: Npt lipid simulation with electrostatics 2"**In reply to:**Carlo Guardiani: "Metadynamics state file"**Next in thread:**Carlo Guardiani: "Re: Metadynamics state file"**Reply:**Carlo Guardiani: "Re: Metadynamics state file"**Messages sorted by:**[ date ] [ thread ] [ subject ] [ author ] [ attachment ]

Hi Carlo, you can refer to the Colvars reference doc here:

http://colvars.github.io/colvars-refman-namd/colvars-refman-namd.html

All derivatives (or forces, if you change the sign) of the Gaussian

functions are computed analytically at each grid point, and those from

different Gaussians are summed together. No numerical derivatives are used

in any case. When the total force of the Gaussian bias is needed at a

point that is not exactly the grid point, the closest grid point is taken.

Not knowing your configuration nor what do you mean by "significantly

higher" I would say that you may have a large grid spacing compared to the

Gaussian width. The default option is set so that the integral of one

Gaussian equals the volume of one grid bin, which is not a fine

interpolation. It works for many small Gaussians:

http://dx.doi.org/10.1080/00268976.2013.813594

but will not be accurate for few, large Gaussians.

I suggest you try changing the width and let the Colvars code recompute the

sum of all Gaussians on the grid, and compare them to your own

interpolation.

Giacomo

On Fri, Feb 17, 2017 at 7:53 AM, Carlo Guardiani <

carlo.guardiani_at_dsf.unica.it> wrote:

*> Dear NAMD experts,
*

*>
*

*> I a novice in running metadynamics simulations
*

*> with NAMD and I would like to ask some clarifications
*

*> on the content of the state file. My well-tempered
*

*> metadynamics simulation is biasing two collective
*

*> variables.
*

*>
*

*> 1) As I understand, the first list of numbers
*

*> represents the biasing potential. I checked
*

*> that the number of values in this list
*

*> corresponds to the number of grid ponts in the
*

*> pmf files. Could you please confirm whether these
*

*> numbers represent the values of the biasing
*

*> potential (sum of gaussians) in correspondence
*

*> of the grid points ?
*

*>
*

*> 2) If this is the case, in which order are these
*

*> grid point values listed ? Is it the same order
*

*> employed in the pmf file, i.e.
*

*>
*

*> [CV1(1) CV2(1)] [CV1(1) CV2(2)] [CV1(1) CV2(3)]
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*> [CV1(1) CV2(n2-2)] [CV1(1) CV2(n2-1)] [CV1(1) CV2(n2)]
*

*>
*

*>
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*>
*

*> [CV1(n1) CV2(1)] [CV1(n1) CV2(2)] [CV1(n1) CV2(3)]
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*> [CV1(n1) CV2(n2-2)] [CV1(n1) CV2(n2-1)] [CV1(n1) CV2(n2)]
*

*>
*

*>
*

*> 3) What is the exact meaning of the second list of numbers ?
*

*> Do they represent the values of the gradients with respect
*

*> to the two CV in corresponcence of the grid points ? In which
*

*> order are they listed ? Is it the following (calling the two
*

*> CV x1 and x2) ?
*

*>
*

*>
*

*> dV[x1(1) x2(1)]/dx1 dV[x1(1) x2(1)]/dx2 dV[x1(1) x2(2)]/dx1
*

*> dV[x1(1) x2(2)]/dx2 dV[x1(1) x2(3)]/dx1 dV[x1(1) x2(3)]/dx2
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*> dV[x1(1) x2(n2-2)]/dx1 dV[x1(1) x2(n2-2)]/dx2 dV[x1(1) x2(n2-1)]/dx1
*

*> dV[x1(1) x2(n2-1)]/dx2 dV[x1(1) x2(n2)]/dx1 dV[x1(1) x2(n2)]/dx2
*

*>
*

*>
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*> dV[x1(n1) x2(1)]/dx1 dV[x1(n1) x2(1)]/dx2 dV[x1(n1) x2(2)]/dx1
*

*> dV[x1(n1) x2(2)]/dx2 dV[x1(n1) x2(3)]/dx1 dV[x1(n1) x2(3)]/dx2
*

*> . . .
*

*> . . .
*

*> . . .
*

*>
*

*> dV[x1(n1) x2(n2-2)]/dx1 dV[x1(n1) x2(n2-2)]/dx2 dV[x1(n1) x2(n2-1)]/dx1
*

*> dV[x1(n1) x2(n2-1)]/dx2 dV[x1(n1) x2(n2)]/dx1 dV[x1(n1) x2(n2)]/dx2
*

*>
*

*>
*

*> 4) How does NAMD compute the gradient ? Does it use analytical or numerical
*

*> derivatives ? And in this latter case, how does it deal with the points
*

*> at
*

*> the boundary of the grid ?
*

*>
*

*>
*

*> 5) Since I am using the keepHills option, I tried to recover the biasing
*

*> potential by summing the guassians whose parameters are listed at the
*

*> end
*

*> of the state file. The PMF that I retrieved is qualitatively similar to
*

*> the one computed by NAMD and shown on the pmf file but my free energy
*

*> values
*

*> are significantly higher (it's not just the matter of an additive
*

*> constant,
*

*> the distance between the deepest minimum and the highest maximum is
*

*> larger).
*

*> The core of the fortran90 program that I wrote is the following:
*

*>
*

*> fact=(Temp + DeltaT)/DeltaT
*

*> open(10,file=inpfile,status='old',form='formatted')
*

*> 11 read (10,*,end=12) istep, weight, x1c, x2c, sig1, sig2
*

*>
*

*>
*

*> do i=1,n1
*

*> do j=1, n2
*

*> x1=x1min+(i-1)*dx1
*

*> x2=x2min+(j-1)*dx2
*

*> exp1=(x1-x1c)**2
*

*> exp1=exp1/(2.d0*sig1**2)
*

*> exp2=(x2-x2c)**2
*

*> exp2=exp2/(2.d0*sig2**2)
*

*> exptot=exp1+exp2
*

*> hills(i,j)=hills(i,j)+ weight*dexp(-exptot)
*

*> end do
*

*> end do
*

*>
*

*> goto 11
*

*> 12 close(10)
*

*>
*

*> open(20,file=outfile1,status='unknown',form='formatted')
*

*> do i=1,n1
*

*> do j=1, n2
*

*> x1=x1min+(i-1)*dx1
*

*> x2=x2min+(j-1)*dx2
*

*> pmf(i,j)=-hills(i,j)*fact
*

*> write(20,*) x1, x2, hills(i,j)
*

*> end do
*

*> write(20,*) " "
*

*> end do
*

*> close(20)
*

*>
*

*>
*

*> Could you please tell me what I did wrong ? Does NAMD employ
*

*> any special trick to compute the biasing potential ?
*

*>
*

*> Thank you very much for your help and kind regards,
*

*>
*

*> Carlo Guardiani
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

-- Giacomo Fiorin Associate Professor of Research, Temple University, Philadelphia, PA Contractor, National Institutes of Health, Bethesda, MD http://goo.gl/Q3TBQU https://github.com/giacomofiorin

**Next message:**Life Sciences Inc: "Re: Npt lipid simulation with electrostatics 2"**Previous message:**Life Sciences Inc: "Re: Npt lipid simulation with electrostatics 2"**In reply to:**Carlo Guardiani: "Metadynamics state file"**Next in thread:**Carlo Guardiani: "Re: Metadynamics state file"**Reply:**Carlo Guardiani: "Re: Metadynamics state file"**Messages sorted by:**[ date ] [ thread ] [ subject ] [ author ] [ attachment ]

*
This archive was generated by hypermail 2.1.6
: Sun Dec 31 2017 - 23:21:04 CST
*