Eigen基础

Eigen基础

目录

• Eigen基础
  • 1 概览
    • 1.1 Eigen是什么
    • 1.2 Eigen的优点
  • 2 Eigen库的安装
    • 2.1 通过源代码
    • 2.2 通过包管理器
  • 3 基础用法
    • 3.1 Matrix类
      • 矩阵类的前三个模板参数
      • 向量
      • 一个特殊值：Dynamic
      • 构造函数
      • 系数访问器
      • 逗号初始化
      • 固定大小的矩阵和动态大小的矩阵
    • 3.2 矩阵和向量运算
      • 加法和减法
      • 矩阵与标量乘法和除法
      • 转置和共轭
      • 矩阵与矩阵、矩阵与向量的乘法
      • 点积和叉积
      • 基本算数运算
  • 4 扩展/定制Eigen

1 概览

1.1 Eigen是什么

Eigen 是C++语言里的一个开源模板库，支持线性代数运算，矩阵和矢量运算，数值分析及其相关的算法。

1.2 Eigen的优点

• 通用性
  • 它支持所有矩阵大小，从小固定大小矩阵到任意大密度矩阵，甚至稀疏矩阵。
  • 它支持所有标准数字类型，包括std::complex、整数，并且很容易扩展到自定义数字类型。
  • 它支持各种矩阵分解和几何特征。
  • 其非支持模块的生态系统由用户提供了许多专用功能，如非线性优化，矩阵函数，多项式求解器，FFT等。
• 快速性
  • 表达式模板允许在适当的时候智能地删除临时变量并启用惰性计算。
  • 对SSE 2/3/4, AVX, AVX2, FMA, AVX512, ARM NEON(32位和64位)，PowerPC AltiVec/VSX(32位和64位)，ZVector (s390x/zEC13) SIMD指令集执行显式矢量化。
  • 固定大小的矩阵是完全优化的:避免了动态内存分配，并且在有意义时展开循环。
  • 对于大尺寸矩阵，特别注意了缓存的友好性。
• 可靠性
  • 为了保证可靠性，算法是经过精心挑选的。对可靠性的权衡有明确的文档记录，并且可以进行非常安全的分解。
  • Eigen通过它自己的测试套件(超过500个可执行文件)、标准BLAS测试套件和LAPACK测试套件的一部分进行了彻底的测试。
• 优雅
  • 由于有了表达式模板，这个API非常干净和富有表现力，而且对c++程序员来说很自然。
  • 在Eigen之上实现算法感觉就像复制伪代码。
• 很好的编译器支持
  • 在许多编译器上运行的测试套件，以保证可靠性并解决任何编译器错误。Eigen 3.4版本支持c++ 03标准，并保持合理的编译时间。Eigen 3.4之后的版本将支持c++ 14标准。

2 Eigen库的安装

2.1 通过源代码

为了使用Eigen，您只需要下载Eigen的源代码。

使用git下载Eigen源代码：

git clone https://gitlab.com/libeigen/eigen.git

目录结构：

.
├── bench
├── blas
├── ci
├── cmake
├── CMakeLists.txt
├── COPYING.APACHE
├── COPYING.BSD
├── COPYING.GPL
├── COPYING.LGPL
├── COPYING.MINPACK
├── COPYING.MPL2
├── COPYING.README
├── CTestConfig.cmake
├── CTestCustom.cmake.in
├── debug
├── demos
├── doc
├── Eigen
├── eigen3.pc.in
├── failtest
├── INSTALL
├── lapack
├── README.md
├── scripts
├── signature_of_eigen3_matrix_library
├── test
└── unsupported
 

Eigen子目录中的头文件是使用Eigen编译程序所需的唯一文件。所有平台的头文件都是相同的。没有必要使用 CMake 或安装任何东西。

在CMake工程中调用Eigen：

在CMakeLists.txt中添加include_directories指定Eigen库路径即可。

或者直接采用相对路径包含头文件。

2.2 通过包管理器

Ubuntu下使用命令sudo apt-get install libeigen3-dev安装Eigen，库文件默认安装在/usr/include/eigen3/Eigen。

在CMake工程中调用Eigen：

# CMakeLists.txt示例
cmake_minimum_required (VERSION 3.0)
project (myproject)
 
find_package (Eigen3 3.3 REQUIRED NO_MODULE)
 
add_executable (example example.cpp)
target_link_libraries (example Eigen3::Eigen)

3 基础用法

3.1 Matrix类

矩阵类的前三个模板参数

Matrix类接受6个模板参数，但现在只要了解前三个参数就足够了。另外三个参数有默认值，在此先不做讨论。

Matrix的三个必选模板参数是：Matrix<typename Scalar, int RowsAtCompileTime, int ColsAtCompileTime>

• Scalar是指标量类型，即矩阵中系数的类型。
• RowsAtCompileTime和ColsAtCompileTime是在编译时已知矩阵的行数和列数。

我们提供了许多方便的类型定义来覆盖通常的情况。例如，Matrix4f是一个4x4的浮点矩阵，它在Eigen中的定义为：typedef Matrix<float, 4, 4> Matrix4f;

向量

在Eigen中，向量只是矩阵的一种特殊情况，有一行或一列。有一列的情况是最常见的，这样的向量称为列向量，通常简称为向量。在另一种情况下它们有一行，它们被称为行向量。下面给出Eigen中一种向量和一种行向量的定义：typedef Matrix<float, 3, 1> Vector3f;、typedef Matrix<int, 1, 2> RowVector2i;。

一个特殊值：Dynamic

Eigen并不局限与在编译时维数已知的矩阵，RowsAtCompileTime和ColsAtCompileTime模板参数可以接受特殊值Dynamic，该值表示编译时大小未知，因此必须作为运行时变量处理。在Eigen的术语中，这样的尺寸被称为动态尺寸；而在编译时已知的大小称为固定大小。例如，MatrixXd表示具有动态大小的双精度矩阵，定义为typedef Matrix<double, Dynamic, Dynamic> MatrixXd;。类似地，VectorXi表示具有动态大小的整型列向量，定义为typedef Matrix<int, Dynamic, 1> VectorXi;。也可以通过Matrix<float, 3, Dynamic>获得一个行数固定、列数不固定的矩阵。

构造函数

使用默认构造函数

Matrix3f a;  // a是一个3 × 3矩阵，带有一个float[9]数组，其中包含未初始化的系数
MatrixXf b;  // b是一个动态大小的矩阵，它的大小目前是0乘0，它的系数数组还没有分配。

接受矩阵大小的构造函数

// 对于矩阵，总是先传递行数。对于向量，只需传递向量的大小。它们按照给定的大小分配系数数组，但不初始化系数
MatrixXf a(10,15); // a是一个10x15的动态大小的矩阵，带有已分配但目前未初始化的系数。
VectorXf b(30);    // b是一个大小为30的动态大小向量，具有已分配但目前未初始化的系数。
 
// 为了在固定大小的矩阵和动态大小的矩阵之间提供统一的API，在固定大小的矩阵上使用这些构造函数是合法的，即使在这种情况下传递大小是没有用的。
Matrix3f a(3,3);   // 合法，但传入尺寸实际无意义。

矩阵和向量也可以由系数列表初始化

// 在c++ 11之前，该特性仅限于大小固定的小型列或大小不超过4的向量
Vector2d a(5.0, 6.0);
Vector3d b(5.0, 6.0, 7.0);
Vector4d c(5.0, 6.0, 7.0, 8.0);
 
// 如果启用了c++ 11，那么可以通过传递任意数量的系数来初始化固定大小的列或行向量
Vector2i a(1, 2);                      // A column vector containing the elements {1, 2}
Matrix<int, 5, 1> b {1, 2, 3, 4, 5};   // A row-vector containing the elements {1, 2, 3, 4, 5}
Matrix<int, 1, 5> c = {1, 2, 3, 4, 5}; // A column vector containing the elements {1, 2, 3, 4, 5}
 
// 系数必须按行分组，并作为initializer list的initializer list进行传递
MatrixXi a {      // construct a 2x2 matrix
      {1, 2},     // first row
      {3, 4}      // second row
};
Matrix<double, 2, 3> b {
      {2, 3, 4},
      {5, 6, 7},
};
 
// 对于列向量或行向量，隐式转置是允许的。这意味着可以从单行初始化列向量
VectorXd a {{1.5, 2.5, 3.5}};             // A column-vector with 3 coefficients
RowVectorXd b {{1.0, 2.0, 3.0, 4.0}};     // A row-vector with 4 coefficients
 

系数访问器

Eigen中基础的系数访问器和修改器是重载的括号运算符。对于矩阵，总是先传递行索引。对于向量，只要传递一个下标。编号从0开始。

// Eg. 矩阵系数的访问和修改
#include <iostream>
#include <Eigen/Dense>
 
using namespace Eigen;
 
int main()
{
  MatrixXd m(2,2);
  m(0,0) = 3;
  m(1,0) = 2.5;
  m(0,1) = -1;
  m(1,1) = m(1,0) + m(0,1);
  std::cout << "Here is the matrix m:\n" << m << std::endl;
  VectorXd v(2);
  v(0) = 4;
  v(1) = v(0) - 1;
  std::cout << "Here is the vector v:\n" << v << std::endl;
}

逗号初始化

矩阵和向量可以通过逗号初始化的方式便捷的设置系数。

// Eg. 逗号初始化
Matrix3f m;
m << 1, 2, 3,
     4, 5, 6,
     7, 8, 9;
std::cout << m;

修改尺寸

矩阵的当前大小可以通过rows()， cols()和size()来获取。这些方法分别返回行数、列数和系数数。resize()方法可以调整动态大小矩阵的大小。

// Eg. 修改尺寸
#include <iostream>
#include <Eigen/Dense>
 
using namespace Eigen;
 
int main()
{
  MatrixXd m(2,5);
  m.resize(4,3);
  std::cout << "The matrix m is of size "
            << m.rows() << "x" << m.cols() << std::endl;
  std::cout << "It has " << m.size() << " coefficients" << std::endl;
  VectorXd v(2);
  v.resize(5);
  std::cout << "The vector v is of size " << v.size() << std::endl;
  std::cout << "As a matrix, v is of size "
            << v.rows() << "x" << v.cols() << std::endl;
}

如果实际矩阵大小不会改变，则resize()方法是无操作；否则它是破坏性的，系数的值可能会改变。如果您想要一个不改变系数的resize()，请使用scentrateversize()。

为了API一致性，所有这些方法仍然可用于固定大小的矩阵。当然，您无法实际调整固定大小的矩阵大小。尝试将固定大小更改为实际不同的值将触发断言失败。

固定大小的矩阵和动态大小的矩阵

什么时候应该使用固定大小（如Matrix4f），什么时候应该选择动态大小（如MatrixXf）？答案是:对于非常小的尺寸，尽可能使用固定尺寸;对于较大的尺寸，或者必须使用动态尺寸。对于较小的内存，特别是小于（大约）16的内存，使用固定大小对性能非常有利，因为它允许Eigen避免动态内存分配和展开循环。在内部，固定大小的特征矩阵只是一个普通的数组

3.2 矩阵和向量运算

Eigen提供了矩阵/向量算术运算，通过重载常见的c++算术运算符，如+，-，*，或通过特殊方法，如dot()， cross()等。对于Matrix类(矩阵和向量)，运算符只被重载支持了现行代数运算。例如matrix1 * matrix2表示矩阵与矩阵的乘积，而vector + scalar是不允许的。

加法和减法

运算符左侧和右侧矩阵的行数和列数必须相同，且具有相同的Scalar类型，因为Eigen不进行自动类型提升。执行此运算的运算符如下所示：

• 二元运算符 + 如 a+b
• 二元运算符 - 如 a-b
• 一元运算符 - 如 -a
• 复合运算符 += 如 a+=b
• 复合运算符 -= 如 a-=b

// Eg. 矩阵加减
#include <iostream>
#include <Eigen/Dense>
 
using namespace Eigen;
 
int main()
{
  Matrix2d a;
  a << 1, 2,
       3, 4;
  MatrixXd b(2,2);
  b << 2, 3,
       1, 4;
  std::cout << "a + b =\n" << a + b << std::endl;
  std::cout << "a - b =\n" << a - b << std::endl;
  std::cout << "Doing a += b;" << std::endl;
  a += b;
  std::cout << "Now a =\n" << a << std::endl;
  Vector3d v(1,2,3);
  Vector3d w(1,0,0);
  std::cout << "-v + w - v =\n" << -v + w - v << std::endl;
}

矩阵与标量乘法和除法

矩阵与标量的乘法和除法也很简单。执行这些运算的运算符是：

• 二元运算符 * 如 matrix*scalar
• 二元运算符 * 如 scalar*matrix
• 二元运算符 / 如 matrix/scalar
• 复合运算符 *= 如 matrix*=scalar
• 复合运算符 /= 如 matrix/=scalar

// Eg. 矩阵的标量乘除法运算
#include <iostream>
#include <Eigen/Dense>
 
using namespace Eigen;
 
int main()
{
  Matrix2d a;
  a << 1, 2,
       3, 4;
  Vector3d v(1,2,3);
  std::cout << "a * 2.5 =\n" << a * 2.5 << std::endl;
  std::cout << "0.1 * v =\n" << 0.1 * v << std::endl;
  std::cout << "Doing v *= 2;" << std::endl;
  v *= 2;
  std::cout << "Now v =\n" << v << std::endl;
}

转置和共轭

矩阵或向量a的转置$a^T$、共轭$\bar a$和共轭转置$a^*$分别由成员函数transpose()、conjugate()和adjoint()得到。

// Eg. 矩阵转置和共轭
MatrixXcf a = MatrixXcf::Random(2,2);
cout << "Here is the matrix a\n" << a << endl;
cout << "Here is the matrix a^T\n" << a.transpose() << endl;
cout << "Here is the conjugate of a\n" << a.conjugate() << endl;
cout << "Here is the matrix a^*\n" << a.adjoint() << endl;

矩阵与矩阵、矩阵与向量的乘法

矩阵与矩阵、矩阵与向量的乘法也是由运算符*完成的。因为向量是矩阵的一种特殊情况，它们在这里也被隐式处理，所以矩阵与向量乘积实际上只是矩阵与矩阵乘积的一种特殊情况，向量与向量的外积也是如此。因此，所有这些情况都由两个操作符处理：

• 二元运算符 * 如 a*b
• 复合运算符 *= 如 a*=b

点积和叉积

对于点积和叉积，需要使用dot()和cross()方法。

// Eg. 矩阵点积和叉积
#include <iostream>
#include <Eigen/Dense>
 
using namespace Eigen;
using namespace std;
int main()
{
  Vector3d v(1,2,3);
  Vector3d w(0,1,2);
 
  cout << "Dot product: " << v.dot(w) << endl;
  double dp = v.adjoint()*w; // automatic conversion of the inner product to a scalar
  cout << "Dot product via a matrix product: " << dp << endl;
  cout << "Cross product:\n" << v.cross(w) << endl;
}

基本算数运算

Eigen还提供一些归约运算来将给定的矩阵或向量简化为单个值，如求和(由sum()计算)、求乘积(prod())或其所有系数的最大值(maxCoeff())和最小值(minCoeff())。

// Eg. 矩阵的基本算数运算
#include <Eigen/Dense>
using namespace std;
int main()
{
  Eigen::Matrix2d mat;
  mat << 1, 2,
         3, 4;
  cout << "Here is mat.sum():       " << mat.sum()       << endl;
  cout << "Here is mat.prod():      " << mat.prod()      << endl;
  cout << "Here is mat.mean():      " << mat.mean()      << endl;
  cout << "Here is mat.minCoeff():  " << mat.minCoeff()  << endl;
  cout << "Here is mat.maxCoeff():  " << mat.maxCoeff()  << endl;
  cout << "Here is mat.trace():     " << mat.trace()     << endl;
}

4 扩展/定制Eigen

Eigen可以通过多种方式进行扩展，例如，通过定义全局方法，通过插件机制在main Eigen的类中插入自定义方法，通过添加对自定义标量类型的支持等。请参阅下面各自的子主题。

• Extending MatrixBase (and other classes)
• Inheriting from Matrix
• Using custom scalar types
• Matrix manipulation via nullary-expressions
• Adding a new expression type