如何用SSE3实现符号函数？

Question

1)是否有一种使用具有以下特点的SSE3 (无SSE4)有效实现符号函数的方法？

• 输入是浮点向量__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)，就像它的输入一样。

更新：似乎科里·尼尔森的答案在这里会奏效：

__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

首先想到的也许是最简单的：

__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 2

如果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时，NaN是真的:参见说明)。

如果没有AVX，这将比科里的版本(在戈德波特编译器浏览器上使用clang3.9)少两个指令。当内联到循环中时，内存源操作数可以是寄存器源操作数。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周期。科里的版本需要并行运行两个cmpps指令来实现延迟，这在Skylake上是可能的。其他CPU将出现资源冲突，导致额外的延迟周期。

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

• _mm_cmpunord_ps(x,x)去拿个面罩.(或相当于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);
}

这可以有效地编译，并且几乎不会延长关键路径延迟时间。它确实需要更多的MOVAPS指令复制寄存器没有AVX，尽管。

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

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

// 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));
}

Answer 3

如果您需要一个用于signum函数向量的float向量，其结果是一个int32_t向量，并且您不关心NaN，那么可以根据以下理论使用整数指令实现一个更有效的版本。

如果您接受一个浮点数，并将这些位重新解释为有符号的两个补整数，则可以得到三个不同的情况( 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在注释中所建议的那样，使用一个浮点比较_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)));
}

在clang3.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进一步改进。

如果需要浮点结果，最好完全停留在浮点域(而不是在浮点/整数域之间交换)，所以使用彼得解决方案之一。