前言
Eigen是一个基于线性代数的C++模板库,主要用于矩阵、向量、数值求解和相关算法,不论是做算法还是C++开发,使用好Eigen库会让你事半功倍。
正文
由于本人是非专业的算法工程师,所以这里也只是做一些简单入门介绍和使用。
类型定义
Eigen库中的核心类就是Matrix
,它表示矩阵,是一个模板类,定义如下:
template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
class Matrix { }
其中前3个模板参数需要我们指定,后面3个参数使用默认即可。
- 第一个参数
typename _Scalar
,表示该参数是一个数据类型,而不是像_Rows
和_Cols
这样的变量。那么现在比如我想定义一个4x4
的float
矩阵,就可以如下定义:
Matrix<float,4,4> matrix44f;
这样定义虽然简单明了,但是有点长,所以在Eigen库中提供了很多开源直接使用的模板类,比如:
Matrix4f
,表示4x4
的float
矩阵;Vector3f
,表示列向量,长度3的float
向量;RowVector2i
,表示行向量,长度2的int
向量;
所以对于一些常见矩阵可以使用内置的模板类。
除了定义已知大小的矩阵,Eigen库还可以定义动态矩阵。那什么是动态矩阵和静态矩阵呢?
静态矩阵就是编译时候就知道大小的矩阵,比如Matrix3d
,在编译时我们就知道这是一个3x3
的double
类型矩阵。
与之对应的是动态矩阵,这种矩阵需要在运行后才知道大小,比如需要从vector<vector<double>>
数据类型中构建矩阵,只有当代码运行时才可以知道构建的矩阵的行列个数。这时可以使用MatrixXd
来表示任意大小的元素类型为double
的矩阵变量,其中X
就表示未知大小,即Dynamic
动态的意思。
对于矩阵类型后面的后缀表示元素类型,Eigen
库定义了4种元素类型:
- d:表示
double
类型; - f:表示
float
类型; - i:表示整数;
- c:表示复数;
比如Matrix2c
,表示一个2x2
元素类型为复数的矩阵。
新建矩阵
对于默认构造函数的矩阵,分配了大小和内存空间,但是没有初始化矩阵元素。而对于基本数据类型,默认初始化的话,其值是随机的,不能使用。比如下面代码:
Eigen::Matrix2f m;
std::cout << "m:" << std::endl << m << std::endl;
Eigen::Vector4i n;
std::cout << "n:" << std::endl << n << std::endl;
//运行结果:
m:
5.88604e-039 0
9.18341e-041 2.8026e-045
n:
-517829874
-2
1952350234
4199631
在Eigen中重载了std::cout
的<<
运算符,所以可以直接使用std::cout
打印矩阵结果。所以定义完矩阵后,我们必须要对进行初始化。对于矩阵的初始化,有如下几种方式:
- 直接赋值,可以使用
<<
进行赋值。
Eigen::Matrix2f m; //2x2的float矩阵
m << 1.9,2.884, //进行初始化,期望是4个值
4.33,5.898;
std::cout << "m:" << std::endl << m << std::endl;
Eigen::Vector4i n; //长度为4的列向量
n << 3,4,2,6;
std::cout << "n:" << std::endl << n << std::endl;
//运行结果:
m:
1.9 2.884
4.33 5.898
n:
3
4
2
6
这种方式只适合于矩阵的大小固定的情况,这种初始化方式很像使用列表初始化的数组,只有知道矩阵的大小,其构造函数才能从左向右、从上向下来进行初始化。比如可以定义动态大小的矩阵MatrixXd
,这时就无法使用<<
进行直接赋值,比如:
Eigen::MatrixXi xi; //动态大小的int矩阵
xi << 1,2,3,4,5,6;
std::cout << "xi:" << std::endl << xi << std::endl;
因为这种定义的变量xi
,并不知道其行列数目,可能是3x2
、2x3
等等,所以用6个int
去赋值,构造函数无法进行构造。
- 对于向量,可以在构造的时候进行初始化。
Eigen::RowVector3d k(1,2,3); //长度为3,类型为double的行向量
std::cout << "k:" << std::endl << k << std::endl;
- 一些特殊函数,比如全部初始化为0,初始化为1,随机数等。
Eigen::MatrixXi zero = MatrixXi::Zero(3,4); //全部初始化为0的矩阵
std::cout << "------ zero ------" << std::endl << zero << std::endl;
Eigen::MatrixXi one = MatrixXi::Ones(4,4); //全部初始化为1的矩阵
std::cout << "------ one ------" << std::endl << one << std::endl;
Eigen::MatrixXi i = MatrixXi::Identity(3,3); //单位矩阵
std::cout << "------ i ------" << std::endl << i << std::endl;
Eigen::Matrix4f random = Matrix4f::Random(); //随机矩阵
std::cout << "------ random ------" << std::endl << random << std::endl;
运行结果:
------ zero ------
0 0 0 0
0 0 0 0
0 0 0 0
------ one ------
1 1 1 1
1 1 1 1
1 1 1 1
1 1 1 1
------ i ------
1 0 0
0 1 0
0 0 1
------ random ------
-0.997497 0.170019 0.64568 0.421003
0.127171 -0.0402539 0.49321 0.0270699
-0.613392 -0.299417 -0.651784 -0.39201
0.617481 0.791925 0.717887 -0.970031
- 使用
Eigen::Map
将已有的内存块数据映射到矩阵,必须要指定需要映射的矩阵类型,是MatrixXd
、VectorXd
等等,在模板参数中指明。
看一下Eigen::Map
函数的定义:
template< e> class Map
: public MapBase<Map<PlainObjectType, MapOptions, StrideType> >
{}
PlanObjectType
:映射数据的等效矩阵类型。MapOptions
:指定指针对齐方式,以字节为单位,默认为Unaligned
,即未对齐。StrideType
:可选择指定步幅。默认情况下,Map
采用内存布局是一个普通的、连续的数组,可以通过指定步幅来定制化。
正常我们只需指明第一个参数即可,使用Map
可以在无任何其他开销的情况下,让Eigen
使用非Eigen
数据结构来对其进行初始化,包括多种数据结构,比如原生数组,比如double[]
、float[]
等,以及STL容器,比如std::vector
、std::array
等。
std::vector<double> vec = {1,2,3,4,5,6}; //std::vector数据类型
Eigen::Map<Eigen::VectorXd> vd(vec.data(),6); //构造成一个长度为6的列向量
std::cout << "------ vd ------" << std::endl << vd << std::endl;
Eigen::Map<Eigen::MatrixXd> xd(vec.data(),3,2); //构造成一个3x2的矩阵
std::cout << "------ xd ------" << std::endl << xd << std::endl;
//运行结果:
------ vd ------
1
2
3
4
5
6
------ xd ------
1 4
2 5
3 6
在上面代码中,我们使用vec.data()
来获取数vector
中的数组数据,同时必须确定矩阵或者向量的大小。
- 对于矩阵来说,还有一种非常常见的初始化方式,就是一行一行或者一列一列进行赋值。比如代码:
Eigen::MatrixXd mat(3,2); //创建一个3x2的double类型矩阵
mat.col(0) << 1,2,3; //对第一列进行初始化
mat.col(1) << 4,5,6; //对第二列进行初始化
这种方式在构建特殊矩阵时经常用到。
矩阵索引
当前矩阵的行数、列数和大小可以分别通过rows()
、cols()
和size()
方法来获取,在遍历Eigen
矩阵时最好获取行数、列数来限制范围。
同时索引下标从0开始,矩阵元素的访问可以通过matrix(i,j)
来访问第i
行第j
列的元素,对于向量还可以使用类似数组取值的vector[i]
来访问第i个元素。比如代码:
std::vector<double> vec = {1,2,3,4,5,6};
Eigen::Map<Eigen::VectorXd> vd(vec.data(),6); //创建一个长度为6的列向量
std::cout << "------ vd ------" << std::endl << vd << std::endl;
Eigen::Map<Eigen::MatrixXd> xd(vec.data(),3,2); //创建一个3x2的矩阵
std::cout << "------ xd ------" << std::endl << xd << std::endl;
for (int i = 0; i < xd.rows(); ++i) {
for (int j = 0; j < xd.cols(); ++j) {
std::cout << "xd.(" << i << "," << j << ")= " << xd(i,j) << std::endl;
}
}
//运行结果:
------ xd ------
1 4
2 5
3 6
xd.(0,0)= 1
xd.(0,1)= 4
xd.(1,0)= 2
xd.(1,1)= 5
xd.(2,0)= 3
xd.(2,1)= 6
还可以利用block()
函数来从Matrix
中取出一个小矩阵来进行处理,语法是matrix:block<p,q>(i,j)
,表示从原矩阵的(i,j)
位置开始,截取出pxq
大小的子矩阵,比如测试代码:
Eigen::MatrixXi randomI = MatrixXi::Random(6,6); //随机值的6x6矩阵
std::cout << "------ randomI ------" << std::endl << randomI << std::endl;
Eigen::MatrixXi smallI = randomI.block<3,3>(1,1); //从(1,1)处截取一个3x3的矩阵
std::cout << "------ smallI ------" << std::endl << smallI << std::endl;
//运行结果:
------ randomI ------
-13389 -12482 3334 -14515 12319 -1243
-4442 -16231 3511 3528 7427 -8673
-11557 -16092 -10937 9283 14938 11869
-10948 -4002 5342 9915 13949 -9516
16007 1037 -1613 651 1289 9163
-1780 2332 -4846 -6490 -11720 11260
------ smallI ------
-16231 3511 3528
-16092 -10937 9283
-4002 5342 9915
数学运算
Eigen
库重载了常用加减乘运算符,而对于矩阵的除法,我们一般是求逆然后再转换为乘法,使用起来比较简单,就不过多介绍了。
还有几个常用的方法,必须要介绍一下,因为在项目中非常常见。
- 矩阵转置。调用
transpose()
方法即可求出一个矩阵的转置矩阵。 - 矩阵求逆。调用
inverse()
方法即可求出一个矩阵的逆矩阵。 - 共轭矩阵。调用
conjugate()
方法可以求出共轭矩阵。 - 求平均值。调用
mean()
方法求平均值。 - 逐元素相乘。调用
cwiseProduct()
函数进行2个矩阵逐元素相乘。 - 最大/小值。调用
maxCoeff()
与minCoeff()
直接求出矩阵的最大和最小值。
总结
当然在实际使用中,使用Eigen的最大场景是求解线性方程组,这个我们后面文章继续探讨。
更多算法相关总结移步专栏: juejin.cn/column/7415…