VMD-L Mailing List

From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Sun Nov 06 2011 - 06:14:50 CST

On Sun, 2011-11-06 at 00:07 -0500, Andrés Morales wrote:
> Dear VMD users:
>
> I am trying to estimate lipid bilayer thickness. I am using a simple
> script to do it, but when I run it in VMD TkConsole the program
> collapses. It is displayed the following message: "A problem caused
> the program to stop working correctly" . I am using a 400 frames dcd

most likely you ran out of memory due sloppy scripting.

> file and the system contains 72 lipids (36 per layer). I used the
> following script:
>
> set outfile [open Espesor.txt
> w];
> puts $outfile "Frame Espesor(P)"
>
> set nf [molinfo top get numframes]
> set sel [atomselect top "lipids"]
> set sel1 [atomselect top all]
> set sel2 [atomselect top "lipids and resid 1 to 36 and name P"]
> set sel3 [atomselect top "lipids and resid 37 to 72 and name P"]
>
> for {set i 0 } {$i < $nf } { incr i } {
> $sel frame $i
> $sel1 frame $i
> $sel2 frame $i
> $sel3 frame $i
>
> $sel1 moveby [vecinvert [measure center $sel weight mass]]
> set sumz 0
> foreach atom [$sel2 get index] {
> set pos [lindex [[atomselect top "index $atom"] get {x y z}] 0]

ouch! this is very ugly, inefficient and causing memory leaks.

> set z1 [lindex $pos 0]
> set sumz [expr abs($sumz)+ abs($z1)]
> }
> set promz1 [expr $sumz / 36]
>
> set sumz2 0
> foreach atom [$sel3 get index] {
> set pos2 [lindex [[atomselect top "index $atom"] get {x y z}] 0]

and another one.

> set z2 [lindex $pos2 0]
> set sumz2 [expr abs($sumz2)+ abs($z2)]
> }
> set promz2 [expr $sumz2 / 36]
>
> set thinckness1 [expr abs($promz1) + abs($promz2)]
>
> puts $outfile "$i $thinckness1"
> }
> close $outfile
>
> Any ideas to solve it?

this is pretty trivial:
you already have the selections, so you don't
have to create a new selection to get access to
the atom coordinates. just replace the for loops
with:

foreach pos [$sel2 get {x y z}]

and

foreach pos2 [$sel3 get {x y z}]

but since you only use the z-value, you can
just "get" the z components and have less useless
data to move around.

axel.

> Thanks a lot
>
> Andres