Why am I getting wrong eigenvectors using eig().

My application is the common one of multimodal vibration: set up a mass M and stiffness K matrix using convenient generalised coordinates, then get the eigenvalues and eigenvectors of K*M^-1.

This works fine normally, but it gives me wrong eigenvectors when used on the standard example of a massive block (usually a car body) mounted on two springs and using the simplest generalised coordinates: vertical displacement of the centre of mass and angle of rotation.

The matrices I've used are: M=[70 0;0 10] K=[300 40;40 22]*1000 Matlab gives correct eigenvalues:[5079 0;0 1406], but the wrong eigenvectors:[1 -1;0.1984 0.7199].

The correct eigenvectors are [0.1984 0.7199;-1 1], which I got using Mathematica and confirmed by using them to give me diagonal generalised mass and stiffness matrices, as they should.

You'll see that the Matlab eigenvectors, which form the modal matrix, have got the correct numbers in them, but they are mixed up, both in position and sign.

Can anybody suggest a reason, and better still a solution?

[Matt J]

Does Mathematica put the eigenvectors in the rows of the matrix or in the columns? Matlab puts them in the columns.

[Matt Kindig]

Also keep in mind that [1 -1] and [-1 1] are the same eigenvector, since eigenvectors are only defined up to some scale factor. Multiplying all of the elements of an eigenvector by a scalar (such as -1) does not change the eigenvector.

Answer 1

[Matt J]

In any case, you can directly verify that MATLAB did not give you "wrong" eigenvectors. The eigenvector matrix, V, which MATLAB gave you, does satisfy the appropriate eigenvalue equations

>> M=[70 0;0 10]; K=[300 40;40 22]*1000;

>> [V,D]=eig(K*inv(M));

>> (K*inv(M))*V- V*D

ans =

     1.0e-12 *

           0    0.2274
     -0.2274    0.1137

[Jack]

Matt... Thanks for the very fast response.

I have thought about this over the weekend. Clearly you are right about getting correct eigenvectors (I should have known better than to doubt Matlab!).

I've now realised that the two forms are the same eigenvectors rotated through pi/2. But of course this rotation changes the mode shapes, so must be physically wrong.

The question now is, why do I get the correctly diagonalised generalised mass and stiffness matrices when I use this rotated modal matrix? I think it may be something to do with my choice of generalised coordinates. But that's not a Matlab problem. I'll have to think some more.

[Heitor Camarini]

Hello Jack,

I am having that exact same problem now, almost 8 years late!! The eigenvectors I am getting from MATLAB do not mach the ones from my HP calculator. Have you been able to solve it? I am trying to calculate multimodal vibrations as well.

Thanks!

[Matt J]

Jack said he solved the problem:

https://www.mathworks.com/matlabcentral/answers/50647-why-am-i-getting-wrong-eigenvectors-using-eig#comment_570410

Answer 2

[Huy Dinh]

Jack sorry I'm 6 years late with this answer lol. I've been having this problem also for the past 2 days and it drove me crazy. I finally found the reason! The eigenvalue problem you were doing is for the matrix

A1 = K*inv(M).

If you do the eigenvalue problem for

A2 = inv(M)*K

instead, you'll get what you wanted! Note that A1 is the transpose of A2, which, in 2D, happens to be the angle transformation using

T =  [cos(90) sin(90)
      sin(90) cos(90)]

where

A1 = A2 * T

For that reason, your eigenvalues are "rotated" 90°. I hope that makes sense

[Jack]

Huy Dinh...

Thanks for your note after all this time. Yes, I found this solution at the time in conversation with a mathematician who took pleasure in pointing out the poverty of my maths: that [A]*[B] is not the same as [B]*[A] !

I should have posted it here then.

[Syed Abdul Rafay]

clear all
clc
syms lambda
syms x
syms i
rhoEpoxy = 1200;
Em = 3e+09;
R=5;
E_L=70e+9;
E_r=200e+9;
rho_L=2702;
rho_r=5700;
nu=3;
A_K=zeros(R+1,R+1);
A_M=zeros(R+1,R+1);
G=zeros(R+1,R+1);
H1_K=zeros(R+1,R+1);
H1_M=zeros(R+1,R+1);
H2_M=zeros(R+1,R+1);
for ri=1:R+1
    r = ri-1;
    for mi=1:R+1
        m = mi-1;
        fun1 = @(ksei1) (ksei1.^(r+m)).*((1-(exp(nu*ksei1)-1)./(exp(nu)-1))+E_r/E_L*(exp(nu*ksei1)-1)./(exp(nu)-1));
        %G(mi,ri)=integral(fun1,0,1);
        fun_h1=@(ksei2,s2) (-2*nu/(exp(nu)-1)*(E_r/E_L-1)*exp(nu*s2)).*(s2.^r).*(ksei2.^m)+...
            (nu^2/(exp(nu)-1)*(E_r/E_L-1)*exp(nu*s2)).*(s2.^r).*(ksei2.^(m+1))-...
            (nu^2/(exp(nu)-1)*(E_r/E_L-1)*exp(nu*s2)).*(s2.^(r+1)).*(ksei2.^m);
        fun_h1_lambda=@(ksei3,s3) (-1/6*((ksei3-s3).^3).*((1-(exp(nu*s3)-1)./(exp(nu)-1))+...
            rho_r/rho_L*(exp(nu*s3)-1)./(exp(nu)-1))).*(s3.^r).*(ksei3.^m);
        fun_h2_lambda=@(ksei4,s4) (1/6*(ksei4.^3).*((1-(exp(nu*s4)-1)./(exp(nu)-1))+...
            rho_r/rho_L*(exp(nu*s4)-1)./(exp(nu)-1))-1/2*(ksei4.^2).*s4.*((1-(exp(nu*s4)-1)./(exp(nu)-1))+...
            rho_r/rho_L*(exp(nu*s4)-1)./(exp(nu)-1))).*(s4.^r).*(ksei4.^m);
        
        options = {'RelTol', 1e-22, 'AbsTol', 1e-24};
        G(mi,ri) = integral(fun1, 0, 1, options{:});
        H1_K(mi,ri) = integral2(fun_h1, 0, 1, 0, @(ksei2) ksei2, options{:});
        H1_M(mi,ri) = integral2(fun_h1_lambda, 0, 1, 0, @(ksei3) ksei3, options{:});
        H2_M(mi,ri) = integral2(fun_h2_lambda, 0, 1, 0, 1, options{:});
        A_K(mi,ri)=G(mi,ri)+H1_K(mi,ri);
        A_M(mi,ri)=H1_M(mi,ri)+H2_M(mi,ri);
    end
end
sort(sqrt(eig(A_K,-A_M)))
% fplot(@(x) rho_L*(1-(exp(nu*x)-1)./(exp(nu)-1))+rho_r*(exp(nu*x)-1)./(exp(nu)-1), [0,1])
% xval = fzero(@(lambda2) det(A_K+lambda2*A_M),0);
% sqrt(xval)

Can you help me also. I have kind of same problem. My first three modes are correct but 4th and 5th mode of the frequency are wrong. There is a big difference

these are the correct values

2.854

21.496

63.679

126.60

As I increase R I start getting complex values don't know why

The equation is correct and even the code. I think there is a problem with eig function