Я ищу более быстрый и хитрый способ умножения двух матриц 4x4 на C. Мое текущее исследование сосредоточено на сборке x86-64 с расширениями SIMD. Пока что яМы создали функцию, которая примерно в 6 раз быстрее, чем простая реализация C, которая превзошла мои ожидания по улучшению производительности. К сожалению, это остается верным, только если для компиляции не используются флаги оптимизации (GCC 4.7). С-O2, C становится быстрее, и мои усилия становятся бессмысленными.

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

Моя функция выглядит следующим образом (AT &Синтаксис T, GNU Assembler):

    .text
    .globl matrixMultiplyASM
    .type matrixMultiplyASM, @function
matrixMultiplyASM:
    movaps   (%rdi), %xmm0    # fetch the first matrix (use four registers)
    movaps 16(%rdi), %xmm1
    movaps 32(%rdi), %xmm2
    movaps 48(%rdi), %xmm3
    xorq %rcx, %rcx           # reset (forward) loop iterator
.ROW:
    movss (%rsi), %xmm4       # Compute four values (one row) in parallel:
    shufps $0x0, %xmm4, %xmm4 # 4x 4FP mul's, 3x 4FP add's 6x mov's per row,
    mulps %xmm0, %xmm4        # expressed in four sequences of 5 instructions,
    movaps %xmm4, %xmm5       # executed 4 times for 1 matrix multiplication.
    addq $0x4, %rsi

    movss (%rsi), %xmm4       # movss + shufps comprise _mm_set1_ps intrinsic
    shufps $0x0, %xmm4, %xmm4 #
    mulps %xmm1, %xmm4
    addps %xmm4, %xmm5
    addq $0x4, %rsi           # manual pointer arithmetic simplifies addressing

    movss (%rsi), %xmm4
    shufps $0x0, %xmm4, %xmm4
    mulps %xmm2, %xmm4        # actual computation happens here
    addps %xmm4, %xmm5        #
    addq $0x4, %rsi

    movss (%rsi), %xmm4       # one mulps operand fetched per sequence
    shufps $0x0, %xmm4, %xmm4 #  |
    mulps %xmm3, %xmm4        # the other is already waiting in %xmm[0-3]
    addps %xmm4, %xmm5
    addq $0x4, %rsi           # 5 preceding comments stride among the 4 blocks

    movaps %xmm5, (%rdx,%rcx) # store the resulting row, actually, a column
    addq $0x10, %rcx          # (matrices are stored in column-major order)
    cmpq $0x40, %rcx
    jne .ROW
    ret
.size matrixMultiplyASM, .-matrixMultiplyASM

Он вычисляет целый столбец результирующей матрицы за одну итерацию, обрабатывая четыре числа с плавающей запятой, упакованные в 128-битные регистры SSE. Полная векторизация возможна с небольшим количеством математики (переупорядочение и агрегирование операций) и /mullpsaddps инструкции по параллельному умножению / сложению пакетов 4xfloat. Код повторно использует регистры, предназначенные для передачи параметров (,,%rdi%rsi%rdx : GNU / Linux ABI), выигрывает от развертывания (внутреннего) цикла и полностью хранит одну матрицу в регистрах XMM, чтобы уменьшить чтение памяти. Как видите, я исследовал эту тему и потратил время на то, чтобы реализовать ее как можно лучше.

Наивное вычисление C, покоряющее мой код, выглядит так:

void matrixMultiplyNormal(mat4_t *mat_a, mat4_t *mat_b, mat4_t *mat_r) {
    for (unsigned int i = 0; i < 16; i += 4)
        for (unsigned int j = 0; j < 4; ++j)
            mat_r->m[i + j] = (mat_b->m[i + 0] * mat_a->m[j +  0])
                            + (mat_b->m[i + 1] * mat_a->m[j +  4])
                            + (mat_b->m[i + 2] * mat_a->m[j +  8])
                            + (mat_b->m[i + 3] * mat_a->m[j + 12]);
}

Я исследовал оптимизированный выход сборкиs C-код, который при сохранении чисел с плавающей запятой в регистрах XMM,не включает никаких параллельных операций - только скалярные вычисления, арифметика указателей и условные переходы. КомпиляторКод кажется менее преднамеренным, но все же он немного более эффективен, чем моя векторизованная версия, которая, как ожидается, будет примерно в 4 раза быстрее. Я'я уверен, что общая идея верна - программисты делают подобные вещи с полезными результатами. Но что здесь не так? Существуют ли какие-либо проблемы с регистрацией или расписанием команд, о которых я не знаю? Знаете ли вы какие-либо инструменты или приемы сборки x86-64, чтобы поддержать мою битву против машины?