From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Thu Nov 26 2015 - 07:27:01 CST

On Thu, Nov 26, 2015 at 7:30 AM, Frank Zack
<frankzack123_at_remove_yahoo.de> wrote:
> Hi all,
>
> i've written a simple,dirty tcl-script to coarse-grain 40 Chains of EPR
> Poylmer (All-Atom,1000-repeating units each chain) to a Martini-like
> representation.
> It perfectly works for the first 4 chains...however during coarsegraining
> the 5th chain vmd gets killed and i have no idea why...

these cases are almost always caused due to memory leaks from creating
atom selections inside loops and not deleting them.

axel.

>
> Do you have any suggestions?
>
>
> proc cg {} {
>
> #package require pbctools
> #package require topotools 1.5
> #topo readlammpsdata stage21.data
>
> set sel [atomselect 0 "residue 40 to 80"]
> set r1 [$sel get residue]
> set firstRes [lindex $r1 0]
> set lastRes [lindex $r1 end]
> set oldcell [pbc get]
> puts "firstRes:$firstRes"
> puts "lastRes:$lastRes"
>
> #always try to free mem to prevent vmd from getting killed
> $sel delete
>
> set counter 0
>
> set cname "1"
> set ppname "13"
>
> #create new mol+atoms before
> set sel_c [atomselect 0 "name '$cname'" ]
> set sel_pp [atomselect 0 "name '$ppname'" ]
> set num_c [$sel_c num]
> set num_pp [$sel_pp num]
> set nrbeads [expr (($num_c/2)+$num_pp)]
>
> $sel_c delete
> $sel_pp delete
>
> set mol [mol new atoms $nrbeads]
> animate dup $mol set logfile [open "cgepr.log" "w"]
> for {set curchain $firstRes} {$curchain <= $lastRes} {incr curchain}
> {
> puts "coarsing chain $curchain"
> set sel [atomselect 0 "residue $curchain"]
> set maxid [lindex [$sel list] end]
> set maxid [expr ($maxid-1)]
> set bbatoms {$cname $ppname}
> set bb [atomselect 0 "residue $curchain and name '$cname'
> '$ppname'"]
> set bblist [$bb list]
>
> $bb delete
> set curid [lindex [$sel list] 0]
> $sel delete
>
> puts "maxid:$maxid"
> puts "curid:$curid"
>
> while {$curid < $maxid} {
> #ADD BEAD
> set cursel [atomselect 0 "index $curid"]
> set nextid [expr ($curid+3)]
> set nextsel [atomselect 0 "index $nextid"]
> set curname [$cursel get name]
> set nextname [$nextsel get name]
> set cg 0
>
> #ignore-terminators
> set test [lsearch $bbatoms $curname]
> if {$test==-1} {incr curid;continue}
> if {$curname==$ppname} {incr curid;continue}
>
> if {$curname==$cname} {
> if {$nextname==$cname} {
> if { $nextid>$maxid } {break}
> set tmp [atomselect top "index
> $curid to $nextid"]
> set cg [measure center $tmp weight
> {12 1 1 12}]
>
> incr counter
> set sel2 [atomselect $mol "index
> [expr ($counter+1)]"]
> $sel2 set name "PE"
> $sel2 set type "1"
> $sel2 set {x y z} $cg
> $sel2 set resid $curchain
>
> incr curid
> incr curid
> incr curid
> $sel2 delete
> $tmp delete
> continue
>
> }
> if {$nextname==$ppname} {
> if { $nextid>$maxid }
> #overnext ist das methyl am pp
> set overnext [expr ($nextid+1)]
> set tmp [atomselect 0 "index $curid
> to $overnext"]
> set cg [measure center $tmp weight
> {12 1 1 12 1 12}]
> #lappend namelist "PP"
> #lappend typelist 2
> #set xyz [list $cg]
> #lappend xyzlist $cg
> #lappend idlist $counter
> #lappend reslist $curchain
> #lappend resiudelist $curchain
> incr counter
>
> set sel2 [atomselect $mol "index
> $counter"]
> $sel2 set name "PP"
> $sel2 set type "2"
> $sel2 set {x y z} $cg
> $sel2 set resid $curchain
>
> incr curid
> incr curid
> incr curid
> $sel2 delete
> $tmp delete
> continue
> }
> }
> }
> }
> #puts "adding $counter atoms"
> #set mol [mol new atoms $counter]
> #animate dup $mol
>
> #set sel [atomselect $mol all]
> #$sel set name $namelist
> #$sel set type $typelist
> #$sel set {x y z} $xyzlist
> #$sel set resid $reslist
>
> #ADDBONDS
>
> puts "adding bonds"
> puts "first:$firstRes, last:$lastRes"
> for {set i $first } {$i <=$last } {incr i} {
> #puts "i:$i"
> set sel [atomselect $mol "all resid $i" ]
> set l [$sel get index]
> set firstIndex [lindex $l 1]
> set lastIndex [lindex $l end]
> puts $logfile "outerloop for chain $i"
> puts $logfile "firstIndex: $firstIndex"
> puts $logfile "lastIndex: $lastIndex"
> for {set k $firstIndex} {$k <=$lastIndex} {incr k} {
> set i1 $k
> set i2 [expr ($k-1)]
> ::TopoTools::addbond $mol $i1 $i2 1 1
> puts $logfile "added bond for chain $i, index1:$i1,
> index2:$i2"
> }
> }
> set sel [atomselect $mol all]
> puts $logfile "guessing angles"
> ::TopoTools::guessangles $sel
> #pbc set cell [lindex $oldcell 0]
> ::PBCTools::pbcset [lindex $oldcell 0]
> topo writelammpsdata -sel $sel cg.data
>

-- 
Dr. Axel Kohlmeyer  akohlmey_at_gmail.com  http://goo.gl/1wk0
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.