Skip to content
GitLab
Commit 1a500878 authored by Pedro Gonnet's avatar Pedro Gonnet
Browse files

still doesn't work, crashes after two steps with undrifted cells.

parent 190a201e
Hide whitespace changes
Inline Side-by-side
src/cell.c
View file @ 1a500878
......@@ -468,7 +468,8 @@ void cell_gunlocktree(struct cell *c) {
* @param buff A buffer with at least max(c->count, c->gcount) entries,
* used for sorting indices.
*/
void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
void cell_split(struct cell *c, ptrdiff_t parts_offset, struct cell_buff *buff,
struct cell_buff *gbuff) {
const int count = c->count, gcount = c->gcount;
struct part *parts = c->parts;
......@@ -480,18 +481,22 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
int bucket_count[8] = {0, 0, 0, 0, 0, 0, 0, 0};
int bucket_offset[9];
/* If the buff is NULL, allocate it, and remember to free it. */
const int allocate_buffer = (buff == NULL);
if (allocate_buffer &&
(buff = (int *)malloc(sizeof(int) * max(count, gcount))) == NULL)
error("Failed to allocate temporary indices.");
#ifdef SWIFT_DEBUG_CHECKS
/* Check that the buffs are OK. */
for (int k = 0; k < count; k++) {
if (buff[k].x[0] != (float)parts[k].x[0] ||
buff[k].x[1] != (float)parts[k].x[1] ||
buff[k].x[2] != (float)parts[k].x[2])
error("Inconsistent buff contents.");
}
#endif /* SWIFT_DEBUG_CHECKS */
/* Fill the buffer with the indices. */
for (int k = 0; k < count; k++) {
const int bid = (parts[k].x[0] > pivot[0]) * 4 +
(parts[k].x[1] > pivot[1]) * 2 + (parts[k].x[2] > pivot[2]);
const int bid = (buff[k].x[0] > pivot[0]) * 4 +
(buff[k].x[1] > pivot[1]) * 2 + (buff[k].x[2] > pivot[2]);
bucket_count[bid]++;
buff[k] = bid;
buff[k].ind = bid;
}
/* Set the buffer offsets. */
......@@ -505,23 +510,25 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
for (int bucket = 0; bucket < 8; bucket++) {
for (int k = bucket_offset[bucket] + bucket_count[bucket];
k < bucket_offset[bucket + 1]; k++) {
int bid = buff[k];
int bid = buff[k].ind;
if (bid != bucket) {
struct part part = parts[k];
struct xpart xpart = xparts[k];
struct cell_buff temp_buff = buff[k];
while (bid != bucket) {
int j = bucket_offset[bid] + bucket_count[bid]++;
while (buff[j] == bid) {
while (buff[j].ind == bid) {
j++;
bucket_count[bid]++;
}
memswap(&parts[j], &part, sizeof(struct part));
memswap(&xparts[j], &xpart, sizeof(struct xpart));
memswap(&buff[j], &bid, sizeof(int));
memswap(&buff[j], &temp_buff, sizeof(struct cell_buff));
bid = temp_buff.ind;
}
parts[k] = part;
xparts[k] = xpart;
buff[k] = bid;
buff[k] = temp_buff;
}
bucket_count[bid]++;
}
......@@ -538,6 +545,14 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
if (count > 0 && gcount > 0) part_relink_gparts(parts, count, parts_offset);
#ifdef SWIFT_DEBUG_CHECKS
/* Check that the buffs are OK. */
for (int k = 0; k < count; k++) {
if (buff[k].x[0] != (float)parts[k].x[0] ||
buff[k].x[1] != (float)parts[k].x[1] ||
buff[k].x[2] != (float)parts[k].x[2])
error("Inconsistent buff contents (k=%i).", k);
}
/* Verify that _all_ the parts have been assigned to a cell. */
for (int k = 1; k < 8; k++)
if (&c->progeny[k - 1]->parts[c->progeny[k - 1]->count] !=
......@@ -564,6 +579,31 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
c->progeny[2]->parts[k].x[1] <= pivot[1] ||
c->progeny[2]->parts[k].x[2] > pivot[2])
error("Sorting failed (progeny=2).");
for (int k = 0; k < c->progeny[3]->count; k++)
if (c->progeny[3]->parts[k].x[0] > pivot[0] ||
c->progeny[3]->parts[k].x[1] <= pivot[1] ||
c->progeny[3]->parts[k].x[2] <= pivot[2])
error("Sorting failed (progeny=3).");
for (int k = 0; k < c->progeny[4]->count; k++)
if (c->progeny[4]->parts[k].x[0] <= pivot[0] ||
c->progeny[4]->parts[k].x[1] > pivot[1] ||
c->progeny[4]->parts[k].x[2] > pivot[2])
error("Sorting failed (progeny=4).");
for (int k = 0; k < c->progeny[5]->count; k++)
if (c->progeny[5]->parts[k].x[0] <= pivot[0] ||
c->progeny[5]->parts[k].x[1] > pivot[1] ||
c->progeny[5]->parts[k].x[2] <= pivot[2])
error("Sorting failed (progeny=5).");
for (int k = 0; k < c->progeny[6]->count; k++)
if (c->progeny[6]->parts[k].x[0] <= pivot[0] ||
c->progeny[6]->parts[k].x[1] <= pivot[1] ||
c->progeny[6]->parts[k].x[2] > pivot[2])
error("Sorting failed (progeny=6).");
for (int k = 0; k < c->progeny[7]->count; k++)
if (c->progeny[7]->parts[k].x[0] <= pivot[0] ||
c->progeny[7]->parts[k].x[1] <= pivot[1] ||
c->progeny[7]->parts[k].x[2] <= pivot[2])
error("Sorting failed (progeny=7).");
#endif
/* Now do the same song and dance for the gparts. */
......@@ -571,11 +611,10 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
/* Fill the buffer with the indices. */
for (int k = 0; k < gcount; k++) {
const int bid = (gparts[k].x[0] > pivot[0]) * 4 +
(gparts[k].x[1] > pivot[1]) * 2 +
(gparts[k].x[2] > pivot[2]);
const int bid = (gbuff[k].x[0] > pivot[0]) * 4 +
(gbuff[k].x[1] > pivot[1]) * 2 + (gbuff[k].x[2] > pivot[2]);
bucket_count[bid]++;
buff[k] = bid;
gbuff[k].ind = bid;
}
/* Set the buffer offsets. */
......@@ -589,20 +628,22 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff) {
for (int bucket = 0; bucket < 8; bucket++) {
for (int k = bucket_offset[bucket] + bucket_count[bucket];
k < bucket_offset[bucket + 1]; k++) {
int bid = buff[k];
int bid = gbuff[k].ind;
if (bid != bucket) {
struct gpart gpart = gparts[k];
struct cell_buff temp_buff = gbuff[k];
while (bid != bucket) {
int j = bucket_offset[bid] + bucket_count[bid]++;
while (buff[j] == bid) {
while (gbuff[j].ind == bid) {
j++;
bucket_count[bid]++;
}
memswap(&gparts[j], &gpart, sizeof(struct gpart));
memswap(&buff[j], &bid, sizeof(int));
memswap(&gbuff[j], &temp_buff, sizeof(struct cell_buff));
bid = temp_buff.ind;
}
gparts[k] = gpart;
buff[k] = bid;
gbuff[k] = temp_buff;
}
bucket_count[bid]++;
}
......
src/cell.h
View file @ 1a500878
......@@ -52,6 +52,12 @@ struct scheduler;
/* Global variables. */
extern int cell_next_tag;
/* Struct to temporarily buffer the particle locations and bin id. */
struct cell_buff {
double x[3];
int ind;
};
/* Mini struct to link cells to tasks. Used as a linked list. */
struct link {
......@@ -272,7 +278,8 @@ struct cell {
((int)(k) + (cdim)[2] * ((int)(j) + (cdim)[1] * (int)(i)))
/* Function prototypes. */
void cell_split(struct cell *c, ptrdiff_t parts_offset, int *buff);
void cell_split(struct cell *c, ptrdiff_t parts_offset, struct cell_buff *buff,
struct cell_buff *gbuff);
void cell_sanitize(struct cell *c);
int cell_locktree(struct cell *c);
void cell_unlocktree(struct cell *c);
......
src/debug.c
View file @ 1a500878
......@@ -227,21 +227,21 @@ int checkCellhdxmax(const struct cell *c, int *depth) {
if (c->h_max != h_max) {
message("%d Inconsistent h_max: cell %f != parts %f", *depth, c->h_max,
h_max);
message("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
error("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
result = 0;
}
if (c->dx_max != dx_max) {
message("%d Inconsistent dx_max: %f != %f", *depth, c->dx_max, dx_max);
message("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
error("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
result = 0;
}
/* Check rebuild criterion. */
if (h_max > c->dmin) {
/* if (h_max > c->dmin) {
message("%d Inconsistent c->dmin: %f > %f", *depth, h_max, c->dmin);
message("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
error("location: %f %f %f", c->loc[0], c->loc[1], c->loc[2]);
result = 0;
}
} */
return result;
}
......
src/space.c
View file @ 1a500878
......@@ -1456,9 +1456,12 @@ void space_map_cells_pre(struct space *s, int full,
* @param s The #space in which the cell lives.
* @param c The #cell to split recursively.
* @param buff A buffer for particle sorting, should be of size at least
* max(c->count, c->gount) or @c NULL.
* c->count or @c NULL.
* @param gbuff A buffer for particle sorting, should be of size at least
* c->gcount or @c NULL.
*/
void space_split_recursive(struct space *s, struct cell *c, int *buff) {
void space_split_recursive(struct space *s, struct cell *c,
struct cell_buff *buff, struct cell_buff *gbuff) {
const int count = c->count;
const int gcount = c->gcount;
......@@ -1473,10 +1476,29 @@ void space_split_recursive(struct space *s, struct cell *c, int *buff) {
struct engine *e = s->e;
/* If the buff is NULL, allocate it, and remember to free it. */
const int allocate_buffer = (buff == NULL);
if (allocate_buffer &&
(buff = (int *)malloc(sizeof(int) * max(count, gcount))) == NULL)
error("Failed to allocate temporary indices.");
const int allocate_buffer = (buff == NULL && gbuff == NULL);
if (allocate_buffer) {
if (count > 0) {
if ((buff = (struct cell_buff *)malloc(sizeof(struct cell_buff) *
count)) == NULL)
error("Failed to allocate temporary indices.");
for (int k = 0; k < count; k++) {
buff[k].x[0] = parts[k].x[0];
buff[k].x[1] = parts[k].x[1];
buff[k].x[2] = parts[k].x[2];
}
}
if (gcount > 0) {
if ((gbuff = (struct cell_buff *)malloc(sizeof(struct cell_buff) *
gcount)) == NULL)
error("Failed to allocate temporary indices.");
for (int k = 0; k < gcount; k++) {
gbuff[k].x[0] = gparts[k].x[0];
gbuff[k].x[1] = gparts[k].x[1];
gbuff[k].x[2] = gparts[k].x[2];
}
}
}
/* Check the depth. */
while (depth > (maxdepth = s->maxdepth)) {
......@@ -1522,15 +1544,18 @@ void space_split_recursive(struct space *s, struct cell *c, int *buff) {
}
/* Split the cell data. */
cell_split(c, c->parts - s->parts, buff);
cell_split(c, c->parts - s->parts, buff, gbuff);
/* Remove any progeny with zero parts. */
struct cell_buff *progeny_buff = buff, *progeny_gbuff = gbuff;
for (int k = 0; k < 8; k++)
if (c->progeny[k]->count == 0 && c->progeny[k]->gcount == 0) {
space_recycle(s, c->progeny[k]);
c->progeny[k] = NULL;
} else {
space_split_recursive(s, c->progeny[k], buff);
space_split_recursive(s, c->progeny[k], progeny_buff, progeny_gbuff);
progeny_buff += c->progeny[k]->count;
progeny_gbuff += c->progeny[k]->gcount;
h_max = max(h_max, c->progeny[k]->h_max);
ti_end_min = min(ti_end_min, c->progeny[k]->ti_end_min);
ti_end_max = max(ti_end_max, c->progeny[k]->ti_end_max);
......@@ -1589,7 +1614,10 @@ void space_split_recursive(struct space *s, struct cell *c, int *buff) {
c->owner = 0; /* Ok, there is really nothing on this rank... */
/* Clean up. */
if (allocate_buffer) free(buff);
if (allocate_buffer) {
if (buff != NULL) free(buff);
if (gbuff != NULL) free(gbuff);
}
}
/**
......@@ -1608,7 +1636,7 @@ void space_split_mapper(void *map_data, int num_cells, void *extra_data) {
for (int ind = 0; ind < num_cells; ind++) {
struct cell *c = &cells_top[ind];
space_split_recursive(s, c, NULL);
space_split_recursive(s, c, NULL, NULL);
}
#ifdef SWIFT_DEBUG_CHECKS
......@@ -1616,7 +1644,7 @@ void space_split_mapper(void *map_data, int num_cells, void *extra_data) {
for (int ind = 0; ind < num_cells; ind++) {
int depth = 0;
if (!checkCellhdxmax(&cells_top[ind], &depth))
message(" at cell depth %d", depth);
error(" at cell depth %d", depth);
}
#endif
}
......
src/threadpool.c
View file @ 1a500878
......@@ -90,7 +90,7 @@ void threadpool_init(struct threadpool *tp, int num_threads) {
/* Initialize the thread counters. */
tp->num_threads = num_threads;
tp->num_threads_waiting = 0;
/* If there is only a single thread, do nothing more as of here as
we will just do work in the (blocked) calling thread. */
if (num_threads == 1) return;
......@@ -147,7 +147,7 @@ void threadpool_init(struct threadpool *tp, int num_threads) {
void threadpool_map(struct threadpool *tp, threadpool_map_function map_function,
void *map_data, size_t N, int stride, int chunk,
void *extra_data) {
/* If we just have a single thread, call the map function directly. */
if (tp->num_threads == 1) {
map_function(map_data, N, extra_data);
......
  • Matthieu Schaller @matthieu ·
    Owner

    Is the crash on the Sedov blast with --enable-debugging-checks ? If so, the problem is unrelated to this branch. Master crashes as well. It's a drift problem, fixed in !292 (merged).

  • Pedro Gonnet @nnrw56 ·
    Author Developer

    It's on EAGLE_12 with --enable-debugging-checks using a single thread.

  • Matthieu Schaller @matthieu ·
    Owner

    Getting

    [00065.7] active.h:cell_is_drifted():53: Cell has not been drifted to the current time c->ti_old=0, e->ti_current=1024

    ?

    Then it's the same issue. We are not drifting correctly.

Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment