VMD-L Mailing List

From: Andr�s Morales (h.andres.m1986_at_gmail.com)
Date: Sun Nov 06 2011 - 00:07:36 CDT

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 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]
  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]
  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?
Thanks a lot

Andres