From 9fe67486cc06a02308a169a33ecbbf841017a27d Mon Sep 17 00:00:00 2001
From: Arthur O'Dwyer <arthur.j.odwyer@gmail.com>
Date: Wed, 1 Dec 2021 19:55:26 -0500
Subject: [PATCH] [libc++] [test] Improve the test for
 `uniform_int_distribution<T>`.

Extracted from https://reviews.llvm.org/D114920
---
 .../rand.dist.uni/rand.dist.uni.int/eval.pass.cpp  | 532 +++++----------------
 1 file changed, 107 insertions(+), 425 deletions(-)

diff --git a/libcxx/test/std/numerics/rand/rand.dis/rand.dist.uni/rand.dist.uni.int/eval.pass.cpp b/libcxx/test/std/numerics/rand/rand.dis/rand.dist.uni/rand.dist.uni.int/eval.pass.cpp
index e3ab0ab9..dd90b79 100644
--- a/libcxx/test/std/numerics/rand/rand.dis/rand.dist.uni/rand.dist.uni.int/eval.pass.cpp
+++ b/libcxx/test/std/numerics/rand/rand.dis/rand.dist.uni/rand.dist.uni.int/eval.pass.cpp
@@ -17,443 +17,125 @@
 
 #include <random>
 #include <cassert>
-#include <vector>
-#include <numeric>
+#include <climits>
 #include <cstddef>
+#include <limits>
+#include <numeric>
+#include <vector>
 
 #include "test_macros.h"
 
 template <class T>
-inline
-T
-sqr(T x)
+T sqr(T x)
 {
     return x * x;
 }
 
-int main(int, char**)
+template <class ResultType, class EngineType>
+void test_statistics(ResultType a, ResultType b)
 {
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::minstd_rand0 G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::minstd_rand G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::mt19937 G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::mt19937_64 G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::ranlux24_base G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::ranlux48_base G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::ranlux24 G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::ranlux48 G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
-    }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::knuth_b G;
-        G g;
-        D d;
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
+    ASSERT_SAME_TYPE(typename std::uniform_int_distribution<ResultType>::result_type, ResultType);
+
+    EngineType g;
+    std::uniform_int_distribution<ResultType> dist(a, b);
+    assert(dist.a() == a);
+    assert(dist.b() == b);
+    std::vector<ResultType> u;
+    for (int i = 0; i < 10000; ++i) {
+        ResultType v = dist(g);
+        assert(a <= v && v <= b);
+        u.push_back(v);
     }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::minstd_rand0 G;
-        G g;
-        D d(-6, 106);
-        for (int i = 0; i < 10000; ++i)
-        {
-            int u = d(g);
-            assert(-6 <= u && u <= 106);
-        }
+
+    // Quick check: The chance of getting *no* hits in any given tenth of the range
+    // is (0.9)^10000, or "ultra-astronomically low."
+    bool bottom_tenth = false;
+    bool top_tenth = false;
+    for (std::size_t i = 0; i < u.size(); ++i) {
+        bottom_tenth = bottom_tenth || (u[i] <= (a + (b / 10) - (a / 10)));
+        top_tenth = top_tenth || (u[i] >= (b - (b / 10) + (a / 10)));
     }
-    {
-        typedef std::uniform_int_distribution<> D;
-        typedef std::minstd_rand G;
-        G g;
-        D d(5, 100);
-        const int N = 100000;
-        std::vector<D::result_type> u;
-        for (int i = 0; i < N; ++i)
-        {
-            D::result_type v = d(g);
-            assert(d.a() <= v && v <= d.b());
-            u.push_back(v);
-        }
-        double mean = std::accumulate(u.begin(), u.end(),
-                                              double(0)) / u.size();
-        double var = 0;
-        double skew = 0;
-        double kurtosis = 0;
-        for (std::size_t i = 0; i < u.size(); ++i)
-        {
-            double dbl = (u[i] - mean);
-            double d2 = sqr(dbl);
-            var += d2;
-            skew += dbl * d2;
-            kurtosis += d2 * d2;
-        }
-        var /= u.size();
-        double dev = std::sqrt(var);
-        skew /= u.size() * dev * var;
-        kurtosis /= u.size() * var * var;
-        kurtosis -= 3;
-        double x_mean = ((double)d.a() + d.b()) / 2;
-        double x_var = (sqr((double)d.b() - d.a() + 1) - 1) / 12;
-        double x_skew = 0;
-        double x_kurtosis = -6. * (sqr((double)d.b() - d.a() + 1) + 1) /
-                            (5. * (sqr((double)d.b() - d.a() + 1) - 1));
-        assert(std::abs((mean - x_mean) / x_mean) < 0.01);
-        assert(std::abs((var - x_var) / x_var) < 0.01);
-        assert(std::abs(skew - x_skew) < 0.01);
-        assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
+    assert(bottom_tenth);  // ...is populated
+    assert(top_tenth);  // ...is populated
+
+    // Now do some more involved statistical math.
+    double mean = std::accumulate(u.begin(), u.end(), 0.0) / u.size();
+    double var = 0;
+    double skew = 0;
+    double kurtosis = 0;
+    for (std::size_t i = 0; i < u.size(); ++i) {
+        double dbl = (u[i] - mean);
+        double d2 = dbl * dbl;
+        var += d2;
+        skew += dbl * d2;
+        kurtosis += d2 * d2;
     }
+    var /= u.size();
+    double dev = std::sqrt(var);
+    skew /= u.size() * dev * var;
+    kurtosis /= u.size() * var * var;
+
+    double expected_mean = double(a) + double(b)/2 - double(a)/2;
+    double expected_var = (sqr(double(b) - double(a) + 1) - 1) / 12;
+
+    double range = double(b) - double(a) + 1.0;
+    assert(range > range / 10);  // i.e., it's not infinity
+
+    assert(std::abs(mean - expected_mean) < range / 100);
+    assert(std::abs(var - expected_var) < expected_var / 50);
+    assert(-0.1 < skew && skew < 0.1);
+    assert(1.6 < kurtosis && kurtosis < 2.0);
+}
+
+template <class ResultType, class EngineType>
+void test_statistics()
+{
+    test_statistics<ResultType, EngineType>(0, std::numeric_limits<ResultType>::max());
+}
+
+int main(int, char**)
+{
+    test_statistics<int, std::minstd_rand0>();
+    test_statistics<int, std::minstd_rand>();
+    test_statistics<int, std::mt19937>();
+    test_statistics<int, std::mt19937_64>();
+    test_statistics<int, std::ranlux24_base>();
+    test_statistics<int, std::ranlux48_base>();
+    test_statistics<int, std::ranlux24>();
+    test_statistics<int, std::ranlux48>();
+    test_statistics<int, std::knuth_b>();
+    test_statistics<int, std::minstd_rand0>(-6, 106);
+    test_statistics<int, std::minstd_rand>(5, 100);
+
+    test_statistics<short, std::minstd_rand0>();
+    test_statistics<int, std::minstd_rand0>();
+    test_statistics<long, std::minstd_rand0>();
+    test_statistics<long long, std::minstd_rand0>();
+
+    test_statistics<unsigned short, std::minstd_rand0>();
+    test_statistics<unsigned int, std::minstd_rand0>();
+    test_statistics<unsigned long, std::minstd_rand0>();
+    test_statistics<unsigned long long, std::minstd_rand0>();
+
+    test_statistics<short, std::minstd_rand0>(SHRT_MIN, SHRT_MAX);
+
+    // http://eel.is/c++draft/rand.req#genl-1.5
+    // The effect of instantiating a template that has a parameter
+    // named IntType is undefined unless the corresponding template
+    // argument is cv-unqualified and is one of short, int, long,
+    // long long, unsigned short, unsigned int, unsigned long,
+    // or unsigned long long.
+    // (We support __int128 as an extension.)
+
+#ifndef _LIBCPP_HAS_NO_INT128
+    test_statistics<__int128_t, std::minstd_rand0>();
+    test_statistics<__uint128_t, std::minstd_rand0>();
+
+    test_statistics<__int128_t, std::minstd_rand0>(-100, 900);
+    test_statistics<__int128_t, std::minstd_rand0>(0, UINT64_MAX);
+    test_statistics<__int128_t, std::minstd_rand0>(std::numeric_limits<__int128_t>::min(), std::numeric_limits<__int128_t>::max());
+    test_statistics<__uint128_t, std::minstd_rand0>(0, UINT64_MAX);
+#endif
 
-  return 0;
+    return 0;
 }
-- 
2.7.4