diff --git a/CHANGELOG.md b/CHANGELOG.md
new file mode 100644
--- /dev/null
+++ b/CHANGELOG.md
@@ -0,0 +1,5 @@
+# accelerate-kullback-liebler
+
+## 0.1.0.0
+
+Initial release
diff --git a/LICENSE b/LICENSE
new file mode 100644
--- /dev/null
+++ b/LICENSE
@@ -0,0 +1,11 @@
+Copyright Vanessa McHale (c) 2019
+
+Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
+
+1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
+
+2. 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.
+
+3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "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 COPYRIGHT HOLDER OR 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/README.md b/README.md
new file mode 100644
--- /dev/null
+++ b/README.md
@@ -0,0 +1,3 @@
+# accelerate-kullback-liebler
+
+Kullback-Liebler divergence on arrays.
diff --git a/accelerate-kullback-liebler.cabal b/accelerate-kullback-liebler.cabal
new file mode 100644
--- /dev/null
+++ b/accelerate-kullback-liebler.cabal
@@ -0,0 +1,137 @@
+cabal-version:   2.0
+name:            accelerate-kullback-liebler
+version:         0.1.0.0
+license:         BSD3
+license-file:    LICENSE
+copyright:       Copyright: (c) 2019 Vanessa McHale
+maintainer:      vamchale@gmail.com
+author:          Vanessa McHale
+synopsis:        Kullback-Liebler divergence
+description:     Kullback-Liebler divergence using accelerate
+category:        Accelerate, Array
+build-type:      Simple
+extra-doc-files:
+    README.md
+    CHANGELOG.md
+
+source-repository head
+    type:     git
+    location: https://github.com/vmchale/accelerate-kullback-liebler
+
+flag development
+    description: Enable `-Werror`
+    default:     False
+    manual:      True
+
+flag gpu
+    description: Enable GPU backend (only affects test/benchmark suite)
+
+library
+    exposed-modules:  Data.Array.Accelerate.KullbackLiebler
+    hs-source-dirs:   src
+    default-language: Haskell2010
+    other-extensions: FlexibleContexts ConstraintKinds
+    ghc-options:      -Wall
+    build-depends:
+        base >=4.3 && <5,
+        accelerate >=1.1.0.0
+
+    if flag(development)
+        ghc-options: -Werror
+
+    if impl(ghc >=8.0)
+        ghc-options:
+            -Wincomplete-uni-patterns -Wincomplete-record-updates
+            -Wredundant-constraints -Widentities
+
+    if impl(ghc >=8.4)
+        ghc-options: -Wmissing-export-lists
+
+library accelerate-test
+    exposed-modules:  Data.Array.Accelerate.Test
+    hs-source-dirs:   common
+    default-language: Haskell2010
+    ghc-options:      -Wall
+    build-depends:
+        base >=4.3 && <5,
+        mwc-random-accelerate -any,
+        accelerate >=1.1.0.0
+
+    if flag(development)
+        ghc-options: -Werror
+
+    if impl(ghc >=8.0)
+        ghc-options:
+            -Wincomplete-uni-patterns -Wincomplete-record-updates
+            -Wredundant-constraints -Widentities
+
+    if impl(ghc >=8.4)
+        ghc-options: -Wmissing-export-lists
+
+test-suite accelerate-kullback-liebler-test
+    type:             exitcode-stdio-1.0
+    main-is:          Spec.hs
+    hs-source-dirs:   test
+    default-language: Haskell2010
+    ghc-options:      -threaded -rtsopts -with-rtsopts=-N -Wall
+    build-depends:
+        base -any,
+        accelerate-kullback-liebler -any,
+        tasty -any,
+        composition-prelude >=2.0.2.0,
+        tasty-hedgehog -any,
+        tasty-hunit -any,
+        accelerate-llvm-native -any,
+        accelerate-test -any,
+        accelerate -any
+
+    if flag(development)
+        ghc-options: -Werror
+
+    if flag(gpu)
+        cpp-options:   -DACCELERATE_GPU
+        build-depends: accelerate-llvm-ptx -any
+
+    if impl(ghc >=8.0)
+        ghc-options:
+            -Wincomplete-uni-patterns -Wincomplete-record-updates
+            -Wredundant-constraints -Widentities
+
+    if impl(ghc >=8.4)
+        ghc-options: -Wmissing-export-lists
+
+    if impl(ghc >=8.2)
+        ghc-options: -Wcpp-undef
+
+benchmark accelerate-kullback-liebler-bench
+    type:               exitcode-stdio-1.0
+    main-is:            Bench.hs
+    build-tool-depends: cpphs:cpphs -any
+    hs-source-dirs:     bench
+    default-language:   Haskell2010
+    ghc-options:        -threaded -rtsopts -with-rtsopts=-N -O3 -Wall
+    build-depends:
+        base -any,
+        accelerate-kullback-liebler -any,
+        accelerate-test -any,
+        accelerate-llvm-native -any,
+        accelerate -any,
+        criterion -any
+
+    if flag(gpu)
+        cpp-options:   -DACCELERATE_GPU
+        build-depends: accelerate-llvm-ptx -any
+
+    if flag(development)
+        ghc-options: -Werror
+
+    if impl(ghc >=8.0)
+        ghc-options:
+            -Wincomplete-uni-patterns -Wincomplete-record-updates
+            -Wredundant-constraints -Widentities
+
+    if impl(ghc >=8.4)
+        ghc-options: -Wmissing-export-lists
+
+    if impl(ghc >=8.2)
+        ghc-options: -Wcpp-undef
diff --git a/bench/Bench.cpphs b/bench/Bench.cpphs
new file mode 100644
--- /dev/null
+++ b/bench/Bench.cpphs
@@ -0,0 +1,47 @@
+{-# LANGUAGE FlexibleContexts #-}
+
+module Main (main) where
+
+import           Criterion.Main
+import qualified Data.Array.Accelerate                 as A
+import           Data.Array.Accelerate.KullbackLiebler
+import qualified Data.Array.Accelerate.LLVM.Native     as CPU
+#ifdef ACCELERATE_GPU
+import qualified Data.Array.Accelerate.LLVM.PTX        as GPU
+#endif
+import           Data.Array.Accelerate.Test
+
+main :: IO ()
+main = do
+    ps <- distEnv
+    qs <- distEnv
+#ifdef ACCELERATE_GPU
+    let res = kullbackLieblerGPU ps qs -- TODO: better way to force this
+#else
+    let res = kullbackLieblerCPU ps qs -- TODO: better way to force this
+#endif
+    putStrLn (res `seq` "Setup complete")
+    defaultMain [
+#ifdef ACCELERATE_GPU
+                  bgroup "GPU"
+                      [ bench "K-L (1e7)" $ nf (kullbackLieblerGPU ps) qs
+                      , bench "entropy (1e7)" $ nf entropyGPU ps
+                      ]
+                , bgroup "CPU"
+#else
+                  bgroup "CPU"
+#endif
+                      [ bench "K-L (1e7)" $ nf (kullbackLieblerCPU ps) qs
+                      , bench "entropy (1e7)" $ nf entropyCPU ps
+                      ]
+                ]
+
+    where distEnv :: IO (A.Acc (A.Vector Double))
+          distEnv = randDist (A.Z A.:. 10000000)
+#ifdef ACCELERATE_GPU
+          kullbackLieblerGPU ps qs = the $ GPU.run $ kullbackLiebler ps qs
+          entropyGPU = the . GPU.run . entropy
+#endif
+          -- TODO: runQ
+          kullbackLieblerCPU ps qs = the $ CPU.run $ kullbackLiebler ps qs
+          entropyCPU = the . CPU.run . entropy
diff --git a/common/Data/Array/Accelerate/Test.hs b/common/Data/Array/Accelerate/Test.hs
new file mode 100644
--- /dev/null
+++ b/common/Data/Array/Accelerate/Test.hs
@@ -0,0 +1,26 @@
+{-# LANGUAGE ConstraintKinds  #-}
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE TypeFamilies     #-}
+
+module Data.Array.Accelerate.Test ( normalize
+                                  , randDist
+                                  , the
+                                  ) where
+
+import qualified Data.Array.Accelerate                   as A
+import           Data.Array.Accelerate.System.Random.MWC (Variate, randomArray,
+                                                          uniformR)
+
+-- | Doesn't actually check the list has one element
+the :: A.Scalar e -> e
+the = head . A.toList
+
+-- | Doesn't check for negative values
+normalize :: A.Floating e => A.Acc (A.Vector e) -> A.Acc (A.Vector e)
+normalize xs =
+    let tot = A.the $ A.sum xs
+    in A.map (/tot) xs
+
+-- | Make a distribution of a given size
+randDist :: (A.Shape sh, Fractional e, Variate e, A.Floating e, sh ~ A.DIM1) => sh -> IO (A.Acc (A.Vector e))
+randDist = fmap (normalize . A.use) . randomArray (uniformR (0.0, 1.0))
diff --git a/src/Data/Array/Accelerate/KullbackLiebler.hs b/src/Data/Array/Accelerate/KullbackLiebler.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Array/Accelerate/KullbackLiebler.hs
@@ -0,0 +1,21 @@
+{-# LANGUAGE ConstraintKinds  #-}
+{-# LANGUAGE FlexibleContexts #-}
+
+module Data.Array.Accelerate.KullbackLiebler ( kullbackLiebler
+                                             , entropy
+                                             , dropZeroes
+                                             ) where
+
+import qualified Data.Array.Accelerate as A
+
+-- | Assumes input is nonzero
+kullbackLiebler :: (A.Floating e) => A.Acc (A.Vector e) -> A.Acc (A.Vector e) -> A.Acc (A.Scalar e)
+kullbackLiebler ps qs = A.sum (A.zipWith (\p q -> p * log (p / q)) ps qs)
+
+-- | Assumes input is nonzero
+entropy :: (A.Floating e) => A.Acc (A.Vector e) -> A.Acc (A.Scalar e)
+entropy = A.sum . A.map (\p -> p * log p)
+
+-- | Drops zero values (for parity with scipy @entropy@ function)
+dropZeroes :: (A.Eq e, Num (A.Exp e)) => A.Acc (A.Vector e) -> A.Acc (A.Vector e)
+dropZeroes = A.afst . A.filter (A./= 0)
diff --git a/test/Spec.cpphs b/test/Spec.cpphs
new file mode 100644
--- /dev/null
+++ b/test/Spec.cpphs
@@ -0,0 +1,51 @@
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE OverloadedLists  #-}
+
+module Main ( main ) where
+
+import           Control.Composition                   ((.$))
+import qualified Data.Array.Accelerate                 as A
+import           Data.Array.Accelerate.KullbackLiebler (kullbackLiebler)
+import qualified Data.Array.Accelerate.LLVM.Native     as CPU
+#ifdef ACCELERATE_GPU
+import qualified Data.Array.Accelerate.LLVM.PTX        as GPU
+#endif
+import           Data.Array.Accelerate.Test            (randDist, the)
+import           Test.Tasty
+import           Test.Tasty.HUnit
+
+main :: IO ()
+main = defaultMain accTests
+
+accTests :: TestTree
+accTests = testGroup "Unit tests"
+    [
+#ifdef ACCELERATE_GPU
+      testCase "Setup (GPU)" $
+        assertBool "simple unit test" =<< selfGPU
+    , testCase "Setup (CPU)" $
+#else
+      testCase "Setup (CPU)" $
+#endif
+        assertBool "simple unit test" =<< selfCPU
+    ]
+
+--TODO: random gen/normalization?
+commonData :: IO (A.Acc (A.Scalar Double))
+commonData = do
+    datums <- randDist (A.Z A.:. 1000000)
+    pure $ kullbackLiebler .$ datums
+
+#ifdef ACCELERATE_GPU
+selfGPU :: IO Bool
+selfGPU = do
+    datums <- commonData
+    let res = the $ GPU.run datums
+    pure $ res == 0
+#endif
+
+selfCPU :: IO Bool
+selfCPU = do
+    datums <- commonData
+    let res = the $ CPU.run datums
+    pure $ res == 0