From 0ad67ed854ad332474d90ab80a7b653cb74cc39a Mon Sep 17 00:00:00 2001
From: Raymond Hettinger <raymond.hettinger@gmail.com>
Date: Sat, 14 Sep 2019 11:01:35 -0700
Subject: [PATCH] bpo-37798: Prevent undefined behavior in direct calls to the
 C helper function.

While the _statistics helper is private and not intended to be called
directly, it is possible to do so.  That would bypass the parameter
validation checks in the pure python calling code.

Accordingly, we need to add parameter checks to avoid undefined
behavior (division by zero) or errors in the underlying functions
(negative square root or non-positive argument for logarithm).
---
 Modules/_statisticsmodule.c | 18 +++++++++++++++++-
 1 file changed, 17 insertions(+), 1 deletion(-)

diff --git a/Modules/_statisticsmodule.c b/Modules/_statisticsmodule.c
index 16a992485237abe..a646e96d0165c08 100644
--- a/Modules/_statisticsmodule.c
+++ b/Modules/_statisticsmodule.c
@@ -32,8 +32,11 @@ _statistics__normal_dist_inv_cdf_impl(PyObject *module, double p, double mu,
 /*[clinic end generated code: output=02fd19ddaab36602 input=24715a74be15296a]*/
 {
     double q, num, den, r, x;
+    if (p <= 0.0 || p >= 1.0 || sigma <= 0.0) {
+        goto error;
+    }
+
     q = p - 0.5;
-    // Algorithm AS 241: The Percentage Points of the Normal Distribution
     if(fabs(q) <= 0.425) {
         r = 0.180625 - q * q;
         // Hash sum-55.8831928806149014439
@@ -53,10 +56,16 @@ _statistics__normal_dist_inv_cdf_impl(PyObject *module, double p, double mu,
                      6.8718700749205790830e+2) * r +
                      4.2313330701600911252e+1) * r +
                      1.0);
+        if (den == 0.0) {
+            goto error;
+        }
         x = num / den;
         return mu + (x * sigma);
     }
     r = (q <= 0.0) ? p : (1.0 - p);
+    if (r <= 0.0 || r >= 1.0) {
+        goto error;
+    }
     r = sqrt(-log(r));
     if (r <= 5.0) {
         r = r - 1.6;
@@ -97,11 +106,18 @@ _statistics__normal_dist_inv_cdf_impl(PyObject *module, double p, double mu,
                      5.99832206555887937690e-1) * r +
                      1.0);
     }
+    if (den == 0.0) {
+        goto error;
+    }
     x = num / den;
     if (q < 0.0) {
         x = -x;
     }
     return mu + (x * sigma);
+
+  error:
+    PyErr_SetString(PyExc_ValueError, "inv_cdf undefined for these parameters");
+    return -1.0;
 }