VMD-L Mailing List

From: Ackad, Edward (eackad_at_siue.edu)
Date: Mon May 23 2022 - 08:47:00 CDT

This worked! Thanks!
It had started with ~40Gb with the dcd loaded but would balloon to beyond 512 Gb before swapping and crashing vmd.

_________________________________________________________
Edward Ackad, Ph.D<https://urldefense.com/v3/__http://www.siue.edu/*7Eeackad__;JQ!!DZ3fjg!7vRIe5g_q_FCIgF4bc7nJEWxtvYBTHPRZk02i0s1lT67mU-I-MApOrmztEcGxvdNEY1MF7NzCKN0egnjzQ$ >
Associate Professor of Physics
Computational Nanophotonics/Biophysics
Southern Illinois University Edwardsville
(618) 650-2390

________________________________
From: Axel Kohlmeyer <akohlmey_at_gmail.com>
Sent: Sunday, May 22, 2022 12:24 PM
To: Ackad, Edward <eackad_at_siue.edu>
Cc: vmd-l_at_ks.uiuc.edu <vmd-l_at_ks.uiuc.edu>
Subject: Re: vmd-l: measure sasa and memory issue

Your problem seems to be the usual one that atom selections created with the atomselect command persist until they go out of scope, but that will not happen within loops.

So *everywhere* you have an [atomselect ....] you will need to do:
# at the beginning of the loop iteration
set sel1 [atomselect ...]
set sel2 [atomselect ...]
..
# now use $sel1 $sel2
..
# at the end of the loop iteration
$sel1 delete
$sel2 delete
..

If you are looping over a long trajectory and are already short on RAM, not doing this can put you over the edge.

On Sun, May 22, 2022 at 12:55 PM Ackad, Edward <eackad_at_siue.edu<mailto:eackad_at_siue.edu>> wrote:
Hi all,
I'm trying to write a script to renormalize the SASA measurement from GaMD runs (the scripts provided only do the reweighting for other measures). So I am trying to calculate the SASA for a small set of residues (using -restrict, which is why I can't use mdtraj's version). I can fit the stride=1 trajectory into memory without issue but once I start to calculate the SASA the memory blows up. I use a set for each frame and output to a file so I don't know why the memory requirements keep escalating. Is there a free or wipe memory option or some issue with my code (shown below)? Or any other ideas to make this not-so-memory intensive?
Thanks!

proc run_sasa {oname last_10ps_frame dcdstride } {
  # create the dictionary of the values for each reside
  set hydrophob [dict create ILE -9.73 ALA -9.57 PHE -9.23 LEU -8.66 MET -7.65 PRO -7.55 VAL -8.26 TRP -8.62 CYS -5.63 GLY -4.57 THR -4.06 SER -3.56 TYR -0.05 GLN 2.2 ASN 2.76 HIS 1.41 >

  set rdis 1.4
  set fname "$oname.dat"
  set outfile [open $fname w]
  set molid top
  set nes "protein and resid 40 to 53"
  set nesphobic "protein and resid 40 to 53 and hydrophobic"
  set bucket "protein and resid 48 51"

  for {set frame 0} {$frame < [molinfo $molid get numframes]} {incr frame} {
    set printline [measure sasa $rdis [atomselect $molid "protein" frame $frame]]
    set printline2 [measure sasa $rdis [atomselect $molid "protein" frame $frame] -restrict [atomselect $molid $bucket frame $frame ] ]
    set printline3 [measure sasa $rdis [atomselect $molid "protein" frame $frame] -restrict [atomselect $molid $nes frame $frame ] ]
    set printline4 [measure sasa $rdis [atomselect $molid "protein" frame $frame] -restrict [atomselect $molid $nesphobic frame $frame ] ]

    # calculate the sasa's for each hydrophobic res
    set totalhydrophobic 0.0
    set totalsasa 0.0
    # loop over all the NES residues
    for { set i 40 } {$i < 54 } {incr i} {
      #use just one atom (alpha == name CA) so we only get 1 name back instead of a list
      set ressel [atomselect top "resid $i and alpha"]
      #determine what the coefficent is using the name of the residue as the index for the dictionary
      set coef [dict get $hydrophob [$ressel get resname]]
      set thisres [measure sasa $rdis [atomselect $molid "protein" frame $frame] -restrict [atomselect $molid "protein and resid $i" frame $frame ] ]
      # sum each residues value scaled by the coef
      set totalsasa [expr $totalhydrophobic + $thisres*$coef]
      # if the residue is hydrophoic (coef <0 then sub it so we can get a positive number)
      if {$coef <0 } {
        set totalhydrophobic [expr $totalhydrophobic - $thisres*$coef]
      }
    }

    #set the actual frame of the frame to match with the output in the log file
    if { $frame <= $last_10ps_frame} {
      set outframe [expr $dcdstride*$frame*10000]
    } else {
      set outframe [expr $dcdstride*$frame*100000]
    }

    puts $outfile "$outframe $printline $printline2 $printline3 $printline4 $totalsasa $totalhydrophobic"
  }
  close $outfile

} 

--
Dr. Axel Kohlmeyer  akohlmey_at_gmail.com<mailto:akohlmey_at_gmail.com>  https://urldefense.com/v3/__http://goo.gl/1wk0__;!!DZ3fjg!7vRIe5g_q_FCIgF4bc7nJEWxtvYBTHPRZk02i0s1lT67mU-I-MApOrmztEcGxvdNEY1MF7NzCKPfPOisBA$ 
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.