I'm trying to convert a piece of MATLAB code to C++, which does:

Code: Select all

  [V,D] = eig(A,B)


According to http://www.mathworks.nl/help/techdoc/ref/eig.html

[V,D] = eig(A,B) produces a diagonal matrix D of generalized eigenvalues and a full matrix V whose columns are the corresponding eigenvectors so that A*V = B*V*D.

In my case, 'A' and 'B' are both symmetric, 6 by 6 matrices. 'B' is a sparse matrix: it only has three non-zero elements. How would you do this in C++, using the Eigen library? Does GeneralizedSelfAdjointEigenSolver do the job? As in the following C++ code?

Code: Select all

  typedef Eigen::Matrix<double, 6, 6> Matrix6d;
  typedef Eigen::Matrix<double, 6, 1> Vector6d;

  // [V,D] = eig(A,B)
  std::pair<Matrix6d, Vector6d> eig(const Matrix6d& A, const Matrix6d& B)
  {
    const Eigen::GeneralizedSelfAdjointEigenSolver<Matrix6d> solver(A, B);
         
    const Matrix6d V = solver.eigenvectors();
    const Vector6d D = solver.eigenvalues();

    return std::pair<Matrix6d, Vector6d>(V, D);
  }


I'm asking, because it looks like my code, using Eigen::GeneralizedSelfAdjointEigenSolve, produces more undefined results (NaN values) than the orginal MATLAB code.

Kind regards,

Niels
--
Niels Dekker
http://www.xs4all.nl/~nd/dekkerware
Scientific programmer at LKEB, Leiden University Medical Center

If B has only 3 non zeros, this means B is not SPD, so you cannot use GeneralizedSelfAdjointEigenSolver.
Actually, I'm wondering whether your problem could be advantageously reduced to a 3x3 eigen value problem.

ggael wrote:If B has only 3 non zeros, this means B is not SPD, so you cannot use GeneralizedSelfAdjointEigenSolver.
Actually, I'm wondering whether your problem could be advantageously reduced to a 3x3 eigen value problem.

Thanks for your reply, Gaël. The MATLAB code I'm trying to convert to C++ is at http://homepages.inf.ed.ac.uk/rbf/CVonl ... fitellip.m (by Robert B. Fisher). It says:

Code: Select all

   % Build 6x6 constraint matrix
   C(6,6) = 0; C(1,3) = -2; C(2,2) = 1; C(3,1) = -2;

   % Solve eigensystem
   [gevec, geval] = eig(S,C);


So the sparse matrix 'B' in my previous mail is actually named 'C' in fitellip.m Note that MATLAB starts index numbering at 1.

Do you think this can be translated to C++, using the Eigen library?

Hi

I have 5x5 matrices P and Q.

P:

92,3135999999999 11,7504000000000 -15,8464000000000 -21,8880000000000 -0,831999999999999
11,7504000000000 15,7696000000000 -4,37760000000000 -0,831999999999999 2,11200000000000
-15,8464000000000 -4,37760000000000 4,24960000000000 2,11200000000000 -1,15200000000000
-21,8880000000000 -0,831999999999999 2,11200000000000 15,0400000000000 -1,44000000000000
-0,831999999999999 2,11200000000000 -1,15200000000000 -1,44000000000000 2,24000000000000

Q:

60,1600000000000 -2,88000000000000 0 0 0
-2,88000000000000 17,2800000000000 -2,88000000000000 0 0
0 -2,88000000000000 8,96000000000000 0 0
0 0 0 1 0
0 0 0 0 1

When I run [V,D] = eig(P,Q); code in MATLAB, I got V matrix like this:

-0,0681677353936515 0,0374198822577883 0,0982046018281875 -0,0205744132768129 0,0254779498548796
-7,37412612061575e-17 -0,208265106287322 0,120417010536798 0,0611738781611838 0,00700447259303914
-0,272670941574606 -0,170780968664988 -0,115928660970322 -0,0272985869476333 -0,0166730304789956
-0,0818012824723817 0,0928465025089220 0,176927171720780 0,0792551550611296 -0,973192413380240
-0,218136753259684 0,237242645867653 -0,0723166534820224 0,938104911994981 0,104219765779281

I wrote and tested Niels_Dekker's code with the same input matrices P and Q then the result is:

-0,0681677 -0,0374199 0,0982046 -0,0205744 0,0254779
5.8993e-16 0,208265 0,120417 0,0611739 0,00700447
-0,272671 0,170781 -0,115929 -0,0272986 -0,016673
-0,0818013 -0,0928465 0,176927 0,0792552 -0,973192
-0,218137 -0,237243 -0,0723167 0,938105 0,10422

As you see that the values of the second column is different from MATLAB's result. Why it has negative values?

Halil

eigenvectors are not signed since they only are supposed to satisfy A * v = l * B * v with v'*B*v=1. If v satisfies these equations, then so does -v.

Sorry my question was wrong: It must be Why the second column of v is different from MATLAB's result? The difference is the second column's multiplication by -1 when compared to MATLAB's result.

That's what I said, both answers (MatLab and Eigen) are equivalent, and there is no reason to prefer one from the other.

Ok, I see. Thanks.