TclForces and getstep

From: Jim Parker (jim.parker_at_prismsciences.com)
Date: Wed Dec 02 2009 - 00:33:37 CST

Hello,
  I'm trying to generate a time variable force (harmonic perturbation)
for use in NAMD. I thought a straight-forward approach would be to
compute the perturbation using the "getstep" function in Tcl. Then
addforce to an atom-id with size =Amplitude*sin(2*pi*freq*t), 0,0

However, executing the code fails with FATAL ERROR: expected
floating-point number. (actual error below)

The problem seems to be that the TCL interpreter is not substituting a
value for "getstep" into the argument for sin() and thus the argument is
not a floating-point number. Can anyone help point me in the right
direction?

I did see that several people mentioned float-point errors in Tcl. I
think this a syntax error, but maybe related to that??

Cheers,
--Jim Parker
UTSA Physics and Astronomy

Error and code follows:
*****************************************
Reason: FATAL ERROR: expected floating-point number but got
"0.1*sin(2*3.1416*1000*getstep/1000000*2)"
    while executing
"vecscale $xvec $compF"
    (procedure "calcforces" line 10)
    invoked from within
"calcforces"

Charm++ fatal error:
FATAL ERROR: expected floating-point number but got
"0.1*sin(2*3.1416*1000*getstep/1000000*2)"
    while executing
"vecscale $xvec $compF"
    (procedure "calcforces" line 10)
    invoked from within
"calcforces"

Aborted
*****************************************
TCL snipet follows

tclForcesScript {
  # save timestep
  set t [getstep] ; #t is in 2fs steps
  print "The current time is $t"
  # The IDs of the atoms defining the oscillation
  set aidleft 1
  set aidright 2
  addatom $aidleft
  addatom $aidright
  set PI 3.1416

  # real frequency = RFfreq * 10^9 Hz
  set RFfreq 1000 ; #1 THz
  # RFampl in angstroms
  set RFampl 0.1
# print "Freq is $RFfreq and amplitude is $RFampl"

  proc calcforces {} {
    global aidleft aidright PI RFfreq RFampl t
    loadcoords p
    # Calculate the harmonic force component
    # note time = t * time/step
    set compF "$RFampl*sin(2*$PI*$RFfreq*$t/1000000*2)"
    # The force to be applied on each atom along x-axis (basevector)
    set basevec {1. 0. 0.}
    addforce $aidleft [vecscale $basevec $compF]
    addforce $aidright [vecscale $basevec -1.0*$compF]
  }
}

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:51:44 CST