MAINT: Restore legacy zipf implementation.

numpy · charris · Jul 26, 2024 · Jul 26, 2024 · Jul 26, 2024 · Jul 26, 2024
commit 0a576b5c93e596a10db676ddabb7069bf637b74a
diff --git a/numpy/random/src/legacy/legacy-distributions.c b/numpy/random/src/legacy/legacy-distributions.c
@@ -388,7 +388,31 @@ int64_t legacy_random_poisson(bitgen_t *bitgen_state, double lam) {
 }
 
 int64_t legacy_random_zipf(bitgen_t *bitgen_state, double a) {
-  return (int64_t)random_zipf(bitgen_state, a);
+  double am1, b;
+
+  am1 = a - 1.0;
+  b = pow(2.0, am1);
+  while (1) {
+    double T, U, V, X;
+
+    U = 1.0 - next_double(bitgen_state);
+    V = next_double(bitgen_state);
+    X = floor(pow(U, -1.0 / am1));
+    /*
+     * The real result may be above what can be represented in a signed
+     * long. Since this is a straightforward rejection algorithm, we can
+     * just reject this value. This function then models a Zipf
+     * distribution truncated to sys.maxint.
+     */
+    if (X > (double)RAND_INT_MAX || X < 1.0) {
+      continue;
+    }
+
+    T = pow(1.0 + 1.0 / X, am1);
+    if (V * X * (T - 1.0) / (b - 1.0) <= T / b) {
+      return (RAND_INT_TYPE)X;
+    }
+  }
 }