10

我有以下功能（来自开源项目“重铸导航”）：

/// Derives the dot product of two vectors on the xz-plane. (@p u . @p v)
///  @param[in]     u       A vector [(x, y, z)]
///  @param[in]     v       A vector [(x, y, z)]
/// @return The dot product on the xz-plane.
///
/// The vectors are projected onto the xz-plane, so the y-values are ignored.
inline float dtVdot2D(const float* u, const float* v)
{
    return u[0]*v[0] + u[2]*v[2];
}

我通过 VS2010 CPU 性能测试运行它，它告诉我，在这个函数中的所有重铸代码库代码行中u[0]*v[0] + u[2]*v[2]CPU 最热。

我如何 CPU 优化（例如通过 SSE 或GLM 之类的 GLSL（如果在这种情况下更容易或更快且合适））这条线？

编辑：调用出现的上下文：

bool dtClosestHeightPointTriangle(const float* p, const float* a, const float* b, const float* c, float& h) {
    float v0[3], v1[3], v2[3];
    dtVsub(v0, c,a);
    dtVsub(v1, b,a);
    dtVsub(v2, p,a);

    const float dot00 = dtVdot2D(v0, v0);
    const float dot01 = dtVdot2D(v0, v1);
    const float dot02 = dtVdot2D(v0, v2);
    const float dot11 = dtVdot2D(v1, v1);
    const float dot12 = dtVdot2D(v1, v2);

    // Compute barycentric coordinates
    const float invDenom = 1.0f / (dot00 * dot11 - dot01 * dot01);
    const float u = (dot11 * dot02 - dot01 * dot12) * invDenom;
    const float v = (dot00 * dot12 - dot01 * dot02) * invDenom;
4

4 回答 4

21

在纸上尝试了一些东西之后，我想出了一些可能对你有用的东西。它是 SSE 中函数的适当并行化/矢量化实现。

然而，它确实需要数据重组，因为我们将同时对 4 个三角形进行并行处理。

我将逐步分解它并在这里和那里使用指令名称，但请使用 C 内在函数（_mm_load_ps()、_mm_sub_ps() 等，它们在 VC 中的 xmmintrin.h 中）——当我谈到寄存器时只是意味着__m128。

步骤1。

我们根本不需要 Y 坐标，因此我们设置了指向 X 和 Z 对的指针。每次调用至少提供 4 对（即总共 4 个三角形）。我将把每个 X 和 Z 对称为一个顶点。

第2步。

使用 MOVAPS（要求指针对齐到 16 位）将每个指针指向的前两个顶点加载到寄存器中。

从a加载的寄存器如下所示： [ a0.x, a0.z, a1.x, a1.z ]

第 3 步。

现在，使用一条减法指令，您可以一次计算 2 个顶点的增量（您的v0、v1、v2）。

不仅要计算前 2 个三角形的v0、v1和v2，还要计算后 2 个三角形！正如我所说，您应该为每个输入提供总共 4 个顶点或 8 个浮点数。只需对该数据重复步骤 2 和 3。

现在我们有 2 对vx寄存器，每对包含 2 个三角形的结果。我将它们称为vx_0（第一对）和vx_1（第二对），其中x是从 0 到 2。

第4步。

点产品。为了并行化重心计算（稍后），我们需要 4 个三角形中每个三角形的每个点积的结果，在 1 个单个寄存器中。

因此，例如，在您计算dot01的地方，我们将做同样的事情，但一次计算 4 个三角形。每个v -register 包含 2 个向量的结果，因此我们首先将它们相乘。

假设u和v - 点积函数中的参数 - 现在是v0_0和v1_0（用于计算dot01）：

将u和v相乘得到： [ (v0_0.x0) * (v1_0.x0), (v0_0.z0) * (v1_0.z0), (v0_0.x1) * (v1_0.x1), (v0_0.z1) * (v1_0.z1)]

由于.x0 / .z0和.x1 / .z1 ，这可能看起来令人困惑，但看看在第 2 步加载的内容：a0，a1。

如果现在这仍然感觉模糊，请拿起一张纸从头开始写。

接下来，仍然对于相同的点积，对v0_1和v1_1（即第二对三角形）进行乘法运算。

现在我们有 2 个寄存器，每个寄存器有 2 个 X 和 Z 对（总共 4 个顶点），相乘并准备加在一起形成 4 个单独的点积。SSE3 有一个指令可以做到这一点，它被称为 HADDPS：

xmm0 = [A、B、C、D] xmm1 = [E、F、G、H]

HADDPS xmm0, xmm1 这样做：

xmm0 = [A+B, C+D, E+F, G+H]

它将从我们的第一个寄存器和第二个寄存器中获取 X 和 Z 对，将它们加在一起并将它们存储在目标寄存器的第一个、第二个、第三个和第四个组件中。Ergo：此时我们已经为所有 4 个三角形得到了这个特殊的点积！

**现在对所有点积重复此过程：dot00等。**

第 5 步。

最后一个计算（据我所提供的代码所见）是重心的东西。这是代码中的 100% 标量计算。但是您现在的输入不是标量点积结果（即单个浮点数），它们是 SSE 向量/寄存器，4 个三角形中的每一个都有一个点积。

因此，如果您通过使用对所有 4 个浮点数进行操作的并行 SSE 操作将其平面化，您最终将得到 1 个寄存器（或结果）携带 4 个高度，每个三角形 1 个高度。

由于我的午休时间已经过去了，我不会把这个拼出来，但考虑到我给出的设置/想法，这是最后一步，应该不难弄清楚。

我知道这个想法有点牵强，需要从上面的代码中得到一些爱，也许还有一些用铅笔和纸的质量时间，但它会很快（如果你愿意，你甚至可以在之后添加 OpenMP ）。

祝你好运 ：）

（请原谅我的模糊解释，如果需要，我可以启动该功能=））

更新

我已经编写了一个实现，但它并没有像我预期的那样进行，主要是因为 Y 组件超出了您最初粘贴在问题中的代码段（我查了一下）。我在这里所做的不仅是要求您将所有点重新排列成 XZ 对并每 4 个输入它们，而且还输入 3 个指针（用于点 A、B 和 C）以及 4 个三角形中每个三角形的 Y 值。从当地的角度来看，这是最快的。我仍然可以对其进行修改，以减少被调用者端的侵入性更改，但请让我知道什么是可取的。

然后是免责声明：这段代码非常简单（我发现它在 SSE 方面与编译器一起工作得很好......他们可以根据需要重新组织，x86/x64 CPU 也参与其中）。还有命名，这不是我的风格，也不是任何人的风格，只要你认为合适就行。

希望它有所帮助，如果没有，我会很乐意再过一遍。如果这是一个商业项目，我猜也可以选择让我加入；）

无论如何，我已经把它放在了 pastebin 上：http: //pastebin.com/20u8fMEb

于 2012-06-27T12:09:56.920 回答
4

您可以使用 SSE 指令实现单点积，但结果不会比现在编写的代码快得多（甚至可能更慢）。您的重写可能会破坏有助于当前版本的编译器优化。

为了从使用 SSE 或 CUDA 重写中获得任何好处，您必须优化调用该点积的循环。对于 CUDA 来说尤其如此，因为执行一个点积的开销会很大。只有将数千个向量发送到 GPU 以计算数千个点积，才能看到加速。同样的想法也适用于 CPU 上的 SSE，但您可能会看到在较少数量的操作上有所改进。不过，它仍将不仅仅是一个点积。

最简单的尝试可能是g++ -ftree-vectorize. GCC 将能够内联您的小函数，然后尝试为您优化循环（实际上它可能已经是，但没有 SSE 指令）。树矢量化器将尝试自动执行您建议手动执行的操作。它并不总是成功的。

于 2012-06-20T18:31:00.963 回答
3

SSE 指令用于优化处理以整数或浮点数表示的大数据块的算法。典型的大小是需要以某种方式处理的数以百万计的数字。优化只处理四个（或二十个）标量的函数是没有意义的。您从 SSE 获得的收益可能会因函数调用开销而损失。一次函数调用处理的合理数字至少为千。通过使用 SSE 内在函数，您可能会获得巨大的性能提升。但很难根据您提供的信息为您提供适合您需求的具体建议。您应该编辑您的问题并提供更高级的问题视图（位于调用堆栈更深处的函数）。例如，每秒调用 dtClosestHeightPointTriangle 方法多少次并不明显？这个数字对于客观判断向 SSE 过渡是否具有实用价值至关重要。数据的组织也很重要。理想情况下，您的数据应存储在尽可能少的线性内存段中，以有效利用 CPU 的缓存子系统。

于 2012-06-27T13:54:23.897 回答
2

您要求提供算法的 SSE 版本，所以它是：

// Copied and modified from xnamathvector.inl
XMFINLINE XMVECTOR XMVector2DotXZ
(
    FXMVECTOR V1, 
    FXMVECTOR V2
)
{
#if defined(_XM_NO_INTRINSICS_)

    XMVECTOR Result;

    Result.vector4_f32[0] =
    Result.vector4_f32[1] =
    Result.vector4_f32[2] =
    Result.vector4_f32[3] = V1.vector4_f32[0] * V2.vector4_f32[0] + V1.vector4_f32[2] * V2.vector4_f32[2];

    return Result;

#elif defined(_XM_SSE_INTRINSICS_)
    // Perform the dot product on x and z
    XMVECTOR vLengthSq = _mm_mul_ps(V1,V2);
    // vTemp has z splatted
    XMVECTOR vTemp = _mm_shuffle_ps(vLengthSq,vLengthSq,_MM_SHUFFLE(2,2,2,2));
    // x+z
    vLengthSq = _mm_add_ss(vLengthSq,vTemp);
    vLengthSq = _mm_shuffle_ps(vLengthSq,vLengthSq,_MM_SHUFFLE(0,0,0,0));
    return vLengthSq;
#else // _XM_VMX128_INTRINSICS_
#endif // _XM_VMX128_INTRINSICS_
}

bool dtClosestHeightPointTriangle(FXMVECTOR p, FXMVECTOR a, FXMVECTOR b, FXMVECTOR c, float& h)
{
    XMVECTOR v0 = XMVectorSubtract(c,a);
    XMVECTOR v1 = XMVectorSubtract(b,a);
    XMVECTOR v2 = XMVectorSubtract(p,a);

    XMVECTOR dot00 = XMVector2DotXZ(v0, v0);
    XMVECTOR dot01 = XMVector2DotXZ(v0, v1);
    XMVECTOR dot02 = XMVector2DotXZ(v0, v2);
    XMVECTOR dot11 = XMVector2DotXZ(v1, v1);
    XMVECTOR dot12 = XMVector2DotXZ(v1, v2);

    // Compute barycentric coordinates
    XMVECTOR invDenom = XMVectorDivide(XMVectorReplicate(1.0f), XMVectorSubtract(XMVectorMultiply(dot00, dot11), XMVectorMultiply(dot01, dot01)));

    XMVECTOR u = XMVectorMultiply(XMVectorSubtract(XMVectorMultiply(dot11, dot02), XMVectorMultiply(dot01, dot12)), invDenom);
    XMVECTOR v = XMVectorMultiply(XMVectorSubtract(XMVectorMultiply(dot00, dot12), XMVectorMultiply(dot01, dot02)), invDenom);
}

XMVector2Dot 取自 xnamathvector.inl，我将其重命名并修改为在 X/Z 坐标上进行操作。

XNAMath 是微软的一个伟大的矢量化跨平台数学库；我也通过导入 sal.h 标头在 OS X 上使用它（我不确定那里的许可问题，所以要小心）。
事实上，任何支持内部 SSE 的平台都应该支持它。

有几点需要注意：

  • 您需要使用 XMLoadFloat3 方法将浮点数加载到 XMVECTOR 中；这会将未对齐的 float3 加载到 __m128 结构中。
  • 您可能不会看到此 SSE 代码（配置文件！！）的任何性能改进，因为将未对齐的浮点数加载到 SSE 寄存器中会降低性能。
  • 这是算法到 SSE 的强力转换，如果你比我更聪明，并且实际尝试理解算法并实现 SSE 友好版本，你会获得更好的运气。
  • 通过将整个应用程序转换为使用 XNA Math/SSE 代码而不是仅仅使用一小部分代码，您将获得更好的运气。至少强制使用对齐的向量类型（XMFLOAT3A 或 struct __declspec(align(16)) myvectortype { };）。
  • 不鼓励直接 SSE 汇编，尤其是在 x64 中，有利于内在函数。
于 2012-06-27T13:29:54.963 回答