скорость умножения матриц зависит от глупых вещей

Вступительное примечание. Большая часть этого поста была переписана, поэтому некоторые из приведенных ниже комментариев больше не будут иметь особого смысла. Пожалуйста, посмотрите за редактированием для деталей, если вам не все равно. Так...

тл;др

• профиль для поиска горячих петель
• используйте постоянные счетчики для них, если это возможно
• попробуйте вручную развернуть их, если нет
• Результаты ОП загадка, никто другой не может воспроизвести что-то настолько экстремальное.

Я согласен с @user4581301 - чем больше компилятор знает заранее, тем больше он может сделать для вас с точки зрения оптимизации вашего кода.

Но вам нужно профилировать этот код — время на настенных часах далеко не уйдет. Я мало что знаю о профилировщиках для gcc (у ​​меня есть хороший профилировщик для MSVC), но вы можете попытать счастья здесь.

Также полезно (как сразу сказал @RetiredNinja) попытаться изучить ассемблер, используя Godbolt в качестве инструмента. , особенно если вы хотите понять столь резкое замедление.

Теперь, сказав все это, ваши тайминги не имеют для меня никакого смысла, поэтому в вашей системе происходит что-то странное. Поэтому я провел несколько собственных тестов, и мои результаты заметно отличаются от ваших. Я провел некоторые из этих тестов на MSVC (потому что у меня есть такие замечательные инструменты профилирования), а некоторые — на gcc на Mac (хотя я думаю, что на самом деле это тайный лязг под капотом). У меня нет линукс бокса, извините.

Во-первых, давайте разберемся с вопросом размещения таких больших объектов в стеке. Это явно неразумно, и я вообще не могу сделать это на MSVC, так как он не поддерживает массивы переменной длины, но мои тесты на Mac показали, что это не повлияло на тайминги после того, как я увеличил ulimit, чтобы получить его. вообще работать (см. здесь). Так что я думаю, что мы можем не учитывать это как фактор, как вы сами говорите в комментариях.

Итак, теперь давайте посмотрим на тайминги, которые я получил на Mac:

VAR=0 USE_STACK=0 N=2000 (known) N2=10 (known)
user    0m10.813s

VAR=1 USE_STACK=0 N=2000 (unknown) N2=10 (known)
user    0m11.008s

VAR=2 USE_STACK=0 N=2000 (known) N2=10 (unknown)
user    0m12.714s

VAR=3 USE_STACK=0 N=2000 (unknown) N2=10 (unknown)
user    0m12.778s

VAR=0 USE_STACK=1 N=2000 (known) N2=10 (known)
user    0m10.617s

VAR=1 USE_STACK=1 N=2000 (unknown) N2=10 (known)
user    0m10.987s

VAR=2 USE_STACK=1 N=2000 (known) N2=10 (unknown)
user    0m12.653s

VAR=3 USE_STACK=1 N=2000 (unknown) N2=10 (unknown)
user    0m12.673s

Там особо не на что смотреть; давайте перейдем к тому, что я наблюдал на MSVC (где я могу профилировать):

VAR=0 USE_STACK=0 N=2000 (known) N2=10 (known)
Elapsed: 0:00:06.89

VAR=1 USE_STACK=0 N=2000 (unknown) N2=10 (known)
Elapsed: 0:00:06.86

VAR=2 USE_STACK=0 N=2000 (known) N2=10 (unknown)
Elapsed: 0:00:10.24

VAR=3 USE_STACK=0 N=2000 (unknown) N2=10 (unknown)
Elapsed: 0:00:10.39

Теперь у нас есть кое-что, во что мы можем вникнуть. Как заметил @geza, код выполняется дольше, когда N2 неизвестно, что полностью соответствует тому, что можно было бы ожидать, поскольку именно здесь будут горячие циклы, и гораздо более вероятно, что компилятор развернет такой цикл, когда он знает, что на самом деле он состоит из небольшого известного числа итераций.

Итак, давайте получим некоторую информацию от профилировщика. Это говорит мне, что горячий цикл (довольно длинный путь) - это внутренний цикл в Mul():

inline void Mul(int N, float* z, float* x, float* y) {
    forall (i, N)
        forall (j, N) {
            double sum = 0;

=> forall (k, N) => sum += x[i*N+k] * y[kN+j]; z[iN+j] = сумма; } }

Опять же, я не могу сказать, что это меня сильно удивляет, и когда я смотрю на код, я вижу, что цикл вообще не разворачивается (код настройки цикла для краткости опущен):

$1:
    movss       xmm0,dword ptr [r9+rsi*4]  
    mulss       xmm0,dword ptr [r8+4]  
    movss       xmm1,dword ptr [r9+r15*4]  
    mulss       xmm1,dword ptr [r8]  
    cvtps2pd    xmm2,xmm0  
    cvtps2pd    xmm0,xmm1  
    movss       xmm1,dword ptr [r8+8]  
    mulss       xmm1,dword ptr [r9]  
    addsd       xmm0,xmm3  
    addsd       xmm2,xmm0  
    cvtps2pd    xmm0,xmm1  
    movss       xmm1,dword ptr [r9+r14*4]  
    movaps      xmm3,xmm2  
    mulss       xmm1,dword ptr [r8+0Ch]  
    add         r9,rbp  
    add         r8,10h  
    addsd       xmm3,xmm0  
    cvtps2pd    xmm0,xmm1  
    addsd       xmm3,xmm0  
    sub         rcx,1  
    jne         $1

Теперь не похоже, что можно было бы сэкономить, развернув этот цикл, поскольку цикл вокруг будет дешевым по сравнению с выполнением всего остального кода там, но если вы посмотрите на дизассемблирование того же самого цикла, когда известно N2, вы получаете сюрприз:

    movss       xmm0,dword ptr [rax-8]  
    mulss       xmm0,dword ptr [rcx-50h]  
    cvtps2pd    xmm2,xmm0  
    movss       xmm0,dword ptr [rcx-28h]  
    mulss       xmm0,dword ptr [rax-4]  
    addsd       xmm2,xmm7  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx]  
    mulss       xmm0,dword ptr [rax]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+28h]  
    mulss       xmm0,dword ptr [rax+4]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+50h]  
    mulss       xmm0,dword ptr [rax+8]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+78h]  
    mulss       xmm0,dword ptr [rax+0Ch]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+0A0h]  
    mulss       xmm0,dword ptr [rax+10h]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+0C8h]  
    mulss       xmm0,dword ptr [rax+14h]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+0F0h]  
    mulss       xmm0,dword ptr [rax+18h]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    movss       xmm0,dword ptr [rcx+118h]  
    mulss       xmm0,dword ptr [rax+1Ch]  
    addsd       xmm2,xmm1  
    cvtps2pd    xmm1,xmm0  
    addsd       xmm2,xmm1  
    cvtpd2ps    xmm0,xmm2  
    movss       dword ptr [rdx+rcx],xmm0  

Теперь нет цикла, и количество инструкций, которые будут выполняться в целом, явно уменьшено. Возможно, MS не такая уж и кучка тупых кудахтанов.

Итак, наконец, в качестве упражнения, давайте просто быстро развернем этот цикл вручную и посмотрим, какие тайминги мы получим (я не смотрел на сгенерированный код):

#define UNROLL_LOOP 1

inline void Mul(int N, float* z, float* x, float* y) {
    forall (i, N)
        forall (j, N) {
            double sum = 0;
#if UNROLL_LOOP
            assert (N == 10);
            sum += x[i*N] * y[0*N+j];
            sum += x[i*N+1] * y[1*N+j];
            sum += x[i*N+2] * y[2*N+j];
            sum += x[i*N+3] * y[3*N+j];
            sum += x[i*N+4] * y[4*N+j];
            sum += x[i*N+5] * y[5*N+j];
            sum += x[i*N+6] * y[6*N+j];
            sum += x[i*N+7] * y[7*N+j];
            sum += x[i*N+8] * y[8*N+j];
            sum += x[i*N+9] * y[9*N+j];
#else
            forall (k, N)
                sum += x[i*N+k] * y[k*N+j];
#endif
            z[i*N+j] = sum;
        }
}

И когда я это сделал, я получил:

VAR=3 USE_STACK=0 N=2000 (unknown) N2=10 (unknown)
Elapsed: 0:00:07.48 (compared with 10.39 / 6.86, not bad, more may be possible).

Так что это процесс, который вам нужно пройти, чтобы проанализировать подобные проблемы с производительностью, и вам нужны хорошие инструменты. Я не знаю, что происходит в вашем случае, потому что я не могу воспроизвести проблему, но развертывание цикла (как и ожидалось) является основным фактором в MSVC, когда неизвестно (небольшое) количество циклов.

Код теста, который я использовал, находится здесь на случай, если кто-то захочет сослаться на него. Я думаю, ты должен мне проголосовать, ОП.

Изменить:

Немного поиграл в Wandbox с gcc 9.0.0. Тайминги (они довольно медленные и немного более неточные, так как мы работаем на общей машине или, что более вероятно, на виртуальной машине):

VAR=0 USE_STACK=0 N=2000 (известно) N2=10 (известно) Истекшее время = ~8 секунд

VAR=3 USE_STACK=0 N=2000 (неизвестно) N2=10 (неизвестно) Истекшее время = ~15,5 с

VAR=3 USE_STACK=0 N=2000 (неизвестно) N2=10 (неизвестно), цикл развернут Истекшее время = ~ 13,5 с

Так что это требует дополнительного изучения — с помощью профилировщика и просмотра сгенерированного кода — и все еще в миллионе миль от того, что получает OP.

Живая демонстрация — вы можете сами поиграть с этим, если хотите попробовать разные вещи, OP.