[gmx-users] Computing vectors normal (perpendicular) to molecules

[Mark Abraham]

On 2/07/2012 11:51 AM, Mark Abraham wrote:
> On 2/07/2012 6:50 AM, Andrew DeYoung wrote:
>> Hi,
>>
>> Suppose that I have a system of benzene molecules (in reality, my 
>> system is
>> more complicated, but for my question, it will be simpler to consider 
>> just
>> an ensemble of benzene molecules).
>>
>> I would like to find a vector normal (i.e., perpendicular) to the 
>> plane of
>> each benzene molecule in my system.  I want to track these "normal 
>> vectors"
>> over time; eventually, I would like to calculate a histogram of 
>> angles theta
>> between the vectors normal to the benzene molecules and some fixed
>> "laboratory" axis.
>>
>> If you have time, can you please tell me if there is any Gromacs utility
>> that already does this -- i.e., computes normal vectors like this?  
>> It does
>> not seem so, even with the functionalities provided by g_angle.
>>
>> This may not be the most rigorous idea, but one way to do this might 
>> be to
>> compute the vector cross product of two vectors.  Suppose that a benzene
>> molecule has carbon atoms named C1, C2, C3, C4, C5, and C6 (arranged 
>> in that
>> order around the ring). Then I could, for example, find two vectors:
>>
>> \vec{C2-C1}
>> \vec{C2-C3}
>>
>> where \vec{C2-C1} is the vector from C2 to C1, and \vec{C2-C3} is the 
>> vector
>> from C2 to C3.
>>
>> The cross product of those two vectors gives a vector that is 
>> perpendicular
>> to both vectors.
>>
>> Is there any way to compute cross products in the Gromacs utilities?  It
>> does not seem so, but I am not certain.
>
> Yes there are tools that can do this, e.g. g_sgangle or perhaps 
> others. Manual 7.4 is a good place to start your search, then to 

...chapter 8, and then to...

> appendix D.
>
>> So, I am thinking that I will need to use g_traj to extract the x, y, 
>> and z
>> coordinates of the atoms of interest (i.e., C1, C2, and C3 in each of 
>> the
>> many benzene molecules in the system).  Then I will need to write a C or
>> Fortran script to find the relevant vectors and compute their cross
>> products, at every timestep in the trajectory.
>>
>> However, to complicate things even more, I may in the future wish to 
>> make
>> dynamic selections of benzene molecules, using g_select.  I may want to
>> consider only the benzene molecules satisfying z<12, for example.  
>> Since the
>> benzene molecules obviously move, the particular benzene molecules 
>> that I
>> consider for the vector/cross product calculation will change over 
>> time. I
>> am able (using the -oi option in g_select in conjunction with a 
>> script of my
>> own) to generate an index file with the indices (whose groups 
>> correspond to
>> the simulation timestep number) of the C1, C2, and C3 atoms in benzene
>> molecules whose centers-of-mass satisfy z<12.  But is there any way I 
>> can
>> feed this index file to g_traj without having to call g_traj as many 
>> times
>> as I have timesteps in my trajectory?
>
> Not until GROMACS 5.
>
> Mark