在numpy数组上使用\uu contains进行非常慢的迭代__

Question

我有一个nx3numpy数组，列出了一个曲面的N个三角形面的3个顶点，还有一个Nx1数组，其中的每个面对应一个值。你知道吗

我想将这些“面”值转换为“顶点”值（尽可能最好），例如通过查找与顶点关联的所有面的平均值。你知道吗

我目前的解决方案适用于较小的N值，但缩放为面数x顶点数，这很快变得不切实际：

def face2vertVal(faces, facesVals, verts):
    # INPUT:
    # faces:     Nx3 array of N vertex IDs for N faces
    # facesVals: Nx1 array of some parameter for each face

    # OUTPUT:
    # vertsVals: Px1 array of the mean value in "facesVals" that corresponds
    #            to each vertex

    import numpy as np

    vertsVals = np.zeros(faces.max()+1)

    for vertex in range(0, faces.max()+1):       
        tmpVals = np.zeros(1)
        for row in range(0, faces.shape[0]):
                if faces[row].__contains__(vertex):
                    tmpVals = np.append(tmpVals, facesVals[row])
        vertsVals[vertex] = tmpVals[1:].mean()
        del tmpVals

    return vertsVals

提前谢谢。你知道吗

编辑：矢量化的方法非常快，但是对于700k个面和350k个顶点需要太多的内存。你知道吗

Answer 1

您的代码有一些对性能有很大影响的问题（例如，不使用np.附加在循环中）。 因此，第一步是改进代码，同时避免不必要的搜索。你知道吗

在下一步中，我们可以使用jit编译器来获得额外的性能。你知道吗

代码

import numpy as np
import numba as nb
import time
from scipy import spatial
N_verts=350000

#This gernerates a mesh for performance Testing
#https://stackoverflow.com/a/50579387/4045774
def make(N):
    verts= np.random.uniform(-10, 10, (N, 3))
    faces = spatial.Delaunay(verts[:, :2]).simplices
    return verts,faces

@nb.njit()
def face2vertVal(faces, facesVals):
  # INPUT:
  # faces:     Nx3 array of N vertex IDs for N faces
  # facesVals: Nx1 array of some parameter for each face

  # OUTPUT:
  # vertsVals: Px1 array of the mean value in "facesVals" that corresponds
  #            to each vertex

  vertsVals = np.zeros(faces.max()+1)
  verts_counts= np.zeros(faces.max()+1)

  for i in range(faces.shape[0]):
    for j in range(faces.shape[1]):
      vertsVals[faces[i,j]]+=facesVals[i]
      verts_counts[faces[i,j]]+= 1.

  vertsVals=vertsVals/verts_counts
  return vertsVals

[verts,faces]=make(N_verts)
facesVals=np.random.rand(faces.shape[0])
res=face2vertVal(faces, facesVals)

性能

N_verts=350 000
Pure-Python:2.12s
With Numba (after the compilation which takes about 0.5s at the first call): 18.7ms

Answer 2

编辑：

如果内存有问题，也可以使用“批处理”版本的算法，一次只处理有限数量的面：

将numpy作为np导入

def face2vertVal(faces, faces_values, batch_size=None):
    batch_size = batch_size or len(faces)
    faces = np.asarray(faces)
    faces_values = np.asarray(faces_values)
    vert_idx = np.arange(faces.max() + 1)
    vertex_values = np.zeros(len(vert_idx))
    vertex_counts = np.zeros(len(vert_idx))
    for i in range(0, len(faces), batch_size):
        faces_batch = faces[i:i + batch_size]
        faces_values_batch = faces_values[i:i + batch_size]
        vertex_faces = np.any(faces_batch == vert_idx[:, np.newaxis, np.newaxis], axis=-1)
        vertex_values += np.sum(vertex_faces * faces_values_batch, axis=1)
        vertex_counts += np.count_nonzero(vertex_faces, axis=1)
    return vertex_values / np.maximum(vertex_counts, 0)

# All at once
vertex_values = face2vertVal([[0, 1, 2], [1, 2, 3], [2, 3, 0]], [1, 2, 3])
print(*vertex_values, sep=', ')
# 2.0, 1.5, 2.0, 2.5

# In batches of two
vertex_values_1b = face2vertVal([[0, 1, 2], [1, 2, 3], [2, 3, 0]], [1, 2, 3], 2)
print(*vertex_values_1b, sep=', ')
# 2.0, 1.5, 2.0, 2.5

您可以手动选择batch_size（或者根据要使用的内存或其他内容使用一些公式）来平衡速度和内存之间的平衡。你知道吗

---

你应该以矢量化的方式进行计算，否则会慢得多。这是一种方法：

import numpy as np

def face2vertVal(faces, faces_values):
    faces = np.asarray(faces)
    faces_values = np.asarray(faces_values)
    vert_idx = np.arange(faces.max() + 1)
    vertex_faces = np.any(faces == vert_idx[:, np.newaxis, np.newaxis], axis=-1)
    vertex_values = np.sum(vertex_faces * faces_values, axis=1) / np.maximum(np.count_nonzero(vertex_faces, axis=1), 0)
    return vertex_values

vertex_values = face2vertVal([[0, 1, 2], [1, 2, 3], [2, 3, 0]], [1, 2, 3])
print(*vertex_values, sep=', ')
# 2.0, 1.5, 2.0, 2.5

请注意，这种方法的缺点是，给定具有P顶点的N面，它需要O(N*P)内存，而非矢量化算法使用O(max(N, P))内存。你知道吗