如何用 SSE3 实现签名功能？

Question

1) 有没有办法使用具有以下特征的 SSE3（无 SSE4）有效地实现sign function？

• 输入是一个浮点向量__m128。
• 输出也应该是__m128，其值为 [-1.0f, 0.0f, 1.0f]

我试过这个，但它没有用（虽然我认为它应该）：

inputVal = _mm_set_ps(-0.5, 0.5, 0.0, 3.0);
comp1 = _mm_cmpgt_ps(_mm_setzero_ps(), inputVal);
comp2 = _mm_cmpgt_ps(inputVal, _mm_setzero_ps());
comp1 = _mm_castsi128_ps(_mm_castps_si128(comp1));
comp2 = _mm_castsi128_ps(_mm_castps_si128(comp2));
signVal = _mm_sub_ps(comp1, comp2);

---

2) 有没有办法创建“标志”功能（我不确定正确的名称）。即，如果A > B，则结果将为1，否则为0。结果应该是浮点数 (__m128)，就像它的输入一样。

更新：看来 Cory Nelson 的回答在这里会起作用：

__m128 greatherThanFlag = _mm_and_ps(_mm_cmpgt_ps(valA, valB), _mm_set1_ps(1.0f));    
__m128 lessThanFlag = _mm_and_ps(_mm_cmplt_ps(valA, valB), _mm_set1_ps(1.0f));

Answer 1

您已经接近了，但您的代码无法正常工作，因为您试图仅使用强制转换将 0/-1 int 转换为浮点数。

试试这个（未经测试）：

inputVal = _mm_set_ps(-0.5, 0.5, 0.0, 3.0);
comp1 = _mm_cmpgt_ps(_mm_setzero_ps(), inputVal);
comp2 = _mm_cmpgt_ps(inputVal, _mm_setzero_ps());
comp1 = _mm_cvtepi32_ps(_mm_castps_si128(comp1)); // 0/-1 => 0.0f/-1.0f
comp2 = _mm_cvtepi32_ps(_mm_castps_si128(comp2));
signVal = _mm_sub_ps(comp1, comp2);

话虽如此，我认为Cory's solution 可能更有效。

Answer 2

首先想到的可能是最简单的：

__m128 sign(__m128 x)
{
    __m128 zero = _mm_setzero_ps();

    __m128 positive = _mm_and_ps(_mm_cmpgt_ps(x, zero), _mm_set1_ps(1.0f));
    __m128 negative = _mm_and_ps(_mm_cmplt_ps(x, zero), _mm_set1_ps(-1.0f));

    return _mm_or_ps(positive, negative);
}

或者，如果你说错了并打算得到一个整数结果：

__m128i sign(__m128 x)
{
    __m128 zero = _mm_setzero_ps();

    __m128 positive = _mm_and_ps(_mm_cmpgt_ps(x, zero),
                                 _mm_castsi128_ps(_mm_set1_epi32(1)));
    __m128 negative = _mm_cmplt_ps(x, zero);

    return _mm_castps_si128(_mm_or_ps(positive, negative));
}

Answer 3

如果您需要signum function 用于float 向量，其中结果是int32_t 向量，并且您不关心NaNs，那么可以使用整数指令实现更高效的版本，基于以下理论。

如果您采用浮点数并将这些位重新解释为有符号 二进制补码 整数，则可以得到 3 种不同的情况（其中 X 是任意的 0 或 1 , 加粗的 MSB 是符号位):

• 0 X X X X X X X X X X X X X X 1，即> 0（或> 0.0f作为浮点数）
• 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0，即== 0（或== 0.0f作为浮点数）
• 1 X X X X X X X X X X X X X X X，即< 0（或<= 0.0f作为浮点数）

最后一种情况不明确，因为它可能是负零-0.0f的特殊浮点情况：

• 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0，即== -0.0f == 0.0f作为浮点数

此时浮点符号函数变成整数函数。

---

使用 SSE3（不是 SSSE3）可用的内在函数，这可以实现为：

inline __m128i _mm_signum_ps(__m128 a)
{
    __m128i x = _mm_castps_si128(a);

    __m128i zero = _mm_setzero_si128();
    __m128i m0 = _mm_cmpgt_epi32(x, zero);
    __m128i m1 = _mm_cmplt_epi32(x, zero);
    __m128i m2 = _mm_cmpeq_epi32(x, _mm_set1_epi32(0x80000000));

    __m128i p = _mm_and_si128(m0, _mm_set1_epi32(+1));

    // Note that since (-1 == 0xFFFFFFFF) in two's complement,
    // n satisfies (n == m1), so the below line is strictly semantic
    // __m128i n = _mm_and_si128(m1, _mm_set1_epi32(-1));
    __m128i n = m1;

    return _mm_andnot_si128(m2, _mm_or_si128(p, n));
}

这个的优化版本是

inline __m128i _mm_signum_ps(__m128 a)
{
    __m128i x = _mm_castps_si128(a);

    __m128i zr = _mm_setzero_si128();
    __m128i m0 = _mm_cmpeq_epi32(x, _mm_set1_epi32(0x80000000));
    __m128i mp = _mm_cmpgt_epi32(x, zr);
    __m128i mn = _mm_cmplt_epi32(x, zr);

    return _mm_or_si128(
      _mm_andnot_si128(m0, mn),
      _mm_and_si128(mp, _mm_set1_epi32(1))
    );
}

正如 Peter 在 cmets 中建议的那样，使用一个浮点比较 _mm_cmplt_ps 而不是两个整数比较 _mm_cmplt_epi32/_mm_cmpeq_epi32 来处理 -0.0f 可以节省 1 个延迟，但它可能会受到绕过延迟延迟的影响，因为在浮点/整数域之间切换，因此最好坚持上面的仅整数实现。或不。由于您需要整数结果，因此您更有可能使用它并无论如何交换到整数域。所以：

inline __m128i _mm_signum_ps(__m128 a)
{
    __m128i x = _mm_castps_si128(a);
    __m128 zerops = _mm_setzero_ps();

    __m128i mn = _mm_castps_si128(_mm_cmplt_ps(a, zerops));
    __m128i mp = _mm_cmpgt_epi32(x, _mm_castps_si128(zerops));

    return _mm_or_si128(mn, _mm_and_si128(mp, _mm_set1_epi32(1)));
}

---

在 clang 3.9 中使用 -march=x86-64 -msse3 -O3 编译为

_mm_signum_ps(float __vector(4)):                # @_mm_signum2_ps(float __vector(4))
        xorps   xmm1, xmm1                       # fp domain
        movaps  xmm2, xmm0                       # fp domain
        cmpltps xmm2, xmm1                       # fp domain
        pcmpgtd xmm0, xmm1                       # int domain
        psrld   xmm0, 31                         # int domain
        por     xmm0, xmm2                       # int domain
        ret

除了cmpltps，这里每条指令的延迟为1，吞吐量为<= 1。我认为这是一个非常有效的解决方案，并且可以通过SSSE3 的_mm_sign_epi32 进一步改进。

---

如果您需要浮点结果，最好完全留在浮点域中（而不是在浮点/整数域之间交换），因此请使用Peter's solutions 之一。

Answer 4

如果sgn(-0.0f) 可以生成-0.0f 而不是+0.0f 的输出，则与@Cory Nelson 的版本相比，您可以保存一两条指令。请参阅下面的版本，它也传播 NaN。

• 根据x != 0.0f 的比较选择 0.0 或 1.0
• 将x 的符号位复制到那个。

// return -0.0 for x=-0.0, otherwise the same as Cory's (except for NaN which neither handle well)
__m128 sgn_fast(__m128 x)
{
    __m128 negzero = _mm_set1_ps(-0.0f);

    // using _mm_setzero_ps() here might actually be better without AVX, since xor-zeroing is as cheap as a copy but starts a new dependency chain
    //__m128 nonzero = _mm_cmpneq_ps(x, negzero);  // -0.0 == 0.0 in IEEE floating point
    __m128 nonzero = _mm_cmpneq_ps(x, _mm_setzero_ps());

    __m128 x_signbit = _mm_and_ps(x, negzero);

    __m128 zeroone = _mm_and_ps(nonzero, _mm_set1_ps(1.0f));
    return _mm_or_ps(zeroone, x_signbit);
}

当输入为 NaN 时，我认为它会根据 NaN 的符号返回 +/-1.0f。 （因为当 x 为 NaN 时 _mm_cmpneq_ps() 为真：参见 the table on the CMPPD instruction）。

没有 AVX，这比 Cory 的版本 (with clang3.9 on the Godbolt compiler explorer) 少了两条指令。当内联到循环中时，内存源操作数可以是寄存器源操作数。 gcc 使用更多指令，执行单独的 MOVAPS 加载并将自身绘制到需要额外 MOVAPS 才能将返回值放入 xmm0 的角落。

    xorps   xmm1, xmm1
    cmpneqps        xmm1, xmm0
    andps   xmm0, xmmword ptr [rip + .LCPI0_0]    # x_signbit
    andps   xmm1, xmmword ptr [rip + .LCPI0_1]    # zeroone
    orps    xmm0, xmm1

关键路径延迟为cmpneqps + andps + orps，例如在 Intel Haswell 上为 3+1+1 个周期。 Cory 的版本需要并行运行两条cmpps 指令来实现这种延迟，这只能在 Skylake 上实现。其他 CPU 将发生资源冲突，从而导致额外的延迟周期。

---

要传播 NaN，因此可能的输出将是 -1.0f、-/+0.0f、1.0f 和 NaN，我们可以利用全 1 位这一事实模式是一个 NaN。

• _mm_cmpunord_ps(x,x) 获取 NaN 掩码。 （或等效的 cmpneqps）
• or 将其添加到结果上，使其保持不变或强制为 NaN。

// return -0.0 for x=-0.0.  Return -NaN for any NaN
__m128 sgn_fast_nanpropagating(__m128 x)
{
    __m128 negzero = _mm_set1_ps(-0.0f);
    __m128 nonzero = _mm_cmpneq_ps(x, _mm_setzero_ps());

    __m128 x_signbit = _mm_and_ps(x, negzero);
    __m128 nanmask   = _mm_cmpunord_ps(x,x);
    __m128 x_sign_or_nan = _mm_or_ps(x_signbit, nanmask);   // apply it here instead of to the final result for better ILP

    __m128 zeroone = _mm_and_ps(nonzero, _mm_set1_ps(1.0f));
    return _mm_or_ps(zeroone, x_sign_or_nan);
}

这样编译效率很高，并且几乎不会延长关键路径延迟。不过，在没有 AVX 的情况下复制寄存器确实需要更多的 MOVAPS 指令。

---

您可以使用SSE4.1 BLENDVPS 做一些有用的事情，但这并不是所有 CPU 上最有效的指令。也很难避免将负零视为非零。

---

如果您想要整数结果，可以使用 SSSE3 _mm_sign_epi32(set1(1), x) 获得 -1、0 或 1 输出。如果-0.0f -> -1 太草率，您可以通过与_mm_cmpneq_ps(x, _mm_setzero_ps()) 的结果进行ANDing 来修复它

// returns -1 for x = -0.0f
__m128i sgn_verysloppy_int_ssse3(__m128 x) {
  __m128i one = _mm_set1_epi32(1);
  __m128i sign = _mm_sign_epi32(one, _mm_castps_si128(x));
  return sign;
}

// correct results for all inputs
// NaN -> -1 or 1 according to its sign bit, never 0
__m128i sgn_int_ssse3(__m128 x) {
  __m128i one = _mm_set1_epi32(1);
  __m128i sign = _mm_sign_epi32(one, _mm_castps_si128(x));

  __m128  nonzero = _mm_cmpneq_ps(x, _mm_setzero_ps());
    return _mm_and_si128(sign, _mm_castps_si128(nonzero));
}