← Back to team overview

yade-dev team mailing list archive

Re: EigenDecomposition wrapper problem

 

The ordering of values will be reflected in the rotation matrix, and the code always uses both simultaneously. So, it should not be a problem to get different ordering. It is just a matter of convention, with no effect on the result.

Bruno

Vincent Richefeu a écrit :
If 'eigen' is used instead of Wm3 to find principal axes of inertia, eigenvalues must NOT be ordered (as Wm3 do).
So don't forgot to disable "ordering" in 'eigen' (if possible)

Le 8 févr. 2010 à 17:05, Anton Gladky a écrit :

Thanks for comments.
You are right, eigendecomposition in yade uses in:
- Tetra.cpp { * The inertia tensor is in global coordinates; by eigendecomposition, we find principal axes}
 - Clump.cpp

We have a next comment there:
http://bazaar.launchpad.net/%7Eyade-dev/yade/trunk/annotate/head%3A/pkg/dem/DataClass/Clump.cpp#L240

/*! @bug eigendecomposition might be wrong. see http://article.gmane.org/gmane.science.physics.yade.devel/99 for message. It is worked around below, however.
    */
Vaclav, it seems, it is your comment.

So, we just let it be?
______________________________

Anton Gladkyy


2010/2/8 Vincent Richefeu <richefeu@xxxxxxxxxxxxxxx <mailto:richefeu@xxxxxxxxxxxxxxx>>

    if numerical recipes, numpy and eigen give the same result, we
    can have some confidence with this result.
    It means that Wm3 result is "less good".
    However, I don't know where eigen values and vectors are needed
    in yade (maybe for the computation of inertia?).

    I don't know if my comments help (sorry), but all I can say is
    that eigenvalues from Wm3 seems correct but unordered, and
    eigenvectors seems not well converged (maybe it is not a problem)


    Le 8 févr. 2010 à 16:21, Anton Gladky a écrit :

    Vincent, what do you think, if we use these "new values" it will
    not work?
    ______________________________

    Anton Gladkyy


    2010/2/8 Vincent Richefeu <richefeu@xxxxxxxxxxxxxxx
    <mailto:richefeu@xxxxxxxxxxxxxxx>>

        Hi,

        At first glance, it seems the eigenvalues (and thus the
        eigenvectors) are not ordered with Wm3.
        The difference can also be due to the fact that the solution
        is not converged.

        If needed, I can send a peace of code (based on "jacobi" in
        numerical recipes) for doing that task

VR

        Le 8 févr. 2010 à 13:46, Anton Gladky a écrit :

        Thanks, Vaclav!
        So, it seems Wm3 is broken for that function?

        I have taken another matrix, which is positive definite:
        1 7 3
        6 2 4
        2 5 3
        ========================
        Eigen gives:

        Rot:
        -0.584754 -0.665682 -0.466054
        -0.617452 0.703421 -0.301128
        -0.526133 -0.266722 0.841805
        Diag:
        11.0907 0 0
        0 -5.19482 0
        0 0 0.104141

        ========================
        Numpy gives (see checkEigen.py):
        Rot:
        array([[-0.58475406, -0.66268234, -0.46225138],
               [-0.61745205,  0.70025076, -0.29867115],
               [-0.52613273, -0.26552022,  0.83493665]])
        Diag:
        array([ 11.09067395,  -5.1948153 ,   0.10414135])

        So, it was similar to Eigen results.
        ========================
        But Wm3:
        Rot:
        0.710682 0.404879 0.57533
        -0.699213 0.316215 0.641178
        0.0776719 -0.857952 0.507825
        Diag:
        -5.55916 0 0
        0 0.10998 0
        0 0 11.4492
        ========================

        Is it the Wm3 library problem?

        ______________________________

        Anton Gladkyy


        2010/2/8 Václav Šmilauer <eudoxos@xxxxxxxx
        <mailto:eudoxos@xxxxxxxx>>


            > I'm trying to create a wrapper for Eigen library.
            Almost all functions
            > are already done, but I have a problem with
            EigenDecomposition.
            > Unfortunately, I never used this functions and do not
            clearly
            > understand what it implements.
            >
            > Eigen library has an EigenSolver for those tasks
            >
            http://eigen.tuxfamily.org/dox/classEigen_1_1EigenSolver.html
            >
            > I tried to get all values from that solver and
            compare with results,
            > what Wm3 gives. It is completely different.

            Using numpy to check:

            >>> from numpy import array
            >>> fron numpy.linalg import eig
            >>>
            a=array([[-26.8141,20.0536,-37.6382],[-17.0536,-4.37217,13.4546],[39.0891,8.34892,-21.6416]])
            >>> eig(a)                                 ##
            http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eig.html

            (array([-26.91393664+42.69101435j,
            -26.91393664-42.69101435j, 1.00000328 +0.j        ]),
array([[ 0.70584131+0.j , 0.70584131+0.j , 0.05979811+0.j ],
                  [-0.03806586+0.30786093j,
            -0.03806586-0.30786093j,  0.89863520+0.j        ],
                  [-0.01840919-0.63657033j,
            -0.01840919+0.63657033j,  0.43460208+0.j        ]]))

            You've picked matrix that has complex eigenvalues, that
            explains the
            difference; eigen handles it just fine ("(real,imag)"
            notation), wm3
            obviously doesn't. (the matrix should be positive
            definite for real
            eigenvalues, iirc).

            Cheers, Vaclav

            > Here is the source matrix:
            >
            > -26.8141 20.0536 -37.6382
            > -17.0536 -4.37217 13.4546
            > 39.0891 8.34892 -21.6416
            >
            >
            > =======================
            > wm3:
            >
            > tRot:
            > -0.698021 -0.0166317 -0.715884
            > 0.339616 -0.887829 -0.310515
            > -0.630419 -0.459872 0.625372
            >
            > tDiag:
            > -70.5641 0 0
            > 0 2.97262 0
            > 0 0 14.7635
            >
            > =======================
            > Eigen:
            > eigenvalues:
            > (-26.91,42.69)
            > (-26.91,-42.69)
            > (1,0)
            >
            > eigenvectors:
            > (-0.117,0.6961) (-0.117,-0.6961) (0.0598,0)
            > (-0.2973,-0.08855) (-0.2973,0.08855) (0.8986,0)
            > (0.6308,0.08732) (0.6308,-0.08732) (0.4346,0)
            >
            > pseudoEigenvalueMatrix:
            > -26.91 42.69 0
            > -42.69 -26.91 0
            > 0 0 1
            >
            > pseudoEigenvectors:
            > -0.1654 0.9844 0.0598
            > -0.4204 -0.1252 0.8986
            > 0.8921 0.1235 0.4346
            > =======================


            _______________________________________________
            Mailing list: https://launchpad.net/~yade-dev
            <https://launchpad.net/%7Eyade-dev>
            Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
            <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx>
            Unsubscribe : https://launchpad.net/~yade-dev
            <https://launchpad.net/%7Eyade-dev>
            More help   : https://help.launchpad.net/ListHelp


        _______________________________________________
        Mailing list: https://launchpad.net/~yade-dev
        <https://launchpad.net/%7Eyade-dev>
        Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
        <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx>
        Unsubscribe : https://launchpad.net/~yade-dev
        <https://launchpad.net/%7Eyade-dev>
        More help   : https://help.launchpad.net/ListHelp


        _______________________________________________
        Mailing list: https://launchpad.net/~yade-dev
        <https://launchpad.net/%7Eyade-dev>
        Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
        <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx>
        Unsubscribe : https://launchpad.net/~yade-dev
        <https://launchpad.net/%7Eyade-dev>
        More help   : https://help.launchpad.net/ListHelp


    _______________________________________________
    Mailing list: https://launchpad.net/~yade-dev
    <https://launchpad.net/%7Eyade-dev>
    Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
    <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx>
    Unsubscribe : https://launchpad.net/~yade-dev
    <https://launchpad.net/%7Eyade-dev>
    More help   : https://help.launchpad.net/ListHelp


    _______________________________________________
    Mailing list: https://launchpad.net/~yade-dev
    <https://launchpad.net/%7Eyade-dev>
    Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
    <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx>
    Unsubscribe : https://launchpad.net/~yade-dev
    <https://launchpad.net/%7Eyade-dev>
    More help   : https://help.launchpad.net/ListHelp


_______________________________________________
Mailing list: https://launchpad.net/~yade-dev <https://launchpad.net/%7Eyade-dev> Post to : yade-dev@xxxxxxxxxxxxxxxxxxx <mailto:yade-dev@xxxxxxxxxxxxxxxxxxx> Unsubscribe : https://launchpad.net/~yade-dev <https://launchpad.net/%7Eyade-dev>
More help   : https://help.launchpad.net/ListHelp

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

_______________________________________________
Mailing list: https://launchpad.net/~yade-dev
Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
Unsubscribe : https://launchpad.net/~yade-dev
More help   : https://help.launchpad.net/ListHelp


--
_______________
Bruno Chareyre
Associate Professor
Grenoble INP
Lab. 3SR
BP 53 - 38041, Grenoble cedex 9 - France
Tél : 33 4 56 52 86 21
Fax : 33 4 76 82 70 43
________________




Follow ups

References