Matlab中重复元素矩阵的优化矩阵运算答案

【问题标题】：Optimized matrix operation of matrix with repeating elements in MatlabMatlab中重复元素矩阵的优化矩阵运算
【发布时间】：2015-09-02 16:57:03
【问题描述】：

我想获得一个大矩阵 (A) 的 exp()，其值在不同的索引处重复。为了加快 exp() 操作，我只对 A 的唯一值执行它，然后重新组装矩阵。然而，矩阵的重新组装非常缓慢。以下代码提供了一个工作示例：

% defintion of a grid
gridSp = 5:5:35*5;

X = repmat(gridSp,35,1);
Z = repmat(gridSp',1,35);

% calculation of the distances
locMat = [X(:) Z(:)];
dist=sqrt(bsxfun(@minus,locMat(:,1),locMat(:,1)').^2 +...
bsxfun(@minus,locMat(:,2),locMat(:,2)').^2); 

sizeDist = size(dist);
uniqueDist = unique(dist,'stable');
[~, Locb] = ismember(dist,uniqueDist);

nn_A = exp(1i*2*pi*rand(sizeDist(1),100));
H_A = zeros(size(nn_A));

freq = linspace(10^-3,10,100);

psdA = 4096*length(freq).*10.*4.*22.6./((1 + 6.*freq*22.6).^(5/3));

for jj=1:100 
    b = exp(-8.8*uniqueDist*sqrt((freq(jj)/15).^2 + 10^-7));
    b = b.*psdA(jj);
    A = b(Locb);
    droptol = max(A(:))*10^-10;
    if min(A(:))<droptol
        A = sparse(A);
        HH_A = ichol(A,struct('type','ict','shape','lower','droptol',droptol));
    else
        HH_A = chol(A,'lower');
    end
    H_A(:,jj) = HH_A*nn_A(:,jj);
end

尤其是矩阵的重组

A = b(Locb);

以及矩阵到稀疏的转换

A = sparse(A);

在最后一个for循环中占用了很多时间。有没有更快的方法来做到这一点？有趣的是：

B = A + A;

比

快很多

A = b(Locb);

我必须比示例中的 100 次迭代更频繁地执行这些操作。

这里是根据要求提供的代码的精简版本（如下）。

% defintion of a grid
gridSp = 5:5:28*5;

X = repmat(gridSp,35,1);
Z = repmat(gridSp',1,35);

% calculation of the distances
locMat = [X(:) Z(:)];
dist=sqrt(bsxfun(@minus,locMat(:,1),locMat(:,1)').^2 +bsxfun(@minus,locMat(:,2),locMat(:,2)').^2); 

uniqueDist = unique(dist,'stable');
[~, Locb] = ismember(dist,uniqueDist);

for jj=1:100 
    b = exp(jj.*uniqueDist);
    A = b(Locb);
end

【问题讨论】：

希望这里有人可以提供一些提示，但我最初的反应是您的 matlab 代码非常好，而您要获得显着性能改进的唯一方法是编写一个 c 程序来完成它并将其编译为 mex 文件
您的问题看起来很有趣。但是，您的代码似乎比您在文本中所说的要多得多。所以我有两个问题：你确定exp-and-reassembly 是瓶颈吗？如果是的话，你能只提供那部分更短的代码吗？
是的，我相信他们的重组是瓶颈。缩短的代码发布在上述问题的编辑中。在示例中 exp 实际上比重新组装更快。在原始代码中，重组速度更快，因为 exp 应用于不同的值范围。
我相信您的算法可能存在数值问题。 dist*sqrt((freq(jj)/15).^2 + 10^-7) 假定值最大为 160，但 exp(-8.8*x) 会出现从 x 低至 85 的值开始的数值下溢。由于这种下溢，结果矩阵只是稀疏的。

标签： matlab matrix sparse-matrix

【解决方案1】：

在您的示例中， dist 的尺寸仅为 980 x 980 在这种情况下，您最好只执行密集矩阵运算，即

for jj=1:100
   A=exp(jj*dist);
end

比

快2倍

for jj=1:100 
    b = exp(jj.*uniqueDist);
    A = b(Locb);
end

对于你给定的例子。

【讨论】：