а вот вам развлекуха на рабочую субботу (у нас, как минимум, рабочая). https://godbolt.org/g/y4G8WK
дано: исходник на сях:
typedef float mat4[16];
typedef float vec4[4];
void matrix4_mul(mat4 M, mat4 A, mat4 B)
{
M[ 0] = A[ 0] * B[ 0] + A[ 1] * B[ 4] + A[ 2] * B[ 8] + A[ 3] * B[12];
M[ 1] = A[ 0] * B[ 1] + A[ 1] * B[ 5] + A[ 2] * B[ 9] + A[ 3] * B[13];
M[ 2] = A[ 0] * B[ 2] + A[ 1] * B[ 6] + A[ 2] * B[10] + A[ 3] * B[14];
M[ 3] = A[ 0] * B[ 3] + A[ 1] * B[ 7] + A[ 2] * B[11] + A[ 3] * B[15];
M[ 4] = A[ 4] * B[ 0] + A[ 5] * B[ 4] + A[ 6] * B[ 8] + A[ 7] * B[12];
M[ 5] = A[ 4] * B[ 1] + A[ 5] * B[ 5] + A[ 6] * B[ 9] + A[ 7] * B[13];
M[ 6] = A[ 4] * B[ 2] + A[ 5] * B[ 6] + A[ 6] * B[10] + A[ 7] * B[14];
M[ 7] = A[ 4] * B[ 3] + A[ 5] * B[ 7] + A[ 6] * B[11] + A[ 7] * B[15];
M[ 8] = A[ 8] * B[ 0] + A[ 9] * B[ 4] + A[10] * B[ 8] + A[11] * B[12];
M[ 9] = A[ 8] * B[ 1] + A[ 9] * B[ 5] + A[10] * B[ 9] + A[11] * B[13];
M[10] = A[ 8] * B[ 2] + A[ 9] * B[ 6] + A[10] * B[10] + A[11] * B[14];
M[11] = A[ 8] * B[ 3] + A[ 9] * B[ 7] + A[10] * B[11] + A[11] * B[15];
M[12] = A[12] * B[ 0] + A[13] * B[ 4] + A[14] * B[ 8] + A[15] * B[12];
M[13] = A[12] * B[ 1] + A[13] * B[ 5] + A[14] * B[ 9] + A[15] * B[13];
M[14] = A[12] * B[ 2] + A[13] * B[ 6] + A[14] * B[10] + A[15] * B[14];
M[15] = A[12] * B[ 3] + A[13] * B[ 7] + A[14] * B[11] + A[15] * B[15];
}
всего лишь перемножение матриц 4х4. ничего сложного
если скормить это свежему gcc (8.1) и дать ключики -O3 -march=skylake -fno-strict-aliasing
то получаем на выхлопе:
matrix4_mul:
vmovss xmm0, DWORD PTR [rsi+4]
vmulss xmm0, xmm0, DWORD PTR [rdx+16]
vmovss xmm1, DWORD PTR [rsi]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx]
vmovss xmm2, DWORD PTR [rsi+8]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+32]
vmovss xmm3, DWORD PTR [rsi+12]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+48]
vmovss DWORD PTR [rdi], xmm0
vmovss xmm0, DWORD PTR [rsi+4]
vmulss xmm0, xmm0, DWORD PTR [rdx+20]
vmovss xmm4, DWORD PTR [rsi]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+4]
vmovss xmm5, DWORD PTR [rsi+8]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+36]
vmovss xmm6, DWORD PTR [rsi+12]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+52]
vmovss DWORD PTR [rdi+4], xmm0
vmovss xmm0, DWORD PTR [rsi+4]
vmulss xmm0, xmm0, DWORD PTR [rdx+24]
vmovss xmm7, DWORD PTR [rsi]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+8]
vmovss xmm1, DWORD PTR [rsi+8]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+40]
vmovss xmm2, DWORD PTR [rsi+12]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+56]
vmovss DWORD PTR [rdi+8], xmm0
vmovss xmm0, DWORD PTR [rsi+4]
vmulss xmm0, xmm0, DWORD PTR [rdx+28]
vmovss xmm3, DWORD PTR [rsi]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+12]
vmovss xmm4, DWORD PTR [rsi+8]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+44]
vmovss xmm5, DWORD PTR [rsi+12]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+60]
vmovss DWORD PTR [rdi+12], xmm0
vmovss xmm0, DWORD PTR [rsi+20]
vmulss xmm0, xmm0, DWORD PTR [rdx+16]
vmovss xmm6, DWORD PTR [rsi+16]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx]
vmovss xmm7, DWORD PTR [rsi+24]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+32]
vmovss xmm1, DWORD PTR [rsi+28]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+48]
vmovss DWORD PTR [rdi+16], xmm0
vmovss xmm0, DWORD PTR [rsi+20]
vmulss xmm0, xmm0, DWORD PTR [rdx+20]
vmovss xmm2, DWORD PTR [rsi+16]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+4]
vmovss xmm3, DWORD PTR [rsi+24]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+36]
vmovss xmm4, DWORD PTR [rsi+28]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+52]
vmovss DWORD PTR [rdi+20], xmm0
vmovss xmm0, DWORD PTR [rsi+20]
vmulss xmm0, xmm0, DWORD PTR [rdx+24]
vmovss xmm5, DWORD PTR [rsi+16]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+8]
vmovss xmm6, DWORD PTR [rsi+24]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+40]
vmovss xmm7, DWORD PTR [rsi+28]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+56]
vmovss DWORD PTR [rdi+24], xmm0
vmovss xmm0, DWORD PTR [rsi+20]
vmulss xmm0, xmm0, DWORD PTR [rdx+28]
vmovss xmm1, DWORD PTR [rsi+16]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+12]
vmovss xmm2, DWORD PTR [rsi+24]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+44]
vmovss xmm3, DWORD PTR [rsi+28]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+60]
vmovss DWORD PTR [rdi+28], xmm0
vmovss xmm0, DWORD PTR [rsi+36]
vmulss xmm0, xmm0, DWORD PTR [rdx+16]
vmovss xmm4, DWORD PTR [rsi+32]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx]
vmovss xmm5, DWORD PTR [rsi+40]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+32]
vmovss xmm6, DWORD PTR [rsi+44]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+48]
vmovss DWORD PTR [rdi+32], xmm0
vmovss xmm0, DWORD PTR [rsi+36]
vmulss xmm0, xmm0, DWORD PTR [rdx+20]
vmovss xmm7, DWORD PTR [rsi+32]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+4]
vmovss xmm1, DWORD PTR [rsi+40]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+36]
vmovss xmm2, DWORD PTR [rsi+44]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+52]
vmovss DWORD PTR [rdi+36], xmm0
vmovss xmm0, DWORD PTR [rsi+36]
vmulss xmm0, xmm0, DWORD PTR [rdx+24]
vmovss xmm3, DWORD PTR [rsi+32]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+8]
vmovss xmm4, DWORD PTR [rsi+40]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+40]
vmovss xmm5, DWORD PTR [rsi+44]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+56]
vmovss DWORD PTR [rdi+40], xmm0
vmovss xmm0, DWORD PTR [rsi+36]
vmulss xmm0, xmm0, DWORD PTR [rdx+28]
vmovss xmm6, DWORD PTR [rsi+32]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+12]
vmovss xmm7, DWORD PTR [rsi+40]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+44]
vmovss xmm1, DWORD PTR [rsi+44]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+60]
vmovss DWORD PTR [rdi+44], xmm0
vmovss xmm0, DWORD PTR [rsi+52]
vmulss xmm0, xmm0, DWORD PTR [rdx+16]
vmovss xmm2, DWORD PTR [rsi+48]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx]
vmovss xmm3, DWORD PTR [rsi+56]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+32]
vmovss xmm4, DWORD PTR [rsi+60]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+48]
vmovss DWORD PTR [rdi+48], xmm0
vmovss xmm0, DWORD PTR [rsi+52]
vmulss xmm0, xmm0, DWORD PTR [rdx+20]
vmovss xmm5, DWORD PTR [rsi+48]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+4]
vmovss xmm6, DWORD PTR [rsi+56]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+36]
vmovss xmm7, DWORD PTR [rsi+60]
vfmadd231ss xmm0, xmm7, DWORD PTR [rdx+52]
vmovss DWORD PTR [rdi+52], xmm0
vmovss xmm0, DWORD PTR [rsi+52]
vmulss xmm0, xmm0, DWORD PTR [rdx+24]
vmovss xmm1, DWORD PTR [rsi+48]
vfmadd231ss xmm0, xmm1, DWORD PTR [rdx+8]
vmovss xmm2, DWORD PTR [rsi+56]
vfmadd231ss xmm0, xmm2, DWORD PTR [rdx+40]
vmovss xmm3, DWORD PTR [rsi+60]
vfmadd231ss xmm0, xmm3, DWORD PTR [rdx+56]
vmovss DWORD PTR [rdi+56], xmm0
vmovss xmm0, DWORD PTR [rsi+52]
vmulss xmm0, xmm0, DWORD PTR [rdx+28]
vmovss xmm4, DWORD PTR [rsi+48]
vfmadd231ss xmm0, xmm4, DWORD PTR [rdx+12]
vmovss xmm5, DWORD PTR [rsi+56]
vfmadd231ss xmm0, xmm5, DWORD PTR [rdx+44]
vmovss xmm6, DWORD PTR [rsi+60]
vfmadd231ss xmm0, xmm6, DWORD PTR [rdx+60]
vmovss DWORD PTR [rdi+60], xmm0
ret
ну прикольно, он умеет в AVX набор. чот там оптимизирует.
но если этот же исходник слегка удобрить прагмами (я их тут не показал), и скормить интеловскому компилятору с хитрыми ключиками, то можно отхватить такое:
matrix4_mul:
vbroadcastss xmm3, DWORD PTR [4+rsi]
vmovups xmm2, XMMWORD PTR [16+rdx]
vmovups xmm1, XMMWORD PTR [32+rdx]
vmovups xmm0, XMMWORD PTR [48+rdx]
vmulps xmm4, xmm3, xmm2
vbroadcastss xmm5, DWORD PTR [rsi]
vbroadcastss xmm6, DWORD PTR [8+rsi]
vbroadcastss xmm8, DWORD PTR [20+rsi]
vbroadcastss xmm13, DWORD PTR [36+rsi]
vbroadcastss xmm7, DWORD PTR [12+rsi]
vbroadcastss xmm10, DWORD PTR [16+rsi]
vbroadcastss xmm15, DWORD PTR [32+rsi]
vbroadcastss xmm11, DWORD PTR [24+rsi]
vbroadcastss xmm3, DWORD PTR [40+rsi]
vbroadcastss xmm12, DWORD PTR [28+rsi]
vfmadd132ps xmm5, xmm4, XMMWORD PTR [rdx]
vmulps xmm9, xmm8, xmm2
vbroadcastss xmm4, DWORD PTR [44+rsi]
vmulps xmm14, xmm13, xmm2
vfmadd213ps xmm6, xmm1, xmm5
vfmadd132ps xmm10, xmm9, XMMWORD PTR [rdx]
vbroadcastss xmm5, DWORD PTR [52+rsi]
vfmadd132ps xmm15, xmm14, XMMWORD PTR [rdx]
vfmadd213ps xmm7, xmm0, xmm6
vfmadd213ps xmm11, xmm1, xmm10
vmulps xmm2, xmm5, xmm2
vbroadcastss xmm6, DWORD PTR [48+rsi]
vfmadd213ps xmm3, xmm1, xmm15
vfmadd213ps xmm12, xmm0, xmm11
vmovups XMMWORD PTR [rdi], xmm7
vbroadcastss xmm7, DWORD PTR [56+rsi]
vfmadd132ps xmm6, xmm2, XMMWORD PTR [rdx]
vfmadd213ps xmm4, xmm0, xmm3
vmovups XMMWORD PTR [16+rdi], xmm12
vfmadd213ps xmm7, xmm1, xmm6
vmovups XMMWORD PTR [32+rdi], xmm4
vbroadcastss xmm1, DWORD PTR [60+rsi]
vfmadd213ps xmm1, xmm0, xmm7
vmovups XMMWORD PTR [48+rdi], xmm1
ret
ну как? разница ощущается? )))))
конечно, предполагаются все условия для аццких оптимизаций: всё выровнено, никаких альясингов, никаких зависимостей и прочее.
развлекаловка в следующем: какие ключи, прагмы и прочее надо скормить в GCC или CLANG, чтобы оно выдало хоть отдалённо что-то подобное?
играться тут (по ссылке): https://godbolt.org/g/y4G8WK
дано: исходник на сях:
-
- Тот и другой код никуда не годится. Надо использовать packed-инструкции, засасывая строки матриц целиком в регистр xmm. И в этом виде их перемножая. Именно поэтому работа с матрицами 4х4 - особый случай, когда эффективность можно резко повысить. - Ксения(07.07.2018 13:24)
- ну ок, давайте лучший код? где он? - Mahagam(07.07.2018 13:32)
- Вы посмотрите по листингу, как компилятор достает данное по одному и тому же индексу: каждый раз вычисляет адрес и производит доступ иои кэширует одинаковые в промежуточной переменной. Может, придется помочь компилятору) - wza5a(07.07.2018 14:11, )
- Базовая операция - нахождение скалярного произведения. У меня она для любой длины вектора - сперва множит остаток не кратный 4 поштучно, а затем запускает цикл умножения четверок. Ксения(1742 знак., 07.07.2018 14:07)
- ну ок, давайте лучший код? где он? - Mahagam(07.07.2018 13:32)
- Тот и другой код никуда не годится. Надо использовать packed-инструкции, засасывая строки матриц целиком в регистр xmm. И в этом виде их перемножая. Именно поэтому работа с матрицами 4х4 - особый случай, когда эффективность можно резко повысить. - Ксения(07.07.2018 13:24)