Avoid hangs and debug asserts on invalid parameters for Zipf

The generalized Zipf distribution is not well defined when its
normalization constant is infinity (when n is inf and s is <= 1).

Also, when s is inf, there is a reasonable limiting distribution,
which can be produced with a small special case in Zipf::new.

Author: mstoeckl

CHANGELOG.md

@@ -19,6 +19,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Fixes
 - Fix `Geometric::new` for small `p > 0` where `1 - p` rounds to 1 (#36)
 - Fix panic in `FisherF::new` on almost zero parameters (#39)
+- Fix hang and debug assertion in `Zipf::new` on invalid parameters (#41)
 
 ## [0.5.2]
 

src/zipf.rs

@@ -66,17 +66,20 @@ where
 /// Error type returned from [`Zipf::new`].
 #[derive(Clone, Copy, Debug, PartialEq, Eq)]
 pub enum Error {
-    /// `s < 0` or `nan`.
+    /// `s < 0` or `s` is `nan`
     STooSmall,
-    /// `n < 1`.
+    /// `n < 1` or `n` is `nan`
     NTooSmall,
+    /// `n = inf` and `s <= 1`
+    IllDefined,
 }
 
 impl fmt::Display for Error {
     fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
         f.write_str(match self {
             Error::STooSmall => "s < 0 or is NaN in Zipf distribution",
-            Error::NTooSmall => "n < 1 in Zipf distribution",
+            Error::NTooSmall => "n < 1 or is NaN in Zipf distribution",
+            Error::IllDefined => "n = inf and s <= 1 in Zipf distribution",
         })
     }
 }
@@ -103,19 +106,25 @@ where
         if n < F::one() {
             return Err(Error::NTooSmall);
         }
+        if n.is_infinite() && s <= F::one() {
+            return Err(Error::IllDefined);
+        }
         let q = if s != F::one() {
             // Make sure to calculate the division only once.
             F::one() / (F::one() - s)
         } else {
             // This value is never used.
             F::zero()
         };
-        let t = if s != F::one() {
-            (n.powf(F::one() - s) - s) * q
+        let t = if s == F::infinity() {
+            F::zero()
+        } else if s != F::one() {
+            let t = (n.powf(F::one() - s) - s) * q;
+            debug_assert!(t > F::zero());
+            t
         } else {
             F::one() + n.ln()
         };
-        debug_assert!(t > F::zero());
         Ok(Zipf { s, t, q })
     }
 
@@ -220,6 +229,16 @@ mod tests {
         // TODO: verify that this is a uniform distribution
     }
 
+    #[test]
+    fn zipf_sample_s_inf() {
+        let d = Zipf::new(10., f64::infinity()).unwrap();
+        let mut rng = crate::test::rng(2);
+        for _ in 0..1000 {
+            let r = d.sample(&mut rng);
+            assert!(r == 1.);
+        }
+    }
+
     #[test]
     fn zipf_sample_large_n() {
         let d = Zipf::new(f64::MAX, 1.5).unwrap();