diff --git a/CHANGELOG.markdown b/CHANGELOG.markdown
--- a/CHANGELOG.markdown
+++ b/CHANGELOG.markdown
@@ -1,3 +1,8 @@
+3.1.4
+-----
+* Added a better "convergence" test for `findZero`
+* Compute `tan` and `tanh` derivatives directly.
+
 3.1.3
 -----
 * Added `conjugateGradientDescent` and `conjugateGradientAscent` to `Numeric.AD.Newton`.
diff --git a/ad.cabal b/ad.cabal
--- a/ad.cabal
+++ b/ad.cabal
@@ -1,5 +1,5 @@
 name:         ad
-version:      3.1.3
+version:      3.1.4
 license:      BSD3
 license-File: LICENSE
 copyright:    (c) Edward Kmett 2010-2012,
@@ -138,10 +138,12 @@
   type:    exitcode-stdio-1.0
   main-is: doctests.hs
   build-depends:
-    base == 4.*,
-    directory >= 1.0 && < 1.2,
-    doctest >= 0.8 && <= 0.9,
-    filepath >= 1.3 && < 1.4,
+    base,
+    directory,
+    doctest >= 0.9.0.1 && <= 0.10,
+    filepath,
     mtl
-  ghc-options: -Wall -Werror -threaded
+  ghc-options: -Wall -threaded
+  if impl(ghc<7.6)
+    ghc-options: -Werror
   hs-source-dirs: tests
diff --git a/src/Numeric/AD/Halley.hs b/src/Numeric/AD/Halley.hs
--- a/src/Numeric/AD/Halley.hs
+++ b/src/Numeric/AD/Halley.hs
@@ -32,7 +32,8 @@
 
 -- | The 'findZero' function finds a zero of a scalar function using
 -- Halley's method; its output is a stream of increasingly accurate
--- results.  (Modulo the usual caveats.)
+-- results.  (Modulo the usual caveats.) If the stream becomes constant
+-- ("it converges"), no further elements are returned.
 --
 -- Examples:
 --
@@ -43,16 +44,16 @@
 -- >>> last $ take 10 $ findZero ((+1).(^2)) (1 :+ 1)
 -- 0.0 :+ 1.0
 findZero :: (Fractional a, Eq a) => (forall s. Mode s => AD s a -> AD s a) -> a -> [a]
-findZero f = go
-    where
-        go x = x : if y == 0 then [] else go (x - 2*y*y'/(2*y'*y'-y*y''))
-            where
-                (y:y':y'':_) = diffs0 f x
+findZero f = go where
+  go x = x : if x == xn then [] else go xn where
+    (y:y':y'':_) = diffs0 f x
+    xn = x - 2*y*y'/(2*y'*y'-y*y'')
 {-# INLINE findZero #-}
 
 -- | The 'inverse' function inverts a scalar function using
 -- Halley's method; its output is a stream of increasingly accurate
--- results.  (Modulo the usual caveats.)
+-- results.  (Modulo the usual caveats.) If the stream becomes constant
+-- ("it converges"), no further elements are returned.
 --
 -- Note: the @take 10 $ inverse sqrt 1 (sqrt 10)@ example that works for Newton's method
 -- fails with Halley's method because the preconditions do not hold!
@@ -64,6 +65,9 @@
 -- function using Halley's method; its output is a stream of
 -- increasingly accurate results.  (Modulo the usual caveats.)
 --
+-- If the stream becomes constant ("it converges"), no further
+-- elements are returned.
+--
 -- >>> last $ take 10 $ fixedPoint cos 1
 -- 0.7390851332151607
 fixedPoint :: (Fractional a, Eq a) => (forall s. Mode s => AD s a -> AD s a) -> a -> [a]
@@ -72,11 +76,11 @@
 
 -- | The 'extremum' function finds an extremum of a scalar
 -- function using Halley's method; produces a stream of increasingly
--- accurate results.  (Modulo the usual caveats.)
+-- accurate results.  (Modulo the usual caveats.) If the stream becomes
+-- constant ("it converges"), no further elements are returned.
 --
 -- >>> take 10 $ extremum cos 1
 -- [1.0,0.29616942658570555,4.59979519460002e-3,1.6220740159042513e-8,0.0]
 extremum :: (Fractional a, Eq a) => (forall s. Mode s => AD s a -> AD s a) -> a -> [a]
 extremum f = findZero (diff (decomposeMode . f . composeMode))
 {-# INLINE extremum #-}
-
diff --git a/src/Numeric/AD/Internal/Classes.hs b/src/Numeric/AD/Internal/Classes.hs
--- a/src/Numeric/AD/Internal/Classes.hs
+++ b/src/Numeric/AD/Internal/Classes.hs
@@ -227,13 +227,13 @@
         --   | otherwise         = lift2_ (**) (\z xi yi -> (yi *! z /! xi, z *! log1 xi)) x y
         sin1      = lift1 sin cos1
         cos1      = lift1 cos $ negate1 . sin1
-        tan1 x    = sin1 x /! cos1 x
+        tan1      = lift1 tan $ recip1 . square1 . cos1
         asin1     = lift1 asin $ \x -> recip1 (sqrt1 (one -! square1 x))
         acos1     = lift1 acos $ \x -> negate1 (recip1 (sqrt1 (one -! square1 x)))
         atan1     = lift1 atan $ \x -> recip1 (one +! square1 x)
         sinh1     = lift1 sinh cosh1
         cosh1     = lift1 cosh sinh1
-        tanh1 x   = sinh1 x /! cosh1 x
+        tanh1     = lift1 tanh $ recip1 . square1 . cosh1
         asinh1    = lift1 asinh $ \x -> recip1 (sqrt1 (one +! square1 x))
         acosh1    = lift1 acosh $ \x -> recip1 (sqrt1 (square1 x -! one))
         atanh1    = lift1 atanh $ \x -> recip1 (one -! square1 x)
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
@@ -36,7 +36,8 @@
 
 -- | The 'findZero' function finds a zero of a scalar function using
 -- Newton's method; its output is a stream of increasingly accurate
--- results.  (Modulo the usual caveats.)
+-- results.  (Modulo the usual caveats.) If the stream becomes constant
+-- ("it converges"), no further elements are returned.
 --
 -- Examples:
 --
@@ -47,16 +48,16 @@
 -- >>> last $ take 10 $ findZero ((+1).(^2)) (1 :+ 1)
 -- 0.0 :+ 1.0
 findZero :: (Fractional a, Eq a) => (forall s. Mode s => AD s a -> AD s a) -> a -> [a]
-findZero f = go
-    where
-        go x = x : if y == 0 then [] else go (x - y/y')
-            where
-                (y,y') = diff' f x
+findZero f = go where
+  go x = x : if x == xn then [] else go xn where
+    (y,y') = diff' f x
+    xn = x - y/y'
 {-# INLINE findZero #-}
 
 -- | The 'inverse' function inverts a scalar function using
 -- Newton's method; its output is a stream of increasingly accurate
--- results.  (Modulo the usual caveats.)
+-- results.  (Modulo the usual caveats.) If the stream becomes
+-- constant ("it converges"), no further elements are returned.
 --
 -- Example:
 --
@@ -70,6 +71,9 @@
 -- function using Newton's method; its output is a stream of
 -- increasingly accurate results.  (Modulo the usual caveats.)
 --
+-- If the stream becomes constant ("it converges"), no further
+-- elements are returned.
+--
 -- >>> last $ take 10 $ fixedPoint cos 1
 -- 0.7390851332151607
 fixedPoint :: (Fractional a, Eq a) => (forall s. Mode s => AD s a -> AD s a) -> a -> [a]
@@ -78,7 +82,8 @@
 
 -- | The 'extremum' function finds an extremum of a scalar
 -- function using Newton's method; produces a stream of increasingly
--- accurate results.  (Modulo the usual caveats.)
+-- accurate results.  (Modulo the usual caveats.) If the stream
+-- becomes constant ("it converges"), no further elements are returned.
 --
 -- >>> last $ take 10 $ extremum cos 1
 -- 0.0
@@ -110,11 +115,12 @@
                 (fx1, xgx1) = gradWith' (,) f x1
 {-# INLINE gradientDescent #-}
 
+-- | Perform a gradient descent using reverse mode automatic differentiation to compute the gradient.
 gradientAscent :: (Traversable f, Fractional a, Ord a) => (forall s. Mode s => f (AD s a) -> AD s a) -> f a -> [f a]
 gradientAscent f = gradientDescent (negate . f)
 {-# INLINE gradientAscent #-}
 
-
+-- | Perform a conjugate gradient descent using reverse mode automatic differentiation to compute the gradient.
 conjugateGradientDescent :: (Traversable f, Fractional a, Ord a) =>
                 (forall s. Mode s => f (AD s a) -> AD s a) -> f a -> [f a]
 conjugateGradientDescent f x0 = go x0 d0 d0
@@ -131,6 +137,7 @@
         di1 = zipWithT (\r d -> r * bi1*d) ri1 di
 {-# INLINE conjugateGradientDescent #-}
 
+-- | Perform a conjugate gradient ascent using reverse mode automatic differentiation to compute the gradient.
 conjugateGradientAscent :: (Traversable f, Fractional a, Ord a) => (forall s. Mode s => f (AD s a) -> AD s a) -> f a -> [f a]
 conjugateGradientAscent f = conjugateGradientDescent (negate . f)
 {-# INLINE conjugateGradientAscent #-}