From: luca (bellucci14_at_unisi.it)
Date: Thu Nov 23 2006 - 02:22:23 CST
Alle 23:50, mercoled¨¬ 22 novembre 2006, hai scritto:
> Hi All,
>
> I need to add a harmonic restraint on two atoms in my protein molecule
> during MD simulation. It seems to me that I should use TclForces. As I
> have never used this feature before, I'd like to know where I can find
> example scripts to follow. Thanks!
>
> Yao Fan
Hi, Yao
There is an example for two atoms :
-------------------------------------------------------
# Scripting
# switch on tclforces
tclforces on
# Call tcl force
tclforcesScript {
# Set atom index
set at1 1300
set at2 4618
addatom $at1
addatom $at2
# Set costant force and eq. dist. V=-0.5*K(r-ro)^2
set k 100.0
set r0 5.907
print " TCLFORCE ACTION "
# Call procedure
proc calcforces {} {
# set change variable
global at1 at2 k r0
# load coordinates
loadcoords c
# set distance for atom 1 and atom 2
set r12 [vecsub $c($at2) $c($at1) ]
set r [veclength $r12 ]
# set unit vector
set n0 [vecnorm $r12]
# Force: F=-dV/dr
set force [expr $k*($r-$r0)]
# set vector force for atom 1 and atom 2
set f1 [vecscale $force $n0 ]
set f2 [vecscale -$force $n0 ]
# add force f1 and f2 to atoms
addforce $at1 $f1
addforce $at2 $f2
}
# Ther is a definition of subroutine veclenght and vecnorm
# Returns: the vector length
proc veclength {v} {
set retval 0
foreach term $v {
set retval [expr $retval + $term * $term ]
}
return [expr sqrt($retval)]
}
# Returns: the normal vector
proc vecnorm {v} {
set sum 0
foreach term $v {
set sum [expr $sum + $term*$term]
}
set sum [expr sqrt($sum)]
set retval {}
foreach term $v {
lappend retval [expr $term / $sum]
}
return $retval
}
}
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:44:12 CST