Что мне нужно сделать, чтобы GCC 4.9 распознал возможность использовать AVX FMA?

Question

Я имею std::vector<double> X,Y оба размера N (с N%16==0) а я хочу посчитать sum(X[i]*Y[i]), Это классический вариант использования Fused Multiply и Add (FMA), который должен быть быстрым на процессорах с поддержкой AVX. Я знаю, что все мои целевые процессоры - Intel, Haswell или новее.

Как мне заставить GCC излучать этот код AVX? -mfma является частью решения, но мне нужны другие переключатели?

И является std::vector<double>::operator[] мешать этому? Я знаю, что могу трансформироваться

size_t N = X.size();
double sum = 0.0;
for (size_t i = 0; i != N; ++i) sum += X[i] * Y[i];

в

size_t N = X.size();
double sum = 0.0;
double const* Xp = &X[0];
double const* Yp = &X[0];
for (size_t i = 0; i != N; ++i) sum += Xp[i] * Yp[i];

так что компилятор может заметить, что &X[0] не меняется в цикле. Но достаточно ли это или даже необходимо?

Текущий компилятор - GCC 4.9.2, Debian 8, но при необходимости может обновиться до GCC 5.

Answer 1

Вы смотрели на сборке? я кладу

double foo(std::vector<double> &X, std::vector<double> &Y) {
    size_t N = X.size();
    double sum = 0.0;
    for (size_t i = 0; i <N; ++i) sum += X[i] * Y[i];
    return sum;
}

в http://gcc.godbolt.org/ и посмотрел на сборку в GCC 4.9.2 с -O3 -mfma и я вижу

.L3:
        vmovsd  (%rcx,%rax,8), %xmm1
        vfmadd231sd     (%rsi,%rax,8), %xmm1, %xmm0
        addq    $1, %rax
        cmpq    %rdx, %rax
        jne     .L3

Так что он использует FMA. Тем не менее, это не векторизация цикла (s в sd означает один (то есть не упакованный) и d означает двойную плавающую точку).

Для векторизации цикла необходимо включить ассоциативную математику, например, с -Ofast, С помощью -Ofast -mavx2 -mfma дает

.L8:
        vmovupd (%rax,%rsi), %xmm2
        addq    $1, %r10
        vinsertf128     $0x1, 16(%rax,%rsi), %ymm2, %ymm2
        vfmadd231pd     (%r12,%rsi), %ymm2, %ymm1
        addq    $32, %rsi
        cmpq    %r10, %rdi
        ja      .L8

Так что теперь это векторизация (pd означает упакованные двойники). Однако, это не развернуто. В настоящее время это ограничение GCC. Вам нужно развернуть несколько раз из-за цепочки зависимостей. Если вы хотите, чтобы компилятор сделал это за вас, тогда подумайте об использовании Clang, который разворачивается четыре раза, в противном случае раскручивайте вручную с помощью встроенных функций.

Обратите внимание, что в отличие от GCC, Clang по умолчанию не использует fma с -mfma, Для того, чтобы использовать FMA с использованием Clang -ffp-contract=fast (например -O3 -mfma -ffp-contract=fast) или же #pragma STDC FP_CONTRACT ON или включить ассоциативную математику, например, -Ofast Вы все равно захотите включить ассоциированную математику, если хотите векторизовать цикл с Clang.

См. Fused multiply add и режимы округления по умолчанию и /questions/746991/kak-ispolzovat-instruktsii-fused-multiply-add-fma-s-sseavx/747005#747005 для получения дополнительной информации о включении fma с различными компиляторами.

---

GCC создает много дополнительного кода для обработки смещения и для N не кратно 8. Вы можете сказать компилятору, что предполагается, что массивы выровнены, используя __builtin_assume_aligned и что N кратно 8, используя N & -8

Следующий код с -Ofast -mavx2 -mfma

double foo2(double * __restrict X, double * __restrict Y, int N) {
    X = (double*)__builtin_assume_aligned(X,32);
    Y = (double*)__builtin_assume_aligned(Y,32);
    double sum = 0.0;
    for (int i = 0; i < (N &-8); ++i) sum += X[i] * Y[i];
    return sum;
}

производит следующую простую сборку

        andl    $-8, %edx
        jle     .L4
        subl    $4, %edx
        vxorpd  %xmm0, %xmm0, %xmm0
        shrl    $2, %edx
        xorl    %ecx, %ecx
        leal    1(%rdx), %eax
        xorl    %edx, %edx
.L3:
        vmovapd (%rsi,%rdx), %ymm2
        addl    $1, %ecx
        vfmadd231pd     (%rdi,%rdx), %ymm2, %ymm0
        addq    $32, %rdx
        cmpl    %eax, %ecx
        jb      .L3
        vhaddpd %ymm0, %ymm0, %ymm0
        vperm2f128      $1, %ymm0, %ymm0, %ymm1
        vaddpd  %ymm1, %ymm0, %ymm0
        vzeroupper
        ret
.L4:
        vxorpd  %xmm0, %xmm0, %xmm0
        ret

Answer 2

Я не уверен, что это поможет вам, но я почти уверен, что большая часть решения.

Вы должны разбить цикл на два: от 0 до N, с шагом M>1. Я бы попробовал с М 16, 8, 4, и посмотреть на АСМ. И внутренний цикл от 0 до М. Не беспокойтесь о математике итератор математики. Gcc достаточно умен с этим.

Gcc должен развернуть внутренний цикл, и он может его SIMD и, возможно, использовать FMA.