3
$\begingroup$

I have a set of weighted points in 3D space (in fact, a molecule) and I'm trying to align the principal axes of this set with the $x$, and $y$ and $z$ axes. To do so, I've first translated my points so their barycenter coincides with the origin. Then, I've calculated the moment of inertia tensor $I$ and its eigenvalues ($\lambda_i$) and eigenvectors ($V_i$).

Then, I need to build the rotation matrix associated with the rotation bringing my vectors $V_i$ onto $(x,y,z)$. I assumed that would be the inverse of the matrix formed by the vectors $V_i$ in columns, but it is not. So, how would you write this matrix from what I have already calculated?

2 Answers 2

2

We have $I V_i = \lambda_i V_i$, where $\lambda_i$ is real. Let $M$ be the matrix whose columns are the normalized eigenvectors of $I$. Then $M$ is orthogonal, $M^T M = M M^T = \mathbb{I}$. Thus, $M^T I M = D = \mathrm{diag}(\lambda_1,\lambda_2,\lambda_3)$ and $M^T M_i = e_i$. Notice that the last equation implies, for example, that $M^T M_1 = e_1 = (1\ 0\ 0)^T$. That is, the transpose of $M$ brings the principal axes to the Cartesian axes.

The simplest way to remember how the various objects transform is to look at the kinetic energy, for example, $$\begin{eqnarray*} \frac{1}{2}\omega^T I \omega &=& \frac{1}{2}\underbrace{\omega^T M}_{\xi^T} \underbrace{M^T I M}_D \underbrace{M^T \omega}_\xi \\ &=& \frac{1}{2} \sum_{i=1}^3 \lambda_i\xi_i^2, \end{eqnarray*}$$ where $\omega$ is the angular velocity in the original frame and $\xi$ is the angular velocity with respect to the principal axes.

If you are using Mathematica, note that the rows of the matrix it spits out are the unnormalized eigenvectors. In another computing environment the convention may be something else, so be careful. Without more information it is impossible to tell where this goes wrong for you.

  • 0
    I don't understand… Are you saying M is the rotation matrix bringing the principal axes to the cartesian axes? Because it doesn't work… and I've tried more than once! PS: My eigenvectors are already normalized. PPS: I don't know what you define as $\omega$ or $\xi$ or $D$2012-05-14
  • 0
    @F'x: I added some clarification. If this doesn't help I'll need to see your code or computations to go further.2012-05-14
1

I know that this this post is really old, but I was just working on the exact same problem, and I thought that maybe at some point in the future someone else might be facing the same issue as well.

I'm writing an algorithm for which I need to compare molecular structures, represented as sets of points in space, and I need these structures to be oriented in exactly the same way so that I can see if the coordinates for each point are the same, and with that if the structures are the same.

I tried the method with multiplying by the transpose of the eigenvector matrix of the inertia tensor. It's correct, the transpose of the matrix containing the eigenvectors orients the molecule properly when applied to it, and the same can be achieved with rotation matrices. With rotation matrices, I could align the inertia tensor that corresponds to the axis of maximum rotation symmetry with the z-axis, by rotating the molecule first around the x- axis to bring it in to the xz plane and then around the y- axis, to bring it at the same position as the z-axis of the cartesian coordinate system, and that works just fine.

But, I found that when this method is used, the resultant orientation doesn't discriminate between the relative directions in space, so the choice of the x- and y- axes is arbitrary, as well as the + and - directions on each axis. I think that maybe that's why the OP had the impression that it doesn't work.

The way that I solved this is by just going through the different choices of the eigenvector directions for the x- and y- axes as well as the + and - directions, which means that for every structure I have to do the checking many times, but for my system it was OK, because it's quite small. I think that for a small system this is a better solution, or perhaps it is possible to also specify the other two eigenvectors to correspond with either the x- or the y- axis, in which case the orientation would be absolute.

Anyway, this is a much simpler approach than any other method avialable.

But, since I am not a mathematician, I'm a chemist, I would like to ask, if my observation that when the transpose of the eigenvector matrix is applied, or when one eigenvector is fixed to be aligned with the z-matrix, the resultant orientation doesn't completely distinguish between the relative directions in space, is correct. And also, if this is correct, would aligning the other two eigenvectors with the x- and y- axes respectively fix this problem? And also, if this is correct, why is this? Is it related to some property of eigenvectors, or to some property of the inertia tensor, or of the inertia eigenvectors, or does it come from some property of matrix multiplication?

If anyone is using this method, I would recommend the diagonalization subroutines available in Lapack for Fortran, rather than Mathematica. The result is always an orthonormal matrix, and calculations are very straightforward.

  • 0
    Hi, I'm a chemist too. What you said is right, remember that if $v$ is an eigenvector of $R$, and the corresponding eigenvalue is $\lambda$: $R v = \lambda v$, but it is true for any $v' = \alpha v$ where $\alpha$ is a scalar. Notice what will happend if $\alpha = -1$... If $M_1$ and $M_2$ are the $N\times 3$ cartesian matrixes of two molecules, I think that a better solution is to find the rotation matrix that best" overlap $R M_1$ and $M_2$. But this method is only suitable when you have paired atoms (when you know exactly a the bijection between atoms in booth molecules).2015-09-02
  • 0
    Oh, thank you :) Sorry, I just saw your reply. It's been a while. You're right, it makes sense. I never thought of that Then, I guess the simplest solution is to consider all possible combinations of directions, when the system is small enough.2016-03-09
  • 0
    Yes, although if you know the bijection above finding the rotation matrix would be simple and fast enough (almost instantaneous) molecules up to hundreds or thousands of atoms. I coded that routing if you are interested I can share it to you.2016-03-09
  • 0
    Yes, I would like to understand your method better. I completed this program last year, and I'm not using it anymore, but I still might use it in the future, it's a program for searching molecular structures. And, the orientation-checking part is still the weak point.2016-03-10
  • 0
    But, I don't understand what you mean exactly when you say the bijection between the atoms in the molecules. I only have them as coordinates in space, they could be anywhere, and I don't have a function that relates the coordinates.2016-03-10
  • 0
    What I did was to move the center of mass of the molecule to the (0,0,0) point. Then, find the eigenvectors, then look at the rotation symmetry around each eigenvector, and then choose the eigenvector with the largest rotational symmetry to be the z-axis. If there is no rotation symmetry, then choose the z-axis eigenvector to be the one corresponding to the largest eigenvalue. Maybe this seems like overkill... It's coming from a chemist's mindset :)2016-03-10
  • 0
    OK, I see what you're saying now. To relate the two eigenvector matrices correct?2016-03-10
  • 0
    I am chemist too. I mean: The nuclei positions of molecule with $N$ atoms can be represented by a matriz of $N \times 3$ (you can hold in other $N$-item list the element specification). Leaving aside the indistinguishableness of particles (for example think about a molecule with all atoms of different elements), if you have two molecules (or the same molecule before and after a modification of its nuclei's coordinates), two alternatives arises:2016-03-10
  • 0
    1 - That you know for each atom $i$ (row $i$) of the first molecule (matrix) which atom (which row) correspond to $i$ in the other molecule. (so you can make a bijective function between rows from two matrix, that is you know exactly which row $j$ of the second matrix correspond to the row $i$ of the first matrix) 2 - That you do not know. The method of finding the best rotation matrix requires the case 1.2016-03-10
  • 0
    It is not always possible. Normally you can't when they come from different sources. But If you are working with a molecule and transform it successively (always using x,y,z coordinates) there shouldn't be a problem.2016-03-10
  • 0
    Yes, I know. Here is how it works. I generate one structure held in an N x 3 matrix. Then another, also in an N x 3 matrix. Now, I want to know if they are the same. Because if they are the same I will discard the second one. Because I'm searching for structures, so if I've already found it before, I don't need it again. I add the elements at the end. They are not important at this stage.2016-03-12
  • 0
    So, my approach to seeing if the new structure is the same as one generated before, was, put them in them in 0,0,0 position, orient them the same way, and check if all of the coordinates are the same. So, for this, I needed absolute orientation. So that they are oriented in exactly the same way. I ran in to the original post above when I was trying to figure out how to orient them properly :)2016-03-12
  • 0
    I don't know which row in matrix1 corresponds to which row in matrix2. It's not the same structure, they're two randomly generated structures.2016-03-12
  • 0
    I got you, I was in the same situation and that is why I arrive to this post. Notice that checking magnitude of eigenvalues of inertia tensor matrix will suffice in some case to state that the molecules are different., and it is very fast. Can I ask you how did you manage rotational symmetry?2016-03-12
  • 0
    The main rotation axis is always along an eigenvector of the inertia tensor. You put your molecule in (0,0,0). Then find the inertia eigenvectors. Then align one eigenvector to be the z-axis. Then rotate the molecule around the z-axis, for an angle of 2pi/8, check if all the coordinates and corresponding elements are the same (definition of molecular symmetry), if not then 2pi/7, check if it's the same... until you get the same molecule after rotation. For a rotation of 2pi/k where the result is the same molecule, k is the order of the rotation axis around that eigenvector.2016-03-14
  • 0
    Then move on to the next eigenvector. Once you've done this for all three eigenvectors, the one with the largest k is the main rotation axis, if such an axis exists.2016-03-14
  • 0
    By the way, about the eigenvalues. Will the eigenvalues be the same, no matter where the molecule is in space? Meaning, if they are not the same, the coordinates are definitely different structures?2016-03-14
  • 0
    If you want to get the full point group symmetry, once you have the main rotation axis, it's easy. There are only so many combinations. You don't need to find all of the symmetry elements.2016-03-14
  • 0
    Oh, and also, for the absolute orientation. After you align the eigenvector corresponding to the rotational axis of maximum rotational symmetry to the z-axis, you still need to align the other two to the y- and x- axis. That way, any two identical structures will be aligned the same way. Except that they may be flipped...2016-03-14