在 C 中使用递归的堆栈溢出
Stack-overflow using recursion in C
我一直在疯狂地寻找答案无济于事；我创建了一个四叉树，它应该对结构声明的超过 1000 个对象的数组进行排序：
typedef struct node
{
char is_leaf;
struct Particle *p;
double m;
double center_x;
double center_y;
double width;
struct node *sw;
struct node *nw;
struct node *se;
struct node *ne;
} node;
使用以下函数生成四叉树：
node* quadtree_insert(node *n, struct Particle *p, double center_x, double center_y, double width)
{
if(n == NULL)
{
n = (node*)malloc(sizeof(node));
n->is_leaf = 1;
n->p = p;
//n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
if(width < 1e-300){
n->width = 1e-300;
}
else
n->width = width;
return n;
}
else{
//n = (node*)malloc(sizeof(node));
double x;
double y;
if(width < 1e-300){
width = 1e-300;
}
if(n->is_leaf == 1) //! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if(x <= center_x && y <= center_y) //! first quadrant
{
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if(x <= center_x && y > center_y) //! second quadrant
{
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if(x > center_x && y <= center_y) //! third quadrant
{
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else //! fourth quadrant
{
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; //! sets branch pointer to nothing...
n->is_leaf = 0;
}
//}
x = (double)p->x_pos;
y = (double)p->y_pos;
if(x <= center_x && y <= center_y) //! first quadrant
{
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if(x <= center_x && y > center_y) //! second quadrant
{
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if(x > center_x && y <= center_y) //! third quadrant
{
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else //! fourth quadrant
{
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
return n;
}
}
这一切都是由：
node *root = NULL;
root = quadtree_insert(root, &particles[0],0.500,0.5,1);
for(i = 1; i < nParticles; i++)
{
quadtree_insert(root, &particles[i],0.5000,0.5,1);
}
其中 "particles" 以 struct Particle* 粒子的形式传递。一个粒子被定义为：
struct Particle {
double mass;
double x_pos;
double y_pos;
double x_vel;
double y_vel;
};
typedef struct Particle * Particle_structure;
每次迭代后，在 for 循环后，根都被释放，代码适用于小于 ~200 的样本，valgrind 不会给出这些错误。添加到这里的是一个中间函数，它对粒子执行一些算术运算：
double quadtree_calculate_forcey(struct Particle *p, node *n, double theta_max, double delta_t, int numP, double epsilon, double G)
{
if(n != NULL)
{
double d_x = (n->center_x - p->x_pos);
double d_y = (n->center_y - p->y_pos);
double r_2 = d_x * d_x + d_y * d_y;
r_2 = sqrt(r_2) + epsilon;
if(theta_max <= (n->width / r_2) && !n->is_leaf){
double a = 0;
if(n->sw != NULL)
a += quadtree_calculate_forcey(p, n->sw, theta_max, delta_t, numP, epsilon, G);
if(n->nw != NULL)
a += quadtree_calculate_forcey(p, n->nw, theta_max, delta_t, numP, epsilon, G);
if(n->se != NULL)
a += quadtree_calculate_forcey(p, n->se, theta_max, delta_t, numP, epsilon, G);
if(n->ne != NULL)
a += quadtree_calculate_forcey(p, n->ne, theta_max, delta_t, numP, epsilon, G);
return a;
}
else
{
double fy;
double mass;
if(d_x == 0 && d_y == 0){ // could be comparing the same star
//printf("RÖÖÖVHATT\n");
return 0;
}
else{
//printf("MASS : %f\n", n->m);
mass = n->m;
//printf("MASS : %f\n", mass);
fy = G * (mass * p->mass/ pow(r_2,3)) * d_y;
//printf("DY:%f DX:%f R_2:%f MASSA:%f\n",d_y, d_x, r_2 - epsilon, mass);
// printf("HIT SKA JAG: %f\n",d_y);
return fy;
}
}
}
return 0.0;
}
铃声是它滚动了一点（清除根，并为新位置重做），所以递归中的变量数量肯定不是问题（？）。当涉及到一些指针声明的分配时，我很确定我会和鱼一起游泳。任何 thoughts/help 都会让你站在 Odin 好的一边！
编辑：我们如何找到质量中心等的示例。节点质量以相同的方式完成。
double set_centerx(node *n)
{
if(n != NULL)
{
if(!(n->is_leaf))
{
double a = set_centerx(n->ne);
double b = set_centerx(n->nw);
double c = set_centerx(n->se);
double d = set_centerx(n->sw);
double m1 = get_mass2(n->ne);
double m2 = get_mass2(n->nw);
double m3 = get_mass2(n->se);
double m4 = get_mass2(n->sw);
n->center_x = (double)(m1*a + m2*b + m3*c + m4*d)/(m1+m2+m3+m4);
return n->center_x;
}
n->center_x =n->p->x_pos;
return n->p->x_pos;
}
return 0;
}
部分分析
我可以容忍地确定问题主要是由于 quadtree_insert() 的主要 else 子句中的大量重复代码（接近重复代码）造成的。我用 Fragment 1A 和 Fragment 1B 的注释标记了开头 — 片段 1B 现在也用 #ifdef DO_REPEAT 和 #endif.
加框
else
{
/* Fragment 1A */
// n = (node*)malloc(sizeof(node));
double x;
double y;
if (width < 1e-300)
{
width = 1e-300;
}
if (n->is_leaf == 1) // ! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; // ! sets branch pointer to nothing...
n->is_leaf = 0;
}
#ifdef DO_REPEAT
/* Fragment 1B */
x = (double)p->x_pos;
y = (double)p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
#endif /* DO_REPEAT */
return n;
}
我拿了你的代码并重新排序了一些片段，并使用了这个 main()。请注意，我不需要这样做；你应该已经制作了一个 MCVE (How to create a Minimal, Complete, and Verifiable Example?) or SSCCE (Short, Self-Contained, Correct Example) — 相同基本思想的两个名称和链接。
static struct Particle particles[] =
{
{ 19.99, 96.07, 62.79, -99.46, 19.70 },
{ 12.94, 1.43, -33.45, 31.80, -66.08 },
{ 6.49, 16.99, -20.83, 92.51, 35.98 },
{ 17.01, -28.85, -94.10, 42.82, -1.30 },
{ 14.27, 85.07, 88.21, 11.22, 16.85 },
{ 15.73, -56.37, 46.85, 27.40, -15.15 },
{ 1.53, -49.44, -64.27, -29.45, -38.25 },
{ 8.03, 92.11, -47.50, 63.77, -29.99 },
{ 8.67, -99.81, 73.19, 18.75, 88.66 },
{ 16.36, 66.33, 14.23, 87.65, 40.01 },
};
enum { nParticles = sizeof(particles) / sizeof(particles[0]) };
int main(void)
{
node *root = NULL;
printf("Particle 0:\n");
root = quadtree_insert(root, &particles[0], 0.500, 0.5, 1);
for (int i = 1; i < nParticles; i++)
{
printf("Particle %d:\n", i);
quadtree_insert(root, &particles[i], 0.5000, 0.5, 1);
}
return 0;
}
我使用 运行dom 数字生成器生成 table 个值：
random -n 10 -T ' { %6:2[1:20]f, %6:2[-100:100]f, %6:2[-100:100]f, %6:2[-100:100]f, %6:2[-100:100]f },'
它在第 6 条上让我崩溃了。概括地说，valgrind 说：
…startup blurb omitted…
Particle 0:
Particle 1:
Particle 2:
Particle 3:
Particle 4:
Particle 5:
Particle 6:
==79528==
==79528== Process terminating with default action of signal 11 (SIGSEGV)
==79528== Access not within mapped region at address 0x104003FD0
==79528== at 0x100008E39: malloc (vg_replace_malloc.c:302)
==79528== If you believe this happened as a result of a stack
==79528== overflow in your program's main thread (unlikely but
==79528== possible), you can try to increase the size of the
==79528== main thread stack using the --main-stacksize= flag.
==79528== The main thread stack size used in this run was 8388608.
==79528==
==79528== HEAP SUMMARY:
==79528== in use at exit: 6,017,700 bytes in 75,080 blocks
==79528== total heap usage: 75,162 allocs, 82 frees, 6,023,876 bytes allocated
==79528==
==79528== LEAK SUMMARY:
==79528== definitely lost: 4,120 bytes in 2 blocks
==79528== indirectly lost: 2,288 bytes in 6 blocks
==79528== possibly lost: 4,880 bytes in 45 blocks
==79528== still reachable: 5,997,720 bytes in 74,904 blocks
==79528== suppressed: 8,692 bytes in 123 blocks
==79528== Rerun with --leak-check=full to see details of leaked memory
即使在我 运行 的 Mac 上，也表明存在问题；压抑的东西还行，其他大部分不行
当没有-DDO_REPEAT编译时（正常编译），示例程序运行s完成。当然，它会泄漏，因为没有释放内存的代码。
Particle 0:
Particle 1:
Particle 2:
Particle 3:
Particle 4:
Particle 5:
Particle 6:
Particle 7:
Particle 8:
Particle 9:
==79683==
==79683== HEAP SUMMARY:
==79683== in use at exit: 26,580 bytes in 191 blocks
==79683== total heap usage: 273 allocs, 82 frees, 32,756 bytes allocated
==79683==
==79683== LEAK SUMMARY:
==79683== definitely lost: 4,200 bytes in 3 blocks
==79683== indirectly lost: 2,368 bytes in 7 blocks
==79683== possibly lost: 4,880 bytes in 45 blocks
==79683== still reachable: 6,440 bytes in 13 blocks
==79683== suppressed: 8,692 bytes in 123 blocks
==79683== Rerun with --leak-check=full to see details of leaked memory
请注意，使用的内存比以前少得多。
如果由于缺少 DO_REPEAT 而删除的代码实际上是关键的，那么错误要么在于该代码，要么在于片段 1A 中完成的设置工作。但是，在我看来，通过递归调用两次插入同一个粒子可能是问题的根源。
我还注意到代码中没有使用 quadtree_calculate_forcey() 函数；它绝对不是 MCVE 的一部分。
需要片段 1B
Do note that the "repeated" code is the same in form, but not in detail. That is, fragment 1B uses n->p where fragment 1A uses p. I think this is purposeful: the idea seems to be to force all data (the Particles) out to leaf nodes.
和radnaskela确认：
It's exactly like John says, every particle should occupy one end-node. My suspicions say its to do with the movement, and planets aligning so that a ridiculous number of iterations have to be done to open two free squares. I'm unable to see a solution to my problem though, any thoughts?
逻辑有问题，虽然我不完全确定是什么问题。我要做的第一步是向 quadtree_insert() 中的代码添加一些检测，如下所示：
static node *quadtree_insert(node *n, struct Particle *p, double center_x, double center_y, double width)
{
printf("Centre (X,Y) = (%6.2f,%6.2f)\n", center_x, center_y);
if (n == NULL)
{
n = (node *)malloc(sizeof(node));
n->is_leaf = 1;
n->p = p;
// n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
if (width < 1e-300)
{
n->width = 1e-300;
}
else
n->width = width;
return n;
}
else
{
// n = (node*)malloc(sizeof(node));
double x;
double y;
if (width < 1e-300)
{
width = 1e-300;
}
if (n->is_leaf == 1) // ! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
printf("Recurse SW 1: ");
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
printf("Recurse NW 1: ");
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
printf("Recurse SE 1: ");
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
printf("Recurse NE 1: ");
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; // ! sets branch pointer to nothing...
n->is_leaf = 0;
}
x = (double)p->x_pos;
y = (double)p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
printf("Recurse SW 2: ");
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
printf("Recurse NW 2: ");
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
printf("Recurse SE 2: ");
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
printf("Recurse NE 2: ");
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
return n;
}
}
当运行10点数据集时，到第5个粒子的输出为：
Particle 0: Centre (X,Y) = ( 0.50, 0.50)
Particle 1: Centre (X,Y) = ( 0.50, 0.50)
Recurse NE 1: Centre (X,Y) = ( 0.75, 0.75)
Recurse SE 2: Centre (X,Y) = ( 0.75, 0.25)
Particle 2: Centre (X,Y) = ( 0.50, 0.50)
Recurse SE 2: Centre (X,Y) = ( 0.75, 0.25)
Recurse SE 1: Centre (X,Y) = ( 1.12, 0.12)
Recurse SE 2: Centre (X,Y) = ( 1.12, 0.12)
Recurse SE 1: Centre (X,Y) = ( 1.69, 0.06)
Recurse SE 2: Centre (X,Y) = ( 1.69, 0.06)
Recurse SW 1: Centre (X,Y) = ( 0.84, 0.03)
Recurse SE 2: Centre (X,Y) = ( 2.53, 0.03)
Particle 3: Centre (X,Y) = ( 0.50, 0.50)
Recurse SW 2: Centre (X,Y) = ( 0.25, 0.25)
Particle 4: Centre (X,Y) = ( 0.50, 0.50)
Recurse NE 2: Centre (X,Y) = ( 0.75, 0.75)
Recurse NE 1: Centre (X,Y) = ( 1.12, 1.12)
Recurse NE 2: Centre (X,Y) = ( 1.12, 1.12)
Recurse NE 1: Centre (X,Y) = ( 1.69, 1.69)
Recurse NE 2: Centre (X,Y) = ( 1.69, 1.69)
Recurse NE 1: Centre (X,Y) = ( 2.53, 2.53)
Recurse NE 2: Centre (X,Y) = ( 2.53, 2.53)
Recurse NE 1: Centre (X,Y) = ( 3.80, 3.80)
Recurse NE 2: Centre (X,Y) = ( 3.80, 3.80)
Recurse NE 1: Centre (X,Y) = ( 5.70, 5.70)
Recurse NE 2: Centre (X,Y) = ( 5.70, 5.70)
Recurse NE 1: Centre (X,Y) = ( 8.54, 8.54)
Recurse NE 2: Centre (X,Y) = ( 8.54, 8.54)
Recurse NE 1: Centre (X,Y) = ( 12.81, 12.81)
Recurse NE 2: Centre (X,Y) = ( 12.81, 12.81)
Recurse NE 1: Centre (X,Y) = ( 19.22, 19.22)
Recurse NE 2: Centre (X,Y) = ( 19.22, 19.22)
Recurse NE 1: Centre (X,Y) = ( 28.83, 28.83)
Recurse NE 2: Centre (X,Y) = ( 28.83, 28.83)
Recurse NE 1: Centre (X,Y) = ( 43.25, 43.25)
Recurse NE 2: Centre (X,Y) = ( 43.25, 43.25)
Recurse NE 1: Centre (X,Y) = ( 64.87, 64.87)
Recurse NE 2: Centre (X,Y) = ( 64.87, 64.87)
Recurse SE 1: Centre (X,Y) = ( 97.31, 32.44)
Recurse NE 2: Centre (X,Y) = ( 97.31, 97.31)
Particle 5: Centre (X,Y) = ( 0.50, 0.50)
Recurse NW 2: Centre (X,Y) = ( 0.25, 0.75)
我对那里的某些处理以及第 4 条的条目数感到有点惊讶。这可能表明存在问题。
然后它处理第 6 个粒子：
Particle 6: Centre (X,Y) = ( 0.50, 0.50)
Recurse SW 2: Centre (X,Y) = ( 0.25, 0.25)
Recurse SW 1: Centre (X,Y) = ( 0.12, 0.12)
Recurse SW 2: Centre (X,Y) = ( 0.12, 0.12)
Recurse SW 1: Centre (X,Y) = ( 0.06, 0.06)
Recurse SW 2: Centre (X,Y) = ( 0.06, 0.06)
Recurse SW 1: Centre (X,Y) = ( 0.03, 0.03)
Recurse SW 2: Centre (X,Y) = ( 0.03, 0.03)
Recurse SW 1: Centre (X,Y) = ( 0.02, 0.02)
Recurse SW 2: Centre (X,Y) = ( 0.02, 0.02)
Recurse SW 1: Centre (X,Y) = ( 0.01, 0.01)
Recurse SW 2: Centre (X,Y) = ( 0.01, 0.01)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
从那里开始变得非常乏味。你要问的问题是"what should be happening with that point"？也有一个转储四叉树结构的函数可能是个好主意。
这基本上表明你还没有充分分析加点的条件。我不清楚为什么当插入位于 SE 或 NE quad运行ts 时将中心坐标乘以 1.5，而当插入位于 SW 或 NW quad运行ts 时乘以 0.5（也不知道为什么要使用加加乘法而不是简单地乘法）。
width 小于 1e-300 的测试有点令人担忧。鉴于您使用值 1.0 调用函数，每次递归时将宽度减半需要一段时间才能变小。
更好的跟踪可能会报告 x 和 y 的值以及中心坐标。
我发现你的代码有几个问题。
首先，您有一个相当大的非平凡代码块的近乎重复。这两个块似乎仅在它们引用 Particle 时有所不同；这需要被分解成一个辅助函数。
其次，考虑这个片段：
else if (x <= center_x && y > center_y) // ! second quadrant
{
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
我将输入 center_x 和 center_y 作为整个区域正方形补丁的中心，width 作为该补丁的边长。在这种情况下，补丁的西北象限的中心是 (center_x - width * 0.25, center_y + width * 0.25)。只有在某些特殊情况下，这才会与递归调用中的计算一致。您实际上可以根据其中心坐标确定补丁的宽度，但不像您尝试的那样直接。
同样适用于所有其他七个递归调用。
第三，考虑一下如果你最终得到两个坐标非常相似——或者更糟，完全相同——的粒子会发生什么。他们可以是模拟中仅有的两个。
第一个粒子最初位于代表整个模拟区域的节点中。当插入第二个时，四叉树遍历到同一个节点。因此，原始粒子被移动到原始节点区域的新创建象限，并且第二个节点的插入过程恢复。但是插入 then again 到达原始粒子所在的相同节点，并重复该过程。也许它会再次重复。再一次。
因此，您的递归没有明确的界限。如果你碰巧得到两个坐标相同的粒子——例如，它们可能被固定在模拟区域的角落——那么你肯定会无休止地递归，最终产生堆栈溢出。然而，即使没有两个粒子具有相同的坐标，也可以想象两个粒子将足够靠近以使递归深入到堆栈溢出。
我怀疑这个问题被上一个问题加重了，但它不依赖于那个问题。
解决这个问题的一种方法是合并足够接近的粒子。这将允许您根据补丁宽度对递归深度进行严格限制。或者，如果递归太深，您可以中止。或者，在这种情况下，您可以让它崩溃，因为它已经崩溃了。
不管怎样，我认为您的 quadtree_insert() 应该是这样的：
// adjust as needed:
#define MIN_PARTICLE_SEPARATION 1e-300
void quadtree_insert_helper(node *n, Particle *p, double center_x,
double center_y, double width) {
width *= 0.5;
double new_center_x = center_x
+ width * ((p->x_pos <= center_x) ? -0.5 : 0.5);
double new_center_y = center_y
+ width * ((p->y_pos <= center_y) ? -0.5 : 0.5);
node **quad;
if (new_center_x <= center_x && new_center_y <= center_y) {
//! first quadrant
quad = &n->sw;
} else if (new_center_x <= center_x && new_center_y > center_y) {
//! second quadrant
quad = &n->nw;
} else if (new_center_x > center_x && new_center_y <= center_y) {
//! third quadrant
quad = &n->se;
} else {
//! fourth quadrant
quad = &n->ne;
}
*quad = quadtree_insert(*quad, p, new_center_x, new_center_y, width);
}
node* quadtree_insert(node *n, struct Particle *p, double center_x,
double center_y, double width) {
if (n == NULL) {
n = malloc(sizeof(*n));
n->is_leaf = 1;
n->p = p;
//n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
n->center_x = center_x;
n->center_y = center_y;
n->width = width;
} else if (n->is_leaf && (with <= MIN_PARTICLE_SEPARATION)) {
// ... merge particles ...
} else {
if (n->is_leaf) { //! that is, if the node is not a branch
quadtree_insert_helper(n, n->p, center_x, center_y, width);
//! the node is now a branch
n->p = NULL;
n->is_leaf = 0;
}
quadtree_insert_helper(n, p, center_x, center_y, width);
}
return n;
}
我一直在疯狂地寻找答案无济于事；我创建了一个四叉树，它应该对结构声明的超过 1000 个对象的数组进行排序：
typedef struct node
{
char is_leaf;
struct Particle *p;
double m;
double center_x;
double center_y;
double width;
struct node *sw;
struct node *nw;
struct node *se;
struct node *ne;
} node;
使用以下函数生成四叉树：
node* quadtree_insert(node *n, struct Particle *p, double center_x, double center_y, double width)
{
if(n == NULL)
{
n = (node*)malloc(sizeof(node));
n->is_leaf = 1;
n->p = p;
//n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
if(width < 1e-300){
n->width = 1e-300;
}
else
n->width = width;
return n;
}
else{
//n = (node*)malloc(sizeof(node));
double x;
double y;
if(width < 1e-300){
width = 1e-300;
}
if(n->is_leaf == 1) //! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if(x <= center_x && y <= center_y) //! first quadrant
{
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if(x <= center_x && y > center_y) //! second quadrant
{
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if(x > center_x && y <= center_y) //! third quadrant
{
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else //! fourth quadrant
{
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; //! sets branch pointer to nothing...
n->is_leaf = 0;
}
//}
x = (double)p->x_pos;
y = (double)p->y_pos;
if(x <= center_x && y <= center_y) //! first quadrant
{
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if(x <= center_x && y > center_y) //! second quadrant
{
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if(x > center_x && y <= center_y) //! third quadrant
{
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else //! fourth quadrant
{
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
return n;
}
}
这一切都是由：
node *root = NULL;
root = quadtree_insert(root, &particles[0],0.500,0.5,1);
for(i = 1; i < nParticles; i++)
{
quadtree_insert(root, &particles[i],0.5000,0.5,1);
}
其中 "particles" 以 struct Particle* 粒子的形式传递。一个粒子被定义为：
struct Particle {
double mass;
double x_pos;
double y_pos;
double x_vel;
double y_vel;
};
typedef struct Particle * Particle_structure;
每次迭代后，在 for 循环后，根都被释放，代码适用于小于 ~200 的样本，valgrind 不会给出这些错误。添加到这里的是一个中间函数，它对粒子执行一些算术运算：
double quadtree_calculate_forcey(struct Particle *p, node *n, double theta_max, double delta_t, int numP, double epsilon, double G)
{
if(n != NULL)
{
double d_x = (n->center_x - p->x_pos);
double d_y = (n->center_y - p->y_pos);
double r_2 = d_x * d_x + d_y * d_y;
r_2 = sqrt(r_2) + epsilon;
if(theta_max <= (n->width / r_2) && !n->is_leaf){
double a = 0;
if(n->sw != NULL)
a += quadtree_calculate_forcey(p, n->sw, theta_max, delta_t, numP, epsilon, G);
if(n->nw != NULL)
a += quadtree_calculate_forcey(p, n->nw, theta_max, delta_t, numP, epsilon, G);
if(n->se != NULL)
a += quadtree_calculate_forcey(p, n->se, theta_max, delta_t, numP, epsilon, G);
if(n->ne != NULL)
a += quadtree_calculate_forcey(p, n->ne, theta_max, delta_t, numP, epsilon, G);
return a;
}
else
{
double fy;
double mass;
if(d_x == 0 && d_y == 0){ // could be comparing the same star
//printf("RÖÖÖVHATT\n");
return 0;
}
else{
//printf("MASS : %f\n", n->m);
mass = n->m;
//printf("MASS : %f\n", mass);
fy = G * (mass * p->mass/ pow(r_2,3)) * d_y;
//printf("DY:%f DX:%f R_2:%f MASSA:%f\n",d_y, d_x, r_2 - epsilon, mass);
// printf("HIT SKA JAG: %f\n",d_y);
return fy;
}
}
}
return 0.0;
}
铃声是它滚动了一点（清除根，并为新位置重做），所以递归中的变量数量肯定不是问题（？）。当涉及到一些指针声明的分配时，我很确定我会和鱼一起游泳。任何 thoughts/help 都会让你站在 Odin 好的一边！
编辑：我们如何找到质量中心等的示例。节点质量以相同的方式完成。
double set_centerx(node *n)
{
if(n != NULL)
{
if(!(n->is_leaf))
{
double a = set_centerx(n->ne);
double b = set_centerx(n->nw);
double c = set_centerx(n->se);
double d = set_centerx(n->sw);
double m1 = get_mass2(n->ne);
double m2 = get_mass2(n->nw);
double m3 = get_mass2(n->se);
double m4 = get_mass2(n->sw);
n->center_x = (double)(m1*a + m2*b + m3*c + m4*d)/(m1+m2+m3+m4);
return n->center_x;
}
n->center_x =n->p->x_pos;
return n->p->x_pos;
}
return 0;
}
部分分析
我可以容忍地确定问题主要是由于 quadtree_insert() 的主要 else 子句中的大量重复代码（接近重复代码）造成的。我用 Fragment 1A 和 Fragment 1B 的注释标记了开头 — 片段 1B 现在也用 #ifdef DO_REPEAT 和 #endif.
else
{
/* Fragment 1A */
// n = (node*)malloc(sizeof(node));
double x;
double y;
if (width < 1e-300)
{
width = 1e-300;
}
if (n->is_leaf == 1) // ! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; // ! sets branch pointer to nothing...
n->is_leaf = 0;
}
#ifdef DO_REPEAT
/* Fragment 1B */
x = (double)p->x_pos;
y = (double)p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
#endif /* DO_REPEAT */
return n;
}
我拿了你的代码并重新排序了一些片段，并使用了这个 main()。请注意，我不需要这样做；你应该已经制作了一个 MCVE (How to create a Minimal, Complete, and Verifiable Example?) or SSCCE (Short, Self-Contained, Correct Example) — 相同基本思想的两个名称和链接。
static struct Particle particles[] =
{
{ 19.99, 96.07, 62.79, -99.46, 19.70 },
{ 12.94, 1.43, -33.45, 31.80, -66.08 },
{ 6.49, 16.99, -20.83, 92.51, 35.98 },
{ 17.01, -28.85, -94.10, 42.82, -1.30 },
{ 14.27, 85.07, 88.21, 11.22, 16.85 },
{ 15.73, -56.37, 46.85, 27.40, -15.15 },
{ 1.53, -49.44, -64.27, -29.45, -38.25 },
{ 8.03, 92.11, -47.50, 63.77, -29.99 },
{ 8.67, -99.81, 73.19, 18.75, 88.66 },
{ 16.36, 66.33, 14.23, 87.65, 40.01 },
};
enum { nParticles = sizeof(particles) / sizeof(particles[0]) };
int main(void)
{
node *root = NULL;
printf("Particle 0:\n");
root = quadtree_insert(root, &particles[0], 0.500, 0.5, 1);
for (int i = 1; i < nParticles; i++)
{
printf("Particle %d:\n", i);
quadtree_insert(root, &particles[i], 0.5000, 0.5, 1);
}
return 0;
}
我使用 运行dom 数字生成器生成 table 个值：
random -n 10 -T ' { %6:2[1:20]f, %6:2[-100:100]f, %6:2[-100:100]f, %6:2[-100:100]f, %6:2[-100:100]f },'
它在第 6 条上让我崩溃了。概括地说，valgrind 说：
…startup blurb omitted…
Particle 0:
Particle 1:
Particle 2:
Particle 3:
Particle 4:
Particle 5:
Particle 6:
==79528==
==79528== Process terminating with default action of signal 11 (SIGSEGV)
==79528== Access not within mapped region at address 0x104003FD0
==79528== at 0x100008E39: malloc (vg_replace_malloc.c:302)
==79528== If you believe this happened as a result of a stack
==79528== overflow in your program's main thread (unlikely but
==79528== possible), you can try to increase the size of the
==79528== main thread stack using the --main-stacksize= flag.
==79528== The main thread stack size used in this run was 8388608.
==79528==
==79528== HEAP SUMMARY:
==79528== in use at exit: 6,017,700 bytes in 75,080 blocks
==79528== total heap usage: 75,162 allocs, 82 frees, 6,023,876 bytes allocated
==79528==
==79528== LEAK SUMMARY:
==79528== definitely lost: 4,120 bytes in 2 blocks
==79528== indirectly lost: 2,288 bytes in 6 blocks
==79528== possibly lost: 4,880 bytes in 45 blocks
==79528== still reachable: 5,997,720 bytes in 74,904 blocks
==79528== suppressed: 8,692 bytes in 123 blocks
==79528== Rerun with --leak-check=full to see details of leaked memory
即使在我 运行 的 Mac 上，也表明存在问题；压抑的东西还行，其他大部分不行
当没有-DDO_REPEAT编译时（正常编译），示例程序运行s完成。当然，它会泄漏，因为没有释放内存的代码。
Particle 0:
Particle 1:
Particle 2:
Particle 3:
Particle 4:
Particle 5:
Particle 6:
Particle 7:
Particle 8:
Particle 9:
==79683==
==79683== HEAP SUMMARY:
==79683== in use at exit: 26,580 bytes in 191 blocks
==79683== total heap usage: 273 allocs, 82 frees, 32,756 bytes allocated
==79683==
==79683== LEAK SUMMARY:
==79683== definitely lost: 4,200 bytes in 3 blocks
==79683== indirectly lost: 2,368 bytes in 7 blocks
==79683== possibly lost: 4,880 bytes in 45 blocks
==79683== still reachable: 6,440 bytes in 13 blocks
==79683== suppressed: 8,692 bytes in 123 blocks
==79683== Rerun with --leak-check=full to see details of leaked memory
请注意，使用的内存比以前少得多。
如果由于缺少 DO_REPEAT 而删除的代码实际上是关键的，那么错误要么在于该代码，要么在于片段 1A 中完成的设置工作。但是，在我看来，通过递归调用两次插入同一个粒子可能是问题的根源。
我还注意到代码中没有使用 quadtree_calculate_forcey() 函数；它绝对不是 MCVE 的一部分。
需要片段 1B
Do note that the "repeated" code is the same in form, but not in detail. That is, fragment 1B uses
n->pwhere fragment 1A usesp. I think this is purposeful: the idea seems to be to force all data (the Particles) out to leaf nodes.
和radnaskela确认：
It's exactly like John says, every particle should occupy one end-node. My suspicions say its to do with the movement, and planets aligning so that a ridiculous number of iterations have to be done to open two free squares. I'm unable to see a solution to my problem though, any thoughts?
逻辑有问题，虽然我不完全确定是什么问题。我要做的第一步是向 quadtree_insert() 中的代码添加一些检测，如下所示：
static node *quadtree_insert(node *n, struct Particle *p, double center_x, double center_y, double width)
{
printf("Centre (X,Y) = (%6.2f,%6.2f)\n", center_x, center_y);
if (n == NULL)
{
n = (node *)malloc(sizeof(node));
n->is_leaf = 1;
n->p = p;
// n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
if (width < 1e-300)
{
n->width = 1e-300;
}
else
n->width = width;
return n;
}
else
{
// n = (node*)malloc(sizeof(node));
double x;
double y;
if (width < 1e-300)
{
width = 1e-300;
}
if (n->is_leaf == 1) // ! that is, if the node is not a branch
{
x = (double)n->p->x_pos;
y = (double)n->p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
printf("Recurse SW 1: ");
n->sw = quadtree_insert(n->sw, n->p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
printf("Recurse NW 1: ");
n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
printf("Recurse SE 1: ");
n->se = quadtree_insert(n->se, n->p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
printf("Recurse NE 1: ");
n->ne = quadtree_insert(n->ne, n->p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
n->p = NULL; // ! sets branch pointer to nothing...
n->is_leaf = 0;
}
x = (double)p->x_pos;
y = (double)p->y_pos;
if (x <= center_x && y <= center_y) // ! first quadrant
{
printf("Recurse SW 2: ");
n->sw = quadtree_insert(n->sw, p, center_x * 0.5, center_y * 0.5, width * 0.5);
}
else if (x <= center_x && y > center_y) // ! second quadrant
{
printf("Recurse NW 2: ");
n->nw = quadtree_insert(n->nw, p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
else if (x > center_x && y <= center_y) // ! third quadrant
{
printf("Recurse SE 2: ");
n->se = quadtree_insert(n->se, p, center_x + center_x * 0.5, center_y * 0.5, width * 0.5);
}
else // ! fourth quadrant
{
printf("Recurse NE 2: ");
n->ne = quadtree_insert(n->ne, p, center_x + center_x * 0.5, center_y + center_y * 0.5, width * 0.5);
}
return n;
}
}
当运行10点数据集时，到第5个粒子的输出为：
Particle 0: Centre (X,Y) = ( 0.50, 0.50)
Particle 1: Centre (X,Y) = ( 0.50, 0.50)
Recurse NE 1: Centre (X,Y) = ( 0.75, 0.75)
Recurse SE 2: Centre (X,Y) = ( 0.75, 0.25)
Particle 2: Centre (X,Y) = ( 0.50, 0.50)
Recurse SE 2: Centre (X,Y) = ( 0.75, 0.25)
Recurse SE 1: Centre (X,Y) = ( 1.12, 0.12)
Recurse SE 2: Centre (X,Y) = ( 1.12, 0.12)
Recurse SE 1: Centre (X,Y) = ( 1.69, 0.06)
Recurse SE 2: Centre (X,Y) = ( 1.69, 0.06)
Recurse SW 1: Centre (X,Y) = ( 0.84, 0.03)
Recurse SE 2: Centre (X,Y) = ( 2.53, 0.03)
Particle 3: Centre (X,Y) = ( 0.50, 0.50)
Recurse SW 2: Centre (X,Y) = ( 0.25, 0.25)
Particle 4: Centre (X,Y) = ( 0.50, 0.50)
Recurse NE 2: Centre (X,Y) = ( 0.75, 0.75)
Recurse NE 1: Centre (X,Y) = ( 1.12, 1.12)
Recurse NE 2: Centre (X,Y) = ( 1.12, 1.12)
Recurse NE 1: Centre (X,Y) = ( 1.69, 1.69)
Recurse NE 2: Centre (X,Y) = ( 1.69, 1.69)
Recurse NE 1: Centre (X,Y) = ( 2.53, 2.53)
Recurse NE 2: Centre (X,Y) = ( 2.53, 2.53)
Recurse NE 1: Centre (X,Y) = ( 3.80, 3.80)
Recurse NE 2: Centre (X,Y) = ( 3.80, 3.80)
Recurse NE 1: Centre (X,Y) = ( 5.70, 5.70)
Recurse NE 2: Centre (X,Y) = ( 5.70, 5.70)
Recurse NE 1: Centre (X,Y) = ( 8.54, 8.54)
Recurse NE 2: Centre (X,Y) = ( 8.54, 8.54)
Recurse NE 1: Centre (X,Y) = ( 12.81, 12.81)
Recurse NE 2: Centre (X,Y) = ( 12.81, 12.81)
Recurse NE 1: Centre (X,Y) = ( 19.22, 19.22)
Recurse NE 2: Centre (X,Y) = ( 19.22, 19.22)
Recurse NE 1: Centre (X,Y) = ( 28.83, 28.83)
Recurse NE 2: Centre (X,Y) = ( 28.83, 28.83)
Recurse NE 1: Centre (X,Y) = ( 43.25, 43.25)
Recurse NE 2: Centre (X,Y) = ( 43.25, 43.25)
Recurse NE 1: Centre (X,Y) = ( 64.87, 64.87)
Recurse NE 2: Centre (X,Y) = ( 64.87, 64.87)
Recurse SE 1: Centre (X,Y) = ( 97.31, 32.44)
Recurse NE 2: Centre (X,Y) = ( 97.31, 97.31)
Particle 5: Centre (X,Y) = ( 0.50, 0.50)
Recurse NW 2: Centre (X,Y) = ( 0.25, 0.75)
我对那里的某些处理以及第 4 条的条目数感到有点惊讶。这可能表明存在问题。
然后它处理第 6 个粒子：
Particle 6: Centre (X,Y) = ( 0.50, 0.50)
Recurse SW 2: Centre (X,Y) = ( 0.25, 0.25)
Recurse SW 1: Centre (X,Y) = ( 0.12, 0.12)
Recurse SW 2: Centre (X,Y) = ( 0.12, 0.12)
Recurse SW 1: Centre (X,Y) = ( 0.06, 0.06)
Recurse SW 2: Centre (X,Y) = ( 0.06, 0.06)
Recurse SW 1: Centre (X,Y) = ( 0.03, 0.03)
Recurse SW 2: Centre (X,Y) = ( 0.03, 0.03)
Recurse SW 1: Centre (X,Y) = ( 0.02, 0.02)
Recurse SW 2: Centre (X,Y) = ( 0.02, 0.02)
Recurse SW 1: Centre (X,Y) = ( 0.01, 0.01)
Recurse SW 2: Centre (X,Y) = ( 0.01, 0.01)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 1: Centre (X,Y) = ( 0.00, 0.00)
Recurse SW 2: Centre (X,Y) = ( 0.00, 0.00)
从那里开始变得非常乏味。你要问的问题是"what should be happening with that point"？也有一个转储四叉树结构的函数可能是个好主意。
这基本上表明你还没有充分分析加点的条件。我不清楚为什么当插入位于 SE 或 NE quad运行ts 时将中心坐标乘以 1.5，而当插入位于 SW 或 NW quad运行ts 时乘以 0.5（也不知道为什么要使用加加乘法而不是简单地乘法）。
width 小于 1e-300 的测试有点令人担忧。鉴于您使用值 1.0 调用函数，每次递归时将宽度减半需要一段时间才能变小。
更好的跟踪可能会报告 x 和 y 的值以及中心坐标。
我发现你的代码有几个问题。
首先，您有一个相当大的非平凡代码块的近乎重复。这两个块似乎仅在它们引用 Particle 时有所不同；这需要被分解成一个辅助函数。
其次，考虑这个片段：
else if (x <= center_x && y > center_y) // ! second quadrant { n->nw = quadtree_insert(n->nw, n->p, center_x * 0.5, center_y + center_y * 0.5, width * 0.5); }
我将输入 center_x 和 center_y 作为整个区域正方形补丁的中心，width 作为该补丁的边长。在这种情况下，补丁的西北象限的中心是 (center_x - width * 0.25, center_y + width * 0.25)。只有在某些特殊情况下，这才会与递归调用中的计算一致。您实际上可以根据其中心坐标确定补丁的宽度，但不像您尝试的那样直接。
同样适用于所有其他七个递归调用。
第三，考虑一下如果你最终得到两个坐标非常相似——或者更糟，完全相同——的粒子会发生什么。他们可以是模拟中仅有的两个。
第一个粒子最初位于代表整个模拟区域的节点中。当插入第二个时，四叉树遍历到同一个节点。因此，原始粒子被移动到原始节点区域的新创建象限，并且第二个节点的插入过程恢复。但是插入 then again 到达原始粒子所在的相同节点，并重复该过程。也许它会再次重复。再一次。
因此，您的递归没有明确的界限。如果你碰巧得到两个坐标相同的粒子——例如，它们可能被固定在模拟区域的角落——那么你肯定会无休止地递归，最终产生堆栈溢出。然而，即使没有两个粒子具有相同的坐标，也可以想象两个粒子将足够靠近以使递归深入到堆栈溢出。
我怀疑这个问题被上一个问题加重了，但它不依赖于那个问题。
解决这个问题的一种方法是合并足够接近的粒子。这将允许您根据补丁宽度对递归深度进行严格限制。或者，如果递归太深，您可以中止。或者，在这种情况下，您可以让它崩溃，因为它已经崩溃了。
不管怎样，我认为您的 quadtree_insert() 应该是这样的：
// adjust as needed:
#define MIN_PARTICLE_SEPARATION 1e-300
void quadtree_insert_helper(node *n, Particle *p, double center_x,
double center_y, double width) {
width *= 0.5;
double new_center_x = center_x
+ width * ((p->x_pos <= center_x) ? -0.5 : 0.5);
double new_center_y = center_y
+ width * ((p->y_pos <= center_y) ? -0.5 : 0.5);
node **quad;
if (new_center_x <= center_x && new_center_y <= center_y) {
//! first quadrant
quad = &n->sw;
} else if (new_center_x <= center_x && new_center_y > center_y) {
//! second quadrant
quad = &n->nw;
} else if (new_center_x > center_x && new_center_y <= center_y) {
//! third quadrant
quad = &n->se;
} else {
//! fourth quadrant
quad = &n->ne;
}
*quad = quadtree_insert(*quad, p, new_center_x, new_center_y, width);
}
node* quadtree_insert(node *n, struct Particle *p, double center_x,
double center_y, double width) {
if (n == NULL) {
n = malloc(sizeof(*n));
n->is_leaf = 1;
n->p = p;
//n->m = 0.1;
n->sw = NULL;
n->se = NULL;
n->nw = NULL;
n->ne = NULL;
n->center_x = center_x;
n->center_y = center_y;
n->width = width;
} else if (n->is_leaf && (with <= MIN_PARTICLE_SEPARATION)) {
// ... merge particles ...
} else {
if (n->is_leaf) { //! that is, if the node is not a branch
quadtree_insert_helper(n, n->p, center_x, center_y, width);
//! the node is now a branch
n->p = NULL;
n->is_leaf = 0;
}
quadtree_insert_helper(n, p, center_x, center_y, width);
}
return n;
}