diff --git a/src/cell.h b/src/cell.h
index fd836206241c55cd6b9da2e29157c11a14c5a892..8e0494fa9db89b451a776140c18b57b54547120f 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 40bf2e919d84c278b1860dcc3360e887aa3be71c..b25980ff2269ca9ea176edcc2a3c771647819133 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 1355ae9fa4bc04f619e475fe6f9205f547bc76bc..cc00b001cf6b576266de02dac885f87d089bd8e4 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 6d3f904d15ac8eb693175037294fccf73bbf1b9b..11cf2cae51f1ab2646d1391d2164c399c77a7bba 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 5792e6be32e83db8905604a38dc2a58bbaf6e58c..27c5df16bffbe7d165237d201ca68ea4ba89dd73 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 172781f0668cf142cd4966c9f952dc7fabbe049e..0461100dc11e7ffbb4616766923142442b4ac943 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 0b4aa4f9197dc5a0cbdc00f15c19efaa58d3be2c..03e95b663a226ce5b5904513a70968390b3c8055 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 eb2605915c13622893cd06d2f8806254d784afba..76dc9fd84db2dab9bc3562693be3cc1413a1f051 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];