From: Vlad Cojocaru (Vlad.Cojocaru_at_eml-r.villa-bosch.de)
Date: Mon Oct 29 2007 - 07:24:02 CDT
Dear NAMD users,
I would like to use Tcl forces to move an atom with constant velocity
along a directon defined by two other atoms in the protein (an internal
direction which would change between time steps in function of how the
protein translates and/or rotates). Below I attach an attempt of script
to do that. However, I am not sure whether this is a correct script and
NAMD complains about the "vecnorm" command.
First I would like to ask whether somebody has performed such a task in
NAMD and if the script below seems reasonablefor this? If yes, could
somebody indicate me how to calculate the direction in each time step
without using "vecnorm" (using only routines implemented in namd)?
Thanks in advance
Best wishes
vlad
----tcl script ---------
set id1 7614
set id2 1193
set id3 3101
set atoms {}
lappend atoms $id1 $id2 $id3
foreach atom $atoms {
addatom $atom
}
# set the SMD output frequency, initialize the time counter,
set Tclfreq 50
set t 0
#set force constant in kcal/mol AA^2 and pulling vel. in AA/ts;
set k 7.2
#set velocity to be 35 A/100 ps
set v 0.00035
# initial coordinates of pulled atom (t = 0)
set x0 52.842
set y0 50.161
set z0 51.822
set outfilename smd_tcl.out
open $outfilename w
proc calcforces {} {
global atoms id1 id2 id3 Tclfreq t k v x0 y0 z0 outfilename
# form local array coords
loadcoords coords
# calculate direction of pulling
set n [vecnorm [expr 0.5*[vecsub $coords($id2) $coords($id3)]]]
# get coords of pulled atom
set r $coords($id1)
# get the new x,y,z components and calculate the force
set x1 [lindex $r 0]
set y1 [lindex $r 1]
set z1 [lindex $r 2]
set nx [lindex $n 0]
set ny [lindex $n 1]
set nz [lindex $n 2]
set vx [expr $v*$nx]
set vy [expr $v*$ny]
set vz [expr $v*$nz]
# calculate forces
set fx [expr $k*($x0+$vx*$t-$x1)]
set fy [expr $k*($y0+$vy*$t-$y1)]
set fz [expr $k*($z0+$vz*$t-$z1)]
lappend f $fx $fy $fz
# add forces
addforce $id1 $f
# set the output
set foo [expr $t % $Tclfreq]
if { $foo == 0 } {
set outfile [open $outfilename a]
set time [expr $t/1000.0]
puts $outfile "$time $r $f"
close $outfile
}
incr t
return
}
-- ---------------------------------------------------------------------------- Dr. Vlad Cojocaru EML Research gGmbH Schloss-Wolfsbrunnenweg 33 69118 Heidelberg Tel: ++49-6221-533266 Fax: ++49-6221-533298 e-mail:Vlad.Cojocaru[at]eml-r.villa-bosch.de http://projects.villa-bosch.de/mcm/people/cojocaru/ ---------------------------------------------------------------------------- EML Research gGmbH Amtgericht Mannheim / HRB 337446 Managing Partner: Dr. h.c. Klaus Tschira Scientific and Managing Director: Prof. Dr.-Ing. Andreas Reuter http://www.eml-r.org ----------------------------------------------------------------------------
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:45:26 CST