从 3 个向量及其原始坐标生成所有可能的三元组对组合

Generate all possible combinations of pairs of triplets from 3 vectors together with their original coordinates

给定 Matlab 中的 3 个行向量 A,B,C,我想生成矩阵 D 报告来自 A,B,C 的三元组对的所有可能组合及其在 A,B,C

我已经编写了一个代码来完成我想要的。由于我正在尝试尽可能多地优化我的代码(代码应该重复数百万次),我想知道您是否可以想出更有效的解决方案。例如,在我的代码中,我没有预先分配矩阵 D,因为我不知道如何获取每对三元组的索引,这样效率不高。

下面的代码会更好地解释它:

clear 
A=[1 2];
B=[-4 -2 5];
C=[8 9 -3 0];

sA=size(A,2);
sB=size(B,2);
sC=size(C,2);
sT=sA*sB*sC;

%Generate the matrix D of dimension [sT*(sT-1)/2]x[12]
%reporting all the possible combinations of pairs of triplets from A,B,C
%together with their original coordinates in A,B,C

[ca, cb, cc] = ndgrid(A, B, C);
T = [ca(:), cb(:), cc(:)];  %matrix of dimension sTx3 reporting all the possible triplets 
                            %from A,B,C

[ca, cb, cc] = ndgrid(1:sA, 1:sB, 1:sC);
Tcoord = [ca(:), cb(:), cc(:)];  %matrix of dimension sTx3 reporting the coordinates of all 
                                 %the possible triplets from A,B,C

D=[];
for w=1:sA*sB*sC
    for r=w+1:sA*sB*sC 
        D=[D; T(w,:) T(r,:) Tcoord(w,:) Tcoord(r,:)];
    end
end

填充矩阵 D 的最后一个嵌套 for 循环效率更高。 OP 在他们的声明中是正确的:

"For example, in my code I am not pre-allocating the matrix D because I don't know how to get the index of each pair of triplets and this is not efficient."

我们可以对这些循环中的大部分工作进行向量化,方法是注意到 OP 在他们关于矩阵最终大小的评论中提到了一种模式 D(即 Generate the matrix D of dimension [sT*(sT-1)/2]x[12]) .任何熟悉系列和序列的人都会熟悉第一个维度。就是Triangle Numbers.

的公式

考虑到这一点,我们可以看到在最终结果中,前 3 列(以及第 7 列到第 9 列)重复了 23 次,然后是 22 次,依此类推,而第 4 列到第 6 列(以及第 10 列)到 12) 是 T/Tcoord 的最后 23 行,T/Tcoord 的最后 22 行等等。在代码中我们有:

D1 = zeros(sT * (sT - 1) / 2, 12);
s = 1;
e = sT - 1;

for w = 1:(sT - 1)
    D1(s:e,[1:3,7:9]) = repmat([T(w,:),Tcoord(w,:)], sT - w, 1);
    D1(s:e,[4:6,10:12]) = [T((w+1):sT,:),Tcoord((w+1):sT,:)];
    s = e + 1;
    e = e + (sT - (w + 1));
end

和运行每个方法使用tictoc 200次,我们发现我们的效率几乎提高了35%

% OP's setup code goes here

tic
for i=1:200
    D=[];
    for w=1:sA*sB*sC
        for r=w+1:sA*sB*sC
            D=[D; T(w,:) T(r,:) Tcoord(w,:) Tcoord(r,:)];
        end
    end
end
toc

tic
for i = 1:200
    D1 = zeros(sT * (sT - 1) / 2, 12);
    s = 1;
    e = sT - 1;

    for w = 1:(sT - 1)
        D1(s:e,[1:3,7:9]) = repmat([T(w,:),Tcoord(w,:)], sT - w, 1);
        D1(s:e,[4:6,10:12]) = [T((w+1):sT,:),Tcoord((w+1):sT,:)];
        s = e + 1;
        e = e + (sT - (w + 1));
    end
end
toc

% Gives same result
isequal(D, D1)

% Timing for 200 runs on 24 total combinations
Elapsed time is 2.09613 seconds.
Elapsed time is 1.35988 seconds.
ans = 1

如果我们使输入向量更大,我们会看到效率有更大的提高。以下是 运行 在以下向量上对每种方法进行 15 次的结果:

A=[1 2 3 4 23];
B=[-4 -2 5 74];
C=[8 9 -3 0];

% Timing for 15 run on 80 total combinations
Elapsed time is 4.00448 seconds.
Elapsed time is 0.379919 seconds.
ans = 1

快了 10 倍多。随着输入向量的大小变大,差距呈指数增长。

A=[1 2 3 4 23];
B=[-4 -2 5 74 28];
C=[8 9 -3 0 -100 -5];

% Timing for 1 run on 150 total combinations
Elapsed time is 3.63065 seconds.
Elapsed time is 0.0481789 seconds.
ans = 1

大约快了 75 倍!!!


更新

OP 在评论中给出了更好的答案:

indices=nchoosek((1:1:sT),2);
D=[T(indices(:,1),:) T(indices(:,2),:) Tcoord(indices(:,1),:) Tcoord(indices(:,2),:)];

这是我用来做基准测试的代码:

clear 
A=[1 2 3 4 23 24 25 26];
B=[-4 -2 5 74 28 10 11 12 13];
C=[8 9 -3 0 -100 -5 60 120];

sA=size(A,2);
sB=size(B,2);
sC=size(C,2);
sT=sA*sB*sC;

tic
for i = 1:10
    [ca, cb, cc] = ndgrid(A, B, C);
    T = [ca(:), cb(:), cc(:)];
    [ca, cb, cc] = ndgrid(1:sA, 1:sB, 1:sC);
    Tcoord = [ca(:), cb(:), cc(:)];

    D1 = zeros(sT * (sT - 1) / 2, 12);
    s = 1;
    e = sT - 1;

    for w = 1:(sT - 1)
        D1(s:e,[1:3,7:9]) = repmat([T(w,:),Tcoord(w,:)], sT - w, 1);
        D1(s:e,[4:6,10:12]) = [T((w+1):sT,:),Tcoord((w+1):sT,:)];
        s = e + 1;
        e = e + (sT - (w + 1));
    end
end
toc

tic
for i = 1:10
    indices=nchoosek((1:1:sT),2);
    D=[T(indices(:,1),:) T(indices(:,2),:) Tcoord(indices(:,1),:) Tcoord(indices(:,2),:)];
end
toc

isequal(D, D1)

结果如下:

% Timing for 10 runs on 576 total combinations
Elapsed time is 1.9834 seconds.
Elapsed time is 0.13818 seconds.
ans = 1

我提供的改进解决方案比原来的解决方案好很多,但与 OP 更新的解决方案不匹配。它的订单速度更快,非常优雅,我可能会补充。