Delaunay 三角剖分代码 (Cython) 中不断出现内核崩溃

Question

编辑：

经过多次测试，发现问题出在以下循环中。我将循环放在它自己的函数中。它的用途是接受指向 Vertex 结构数组的指针并删除所有重复的 Vertices。这是通过检查它们的 ID 而不是比较它们的坐标来完成的。因为每个 Vertex 一开始就被赋予了一个唯一的 id，所以它变得更加简单。

要修复的代码将简化为这样

基本结构：

ctypedef struct Vertex:
    double x
    double y
    int vid_nr

    


cdef Vertex create_vertex(double x_coord, double y_coord, int vid_number):
    cdef Vertex v
    v.x = x_coord
    v.y = y_coord
    v.vid_nr = vid_number
    return v


cdef Vertex* alloc_varr(int size, double var):
    cdef Vertex* varray = <Vertex*>malloc(size * sizeof(Vertex))
    cdef size_t i
    for i in range(size):
        varray[i] = create_vertex(var, var, i)
    return varray

基本的释放方法。已测试且功能符合预期

cdef Vertex* dealloc_from_varray(Vertex* original_array, int current_size, int remove_index):
    #careful remove_index is actual position - 1

    cdef Vertex* new_array = <Vertex*>malloc((current_size-1) * sizeof(Vertex))
    cdef size_t i
    cdef size_t j
    cdef int upper_limit = remove_index
    for i in range(upper_limit):
        new_array[i] = original_array[i]
    for j in range(upper_limit+1, current_size):
        new_array[j-1] = original_array[j]
        
    free(original_array)
    return new_array

现在通过 vid_nr 取消分配重复顶点的方法

cdef void dealloc_non_uniques(Vertex* free_varray, int size):
    cdef int k = 0 # Initialize k
    cdef int l = 1 #initialized at first iteration but to be safe
    cdef int vsize = size
    #dealloc all non-unique verrtices from the free arr
    while k < vsize - 1:
        l = k + 1
        while l < vsize:
            if free_varray[l].vid_nr == free_varray[k].vid_nr:
                free_varray = dealloc_from_varray(free_varray, vsize, l)
                vsize = vsize - 1
            else:
                l = l + 1
        k = k + 1
    size = vsize

下面简单测试一下

cdef int varr_size = 7
cdef Vertex* og_array = alloc_varr(varr_size, 0)

og_array[0].x = 0.2
og_array[0].y = 2
og_array[1].x = -1.8
og_array[1].y = 0.15
og_array[2].x = -0.1
og_array[2].y = -2.3
og_array[3].x = 2
og_array[3].y = -0.3
og_array[4].x = -1.8
og_array[4].y = 6
og_array[5].x = -3.9
og_array[5].y = -2
og_array[6].x = 1.7
og_array[6].y = -3.4
og_array[6].vid_nr = og_array[5].vid_nr
og_array[2].vid_nr = og_array[4].vid_nr
og_array[3].vid_nr = og_array[4].vid_nr




print(varr_size)
for i in range(varr_size):
    print(og_array[i])


#print("\ntest dealloc_non_uniques:\nwe expect varr_size to be 4")

#dealloc_non_uniques(og_array, varr_size)

一旦我删除注释 # 并运行

dealloc_non_uniques

 方法，

将使内核崩溃。为什么？

Answer 1

经过很长一段时间我终于发现问题所在了。 首先，声明的变量是不可变的。因此，当我将整个网格的大小声明为全局大小并稍后出于绘图原因调用它时，我无法使用

triangulate

函数即时更新它。这很容易解决。我可以在内存中分配适当的大小，例如

<int*>global_mesh_size = ....

，然后通过将

global_mesh_size[0]

等于我在函数本身内部创建的临时 mesh_size

int

来在三角函数内修改它。

我在这里发布的第一个版本中的第二个错误是循环的编写方式。由于最终的网格大小未知，因此当我迭代其元素并打算删除时，我必须确保动态调整循环限制。 这无法在

for range()

循环中完成，因为

range()

接受初始输入并创建在循环结束之前保持静态的文字范围。 因此，有必要在函数中手动更新范围，并使用带有简单中断条件的

while

循环。这使我能够从 Final_mesh 中释放坏三角形，而不会触及越界元素。

第三个错误是

sort_by_angle()

函数中的一个非常简单的语法错误，其中我不小心使用 .x 而不是 .y 坐标来计算临时向量幅度，然后再进行范数和 x 轴上的投影。 这导致大多数时候顶点没有被正确排序。这是发现的最耗时的错误，因为有很多行代码，而且我最终添加了很多

print()

语句来查看代数是否正确完成。

以下测试单元的总体情况

cdef int varr_ssize = 4
cdef Vertex* og_array = alloc_varr(varr_ssize, 0)
cdef int* global_mesh_size = <int*>malloc(sizeof(int))
global_mesh_size[0] = 0 #dynamically changed through triangulate fx

og_array[0].x = 0.2
og_array[0].y = 2
og_array[1].x = -1.8
og_array[1].y = 0.15
og_array[2].x = -0.1
og_array[2].y = -2.3
og_array[3].x = 2
og_array[3].y = -0.3


cdef Face* testmesh = triangulate(og_array, varr_ssize, global_mesh_size)

print("\nglobal mesh size is:", global_mesh_size[0])
for f in range(global_mesh_size[0]):
    print(testmesh[f])

free(testmesh)
free(og_array)
free(global_mesh_size)

给出了

的正确三角测量

global mesh size is: 2

{'ver1': {'x': 2.0, 'y': -0.3, 'vid_nr': 3}, 'ver2': {'x': 0.2, 'y': 2.0, 'vid_nr': 0}, 'ver3': {'x': -1.8, 'y': 0.15, 'vid_nr': 1}, 'fid_nr': 0}
{'ver1': {'x': 2.0, 'y': -0.3, 'vid_nr': 3}, 'ver2': {'x': -1.8, 'y': 0.15, 'vid_nr': 1}, 'ver3': {'x': -0.1, 'y': -2.3, 'vid_nr': 2}, 'fid_nr': 1}

为了避免该线程混乱，我只会将 Github 页面链接到代码，因为它太长，无法在论坛布局中舒适地阅读。 像往常一样，只有笔记本中的前几个单元是相关的。第一个输出以下的任何内容仅用于测试目的

通过 Bowyer-Watson 算法进行 Delaunay 三角剖分

感谢所有帮助过的人<3