From 6466e76c7a5fb8b0cee4e9e190f53f7b313fd564 Mon Sep 17 00:00:00 2001
From: Bert Vandenbroucke <bert.vandenbroucke@gmail.com>
Date: Fri, 8 Jun 2018 19:21:42 +0100
Subject: [PATCH] Cleaned up some more GizmoMFM code.
---
src/hydro/GizmoMFM/hydro_iact.h | 42 ++++++++-----------
.../GizmoMFM/hydro_slope_limiters_cell.h | 10 ++---
.../GizmoMFM/hydro_slope_limiters_face.h | 10 +++--
3 files changed, 29 insertions(+), 33 deletions(-)
diff --git a/src/hydro/GizmoMFM/hydro_iact.h b/src/hydro/GizmoMFM/hydro_iact.h
index 947dd3de01..d30aaa7252 100644
--- a/src/hydro/GizmoMFM/hydro_iact.h
+++ b/src/hydro/GizmoMFM/hydro_iact.h
@@ -57,7 +57,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_density(
const float r = sqrtf(r2);
/* Compute density of pi. */
- const float hi_inv = 1.f / hi;
+ const float hi_inv = 1.0f / hi;
const float xi = r * hi_inv;
kernel_deval(xi, &wi, &wi_dx);
@@ -75,7 +75,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_density(
pi->geometry.centroid[2] -= dx[2] * wi;
/* Compute density of pj. */
- const float hj_inv = 1.f / hj;
+ const float hj_inv = 1.0f / hj;
const float xj = r * hj_inv;
kernel_deval(xj, &wj, &wj_dx);
@@ -123,7 +123,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_density(
/* Get r and h inverse. */
const float r = sqrtf(r2);
- const float hi_inv = 1.f / hi;
+ const float hi_inv = 1.0f / hi;
const float xi = r * hi_inv;
kernel_deval(xi, &wi, &wi_dx);
@@ -220,7 +220,7 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
float r2, const float *dx, float hi, float hj, struct part *restrict pi,
struct part *restrict pj, int mode, float a, float H) {
- const float r_inv = 1.f / sqrtf(r2);
+ const float r_inv = 1.0f / sqrtf(r2);
const float r = r2 * r_inv;
/* Initialize local variables */
@@ -258,21 +258,17 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
} else
vmax = 0.0f;
- float dvdr = (pi->v[0] - pj->v[0]) * dx[0] + (pi->v[1] - pj->v[1]) * dx[1] +
- (pi->v[2] - pj->v[2]) * dx[2];
-
/* Velocity on the axis linking the particles */
- float dvdotdx = (Wi[1] - Wj[1]) * dx[0] + (Wi[2] - Wj[2]) * dx[1] +
- (Wi[3] - Wj[3]) * dx[2];
- dvdotdx = min(dvdotdx, dvdr);
+ float dvdr = (Wi[1] - Wj[1]) * dx[0] + (Wi[2] - Wj[2]) * dx[1] +
+ (Wi[3] - Wj[3]) * dx[2];
- /* We only care about this velocity for particles moving towards each others
+ /* We only care about this velocity for particles moving towards each other
*/
- dvdotdx = min(dvdotdx, 0.f);
+ const float dvdotdx = min(dvdr, 0.0f);
/* Get the signal velocity */
/* the magical factor 3 also appears in Gadget2 */
- vmax -= 3.f * dvdotdx * r_inv;
+ vmax -= 3.0f * dvdotdx * r_inv;
/* Store the signal velocity */
pi->timestepvars.vmax = max(pi->timestepvars.vmax, vmax);
@@ -280,14 +276,14 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
/* Compute kernel of pi. */
float wi, wi_dx;
- const float hi_inv = 1.f / hi;
+ const float hi_inv = 1.0f / hi;
const float hi_inv_dim = pow_dimension(hi_inv);
const float xi = r * hi_inv;
kernel_deval(xi, &wi, &wi_dx);
/* Compute kernel of pj. */
float wj, wj_dx;
- const float hj_inv = 1.f / hj;
+ const float hj_inv = 1.0f / hj;
const float hj_inv_dim = pow_dimension(hj_inv);
const float xj = r * hj_inv;
kernel_deval(xj, &wj, &wj_dx);
@@ -298,9 +294,9 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
const float wi_dr = hidp1 * wi_dx;
const float wj_dr = hjdp1 * wj_dx;
dvdr *= r_inv;
- if (pj->rho > 0.)
+ if (pj->rho > 0.0f)
pi->force.h_dt -= pj->conserved.mass * dvdr / pj->rho * wi_dr;
- if (mode == 1 && pi->rho > 0.)
+ if (mode == 1 && pi->rho > 0.0f)
pj->force.h_dt -= pi->conserved.mass * dvdr / pi->rho * wj_dr;
/* Compute (square of) area */
@@ -342,10 +338,10 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
/* if the interface has no area, nothing happens and we return */
/* continuing results in dividing by zero and NaN's... */
- if (Anorm2 == 0.f) return;
+ if (Anorm2 == 0.0f) return;
/* Compute the area */
- const float Anorm_inv = 1. / sqrtf(Anorm2);
+ const float Anorm_inv = 1.0f / sqrtf(Anorm2);
const float Anorm = Anorm2 * Anorm_inv;
#ifdef SWIFT_DEBUG_CHECKS
@@ -375,11 +371,9 @@ __attribute__((always_inline)) INLINE static void runner_iact_fluxes_common(
/* Compute interface velocity */
/* eqn. (9) */
- float xijdotdx = xij_i[0] * dx[0] + xij_i[1] * dx[1] + xij_i[2] * dx[2];
- xijdotdx *= r_inv * r_inv;
- const float vij[3] = {vi[0] + (vi[0] - vj[0]) * xijdotdx,
- vi[1] + (vi[1] - vj[1]) * xijdotdx,
- vi[2] + (vi[2] - vj[2]) * xijdotdx};
+ const float vij[3] = {vi[0] + (vi[0] - vj[0]) * xfac,
+ vi[1] + (vi[1] - vj[1]) * xfac,
+ vi[2] + (vi[2] - vj[2]) * xfac};
/* complete calculation of position of interface */
/* NOTE: dx is not necessarily just pi->x - pj->x but can also contain
diff --git a/src/hydro/GizmoMFM/hydro_slope_limiters_cell.h b/src/hydro/GizmoMFM/hydro_slope_limiters_cell.h
index a33d03dee8..329ca5fda4 100644
--- a/src/hydro/GizmoMFM/hydro_slope_limiters_cell.h
+++ b/src/hydro/GizmoMFM/hydro_slope_limiters_cell.h
@@ -108,7 +108,7 @@ __attribute__((always_inline)) INLINE static void hydro_slope_limit_cell(
gradtrue *= p->limiter.maxr;
const float gradmax = p->limiter.rho[1] - p->rho;
const float gradmin = p->rho - p->limiter.rho[0];
- const float gradtrue_inv = 1.f / gradtrue;
+ const float gradtrue_inv = 1.0f / gradtrue;
const float alpha =
min3(1.0f, gradmax * gradtrue_inv, gradmin * gradtrue_inv);
p->gradients.rho[0] *= alpha;
@@ -122,7 +122,7 @@ __attribute__((always_inline)) INLINE static void hydro_slope_limit_cell(
gradtrue *= p->limiter.maxr;
const float gradmax = p->limiter.v[0][1] - p->v[0];
const float gradmin = p->v[0] - p->limiter.v[0][0];
- const float gradtrue_inv = 1.f / gradtrue;
+ const float gradtrue_inv = 1.0f / gradtrue;
const float alpha =
min3(1.0f, gradmax * gradtrue_inv, gradmin * gradtrue_inv);
p->gradients.v[0][0] *= alpha;
@@ -136,7 +136,7 @@ __attribute__((always_inline)) INLINE static void hydro_slope_limit_cell(
gradtrue *= p->limiter.maxr;
const float gradmax = p->limiter.v[1][1] - p->v[1];
const float gradmin = p->v[1] - p->limiter.v[1][0];
- const float gradtrue_inv = 1.f / gradtrue;
+ const float gradtrue_inv = 1.0f / gradtrue;
const float alpha =
min3(1.0f, gradmax * gradtrue_inv, gradmin * gradtrue_inv);
p->gradients.v[1][0] *= alpha;
@@ -150,7 +150,7 @@ __attribute__((always_inline)) INLINE static void hydro_slope_limit_cell(
gradtrue *= p->limiter.maxr;
const float gradmax = p->limiter.v[2][1] - p->v[2];
const float gradmin = p->v[2] - p->limiter.v[2][0];
- const float gradtrue_inv = 1.f / gradtrue;
+ const float gradtrue_inv = 1.0f / gradtrue;
const float alpha =
min3(1.0f, gradmax * gradtrue_inv, gradmin * gradtrue_inv);
p->gradients.v[2][0] *= alpha;
@@ -164,7 +164,7 @@ __attribute__((always_inline)) INLINE static void hydro_slope_limit_cell(
gradtrue *= p->limiter.maxr;
const float gradmax = p->limiter.P[1] - p->P;
const float gradmin = p->P - p->limiter.P[0];
- const float gradtrue_inv = 1.f / gradtrue;
+ const float gradtrue_inv = 1.0f / gradtrue;
const float alpha =
min3(1.0f, gradmax * gradtrue_inv, gradmin * gradtrue_inv);
p->gradients.P[0] *= alpha;
diff --git a/src/hydro/GizmoMFM/hydro_slope_limiters_face.h b/src/hydro/GizmoMFM/hydro_slope_limiters_face.h
index 8f0636c992..a7c0e37d39 100644
--- a/src/hydro/GizmoMFM/hydro_slope_limiters_face.h
+++ b/src/hydro/GizmoMFM/hydro_slope_limiters_face.h
@@ -56,15 +56,17 @@ hydro_slope_limit_face_quantity(float phi_i, float phi_j, float phi_mid0,
float phiplus, phiminus, phi_mid;
- if (same_signf(phimax + delta1, phimax))
+ if (same_signf(phimax + delta1, phimax)) {
phiplus = phimax + delta1;
- else
+ } else {
phiplus = phimax / (1.0f + delta1 / (fabsf(phimax) + FLT_MIN));
+ }
- if (same_signf(phimin - delta1, phimin))
+ if (same_signf(phimin - delta1, phimin)) {
phiminus = phimin - delta1;
- else
+ } else {
phiminus = phimin / (1.0f + delta1 / (fabsf(phimin) + FLT_MIN));
+ }
if (phi_i < phi_j) {
const float temp = min(phibar + delta2, phi_mid0);
--
GitLab