矢量化矩阵乘法

我有一个关于在Octave上运行SOR算法的学校项目，但是我的效率很低。所以我有这段代码：矢量化矩阵乘法

for ii=1:n 
    r = 1/A(ii,ii); 
    for jj=1:n 
     if (ii!=jj) 
      A(ii,jj) = A(ii,jj)*r; 
     end; 
    end; 
    b(ii,1) = b(ii,1)*r; 
    x(ii,1) = b(ii,1); 
end;

我该如何实现这一点？我的第一次尝试是这样的：

for ii=1:n 
    r = 1/A(ii,ii); 
    A(find(eye(length(A))!=1)) = A(find(eye(length(A))!=1))*r; 
    b(ii,1) = b(ii,1)*r; 
    x(ii,1) = b(ii,1); 
end;

但我不确定它有多大的帮助。有没有更好的和/或更有效的方法来做到这一点？

谢谢！

来源

2013-04-15 dccarmo

你能解释代码实际在做什么吗？（一些矩阵分解可能......） – tashuhka

它创建了SOR算法的迭代矩阵。 – dccarmo

发布了一个基于bsxfun的解决方案[here]（http：// stackoverflow。com/a/26128876/3293881）如果你想看看！ – Divakar

你完全可以避免循环我相信。你必须看到你正在采取A的对角元素的倒数，那么你为什么要使用循环。直接做。这是第一步。删除Inf，现在您想将r乘以相应行的相应非对角元素，对吗？因为您将与（1,2），（1,3），...，（1，n）相同的r相乘，因此使用repmat来构造这样的矩阵，其将沿着列具有复制的元素。但是R具有非零对角线元素。因此使它们为零。现在你会得到你的A，但对角元素将为零。因此，您只需将其从原始A中添加回来即可。这可以通过A=A.*R+A.*eye(size(A,1))完成。

矢量化来自经验，最重要的是分析你的代码。想想在每一步是否要使用循环，如果不是用等效命令替换该步骤，其他代码将会跟随（例如，我构建了一个矩阵R，而您正在构建单个元素r。所以我只考虑转换r - >R然后剩下的代码就会落到位）。

的代码如下：

R=1./(A.*eye(size(A,1))); %assuming matrix A is square and it does not contain 0 on the main diagonal 
R=R(~isinf(R)); 
R=R(:); 
R1=R; 
R=repmat(R,[1,size(A,2)]); 
R=R.*(true(size(A,1))-eye(size(A,1))); 
A=A.*R+A.*eye(size(A,1)); %code verified till here since A comes out to be the same 

b = b.*R1; 
x=b;

来源

2013-04-16 00:35:44

太棒了，Parag！它确实有用！你有什么资源可以学习更多关于矢量化的知识吗？我很难说明你做了什么，并且还想将我在代码中的其他循环向量化！ – dccarmo

@dccarmo我刚查过。 'b'不出来是正确的。你能明白为什么，因为你比我更了解代码。或者我犯了一些错误。我会尝试通过编辑我的答案来告诉逻辑 –

我刚刚在这里检查，而且b和原始代码一样，我检查它是否错误？ – dccarmo

我假设有矩阵：

A (NxN) 
b (Nx1)

的代码：

d = diag(A); 
A = diag(1 ./ d) * A + diag(d - 1); 
b = b ./ d; 
x = b;

来源

2013-04-16 17:58:27 Serg

随机碰着到这个问题和第一眼看起来考虑到这个问题被标记为vectorization问题，所以感兴趣。

我能够想出一个bsxfun的矢量化解决方案，也使用diagonal indexing。这个解决方案似乎给了我一个3-4x超过循环代码与体面大小的输入。

假设你对于加速这个问题的速度有所提高仍然有点感兴趣，那么我会非常想知道你将会获得哪些加速。这里的代码 -

diag_ind = 1:size(A,1)+1:numel(A); 
diag_A = A(diag_ind(:)); 
A = bsxfun(@rdivide,A,diag_A); 
A(diag_ind) = diag_A; 
b(:,1) = b(:,1)./diag_A; 
x(:,1) = b(:,1);

让我知道！

来源

2014-09-30 19:47:04 Divakar

矢量化矩阵乘法

回答

相关问题