Почему моя Штрассенская матрица умножается?

Я написал две программы Matrix Multiplications в С++: Regular MM (источник) и Strassen MM (источник), оба из которых работают на квадратных матрицах размеров 2 ^ kx 2 ^ k (другими словами, квадратные матрицы четного размера).

Результаты просто ужасны. Для матрицы 1024 x 1024, Regular MM принимает 46.381 sec, а Strassen MM принимает 1484.303 sec (25 minutes!!!!).

Я попытался сохранить код как можно более простым. Другие примеры Strassen MM, найденные в Интернете, не так сильно отличаются от моего кода. Одна проблема с кодом Штрассена очевидна - у меня нет точки отсечки, которая переключается на обычный MM.

Какие еще проблемы имеет мой код Штрассена MM?

Спасибо!

Прямые ссылки на источники
http://pastebin.com/HqHtFpq9
http://pastebin.com/USRQ5tuy

Edit1. Кулак, много замечательных советов. Спасибо, что нашли время и обмен знаниями.

Я внедрил изменения (сохранил весь свой код), добавил точку отсечения. MM матрицы 2048x2048, с отсечкой 512 уже дает хорошие результаты. Обычный MM: 191.49s Strassen MM: 112.179s Значительное улучшение. Результаты были получены на доисторическом Lenovo X61 TabletPC с процессором Intel Centrino, используя Visual Studio 2012. Я сделаю больше проверок (чтобы убедиться, что у меня есть правильные результаты), и опубликую результаты.

Ответ 1

Одна проблема с кодом Штрассена очевидна - у меня нет точки отсечения, который переключается на обычный MM.

Справедливости ради следует сказать, что рекурсия до 1 пункта является основной проблемой (если не всей). Попытка угадать на других узких местах производительности, не обращаясь к этому, почти спорна из-за массивного удара производительности, который он приносит. (Другими словами, вы сравниваете Яблоки с апельсинами.)

Как обсуждалось в комментариях, выравнивание кеша может иметь эффект, но не этот масштаб. Furthemore, выравнивание кеша, скорее всего, повредит регулярному алгоритму больше, чем алгоритм Штрассена, поскольку последний не обращает внимания на кеш.

void strassen(int **a, int **b, int **c, int tam) {

    // trivial case: when the matrix is 1 X 1:
    if (tam == 1) {
            c[0][0] = a[0][0] * b[0][0];
            return;
    }

Это слишком мало. Хотя алгоритм Штрассена имеет меньшую сложность, он имеет гораздо большую константу Big-O. Во-первых, у вас есть служебный вызов функции вплоть до 1 элемента.

Это аналогично использованию слияния или быстрой сортировки и рекурсии вплоть до одного элемента. Чтобы быть эффективными, вам нужно остановить рекурсию, когда размер станет небольшим и вернуться к классическому алгоритму.

В режиме быстрой сортировки/слияния вы вернетесь к сортировке вставки или выбора с низкой загрузкой O(n^2). Здесь вы возвращаетесь к нормальной матрице O(n^3).

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

Для чего-то вроде умножения Штрассена, где преимущество только O(2.8074) над классическим O(n^3), не удивляйтесь, если этот порог окажется очень высоким. (тысячи элементов?)

В некоторых приложениях может быть много алгоритмов с уменьшением сложности, но увеличение Big-O. В результате несколько алгоритмов становятся оптимальными при разных размерах.

Большое целочисленное умножение - это печально известный пример:

Умножение школы: O (N ^ 2), оптимальное для < ~ 100 цифр *
Умножение Карацубы: O (N ^ 1.585) быстрее, чем указано выше, на ~ 100 цифр *
Toom-Cook 3-way: O (N ^ 1.465) быстрее, чем Карацуба на 3000 знаков *
FFT с плавающей запятой: O ( > N log (N)) быстрее, чем Karatsuba/Toom-3 на ~ 700 цифр *
Алгоритм Schönhage-Strassen (SSA): O (N log (n) loglog (n)) быстрее, чем FFT при ~ a миллиард цифр *
Теоретико-множественное преобразование с фиксированной шириной: O (N log (n) быстрее, чем SSA на ~ несколько миллиардов цифр? *

_{* Обратите внимание, что эти примерные пороговые значения являются приблизительными и могут сильно варьироваться - часто более чем в 10 раз.}

Ответ 2

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

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

Изменить: Вероятно, это лишь незначительно влияет на проблему. Вероятно, проблема заключается в том, что Luchian Grigore ссылается на участие проблем с контентом в строке с двумя аргументами.

Я подтвердил, что моя озабоченность верна для наивного алгоритма. Время для наивного алгоритма уменьшается почти на 50%, если массив смежный. Здесь код для этого (с использованием класса SquareMatrix, который зависит от С++ 11) на pastebin.