```{-# OPTIONS_GHC -fno-warn-name-shadowing #-}
module Math.Root.Bracket where

-- |Predicate that returns 'True' whenever the given pair of points brackets
-- a root of the given function.
brackets :: (Eq a, Eq b, Num b) => (a -> b) -> (a,a) -> Bool
brackets f (x1,x2)
| x1 == x2  = f x1 == 0
| otherwise = signum (f x1) /= signum (f x2)

-- |@bracket f x1 x2@: Given a function and an initial guessed range x1 to x2,
-- this function expands the range geometrically until a root is bracketed by
-- the returned values, returning a list of the successively expanded ranges.
-- The list will be finite if and only if the sequence yields a bracketing pair.
bracket :: (Fractional a, Eq a, Num b, Ord b) =>
(a -> b) -> a -> a -> [(a, a)]
bracket f x1 x2
| x1 == x2  = error "bracket: empty range"
| otherwise = go x1 (f x1) x2 (f x2)
where
factor = 1.618 -- golden ratio (close enough to it, anyway)
go x1 f1 x2 f2
| signum f1 /= signum f2    = [(x1, x2)]
| abs f1 < abs f2           = (x1, x2) : go x1' (f x1') x2 f2
| otherwise                 = (x1, x2) : go x1 f1 x2' (f x2')
where
x1' = x1 - factor * w
x2' = x2 + factor * w
w = x2 - x1

-- |@subdivideAndBracket f x1 x2 n@: Given a function defined on the interval
-- [x1,x2], subdivide the interval into n equally spaced segments and search
-- for zero crossings of the function.  The returned list will contain all
-- bracketing pairs found.
subdivideAndBracket :: (Num b, Eq b, Fractional a, Integral c) =>
(a -> b) -> a -> a -> c -> [(a, a)]
subdivideAndBracket f x1 x2 n =
[ (x1, x2)
| ((x1, y1), (x2, y2)) <- zip xys (tail xys)
, signum y1 /= signum y2
]
where
dx = (x2 - x1) / fromIntegral n
xs = x1 : [x1 + dx * fromIntegral i | i <- [1..n]]
xys = map (\x -> (x, f x)) xs

```