I've been looking frantically for an answer to no avail; I've created a quadtree which is supposed to sort an array of more than 1000 objects declared by the struct:
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;
into a quadtree with the function:
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;
}
}
This is all done by:
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);
}
where "particles" is passed on on the form struct Particle* particles. And a particle is defined such that:
struct Particle {
double mass;
double x_pos;
double y_pos;
double x_vel;
double y_vel;
};
typedef struct Particle * Particle_structure;
The root is free'd after every iteration, after the for-loop, and the code works for samples smaller than ~200 and valgrind gives no errors for these. Added to this is a middle function, which performs some arithmetic on the particles:
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;
}
The ringer is that its rolls for a bit (clearing the root, and redoing it for new positions), so the number of variables in the recursion surely cannot be the problem(?). I'm pretty sure I'm out swimming with the fishes when it comes to some allocation of pointer declaration. Any thoughts/help would put you on the good side of Odin!
Edit: An example from how we find mass-center etc. The node-mass is done the same way.
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;
}