by
In this scientific computing, numerical analysis, and control engineering tutorial, we will learn how to compute eigenvalues and eigenvectors of matrices in C++ by using the Eigen matrix library. The YouTube video accompanying this post is given below.
Eigenvector and Eigenvalues in Eigen C++ Matrix Library
The first step is to install the Eigen matrix library. If you are using the VS Code editor, you can follow this tutorial to do that. On the other hand, if you are using Microsoft Visual Studio, then follow this tutorial. The tutorial presented here is based on the VS Code editor. The GitHub page with all the code files is given here. The first step is to include the necessary libraries and namespaces:
#include<iostream>
#include<vector>
#include <iterator>
#include<Eigen/Dense>
#include <Eigen/Eigenvalues>
using namespace std;
using namespace Eigen;
The include statement
#include<Eigen/Dense>
is used to include a set of header files that are used for dense matrix computations. On the other hand, the include statement
#include <Eigen/Eigenvalues>
is used to include the header files for algorithms that are used for eigenvalue computations. The code presented below should be placed inside the “main()” C++ function.
First, define a random matrix whose eigenvalues and eigenvectors we will compute. Also, we define complex matrices for storing the eigenvalues and eigenvectors.
// create a random matrix
MatrixXd randomMatrix;
randomMatrix.setRandom(3,3);
// this vector/matrix will be used to store the eigenvalues
// The matrix has to be complex
// This is the typedef for MatrixXcd
// typedef Matrix< std::complex<double> , Dynamic , Dynamic > MatrixXcd
// That is, this matrix can store std::complex<double> data type
MatrixXcd eigenvalueMatrix;
// this matrix will be used to store the eigenvectors
// the matrix needs to be complex since the eigenvalues can be complex
MatrixXcd eigenvectorMatrix;
// print the random matrix
cout <<endl<< randomMatrix << endl;
The data type “MatrixXd” is a type definition for a matrix of double real entries whose dimensions can be changed during the compile time. On the other hand, “MatrixXcd” is a type definition for a matrix of double complex entries whose dimensions can be changed during the compile time. The complex matrix “eigenvalueMatrix” is used to store the complex eigenvalues. The matrix “eigenvectorMatrix” is used to store the eigenvectors of the matrix “randomMatrix” column-wise.
Next, we compute the eigenvalues and eigenvectors
EigenSolver<MatrixXd> eigenValueSolver(randomMatrix);
eigenvalueMatrix=eigenValueSolver.eigenvalues();
// the eigenvectors are stored column wise
eigenvectorMatrix=eigenValueSolver.eigenvectors();
// print the eigenvalues
cout <<endl<<"Computed eigenvalues:"<<endl<<eigenvalueMatrix<<endl;
// print the eigenvectors
cout <<endl<<"Computed eigenvectors:"<<endl<<eigenvectorMatrix<<endl;
// extract real and imaginary parts of the eigenvalues
MatrixXd realPartEigenValues;
MatrixXd imagPartEigenValues;
realPartEigenValues=eigenvalueMatrix.real();
imagPartEigenValues=eigenvalueMatrix.imag();
// print the real parts of the eigenvalues
cout<<endl<<realPartEigenValues<<endl;
// print the complex parts of the eigenvalues
cout<<endl<<imagPartEigenValues<<endl;
The following code line
EigenSolver<MatrixXd> eigenValueSolver(randomMatrix);
is used to create the EigenSolver object. We pass to the constructor the matrix “randomMatrix”. The eigenvalues and eigenvectors are computed as follows
eigenvalueMatrix=eigenValueSolver.eigenvalues();
eigenvectorMatrix=eigenValueSolver.eigenvectors();
The eigenvectors are stored column-wise in the matrix “eigenvectorMatrix”. We extract the real and imaginary parts of the computed eigenvalues as follows
realPartEigenValues=eigenvalueMatrix.real();
imagPartEigenValues=eigenvalueMatrix.imag();
The following piece of code will double-check the result. That is, it will pick one computed eigenvalue and eigenvector and it will check if the defining equation is satisfied.
// let us double check one eigenvalue and eigenvector
MatrixXcd testLeft;
MatrixXcd testRight;
// definition of the eigenvalue and the eigenvector
// test1Right = A*v
// test2Right = lambda*v
// A*v=lambda*v
// if everything is computed correctly
// testLeft = testRight
testLeft=randomMatrix*eigenvectorMatrix.col(0);
testRight=eigenvalueMatrix(0,0)*eigenvectorMatrix.col(0);
cout<<endl<<"Left side vector:"<<endl<<testLeft<<endl;
cout<<endl<<"Right side vector:"<<endl<<testRight<<endl;
cout<<endl<<"Left minus right side vector:"<<endl<<testLeft-testRight<<endl;
This code will test if this equation
(1) 
is satisfied, where 
is the random matrix, 
is the computed eigenvector, and 
is the computed eigenvalue. We test is the left-hand side of this equation is equal to the right-hand side. That is if 
equals 
. We also test if the left-hand side of this expression:
(2) 
is equal to zero. By running this piece of code, if we can see that the eigenvalues and eigenvectors are correctly computed. Next, we explain how to perform basic operations on computed eigenvalues. The following piece of code will double the computed eigenvalues and it will store them in C++ vector STL object.
// let us do some manipulations with eigenvalues
// for example let us define a vector structure of eigenvalues
// we need std:: complex<double>
std::vector<std::complex<double>> doubledEigenValueVector;
// eigenvalues are of the type complex <double>
complex<double> eigenValue;
for (int i=0; i<eigenvalueMatrix.rows(); i++)
{
eigenValue=eigenvalueMatrix(i,0);
cout<<endl<<" Real part of the eigenvalue "<<i<<" is "<<eigenValue.real()<<endl;
cout<<" Imaginary part of the eigenvalue "<<i<<" is "<<eigenValue.imag()<<endl;
// double the eigenvalues
doubledEigenValueVector.push_back(2.0*eigenValue);
}
// print the results
for (auto it=doubledEigenValueVector.begin();
it !=doubledEigenValueVector.end();
it++)
{
cout<<endl<<"Original eigenvalue: "<<eigenvalueMatrix(distance(doubledEigenValueVector.begin(), it),0)<<endl;
cout<<"Doubled eigenvalue: "<< *it <<endl;
}
The above code is self explanatory. Here, we should just say that eigenvalues as complex numbers are stored by using the STL “complex <double>” data type.
