matlab - 如何加速包含 kronecker 产品的 for 循环？

Question

我想加速以下 for 循环。

% Use random matrices for testing. 
% Elapsed time of the following code is around 19 seconds. 
% So, when N is large, it's very slow.

n= 800; N =100; k =10;
x = rand(n,N); S = rand(k,N); H = 0;
for i = 1: size(x,2)
  X = x(:,i)*x(:,i)' ;
  DW = diag( S(:,i) ) - S(:,i)*S(:,i)' ;  
  H = H + kron(X,DW);
end

我的尝试:

kron(X, DW) = kron(x(:,i)*x(:,i)' ,diag(S(:,i))) - 
              kron(x(:,i)*x(:,i)', S(:,i)*S(:,i)'); 

我们可以使用 和 来重写上面的等式。

kron(x(:,i)*x(:,i)',diag(S(:,i))) = 
  kron(x(:,i), sqrt( diag(S(:,i))))* 
  kron(x(:,i)*x(:,i)',diag(S(:,i)))' ; 

(因为 S 是非负的，所以我们可以取 sqrt )

kron(x(:,i)*x(:,i)', S(:,i)*S(:,i)') = 
  kron( x(:,i), S(:,i))*
  kron( x(:,i), S(:,i))'

因此，我们只需要计算 kron(x(:,i), sqrt( diag(S(:,i))))和 kron(x(:,i), S(:,i)) .
以下是代码:

x = rand(n,N); S = rand(k,N);
H1= 0; K_D= zeros(n*k, k*1, N); K_S = zeros(n*k,N); 

%K_D records kron( x(:,i), sqrt (diag(S(:,i)) ) ), K_S records kron(x(:,i), S(:,i));

for i = 1:N
      D_half = sqrt( diag(S(:,i)));
      K_D(:,:,i) = kron( x(:,i),D_half);
      K_S(:,i) =  reshape (S(:,i)*x(:,i)',[],1);
end

K_D = reshape(K_D,n*k,[]);
H = K_D*K_D' - K_S*K_S';

新代码节省了大量时间。已用时间约为 1第二。但我还是想加快速度。
有人可以帮我加快上面的代码(我的尝试)吗？或者有人有新的想法/方法来加速我原来的问题吗？
非常感谢!

Answer 1

表达式 H1 = A(:) .* B(:).'产生与 H = kron(A,B) 相同的结果除了元素的顺序改变了。因此可以在不使用 kron 的情况下计算预期结果:

X = reshape(reshape(x, n, 1, N) .* reshape(x, 1, n, N), [], N);
S1 = -reshape(reshape(S, k, 1, N) .* reshape(S, 1, k, N), [], N);
S1(1:k+1:end, :) = reshape(S1 (1:k+1:end, :), size(S)) + S;
S1 = reshape(S1,  [], N);
H1 = X * S1.';
H = reshape(permute(reshape(permute(reshape(H1,[],k,k),[3 2 1]),k,n*k,[]),[2,1,3]),n*k,n*k);

如果您对求和的结果感兴趣，而不管元素计算的顺序如何 H1就足够了，结果超过 3X与您的第二种方法相比加速。但是您可以使用 reshape和 permute恢复带来近乎 0.3X的元素的顺序加速。
如 H在循环中计算并且它的大小在循环中是恒定的，您可以使用预先计算的索引来重新排列元素:

idx = reshape(permute(reshape(permute(reshape(1:(n*k)^2,[],k,k),[3 2 1]),k,n*k,[]),[2,1,3]),n*k,n*k);
for ....
    X = ...
    ....
    H = H1(idx);
end