diff --git a/src/hydro/Default/hydro.h b/src/hydro/Default/hydro.h index fd18233b122337036631d926d75b16ba6aa2e209..3215390e83db1e5f134768211c52d2dec32cecad 100644 --- a/src/hydro/Default/hydro.h +++ b/src/hydro/Default/hydro.h @@ -70,9 +70,10 @@ __attribute__((always_inline)) * and add the self-contribution term. * * @param p The particle to act upon + * @param time The current time */ __attribute__((always_inline)) - INLINE static void hydro_end_density(struct part* p) { + INLINE static void hydro_end_density(struct part* p, float time) { /* Some smoothing length multiples. */ const float h = p->h; @@ -166,12 +167,15 @@ __attribute__((always_inline)) * * @param p The particle * @param xp The extended data of the particle - * @param dt The time-step over which to drift + * @param t0 The time at the start of the drift + * @param t1 The time at the end of the drift */ __attribute__((always_inline)) INLINE static void hydro_predict_extra( - struct part* p, struct xpart* xp, float dt) { + struct part* p, struct xpart* xp, float t0, float t1) { float u, w; + const float dt = t1 - t0; + /* Predict internal energy */ w = p->force.u_dt / p->u * dt; if (fabsf(w) < 0.01f) /* 1st order expansion of exp(w) */ @@ -196,3 +200,22 @@ __attribute__((always_inline)) p->force.h_dt *= p->h * 0.333333333f; } + +/** + * @brief Kick the additional variables + * + * @param p The particle to act upon + * @param dt The time-step for this kick + */ +__attribute__((always_inline)) + INLINE static void hydro_kick_extra(struct part* p, float dt) {} + +/** + * @brief Converts hydro quantity of a particle + * + * Requires the density to be known + * + * @param p The particle to act upon + */ +__attribute__((always_inline)) + INLINE static void hydro_convert_quantities(struct part* p) {} diff --git a/src/hydro/Gadget2/hydro.h b/src/hydro/Gadget2/hydro.h index 53366bd84d6b102731f87a64b2af4191f99080bb..30d714e25392102d6b949bf5035c015a3005b038 100644 --- a/src/hydro/Gadget2/hydro.h +++ b/src/hydro/Gadget2/hydro.h @@ -203,7 +203,7 @@ __attribute__((always_inline)) INLINE static void hydro_predict_extra( */ __attribute__((always_inline)) INLINE static void hydro_end_force(struct part* p) { - + p->entropy_dt *= (const_hydro_gamma - 1.f) * powf(p->rho, -(const_hydro_gamma - 1.f)); } diff --git a/src/hydro/Gadget2/hydro_iact.h b/src/hydro/Gadget2/hydro_iact.h index 31528b63f5b3dc89ce533a3a7acafabee3c91332..260f78ba0ee6ecb1e19b169520f89808a01514db 100644 --- a/src/hydro/Gadget2/hydro_iact.h +++ b/src/hydro/Gadget2/hydro_iact.h @@ -197,7 +197,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_force( const float r_inv = 1.0f / r; /* Get some values in local variables. */ - // const float mi = pi->mass; + const float mi = pi->mass; const float mj = pj->mass; const float rhoi = pi->rho; const float rhoj = pj->rho; @@ -280,6 +280,10 @@ __attribute__((always_inline)) INLINE static void runner_iact_force( pj->a[1] += acc * dx[1]; pj->a[2] += acc * dx[2]; + /* Get the time derivative for h. */ + pi->force.h_dt -= mj * dvdr / rhoj * wi_dr; + pj->force.h_dt -= mi * dvdr / rhoi * wj_dr; + /* Update the signal velocity. */ pi->v_sig = fmaxf(pi->v_sig, v_sig); pj->v_sig = fmaxf(pj->v_sig, v_sig); @@ -366,6 +370,9 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_force( pi->a[1] -= acc * dx[1]; pi->a[2] -= acc * dx[2]; + /* Get the time derivative for h. */ + pi->force.h_dt -= mj * dvdr / rhoj * wi_dr; + /* Update the signal velocity. */ pi->v_sig = fmaxf(pi->v_sig, v_sig); diff --git a/src/runner.c b/src/runner.c index 66c1e02debee3394373e98a6e7ed41b639b5aba4..eeee039344b51dec5e4ada5f5b8d2f3e3616771f 100644 --- a/src/runner.c +++ b/src/runner.c @@ -851,8 +851,13 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) { if (is_fixdt || p->t_end <= t_current) { /* First, finish the force loop */ + p->force.h_dt *= p->h * 0.333333333f; + + /* And do the same of the extra variable */ hydro_end_force(p); + /* Now we are ready to compute the next time-step size */ + if (is_fixdt) { /* Now we have a time step, proceed with the kick */