diff --git a/src/engine.c b/src/engine.c
index b3f1b77e6116b87aaace357921d2bad19f638912..c0591a6f945dc18cc5feddac332bf53e56f79337 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -2503,19 +2503,18 @@ void engine_print_stats(struct engine *e) {
/* } */
#endif
- // const double e_tot = e_kin + e_int + e_pot;
+ const double E_tot = stats.E_kin + stats.E_int + stats.E_pot;
/* Print info */
- /* if (e->nodeID == 0) { */
- /* fprintf(e->file_stats, */
- /* " %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e
- * %14e " */
- /* "%14e\n", */
- /* e->time, mass, e_tot, e_kin, e_int, e_pot, e_rad, entropy,
- * mom[0], */
- /* mom[1], mom[2], ang_mom[0], ang_mom[1], ang_mom[2]); */
- /* fflush(e->file_stats); */
- /* } */
+ if (e->nodeID == 0) {
+ fprintf(e->file_stats,
+ " %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e %14e "
+ "%14e\n",
+ e->time, stats.mass, E_tot, stats.E_kin, stats.E_int, stats.E_pot,
+ stats.E_rad, stats.entropy, stats.mom[0], stats.mom[1],
+ stats.mom[2], stats.ang_mom[0], stats.ang_mom[1], stats.ang_mom[2]);
+ fflush(e->file_stats);
+ }
if (e->verbose)
message("took %.3f %s.", clocks_from_ticks(getticks() - tic),
@@ -2726,12 +2725,6 @@ void engine_step(struct engine *e) {
fflush(e->file_timesteps);
}
- /* Save some statistics */
- if (e->time - e->timeLastStatistics >= e->deltaTimeStatistics) {
- engine_print_stats(e);
- e->timeLastStatistics += e->deltaTimeStatistics;
- }
-
/* Drift only the necessary particles, that means all particles
* if we are about to repartition. */
int repart = (e->forcerepart != REPART_NONE);
@@ -2828,6 +2821,12 @@ void engine_step(struct engine *e) {
engine_launch(e, e->nr_threads, mask, submask);
TIMER_TOC(timer_runners);
+ /* Save some statistics */
+ if (e->time - e->timeLastStatistics >= e->deltaTimeStatistics) {
+ engine_print_stats(e);
+ e->timeLastStatistics += e->deltaTimeStatistics;
+ }
+
TIMER_TOC2(timer_step);
clocks_gettime(&time2);
diff --git a/src/statistics.c b/src/statistics.c
index dfa65516c2d1441b82f0f8ab9d7781d321ef6a32..ba36e0c44470fa8cdbdb51cf7e4e731620253c78 100644
--- a/src/statistics.c
+++ b/src/statistics.c
@@ -27,6 +27,7 @@
#include "statistics.h"
/* Local headers. */
+#include "atomic.h"
#include "cooling.h"
#include "engine.h"
#include "error.h"
@@ -38,8 +39,6 @@
*/
struct index_data {
const struct space *s;
- int chunk_size;
- int num_threads;
struct statistics *stats;
};
@@ -50,10 +49,9 @@ void stats_collect_part_mapper(void *map_data, int nr_parts, void *extra_data) {
const struct space *s = data->s;
struct part *restrict parts = (struct part *)map_data;
struct xpart *restrict xparts = s->xparts + (ptrdiff_t)(parts - s->parts);
- const int chunk_size = data->chunk_size;
- const int num_threads = data->num_threads;
const int ti_current = s->e->ti_current;
const double timeBase = s->e->timeBase;
+ struct statistics *const global_stats = data->stats;
/* Local accumulator */
struct statistics stats;
@@ -97,61 +95,35 @@ void stats_collect_part_mapper(void *map_data, int nr_parts, void *extra_data) {
stats.entropy += m * hydro_get_entropy(p, dt);
}
- /* Where do we write back ? */
- const int thread_id =
- ((ptrdiff_t)(parts - s->parts) / chunk_size) % num_threads;
-
/* Now write back to memory */
- data->stats[thread_id].E_kin += stats.E_kin;
- data->stats[thread_id].E_int += stats.E_int;
- data->stats[thread_id].E_pot += stats.E_pot;
- data->stats[thread_id].E_rad += stats.E_rad;
- data->stats[thread_id].entropy += stats.entropy;
- data->stats[thread_id].mass += stats.mass;
- data->stats[thread_id].mom[0] += stats.mom[0];
- data->stats[thread_id].mom[1] += stats.mom[1];
- data->stats[thread_id].mom[2] += stats.mom[2];
- data->stats[thread_id].ang_mom[0] += stats.ang_mom[0];
- data->stats[thread_id].ang_mom[1] += stats.ang_mom[1];
- data->stats[thread_id].ang_mom[2] += stats.ang_mom[2];
+ if (lock_lock(&global_stats->lock) == 0) {
+
+ global_stats->E_int += stats.E_int;
+ global_stats->E_pot += stats.E_pot;
+ global_stats->E_rad += stats.E_rad;
+ global_stats->entropy += stats.entropy;
+ global_stats->mass += stats.mass;
+ global_stats->mom[0] += stats.mom[0];
+ global_stats->mom[1] += stats.mom[1];
+ global_stats->mom[2] += stats.mom[2];
+ global_stats->ang_mom[0] += stats.ang_mom[0];
+ global_stats->ang_mom[1] += stats.ang_mom[1];
+ global_stats->ang_mom[2] += stats.ang_mom[2];
+ }
+ if (lock_unlock(&global_stats->lock) != 0)
+ error("Failed to unlock stats holder.");
}
void stats_collect(const struct space *s, struct statistics *stats) {
- const int num_threads = s->e->threadpool.num_threads;
- const int chunk_size = s->nr_parts / num_threads;
+ const int chunk_size = 1000;
/* Prepare the data */
struct index_data extra_data;
extra_data.s = s;
- extra_data.chunk_size = chunk_size;
- extra_data.num_threads = num_threads;
- extra_data.stats = malloc(sizeof(struct statistics) * num_threads);
- if (extra_data.stats == NULL)
- error("Impossible to allocate memory for statistics");
- bzero(&extra_data.stats, sizeof(struct statistics) * num_threads);
+ extra_data.stats = stats;
/* Run parallel collection of statistics */
threadpool_map(&s->e->threadpool, stats_collect_part_mapper, s->parts,
s->nr_parts, sizeof(struct part), chunk_size, &extra_data);
-
- /* Collect data from all the threads in the pool */
- for (int i = 0; i < num_threads; ++i) {
-
- stats->E_kin += extra_data.stats[i].E_kin;
- stats->E_int += extra_data.stats[i].E_int;
- stats->E_pot += extra_data.stats[i].E_pot;
- stats->E_rad += extra_data.stats[i].E_rad;
- stats->entropy += extra_data.stats[i].entropy;
- stats->mass += extra_data.stats[i].mass;
- stats->mom[0] += extra_data.stats[i].mom[0];
- stats->mom[1] += extra_data.stats[i].mom[1];
- stats->mom[2] += extra_data.stats[i].mom[2];
- stats->ang_mom[0] += extra_data.stats[i].ang_mom[0];
- stats->ang_mom[1] += extra_data.stats[i].ang_mom[1];
- stats->ang_mom[2] += extra_data.stats[i].ang_mom[2];
- }
-
- /* Free stuff */
- free(extra_data.stats);
}
diff --git a/src/statistics.h b/src/statistics.h
index 4931d33cb351d394237621a8b3ac79bdb5cc100e..1910e3e2becd79e482508c6d56608a848fa938b7 100644
--- a/src/statistics.h
+++ b/src/statistics.h
@@ -23,6 +23,7 @@
#include "../config.h"
/* Local headers. */
+#include "lock.h"
#include "space.h"
/**
@@ -53,6 +54,9 @@ struct statistics {
/*! Angular momentum */
double ang_mom[3];
+
+ /*! Lock for threaded access */
+ swift_lock_type lock;
};
void stats_collect(const struct space* s, struct statistics* stats);