MATLAB中矩阵二次方程的向量化

Question

我正在尝试在 Matlab 中“矢量化”这个循环以提高计算效率

for t=1:T
   j=1;
   for m=1:M
      for n=1:N
         y(t,j) = v{m,n} + data(t,:)*b{m,n} +  data(t,:)*f{m,n}*data(t,:)';
         j=j+1;
      end
   end
end

其中v 是一个 (M x N) 标量单元。 b 是 (K x 1) 个向量的 (M x N) 单元。 f 是 (K x K) 矩阵的 (M x N) 单元格。 data 是一个 (T x K) 数组。

举个例子来说明我用来向量化同一个循环而没有二次项的代码是：

B = [reshape(cell2mat(v)',1,N*M);cell2mat(reshape(b'),1,M*N)];
X = [ones(T,1),data];
y = X*B;

谢谢！

Answer 1

这是针对性能的最矢量化形式 -

% Extract as multi-dim arrays
vA = reshape([v{:}],M,N);
bA = reshape([b{:}],K,M,N);
fA = reshape([f{:}],K,K,M,N);

% Perform : data(t,:)*f{m,n} for all iterations
data_f_mult = reshape(data*reshape(fA,K,[]),T,K,M,N);

% Now there are three parts :
% v{m,n}
% data(t,:)*b{m,n}
% data(t,:)*f{m,n}*data(t,:)';

% Compute those parts one by one
parte1 = vA(:).';
parte2 = data*reshape(bA,[],M*N);

parte3 = zeros(T,M*N);
for t = 1:T
    parte3(t,:) = data(t,:)*reshape(data_f_mult(t,:,:),K,[]);
end

% Finally sum those up and to present in desired format permute dims
sums = bsxfun(@plus, parte1, parte2 + parte3);
out = reshape(permute(reshape(sums,T,M,N),[1,3,2]),[],M*N);

Answer 2

对于那些感兴趣的人，我找到了解决方案

f = f';
tMat = blkdiag(f{:})+(blkdiag(f{:}))';
y2BB = [reshape(cell2mat(v)',1,N*M);...
        cell2mat(reshape(b',1,M*N));...
        reshape(diag(blkdiag(f{:})),K,N*M);...
        reshape(tMat((tril(tMat,-1)~=0)),sum(1:K-1),M*N)]; 
y2YBar = [ones(T,1),data,data.^2];

jj=1;
kk=1;
ll=1;
for k=1:sum(1:K-1)
    y2YBar = [y2YBar,data(:,jj).*data(:,kk+jj)];
    if kk<(K-ll)
        kk=kk+1;
    else
        kk=1;
        jj=jj+1;
        ll=ll+1;
    end
end
y = y2YBar*y2BB;