[Eigen中文文档] 深入了解 Eigen - Eigen内部发生了什么(一)

简介: 本页的目标是了解 Eigen 如何编译。

文档总目录

英文原文(What happens inside Eigen, on a simple example)

考虑以下示例程序:

#include<Eigen/Core>

int main()
{
   
    int size = 50;
    // VectorXf is a vector of floats, with dynamic size.
    Eigen::VectorXf u(size), v(size), w(size);
    u = v + w;
}

本页的目标是了解 Eigen 如何编译,假设启用了 SSE2 矢量化(GCC 选项 -msse2)。

为什么讨论这个问题

也许你认为上面的示例程序很简单,编译它应该不涉及任何非常有趣的事情。在开始之前,让我们解释一下编译它时的一些不普通的部分——也就是生成优化代码——以便我们在这里解释Eigen的复杂性真正有用。

看一下这行代码:

u = v + w;   //   (*)

编译它的第一件重要的事情是数组应该只遍历一次,就像:

for(int i = 0; i < size; i++) u[i] = v[i] + w[i];

问题是,如果我们创建一个简单的 C++ 库,其中 VectorXf 类有一个返回 VectorXf 的运算符+,那么代码行 (*) 将相当于:

VectorXf tmp = v + w;
VectorXf u = tmp;

显然,这里引入临时tmp是没有用的。它对性能有非常糟糕的影响,首先是因为 tmp 的创建需要在此上下文中进行动态内存分配,其次是因为现在有两个 for 循环:

for(int i = 0; i < size; i++) tmp[i] = v[i] + w[i];
for(int i = 0; i < size; i++) u[i] = tmp[i];

如果我们需要进行两次数组遍历而不是一次,这对性能来说非常糟糕,因为这意味着我们需要做很多冗余的内存访问。

编译上述程序的第二个重要方面是正确使用SSE2指令。请注意,Eigen还支持AltiVec,我们在这里讨论的所有内容也适用于AltiVec

AltiVec一样,SSE2是一组指令,允许一次处理128位数据。由于float占32位,这意味着SSE2指令可以一次处理4个float。这意味着,如果正确使用,它们可以使我们的计算速度提高4倍。

然而,在上面的程序中,我们选择了size=50,因此我们的向量由50个float组成,而50不是4的倍数。这意味着我们不能指望使用SSE2指令来执行所有的计算。我们应该着眼于次优解,即使用SSE2指令处理前48个系数,因为48是小于50的最大4的倍数,然后单独处理第49个和第50个系数。类似这样:

for(int i = 0; i < 4*(size/4); i+=4) u.packet(i)  = v.packet(i) + w.packet(i);
for(int i = 4*(size/4); i < size; i++) u[i] = v[i] + w[i];

因此,让我们逐行查看示例程序,然后跟随 Eigen 进行编译。

构建向量

我们分析一下第一行代码:

Eigen::VectorXf u(size), v(size), w(size);

首先,VectorXftypedef 定义类型:

typedef Matrix<float, Dynamic, 1> VectorXf;

类模板Matrixsrc/Core/util/ForwardDeclarations.h中声明,有6个模板参数,但最后3个参数由前3个自动确定,所以暂时不解释它们。在这里,Matrix<float,Dynamic,1>表示一个浮点数矩阵,具有1列和动态的行数。

Matrix 类继承了基类 MatrixBase。可以说 MatrixBase 统一了矩阵/向量和所有表达式类型,下面将详细介绍。

当只执行 Eigen::VectorXf u(size); 时,调用的构造函数是 Matrix::Matrix(int),位于 src/Core/Matrix.h 中。除了一些断言之外,它所做的只是构造 m_storage 成员,其类型为 DenseStorage<float, Dynamic, Dynamic, 1>

你可能会想,把存储放在一个单独的类中是不是设计过度了?原因是Matrix类模板涵盖了各种矩阵和向量:固定大小和动态大小。这两种情况下的存储方法不同。对于固定大小的矩阵,矩阵系数被存储为一个普通的成员数组。对于动态大小的矩阵,系数将被存储为指向动态分配数组的指针。因此,需要将存储从Matrix类中抽象出来。这就是DenseStorage的作用。

让我们看看在 src/Core/DenseStorage.h 中的这个构造函数。这里有许多DenseStorages的部分模板特化,单独处理维度为Dynamic或在编译时固定的情况。我们要看的部分特化是:

template<typename T, int Cols_> class DenseStorage<T, Dynamic, Dynamic, Cols_>

这里,调用的构造函数是 DenseStorage::DenseStorage(int size, int rows, int columns),其中 size=50rows=50columns=1

这个构造函数如下:

inline DenseStorage(int size, int rows, int) : m_data(internal::aligned_new<T>(size)), m_rows(rows) {
   }

这里,m_data成员是矩阵系数的实际数组,它是动态分配的。与其调用new[]malloc(),Eigen 定义了自己的internal::aligned_new,它位于 src/Core/util/Memory.h中。它的作用是,如果启用了向量化,那么它会使用平台特定的调用来分配一个128位对齐的数组,因为对于使用SSE2AltiVec进行向量化非常有用。如果未启用向量化,则相当于标准的new[]

构造函数还将m_rows成员设置为size。请注意,这里没有m_columns成员:实际上,在DenseStorage的这个部分特化中,列数是编译时确定的,因为Cols_模板参数与Dynamic不同。在这种情况下,Cols_1,这意味着我们的向量只是一个只有1列的矩阵。因此,没有必要将列数存储为运行时变量。

当调用VectorXf::data()获取系数数组的指针时,它会返回DenseStorage::data(),后者返回m_data成员。

当调用VectorXf::size()获取向量的大小时,实际上是调用了基类MatrixBase中的一个方法。它确定向量是一个列向量,因为ColsAtCompileTime==1(这来自于typedef VectorXf中的模板参数)。它推断出大小就是行数,所以返回VectorXf::rows(),后者返回DenseStorage::rows(),后者返回构造函数设置的m_rows成员。

相关文章
|
8月前
|
存储 编译器
[Eigen中文文档] 深入了解 Eigen - 类层次结构
本页面介绍了Eigen类层次结构中 Core 类的设计及其相互关系。一般用户可能不需要关注这些细节,但对于高级用户和Eigen开发人员可能会有用。
134 0
|
8月前
|
C++
[Eigen中文文档] 按值将Eigen对象传递给函数
对于 Eigen,这一点更为重要:按值传递固定大小的可向量化 Eigen 对象不仅效率低下,而且可能是非法的或使程序崩溃! 原因是这些 Eigen 对象具有对齐修饰符,在按值传递时会不遵守这些修饰符。
77 0
|
8月前
|
存储 编译器
|
8月前
|
存储 算法 NoSQL
[Eigen中文文档] 稀疏矩阵操作
在许多应用中(例如,有限元方法),通常要处理非常大的矩阵,其中只有少数系数不为零。在这种情况下,可以通过使用仅存储非零系数的特殊表示来减少内存消耗并提高性能。这样的矩阵称为稀疏矩阵。
236 0
|
8月前
|
存储
[Eigen中文文档] 就地矩阵分解
从 Eigen 3.3 开始,LU、Cholesky 和 QR 分解可以就地操作,即直接在给定的输入矩阵内操作。当处理大矩阵时,或者当可用内存非常有限(嵌入式系统)时,此功能特别有用。
49 0
|
8月前
|
存储 缓存
[Eigen中文文档] 深入了解 Eigen - 惰性求值与混叠(Aliasing)
Eigen具有智能的编译时机制,可以实现惰性求值并在适当的情况下删除临时变量。它会自动处理大多数情况下的混叠问题,例如矩阵乘积。自动行为可以通过使用MatrixBase::eval()和MatrixBase::noalias()方法手动覆盖。
98 0
|
8月前
|
存储 编译器 调度
|
8月前
|
并行计算 算法 安全
[Eigen中文文档] Eigen 和多线程
某些 Eigen 算法可以利用硬件中存在的多个内核。
220 0
|
8月前
|
存储 NoSQL API
[Eigen中文文档] Matrix类
在Eigen中,所有矩阵和向量都是Matrix模板类的对象。向量只是行数或者列数为1的特殊矩阵。
224 1
|
8月前
|
测试技术 API C++
[Eigen中文文档] 扩展/自定义Eigen(一)
在本节中,将介绍如何向MatrixBase添加自定义方法。由于所有表达式和矩阵类型都继承自MatrixBase,因此将方法添加到MatrixBase会立即使其对所有表达式可用!一个典型的用例是,使Eigen与另一个API兼容。
118 0