我想在SSE中实现一个简单的功能(像Izhikevich spiking neuron model这样的程序)。它应该使用16位带符号整数(8.8个固定点),并且需要在某些积分步骤中检查溢出条件,并设置SSE掩码(如果发生溢出):
// initialized like following:
short I = 0x1BAD; // current injected to neuron
short vR = 0xF00D; // some reset threshold when spiked (negative)
// step to be vectorized:
short v0 = vReset;
for(;;) {
// v0*v0/16 likely overflows => use 32 bit (16.16)
short v0_sqr = ((int)v0)*((int)v0) / (1<<(8+4)); // not sure how "(v0*v0)>>(8+4)" would affect sign..
// or ((int)v0)*((int)v0) >> (8+4); // arithmetic right shift
// original paper used v' = (v0^2)/25 + ...
short v1 = v0_sqr + v0 + I;
int m; // mask is set when neuron fires
if(v1_overflows_during_this_operation()) { // "v1 > 0x7FFF" - way to detect?
m=0xFFFFFFFF;
else
m=0;
v0 = ( v1 & ~m ) | (vR & m );
}
但是我还没有找到
_mm_mul_epi16()
指令来检查乘法的高位字。为什么以及应该如何在SSE中实现此类任务v1_overflows_during_this_operation()
? 最佳答案
与32x32 => 64不同,没有扩展16x16-> 32 SSE乘法指令。
取而代之的是 _mm_mulhi_epi16
and _mm_mulhi_epu16
,它为您提供整个结果的带符号或无符号的上半部分。
(还有_mm_mullo_epi16
,它打包16x16 => 16位低位半截尾乘法,这对于有符号或无符号是相同的)。
您可以使用_mm_unpacklo/hi_epi16
将低半部分/高半部分交织成一对具有32位元素的 vector ,但这会非常慢。但是,是的,您可以将_mm_srai_epi32(v, 8+4)
算术右移12,然后重新打包,也许使用_mm_packs_epi32
(将符号饱和度恢复为16位)。那我想检查一下饱和度吗?
您的用例很不寻常。有 _mm_mulhrs_epi16
,它为您提供高17位,四舍五入然后被截断为16位。 (请参阅说明)。这对于某些定点算法很有用,在定点算法中,按比例缩放输入以将结果放在上半部分,并且您想舍入包括下半部分而不是截断。
实际上,您可能会使用_mm_mulhrs_epi16
或_mm_mulhi_epi16
作为保持最高精度的最佳选择,也许是在平方前左移v0
到恰好高半数会为您提供(v0*v0) >> (8+4)
的程度。
当然好!获得另一或第二个精度将使您的性能损失2倍,因为您必须计算另一个乘法结果以检查溢出,或有效地扩展到32位(将每个 vector 的元素数量减少一半) ), 如上所述。
通过比较结果,v0 = ( v1 & ~m ) | (vR & m );
成为SSE4.1混合: _mm_blendv_epi8
。
如果您的vThreshold
的顶部有2个未设置的位,则您有向左移动的空间,而不会丢失任何最重要的位。由于mulhi
提供了(v0*v0) >> 16
,因此您可以执行以下操作:
// losing the high 2 bits of v0
__m128i v0_lshift2 = _mm_slli_epi16(v0, 2); // left by 2 before squaring
__m128i v0_sqr_asr12 = _mm_mulhi_epi16(v0_lshift2, v0_lshift2);
__m128i v1 = _mm_add_epi16(v0, I);
v1 = _mm_add_epi16(v1, v0_sqr_asr12);
// v1 = ((v0<<2)* (int)(v0<<2))) >> 16) + v0 + I
// v1 = ((v0*(int)v0) >> 12) + v0 + I
平方前左移2与平方后左移4相同(完整的32位结果)。它将我们想要的16位准确地放入高16位。
但是,如果您的
v0
太接近全范围,以至于在向左移动时可能会溢出,那么这将是不可用的。否则,您可能会失去
v0
的6个低位,然后再乘以通过算术右移舍入到-Infinity会失去6位精度,但是不可能发生溢出。
// losing the low 6 bits of v0
__m128i v0_asr6 = _mm_srai_epi16(v0, 6);
__m128i v0_sqr_asr12 = _mm_mullo_epi16(v0_asr6, v0_asr6);
__m128i v1 = _mm_add_epi16(v0, I);
v1 = _mm_add_epi16(v1, v0_sqr_asr12);
// v1 = (v0>>6) * (int)(v0>>6)) + v0 + I
// v1 ~= ((v0*(int)v0) >> 12) + v0 + I
我认为您会以这种方式失去更多的精度,因此最好将
vThreshold
设置得足够小,以使您有足够的开销来使用上半倍乘法。这种方式包括可能更差的舍入。如果我们可以有效地设置
pmulhrsw
而不是截断,它甚至会更好。但是我不认为可以,因为右移1是一个奇数。我认为我们需要进行2个单独的输入,一个v0_lshift2
和一个仅左移1个。关于c - 如何检查SSE中16位整数乘法的溢出?,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/52598358/