From: John Stone (johns_at_ks.uiuc.edu)
Date: Fri Jan 27 2012 - 17:17:43 CST

Can you be more specific about what the symptom of your "memory issue" is?
Also, I would suggest that you replace the "mol load" command with:
    mol new MD_protein_in_water/build/solvated.psf waitfor all
    mol addfile MD_protein_in_water/prod/equil.dcd waitfor all
 
How big is your trajectory, and how much memory does your machine have?
You won't be able to load a trajectory larger than the available RAM
without encountering performance problems. Since you are only looping
over 5 frames in your current script, you should probably change the
"mol addfile" command example I gave above so that it only loads the same
frame range that you plan to process.

Cheers,
  John

On Fri, Jan 27, 2012 at 05:26:29PM -0500, Mustafa Tekpinar wrote:
> Hi,
> I have a molecular dynamics trajectory of (protein+ligand+water+ions) box.
> I want to cut a sphere whose center is center of mass of (protein+ligand).
> After determining the vector from center of mass to the farthest atom, I
> want add some buffer distance to that radius. Then, I want to save all
> atoms (except ions) in that radius. I have a script but it seems like
> there is a memory issue. I'd be glad if some experts points me my
> mistake(s).
> Here is my script:
> //=====================================================================
> set hydra_shell 7.0
> set frm_beg 0
> set frm_end 5
> #Load protein box simulation
>
>
> mol load psf MD_protein_in_water/build/solvated.psf dcd
> MD_protein_in_water/prod/equil.dcd
> #Select protein and (if exists) ligands
>
>
> set protAndLigand [atomselect top "not water and not ions"]
> for {set frame $frm_beg} {$frame < $frm_end} {incr frame} {
> # get the correct frame
>
>
> $protAndLigand frame $frame
> ## Determine the center of mass of the (protein+ligands) and store the
> coordinates
>
> set CoM [measure center $protAndLigand weight mass]
> set x1 [lindex $CoM 0]
> set y1 [lindex $CoM 1]
> set z1 [lindex $CoM 2]
> ## Determine the distance of the farthest atom from the center of mass
>
>
> foreach atom [ $protAndLigand get index] {
> set pos [lindex [[atomselect top "index $atom" frame $frame ] get
> {x y z}] 0]
> set x2 [lindex $pos 0]
> set y2 [lindex $pos 1]
> set z2 [lindex $pos 2]
> set dist [expr (($x2-$x1)*($x2-$x1) + ($y2-$y1)*($y2-$y1) +
> ($z2-$z1)*($z2-$z1))]
> if {$dist > $max} {
> set max $dist
> }
> }
> set max [expr sqrt($max)]
> #Add some buffer of water to max distance
> set max [expr $max + $hydra_shell]
> set maxSqrd [expr $max*$max]
> puts "Radius squared is $maxSqrd]"
> #Select protein and (if exists) ligand +plus water within the radius!
>
>
> set pwSphr [atomselect top "(not water and not ions) or (same residue
> as {water and ((x-$x1)*(x-$x1) + (y-$y1)*(y-$y1) +
> (z-$z1)*(z-$z1))<($maxSqrd)})" frame $frame]
> #write coordinates.
>
>
> $pwSphr writepdb sphere${frame}.pdb
> $pwSphr delete
> }
> $protAndLigand delete
> exit
> //=====================================================================
> Thanks in advance,
> mtekpinar

-- 
NIH Resource for Macromolecular Modeling and Bioinformatics
Beckman Institute for Advanced Science and Technology
University of Illinois, 405 N. Mathews Ave, Urbana, IL 61801
http://www.ks.uiuc.edu/~johns/           Phone: 217-244-3349
http://www.ks.uiuc.edu/Research/vmd/       Fax: 217-244-6078