diff --git a/src/kernel_long_gravity.h b/src/kernel_long_gravity.h
index ec31c2743079da22d1f3dd0c8683adf674aca1e3..2e0097ded217769af832af049b6c00b1305ac339 100644
--- a/src/kernel_long_gravity.h
+++ b/src/kernel_long_gravity.h
@@ -78,7 +78,7 @@ __attribute__((always_inline)) INLINE static void fourier_kernel_long_grav_eval(
 #else
   const double u = sqrt(u2);
   const double arg = M_PI_2 * u;
-  *W = arg / sinh(arg);
+  *W = arg / (sinh(arg) + FLT_MIN);
 #endif
 }
 
diff --git a/src/runner_doiact_fft.c b/src/runner_doiact_fft.c
index 26b59f9f6b864445df9190c6041ee684c456ba22..a8a882e4e765c2679b02347c96c2c7f1f49339f7 100644
--- a/src/runner_doiact_fft.c
+++ b/src/runner_doiact_fft.c
@@ -231,7 +231,7 @@ void runner_do_grav_fft(struct runner* r, int timer) {
         const int kz = (k > N_half ? k - N : k);
         const double kz_d = (double)kz;
         const double fz = k_fac * kz_d;
-        const double sinc_kz_inv = (kz != 0) ? fz / sin(fz) : 1.;
+        const double sinc_kz_inv = (kz != 0) ? fz / (sin(fz) + FLT_MIN) : 1.;
 
         /* Norm of vector in Fourier space */
         const double k2 = (kx_d * kx_d + ky_d * ky_d + kz_d * kz_d);
@@ -242,7 +242,7 @@ void runner_do_grav_fft(struct runner* r, int timer) {
         /* Green function */
         double W;
         fourier_kernel_long_grav_eval(k2 * a_smooth2, &W);
-        const double green_cor = green_fac * W / k2;
+        const double green_cor = green_fac * W / (k2 + FLT_MIN);
 
         /* Deconvolution of CIC */
         const double CIC_cor = sinc_kx_inv * sinc_ky_inv * sinc_kz_inv;