From: Carlo Guardiani (carlo.guardiani_at_dsf.unica.it)
Date: Fri Feb 17 2017 - 06:53:31 CST

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

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