векторизация и Уплотненное Умножение матриц

Question

Вот исходный код:

K = zeros(N*N) 
for a=1:N 
    for i=1:I 
     for j=1:J 
      M = kron(X(:,:,a).',Y(:,:,a,i,j)); 

      %A function that essentially adds M to K. 
     end 
    end 
end 

Цель состоит в том, чтобы векторизации вызовов kroniker умножения. Моя интуиция состоит в том, чтобы рассматривать X и Y как контейнеры матриц (для справки, ломтики X и Y, подаваемые в kron, являются квадратными матрицами порядка 7x7). В этой контейнерной схеме X появляется в виде контейнера с 1-D и Y в виде 3-D контейнера. Следующей догадкой я должен был преобразовать Y в двумерный контейнер или, еще лучше, в одномерный контейнер, а затем сделать элементное умножение X и Y. Вопросы: как бы это изменить, таким образом, чтобы сохранить след M и может ли Matlab справиться с этой идеей в этой идее контейнера или сделать контейнеры еще более переделанными, чтобы разоблачить внутренние матричные элементы?

Answer 1

Подход № 1: Умножение матриц с 6D перестановочны

% Get sizes 
[m1,m2,~] = size(X); 
[n1,n2,N,n4,n5] = size(Y); 

% Lose the third dim from X and Y with matrix-multiplication 
parte1 = reshape(permute(Y,[1,2,4,5,3]),[],N)*reshape(X,[],N).'; 

% Rearrange the leftover dims to bring kron format 
parte2 = reshape(parte1,[n1,n2,I,J,m1,m2]); 

% Lose dims correspinding to last two dims coming in from Y corresponding 
% to the iterative summation as suggested in the question 
out = reshape(permute(sum(sum(parte2,3),4),[1,6,2,5,3,4]),m1*n1,m2*n2) 

Подход № 2: Простой 7D переставлять

% Get sizes 
[m1,m2,~] = size(X); 
[n1,n2,N,n4,n5] = size(Y); 

% Perform kron format elementwise multiplication betwen the first two dims 
% of X and Y, keeping the third dim aligned and "pushing out" leftover dims 
% from Y to the back 
mults = bsxfun(@times,permute(X,[4,2,5,1,3]),permute(Y,[1,6,2,7,3,4,5])); 

% Lose the two dims with summation reduction for final output 
out = sum(reshape(mults,m1*n1,m2*n2,[]),3); 

---

Проверка

Вот настройки для запуска оригинала и предлагаемые подходы -

% Setup inputs 
X = rand(10,10,10); 
Y = rand(10,10,10,10,10); 

% Original approach 
[n1,n2,N,I,J] = size(Y); 
K = zeros(100); 
for a=1:N 
    for i=1:I 
     for j=1:J 
      M = kron(X(:,:,a).',Y(:,:,a,i,j)); 
      K = K + M; 
     end 
    end 
end 

% Approach #1 
[m1,m2,~] = size(X); 
[n1,n2,N,n4,n5] = size(Y); 
mults = bsxfun(@times,permute(X,[4,2,5,1,3]),permute(Y,[1,6,2,7,3,4,5])); 
out1 = sum(reshape(mults,m1*n1,m2*n2,[]),3); 

% Approach #2 
[m1,m2,~] = size(X); 
[n1,n2,N,n4,n5] = size(Y); 
parte1 = reshape(permute(Y,[1,2,4,5,3]),[],N)*reshape(X,[],N).'; 
parte2 = reshape(parte1,[n1,n2,I,J,m1,m2]); 
out2 = reshape(permute(sum(sum(parte2,3),4),[1,6,2,5,3,4]),m1*n1,m2*n2); 

После запуска мы видим макс. абсолютное отклонение с предлагаемыми подходами к первоначальному -

>> error_app1 = max(abs(K(:)-out1(:))) 
error_app1 = 
    1.1369e-12 
>> error_app2 = max(abs(K(:)-out2(:))) 
error_app2 = 
    1.1937e-12 

Ценности выглядят хорошо для меня!

---

Бенчмаркинг

Timing эти три подхода, используя тот же большой набор данных, который используется для проверки, мы получаем что-то вроде этого -

----------------------------- With Loop 
Elapsed time is 1.541443 seconds. 
----------------------------- With BSXFUN 
Elapsed time is 1.283935 seconds. 
----------------------------- With MATRIX-MULTIPLICATION 
Elapsed time is 0.164312 seconds. 

Похоже матрицы-умножения делает довольно хорошо для набора данных этих размеров!