在C++中实现矩阵运算（Eigen基础使用）

很多使用过matlab的人肯定对于其中方便的矩阵运算有很深的印象，那么如何在C++上同样实现像matlab一样的便捷计算呢？最近新学到了一个名为Eigen的库，这里记录一下希望加深一下印象也希望可以帮助到他人。

什么是Eigen

其实说白了Eigen就是一个库，在里面为矩阵以及向量进行了定义，并且定义了很多矩阵运算，让你像使用matlab一样在C++中实现矩阵运算（或者说就是很直白的运算），我在这里介绍的也只是其中矩阵的部分，其他部分功能如图所示：

Eigen在windows下的安装大家可以参考：http://blog.csdn.net/abcjennifer/article/details/7781936；

Linux下的安装可以参考：https://zhuanlan.zhihu.com/p/36706885

矩阵

矩阵模板

矩阵模板函数中一共包含六个模板参数，前三个是比较常用的，分别表示矩阵元素的类型、行数、列数。在矩阵定义的时候可以使用 Dynamic 来表示行或者列数未知。

template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
class Eigen::Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols >Eigen::Matrix<int, 3, 4> mat1;              //  3x4 的 int 类型的矩阵 mat1
Eigen::Matrix<double, 3, Dynamic> mat2;     //  3x? 的 double 类型的矩阵 mat2
Eigen::Matrix<float, Dynamic, 4> mat3;      //  ?x4 的 float 类型的矩阵 mat3
Eigen::Matrix<long, Dynamic, Dynamic> mat4; //  ?x? 的 long 类型的矩阵 mat4

在 Eigen 中 typedef 了很多矩阵的类型，通常命名为 Matrix 前缀加一个长度为 1∼4的字符串 S 的命名——MatrixS。其中 S 可以用来判断该矩阵类型，数字 n 表示 n ∗ n，n 的范围是2∼4，字母 d、f、i、c 表示 double、float、int、complex,另外 X 表示行或者列数未知的矩阵。

typedef Matrix<std::complex<double>, 2, 2> Eigen::Matrix2cd;            //  2x2 的 cd 类型的矩阵
typedef Matrix<double, 2, 2> Eigen::Matrix2d;                           //  2x2 的 d 类型的矩阵
typedef Matrix<std::complex<double>, 2, Dynamic> Eigen::Matrix2Xcd;     //  2x? 的 cd 类型的矩阵
typedef Matrix<std::complex<float>, Dynamic, 2> Eigen::MatrixX2cf;      //  ?x2 的 cf 类型的矩阵
typedef Matrix<std::complex<double>, Dynamic, Dynamic> Eigen::MatrixXcd;//  ?x? 的 cd 类型的矩阵
typedef Matrix<int, Dynamic, Dynamic> Eigen::MatrixXi;                  //  ?x? 的 i 类型的矩阵

行优先/列优先

在Eigen库中的矩阵默认存储是列优先的，这一点和matlab很相似，而我们知道在C++中的存储方式是行有限，所以如果我们为了自己不搞混，我们可以使用矩阵模板中的第四个参数（默认是 ColMajor）

Matrix<int, 5, 6, RowMajor> matRow; //  行优先的 5x6 的 int 类型矩阵

静态矩阵和动态矩阵

在 Eigen 中可以用 Dynamic 表示行或者列数未知，所以在定义一个矩阵时并不能确定矩阵的大小，只有在运行时才可以确定大小，然后进行动态分配，而静态矩阵则是在定义时便明确给定了行数以及列数，在编译时就可以分配好内存。

MatrixXd m = MatrixXd::Random(3,3);
Matrix3d m = Matrix3d::Random();

MatrixXd表示任意大小的元素类型为double的矩阵变量，其大小只有在运行时被赋值之后才能知道；MatrixXd::Random(3,3)表示产生一个元素类型为double的3*3的临时矩阵对象。

Matrix3d表示元素类型为double大小为3*3的矩阵变量，其大小在编译时就知道。

访问矩阵元素/矩阵元素赋值

与 C++ 数组的操作不同的是，Eigen::Matrix 是不能通过 [] 来访问赋值数据的，而是需要通过 ()。矩阵之间也可以通过 = 来进行赋值（拷贝）。（又感觉像matlab一样了= =，其中必有猫腻）

x = mat(a, b);  //  获取到矩阵 mat 的 a 行 b 列的元素并赋值给 x
mat(b, a) = x;  //  将 x 赋值给矩阵 mat 的 b 行 a 列
mat1 = mat2;    //  将矩阵 mat2 赋值（拷贝）给矩阵 mat1

通过 = 进行矩阵之间的拷贝时，如果左右两侧矩阵尺寸不一样并且左侧矩阵为动态矩阵，那么会将左侧矩阵的尺寸修改为与右侧一致。

在Eigen中重载了**"<<"**操作符，通过该操作符即可以一个一个元素的进行赋值，也可以一块一块的赋值。

Matrix3f m;
m << 1, 2, 3,
4, 5, 6,
7, 8, 9;
std::cout << m;

Eigen矩阵还可以进行分块操作，通过成员函数 block() 获取某一部分矩阵。

mat = mat1.block(i, j, p, q);   //  从矩阵 mat1 的 i 行 j 列开始获取一个 p 行 q 列的子矩阵
mat = mat1.block<p, q>(i, j);   //  从矩阵 mat1 的 i 行 j 列开始获取一个 p 行 q 列的子矩阵（动态矩阵）

Eigen矩阵可以使用成员函数 row()、col() 来获取某一行或者某一列。

mat = mat1.row(i);  //  获取 mat1 的第 i 行
mat = mat1.col(j);  //  获取 mat1 的第 j 列

Eigen矩阵还可以使用成员函数 fill()进行统一赋值。

mat.fill(n);    //  将 mat 的所有元素均赋值为 n

简单的示例：

// ==========================矩阵初值定义==============================================MatrixXd m(2, 2);//MatrixXd表示是任意尺寸的矩阵ixj, m(2,2)代表一个2x2的方块矩阵m(0, 0) = 3;//代表矩阵元素a11m(1, 0) = 2.5;//a21m(0, 1) = -1;//a12m(1, 1) = m(1, 0) + m(0, 1);//a22=a21+a12cout << "m="<<endl<<m << endl;//输出矩阵mcout << "m(2)="<<m(2) << endl; //这里注意，这里的二维数组是列有限，m(2)取得是第二列第一行的数

结果：

修改矩阵大小

当前矩阵的行数、列数、大小可以通过rows()，cols()和size()来获取，上面我们介绍过用“=”改变矩阵大小的方式除此之外，对于动态矩阵可以通过resize()函数来动态修改矩阵的大小.

需注意：

(1) 固定大小的矩阵是不能使用resize（）来修改矩阵的大小；

(2) resize（）函数会析构掉原来的数据，因此调用resize（）函数之后将不能保证元素的值不改变。

//==================================== 矩阵改变大小============================================
MatrixXd m6(3,3);cout << " the number of rows " << m6.rows() << endl;cout << " the number of cols " << m6.cols() << endl;cout << " the number of size " << m6.size() << endl;m6.resize(4, 5);cout << "matrix change!!!" << endl;cout << " the number of rows " << m6.rows() << endl;cout << " the number of cols " << m6.cols() << endl;cout << " the number of size " << m6.size() << endl;

输出结果：

矩阵计算

在矩阵计算方面，Eigen提供了许多非常方便的重载后的计算方式供我们使用，我们这里主要介绍平时经常用到的一些矩阵计算方式：

四则运算（+，-，*）

废话不说，直接上代码和结果

//=====================================加减法=============================================Matrix2d a;//这里直接定义的a是个2*2的方阵a << 1, 2,3, 4;MatrixXd b(2, 2);b << 2, 3,1, 4;cout << "a + b =\n" << a + b << endl;//矩阵加法cout << "a - b =\n" << a - b << endl;//矩阵减法cout << "Doing a += b;" << endl;a += b;//a = a + b，同时重新赋值acout << "Now a =\n" << a << endl;Vector3d v4(1, 2, 3);Vector3d w(1, 0, 0);cout << "v4 + w  =\n" <<  w + v4 << endl;//向量加减法cout << "===============next section===============================================" << endl << endl;//=================================矩阵数乘===============================================Matrix2d a2;a2 << 1, 2,3, 4;Vector3d v5(1, 2, 3);cout << "a2 * 2.5 =\n" << a2 * 2.5 << endl;//矩阵数乘cout << "0.1 * v5 =\n" << 0.1 * v5 << endl;//向量数乘

其中存在向量vector这一种数据类型，其实可以简单理解为列数为1的一个矩阵。

//====================================矩阵/向量乘法========================================Matrix2d m4;m4 << 1, 2,3, 4;Vector2d u6(-1, 1), v6(2, 0);cout << "Here is m4*m4:\n" << m4 * m4 << endl;//矩阵-矩阵cout << "Here is m4*u6:\n" << m4 * u6 << endl;//矩阵-向量cout << "Here is u6^T*m4:\n" << u6.transpose()*m4 << endl;//向量-矩阵cout << "Here is u6^T*v6:\n" << u6.transpose()*v6 << endl;//向量-向量cout << "Here is u6*v6^T:\n" << u6 * v6.transpose() << endl;//向量-向量cout << "===============next section===============================================" << endl << endl;

矩阵转置、共轭、伴随以及逆矩阵

//==============================转置，共轭,伴随以及逆矩阵=======================================MatrixXcf a3 = MatrixXcf::Random(2, 2);//定义2-by-2随机矩阵cout << "Here is the matrix a3\n" << a3 << endl;//矩阵a3cout << "Here is the matrix a3^T\n" << a3.transpose() << endl;//a3的转置cout << "Here is the matrix a3^H\n" << a3.conjugate() << endl;//a3的共轭cout << "Here is the conjugate of a3\n" << a3.conjugate() << endl;//a3的伴随矩阵cout << "Here is the matrix a3^{-1}\n" << a3.inverse() << endl;//a3的逆cout << "===============next section===============================================" << endl << endl;

矩阵（向量）点乘/叉乘

//=====================================点乘/叉乘============================================Vector3d v7(1, 2, 3);Vector3d w7(0, 1, 2);cout << "Dot product: " << v7.dot(w7) << endl;//向量点乘cout << "Cross product:\n" << v7.cross(w7) << endl;//向量叉乘cout << "===============next section===============================================" << endl << endl;

其他常见矩阵操作

//==================================== 其他功能==============================================Matrix3d m5;m5 << 1, 2, 3,1, 2, 1,0, 2, 4;cout << "Here is m5.determinant():   " << m5.determinant() << endl;// 行列式cout << "Here is m5.sum():           " << m5.sum() << endl;//所有元素之和cout << "Here is m5.prod():          " << m5.prod() << endl;//所有元素之积cout << "Here is m5.mean():          " << m5.mean() << endl;//元素的平均数cout << "Here is m5.minCoeff():      " << m5.minCoeff() << endl;//最小元素cout << "Here is m5.maxCoeff():      " << m5.maxCoeff() << endl;//最大元素cout << "Here is m5.trace():         " << m5.trace() << endl;//迹（对角元素之和)cout << "===============next section===============================================" << endl << endl;

参考博客：
Eigen 学习笔记（一）
Eigen介绍及简单使用
鸽子的学习笔记