From 4a20260bb26a261312e36a25f6ecb503e8d3fbe0 Mon Sep 17 00:00:00 2001
From: Matthieu Schaller <matthieu.schaller@durham.ac.uk>
Date: Wed, 16 Aug 2017 15:00:41 +0100
Subject: [PATCH] Use floating-point precision throughout the M2L kernel
evaluation.
---
src/gravity_derivatives.h | 14 ++++++++++++++
src/gravity_properties.c | 2 +-
src/gravity_properties.h | 6 +++---
src/multipole.h | 16 ++++++++--------
4 files changed, 26 insertions(+), 12 deletions(-)
diff --git a/src/gravity_derivatives.h b/src/gravity_derivatives.h
index f6eb888cf3..e808f53454 100644
--- a/src/gravity_derivatives.h
+++ b/src/gravity_derivatives.h
@@ -85,6 +85,20 @@ struct potential_derivatives {
#endif
};
+/**
+ * @brief Compute all the relevent derivatives of the softened and truncated
+ * gravitational potential.
+ *
+ * @param r_x x-component of distance vector
+ * @param r_y y-component of distance vector
+ * @param r_z z-component of distance vector
+ * @param r2 Square norm of distance vector
+ * @param r_inv Inverse norm of distance vector
+ * @param eps Softening length.
+ * @param eps2 Square of softening length.
+ * @param eps_inv Inverse of softening length.
+ * @param pot (return) The structure containing all the derivatives.
+ */
__attribute__((always_inline)) INLINE static void compute_potential_derivatives(
float r_x, float r_y, float r_z, float r2, float r_inv, float eps,
float eps2, float eps_inv, struct potential_derivatives *pot) {
diff --git a/src/gravity_properties.c b/src/gravity_properties.c
index 8110b218eb..4782875ff9 100644
--- a/src/gravity_properties.c
+++ b/src/gravity_properties.c
@@ -59,7 +59,7 @@ void gravity_props_init(struct gravity_props *p,
/* Softening lengths */
p->epsilon = 3. * parser_get_param_double(params, "Gravity:epsilon");
p->epsilon2 = p->epsilon * p->epsilon;
- p->epsilon_inv = 1. / p->epsilon;
+ p->epsilon_inv = 1.f / p->epsilon;
}
void gravity_props_print(const struct gravity_props *p) {
diff --git a/src/gravity_properties.h b/src/gravity_properties.h
index c663f628b4..886f29ff08 100644
--- a/src/gravity_properties.h
+++ b/src/gravity_properties.h
@@ -58,13 +58,13 @@ struct gravity_props {
double theta_crit_inv;
/*! Softening length */
- double epsilon;
+ float epsilon;
/*! Square of softening length */
- double epsilon2;
+ float epsilon2;
/*! Inverse of softening length */
- double epsilon_inv;
+ float epsilon_inv;
};
void gravity_props_print(const struct gravity_props *p);
diff --git a/src/multipole.h b/src/multipole.h
index 7eeccff73b..fc7866a289 100644
--- a/src/multipole.h
+++ b/src/multipole.h
@@ -1516,14 +1516,14 @@ INLINE static void gravity_M2L(struct grav_tensor *l_b,
const double dim[3]) {
/* Recover some constants */
- const double eps = props->epsilon;
- const double eps2 = props->epsilon2;
- const double eps_inv = props->epsilon_inv;
+ const float eps = props->epsilon;
+ const float eps2 = props->epsilon2;
+ const float eps_inv = props->epsilon_inv;
/* Compute distance vector */
- double dx = pos_b[0] - pos_a[0];
- double dy = pos_b[1] - pos_a[1];
- double dz = pos_b[2] - pos_a[2];
+ float dx = (float)(pos_b[0] - pos_a[0]);
+ float dy = (float)(pos_b[1] - pos_a[1]);
+ float dz = (float)(pos_b[2] - pos_a[2]);
/* Apply BC */
if (periodic) {
@@ -1533,8 +1533,8 @@ INLINE static void gravity_M2L(struct grav_tensor *l_b,
}
/* Compute distance */
- const double r2 = dx * dx + dy * dy + dz * dz;
- const double r_inv = 1. / sqrt(r2);
+ const float r2 = dx * dx + dy * dy + dz * dz;
+ const float r_inv = 1. / sqrtf(r2);
/* Compute all derivatives */
struct potential_derivatives pot;
--
GitLab