From: Jerome Henin (jhenin_at_cmm.chem.upenn.edu)
Date: Mon Oct 29 2007 - 11:20:29 CDT
Dear Vlad,
The ABF code needs the same kind of vector routines that you are
using. Our solution was to steal the file vectors.tcl from VMD, which
implements some useful operations. It lives in the lib/abf directory
of NAMD, you can source it into your script.
As a side note, the "abscissa" variable for ABF has been modified
since NAMD 2.6 and now works in a way similar to what you want to do,
with two groups of atoms defining the direction. The difference is
that ABF won't pull at constant velocity. In case you are interested
anyway, the new file should be publicly available very soon.
Best,
Jerome
On 10/29/07, Vlad Cojocaru <Vlad.Cojocaru_at_eml-r.villa-bosch.de> wrote:
> 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