matlab - 如何向量化依赖的 For 循环

Question

我正在开发一个函数，该函数接受 1xn 向量 x 作为输入并返回 nxn 矩阵 L。
我想通过向量化循环来加快速度，但有一个问题让我困惑:循环索引 b 取决于循环索引 a。任何帮助将不胜感激。

x = x(:); 
n = length(x);
L = zeros(n, n);
for a = 1 : n,
    for b = 1 : a-1,
        c = b+1 : a-1;
        if all(x(c)' < x(b) + (x(a) - x(b)) * ((b - c)/(b-a))),
            L(a,b) = 1;
        end
    end
end

Answer 1

从快速测试来看，您似乎只对下三角形进行了一些操作。您也许可以使用与此类似的 ind2sub 和 arrayfun 等丑陋技巧进行矢量化

tril_lin_idx = find(tril(ones(n), -1));
[A, B] = ind2sub([n,n], tril_lin_idx);
C = arrayfun(@(a,b) b+1 : a-1, A, B, 'uniformoutput', false); %cell array
f = @(a,b,c) all(x(c{:})' < x(b) + (x(a) - x(b)) * ((b - c{:})/(b-a)));
L = zeros(n, n);
L(tril_lin_idx) = arrayfun(f, A, B, C);

我无法测试它，因为我没有 x 并且我不知道预期的结果。我通常喜欢矢量化解决方案，但这可能有点太过了:)。我会坚持使用你的显式 for 循环，这可能会更清晰，并且 Matlab 的 JIT 应该能够轻松加速。您可以将 if 替换为 L(a,b) = all(...)。

编辑1

更新版本，以防止浪费 C 上的 ~ n^3 空间:

tril_lin_idx = find(tril(ones(n), -1));
[A, B] = ind2sub([n,n], tril_lin_idx);
c = @(a,b) b+1 : a-1;
f = @(a,b) all(x(c(a, b))' < x(b) + (x(a) - x(b)) * ((b - c(a, b))/(b-a)));
L = zeros(n, n);
L(tril_lin_idx) = arrayfun(f, A, B);

编辑2

轻微变体，不使用 ind2sub，并且应该更容易修改，以防 b 以更复杂的方式依赖于 a。我内联了 c 以提高速度，似乎特别是调用函数句柄的成本很高。

[A,B] = ndgrid(1:n);
v = B<A; % which elements to evaluate
f = @(a,b) all(x(b+1:a-1)' < x(b) + (x(a) - x(b)) * ((b - (b+1:a-1))/(b-a)));
L = false(n);
L(v) = arrayfun(f, A(v), B(v));