diff --git a/Changelog.md b/Changelog.md
--- a/Changelog.md
+++ b/Changelog.md
@@ -1,3 +1,11 @@
+0.5
+
+* add explicit (**) implementation for AD (@msakai https://github.com/ocramz/ad-delcont/pull/4 )
+
+0.4
+
+* add unit tests
+
 0.3
 
 * add `grad`
diff --git a/README.md b/README.md
--- a/README.md
+++ b/README.md
@@ -1,5 +1,5 @@
 # ad-delcont
 
-Reverse-mode Automatic Differentiation using delimited continuations (`shift` and `reset`) , as shown in [1]
+Reverse-mode Automatic Differentiation using delimited continuations (`shift` and `reset`)
 
 Introductory blog post : http://ocramz.github.io/haskell/automatic-differentiation/2021/07/19/ad-delcont.html
diff --git a/ad-delcont.cabal b/ad-delcont.cabal
--- a/ad-delcont.cabal
+++ b/ad-delcont.cabal
@@ -1,5 +1,5 @@
 name:                ad-delcont
-version:             0.3.0.0
+version:             0.5.0.0
 synopsis:            Reverse-mode automatic differentiation with delimited continuations
 description:         Reverse-mode automatic differentiation using delimited continuations (@shift@/@reset@).
                      The package exposes a small and easily extensible user interface to automatic differentiation combinators. It's also lightweight as a dependency, since it only requires @base@ and @transformers@.
@@ -37,18 +37,23 @@
   other-modules:       Numeric.AD.DelCont.Internal
   build-depends:       base >= 4.7 && < 5
                      , transformers >= 0.5
+                     -- -- DEBUG
+                     -- , ad 
+                     -- , hspec >= 2.7.10
 
--- -- TODO
--- test-suite spec
---   default-language:    Haskell2010
---   ghc-options:         -Wall
---   type:                exitcode-stdio-1.0
---   hs-source-dirs:      test
---   main-is:             Spec.hs
---   build-depends:       base
---                      , ad-delcont
---                      , hspec
---                      , QuickCheck
+-- TODO
+test-suite spec
+  default-language:    Haskell2010
+  ghc-options:         -Wall
+  type:                exitcode-stdio-1.0
+  hs-source-dirs:      test
+  main-is:             Spec.hs
+  other-modules:       Numeric.AD.DelContSpec
+  build-depends:       base
+                     , ad-delcont
+                     , hspec
+                     --
+                     , ad
 
 source-repository head
   type:     git
diff --git a/src/Numeric/AD/DelCont/Internal.hs b/src/Numeric/AD/DelCont/Internal.hs
--- a/src/Numeric/AD/DelCont/Internal.hs
+++ b/src/Numeric/AD/DelCont/Internal.hs
@@ -3,8 +3,6 @@
 {-# LANGUAGE DeriveFunctor #-}
 {-# language GeneralizedNewtypeDeriving #-}
 {-# LANGUAGE LambdaCase #-}
--- {-# LANGUAGE MultiParamTypeClasses #-}
--- {-# language TypeFamilies #-}
 {-# OPTIONS_GHC -Wno-unused-imports -Wno-unused-top-binds #-}
 module Numeric.AD.DelCont.Internal
   (rad1, rad2, grad,
@@ -25,7 +23,6 @@
 
 import Prelude hiding (read)
 
-
 -- | Dual numbers
 data D a da = D { primal :: a, dual :: da } deriving (Show, Functor)
 instance Eq a => Eq (D a da) where
@@ -51,8 +48,9 @@
 --
 -- As one expects from a constant, its value will be used for computing the result, but it will be discarded when computing the sensitivities.
 auto :: a -> AD s a da
-auto x = AD0 $ lift $ var x undefined
-
+auto x = AD0 $ lift $ var x undefined -- NB blows up with -XStrict
+autoStrict :: a -> da -> AD0 s (DVar s a da)
+autoStrict x dx = AD0 $ lift $ var x dx
 
 -- | Mutable references in the continuation monad
 newtype AD0 s a = AD0 { unAD0 :: forall x . ContT x (ST s) a } deriving (Functor)
@@ -182,6 +180,7 @@
   exp = op1Num $ \x -> (exp x, (exp x *))
   log = op1Num $ \x -> (log x, (/x))
   sqrt = op1Num $ \x -> (sqrt x, (/ (2 * sqrt x)))
+  (**) = op2Num $ \x y -> (x ** y, (* (y * x ** (y - 1))), (* (x ** y * log x)))
   logBase = op2Num $ \x y ->
                        let
                          dx = - logBase x y / (log x * x)
@@ -272,6 +271,22 @@
   xs_bar <- traverse readSTRef xrs
   let xs_bar_d = dual <$> xs_bar
   pure (z, xs_bar_d)
+
+
+-- jacg zeroa onea f xs = runST $ do -- -- Jacobian TODO
+--   xrs <- traverse (`var` zeroa) xs
+--   zr' <- evalContT $
+--     resetT $ do
+--       let
+--         zads = f (fmap pure xrs) -- traversable of AD results
+--       for zads $ \zad -> do
+--         zr <- zad
+--         lift $ modifySTRef' zr (withD (const onea))
+--         pure zr
+--   undefined
+
+for :: (Applicative f, Traversable t) => t a -> (a -> f b) -> f (t b)
+for = flip traverse
 
 
 -- | Evaluate (forward mode) and differentiate (reverse mode) a unary function
diff --git a/test/Numeric/AD/DelContSpec.hs b/test/Numeric/AD/DelContSpec.hs
new file mode 100644
--- /dev/null
+++ b/test/Numeric/AD/DelContSpec.hs
@@ -0,0 +1,49 @@
+{-# options_ghc -Wno-type-defaults -Wno-unused-imports -Wno-missing-methods -Wno-orphans #-}
+module Numeric.AD.DelContSpec where
+
+import Numeric.AD.DelCont (grad, rad1)
+import qualified Numeric.AD as AD (grad)
+import Test.Hspec (Spec, hspec, describe, it, shouldBe, shouldSatisfy)
+
+
+spec :: Spec
+spec = do
+  describe "rad1 : Unary functions" $ do
+    let
+      x0 = 1.2
+      x1 = 0.0
+    it "(** 2)" $ do
+      let
+        f x = x ** 2
+        (_, dfdx) = rad1 f x0
+        (_, dfdx1) = rad1 f x1
+      dfdx `shouldBe` 2.4
+      dfdx1 `shouldBe` 0.0
+    it "sqrt" $ do
+      let
+        (_, dfdx) = rad1 sqrt x0
+        xhat = 1 / (2 * (sqrt x0))
+      dfdx `shouldBe` xhat
+    it "reciprocal" $ do
+      let
+        (_, dfdx) = rad1 recip x0
+        xhat = negate $ 1 / (x0 **2)
+      dfdx `shouldBe` xhat
+  describe "grad : Multivariate functions" $ do
+    let
+      x = [1.2, 1.3]
+    it "inverse of the L2 norm" $ do
+      let
+        f z = 1 / norm z
+        (_, gradf) = grad f x -- 'ad-delcont'
+        gradAD = AD.grad f x -- 'ad'
+      norm (gradf - gradAD) `shouldSatisfy` (<= 1e-12)
+
+instance Num a => Num [a] where
+  (-) = zipWith (-)
+
+norm :: Floating a => [a] -> a
+norm = sqrt . dot
+
+dot :: Num c => [c] -> c
+dot xs = sum $ zipWith (*) xs xs
diff --git a/test/Spec.hs b/test/Spec.hs
new file mode 100644
--- /dev/null
+++ b/test/Spec.hs
@@ -0,0 +1,1 @@
+{-# OPTIONS_GHC -F -pgmF hspec-discover #-}