vpmadd52luq

浮動小数の乗算は、整数の乗算よりbit数が少なく済みます。


IEEE754浮動小数だと、

(2^expA * mantA) * (2^expB * mantB) 

は、

(2^(expA+expB) * (mantA * mantB))

になります。つまり、double だと、11bitの加算 + 52bitの乗算、float だと、8bitの加算 + 23bitの乗算があればよいわけです。

これは、long long が64bitの乗算、int が 32bit の乗算を必要とするのに比べると、ハードウェア的にリソースが少なくて済みますね。

世間的に、double の乗算は必要性がかなり高いですが、long long の乗算はそれほど使う機会は無いです。

結果として、x64 + SSE2,AVXのCPUでは乗算命令が、

|                   | simd | scalar |
|------------------|-----|-------|
| 32bit 整数        | 有る | 有る   |
| 64bit 整数        | **無い** | 有る   |
| 32bit 浮動小数    | 有る | 有る   |
| 64bit 浮動小数    | 有る | 有る   |

というようになっていて、64bit 整数乗算は、並列に実行できませんでした。まあでも、52bit乗算が付いているのだから、それを使えばいいんではない？という気がしますね。それを実装したのが、vpmadd52luq, vpmadd52huq です。


IEEE754とか考えなくても、long long の乗算が遅くて、double の乗算が速いCPUがあったとすると、long long を一旦double に変換して、乗算して、long long に戻すテクニックとか考えられますが、それの無駄を省いているとも言えそうです。

vpmadd52luq は、52bit x 52bit で乗算して出る、104bitのうち下位52bit を、vpmadd52huq は 104bit のうち上位52bitを取得して、それと 1個目のオペランドの値を64bit加算します。


```cpp
#include <immintrin.h>
#include <stdio.h>

unsigned long long in0[8] = {0x8000000000000001ULL, 0x0000000000000001ULL, 0};
unsigned long long in1[8] = {0xffffffffffffffffULL, 0x0008000000000000ULL, ~0ULL};
unsigned long long in2[8] = {0x0000000080000000ULL, 0x0008000000000000ULL, ~0ULL};
unsigned long long out[8];


static void
test_l(void)
{
    register __m512i dest __asm__("zmm0") = _mm512_loadu_si512(in0);
    register __m512i src1 __asm__("zmm1") = _mm512_loadu_si512(in1);
    register __m512i src2 __asm__("zmm2") = _mm512_loadu_si512(in2);

    __asm__ __volatile__ (".byte 0x62, 0xf2, 0xf5, 0x48, 0xb4, 0xc2 # vpmadd52luq %%zmm2, %%zmm1, %%zmm0\n\t"
                          :[dest]"+v"(dest)
                          :[src1]"v"(src1),
                           [src2]"v"(src2));

    _mm512_storeu_si512(out, dest);

    int i;
    for (i=0; i<3; i++) {
        long long r64 = in1[i] * in2[i] + in0[i];
        printf("%2d:%016llx, %016llx\n", i, out[i], r64);
    }
}

static void
test_h(void)
{
    register __m512i dest __asm__("zmm0") = _mm512_loadu_si512(in0);
    register __m512i src1 __asm__("zmm1") = _mm512_loadu_si512(in1);
    register __m512i src2 __asm__("zmm2") = _mm512_loadu_si512(in2);

    __asm__ __volatile__ (".byte 0x62, 0xf2, 0xf5, 0x48, 0xb5, 0xc2 # vpmadd52huq %%zmm2, %%zmm1, %%zmm0\n\t"
                          :[dest]"+v"(dest)
                          :[src1]"v"(src1),
                           [src2]"v"(src2));

    _mm512_storeu_si512(out, dest);

    int i;
    for (i=0; i<3; i++) {
        printf("%2d:%016llx, ________________\n", i, out[i]);
    }
}


int
main()
{
    printf("   %-16s, %-16s\n", "madd52", "scalar");
    puts("vpmadd52luq");
    test_l();
    puts("vpmadd52huq");
    test_h();
}

```

```
 $ ./sde -cnl -- ./a.out 
   madd52          , scalar          
vpmadd52luq
 0:800fffff80000001, 7fffffff80000001
 1:0000000000000001, 0000000000000001
 2:0000000000000004, 0000000000000004
vpmadd52huq
 0:8000000080000000, ________________
 1:0004000000000001, ________________
 2:000ffffffffffffc, ________________
```


これはまた別の AVX-512IFMA52 という拡張になっていて、KNLでもSkylakeでも動きません。

"cnl" は、Cannonlake のことで、Skylakeの次のCPUのようです。


明日は @tanakmura がそういえば書いていなかった kor について書きます。