From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Mon Oct 15 2012 - 14:12:33 CDT

i am adding some comments below.

On Mon, Oct 15, 2012 at 5:53 PM, Vikas Varshney <vv0210_at_gmail.com> wrote:
> Dear John, Axel,
>
> I have a trajectory of my system of interest (CNTs) which I want to
> visualize with certain color scheme associated with retrival of bonds
> information.
>
> Basically, I am stretching a CNT and would like to color the atoms based on
> the longest bond distance. I have looked over the forum to get some insights
> on how to color atoms based on user defined parameters and have the
> following TCl script which hopefully will do the job. I need some help with
> how to retrive bond information so that I can calculated bond lengths for
> each atom with its bonding pairs.
>
> I am using Axel's topotools package to load the lammps data file which I
> think is loading bond information as well. Below is the template of the TCl
> script
>
> -------------------
> #/usr/bin/tclsh
>
> set lmpdata [atomselect 0 all] #This contains topology information
> set lmptrj [atomselect 1 all] #This is the trajectory file

this looks as if you are loading the dump into a different (vmd) "molecule".
it would be much better to load the trajectory "into" the same "molecule".
if your dump contains the atom ids, the order should be preserved.

> for {set i 0} {$i < [molinfo 1 get numframes]} {incr i} {
> $lmptrj frame $i
> foreach atomind [$allsel get index] {
> set csel [atomselect 1 "index $atomind" frame $i]

> # Get indices of other atoms that are connected to atom "atomind"

this is hugely inefficient.

since you seem to be using a conventional force field
or at least a topology file that has contains bond information,
you can do the following.

topotools has a command to retrieve a "bond list" this is a
list of lists which contains a list of all bonds with each entry
containing the two atom indices forming the bond and additional
information, as requested. this information will not change
over the trajectory, so it can be retrieved up front and stored.

then you should take the selection for all atoms (also defined
outside the loop) and advance it to the current frame as in your
script fragment and then retrieve all coordinate triples via
$sel get {x y z}. this will also produce a list of lists.

now you can loop over the bond list and for each bond get
the two indices (e.g. using lassign) and then use lindex
on the coordinate list to get the two coordinate vectors and
use vecdist to compute the distance (i.e. bond length).

now the big question is, how do you want to convert this
information into a color code, since color is stored on a
per atom basis and each atom has multiple bonding partners
(3) and you have 4 user fields to store per atom per frame
information. so theoretically it should be possible with
3 representations using one of 3 user fields to handle
this, but you'd have to figure out a way to select the bonds
accordingly, which i currently have no idea how to do.

> # Calculate distance of each bond pair with one atom "atomind" using vector
> operators or some other functions
> # Get the largest distance (ldist)
> # Color atom
> $csel set user $ldist
> $csel delete

creating and deleting selections during a loop and particularly
on a per atom basis, is going to be very slow, and in your
case not really needed.

axel.

> }
> }
>
> $lmpsdata delete
> $lmptrj delete
>
> -------------------
>
> I will highly appreciate any other suggestions as well.
>
> Best Regards,
> Vikas

-- 
Dr. Axel Kohlmeyer  akohlmey_at_gmail.com  http://goo.gl/1wk0
International Centre for Theoretical Physics, Trieste. Italy.