1#![allow(dead_code)]
9
10#[cfg(test)]
11#[macro_use]
12extern crate approx;
13extern crate num_bigint;
14
15use num_bigint::BigUint;
16use num_traits::cast::ToPrimitive;
17use num_traits::{One, Zero};
18use std::convert::TryFrom;
19use std::f64;
20
21pub fn bell(n: usize) -> BigUint {
35 let mut r1: Vec<BigUint> = vec![Zero::zero(); n];
36 let mut r2: Vec<BigUint> = vec![Zero::zero(); n];
37 r1[0] = One::one();
38 for k in 1..n {
39 r2[0] = r1[k - 1].clone();
40 for i in 1..(k + 1) {
41 r2[i] = r1[i - 1].clone() + &r2[i - 1];
42 }
43 let tmp = r1;
44 r1 = r2;
45 r2 = tmp;
46 }
47 r1[n - 1].clone()
48}
49
50pub fn lbell(n: usize) -> f64 {
60 let value = bell(n);
61 let n_bits = value.bits();
62 let threshold = 1022_u64;
63 let log2 = if n_bits > threshold {
64 let n_shifted_bits = value.bits() - threshold;
65 let shifted_value = value >> n_shifted_bits;
66 if shifted_value.bits() > threshold {
67 return f64::INFINITY;
68 }
69 let y: f64 = shifted_value.to_f64().unwrap();
70 (n_shifted_bits as f64) + y.log2()
71 } else {
72 value.to_f64().unwrap().log2()
73 };
74 log2 / f64::consts::LOG2_E
75}
76
77#[doc(hidden)]
80#[no_mangle]
81pub extern "C" fn dahl_bellnumber__bell(n: i32) -> f64 {
82 if n < 0 {
83 return 0.0;
84 }
85 match usize::try_from(n) {
86 Ok(n) => match bell(n).to_f64() {
87 Some(x) => x,
88 None => f64::INFINITY,
89 },
90 Err(_) => f64::INFINITY,
91 }
92}
93
94#[doc(hidden)]
97#[no_mangle]
98pub extern "C" fn dahl_bellnumber__lbell(n: i32) -> f64 {
99 if n < 0 {
100 return 0.0;
101 }
102 match usize::try_from(n) {
103 Ok(n) => lbell(n),
104 Err(_) => f64::INFINITY,
105 }
106}
107
108#[cfg(test)]
109mod tests {
110 use super::*;
111
112 #[test]
113 fn test_lbell() {
114 assert_relative_eq!(lbell(220), 714.4032630589774);
115 assert_relative_eq!(bell(5).to_f64().unwrap(), 52.0);
116 }
117
118}