Commit 1c6b831c authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Internal energy computed correctly

parent 54523c79
......@@ -241,10 +241,13 @@ __attribute__((always_inline))
* @brief Returns the internal energy of a particle
*
* @param p The particle of interest
* @param dt Time since the last kick
*/
__attribute__((always_inline))
INLINE static float hydro_get_internal_energy(struct part* p) {
__attribute__((always_inline)) INLINE static float hydro_get_internal_energy(
const struct part* p, float dt) {
const float entropy = p->entropy + p->entropy_dt * dt;
return p->entropy * powf(p->rho, const_hydro_gamma - 1.f) *
return entropy * powf(p->rho, const_hydro_gamma - 1.f) *
(1.f / (const_hydro_gamma - 1.f));
}
......@@ -851,7 +851,7 @@ void runner_do_drift(struct runner *r, struct cell *c, int timer) {
/* Collect energies. */
e_kin += 0.5 * m * (v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
e_pot += 0.;
e_int += 0.;
e_int += m * hydro_get_internal_energy(p, half_dt);
}
/* Now, get the maximal particle motion from its square */
......@@ -872,16 +872,16 @@ void runner_do_drift(struct runner *r, struct cell *c, int timer) {
/* Collect */
dx_max = fmaxf(dx_max, cp->dx_max);
h_max = fmaxf(h_max, cp->h_max);
mass += cp->mass;
e_kin += cp->e_kin;
e_int += cp->e_int;
e_pot += cp->e_pot;
mom[0] += cp->mom[0];
mom[1] += cp->mom[1];
mom[2] += cp->mom[2];
ang_mom[0] += cp->ang_mom[0];
ang_mom[1] += cp->ang_mom[1];
ang_mom[2] += cp->ang_mom[2];
mass += cp->mass;
e_kin += cp->e_kin;
e_int += cp->e_int;
e_pot += cp->e_pot;
mom[0] += cp->mom[0];
mom[1] += cp->mom[1];
mom[2] += cp->mom[2];
ang_mom[0] += cp->ang_mom[0];
ang_mom[1] += cp->ang_mom[1];
ang_mom[2] += cp->ang_mom[2];
}
}
......@@ -898,7 +898,7 @@ void runner_do_drift(struct runner *r, struct cell *c, int timer) {
c->ang_mom[0] = ang_mom[0];
c->ang_mom[1] = ang_mom[1];
c->ang_mom[2] = ang_mom[2];
if (timer) TIMER_TOC(timer_drift);
}
......
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