玩转 ArrayFire：06 矢量化介绍

---

前言

目前，程序员和数据科学家希望利用快速并行计算设备的优势。为了从当前的并行硬件和科学计算软件中获得最佳性能，有必要使用向量化代码。然而，编写向量化代码可能不是立即直观的。ArrayFire 提供了许多向量化给定代码段的方法。在本篇中，我们将介绍几种使用 ArrayFire 对代码进行向量化的方法，并讨论每种方法的优缺点。

---

一、通用/默认矢量化

    就其本质而言，ArrayFire 是一个矢量化库。大多数函数对 array 进行整体操作——对所有元素进行并行操作。只要可能，应该使用现有的矢量化函数，而不是手工索引到数组中。例如，考虑下面的代码:

	af::array a = af::range(10); // [0,  9]
	for(int i = 0; i < a.dims(0); ++i)
	{
	    a(i) = a(i) + 1;         // [1, 10]
	}

    尽管代码完全有效，但效率非常低，因为它导致内核只对一个数据进行操作。相反，开发人员应该使用 ArrayFire 的 + 操作符重载：

	af::array a = af::range(10);  // [0,  9]
	a = a + 1;                    // [1, 10]

    这段代码将导致一个内核并行地操作 a 的所有 10 个元素。
    大多数 ArrayFire 函数都是矢量化的。其中的一小部分包括：

操作符类型	函数
算术运算	+, -, *, /, %, >>, <<
逻辑运算	&&, ||, <, >, ==, !=
数值函数	abs(), floor(), round(), min(), max()
复数运算	real(), imag(), conj()
指数和对数函数	exp(), log(), expm1(), log1p()
三角函数	sin(), cos(), tan()
双曲函数	sinh(), cosh(), tanh()

    除了元素操作之外，许多其他函数也在 ArrayFire 中矢量化。
    请注意，即使执行某种形式的聚合（如 sum() 或 min() ）、信号处理（如 convolve() ）、甚至图像处理函数（如 rotate() ），ArrayFire 都支持对不同的列或图像进行矢量化。例如，如果我们有宽为WIDTH、高为HEIGHT的NUM个图像，我们可以用如下向量方式卷积每个图像：

	float g_coef[] = { 1, 2, 1,
	                   2, 4, 2,
	                   1, 2, 1 };
	af::array filter = 1.f/16 * af::array(3, 3, f_coef);
	af::array signal = randu(WIDTH, HEIGHT, NUM);
	af::array conv = convolve2(signal, filter);

    类似地，你可以使用如下代码在一次调用中将 100 张图像旋转 45 度：

	// Construct an array of 100 WIDTH x HEIGHT images of random numbers
	af::array imgs = randu(WIDTH, HEIGHT, 100);
	// Rotate all of the images in a single command
	af::array rot_imgs = rotate(imgs, 45);

    虽然 ArrayFire 中的大多数函数都支持矢量化，但也有一些不支持。最明显的是，所有的线性代数函数。即使它们不是矢量化的，线性代数操作仍然在你的硬件上并行执行。
    想要矢量化那些使用 ArrayFire 编写的任何代码，使用内置的矢量化操作是最佳首选方法。

二、GFOR：并行的 for 循环

     ArrayFire 中提出的另一种新的矢量化方法是 GFOR 循环替换构造。GFOR 允许在 GPU 或设备上并行地启动循环的所有迭代，只要迭代是独立的。标准的 for 循环按顺序执行每个迭代，而 ArrayFire 的 gfor 循环则同时(并行地)执行每个迭代。ArrayFire 通过“平铺”所有循环迭代的值来实现这一点，然后在这些“平铺”的值上一次执行计算。你可以把 gfor 看作是对你的代码执行自动矢量化，例如，你编写了一个 gfor 循环，对 vector 的每个元素递增，但在幕后，ArrayFire 会重写它，以并行地对整个 vector 进行操作。
    本篇开头的 for 循环示例可以使用 GFOR 重写，如下所示：

	af::array a = af::range(10);
	gfor(seq i, n)
	    a(i) = a(i) + 1;

    在这种情况下，gfor 循环的每个实例都是独立的，因此 ArrayFire 将自动平铺设备内存中的 a 数组，并且并行执行增量内核。
    看看另一个例子，你可以在 for 循环中的每个矩阵切片上运行 accum() ，或者你可以“矢量化”并简单地在 gfor 循环操作中完成所有操作：

	   // runs each accum() in sequence
	for (int i = 0; i < N; ++i)
	   B(span,i) = accum(A(span,i));
	// runs N accums in parallel
	gfor (seq i, N)
	   B(span,i) = accum(A(span,i));

    然而，回到我们前面的矢量化技术，accum() 已经矢量化了，只需用：

	B = accum(A);

    最好尽可能地运用向量化计算，以避免 for 循环和 gfor 循环中的开销。然而，gfor循环结构在广播样式操作的狭窄情况下是最有效的。考虑这样一种情况：我们有一个常量向量，我们希望将其应用于变量集合，例如表示多个向量的线性组合的值。将一组常量广播到多个向量在 gfor 循环中可以很好地工作：

	const static int p=4, n=1000;
	af::array consts = af::randu(p);
	af::array var_terms = randn(p, n);
	gfor(seq i, n)
	    combination(span, i) = consts * var_terms(span, i);

    使用 GFOR 需要遵循几条规则和多条指导方针以获得最佳性能。这个向量化方法的详细信息将在后面篇节中进行介绍。

三、批处理

     batchFunc() 函数允许将现有的 ArrayFire 函数广泛应用于多个数据集。实际上， batchFunc() 允许 ArrayFire 函数以“批处理”模式执行。在这种模式下，函数将找到一个包含要处理的“批”数据的维度，并将该过程并行化。
    考虑下面的例子。这里，我们创建一个滤波器，并将其应用到每个权向量。简单的解决方案是使用 for 循环，就像我们之前看到的那样：

	// Create the filter and the weight vectors
	af::array filter = randn(5, 1);
	af::array weights = randu(5, 5);
	// Apply the filter using a for-loop
	af::array filtered_weights = constant(0, 5, 5);
	for(int i=0; i<weights.dims(1); ++i){
	    filtered_weights.col(i) = filter * weights.col(i);
	}

    然而，正如我们上面所讨论的，这个解决方案将非常低效。一个人可能会试图实现一个矢量化的解决方案如下：

	// Create the filter and the weight vectors
	af::array filter = randn(1, 5);
	af::array weights = randu(5, 5);
	af::array filtered_weights = filter * weights; // fails due to dimension mismatch

    然而，滤波器的尺寸和权重不匹配，因此 ArrayFire 将产生一个运行时错误。batchfunc() 就是为了解决这一特定的问题创建的。函数声明如下：

	array batchFunc(const array &lhs, const array &rhs, batchFunc_t func);

    其中 batchFunc_t 是一个函数指针的形式：

	typedef array (*batchFunc_t) (const array &lhs, const array &rhs);

    因此，要使用 batchFunc() ，我们需要提供希望作为批处理操作应用的函数。为了便于说明，让我们按照下面的格式“实现”一个乘法函数。

af::array my_mult (const af::array &lhs, const af::array &rhs){
    return lhs * rhs;
}

    最后的批处理调用并不比我们想象的理想语法困难多少。

	// Create the filter and the weight vectors
	af::array filter = randn(1, 5);
	af::array weights = randu(5, 5);
	// Apply the batch function
	af::array filtered_weights = batchFunc( filter, weights, my_mult );

    批处理函数可以与前面提到的许多向量化的 ArrayFire 函数一起使用。如果将这些函数包装在一个匹配 batchFunc_t 声明的 helper 函数中，它甚至可以使用这些函数的一个组合。batchfunc() 的一个限制是，目前它不能在 gfor() 循环中使用。

四、高级矢量化

    我们已经看到了 ArrayFire 为矢量化代码提供的不同方法。将它们组合在一起是一个稍微复杂一些的过程，需要考虑数据维度和布局、内存使用、嵌套顺序等。官方博客上有一个关于这些因素的很好的例子和讨论：

    http://arrayfire.com/how-to-write-vectorized-code/

    值得注意的是，博客中讨论的内容已经转换为一个方便的 af::nearestNeighbour() 函数。在从头开始编写之前，检查 ArrayFire 是否已经有了实现。ArrayFire 默认的矢量化特性和大量的函数集合除了替换几十行代码之外，还可以加快速度！

---