diff --git a/CHANGELOG.markdown b/CHANGELOG.markdown
--- a/CHANGELOG.markdown
+++ b/CHANGELOG.markdown
@@ -1,3 +1,11 @@
+4.1
+---
+* Fixed a bug in the type of `conjugateGradientAscent` and `conjugateGradientDescent` that prevent users from being able to ever call it.
+
+4.0.0.1
+-------
+* Added the missing `instances.h` header file to `extra-source-files`.
+
 4.0
 ---
 * An overhaul permitting monomorphic modes was completed by @alang9.
diff --git a/README.markdown b/README.markdown
--- a/README.markdown
+++ b/README.markdown
@@ -60,7 +60,7 @@
     [0.8414709848078965,0.5403023058681398,-0.8414709848078965,-0.5403023058681398,0.8414709848078965,0.5403023058681398,-0.8414709848078965,-0.5403023058681398,0.8414709848078965,0.5403023058681398]
 
 or if your function takes multiple inputs, you can use grads, which returns an 'f-branching stream' of derivatives, that you can
-inspect lazily. Somewhat more intuitive answers can be obtained by converting the stream into the polymorphically recursive 
+inspect lazily. Somewhat more intuitive answers can be obtained by converting the stream into the polymorphically recursive
 `Jet` data type. With that we can look at a single "layer" of the answer at a time:
 
 The answer:
@@ -78,12 +78,12 @@
     Prelude Numeric.AD Numeric.AD.Types> headJet $ tailJet $ tailJet $ jet $  grads (\[x,y] -> exp (x * y)) [1,2]
     [[29.5562243957226,22.16716829679195],[22.16716829679195,7.38905609893065]]
 
-Or even higher order tensors of derivatives.
+Or even higher order tensors of derivatives as a jet.
 
     Prelude Numeric.AD Numeric.AD.Types> headJet $ tailJet $ tailJet $ tailJet $ jet $  grads (\[x,y] -> exp (x * y)) [1,2]
     [[[59.1124487914452,44.3343365935839],[44.3343365935839,14.7781121978613]],[[44.3343365935839,14.7781121978613],[14.7781121978613,7.38905609893065]]]
 
-Note the redundant values caused by the various symmetries in the tensors. The `ad` library is careful to compute 
+Note the redundant values caused by the various symmetries in the tensors. The `ad` library is careful to compute
 each distinct derivative only once, lazily and to share the resulting computation.
 
 Overview
diff --git a/ad.cabal b/ad.cabal
--- a/ad.cabal
+++ b/ad.cabal
@@ -1,5 +1,5 @@
 name:         ad
-version:      4.0.0.1
+version:      4.1
 license:      BSD3
 license-File: LICENSE
 copyright:    (c) Edward Kmett 2010-2014,
@@ -11,7 +11,7 @@
 homepage:     http://github.com/ekmett/ad
 bug-reports:  http://github.com/ekmett/ad/issues
 build-type:   Custom
-cabal-version: >= 1.8
+cabal-version: >= 1.10
 extra-source-files:
   .ghci
   .gitignore
@@ -80,9 +80,10 @@
   location: git://github.com/ekmett/ad.git
 
 library
-  extensions: CPP
+  default-extensions: CPP
   hs-source-dirs: src
   include-dirs: include
+  default-language: Haskell2010
 
   other-extensions:
     BangPatterns
@@ -137,6 +138,7 @@
     Numeric.AD.Internal.Identity
     Numeric.AD.Internal.Kahn
     Numeric.AD.Internal.On
+    Numeric.AD.Internal.Or
     Numeric.AD.Internal.Reverse
     Numeric.AD.Internal.Sparse
     Numeric.AD.Internal.Tower
diff --git a/src/Numeric/AD.hs b/src/Numeric/AD.hs
--- a/src/Numeric/AD.hs
+++ b/src/Numeric/AD.hs
@@ -127,7 +127,6 @@
   , gradientAscent
   , conjugateGradientDescent
   , conjugateGradientAscent
-
   ) where
 
 import Control.Applicative
diff --git a/src/Numeric/AD/Internal/Or.hs b/src/Numeric/AD/Internal/Or.hs
new file mode 100644
--- /dev/null
+++ b/src/Numeric/AD/Internal/Or.hs
@@ -0,0 +1,203 @@
+{-# LANGUAGE CPP #-}
+{-# LANGUAGE GADTs #-}
+{-# LANGUAGE TypeFamilies #-}
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE UndecidableInstances #-}
+#if __GLASGOW_HASKELL__ >= 707
+{-# LANGUAGE DeriveDataTypeable #-}
+#endif
+{-# OPTIONS_HADDOCK not-home #-}
+
+-----------------------------------------------------------------------------
+-- |
+-- Copyright   :  (c) Edward Kmett 2014
+-- License     :  BSD3
+-- Maintainer  :  ekmett@gmail.com
+-- Stability   :  experimental
+-- Portability :  GHC only
+--
+-----------------------------------------------------------------------------
+
+module Numeric.AD.Internal.Or
+  ( Or(..)
+  , F, T
+  , runL, runR
+  , Chosen(..)
+  , chosen
+  , unary
+  , binary
+  ) where
+
+import Control.Applicative
+import Data.Number.Erf
+#if __GLASGOW_HASKELL__ >= 707
+import Data.Typeable
+#endif
+import Numeric.AD.Mode
+
+runL :: Or a b F -> a
+runL (L a) = a
+
+runR :: Or a b T -> b
+runR (R b) = b
+
+------------------------------------------------------------------------------
+-- On
+------------------------------------------------------------------------------
+
+chosen :: (a -> r) -> (b -> r) -> Or a b s -> r
+chosen f _ (L a) = f a
+chosen _ g (R b) = g b
+
+unary :: (a -> a) -> (b -> b) -> Or a b s -> Or a b s
+unary f _ (L a) = L (f a)
+unary _ g (R a) = R (g a)
+
+binary :: (a -> a -> a) -> (b -> b -> b) -> Or a b s -> Or a b s -> Or a b s
+binary f _ (L a) (L b) = L (f a b)
+binary _ g (R a) (R b) = R (g a b)
+binary _ _ _ _ = impossible
+
+data F
+data T
+
+class Chosen s where
+  choose :: a -> b -> Or a b s
+
+instance Chosen F where
+  choose x _ = L x
+
+instance Chosen T where
+  choose _ x = R x
+
+#ifndef HLINT
+-- | The choice between two AD modes is an AD mode in its own right
+data Or a b s where
+  L :: a -> Or a b F
+  R :: b -> Or a b T
+#if __GLASGOW_HASKELL__ >= 707
+  deriving Typeable
+#endif
+#endif
+
+impossible :: a
+impossible = error "Numeric.AD.Internal.Or: impossible case"
+
+instance (Eq a, Eq b) => Eq (Or a b s) where
+  L a == L b = a == b
+  R a == R b = a == b
+  _ == _ = impossible
+
+instance (Ord a, Ord b) => Ord (Or a b s) where
+  L a `compare` L b = compare a b
+  R a `compare` R b = compare a b
+  _ `compare` _ = impossible
+
+instance (Enum a, Enum b, Chosen s) => Enum (Or a b s) where
+  pred = unary pred pred
+  succ = unary succ succ
+  toEnum i = choose (toEnum i) (toEnum i)
+  fromEnum = chosen fromEnum fromEnum
+  enumFrom (L a) = L <$> enumFrom a
+  enumFrom (R a) = R <$> enumFrom a
+  enumFromThen (L a) (L b) = L <$> enumFromThen a b
+  enumFromThen (R a) (R b) = R <$> enumFromThen a b
+  enumFromThen _     _     = impossible
+  enumFromTo (L a) (L b) = L <$> enumFromTo a b
+  enumFromTo (R a) (R b) = R <$> enumFromTo a b
+  enumFromTo _     _     = impossible
+  enumFromThenTo (L a) (L b) (L c) = L <$> enumFromThenTo a b c
+  enumFromThenTo (R a) (R b) (R c) = R <$> enumFromThenTo a b c
+  enumFromThenTo _     _     _     = impossible
+
+instance (Bounded a, Bounded b, Chosen s) => Bounded (Or a b s) where
+  maxBound = choose maxBound maxBound
+  minBound = choose minBound minBound
+
+instance (Num a, Num b, Chosen s) => Num (Or a b s) where
+  (+) = binary (+) (+)
+  (-) = binary (-) (-)
+  (*) = binary (*) (*)
+  negate = unary negate negate
+  abs = unary abs abs
+  signum = unary signum signum
+  fromInteger = choose <$> fromInteger <*> fromInteger
+
+instance (Real a, Real b, Chosen s) => Real (Or a b s) where
+  toRational = chosen toRational toRational
+
+instance (Fractional a, Fractional b, Chosen s) => Fractional (Or a b s) where
+  (/) = binary (/) (/)
+  recip = unary recip recip
+  fromRational = choose <$> fromRational <*> fromRational
+
+instance (RealFrac a, RealFrac b, Chosen s) => RealFrac (Or a b s) where
+  properFraction (L a) = case properFraction a of
+    (b, c) -> (b, L c)
+  properFraction (R a) = case properFraction a of
+    (b, c) -> (b, R c)
+  truncate = chosen truncate truncate
+  round = chosen round round
+  ceiling = chosen ceiling ceiling
+  floor = chosen floor floor
+
+instance (Floating a, Floating b, Chosen s) => Floating (Or a b s) where
+  pi = choose pi pi
+  exp = unary exp exp
+  sqrt = unary sqrt sqrt
+  log = unary log log
+  (**) = binary (**) (**)
+  logBase = binary logBase logBase
+  sin = unary sin sin
+  tan = unary tan tan
+  cos = unary cos cos
+  asin = unary asin asin
+  atan = unary atan atan
+  acos = unary acos acos
+  sinh = unary sinh sinh
+  tanh = unary tanh tanh
+  cosh = unary cosh cosh
+  asinh = unary asinh asinh
+  atanh = unary atanh atanh
+  acosh = unary acosh acosh
+
+instance (Erf a, Erf b, Chosen s) => Erf (Or a b s) where
+  erf = unary erf erf
+  erfc = unary erfc erfc
+  erfcx = unary erfcx erfcx
+  normcdf = unary normcdf normcdf
+
+instance (InvErf a, InvErf b, Chosen s) => InvErf (Or a b s) where
+  inverf = unary inverf inverf
+  inverfc = unary inverfc inverfc
+  invnormcdf = unary invnormcdf invnormcdf
+
+instance (RealFloat a, RealFloat b, Chosen s) => RealFloat (Or a b s) where
+  floatRadix = chosen floatRadix floatRadix
+  floatDigits = chosen floatDigits floatDigits
+  floatRange = chosen floatRange floatRange
+  decodeFloat = chosen decodeFloat decodeFloat
+  encodeFloat i j = choose (encodeFloat i j) (encodeFloat i j)
+  exponent = chosen exponent exponent
+  significand = unary significand significand
+  scaleFloat = unary <$> scaleFloat <*> scaleFloat
+  isNaN = chosen isNaN isNaN
+  isInfinite = chosen isInfinite isInfinite
+  isDenormalized = chosen isDenormalized isDenormalized
+  isNegativeZero = chosen isNegativeZero isNegativeZero
+  isIEEE = chosen isIEEE isIEEE
+  atan2 = binary atan2 atan2
+
+type instance Scalar (Or a b s) = Scalar a
+
+instance (Mode a, Mode b, Chosen s, Scalar a ~ Scalar b) => Mode (Or a b s) where
+  auto = choose <$> auto <*> auto
+  isKnownConstant = chosen isKnownConstant isKnownConstant
+  isKnownZero = chosen isKnownZero isKnownZero
+  x *^ L a = L (x *^ a)
+  x *^ R a = R (x *^ a)
+  L a ^* x = L (a ^* x)
+  R a ^* x = R (a ^* x)
+  L a ^/ x = L (a ^/ x)
+  R a ^/ x = R (a ^/ x)
+  zero = choose zero zero
diff --git a/src/Numeric/AD/Newton.hs b/src/Numeric/AD/Newton.hs
--- a/src/Numeric/AD/Newton.hs
+++ b/src/Numeric/AD/Newton.hs
@@ -33,9 +33,11 @@
 import Data.Traversable
 import Numeric.AD.Mode
 import Numeric.AD.Mode.Forward (diff, diff')
-import Numeric.AD.Mode.Reverse (grad, gradWith')
+import Numeric.AD.Mode.Reverse as Reverse (gradWith')
+import Numeric.AD.Mode.Kahn as Kahn (Kahn, grad)
 import Numeric.AD.Internal.Combinators
 import Numeric.AD.Internal.Forward (Forward)
+import Numeric.AD.Internal.Or
 import Numeric.AD.Internal.On
 import Numeric.AD.Internal.Reverse (Reverse, Tape)
 
@@ -108,7 +110,7 @@
 gradientDescent :: (Traversable f, Fractional a, Ord a) => (forall s. Reifies s Tape => f (Reverse a s) -> Reverse a s) -> f a -> [f a]
 gradientDescent f x0 = go x0 fx0 xgx0 0.1 (0 :: Int)
   where
-    (fx0, xgx0) = gradWith' (,) f x0
+    (fx0, xgx0) = Reverse.gradWith' (,) f x0
     go x fx xgx !eta !i
       | eta == 0     = [] -- step size is 0
       | fx1 > fx     = go x fx xgx (eta/2) 0 -- we stepped too far
@@ -119,7 +121,7 @@
       where
         zeroGrad = all (\(_,g) -> g == 0)
         x1 = fmap (\(xi,gxi) -> xi - eta * gxi) xgx
-        (fx1, xgx1) = gradWith' (,) f x1
+        (fx1, xgx1) = Reverse.gradWith' (,) f x1
 {-# INLINE gradientDescent #-}
 
 -- | Perform a gradient descent using reverse mode automatic differentiation to compute the gradient.
@@ -135,22 +137,34 @@
 -- 1
 -- >>> rosenbrock (conjugateGradientDescent rosenbrock [0, 0] !! 5) < 0.1
 -- True
-conjugateGradientDescent :: (Traversable f, Ord a, Fractional a) => (forall t. (Mode t, a ~ Scalar t, Num t) => f t -> t) -> f a -> [f a]
+conjugateGradientDescent
+  :: (Traversable f, Ord a, Fractional a)
+  => (forall s1 s2 s3 s4. Chosen s4 => f (Or (On (Forward (Forward a s1) s2)) (Kahn a s3) s4) -> Or (On (Forward (Forward a s1) s2)) (Kahn a s3) s4)
+  -> f a -> [f a]
 conjugateGradientDescent f = conjugateGradientAscent (negate . f)
 {-# INLINE conjugateGradientDescent #-}
 
+lfu :: Functor f => (f (Or a b F) -> Or a b F) -> f a -> a
+lfu f = runL . f . fmap L
+
+rfu :: Functor f => (f (Or a b T) -> Or a b T) -> f b -> b
+rfu f = runR . f . fmap R
+
 -- | Perform a conjugate gradient ascent using reverse mode automatic differentiation to compute the gradient.
-conjugateGradientAscent :: (Traversable f, Ord a, Fractional a) => (forall t. (Mode t, a ~ Scalar t, Num t) => f t -> t) -> f a -> [f a]
+conjugateGradientAscent
+  :: (Traversable f, Ord a, Fractional a)
+  => (forall s1 s2 s3 s4. Chosen s4 => f (Or (On (Forward (Forward a s1) s2)) (Kahn a s3) s4) -> Or (On (Forward (Forward a s1) s2)) (Kahn a s3) s4)
+  -> f a -> [f a]
 conjugateGradientAscent f x0 = takeWhile (all (\a -> a == a)) (go x0 d0 d0 delta0)
   where
     dot x y = sum $ zipWithT (*) x y
-    d0 = grad f x0
+    d0 = Kahn.grad (rfu f) x0
     delta0 = dot d0 d0
     go xi _ri di deltai = xi : go xi1 ri1 di1 deltai1
       where
-        ai = last $ take 20 $ extremum (\a -> f $ zipWithT (\x d -> auto x + a * auto d) xi di) 0
+        ai = last $ take 20 $ extremum (\a -> lfu f $ zipWithT (\x d -> auto x + a * auto d) xi di) 0
         xi1 = zipWithT (\x d -> x + ai*d) xi di
-        ri1 = grad f xi1
+        ri1 = Kahn.grad (rfu f) xi1
         deltai1 = dot ri1 ri1
         bi1 = deltai1 / deltai
         di1 = zipWithT (\r d -> r + bi1 * d) ri1 di