Re: Distance Restraint

From: Jérôme Hénin (jhenin_at_ifr88.cnrs-mrs.fr)
Date: Tue Nov 23 2010 - 11:27:03 CST

Dear Jacopo,

The best way to apply restraints is through the "collective variables"
feature (colvars for short). See the doc there:
http://www.ks.uiuc.edu/Research/namd/2.7/ug/node47.html

For a distance restraint, you should include the following in the NAMD
config file:

colvars on
colvarsConfig distance_restraint.in

and provide a file distance_restraint.in with contents like this:

colvarsTrajFrequency 100
colvarsRestartFrequency 10000

colvar {
 name d

 distance {
   group1 {
      atomNumbers 25371
   }
    group2 {
     atomNumbers 8186
    }
 }
}

harmonic {
 colvars d
  centers 6.0
 forceConstant 50.0
}

In addition, this will output a trajectory of the distance value to the file
<outputName>.colvars.traj.

Best,
Jerome

On 23 November 2010 17:23, Jacopo Sgrignani <sgrigna_at_sissa.it> wrote:
> Dear Namd users
> i'm new to NAMD and i would like to use an harmonic restraint between to
> atoms.
> I found a script in this list:
>
> # Scripting
> # switch on tclforces
> tclforces on
> # Call tcl force
>
> tclforcesScript {
>
> # Set atom index
> set at1 25371
> set at2 8186
> addatom $at1
> addatom $at2
>
> # Set costant force and eq. dist. V=-0.5*K(r-ro)^2
> set k 050.0
> set r0 6.000
> 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
> }
> }
>
> But now i get this error:
>
> TCL: Running for 700 steps
> TCL: expected floating-point number but got "--17.996288924"
> FATAL ERROR: expected floating-point number but got "--17.996288924"
> while executing
> "vecscale -$force $n0 "
> (procedure "calcforces" line 16)
> invoked from within
> "calcforces"
> ------------- Processor 0 Exiting: Called CmiAbort ------------
> Reason: FATAL ERROR: expected floating-point number but got
"--17.996288924"
> while executing
> "vecscale -$force $n0 "
> (procedure "calcforces" line 16)
> invoked from within
> "calcforces"
>
>
> Could anybody help me?
>
> Thanks
>
> Jacopo
>
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 05:23:24 CST