From 55fd3914ea05d64a6dd413fad0f5d933ed199ced Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Wed, 7 Sep 2016 19:07:13 +0100
Subject: [PATCH] Collect the radiated energy of all particles when computing
statistics and output the result in the log file.
---
src/cell.h | 2 +-
src/cooling/const_du/cooling.h | 16 +++++-
src/cooling/const_du/cooling_struct.h | 3 +-
src/cooling/const_lambda/cooling.h | 13 ++++-
src/cooling/const_lambda/cooling_struct.h | 3 +-
src/cooling/none/cooling.h | 13 +++++
src/engine.c | 62 ++++++++++++-----------
src/runner.c | 6 ++-
8 files changed, 81 insertions(+), 37 deletions(-)
diff --git a/src/cell.h b/src/cell.h
index fd83620624..8e0494fa9d 100644
--- a/src/cell.h
+++ b/src/cell.h
@@ -182,7 +182,7 @@ struct cell {
double mom[3], ang_mom[3];
/* Mass, potential, internal and kinetic energy of particles in this cell. */
- double mass, e_pot, e_int, e_kin, entropy;
+ double mass, e_pot, e_int, e_kin, e_rad, entropy;
/* Number of particles updated in this cell. */
int updated, g_updated;
diff --git a/src/cooling/const_du/cooling.h b/src/cooling/const_du/cooling.h
index 40bf2e919d..b25980ff22 100644
--- a/src/cooling/const_du/cooling.h
+++ b/src/cooling/const_du/cooling.h
@@ -81,7 +81,7 @@ __attribute__((always_inline)) INLINE static void cooling_cool_part(
hydro_set_internal_energy(p, u_new);
/* Store the radiated energy */
- xp->cooling_data.radiated_energy += u_new - u_old;
+ xp->cooling_data.radiated_energy += hydro_get_mass(p) * (u_old - u_new);
}
/**
@@ -123,6 +123,20 @@ __attribute__((always_inline)) INLINE static void cooling_init_part(
xp->cooling_data.radiated_energy = 0.f;
}
+/**
+ * @brief Returns the total radiated energy by this particle.
+ *
+ * In this simple example we jsut return the quantity accumulated in the
+ * #cooling_xpart_data of this particle.
+ *
+ * @param xp The extended particle data
+ */
+__attribute__((always_inline)) INLINE static float cooling_get_radiated_energy(
+ const struct xpart* restrict xp) {
+
+ return xp->cooling_data.radiated_energy;
+}
+
/**
* @brief Initialises the cooling function properties from the parameter file
*
diff --git a/src/cooling/const_du/cooling_struct.h b/src/cooling/const_du/cooling_struct.h
index 1355ae9fa4..cc00b001cf 100644
--- a/src/cooling/const_du/cooling_struct.h
+++ b/src/cooling/const_du/cooling_struct.h
@@ -53,8 +53,7 @@ struct cooling_function_data {
*/
struct cooling_xpart_data {
- /*! Amount of energy radiated away by this particle since the start of the run
- */
+ /*! Energy radiated away by this particle since the start of the run */
float radiated_energy;
};
diff --git a/src/cooling/const_lambda/cooling.h b/src/cooling/const_lambda/cooling.h
index 6d3f904d15..11cf2cae51 100644
--- a/src/cooling/const_lambda/cooling.h
+++ b/src/cooling/const_lambda/cooling.h
@@ -116,7 +116,7 @@ __attribute__((always_inline)) INLINE static void cooling_cool_part(
/* u_old, du_dt, dt, du_dt * dt, u_new, u_new_test); */
/* Store the radiated energy */
- xp->cooling_data.radiated_energy += u_new - u_old;
+ xp->cooling_data.radiated_energy += hydro_get_mass(p) * (u_old - u_new);
}
/**
@@ -154,6 +154,17 @@ __attribute__((always_inline)) INLINE static void cooling_init_part(
xp->cooling_data.radiated_energy = 0.f;
}
+/**
+ * @brief Returns the total radiated energy by this particle.
+ *
+ * @param xp The extended particle data
+ */
+__attribute__((always_inline)) INLINE static float cooling_get_radiated_energy(
+ const struct xpart* restrict xp) {
+
+ return xp->cooling_data.radiated_energy;
+}
+
/**
* @brief Initialises the cooling properties.
*
diff --git a/src/cooling/const_lambda/cooling_struct.h b/src/cooling/const_lambda/cooling_struct.h
index 5792e6be32..27c5df16bf 100644
--- a/src/cooling/const_lambda/cooling_struct.h
+++ b/src/cooling/const_lambda/cooling_struct.h
@@ -48,13 +48,12 @@ struct cooling_function_data {
float cooling_tstep_mult;
};
-
/**
* @brief Properties of the cooling stored in the particle data.
*/
struct cooling_xpart_data {
- /*! Amount of energy radiated away by this particle since the start of the run */
+ /*! Energy radiated away by this particle since the start of the run */
float radiated_energy;
};
diff --git a/src/cooling/none/cooling.h b/src/cooling/none/cooling.h
index 172781f066..0461100dc1 100644
--- a/src/cooling/none/cooling.h
+++ b/src/cooling/none/cooling.h
@@ -83,6 +83,19 @@ __attribute__((always_inline)) INLINE static float cooling_timestep(
__attribute__((always_inline)) INLINE static void cooling_init_part(
const struct part* restrict p, struct xpart* restrict xp) {}
+/**
+ * @brief Returns the total radiated energy by this particle.
+ *
+ * No cooling, so return 0.
+ *
+ * @param xp The extended particle data
+ */
+__attribute__((always_inline)) INLINE static float cooling_get_radiated_energy(
+ const struct xpart* restrict xp) {
+
+ return 0.f;
+}
+
/**
* @brief Initialises the cooling properties.
*
diff --git a/src/engine.c b/src/engine.c
index 0b4aa4f919..03e95b663a 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -2465,7 +2465,8 @@ void engine_print_stats(struct engine *e) {
const ticks tic = getticks();
const struct space *s = e->s;
- double e_kin = 0.0, e_int = 0.0, e_pot = 0.0, entropy = 0.0, mass = 0.0;
+ double e_kin = 0.0, e_int = 0.0, e_pot = 0.0, e_rad = 0.0;
+ double entropy = 0.0, mass = 0.0;
double mom[3] = {0.0, 0.0, 0.0}, ang_mom[3] = {0.0, 0.0, 0.0};
/* Collect the cell data. */
@@ -2476,6 +2477,7 @@ void engine_print_stats(struct engine *e) {
e_kin += c->e_kin;
e_int += c->e_int;
e_pot += c->e_pot;
+ e_rad += c->e_rad;
entropy += c->entropy;
mom[0] += c->mom[0];
mom[1] += c->mom[1];
@@ -2488,33 +2490,35 @@ void engine_print_stats(struct engine *e) {
/* Aggregate the data from the different nodes. */
#ifdef WITH_MPI
{
- double in[11] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
- double out[11];
+ double in[12] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
+ double out[12];
out[0] = e_kin;
out[1] = e_int;
out[2] = e_pot;
- out[3] = mom[0];
- out[4] = mom[1];
- out[5] = mom[2];
- out[6] = ang_mom[0];
- out[7] = ang_mom[1];
- out[8] = ang_mom[2];
- out[9] = mass;
- out[10] = entropy;
+ out[3] = e_rad;
+ out[4] = mom[0];
+ out[5] = mom[1];
+ out[6] = mom[2];
+ out[7] = ang_mom[0];
+ out[8] = ang_mom[1];
+ out[9] = ang_mom[2];
+ out[10] = mass;
+ out[11] = entropy;
if (MPI_Reduce(out, in, 11, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD) !=
MPI_SUCCESS)
error("Failed to aggregate stats.");
e_kin = out[0];
e_int = out[1];
e_pot = out[2];
- mom[0] = out[3];
- mom[1] = out[4];
- mom[2] = out[5];
- ang_mom[0] = out[6];
- ang_mom[1] = out[7];
- ang_mom[2] = out[8];
- mass = out[9];
- entropy = out[10];
+ e_rad = out[3];
+ mom[0] = out[4];
+ mom[1] = out[5];
+ mom[2] = out[6];
+ ang_mom[0] = out[7];
+ ang_mom[1] = out[8];
+ ang_mom[2] = out[9];
+ mass = out[10];
+ entropy = out[11];
}
#endif
@@ -2522,11 +2526,11 @@ void engine_print_stats(struct engine *e) {
/* Print info */
if (e->nodeID == 0) {
- fprintf(
- e->file_stats,
- " %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, entropy, mom[0], mom[1],
- mom[2], ang_mom[0], ang_mom[1], ang_mom[2]);
+ 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);
}
@@ -3358,11 +3362,11 @@ void engine_init(struct engine *e, struct space *s,
engine_default_energy_file_name);
sprintf(energyfileName + strlen(energyfileName), ".txt");
e->file_stats = fopen(energyfileName, "w");
- fprintf(
- e->file_stats,
- "#%14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s\n",
- "Time", "Mass", "E_tot", "E_kin", "E_int", "E_pot", "Entropy", "p_x",
- "p_y", "p_z", "ang_x", "ang_y", "ang_z");
+ fprintf(e->file_stats,
+ "#%14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s %14s "
+ "%14s\n",
+ "Time", "Mass", "E_tot", "E_kin", "E_int", "E_pot", "E_rad",
+ "Entropy", "p_x", "p_y", "p_z", "ang_x", "ang_y", "ang_z");
fflush(e->file_stats);
char timestepsfileName[200] = "";
diff --git a/src/runner.c b/src/runner.c
index eb2605915c..76dc9fd84d 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -734,7 +734,8 @@ static void runner_do_drift(struct cell *c, struct engine *e) {
const double dt = (ti_current - ti_old) * timeBase;
float dx_max = 0.f, dx2_max = 0.f, h_max = 0.f;
- double e_kin = 0.0, e_int = 0.0, e_pot = 0.0, entropy = 0.0, mass = 0.0;
+ double e_kin = 0.0, e_int = 0.0, e_pot = 0.0, e_rad = 0.0;
+ double entropy = 0.0, mass = 0.0;
double mom[3] = {0.0, 0.0, 0.0};
double ang_mom[3] = {0.0, 0.0, 0.0};
@@ -806,6 +807,7 @@ static void runner_do_drift(struct cell *c, struct engine *e) {
e_kin += 0.5 * m * (v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
e_pot += 0.;
e_int += m * hydro_get_internal_energy(p, half_dt);
+ e_rad += cooling_get_radiated_energy(xp);
/* Collect entropy */
entropy += m * hydro_get_entropy(p, half_dt);
@@ -832,6 +834,7 @@ static void runner_do_drift(struct cell *c, struct engine *e) {
e_kin += cp->e_kin;
e_int += cp->e_int;
e_pot += cp->e_pot;
+ e_rad += cp->e_rad;
entropy += cp->entropy;
mom[0] += cp->mom[0];
mom[1] += cp->mom[1];
@@ -849,6 +852,7 @@ static void runner_do_drift(struct cell *c, struct engine *e) {
c->e_kin = e_kin;
c->e_int = e_int;
c->e_pot = e_pot;
+ c->e_rad = e_rad;
c->entropy = entropy;
c->mom[0] = mom[0];
c->mom[1] = mom[1];
--
GitLab