From b9b5cc5ccc4d027a9b44be7f02b64d1d00624457 Mon Sep 17 00:00:00 2001
From: Martin Larralde <martin.larralde@ens-cachan.fr>
Date: Wed, 19 Dec 2018 20:43:41 +0100
Subject: [PATCH 1/4] Add `Average` trait with fast average computation

---
 src/average.rs | 77 ++++++++++++++++++++++++++++++++++++++++++++++++++
 src/lib.rs     |  4 +++
 2 files changed, 81 insertions(+)
 create mode 100644 src/average.rs

diff --git a/src/average.rs b/src/average.rs
new file mode 100644
index 0000000..42cba41
--- /dev/null
+++ b/src/average.rs
@@ -0,0 +1,77 @@
+use core::ops::{Add, BitAnd, BitOr, BitXor, Shr, Sub};
+use Integer;
+
+pub trait Average: Integer {
+    /// Returns the floor value of the average of `self` and `other`.
+    /// -- `⌊(self + other)/2⌋
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// # extern crate num_traits;
+    /// # use num_traits::bounds::Bounded;
+    /// use num_integer::Average;
+    ///
+    /// assert_eq!(( 3).average_ceil(&10),  7);
+    /// assert_eq!((-2).average_ceil(&-5), -3);
+    /// assert_eq!(( 4).average_ceil(& 4),  4);
+    ///
+    /// assert_eq!(u8::max_value().average_ceil(&2), 129);
+    /// ```
+    ///
+    fn average_ceil(&self, other: &Self) -> Self;
+
+    /// Returns the ceil value of the average of `self` and `other`.
+    /// -- `⌈(self + other)/2⌉
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// # extern crate num_traits;
+    /// # use num_traits::bounds::Bounded;
+    /// use num_integer::Average;
+    ///
+    /// assert_eq!(( 3).average_floor(&10),  6);
+    /// assert_eq!((-2).average_floor(&-5), -4);
+    /// assert_eq!(( 4).average_floor(& 4),  4);
+    ///
+    /// assert_eq!(u8::max_value().average_floor(&2), 128);
+    /// ```
+    ///
+    fn average_floor(&self, other: &Self) -> Self;
+}
+
+impl<I> Average for I
+where
+    I: Integer + Add<I, Output = I> + Sub<I, Output = I> + Shr<usize, Output = I>,
+    for<'a, 'b> &'a I:
+        BitAnd<&'b I, Output = I> + BitOr<&'b I, Output = I> + BitXor<&'b I, Output = I>,
+{
+    // The Henry Gordon Dietz, implementation as shown in the Hacker's Delight,
+    // see http://aggregate.org/MAGIC/#Average%20of%20Integers
+
+    /// Returns the floor value of the average of `self` and `other`.
+    #[inline]
+    fn average_floor(&self, other: &I) -> I {
+        (self & other) + ((self ^ other) >> 1)
+    }
+
+    /// Returns the ceil value of the average of `self` and `other`.
+    #[inline]
+    fn average_ceil(&self, other: &I) -> I {
+        (self | other) - ((self ^ other) >> 1)
+    }
+}
+
+/// Returns the floor value of the average of `x` and `y` --
+/// see [Average::average_floor](trait.Average.html#tymethod.average_floor).
+#[inline]
+pub fn average_floor<T: Average>(x: &T, y: &T) -> T {
+    x.average_floor(y)
+}
+/// Returns the ceil value of the average of `x` and `y` --
+/// see [Average::average_floor](trait.Average.html#tymethod.average_floor).
+#[inline]
+pub fn average_ceil<T: Average>(x: &T, y: &T) -> T {
+    x.average_ceil(y)
+}
diff --git a/src/lib.rs b/src/lib.rs
index 5814f65..a938029 100644
--- a/src/lib.rs
+++ b/src/lib.rs
@@ -29,6 +29,10 @@ mod roots;
 pub use roots::Roots;
 pub use roots::{cbrt, nth_root, sqrt};
 
+mod average;
+pub use average::Average;
+pub use average::{average_ceil, average_floor};
+
 pub trait Integer: Sized + Num + PartialOrd + Ord + Eq {
     /// Floored integer division.
     ///

From 3bb1ba03d9593b841593748516c3fd5c03ee419c Mon Sep 17 00:00:00 2001
From: Martin Larralde <martin.larralde@ens-cachan.fr>
Date: Wed, 19 Dec 2018 20:45:17 +0100
Subject: [PATCH 2/4] Add benchmarks using two references (safe and unsafe)

---
 benches/average.rs | 185 +++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 185 insertions(+)
 create mode 100644 benches/average.rs

diff --git a/benches/average.rs b/benches/average.rs
new file mode 100644
index 0000000..4b20dd2
--- /dev/null
+++ b/benches/average.rs
@@ -0,0 +1,185 @@
+//! Benchmark sqrt and cbrt
+
+#![feature(test)]
+
+extern crate num_integer;
+extern crate num_traits;
+extern crate test;
+
+use num_integer::{Average, Integer};
+use num_traits::checked_pow;
+use num_traits::{AsPrimitive, PrimInt, WrappingAdd, WrappingMul};
+use test::{black_box, Bencher};
+
+trait BenchInteger: Integer + PrimInt + WrappingAdd + WrappingMul + 'static {}
+
+impl<T> BenchInteger for T where T: Integer + PrimInt + WrappingAdd + WrappingMul + 'static {}
+
+trait NaiveAverage {
+    fn naive_average_floor(&self, other: &Self) -> Self;
+    fn naive_average_ceil(&self, other: &Self) -> Self;
+}
+
+trait UncheckedAverage {
+    fn unchecked_average_floor(&self, other: &Self) -> Self;
+    fn unchecked_average_ceil(&self, other: &Self) -> Self;
+}
+
+fn bench<T, F>(b: &mut Bencher, v: &[(T, T)], f: F)
+where
+    T: Integer,
+    F: Fn(&T, &T) -> T,
+{
+    b.iter(|| {
+        for (x, y) in v {
+            black_box(f(x, y));
+        }
+    });
+}
+
+// Simple PRNG so we don't have to worry about rand compatibility
+fn lcg<T>(x: T) -> T
+where
+    u32: AsPrimitive<T>,
+    T: BenchInteger,
+{
+    // LCG parameters from Numerical Recipes
+    // (but we're applying it to arbitrary sizes)
+    const LCG_A: u32 = 1664525;
+    const LCG_C: u32 = 1013904223;
+    x.wrapping_mul(&LCG_A.as_()).wrapping_add(&LCG_C.as_())
+}
+
+macro_rules! naive_average {
+    ($T:ident) => {
+        impl super::NaiveAverage for $T {
+            fn naive_average_floor(&self, other: &$T) -> $T {
+                let (x, y) = if self > other {
+                    (self, other)
+                } else {
+                    (other, self)
+                };
+                let diff = x - y;
+                y + (diff / 2)
+            }
+            fn naive_average_ceil(&self, other: &$T) -> $T {
+                let (x, y) = if self > other {
+                    (self, other)
+                } else {
+                    (other, self)
+                };
+                let diff = x - y;
+                x - (diff / 2)
+            }
+        }
+    };
+}
+
+macro_rules! unchecked_average {
+    ($T:ident) => {
+        impl super::UncheckedAverage for $T {
+            fn unchecked_average_floor(&self, other: &$T) -> $T {
+                self.saturating_add(*other) / 2
+            }
+            fn unchecked_average_ceil(&self, other: &$T) -> $T {
+                self.saturating_add(*other) / 2
+            }
+        }
+    };
+}
+
+macro_rules! bench_average {
+    ($($T:ident),*) => {$(
+        mod $T {
+            use test::Bencher;
+            use num_integer::Average;
+            use UncheckedAverage;
+            use NaiveAverage;
+
+            naive_average!($T);
+            unchecked_average!($T);
+
+            const SIZE: $T = 30;
+
+            fn overflowing() -> Vec<($T, $T)> {
+                (($T::max_value()-SIZE)..$T::max_value())
+                    .flat_map(|x| -> Vec<_> {
+                        (($T::max_value()-100)..($T::max_value()-100+SIZE))
+                            .map(|y| (x, y))
+                            .collect()
+                    })
+                    .collect()
+            }
+
+            fn small() -> Vec<($T, $T)> {
+                (0..SIZE)
+                   .flat_map(|x| -> Vec<_> {(0..SIZE).map(|y| (x, y)).collect()})
+                   .collect()
+            }
+
+            fn rand() -> Vec<($T, $T)> {
+                small()
+                    .into_iter()
+                    .map(|(x, y)| (super::lcg(x), super::lcg(y)))
+                    .collect()
+            }
+
+            #[bench]
+            fn average_floor_small(b: &mut Bencher) {
+                let v = small();
+                super::bench(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_small_naive(b: &mut Bencher) {
+                let v = small();
+                super::bench(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_small_unchecked(b: &mut Bencher) {
+                let v = small();
+                super::bench(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_overflowing(b: &mut Bencher) {
+                let v = overflowing();
+                super::bench(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_overflowing_naive(b: &mut Bencher) {
+                let v = overflowing();
+                super::bench(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_overflowing_unchecked(b: &mut Bencher) {
+                let v = overflowing();
+                super::bench(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_rand(b: &mut Bencher) {
+                let v = rand();
+                super::bench(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_rand_naive(b: &mut Bencher) {
+                let v = rand();
+                super::bench(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+            }
+
+            #[bench]
+            fn average_floor_rand_unchecked(b: &mut Bencher) {
+                let v = rand();
+                super::bench(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+            }
+        }
+    )*}
+}
+
+bench_average!(i8, i16, i32, i64, i128, isize);
+bench_average!(u8, u16, u32, u64, u128, usize);

From b470de6504b965cdc6f9a0dcddd2a4e2f64ace03 Mon Sep 17 00:00:00 2001
From: Martin Larralde <martin.larralde@ens-cachan.fr>
Date: Wed, 19 Dec 2018 21:20:34 +0100
Subject: [PATCH 3/4] Fix doctests in `src/average.rs`

---
 benches/average.rs |  3 +--
 src/average.rs     | 13 +++++--------
 2 files changed, 6 insertions(+), 10 deletions(-)

diff --git a/benches/average.rs b/benches/average.rs
index 4b20dd2..480d669 100644
--- a/benches/average.rs
+++ b/benches/average.rs
@@ -6,8 +6,7 @@ extern crate num_integer;
 extern crate num_traits;
 extern crate test;
 
-use num_integer::{Average, Integer};
-use num_traits::checked_pow;
+use num_integer::Integer;
 use num_traits::{AsPrimitive, PrimInt, WrappingAdd, WrappingMul};
 use test::{black_box, Bencher};
 
diff --git a/src/average.rs b/src/average.rs
index 42cba41..eca7ba2 100644
--- a/src/average.rs
+++ b/src/average.rs
@@ -1,15 +1,14 @@
 use core::ops::{Add, BitAnd, BitOr, BitXor, Shr, Sub};
 use Integer;
 
+/// Provides methods to compute the average of two integers, without overflows.
 pub trait Average: Integer {
-    /// Returns the floor value of the average of `self` and `other`.
-    /// -- `⌊(self + other)/2⌋
+    /// Returns the ceiling value of the average of `self` and `other`.
+    /// -- `⌈(self + other)/2⌉`
     ///
     /// # Examples
     ///
     /// ```
-    /// # extern crate num_traits;
-    /// # use num_traits::bounds::Bounded;
     /// use num_integer::Average;
     ///
     /// assert_eq!(( 3).average_ceil(&10),  7);
@@ -21,14 +20,12 @@ pub trait Average: Integer {
     ///
     fn average_ceil(&self, other: &Self) -> Self;
 
-    /// Returns the ceil value of the average of `self` and `other`.
-    /// -- `⌈(self + other)/2⌉
+    /// Returns the floor value of the average of `self` and `other`.
+    /// -- `⌊(self + other)/2⌋`
     ///
     /// # Examples
     ///
     /// ```
-    /// # extern crate num_traits;
-    /// # use num_traits::bounds::Bounded;
     /// use num_integer::Average;
     ///
     /// assert_eq!(( 3).average_floor(&10),  6);

From 93c0f7664ada459a2f1ef5937323c59dd270ca18 Mon Sep 17 00:00:00 2001
From: Martin Larralde <martin.larralde@ens-cachan.fr>
Date: Wed, 19 Dec 2018 22:24:52 +0100
Subject: [PATCH 4/4] Fix naive average implementation in tests

---
 benches/average.rs | 34 ++++++++++++++++++++--------------
 src/average.rs     |  6 +++---
 2 files changed, 23 insertions(+), 17 deletions(-)

diff --git a/benches/average.rs b/benches/average.rs
index 480d669..a53aedc 100644
--- a/benches/average.rs
+++ b/benches/average.rs
@@ -53,22 +53,28 @@ macro_rules! naive_average {
     ($T:ident) => {
         impl super::NaiveAverage for $T {
             fn naive_average_floor(&self, other: &$T) -> $T {
-                let (x, y) = if self > other {
-                    (self, other)
-                } else {
-                    (other, self)
-                };
-                let diff = x - y;
-                y + (diff / 2)
+                match self.checked_add(*other) {
+                    Some(z) => z/2,
+                    None => if self > other {
+                        let diff = self - other;
+                        other + diff/2
+                    } else {
+                        let diff = other - self;
+                        self + diff/2
+                    }
+                }
             }
             fn naive_average_ceil(&self, other: &$T) -> $T {
-                let (x, y) = if self > other {
-                    (self, other)
-                } else {
-                    (other, self)
-                };
-                let diff = x - y;
-                x - (diff / 2)
+                match self.checked_add(*other).and_then(|x| x.checked_add(1)) {
+                    Some(z) => z/2,
+                    None => if self > other {
+                        let diff = self - other;
+                        self - diff/2
+                    } else {
+                        let diff = other - self;
+                        other - diff/2
+                    }
+                }
             }
         }
     };
diff --git a/src/average.rs b/src/average.rs
index eca7ba2..1f95fa5 100644
--- a/src/average.rs
+++ b/src/average.rs
@@ -44,7 +44,7 @@ where
     for<'a, 'b> &'a I:
         BitAnd<&'b I, Output = I> + BitOr<&'b I, Output = I> + BitXor<&'b I, Output = I>,
 {
-    // The Henry Gordon Dietz, implementation as shown in the Hacker's Delight,
+    // The Henry Gordon Dietz implementation as shown in the Hacker's Delight,
     // see http://aggregate.org/MAGIC/#Average%20of%20Integers
 
     /// Returns the floor value of the average of `self` and `other`.
@@ -66,8 +66,8 @@ where
 pub fn average_floor<T: Average>(x: &T, y: &T) -> T {
     x.average_floor(y)
 }
-/// Returns the ceil value of the average of `x` and `y` --
-/// see [Average::average_floor](trait.Average.html#tymethod.average_floor).
+/// Returns the ceiling value of the average of `x` and `y` --
+/// see [Average::average_ceil](trait.Average.html#tymethod.average_ceil).
 #[inline]
 pub fn average_ceil<T: Average>(x: &T, y: &T) -> T {
     x.average_ceil(y)