Eigen common operations

        Eigen It's a high-level C ++ library , Support linear algebra , Matrix and vector operations , Numerical analysis and its related algorithms . This routine contains header files ：

#include <iostream>

using namespace std;

#include <ctime>
// Eigen  The core part of the 
#include <Eigen/Core>
#include <Eigen/Dense>
using namespace Eigen;

#define MATRIX_SIZE 5

1、 data structure

        Eigen The data structure is as follows , Single precision and double precision conversion d Replace with f.

         Rotation matrix （3*3）：Eigen::Matrix3d.

         Rotating vector （3*1）：Eigen::AngleAxisd.

         Euler Angle （3*1）：Eigen::Vector3d.

         Four elements （4*1）：Eigen::Quaterniond.

         Euclidean transformation matrix （4*4）：Eigen::Isometry3d.

         Affine transformation （4*4）：Eigen::Affine3d.

         Projective transformation （4*4）：Eigen::Projective3d.

2、 Variable declarations

        Eigen All vectors and matrices in are Eigen::Matrix, It's a template class . Its first three parameters are ： data type , That's ok , Column .

        1） Make a statement 2*3 Of float matrix

Matrix<float, 2, 3> matrix_23;

        2） Built in type

        Eigen adopt typedef There are many built-in types , But the bottom is still Eigen::Matrix.

Vector3d v_3d; //  Equivalent to Eigen::Matrix<double, 3, 1>, Three dimensional vector 
Matrix3d matrix_33; //  Equivalent to Eigen::Matrix<double, 3, 3>, Three times three matrix 

         If you're not sure about the size of the matrix , Dynamic size matrices can be used ：

Matrix<double, Dynamic, Dynamic> matrix_dynamic;
MatrixXd matrix_x; //  More simple point 

3、 initialization

        Eigen All vectors and matrices in are Eigen::Matrix, It's a template class . Its first three parameters are ： data type , That's ok , Column .

Matrix3d::Zero(); // doubel type ,3*3 All zero matrix 
Matrix3d::Random(); // doubel type ,3*3 Random number matrix 
Matrix3d::Identity(); // doubel type ,3*3 Unit matrix E
Vector3d(0, 0, 1);  // doubel type ,3*1 Column vector 

         The matrix is given initial value ：

Matrix3d matrix_33 = Matrix3d::Zero(); // Initialize to 3*3 All zero matrix 

         assignment , Show ：

//  input data （ initialization ）
Matrix<float, 2, 3> matrix_23;
matrix_23 << 1, 2, 3, 4, 5, 6;


//  Output 
cout << "matrix 2x3 from 1 to 6: \n" << matrix_23 << endl;

 //  use () Access the elements in the matrix 
cout << "print matrix 2x3: " << endl;
for (int i = 0; i < 2; i++)
{
    for (int j = 0; j < 3; j++) cout << matrix_23(i, j) << "\t";
    cout << endl;
}

         result ：

4、 Matrix operations

         Defining variables ：

Matrix<float, 2, 3> matrix_23; // 2*3 matrix ,float
matrix_23 << 1, 2, 3, 4, 5, 6;

Vector3d v_3d; // 3*1 vector ,double
v_3d << 3, 2, 1;

Matrix<float, 3, 1> vd_3d; // 3*1 vector ,float
vd_3d << 4, 5, 6;

         stay Eigen You can't mix two different types of matrices , This syntax is wrong .

Matrix<double, 2, 1> result_wrong_type = matrix_23 * v_3d;

         Again, you can't mistake the dimension of the matrix , This syntax is wrong .        

Eigen::Matrix<double, 2, 3> result_wrong_dimension = matrix_23.cast<double>() * v_3d;

         Matrix display call ：

// 2_3*3_1 = 2_1
Matrix<double, 2, 1> result = matrix_23.cast<double>() * v_3d;
cout << "[1,2,3;4,5,6]*[3,2,1]=" << result.transpose() << endl;

Matrix<float, 2, 1> result2 = matrix_23 * vd_3d;
cout << "[1,2,3;4,5,6]*[4,5,6]: " << result2.transpose() << endl;

5、 Matrix operations

         Some matrix operations , Four operations will not be demonstrated , Direct use +-*/ that will do .

Matrix<float, 3, 3> matrix_33;      //  Random number matrix 
matrix_33 << 1, 5, 2, 2, 4, 1, 3, 6, 2;
	
cout << "random matrix: \n" << matrix_33 << endl;
cout << "transpose: \n" << matrix_33.transpose() << endl;      //  Transposition 
cout << "sum: " << matrix_33.sum() << endl;            //  The elements and 
cout << "trace: " << matrix_33.trace() << endl;          //  trace 
cout << "times 10: \n" << 10 * matrix_33 << endl;               //  Number multiplication 
cout << "inverse: \n" << matrix_33.inverse() << endl;        //  The inverse 
cout << "det: " << matrix_33.determinant() << endl;    //  determinant 

         result ：

6、 The eigenvalue , Diagonalization

         Diagonalization , Solve eigenvalues .

Matrix<double, 3, 3> matrix_33;
matrix_33 <<1, 5, 2, 2, 4, 1, 3, 6, 2;

//  Diagonalization 
SelfAdjointEigenSolver<Matrix3d> eigen_solver(matrix_33.transpose() * matrix_33);

//  Solve eigenvalues 
cout << "Eigen values = \n" << eigen_solver.eigenvalues() << endl;
cout << "Eigen vectors = \n" << eigen_solver.eigenvectors() << endl;

         result ：

7、 solve equations

         solve matrix_NN * x = v_Nd This equation ,N The size of is defined in the previous macro , here N=50, It is generated from random numbers .

Matrix<double, MATRIX_SIZE, MATRIX_SIZE> matrix_NN
 	= MatrixXd::Random(MATRIX_SIZE, MATRIX_SIZE);
matrix_NN = matrix_NN * matrix_NN.transpose();  //  Guaranteed positive semidefinite 
Matrix<double, MATRIX_SIZE, 1> v_Nd = MatrixXd::Random(MATRIX_SIZE, 1);

        1） Direct inversion is naturally the most direct , But the amount of inversion is large .

clock_t time_stt = clock(); //  timing 
//  Direct inversion 
Matrix<double, MATRIX_SIZE, 1> x = matrix_NN.inverse() * v_Nd; 
cout << "time of normal inverse is "
	<< 1000 * (clock() - time_stt) / (double)CLOCKS_PER_SEC << "ms" << endl;
cout << "x = " << x.transpose() << endl;

        2） Matrix decomposition is usually used to find , for example QR decompose , It's going to be a lot faster .

time_stt = clock();
x = matrix_NN.colPivHouseholderQr().solve(v_Nd);
cout << "time of Qr decomposition is "
	<< 1000 * (clock() - time_stt) / (double)CLOCKS_PER_SEC << "ms" << endl;
cout << "x = " << x.transpose() << endl;

        3） For positive definite matrices , You can also use cholesky Decompose to solve the equation .

time_stt = clock();
x = matrix_NN.ldlt().solve(v_Nd);
cout << "time of ldlt decomposition is "
<< 1000 * (clock() - time_stt) / (double)CLOCKS_PER_SEC << "ms" << endl;
cout << "x = " << x.transpose() << endl;

         result ：

x = -0.177629 0.179189 7.66884e-05 -0.493822 0.568221