Handle edge cases in laplace functions

Avoid division by zero.
Avoid taking the log of zero.
Silence clang type conversion warnings using round and trunc.
The existing values returned by the laplace functions do not change.

Add tests for laplace edge cases.
These changes pass the existing unit tests without modification.

Related to HS stats in #13192.

1 files changed, 24 insertions, 9 deletions

diff --git a/src/common/util.c b/src/common/util.c
index f8d1b7be4b..52b3e04946 100644
--- a/src/common/util.c
+++ b/src/common/util.c
@@ -536,13 +536,18 @@ int64_t
 sample_laplace_distribution(double mu, double b, double p)
 {
   double result;
-
   tor_assert(p >= 0.0 && p < 1.0);
+
   /* This is the "inverse cumulative distribution function" from:
    * http://en.wikipedia.org/wiki/Laplace_distribution */
+  if (p == 0.0) {
+    /* Avoid taking log(0.0) == -INFINITY, as some processors or compiler
+     * options can cause the program to trap. */
+    return INT64_MIN;
+  }
+
   result =  mu - b * (p > 0.5 ? 1.0 : -1.0)
                    * tor_mathlog(1.0 - 2.0 * fabs(p - 0.5));
-
   if (result >= INT64_MAX)
     return INT64_MAX;
   else if (result <= INT64_MIN)
@@ -551,18 +556,28 @@ sample_laplace_distribution(double mu, double b, double p)
     return (int64_t) result;
 }
 
-/** Add random noise between INT64_MIN and INT64_MAX coming from a
- * Laplace distribution with mu = 0 and b = <b>delta_f</b>/<b>epsilon</b>
- * to <b>signal</b> based on the provided <b>random</b> value in
- * [0.0, 1.0[. */
+/** Add random noise between INT64_MIN and INT64_MAX coming from a Laplace
+ * distribution with mu = 0 and b = <b>delta_f</b>/<b>epsilon</b> to
+ * <b>signal</b> based on the provided <b>random</b> value in [0.0, 1.0[.
+ * The epislon value must be between ]0.0, 1.0]. delta_f must be greater
+ * than 0. */
 int64_t
 add_laplace_noise(int64_t signal, double random, double delta_f,
                   double epsilon)
 {
-  int64_t noise = sample_laplace_distribution(
-               0.0, /* just add noise, no further signal */
-               delta_f / epsilon, random);
+  int64_t noise;
+
+  /* epsilon MUST be between ]0.0, 1.0] */
+  tor_assert(epsilon > 0.0 && epsilon <= 1.0);
+  /* delta_f MUST be greater than 0. */
+  tor_assert(delta_f > 0.0);
+
+  /* Just add noise, no further signal */
+  noise = sample_laplace_distribution(0.0,
+                                      delta_f / epsilon,
+                                      random);
 
+  /* Clip (signal + noise) to [INT64_MIN, INT64_MAX] */
   if (noise > 0 && INT64_MAX - noise < signal)
     return INT64_MAX;
   else if (noise < 0 && INT64_MIN - noise > signal)