You can do this nicely and easily using the NAMDenergy plugin of VMD.
Or, you can do it directly with NAMD as follows:
REMARK original generated coordinate pdb file ATOM 1 N ARG 1 -9.063 0.480 -0.011 1.00 1.00 A ATOM 2 HT1 ARG 1 -9.839 -0.024 0.368 1.00 1.00 A ATOM 3 HT2 ARG 1 -9.096 0.458 -1.010 1.00 1.00 A ATOM 4 HT3 ARG 1 -9.082 1.427 0.310 1.00 1.00 A ATOM 5 CA ARG 1 -7.825 -0.148 0.437 1.00 1.00 A ....
structure my.psf
paraTypeCharmm on
parameters par_all27_prot_na.inp
outputname junk
numsteps 1
exclude scaled1-4
temperature 0
COMmotion yes
cutoff 16 # <====
switchdist 14 # <====
dielectric 1.0 # <====
pairInteraction on
pairInteractionGroup1 1
pairInteractionFile my.pdb
pairInteractionSelf on
coordinates my.pdb
set ts 0
coorfile open dcd my.dcd
while { ![coorfile read] } {
firstTimestep $ts
run 0
incr ts 1
}
coorfile close
/usr/local/namd/namd2 script | grep '^ENERGY:' > LOG &awk '{print $12}' LOG | plot to plot TOTALawk '{print $7}' LOG | plot to plot ELECTawk '{print $8}' LOG | plot to plot VDWawk '{print $7+$8}' LOG | plot to plot non-bonded.If you want to calculate interaction energies between two species, you need to modify the B-factor column of my.pdb so that the first species will have the number 1.00, the second the number 2.00. In addition you need to tell NAMD what you have done by adding a couple of lines in the script. On the whole, I think it will be safer (and surely easier) to avoid all this by using VMD's NAMDenergy plugin for this type of calculation.
Take a good look at the dielectric constant used in the NAMD script. The value of 1.0 (corresponding to vacuum) immediately tells you that you should consider the magnitude of the electrostatics contribution to the total energy a very rough (if not very wrong) approximation.






