VMD-L Mailing List

From: Peter Freddolino (petefred_at_ks.uiuc.edu)
Date: Mon Apr 23 2007 - 10:10:20 CDT

Hi Jianhui,
while we'd need a few more details (such as where it dies) to figure out
exactly what is going on, I can recommend a couple of changes that
should make this easier and more memory efficient:

-Use lsort -unique instead of bothering to define a new procedure for
removing duplicate elements from a list

-Write your main analysis as nested loops instead of first creating all
the selections you could possibly need and iterating over them, eg:

for {set i 1} {$i <= 1536} {incr i} {
    set thiswatsel [atomselect top "water and within 5.5 of resid $i"]
   
    for {set j 0} {$j < $numfrm} {incr j} {
       $thiswatsel frame $j
       $thiswatsel update
       set closeresids [lsort -unique -integer [$thiswatsel get resid]]
       .....

    }

     $thiswatsel delete
}

You can use a list or array to store the data and unpack it at the end
as desired.

Peter

Jianhui Tian wrote:
> Hi all,
>
> I have a water box of 1536 molecules. Now what I want to do is to make
> atomselect of each of the molecule and then get the residue list of
> other molecules within some distance of each water molecule. The
> problem is that the tcl code will get killed. I dont know the reason.
> Is it due to the memory limitation. I only have 512 memory in my PC
> now. How can I do what I want to do here? I also posted the code I
> write about this part.
>
>
> #########################################
> #########################################
> proc simplifylist {slist} {
> set newlist { }
> foreach elem $slist {
> if {[lsearch $newlist $elem] == -1} {
> lappend newlist $elem
> }
> }
> set newlist [lsort $newlist]
> return $newlist
> }
>
> ########################################
> # Load the trajectoty files.
> mol new prmtop type parm7
> mol addfile dcd type dcd waitfor all
>
> set numfrm [molinfo top get numframes]
> set systm [atomselect top all]
>
> #######################################
> #Get all the atomselect.
> for {set i 1} {$i <= 1536} {incr i} {
> set atom($i) [atomselect top "water within 5.5 of resid $i"]
> }
>
> for {set i 0} {$i < $numfrm} {incr i} {
>
> ######### Get all the connected list for each of the water molecule.
> for {set j 1} {$j <= 1536} {incr j} {
> $atom($j) frame $i
> $atom($j) update
> set list($j) [$atom($j) get resid]
> set list($j) [simplifylist $list($j)]
> puts "frame $i residue $j: $list($j)\n"
> }
> }
>
> exit
>
>