EigenDecompositionImpl.getV() returns eigen matrix with indeterminate 
determinant
---------------------------------------------------------------------------------

                 Key: MATH-416
                 URL: https://issues.apache.org/jira/browse/MATH-416
             Project: Commons Math
          Issue Type: Improvement
    Affects Versions: 2.1
         Environment: Mac OS X 10.6.4
            Reporter: Tom Milac



A call to EigenDecompositionImpl.getV() returns a RealMatrix the columns of 
which are the eigenvectors of the matrix with which EigenDecompositionImpl is 
constructed. Because EigenDecompositionImpl works only with real, symmetric 
matrices, the eigenvectors (columns) returned are orthogonal.  In addition, the 
eigenvectors are normalized to have 2-norm = 1. Unfortunately, for 3x3 input 
matrices, the determinant of the eigenvector matrix is indeterminate, sometimes 
+1 and other times -1.  The -1 output can be 
'repaired' simply by multiplying the matrix by -1.  Example code is included 
below.

Because the columns are eigenvectors, the result with either determinant is 
correct.  However, in the case that the matrix returned is to be interpreted as 
specifying a coordinate system, the principal axes of a body in my case,  the 
+1 result specifies a right-handed coordinate for the principal coordinate 
system of the body, and the -1 result specifies a left-handed coordinate 
system.  Once discovered, this indeterminacy is easy to deal with, but 
nevertheless an inconvenience.

I believe it would improve EigenDecompositionImpl.getV() to return an 
eigenvector matrix with a consistent determinant = +1.

Tom Milac

---------------------------------------------------------

import org.apache.commons.math.geometry.NotARotationMatrixException;
import org.apache.commons.math.geometry.Rotation;
import org.apache.commons.math.linear.Array2DRowRealMatrix;
import org.apache.commons.math.linear.EigenDecompositionImpl;
import org.apache.commons.math.linear.InvalidMatrixException;
import org.apache.commons.math.linear.LUDecompositionImpl;
import org.apache.commons.math.linear.RealMatrix;

/**
 *
 * @author Tom Milac
 */
public class BugReport {

    /**
     * Moment of inertia tensor #1.
     */
    public static final double[][] MOI1 =
            {{128.52722633757742, -29.11849805467669, 8.577081342861376},
             {-29.11849805467669, 521.3276639228706,  35.512665035385666},
             {8.577081342861376,  35.512665035385666, 490.2479495932442}};

    /**
     * Moment of inertia tensor #2.
     */
    public static final double[][] MOI2 =
            {{440.09350934414175, 44.23154125186637, -9.41455073681743},
             {44.23154125186637,  387.1291457565648, -38.07596950448303},
             {-9.41455073681743, -38.07596950448303, 762.0451513430822}};

    /**
     * Constructor.
     */
    public BugReport() {
    }

    /**
     * Main.
     */
    public static void main(String[] args) {

        // Compute the principal axes (eigenvectors) of the #1 moment
        // of inertia tensor.
        RealMatrix moi1  = new Array2DRowRealMatrix(MOI1);
        RealMatrix axes1 = null;
        EigenDecompositionImpl eigenDecompositionImpl = null;
        try {
            eigenDecompositionImpl = new EigenDecompositionImpl(moi1, 0.0d);

            axes1 = eigenDecompositionImpl.getV();
        } catch (InvalidMatrixException ex) {
            System.err.println("MOI1: InvalidMatrixException thrown.");
            System.err.println("Exiting ...");
            System.exit(-1);
        }

        // Compute the principal axes (eigenvectors) of the #2 moment
        // of inertia tensor.
        RealMatrix moi2  = new Array2DRowRealMatrix(MOI2);
        RealMatrix axes2 = null;
        try {
            eigenDecompositionImpl = new EigenDecompositionImpl(moi2, 0.0d);

            axes2 = eigenDecompositionImpl.getV();
        } catch (InvalidMatrixException ex) {
            System.err.println("MOI2: InvalidMatrixException thrown.");
            System.err.println("Exiting ...");
            System.exit(-1);
        }

        // Determinant of axes 1 eigenvector matrix = -1.  If the matrix
        // is interpreted as a Rotation, throws and Exception.
        System.out.print("Determinant of the #1 moment of inertia tensor = ");
        System.out.println(new LUDecompositionImpl(axes1).getDeterminant());

        try {
            Rotation axes1_rotation = new Rotation(axes1.getData(), 1.0E-7);
        } catch (NotARotationMatrixException ex) {
            System.out.println("NotARotationMatrixException thrown for 
'axes1'.");
        }

        System.out.println();

        // Determinant of axes 2 eigenvector matrix = +1.
        System.out.print("Determinant of the #2 moment of inertia tensor = ");
        System.out.println(new LUDecompositionImpl(axes2).getDeterminant());

        try {
            Rotation axes2_rotation = new Rotation(axes2.getData(), 1.0E-7);
        } catch (NotARotationMatrixException ex) {
            System.out.println("NotARotationMatrixException thrown for 
'axes2'.");
        }
    }
}

-- 
This message is automatically generated by JIRA.
-
You can reply to this email to add a comment to the issue online.

Reply via email to