Алгоритм Штрассена

Разберем тут же еще один вариант. Этот вариант не предполагает объявления семи матриц M1 M7, и тогда четыре блока матрицы C можно найти следующим образом:

Очевидно, что такой вариант позволяет сэкономить оперативную память, что имеет большое преимущество при умножении матриц размера более 4000 на 4000, но он, скорее всего, проиграет в скорости потому, что тут получается 12 вызовов функции, а не семь, как в предыдущем варианте. В следующий таблице приведён результат работы этих двух программ.

Из таблицы очень хорошо видно, что вариант без M1 M7 по времени отстает очень сильно. Заметно, что сейчас алгоритм Штрассена работает медленнее, чем умножение с двумя циклами. Так получается потому, что выход из рекурсии происходит очень поздно. Оставим пока этот вопрос, потому что на данном этапе возникла более важная проблема. Если попробовать умножить с помощью этой программы матрицы размера 2048 2048, оставив тот же выход из рекурсии (это 8), то возникает следующая ошибка:

Operating system error : Cannot allocate memory

Allocation would exceed memory limit

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

11

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

real , intent ( inout ), allocatable :: A(: , :) , B(: , :) real , intent ( inout ), allocatable :: C(: , :)

end subroutine end interface

integer size1 , size2 , size3 , m , N , i , j , k real x , start_strassen , finish_strassen

size1 =100 size2 =120 size3 =102

m= max ( size1 , size2 , size3 )

end do

11 allocate ( first (N , N )) allocate ( second (N , N )) do i=1 , size1

do j=1 , size2

call random_number (x) first (i ,j )= x *10

end do

end do

do i=1 , size2

do j=1 , size3

call random_number (x) second (i , j )= x *10

end do

end do

first ( size1 :N , size2 :N) = 0. second ( size2 :N , size3 :N) = 0. call cpu_time ( start_strassen )

call strassen ( first , second , third , N) call cpu_time ( finish_strassen )

print *, ’stassen calc time is ’, finish_strassen - start_strassen deallocate ( third )

end

Матрицы first и second играют роль A и B, память под эти матрицы выделяется в главной программе, потом они заполняются случайными значениями и дополняются нулевыми строками и столбцами, освобождение памяти от этих матриц происходит во время выполнения алгоритма во внешней процедуре. Для результирующей матрицы third происходит наоборот: память для нее выделяется в алгоритме(подпрограмме), а освобождается после вызова внешней подпрограммы из главной. Выделение и освобождение памяти должно происходить вовремя, и, самое

12

главное, один раз. Далее приведён текст подпрограммы.

deallocate (B)

13

call strassen (A1 , B1 , M1 , P) call strassen (A2 , B2 , M2 , P) call strassen (A3 , B3 , M3 , P) call strassen (A4 , B4 , M4 , P) call strassen (A5 , B5 , M5 , P) call strassen (A6 , B6 , M6 , P) call strassen (A7 , B7 , M7 , P)

allocate (C(S , S ))

В этой подпрограмме алгоритм реализован немного по-другому, нежели в предыдущей. Здесь явно вычисляются матрицы B1 B7 и A1 A7, которые потом используются для рекурсивной подпрограммы внутри нее самой. А M1 M7 - неявные результаты внутренних вызовов, из которых потом складывается результирующая матрица C. Эта программа позволяет умножать матрицы размера не более 8192 8192. Почему? Одна матрица размера 16384 16384 (это 214) занимает 1024 Мб, т.к. в такой матрице 228 элементов, каждый по 4 байта, всего выходит 230 байт занимает одна матрица. А 230 байт - это 210 мегабайт, а это и есть 1024. Размер оперативной памяти моего комьютера - 4 Гб, чего будет не достаточно. При вызове подпрограммы требуется две матрицы, т.е. больше 2 Гб будет занято. В ходе алгоритма нужно ввести семь матриц в два раза меньше - это чуть менее 2 Гб. На этом ресурсы исчерпаны, так что придется ограничиться матрицами 8192 8192. Вообще-то, это, конечно, не предел. Увеличить оперативную память системы возможно, выделив для нее часть жесткого диска. Это, конечно, замедлит работу программы.

Пришла пора выяснить, как быстро работает алгоритм и насколько эффективна написанная программа. Для этого начнем с матриц 128 128, будем менять выход из рекурсии и сравнивать результат с самым быстрым вариантом простого умножения с двумя циклами. В строках выход будет дан размер матриц, с которого начитает работать нерекурсивная ветвь, в каждой четвертой строке - время, за которое умножает алгоритм с двумя циклами, а в каждой третьей - время работы алгоритма Штрассена.

14

До 1024 лучший результат алгоритма достагается при самом позднем выходе из рекурсии, для 1024 1024 получилось, что оптимальный выход в середине, при размере 128. Разница во времени хоть и есть, но она все равно очень мала. Кроме этого, я хочу обратить внимание на один очень важный фактор. Время работы умножения с двумя циклами в предыдущих таблицах дано для матриц наибольшего возможного размера - 128, 256 и т.д. Но на практике матрицы размера степени двойки встречаются не очень часто, а значит, чтобы умножить, к примеру, две матрицы 700 700 по алгоритму, придется дополнять их нулями, в то время как умножение двумя циклами справится с этой задачей быстрее. Вывод напрашивается сам собой - пока алгоритм Штрассена не эффективен. Посмотрим, что будет дальше.

На этот раз алгоритм Штрассена обгоняет two loops mul на 1,45 секунды, что уже можно назвать существенным. Но стоит заметить, что использовать алгоритм для матриц меньше 1900 (two loops mul - 4.6042881) не рационально, т.к. умножение с двумя циклами справится с этой задачей быстрее.

15

Стоит заметить, что в данном случае оптимальный выход из рекурсии - 256. Посмотрим теперь, как ведет себя простое умножение. Будем брать средний размер матриц между текущим и предыдущим, потом, в зависимости от результата, еще в два раза меньше или больше.

Дополнять матрицы до 4096 стоит начиная, примерно, с 3600, а до 8196 - с 6900, что хорошо видно из таблицы.

Итак, подведем итоги. Нельзя назвать точный размер матриц, с которого нужно начинать использовать алгоритм Штрассена. Зато получились некоторые промежутки или, точнее, три множества матриц, для которых алгоритм выгоден, а именно:

1900 2048;

3600 4096;

6900 8192.

Для матриц, размер которых укладывается в эти промежутки, использование алгоритма Штрассена позволит ускорить процесс умножения. Конечно, умножать такие матрицы приходится не часто, но прикладные задачи, в которых необходимо переумножить, к примеру, матриц сорок 4000 на 4000, для которых полезен данный алгоритм, не мало. К сожалению, я не могу ничего сказать о матрицах, размер которых превосходит 8192, можно лишь заметить тенденцию увеличения эффективности алгоритма и предположить, что наступит момент, когда умножение с двумя циклами уже не будет уступать Штрассену.

16

Список литературы

[1]Дж. Макконнел. Основы современных алгоритмов.

[2]А. В. Левитин. Алгоритмы: введение в разработку и анализ.

[3]Кормен, Лейзерсон, Ривест, Штайн. Алгоритмы: построение и анализ.

[4]Metcalf M., Reid J., Cohen M. Fortran 95 2003 explained

[5]GNU Compiler Collection: Optimize Options http://gcc.gnu.org/onlinedocs/gcc/Optimize-Options.html

[6]GNU Compiler Collection: Intel 386 and AMD x86-64 Options http://gcc.gnu.org/onlinedocs/gcc/i386-and-x86002d64-Options.html

Есть вопросы - пишите ge-alena@yandex.ru. Герасимова Алёна

17