diff --git a/LICENSE b/LICENSE
new file mode 100644
--- /dev/null
+++ b/LICENSE
@@ -0,0 +1,27 @@
+Copyright 2008, Ross Paterson.
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are met:
+
+- Redistributions of source code must retain the above copyright notice,
+this list of conditions and the following disclaimer.
+ 
+- Redistributions in binary form must reproduce the above copyright notice,
+this list of conditions and the following disclaimer in the documentation
+and/or other materials provided with the distribution.
+ 
+- The names of the contributors may not be used to endorse or promote
+products derived from this software without specific prior written
+permission.
+
+THIS SOFTWARE IS PROVIDED "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES,
+INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY
+AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
+THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
diff --git a/Numeric/Search/Bounded.hs b/Numeric/Search/Bounded.hs
new file mode 100644
--- /dev/null
+++ b/Numeric/Search/Bounded.hs
@@ -0,0 +1,74 @@
+-----------------------------------------------------------------------------
+-- |
+-- Module      :  Numeric.Search.Bounded
+-- Copyright   :  (c) Ross Paterson 2008
+-- License     :  BSD-style
+-- Maintainer  :  ross@soi.city.ac.uk
+-- Stability   :  experimental
+-- Portability :  portable
+--
+-- Searching unbounded intervals within bounded integral types for the
+-- boundary of an upward-closed set, using a combination of exponential
+-- and binary search.
+--
+-----------------------------------------------------------------------------
+--
+module Numeric.Search.Bounded (search, searchFrom, searchTo) where
+
+import Numeric.Search.Range
+
+-- | /O(log(abs n))/.
+-- Search a bounded integral type.
+--
+-- If @p@ is an upward-closed predicate, @search p@ returns
+-- @Just n@ if and only if @n@ is the least such satisfying @p@.
+search :: (Bounded a, Integral a) => (a -> Bool) -> Maybe a
+search p
+  | p 0 = Just (searchDown p minBound 0)
+  | otherwise = searchUp p 1 maxBound
+
+-- | /O(log(abs n))/.
+-- Search the part of a bounded integral type from a given point.
+--
+-- If @p@ is an upward-closed predicate, @searchFrom p l@ returns
+-- @Just n@ if and only if @n@ is the least @n >= l@ satisfying @p@.
+searchFrom :: (Bounded a, Integral a) => (a -> Bool) -> a -> Maybe a
+searchFrom p l
+  | l <= 0 && p 0 = Just (searchDown p l 0)
+  | otherwise = searchUp p (max 1 l) maxBound
+
+-- | /O(log(abs n))/.
+-- Search the part of a bounded integral type up to a given point.
+--
+-- If @p@ is an upward-closed predicate, @searchTo p h@ returns
+-- @Just n@ if and only if @n@ is the least @n <= h@ satisfying @p@.
+searchTo :: (Bounded a, Integral a) => (a -> Bool) -> a -> Maybe a
+searchTo p h
+  | p h' = Just (searchDown p minBound h')
+  | otherwise = searchUp p 1 h
+  where h' = min 0 h
+
+-- @h <= 0 && p h@
+searchDown :: (Integral a) => (a -> Bool) -> a -> a -> a
+searchDown p l h
+  | l `quot` 2 >= h = searchSafeRange p l h
+  | p h' = searchDown p l h'
+  | otherwise = searchSafeRange p (h'+1) h
+  where h' = h*2 - 1
+
+-- @0 < l@
+searchUp :: (Integral a) => (a -> Bool) -> a -> a -> Maybe a
+searchUp p l h
+  | h `div` 2 <= l = searchFromTo p l h
+  | p l' = Just (searchSafeRange p l l')
+  | otherwise = searchUp p (l'+1) h
+  where l' = l*2 + 1
+
+-- | Like 'search', but assuming @l <= h && p h@.
+searchSafeRange :: Integral a => (a -> Bool) -> a -> a -> a
+searchSafeRange p l h
+  | l == h = l
+  | p m = searchSafeRange p l m
+  | otherwise = searchSafeRange p (m+1) h
+  -- Stay within @min 0 l .. max 0 h@ to avoid overflow.
+  where m = l `div` 2 + h `div` 2	-- If l < h, then l <= m < h
diff --git a/Numeric/Search/Integer.hs b/Numeric/Search/Integer.hs
new file mode 100644
--- /dev/null
+++ b/Numeric/Search/Integer.hs
@@ -0,0 +1,125 @@
+-----------------------------------------------------------------------------
+-- |
+-- Module      :  Numeric.Search.Integer
+-- Copyright   :  (c) Ross Paterson 2008
+-- License     :  BSD-style
+-- Maintainer  :  ross@soi.city.ac.uk
+-- Stability   :  experimental
+-- Portability :  portable
+--
+-- Searching unbounded intervals of integers for the boundary of an
+-- upward-closed set, using a combination of exponential and binary
+-- search.
+--
+-----------------------------------------------------------------------------
+
+module Numeric.Search.Integer (
+	-- * One-dimensional searches
+	search, searchFrom, searchTo,
+	-- * Two-dimensional searches
+	search2) where
+
+import Data.Maybe (fromMaybe)
+
+-- | /O(log(abs n))/.
+-- Search the integers.
+--
+-- If @p@ is an upward-closed predicate, @search p@ returns the least
+-- @n@ satisfying @p@.  If no such @n@ exists, either because no integer
+-- satisfies @p@ or all do, @search p@ does not terminate.
+--
+-- For example, the following function computes discrete logarithms (base 2):
+--
+-- > discreteLog :: Integer -> Integer
+-- > discreteLog n = search (\ k -> 2^k <= n)
+--
+search :: (Integer -> Bool) -> Integer
+search p = fromMaybe (searchFrom p 1) (searchTo p 0)
+
+-- | /O(log(n-l))/.
+-- Search the integers from a given lower bound.
+--
+-- If @p@ is an upward-closed predicate,
+-- @searchFrom p l = 'search' (\\ i -> i >= l && p i)@.
+-- If no such @n@ exists (because no integer satisfies @p@),
+-- @searchFrom p@ does not terminate.
+searchFrom :: (Integer -> Bool) -> Integer -> Integer
+searchFrom p = search_from 1
+  where search_from step l
+	  | p l' = searchIntegerRange p l (l'-1)
+	  | otherwise = search_from (2*step) (l'+1)
+	  where l' = l + step
+
+-- | /O(log(h-n))/.
+-- Search the integers up to a given upper bound.
+--
+-- If @p@ is an upward-closed predicate, @searchTo p h == 'Just' n@
+-- if and only if @n@ is the least number @n <= h@ satisfying @p@.
+searchTo :: (Integer -> Bool) -> Integer -> Maybe Integer
+searchTo p h0
+  | p h0 = Just (search_to 1 h0)
+  | otherwise = Nothing
+  where search_to step h		-- @step >= 1 && p h@
+	  | p h' = search_to (2*step) h'
+	  | otherwise = searchSafeRange p (h'+1) h
+	  where h' = h - step
+
+-- | /O(m log(n\/m))/.
+-- Two-dimensional search, using an algorithm due described in
+--
+-- * Richard Bird, /Saddleback search: a lesson in algorithm design/,
+--   in /Mathematics of Program Construction/ 2006,
+--   Springer LNCS vol. 4014, pp82-89.
+--
+-- If @p@ is closed upwards in each argument on non-negative integers,
+-- @search2 p@ returns the minimal non-negative pairs satisfying @p@,
+-- listed in order of increasing x-coordinate.
+--
+-- /m/ and /n/ refer to the smaller and larger dimensions of the
+-- rectangle containing the boundary.
+--
+-- For example,
+--
+-- > search2 (\ x y -> x^2 + y^2 >= 25)  ==  [(0,5),(3,4),(4,3),(5,0)]
+--
+search2 :: (Integer -> Integer -> Bool) -> [(Integer,Integer)]
+search2 p = search2Rect p 0 0 hx hy []
+  where	hx = searchFrom (\ x -> p x 0) 0
+	hy = searchFrom (\ y -> p 0 y) 0
+
+search2Rect :: (Integer -> Integer -> Bool) ->
+	Integer -> Integer -> Integer -> Integer ->
+	[(Integer,Integer)] -> [(Integer,Integer)]
+search2Rect p lx ly hx hy
+  | lx > hx || ly > hy = id
+  | lx == hx && ly == hy = if p lx ly then ((lx, ly) :) else id
+  | hx-lx > hy-ly =
+	let	mx = (lx+hx) `div` 2
+		my = searchIntegerRange (\ y -> p mx y) ly hy
+	in search2Rect p lx my mx hy . search2Rect p (mx+1) ly hx (my-1)
+  | otherwise =
+	let	mx = searchIntegerRange (\ x -> p x my) lx hx
+		my = (ly+hy) `div` 2
+	in search2Rect p lx (my+1) (mx-1) hy . search2Rect p mx ly hx my
+
+-- | Search a bounded interval of integers.
+--
+-- If @p@ is an upward-closed predicate,
+--
+-- > searchIntegerRange p l h = 'search' (\ i -> i >= l && p i || i > h)
+--
+-- Cost: /O(log(h-l))/ calls to @p@.
+searchIntegerRange :: (Integer -> Bool) -> Integer -> Integer -> Integer
+searchIntegerRange p l h
+  | h < l = h+1
+  | p m = searchIntegerRange p l (m-1)
+  | otherwise = searchIntegerRange p (m+1) h
+  where m = (l+h) `div` 2
+
+-- | Like 'search', but assuming @l <= h && p h@.
+searchSafeRange :: (Integer -> Bool) -> Integer -> Integer -> Integer
+searchSafeRange p l h
+  | l == h = l
+  | p m = searchSafeRange p l m
+  | otherwise = searchSafeRange p (m+1) h
+  where m = (l + h) `div` 2	-- If l < h, then l <= m < h
diff --git a/Numeric/Search/Range.hs b/Numeric/Search/Range.hs
new file mode 100644
--- /dev/null
+++ b/Numeric/Search/Range.hs
@@ -0,0 +1,47 @@
+-----------------------------------------------------------------------------
+-- |
+-- Module      :  Numeric.Search.Range
+-- Copyright   :  (c) Ross Paterson 2008
+-- License     :  BSD-style
+-- Maintainer  :  ross@soi.city.ac.uk
+-- Stability   :  experimental
+-- Portability :  portable
+--
+-- Binary search of a bounded interval of an integral type for the
+-- boundary of an upward-closed set, using a combination of exponential
+-- and binary search.
+--
+-----------------------------------------------------------------------------
+
+module Numeric.Search.Range (searchFromTo) where
+
+-- | /O(log(h-l))/.
+-- Search a bounded interval of some integral type.
+--
+-- If @p@ is an upward-closed predicate, @searchFromTo p l h == Just n@
+-- if and only if @n@ is the least number @l <= n <= h@ satisfying @p@.
+--
+-- For example, the following function determines the first index (if any)
+-- at which a value occurs in an ordered array:
+--
+-- > searchArray :: Ord a => a -> Array Int a -> Maybe Int
+-- > searchArray x array = do
+-- >   let (lo, hi) = bounds array
+-- >   k <- searchFromTo (\ i -> array!i >= x) lo hi
+-- >   guard (array!k == x)
+-- >   return k
+--
+searchFromTo :: Integral a => (a -> Bool) -> a -> a -> Maybe a
+searchFromTo p l h
+  | l > h = Nothing
+  | p h = Just (searchSafeRange p l h)
+  | otherwise = Nothing
+
+-- | Like 'searchFromTo', but assuming @l <= h && p h@.
+searchSafeRange :: Integral a => (a -> Bool) -> a -> a -> a
+searchSafeRange p l h
+  | l == h = l
+  | p m = searchSafeRange p l m
+  | otherwise = searchSafeRange p (m+1) h
+  -- Stay within @min 0 l .. max 0 h@ to avoid overflow.
+  where m = l `div` 2 + h `div` 2	-- If l < h, then l <= m < h
diff --git a/Setup.hs b/Setup.hs
new file mode 100644
--- /dev/null
+++ b/Setup.hs
@@ -0,0 +1,2 @@
+import Distribution.Simple
+main = defaultMainWithHooks defaultUserHooks
diff --git a/binary-search.cabal b/binary-search.cabal
new file mode 100644
--- /dev/null
+++ b/binary-search.cabal
@@ -0,0 +1,21 @@
+Name:           binary-search
+Version:        0.0
+Build-Depends:  base
+Build-Type:     Simple
+License:        BSD3
+license-file:   LICENSE
+Author:         Ross Paterson <ross@soi.city.ac.uk>
+Maintainer:     Ross Paterson <ross@soi.city.ac.uk>
+Category:       Algorithms
+Synopsis:       Binary and exponential searches
+Description:    These modules address the problem of finding the boundary
+                of an upward-closed set of integers, using a combination
+                of exponential and binary searches.  Variants are provided
+                for searching within bounded and unbounded intervals of
+                both 'Integer' and bounded integral types.
+Exposed-Modules:
+                Numeric.Search.Bounded
+                Numeric.Search.Integer
+                Numeric.Search.Range
+Extra-Source-Files:
+                search-test.hs
diff --git a/search-test.hs b/search-test.hs
new file mode 100644
--- /dev/null
+++ b/search-test.hs
@@ -0,0 +1,95 @@
+module Main (main) where
+
+import Test.QuickCheck
+import Numeric.Search.Bounded as B
+import Numeric.Search.Integer as I
+import Numeric.Search.Range
+
+main :: IO ()
+main = flip mapM_ tests $ \ (Test n t) -> do
+	putStrLn $ "Testing: " ++ n
+	t
+
+data Test = Test String (IO ())
+
+mkTest :: Testable a => String -> a -> Test
+mkTest n t = Test n (test t)
+
+tests :: [Test]
+tests = [
+	mkTest "searchIntegers" prop_searchIntegers,
+	mkTest "searchIntegersFrom" prop_searchIntegersFrom,
+	mkTest "searchIntegersTo" prop_searchIntegersTo,
+	mkTest "searchIntegersTo (const False)" prop_searchIntegersToF,
+	mkTest "searchFromTo" prop_searchFromTo,
+	mkTest "searchFromTo (const False)" prop_searchFromToF,
+	mkTest "searchBounded" prop_searchBounded,
+	mkTest "searchBounded (const False)" prop_searchBoundedF,
+	mkTest "searchBoundedFrom" prop_searchBoundedFrom,
+	mkTest "searchBoundedFrom (const False)" prop_searchBoundedFromF,
+	mkTest "searchBoundedTo" prop_searchBoundedTo,
+	mkTest "searchBoundedTo (const False)" prop_searchBoundedToF]
+
+-- Every upward closed predicate is equivalent to either (const False),
+-- or (>= n) for some n.
+
+prop_searchIntegers :: Integer -> Bool
+prop_searchIntegers n =
+	I.search (>= n)  ==  n
+
+-- I.search (const False) does not terminate
+
+--	I.searchFrom p l  ==  I.search (\ i -> i >= l && p i)
+
+prop_searchIntegersFrom :: Integer -> Integer -> Bool
+prop_searchIntegersFrom n l =
+	I.searchFrom (>= n) l  ==  max l n
+
+-- I.searchFrom (const False) l does not terminate
+
+--	I.searchTo p h  ==  if n > h then Nothing else Just n
+--		let k = I.search (\ i -> i > h || p i)
+--		in if k <= h then Just k else Nothing
+
+prop_searchIntegersTo :: Integer -> Integer -> Bool
+prop_searchIntegersTo n h =
+	I.searchTo (>= n) h  ==  if n <= h then Just n else Nothing
+
+prop_searchIntegersToF :: Integer -> Bool
+prop_searchIntegersToF h =
+	I.searchTo (const False) h  ==  Nothing
+
+--	searchFromTo p l h  ==  I.search (\ i -> i < l || i <= h && p i)
+
+prop_searchFromTo :: Int -> Int -> Int -> Bool
+prop_searchFromTo n l h =
+	searchFromTo (>= n) l h  ==  if k <= h then Just k else Nothing
+  where k = max n l
+
+prop_searchFromToF :: Int -> Int -> Bool
+prop_searchFromToF l h =
+	searchFromTo (const False) l h  ==  Nothing
+
+prop_searchBounded :: Int -> Bool
+prop_searchBounded n =
+	B.search (>= n)  ==  Just n
+
+prop_searchBoundedF :: Bool
+prop_searchBoundedF =
+	B.search (const False :: Int -> Bool)  ==  Nothing
+
+prop_searchBoundedFrom :: Int -> Int -> Bool
+prop_searchBoundedFrom n l =
+	B.searchFrom (>= n) l  ==  Just (max l n) 
+
+prop_searchBoundedFromF :: Int -> Bool
+prop_searchBoundedFromF l =
+	B.searchFrom (const False) l  ==  Nothing
+
+prop_searchBoundedTo :: Int -> Int -> Bool
+prop_searchBoundedTo n h =
+	B.searchTo (>= n) h  ==  if n <= h then Just n else Nothing
+
+prop_searchBoundedToF :: Int -> Bool
+prop_searchBoundedToF h =
+	B.searchTo (const False) h  ==  Nothing
