Почему в matlab сильно замедляется выполнение алгоритма Винограда-Штрассена?
По учебе надо сравнить классический алгоритм умножения матриц и алгоритм Винограда-Штрассена. Мне непонятно, почему алгоритм Винограда-Штрассена прям сильно замедляется по сравнению с классическим. Я понимаю, что MatLab оптимизированнее выполняет именно классический алгоритм, особенно не с супер супер большими матрицами. И что на использование алгоритма Винограда-Штрассена тратится время на выделение в памяти места для промежуточных матриц и все такое. Но просто замедление при размере матрицы очень большое - приложил на картинке. На картинке классический алгоритм при порядке матрицы 2^10 выполняется пару секунд, а алгоритм Винограда-Штрассена 7 минут.
Вот код, мб не так алгоритм понял. Параметр type нужен просто для вызова одной функции с алгоритмами Штрассена и Винограда-Штрассена. Больше вродь нет особенностей
clc;
function AA = split_to_2x2_blocks(matrix) % разделение матрицы на 4 части
buf = width(matrix);
A11 = matrix(1:buf/2, 1:buf/2);
A12 = matrix(1:buf/2, buf/2+1:buf);
A21 = matrix(buf/2+1:buf, 1:buf/2);
A22 = matrix(buf/2+1:buf, buf/2+1:buf);
AA = {A11, A12
A21, A22};
end
function CC = strassen_mul_2x2(lb, rb, type) % умножение матрицы после разделения, рекурсивная чушь тут и все такое
d = strassen_mul(cell2mat(lb(1,1)) + cell2mat(lb(2,2)), cell2mat(rb(1,1)) + cell2mat(rb(2,2)), type);
d_1 = strassen_mul(cell2mat(lb(1,2)) - cell2mat(lb(2,2)), cell2mat(rb(2,1)) + cell2mat(rb(2,2)), type);
d_2 = strassen_mul(cell2mat(lb(2,1)) - cell2mat(lb(1,1)), cell2mat(rb(1,1)) + cell2mat(rb(1,2)), type);
left = strassen_mul(cell2mat(lb(2,2)), cell2mat(rb(2,1)) - cell2mat(rb(1,1)), type);
right = strassen_mul(cell2mat(lb(1,1)), cell2mat(rb(1,2)) - cell2mat(rb(2,2)), type);
top = strassen_mul(cell2mat(lb(1,1)) + cell2mat(lb(1,2)), cell2mat(rb(2,2)), type);
bottom = strassen_mul(cell2mat(lb(2,1)) + cell2mat(lb(2,2)), cell2mat(rb(1,1)), type);
c1 = d + d_1 + left - top;
c2 = right + top;
c3 = left + bottom;
c4 = d + d_2 + right - bottom;
CC = [c1, c2
c3, c4];
end
function CC = vinograd_strassen_mul_2x2(lb, rb, type)
s1 = cell2mat(lb(2,1)) + cell2mat(lb(2,2));
s2 = s1 - cell2mat(lb(1,1));
s3 = cell2mat(lb(1,1)) - cell2mat(lb(2,1));
s4 = cell2mat(lb(1,2)) - s2;
s5 = cell2mat(rb(1,2)) - cell2mat(rb(1,1));
s6 = cell2mat(rb(2,2)) - s5;
s7 = cell2mat(rb(2,2)) - cell2mat(rb(1,2));
s8 = s6 - cell2mat(rb(2,1));
p1 = strassen_mul(s2, s6, type);
p2 = strassen_mul(cell2mat(lb(1,1)), cell2mat(rb(1,1)), type);
p3 = strassen_mul(cell2mat(lb(1,2)), cell2mat(rb(2,1)), type);
p4 = strassen_mul(s3, s7, type);
p5 = strassen_mul(s1, s5, type);
p6 = strassen_mul(s4, cell2mat(rb(2,2)), type);
p7 = strassen_mul(cell2mat(lb(2,2)), s8, type);
t1 = p1 + p2;
t2 = t1 + p4;
CC = [p2 + p3, t1 + p5 + p6
t2 - p7, t2 + p5];
end
function c = default_mul(left, right) % классический метод умножения матриц
c = zeros(length(left));
for ii = 1:length(left)
for jj = 1:length(right)
for rr = 1:length(left)
c(ii, jj) = c(ii, jj) + left(ii, rr)* right(rr, jj);
end
end
end
end
function c = strassen_mul(left, right, type ) % типо сюда пихаем матрицу и её будущие части, а потом решаем
if length(left) == 2
c = default_mul(left, right);
else
%if type == 0
% c = strassen_mul_2x2(split_to_2x2_blocks(left), split_to_2x2_blocks(right), type);
%elseif type == 1
c = vinograd_strassen_mul_2x2(split_to_2x2_blocks(left), split_to_2x2_blocks(right), type);
%else
% disp('Не выбран вариант')
%end
end
end
%
% a = [1 2 3 4
% 5 6 7 8
% 9 10 11 12
% 13 14 15 16];
% b = [1 2 3 4
% 5 6 7 8
% 9 10 11 12
% 13 14 15 16];
% a = randi([1, 10], 8);
% b = randi([1, 10], 8);
% count_mas = [0, 0, 0]; % сложение, умножение, умножение матриц
% disp(a)
% disp(b)
% c = strassen_mul(a,b);
% disp(c)
% disp(sum)
% disp(umnoj)
% disp(umnoj_mat)
max_stepen = 6;
all_time = zeros(max_stepen, 2);
for t = 1:max_stepen
step = 2^t;
a = randi([1,10], step);
b = randi([1,10], step);
tic
c = default_mul(a, b);
all_time(t, 1) = toc;
tic
c = strassen_mul(a,b, 1);
all_time(t, 2) = toc;
end
t = 1:max_stepen;
hold on
grid on
title("График времени работы без фоновых процессов")
plot(t, all_time(:, 1),'-gs', 'LineWidth',2, 'MarkerSize',10,'Color', [0 0.4470 0.7410])
plot(t,all_time(:, 2),':gs','LineWidth',2,'MarkerSize',10, 'Color', [0.9290 0.6940 0.1250])
xlabel('Порядок матриц')
ylabel('Время работы')
legend("Классический метод", "Метод Винограда-Штрассена")
Ответы (1 шт):
Вот вариант кода. Я не тестировал его на скорость исполнения.
ma=size(A);
mb=size(B);
k=2^(ceil(log2(max([ma mb]))));
Z=zeros(k);
Ar=Z;
Ar(1:ma(1),1:ma(2))=A;
Br=Z;
Br(1:mb(1),1:mb(2))=B;
A11=Ar(1:k/2,1:k/2);
A21=Ar(k/2+1:end,1:k/2);
A12=Ar(1:k/2,k/2+1:end);
A22=Ar(k/2+1:end,k/2+1:end);
B11=Br(1:k/2,1:k/2);
B21=Br(k/2+1:end,1:k/2);
B12=Br(1:k/2,k/2+1:end);
B22=Br(k/2+1:end,k/2+1:end);
S1=A21+A22;
S2=S1-A11;
S3=A11-A21;
S4=A12-S2;
S5=B12-B11;
S6=B22-S5;
S7=B22-B12;
S8=S6-B21;
P1=S2*S6;
P2=A11*B11;
P3=A12*B21;
P4=S3*S7;
P5=S1*S5;
P6=S4*B22;
P7=A22*S8;
T1=P1+P2;
T2=T1+P4;
C=[P2+P3,T1+P5+P6;T2-P7,T2+P5];
C=C(1:ma(1),1:mb(2));