使用Numpy快速生成网格三角形网格

Question

考虑一个常规矩阵，代表如图所示编号的节点：

我想要列出图中所代表的所有三角形，并将其结果保存在以下二维列表中：[[0,1,4],[1,5,4],[1,2,5],[2,6,5],...,[11,15,14]] 假设矩阵的尺寸是（NrXNc）（在这种情况下为4X4），我可以通过以下代码实现此结果：

def MakeFaces(Nr,Nc):
    Nfaces=(Nr-1)*(Nc-1)*2
    Faces=np.zeros((Nfaces,3),dtype=np.int32)
    for r in range(Nr-1):
        for c in range(Nc-1):
            fi=(r*(Nc-1)+c)*2
            l1=r*Nc+c
            l2=l1+1
            l3=l1+Nc
            l4=l3+1
            Faces[fi]=[l1,l2,l3]
            Faces[fi+1]=[l2,l4,l3]
    return Faces

然而，双层循环操作使得这种方法非常缓慢。有没有一种聪明的方式利用numpy来更快地完成这个任务？

- Miguel

Answer 1

我们可以在NumPy环境中使用完美的slicing和multi-dim assignment来玩一个基于multi-dimensional的游戏 -

def MakeFacesVectorized1(Nr,Nc):

    out = np.empty((Nr-1,Nc-1,2,3),dtype=int)

    r = np.arange(Nr*Nc).reshape(Nr,Nc)

    out[:,:, 0,0] = r[:-1,:-1]
    out[:,:, 1,0] = r[:-1,1:]
    out[:,:, 0,1] = r[:-1,1:]

    out[:,:, 1,1] = r[1:,1:]
    out[:,:, :,2] = r[1:,:-1,None]

    out.shape =(-1,3)
    return out

运行时测试和验证 -

In [226]: Nr,Nc = 100, 100

In [227]: np.allclose(MakeFaces(Nr, Nc), MakeFacesVectorized1(Nr, Nc))
Out[227]: True

In [228]: %timeit MakeFaces(Nr, Nc)
100 loops, best of 3: 11.9 ms per loop

In [229]: %timeit MakeFacesVectorized1(Nr, Nc)
10000 loops, best of 3: 133 µs per loop

In [230]: 11900/133.0
Out[230]: 89.47368421052632

对于 Nr, Nc = 100, 100，速度提高了90x！

- Divakar

Answer 2

如果你正确转换问题，就可以在没有任何显式循环的情况下实现类似的结果。一种方法是将结果想象成三个数组，每个数组包含一个顶点：第一、第二和第三个。然后，你可以通过zip或其他方式将数组转换为任何格式，在相对便宜的操作中完成。

你从实际矩阵开始。这将使索引和选择元素变得更加容易：

m = np.arange(Nr * Nc).reshape(Nr, Nc)

第一个数组将包含所有90度的拐角：

c1 = np.concatenate((m[:-1, :-1].ravel(), m[1:, 1:].ravel()))

m[:-1, :-1] 是位于顶部的角落，m[1:, 1:] 是位于底部的角落。

第二个数组将包含相应的顶部锐角：

c2 = np.concatenate((m[:-1, 1:].ravel(), m[:-1, 1:].ravel()))

第三个数组将包含底部角落：

c2 = np.concatenate((m[1:, :-1].ravel(), m[1:, :-1].ravel()))

现在您可以通过压缩来获得与原始数组类似的数组：

faces = list(zip(c1, c2, c3))

我相信你可以找到改进这个算法的方法，但这只是一个开始。

- Mad Physicist