diff --git a/benches/average.rs b/benches/average.rs
new file mode 100644
index 0000000..05d824c
--- /dev/null
+++ b/benches/average.rs
@@ -0,0 +1,414 @@
+//! Benchmark sqrt and cbrt
+
+#![feature(test)]
+
+extern crate num_integer;
+extern crate num_traits;
+extern crate test;
+
+use num_integer::Integer;
+use num_traits::{AsPrimitive, PrimInt, WrappingAdd, WrappingMul};
+use std::cmp::{max, min};
+use std::fmt::Debug;
+use test::{black_box, Bencher};
+
+// --- Utilities for RNG ----------------------------------------------------
+
+trait BenchInteger: Integer + PrimInt + WrappingAdd + WrappingMul + 'static {}
+
+impl<T> BenchInteger for T where T: Integer + PrimInt + WrappingAdd + WrappingMul + 'static {}
+
+// 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_())
+}
+
+// --- Alt. Implementations -------------------------------------------------
+
+trait NaiveAverage {
+    fn naive_average_ceil(&self, other: &Self) -> Self;
+    fn naive_average_floor(&self, other: &Self) -> Self;
+}
+
+trait UncheckedAverage {
+    fn unchecked_average_ceil(&self, other: &Self) -> Self;
+    fn unchecked_average_floor(&self, other: &Self) -> Self;
+}
+
+trait ModuloAverage {
+    fn modulo_average_ceil(&self, other: &Self) -> Self;
+    fn modulo_average_floor(&self, other: &Self) -> Self;
+}
+
+macro_rules! naive_average {
+    ($T:ident) => {
+        impl super::NaiveAverage for $T {
+            fn naive_average_floor(&self, other: &$T) -> $T {
+                match self.checked_add(*other) {
+                    Some(z) => z.div_floor(&2),
+                    None => {
+                        if self > other {
+                            let diff = self - other;
+                            other + diff.div_floor(&2)
+                        } else {
+                            let diff = other - self;
+                            self + diff.div_floor(&2)
+                        }
+                    }
+                }
+            }
+            fn naive_average_ceil(&self, other: &$T) -> $T {
+                match self.checked_add(*other) {
+                    Some(z) => z.div_ceil(&2),
+                    None => {
+                        if self > other {
+                            let diff = self - other;
+                            self - diff.div_floor(&2)
+                        } else {
+                            let diff = other - self;
+                            other - diff.div_floor(&2)
+                        }
+                    }
+                }
+            }
+        }
+    };
+}
+
+macro_rules! unchecked_average {
+    ($T:ident) => {
+        impl super::UncheckedAverage for $T {
+            fn unchecked_average_floor(&self, other: &$T) -> $T {
+                self.wrapping_add(*other) / 2
+            }
+            fn unchecked_average_ceil(&self, other: &$T) -> $T {
+                (self.wrapping_add(*other) / 2).wrapping_add(1)
+            }
+        }
+    };
+}
+
+macro_rules! modulo_average {
+    ($T:ident) => {
+        impl super::ModuloAverage for $T {
+            fn modulo_average_ceil(&self, other: &$T) -> $T {
+                let (q1, r1) = self.div_mod_floor(&2);
+                let (q2, r2) = other.div_mod_floor(&2);
+                q1 + q2 + (r1 | r2)
+            }
+            fn modulo_average_floor(&self, other: &$T) -> $T {
+                let (q1, r1) = self.div_mod_floor(&2);
+                let (q2, r2) = other.div_mod_floor(&2);
+                q1 + q2 + (r1 * r2)
+            }
+        }
+    };
+}
+
+// --- Bench functions ------------------------------------------------------
+
+fn bench_unchecked<T, F>(b: &mut Bencher, v: &[(T, T)], f: F)
+where
+    T: Integer + Debug + Copy,
+    F: Fn(&T, &T) -> T,
+{
+    b.iter(|| {
+        for (x, y) in v {
+            black_box(f(x, y));
+        }
+    });
+}
+
+fn bench_ceil<T, F>(b: &mut Bencher, v: &[(T, T)], f: F)
+where
+    T: Integer + Debug + Copy,
+    F: Fn(&T, &T) -> T,
+{
+    for &(i, j) in v {
+        let rt = f(&i, &j);
+        let (a, b) = (min(i, j), max(i, j));
+        // if both number are the same sign, check rt is in the middle
+        if (a < T::zero()) == (b < T::zero()) {
+            if (b - a).is_even() {
+                assert_eq!(rt - a, b - rt);
+            } else {
+                assert_eq!(rt - a, b - rt + T::one());
+            }
+        // if both number have a different sign,
+        } else {
+            if (a + b).is_even() {
+                assert_eq!(rt, (a + b) / (T::one() + T::one()))
+            } else {
+                assert_eq!(rt, (a + b + T::one()) / (T::one() + T::one()))
+            }
+        }
+    }
+    bench_unchecked(b, v, f);
+}
+
+fn bench_floor<T, F>(b: &mut Bencher, v: &[(T, T)], f: F)
+where
+    T: Integer + Debug + Copy,
+    F: Fn(&T, &T) -> T,
+{
+    for &(i, j) in v {
+        let rt = f(&i, &j);
+        let (a, b) = (min(i, j), max(i, j));
+        // if both number are the same sign, check rt is in the middle
+        if (a < T::zero()) == (b < T::zero()) {
+            if (b - a).is_even() {
+                assert_eq!(rt - a, b - rt);
+            } else {
+                assert_eq!(rt - a + T::one(), b - rt);
+            }
+        // if both number have a different sign,
+        } else {
+            if (a + b).is_even() {
+                assert_eq!(rt, (a + b) / (T::one() + T::one()))
+            } else {
+                assert_eq!(rt, (a + b - T::one()) / (T::one() + T::one()))
+            }
+        }
+    }
+    bench_unchecked(b, v, f);
+}
+
+// --- Bench implementation -------------------------------------------------
+
+macro_rules! bench_average {
+    ($($T:ident),*) => {$(
+        mod $T {
+            use test::Bencher;
+            use num_integer::{Average, Integer};
+            use super::{UncheckedAverage, NaiveAverage, ModuloAverage};
+            use super::{bench_ceil, bench_floor, bench_unchecked};
+
+            naive_average!($T);
+            unchecked_average!($T);
+            modulo_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()
+            }
+
+            mod ceil {
+
+                use super::*;
+
+                mod small {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = small();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = small();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.naive_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = small();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = small();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.modulo_average_ceil(y));
+                    }
+                }
+
+                mod overflowing {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.naive_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.modulo_average_ceil(y));
+                    }
+                }
+
+                mod rand {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = rand();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = rand();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.naive_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = rand();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_ceil(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = rand();
+                        bench_ceil(b, &v, |x: &$T, y: &$T| x.modulo_average_ceil(y));
+                    }
+                }
+
+            }
+
+            mod floor {
+
+                use super::*;
+
+                mod small {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = small();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = small();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = small();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = small();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.modulo_average_floor(y));
+                    }
+                }
+
+                mod overflowing {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = overflowing();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.modulo_average_floor(y));
+                    }
+                }
+
+                mod rand {
+
+                    use super::*;
+
+                    #[bench]
+                    fn optimized(b: &mut Bencher) {
+                        let v = rand();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.average_floor(y));
+                    }
+
+                    #[bench]
+                    fn naive(b: &mut Bencher) {
+                        let v = rand();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.naive_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn unchecked(b: &mut Bencher) {
+                        let v = rand();
+                        bench_unchecked(b, &v, |x: &$T, y: &$T| x.unchecked_average_floor(y));
+                    }
+
+                    #[bench]
+                    fn modulo(b: &mut Bencher) {
+                        let v = rand();
+                        bench_floor(b, &v, |x: &$T, y: &$T| x.modulo_average_floor(y));
+                    }
+                }
+
+            }
+
+        }
+    )*}
+}
+
+bench_average!(i8, i16, i32, i64, i128, isize);
+bench_average!(u8, u16, u32, u64, u128, usize);
diff --git a/src/average.rs b/src/average.rs
new file mode 100644
index 0000000..29cd11e
--- /dev/null
+++ b/src/average.rs
@@ -0,0 +1,78 @@
+use core::ops::{BitAnd, BitOr, BitXor, Shr};
+use Integer;
+
+/// Provides methods to compute the average of two integers, without overflows.
+pub trait Average: Integer {
+    /// Returns the ceiling value of the average of `self` and `other`.
+    /// -- `⌈(self + other)/2⌉`
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// 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);
+    /// assert_eq!(i8::min_value().average_ceil(&-1), -64);
+    /// assert_eq!(i8::min_value().average_ceil(&i8::max_value()), 0);
+    /// ```
+    ///
+    fn average_ceil(&self, other: &Self) -> Self;
+
+    /// Returns the floor value of the average of `self` and `other`.
+    /// -- `⌊(self + other)/2⌋`
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// 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);
+    /// assert_eq!(i8::min_value().average_floor(&-1), -65);
+    /// assert_eq!(i8::min_value().average_floor(&i8::max_value()), -1);
+    /// ```
+    ///
+    fn average_floor(&self, other: &Self) -> Self;
+}
+
+impl<I> Average for I
+where
+    I: Integer + 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 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)
+}
diff --git a/src/lib.rs b/src/lib.rs
index bda2981..0281954 100644
--- a/src/lib.rs
+++ b/src/lib.rs
@@ -30,6 +30,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.
     ///
diff --git a/tests/average.rs b/tests/average.rs
new file mode 100644
index 0000000..9fd8cf1
--- /dev/null
+++ b/tests/average.rs
@@ -0,0 +1,100 @@
+extern crate num_integer;
+extern crate num_traits;
+
+macro_rules! test_average {
+    ($I:ident, $U:ident) => {
+        mod $I {
+            mod ceil {
+                use num_integer::Average;
+
+                #[test]
+                fn same_sign() {
+                    assert_eq!((14 as $I).average_ceil(&16), 15 as $I);
+                    assert_eq!((14 as $I).average_ceil(&17), 16 as $I);
+
+                    let max = $crate::std::$I::MAX;
+                    assert_eq!((max - 3).average_ceil(&(max - 1)), max - 2);
+                    assert_eq!((max - 3).average_ceil(&(max - 2)), max - 2);
+                }
+
+                #[test]
+                fn different_sign() {
+                    assert_eq!((14 as $I).average_ceil(&-4), 5 as $I);
+                    assert_eq!((14 as $I).average_ceil(&-5), 5 as $I);
+
+                    let min = $crate::std::$I::MIN;
+                    let max = $crate::std::$I::MAX;
+                    assert_eq!(min.average_ceil(&max), 0 as $I);
+                }
+            }
+
+            mod floor {
+                use num_integer::Average;
+
+                #[test]
+                fn same_sign() {
+                    assert_eq!((14 as $I).average_floor(&16), 15 as $I);
+                    assert_eq!((14 as $I).average_floor(&17), 15 as $I);
+
+                    let max = $crate::std::$I::MAX;
+                    assert_eq!((max - 3).average_floor(&(max - 1)), max - 2);
+                    assert_eq!((max - 3).average_floor(&(max - 2)), max - 3);
+                }
+
+                #[test]
+                fn different_sign() {
+                    assert_eq!((14 as $I).average_floor(&-4), 5 as $I);
+                    assert_eq!((14 as $I).average_floor(&-5), 4 as $I);
+
+                    let min = $crate::std::$I::MIN;
+                    let max = $crate::std::$I::MAX;
+                    assert_eq!(min.average_floor(&max), -1 as $I);
+                }
+            }
+        }
+
+        mod $U {
+            mod ceil {
+                use num_integer::Average;
+
+                #[test]
+                fn bounded() {
+                    assert_eq!((14 as $U).average_ceil(&16), 15 as $U);
+                    assert_eq!((14 as $U).average_ceil(&17), 16 as $U);
+                }
+
+                #[test]
+                fn overflow() {
+                    let max = $crate::std::$U::MAX;
+                    assert_eq!((max - 3).average_ceil(&(max - 1)), max - 2);
+                    assert_eq!((max - 3).average_ceil(&(max - 2)), max - 2);
+                }
+            }
+
+            mod floor {
+                use num_integer::Average;
+
+                #[test]
+                fn bounded() {
+                    assert_eq!((14 as $U).average_floor(&16), 15 as $U);
+                    assert_eq!((14 as $U).average_floor(&17), 15 as $U);
+                }
+
+                #[test]
+                fn overflow() {
+                    let max = $crate::std::$U::MAX;
+                    assert_eq!((max - 3).average_floor(&(max - 1)), max - 2);
+                    assert_eq!((max - 3).average_floor(&(max - 2)), max - 3);
+                }
+            }
+        }
+    };
+}
+
+test_average!(i8, u8);
+test_average!(i16, u16);
+test_average!(i32, u32);
+test_average!(i64, u64);
+#[cfg(has_i128)]
+test_average!(i128, u128);
+test_average!(isize, usize);