From 35fd53ea3d0cd2af2d1ee0964e23a6627f64952a Mon Sep 17 00:00:00 2001
From: Bert Vandenbroucke <bert.vandenbroucke@ugent.be>
Date: Thu, 18 Aug 2016 09:47:42 +0100
Subject: [PATCH] Made a small but important change in computation of hydro
 acceleration. Added total energy to snapshots instead of thermal energy.

---
 src/hydro/Gizmo/hydro.h    |  6 +++---
 src/hydro/Gizmo/hydro_io.h | 36 ++++++++++++++++++++++++++++++++++--
 src/task.c                 |  8 ++++----
 3 files changed, 41 insertions(+), 9 deletions(-)

diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h
index a97f684e49..ae4bbfde86 100644
--- a/src/hydro/Gizmo/hydro.h
+++ b/src/hydro/Gizmo/hydro.h
@@ -343,13 +343,13 @@ __attribute__((always_inline)) INLINE static void hydro_end_force(
            this is indeed the case. */
   if (p->force.dt) {
     p->a_hydro[0] =
-        (p->conserved.momentum[0] / p->conserved.mass - p->primitives.v[0]) /
+        (p->conserved.momentum[0] / p->conserved.mass - p->force.v_full[0]) /
         p->force.dt;
     p->a_hydro[1] =
-        (p->conserved.momentum[1] / p->conserved.mass - p->primitives.v[1]) /
+        (p->conserved.momentum[1] / p->conserved.mass - p->force.v_full[1]) /
         p->force.dt;
     p->a_hydro[2] =
-        (p->conserved.momentum[2] / p->conserved.mass - p->primitives.v[2]) /
+        (p->conserved.momentum[2] / p->conserved.mass - p->force.v_full[2]) /
         p->force.dt;
 
     p->du_dt = p->conserved.flux.energy / p->force.dt;
diff --git a/src/hydro/Gizmo/hydro_io.h b/src/hydro/Gizmo/hydro_io.h
index aaae075fe8..8bd8f43626 100644
--- a/src/hydro/Gizmo/hydro_io.h
+++ b/src/hydro/Gizmo/hydro_io.h
@@ -55,14 +55,45 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
                                 UNIT_CONV_DENSITY, parts, primitives.rho);
 }
 
+/**
+ * @brief Get the internal energy of a particle
+ *
+ * @param e #engine.
+ * @param p Particle.
+ * @return Internal energy of the particle
+ */
 float convert_u(struct engine* e, struct part* p) {
   return p->primitives.P / hydro_gamma_minus_one / p->primitives.rho;
 }
 
+/**
+ * @brief Get the entropic function of a particle
+ *
+ * @param e #engine.
+ * @param p Particle.
+ * @return Entropic function of the particle
+ */
 float convert_A(struct engine* e, struct part* p) {
   return p->primitives.P / pow_gamma(p->primitives.rho);
 }
 
+/**
+ * @brief Get the total energy of a particle
+ *
+ * @param e #engine.
+ * @param p Particle.
+ * @return Total energy of the particle
+ */
+float convert_Etot(struct engine* e, struct part* p) {
+  float momentum2;
+
+  momentum2 = p->conserved.momentum[0] * p->conserved.momentum[0] +
+              p->conserved.momentum[1] * p->conserved.momentum[1] +
+              p->conserved.momentum[2] * p->conserved.momentum[2];
+
+  return p->conserved.energy + 0.5f * momentum2 / p->conserved.mass;
+}
+
 /**
  * @brief Specifies which particle fields to write to a dataset
  *
@@ -101,8 +132,9 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
       "Entropy", FLOAT, 1, UNIT_CONV_ENTROPY, parts, primitives.P, convert_A);
   list[11] = io_make_output_field("Pressure", FLOAT, 1, UNIT_CONV_PRESSURE,
                                   parts, primitives.P);
-  list[12] = io_make_output_field("TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY,
-                                  parts, conserved.energy);
+  list[12] =
+      io_make_output_field_convert_part("TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY,
+                                        parts, conserved.energy, convert_Etot);
 }
 
 /**
diff --git a/src/task.c b/src/task.c
index 37d6623e81..f80068d465 100644
--- a/src/task.c
+++ b/src/task.c
@@ -48,10 +48,10 @@
 
 /* Task type names. */
 const char *taskID_names[task_type_count] = {
-    "none",       "sort",       "self",    "pair",          "sub_self",
-    "sub_pair",   "init",       "ghost",   "extra_ghost",
-    "kick",       "kick_fixdt", "send",    "recv",          "grav_gather_m",
-    "grav_fft",   "grav_mm",    "grav_up", "grav_external"};
+    "none",       "sort",    "self",         "pair",          "sub_self",
+    "sub_pair",   "init",    "ghost",        "extra_ghost",   "kick",
+    "kick_fixdt", "send",    "recv",         "grav_gather_m", "grav_fft",
+    "grav_mm",    "grav_up", "grav_external"};
 
 const char *subtaskID_names[task_subtype_count] = {
     "none", "density", "gradient", "force", "grav", "tend"};
-- 
GitLab