From: Panel Nicolas (M.) (nicolas.panel_at_polytechnique.edu)
Date: Wed Jun 03 2015 - 10:00:53 CDT
Dear namd user,
This is a question about how is calculated the "half-harmonic" potential in the colvars module.
I'm doing alchemical transformations with an hybrid residue (one backbone and two sidechains).
To reduce the fluctuations of the sidechain for low values of lambda I restraint the chi1 angle with
a lower and an upper wall using the colvars module of NAMD 2.10.
To estimate the value of the restraint at each frame I defined the outputAppliedForce in my input file
but the force printed in the output seems to be only the difference between the chi1 and the reference
multiplied by the force constant. So this potential is not harmonic but linear.
After these observations I have two questions:
1) Does the applied force printed in the output reflect the value of the restraint?
2) If not how is calculated this force? Is it with the formula U(x)=k(1-cos(n(x-xref))), with n the dihedral periodicity?
My input file is as following:
#########################################################################
colvarsTrajFrequency 100
colvar {
name chi1P0
outputAppliedForce on
outputSystemForce off
dihedral {
group1 {
psfSegID PROB
atomNameResidueRange N 8-8
}
group2 {
psfSegID PROB
atomNameResidueRange CA 8-8
}
group3 {
psfSegID PROB
atomNameResidueRange CBF 8-8
}
group4 {
psfSegID PROB
atomNameResidueRange CG 8-8
}
}
lowerWall -90.0
LowerWallConstant 0.01
upperWall -50.0
upperWallConstant 0.01
}
#########################################################################
An example of output:
#########################################################################
# step Protcenterdist Distpetprot Distpetprot2 chi1P0 fa_chi1P0
0 5.95219485967896e-01 2.94944522997649e+00 2.84959860158793e+00 -5.13403693303894e+01 0.00000000000000e+00
100 5.66996723173063e-01 2.99280703876506e+00 2.80940631466381e+00 -6.06827217894346e+01 0.00000000000000e+00
200 5.45366266152213e-01 3.02910909719973e+00 2.80181010781602e+00 -6.20402736575662e+01 0.00000000000000e+00
300 6.00054749415614e-01 2.98001984780081e+00 2.81790301610687e+00 -4.52150725810083e+01 -4.78492741899169e-02
400 6.14995241200626e-01 2.96781468532476e+00 2.82604151332475e+00 -5.84744993536359e+01 0.00000000000000e+00
500 6.73079647049610e-01 2.98836913063312e+00 2.71093661129984e+00 -3.97651865050939e+01 -1.02348134949061e-01
600 6.76784683432109e-01 3.11007557140919e+00 2.78390677338792e+00 -5.32106437211886e+01 0.00000000000000e+00
700 6.45534234285689e-01 2.89431851457446e+00 2.78731307706688e+00 -5.06243344344155e+01 0.00000000000000e+00
800 6.01275984753345e-01 2.93146906379729e+00 2.78613058748865e+00 -5.56683687674018e+01 0.00000000000000e+00
#########################################################################
Thanks a lot,
-- Nicolas Panel PhD student École Polytechnique Palaiseau Laboratoire de Biochimie, BIOC Tel: 01.69.33.48.59 ----------
This archive was generated by hypermail 2.1.6 : Thu Dec 31 2015 - 23:21:54 CST