Brief intro to Metadynamics with NAMD

If you want to run metadynamics with NAMD, all you have to do is:

1) include these lines in your MD configuration file:

colvars on
colvarsConfig colvars_dihedral.conf        # colective variable configuration file
#colvarsInput meta_3.restart.colvars.state # COMMENT OUT UNLEES YOU RESTART A JOB

2) Create a colvars configuration file.

This is an example for distance constraint:


colvarsTrajFrequency 1000
colvar {
name r1

lowerBoundary 4.0
lowerWallConstant 50.0
upperBoundary 12.0
upperWallConstant 50.0
width 0.05

outputVelocity on
outputAppliedForce on
distance {
group1 { atomNumbers 1 10 13 15 21 23 25 30 31 32 35 36 }
group2 { atomNumbers 122 123 124 125 126 127 }
}
}
metadynamics {
name metadyn_dist1
colvars r1
hillWeight 0.01 # Default is 0.01 - higher numbers = faster
dumpFreeEnergyFile yes
newHillFrequency 100 # Default is 100 - lower numbers = faster
}


The second example uses two collective variables (dihedrals):


colvarsTrajFrequency 1000
colvar {
name d1
lowerBoundary -180.0
upperBoundary 180.0
width 1.0
outputVelocity on
outputAppliedForce on
dihedral {
group1 { atomNumbers 35 }
group2 { atomNumbers 36 }
group3 { atomNumbers 1 }
group4 { atomNumbers 2 }
}
}
colvar {
name d2
lowerBoundary -180.0
upperBoundary 180.0
width 1.0
outputVelocity on
outputAppliedForce on
dihedral {
group1 { atomNumbers 21 }
group2 { atomNumbers 23 }
group3 { atomNumbers 25 }
group4 { atomNumbers 8 }
}
}
metadynamics {
name metadyn_dih1
colvars d1 d2
hillWeight 0.01 # Default is 0.01 - higher numbers = faster
dumpFreeEnergyFile yes
newHillFrequency 100 # Default is 100 - lower numbers = faster
}

 

Collective variables can be fairly sophisticated. Refer to the NAMD manual for furher details.

The free energy profile is saved in a file with the "pmf" termination. If it's a 1D colective variable, visualization is trivial (e.g. use xmgrace). If you are using 2 collective variables, R is a good tool to visualize the free energy surface as a heat map:

First delete comments and empty lines from the pmf file:

awk '{if (NF==3){print}}' metadyn.pmf > input.txt

Then launch R and execute the following:

r=read.table("input.txt")
x=matrix(unlist(r[,3]),ncol=360,byrow=T)
image((1:nrow(x))-180,(1:ncol(x))-180,x, main="2D METADYNAMICS", xlab="dihedral 1", ylab="dihedral 2")
contour((1:nrow(x))-180,(1:ncol(x))-180,x, add=TRUE)For 2D plots R provides

The result should be something like this: