Re: tcl force

From: Seungho Choe (sec53_at_pitt.edu)
Date: Sat Oct 17 2009 - 13:08:11 CDT

Hi JC and NAMD users,

I'm trying to compare both methods: tcl force & SMD.
Could you please look at my tcl script again and tell me whether it is
correct ?
I wonder whether I'm missing some factors in the script. Thanks a lot.

Best,
Seungho

JC Gumbart wrote:
> Yes, you can. But it will only work in 1D as you've proposed, e.g., not n =
> (1,1,1).
>
> Regarding your initial question, see the addenergy term in the manual:
> http://www.ks.uiuc.edu/Research/namd/2.7b1/ug/node43.html
>
>
> -----Original Message-----
> From: owner-namd-l_at_ks.uiuc.edu [mailto:owner-namd-l_at_ks.uiuc.edu] On Behalf
> Of Seungho Choe
> Sent: Saturday, October 17, 2009 12:21 AM
> To: Seungho Choe
> Cc: namd-l_at_ks.uiuc.edu
> Subject: Re: namd-l: tcl force
>
> Hi,
>
> I wonder whether I can use SMD with v=0, n= (0,0,1) to restrain the
> center of mass of a molecule in z direction ?
> It seems that in this case I don't need to worry about the energy term
> that I mentioned in the previous my post. Could someone clarify this ?
> Apology if this issue was already discussed before.
>
> Best,
> Seungho
>
>
>
>
> Seungho Choe wrote:
>
>> Dear all,
>>
>> I'm trying to restrain the center of mass (only z position) of a
>> molecule. Below is my tcl script. What I'm wondering is whether a new
>> energy term from the added force is applied correctly to the total
>> energy. Do I need to explicitly add the energy term in the script ? Is
>> there any way to check this ? Any comments are welcome. Thank you in
>> advance.
>>
>> Best,
>> Seungho
>>
>> ************************************************************************
>> set numatoms 63230
>>
>> set atoms {}
>> for {set i 63207} { $i <= $numatoms} { incr i} {
>> lappend atoms $i
>> }
>>
>> print "atoms:$atoms"
>>
>> foreach atom $atoms {
>> addatom $atom
>> }
>>
>> set groupid [addgroup $atoms]
>>
>>
>> # Set spring constant
>> set k 7.0
>> set z0 -3
>>
>> proc calcforces {} {
>>
>> global atoms groupid numatoms k z0
>> loadcoords p
>>
>> #print "$p($groupid)"
>>
>> set z1 [lindex $p($groupid) 2]
>>
>> print "z comp:$z1"
>>
>> set r01 [expr $z1-$z0]
>> #print "r01:$r01"
>>
>> if {$r01 > 0} {
>> set n01 {0 0 -1}
>> } elseif {$r01 < 0} {
>> set n01 {0 0 +1}
>> }
>>
>> # Calculate the "force"
>> set force [expr abs($k*$r01)]
>>
>> set forcez [vecscale $force $n01]
>>
>> print "forcez:$forcez"
>>
>> addforce $groupid $forcez
>>
>> }
>>
>>
>>
>
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:53:23 CST