Introduction to Eigen C++ Matrix Library

Eigen linear algebra library is a powerful C++ library for performing matrix-vector and linear algebra computations. This library can be used for the design and implementation of model-based controllers, as well as other algorithms, such as machine learning and signal processing algorithms. In this tutorial, I give an intro to the Eigen library. A quick list of Eigen commands that are equivalent to MATLAB commands can be found here. The video lecture accompanying this post is given below.

1.Installation

The installation of Eigen is relatively simple. Eigen is a template library defined using header source files. Download the source files and add these files to the C++ compiler path. In Microsoft Visual Studio this is relatively simple. To add Eigen to the path of your project in MS Visual Studio follow these steps. Click on “Project” menu and under this menu, click on “ConsolApplicationX Properties ” (in my case, the name of the project is “ConsolApplicationX”, however, in your case, this might be different). Then, expand “Configuration Properties”, and then expand “C/C++” option. Under the “C/C++” option, click on “General” and click on “Additional Include Directories”. You can add the path of your Eigen (do not forget to “unzip” your folder) folder to the “Additional Include Directories”. That is all. You are ready to use Eigen. Also, in the codes we are going to explain in this tutorials, the following code needs to be added before the main() function:

#include <iostream>
#include<Eigen/Dense>

using namespace std;
using namespace Eigen;

We include the Eigen library header files by “#include<Eigen/Dense>” and we use “using namespace Eigen” to simplify the notation, otherwise in the front of every command, you will have to use “Eigen::”

2. Defining/Declaring Matrices

Let us start with matrix definitions/declarations. There are many ways for defining matrices in Eigen. The following code line defines a 3×3 matrix of floats, sets its entries to zero and prints out the matrix:

	// define 3x3 matrix of floats and set its entries to zero -explicit declaration
	Matrix <float, 3, 3> matrixA;
	matrixA.setZero();
	cout << matrixA <<endl;

In code line 2, the first parameter is the variable type. Change this to double, int, or any other standard variable type. The second two parameters are matrix dimensions. The code line 3 is self-explanatory. Finally, in code line 4, we print the matrix. Eigen also provides built-in type definitions. Using the built-in typedef variables, we can define the same matrix as follows

	//define 3x3 matrix of floats and set its entries to zero -typedef declaration
	Matrix3f matrixA1;
	matrixA1.setZero();
	cout <<"\n"<<matrixA1<<endl;

The extension “3f” in “Matrix3f” stands for “3 x 3” matrix of floats. If we would like to define a 4×4 matrix of integers, the typedef notation would be Matrix4i.
This type of matrix definition might work well when the matrix dimensions are known at compile time. However, in many cases, we would like to have dynamic matrix dimensions or to set them during the compile time. Also, the fixed-size declaration at the compile time might be efficient for small-scale matrices, however, for larger problems and matrices we need to use “dynamic matrices”. The following code line defines a dynamic matrix, whose dimensions can be set during the run time:

	// define a dynamic matrix, explicit declaration
	Matrix <float, Dynamic, Dynamic> matrixB; 

This code line defines a dynamic matrix of size 0, whose array entries have not been allocated at the compile time. We can also use the built-in Eigen type declaration for dynamics matrices:

//define a dynamic matrix, typedef declaration
MatrixXf matrixB1

where MatrixXf is a typedef declaration in which character “X” stands for unknown size and $f$ stands for float. Constructor that allocates the memory for matrix entries but do not initialize them is given by

// constructor
	MatrixXd matrixC(10, 10);

2. Assigning/accessing and printing matrix/vector entries

The entries of a matrix can be accessed and printed as follows

	Matrix4f matrixD;
	matrixD << 1,   2,  3,  4,
 	    	5,   6,  7,  8,
		    9,  10, 11, 12,
		    13, 14, 15, 16;
	cout << endl<<  matrixD << endl;

Notice that matrix entries are indexed from 0. That is, the entry on the first row and first column of the matrix “matrixC1” is denoted by “matrixC1(0,0)”. The matrix is entered using comma-separated values. The matrix is printed using the classical “cout<<” framework. The following code defines, a dynamic matrix, resizes it to a 3×3 matrix and assigns its entries to zero

	// define a dynamic matrix, resize it to a 3x3 matrix, 
	// and set its entries to zero, and print the matrix
	Matrix3f matrixD1;
	matrixD1.resize(3, 3);
	matrixD1.setZero(3, 3);
	cout <<endl<< matrixD1;

Defining/declaring zero matrices, matrices of ones, and matrices consisting of constants:

// setting matrix entries - two approaches
	int rowsNumber = 5;
	int columnNumber= 5;

	// matrix of zeros
	MatrixXf matrix1zeros;
	matrix1zeros = MatrixXf::Zero(rowsNumber, columnNumber);
	cout << "\n \n"<< matrix1zeros<<endl;
	// another option:
	MatrixXf matrix1zeros1;
	matrix1zeros1.setZero(rowsNumber, columnNumber);
	cout << "\n \n" << matrix1zeros1 << endl;

	//matrix of ones
	MatrixXf matrix1ones;
	matrix1ones = MatrixXf::Ones(rowsNumber, columnNumber);
	cout << "\n \n" << matrix1ones << endl;
	//another option
	MatrixXf matrix1ones1;
	matrix1ones1.setOnes(rowsNumber, columnNumber);
	cout << "\n \n" << matrix1ones1 << endl;

	//matrix of constants
	float value = 1.1;
	MatrixXf matrix1const;
	matrix1const = MatrixXf::Constant(rowsNumber, columnNumber,value);
	cout << "\n \n" << matrix1const << endl;
	//another option
	MatrixXf matrix1const1;
	matrix1const1.setConstant(rowsNumber, columnNumber, value);
	cout << "\n \n" << matrix1const1 << endl;

Definition/declaration of identity matrices:

//identity matrix, two approaches

	rowNumber = 10;
	columnNumber = 10;
	
	MatrixXd matrixIdentity;
	matrixIdentity = MatrixXd::Identity(rowNumber,columnNumber);
	cout << "\n \n" << matrixIdentity << endl;
	
	MatrixXd matrixIdentity1;
	matrixIdentity1.setIdentity(rowNumber, columnNumber);
	cout << "\n \n" << matrixIdentity1 << endl;

Matrix blocks, can be accessed using the following notation: “matrix_object_name”.block(i,j,rows,cols) or “matrix_object_name”.block<rows,cols>(i,j). This notation selects the specified numbers of rows and columns of the matrix “matrix_object_name” starting from the entry (i,j). Note that in contrast to the MATLAB notation, the entries counted from 0. That is, the (0,0) corresponds to the (1,1) entry of a matrix.

//accessing matrix blocks
	MatrixXd matrixV(4,4);
	matrixV << 101, 102, 103, 104,
		105, 106, 107, 108,
		109, 110, 111, 112,
		113, 114, 115, 116;
	//access the matrix composed of 1:2 rows and 1:2 columns of matrixV
	MatrixXd matrixVpartition = matrixV.block(0, 0, 2, 2);
	cout << "\n \n" << matrixVpartition << endl;
	
	MatrixXd matrixVpartition2 = matrixV.block(1,1, 2, 2);
	cout << "\n \n" << matrixVpartition2 << endl;

Columns and rows of a matrix can be accessed using the following notation:

	//accessing columns and rows of a matrix

	cout<<"\n\n"<<"Row 1 of matrixV is \n "<<matrixV.row(0);
	cout << "\n\n" << "Column 1 of matrixV is \n" << matrixV.col(0);

We construct a diagonal matrix using this approach

//create a diagonal matrix out of a vector
	Matrix <double, 3, 1> vector1;
	vector1 << 1, 2, 3;
	MatrixXd diagonalMatrix;
	diagonalMatrix = vector1.asDiagonal();
	cout << " Diagonal matrix is \n\n" << diagonalMatrix;

3. Basic matrix operations

Here we explain how to perform matrix addition/subtraction, matrix multiplication, matrix transpose, and matrix inverse. Matrix addition can be performed using the following approach:

// basic matrix operations
	//matrix addition
	MatrixXd A1(2, 2);
	MatrixXd B1(2, 2);

	A1 << 1, 2,
		3, 4;
	B1 << 3, 4,
		5, 6;
	MatrixXd C1 = A1 + B1;
	cout << " \n\n The sum of A1 and B1 is\n\n" << C1 << endl;
	// similarly you can subtract A1 and B1

Matrix multiplication and matrix-scalar multiplication can be performed using the following approach:

//matrix multiplication
	MatrixXd D1 = A1*B1;
	cout << " \n\n The matrix product of A1 and B1 is\n\n" << D1 << endl;

	//multiplication by a scalar
	int s1 = 2;
	MatrixXd F1;
	F1 = s1 * A1;
	cout << " \n\n The matrix product of the scalar 2 and the matrix A1 is\n\n" << F1 << endl;

Matrix transposes can be computed using the following approaches

//matrix transpose
	MatrixXd At;
	MatrixXd R1;
	// we can obain a transpose of A1 as follows
	At = A1.transpose();
	cout << "\n\n Original matrix A1\n\n" << A1;
	cout << "\n\n Its transpose\n\n " << At;

	// we can use a transpose operator in expressions
	R1 = A1.transpose() + B1;
	cout << "\n\n Matrix R1=A1^{T}+B1 is\n\n" << R1;

	// here we should be careful, the operation .transpose() might lead to unexpected results in this scenarios
    //	A1 = A1.transpose();
	// cout << " \n\n This should be a transpose of the matrix A1\n\n" << A1 << endl;

	// the correct and safe way to do the matrix transpose is the following
	 A1.transposeInPlace();
	cout << " \n\n This should be a transpose of the matrix A1\n\n" << A1 << endl;
	//restore A1 matrix to its original state
	A1.transposeInPlace();

The expression A1.transpose() returns a matrix object without performing the actual transpose on A1. This is perfectly fine as long as we do not write A1=A1.transpose(). This expression will not transpose matrix A, instead, the in-place matrix transpose A1.transposeInPlace() should be used to transpose the matrix A directly.

Matrix inverses can be computed using the following approach

	MatrixXd G1;
	G1 = A1.inverse();
	cout << "\n\n The inverse of the matrix A1 is\n\n" << G1;
	cout << "\n\n Double check A1^{-1}*A1=\n\n" << G1*A1;

That will be all in this introduction. In next posts we will explain how to use Eigen to simulate a state-space model and how to compute the Kalman filter.