// 来自高翔SLAM十四讲
#include<iostream>
#include<ctime>
#include<eigen3/Eigen/Core>
// 稠密矩阵的代数运算(逆、特征值等)
#include<eigen3/Eigen/Dense>using namespace std;#define MATRIX_SIZE 5int main(){// eigen以矩阵为基本数据元素。它是一个模板类。它的前三个参数为：数据类型，行，列// 声明一个2*3的float矩阵Eigen::Matrix<float,2,3> matrix23;// 同时，eigen通过typedef提供了许多内置类型，不过底层仍旧是eigen::Matrix// 例如：Vector3d实质上是Eigen::Matrix<double,3,1>Eigen::Vector3d v_3d;// Matrix3d实质上是Eigen::Matrix<double,3,3>Eigen::Matrix3d matrix33 = Eigen::Matrix3d::Zero(); // 初始化为0// 如果大小不确定的矩阵，可以使用动态大小的矩阵Eigen::Matrix<double, Eigen::Dynamic, -1> matrix_dynamic;// 更简单的Eigen::MatrixXd matrix_x;// 还有很多类似的类型，自己看吧// 下面是对矩阵的操作// 输入数据matrix23 << 1,2,3,4,5,6; // 重载了 <<// 输出cout<<matrix23<<endl;// 用()访问矩阵中的元素for(int i=0; i<1; i++){for(int j=0; j<2; j++){cout<<matrix23(i,j)<<endl;}}v_3d << 3,2,1;// 矩阵与向量相乘(实际上仍然是矩阵和矩阵)// 不同类型的矩阵不能相乘// Eigen::Matrix<double, 2, 1> result = matrix23 * v_3d;// 应该显式转换Eigen::Matrix<double, 2, 1> result = matrix23.cast<double>() * v_3d;cout<< result <<endl;// 生成随机的矩阵matrix33 = Eigen::Matrix3d::Random();cout<< matrix33 <<endl;// 转置cout<< matrix33.transpose() <<endl;// 各个元素的和cout<< matrix33.sum() <<endl;// 迹cout<< matrix33.trace() <<endl;// 数乘cout<< 10*matrix33 <<endl;// 逆cout<< matrix33.inverse() <<endl;// 行列式cout<< matrix33.determinant() <<endl;// 特征值// 是对称矩阵可以保证对角化成功Eigen::SelfAdjointEigenSolver<Eigen::Matrix3d> eigen_solver(matrix33.transpose() * matrix33);cout<<" eigen value = "<<eigen_solver.eigenvalues() <<endl; // 特征值cout<<"eigen vector = "<<eigen_solver.eigenvectors() <<endl; // 特征向量// 解方程// 我们求解matrix_NN * x = v_Nd 这个方程// 直接求逆当然是最简单的，但是直接求逆的计算量特别大Eigen::Matrix<double, MATRIX_SIZE, MATRIX_SIZE>matrix_NN;matrix_NN = Eigen::MatrixXd::Random(MATRIX_SIZE,MATRIX_SIZE);Eigen::Matrix<double, MATRIX_SIZE, 1> v_Nd;v_Nd = Eigen::MatrixXd::Random(MATRIX_SIZE, 1);// 计时clock_t time_stt = clock();// 直接求逆Eigen::Matrix<double, MATRIX_SIZE, 1> x = matrix_NN.reverse() * v_Nd;cout<<" time use in normal inverse is "<< 1000 * (clock() - time_stt)/(double)CLOCKS_PER_SEC << "ms"<<endl;// 通常用矩阵分解来求，例如QR分解，速度会快很多,然而这边显示时间却慢很多time_stt = clock();x = matrix_NN.colPivHouseholderQr().solve(v_Nd);cout<<"time use in Qr composition is "<<1000 * (clock() - time_stt)/(double)CLOCKS_PER_SEC<<"ms"<<endl;cout<<endl;return 0;
}