diff --git a/ChangeLog.md b/ChangeLog.md
--- a/ChangeLog.md
+++ b/ChangeLog.md
@@ -2,5 +2,12 @@
 # Changelog for circular
 
 
+## 0.1.1
+
+-   Remove `mean`.
+-   Add benchmark.
+-   Many small improvements.
+
+
 ## Unreleased changes
 
diff --git a/README.md b/README.md
--- a/README.md
+++ b/README.md
@@ -10,15 +10,16 @@
 
 -   memory usage is constant
 -   they are fast, especially when summary statistics need to be
-    computed across the stack
+    computed across the stack (use `unsafePush`, if possible)
 -   they can be saved, and restored using JSON format
 
-When the stack is full, new elements pushed on the stack replace the oldest
-(deepest) elements on the stack. Complex circular behavior can arise when pushes
-and pops are mixed. QuickCheck and unit tests with HSpec give promising results
-— have a look yourself.
+When the stack is full, new, pushed elements replace the oldest (deepest)
+elements on the stack. Complex circular behavior can arise when pushes and pops
+are mixed. QuickCheck and unit tests with HSpec give promising results — have
+a look yourself.
 
-I use them, for example, as the data type for traces of Markov chains.
+I use circular stacks, for example, as the data type for traces of Markov
+chains.
 
 `Circular` is actively developed and functions may be removed, renamed, or
 changed. New ideas are welcome!
diff --git a/bench/Bench.hs b/bench/Bench.hs
new file mode 100644
--- /dev/null
+++ b/bench/Bench.hs
@@ -0,0 +1,54 @@
+-- |
+-- Module      :  Main
+-- Description :  Benchmark circular stacks
+-- Copyright   :  (c) Dominik Schrempf, 2020
+-- License     :  GPL-3.0-or-later
+--
+-- Maintainer  :  dominik.schrempf@gmail.com
+-- Stability   :  unstable
+-- Portability :  portable
+--
+-- Creation date: Sat Jun 20 21:12:38 2020.
+module Main
+  ( main
+  )
+where
+
+import Criterion.Main
+import qualified Data.Stack.Circular as S
+import qualified Data.Vector.Unboxed as U
+import qualified Data.Vector as V
+
+-- When using foldl or foldl', list is much slower than cstack.
+
+list ::  Int -> Int
+list l = sum $ take 1000 $ foldl (flip (:)) [] [0..l]
+
+cstackV ::  Int -> Int
+cstackV l = S.sum $ foldl (flip S.unsafePush) (S.unsafeEmpty 1000 :: S.CStack V.Vector Int) [0..l]
+
+cstackU ::  Int -> Int
+cstackU l = S.sum $ foldl (flip S.unsafePush) (S.unsafeEmpty 1000 :: S.CStack U.Vector Int) [0..l]
+
+-- When using foldr, cstack is slower by far. This is because of the
+-- lazyness. However, for stacks, by definition, the last added elements are
+-- of interest.
+
+-- -- The safe operations are very slow.
+
+-- cstackVSafe ::  Int -> Int
+-- cstackVSafe l = S.sum $ foldl (flip S.push) (S.empty 1000 :: S.CStack V.Vector Int) [0..l]
+
+-- cstackUSafe ::  Int -> Int
+-- cstackUSafe l = S.sum $ foldl (flip S.push) (S.empty 1000 :: S.CStack U.Vector Int) [0..l]
+
+
+main :: IO ()
+main = do
+  let l = 1000000 :: Int
+  print $ list l
+  print $ cstackU l
+  defaultMain
+    [ bench "list, foldl" $ whnf list l
+    , bench "cstack, foldl" $ whnf cstackV l
+    , bench "cstack unboxed, foldl" $ whnf cstackU l ]
diff --git a/circular.cabal b/circular.cabal
--- a/circular.cabal
+++ b/circular.cabal
@@ -5,7 +5,7 @@
 -- see: https://github.com/sol/hpack
 
 name:           circular
-version:        0.1.0
+version:        0.1.1
 synopsis:       Circular fixed-sized mutable vectors
 description:    Please see the README on GitHub at <https://github.com/dschrempf/circular#readme>
 category:       Math, Data Structures
@@ -56,5 +56,21 @@
     , hspec
     , hspec-discover
     , quickcheck-instances
+    , vector
+  default-language: Haskell2010
+
+benchmark circular-bench
+  type: exitcode-stdio-1.0
+  main-is: Bench.hs
+  other-modules:
+      Paths_circular
+  hs-source-dirs:
+      bench
+  ghc-options: -Wall -threaded -rtsopts -with-rtsopts=-N
+  build-depends:
+      aeson
+    , base >=4.7 && <5
+    , circular
+    , criterion
     , vector
   default-language: Haskell2010
diff --git a/src/Data/Stack/Circular.hs b/src/Data/Stack/Circular.hs
--- a/src/Data/Stack/Circular.hs
+++ b/src/Data/Stack/Circular.hs
@@ -19,17 +19,21 @@
 
     -- * Construction
     empty,
+    unsafeEmpty,
 
     -- * Conversion
     toVector,
     toVectorN,
+    unsafeToVectorN,
     fromVector,
+    unsafeFromVector,
 
     -- * Accessors
     get,
     pop,
     push,
     unsafePush,
+    reset,
 
     -- * Queries
     isFull,
@@ -43,9 +47,10 @@
     -- For reasons of efficiency, __commutativity__ of the combining function is
     -- __assumed__ for fold-like functions provided in this section! That is,
     -- the order of elements of the stack must not matter.
+    foldl,
+    foldl',
     foldl1',
     sum,
-    mean,
     product,
   )
 where
@@ -56,7 +61,7 @@
 import qualified Data.Vector.Generic as V
 import Data.Vector.Generic (Vector)
 import qualified Data.Vector.Generic.Mutable as M
-import Prelude hiding (product, sum)
+import Prelude hiding (foldl, product, sum)
 
 -- | Circular stacks with fxed maximum size are just normal vectors with a
 -- pointer to the last element.
@@ -109,6 +114,12 @@
   | m <= i + 1 = i + 1 - m
   | otherwise = i + 1 - m + n
 
+-- -- Do not check for empty stack.
+-- unsafeStartIndex :: Int -> Int -> Int -> Int
+-- unsafeStartIndex i m n
+--   | m <= i + 1 = i + 1 - m
+--   | otherwise = i + 1 - m + n
+
 -- | A circular stack without an element but of a given maximum size. At this
 -- state, it is not very useful :). O(n).
 empty :: Vector v a => Int -> CStack v a
@@ -116,6 +127,10 @@
   | n <= 0 = error "empty: maximum size must be 1 or larger"
   | otherwise = CStack (V.create $ M.unsafeNew n) 0 0
 
+-- | See 'empty'; do no check that length is strictly positive.
+unsafeEmpty :: Vector v a => Int -> CStack v a
+unsafeEmpty n = CStack (V.create $ M.unsafeNew n) 0 0
+
 -- | Convert a circular stack to a vector. The first element of the returned
 -- vector is the deepest (oldest) element of the stack, the last element of the
 -- returned vector is the current (newest) element of the stack.
@@ -140,7 +155,7 @@
 -- This is a relatively expensive operation. O(N).
 toVectorN :: Vector v a => Int -> CStack v a -> v a
 toVectorN k (CStack v i m)
-  | k < 0 = error "toVectorN: negative n"
+  | k < 0 = error "toVectorN: negative N"
   | k > m = error "toVectorN: stack too small"
   | k == 0 = V.empty
   | i' + k <= n = V.unsafeSlice i' k v
@@ -149,6 +164,17 @@
     n = V.length v
     i' = startIndex i k n
 
+-- | See 'toVectorN' but do not check that N is positive.
+unsafeToVectorN :: Vector v a => Int -> CStack v a -> v a
+unsafeToVectorN k (CStack v i m)
+  | k > m = error "toVectorN: stack too small"
+  | k == 0 = V.empty
+  | i' + k <= n = V.unsafeSlice i' k v
+  | otherwise = V.unsafeDrop i' v V.++ V.unsafeTake (i + 1) v
+  where
+    n = V.length v
+    i' = startIndex i k n
+
 -- | Convert a vector to a circular stack. The first element of the vector is
 -- the deepest (oldest) element of the stack, the last element of the vector is
 -- the current (newest) element of the stack. O(n).
@@ -161,9 +187,16 @@
   where
     n = V.length v
 
+-- | See 'fromVector' but do not check for empty vector.
+unsafeFromVector :: Vector v a => v a -> CStack v a
+unsafeFromVector v = CStack v (n - 1) n
+  where
+    n = V.length v
+
 -- | Get the last element without changing the stack. O(1).
 get :: Vector v a => CStack v a -> a
 get (CStack v i _) = V.unsafeIndex v i
+{-# INLINE get #-}
 
 -- Select the previous element without changing the stack.
 previous :: Vector v a => CStack v a -> CStack v a
@@ -179,6 +212,7 @@
 -- The stack must be non-empty.
 pop :: Vector v a => CStack v a -> (a, CStack v a)
 pop c = (get c, previous c)
+{-# INLINE pop #-}
 
 -- Replace an element in a vector.
 set :: Vector v a => Int -> a -> v a -> v a
@@ -197,7 +231,7 @@
   where
     n = V.length v
 
--- | Push an element on the stack. O(n).
+-- | Push an element on the stack. Slow! If possible, use 'unsafePush'. O(n).
 push :: Vector v a => a -> CStack v a -> CStack v a
 push x c = put x $ next c
 
@@ -216,14 +250,37 @@
 -- Be careful; the internal vector is mutated! The immutable circular stack may
 -- not be used after this operation.
 unsafePush :: Vector v a => a -> CStack v a -> CStack v a
-unsafePush x c = unsafePut x $ next c
+unsafePush x = unsafePut x . next
 
+-- | Reset the stack. O(1).
+reset :: CStack v a -> CStack v a
+reset (CStack v _ _) = CStack v 0 0
+
 -- | Check if the stack is full.
 isFull :: Vector v a => CStack v a -> Bool
 isFull (CStack v _ m) = V.length v == m
 
--- | Compute summary statistics of the elements on the stack using a custom
--- commutative `mappend` function.
+-- | Left fold. O(m).
+foldl :: (Vector v a, Vector v b) => (a -> b -> a) -> a -> CStack v b -> a
+foldl f x (CStack v i m)
+  | m == n = V.foldl f x v
+  | i' + m <= n = V.foldl f x $ V.unsafeSlice i' m v
+  | otherwise = V.foldl f (V.foldl f x (V.unsafeDrop i' v)) (V.unsafeTake (i + 1) v)
+  where
+    n = V.length v
+    i' = startIndex i m n
+
+-- | Left fold with strict accumulator. O(m).
+foldl' :: (Vector v a, Vector v b) => (a -> b -> a) -> a -> CStack v b -> a
+foldl' f x (CStack v i m)
+  | m == n = V.foldl' f x v
+  | i' + m <= n = V.foldl' f x $ V.unsafeSlice i' m v
+  | otherwise = V.foldl' f (V.foldl' f x (V.unsafeDrop i' v)) (V.unsafeTake (i + 1) v)
+  where
+    n = V.length v
+    i' = startIndex i m n
+
+-- | Left fold on non-empty vectors with strict accumulator. O(m).
 foldl1' :: Vector v a => (a -> a -> a) -> CStack v a -> a
 foldl1' f (CStack v i m)
   | m == n = V.foldl1' f v
@@ -233,7 +290,7 @@
     n = V.length v
     i' = startIndex i m n
 
--- | Compute the sum of the elements on the stack.
+-- | Compute the sum of the elements on the stack. O(m).
 sum :: (Num a, Vector v a) => CStack v a -> a
 sum (CStack v i m)
   | m == n = V.sum v
@@ -243,14 +300,7 @@
     n = V.length v
     i' = startIndex i m n
 
--- | Compute the mean of the elements on the stack.
-mean :: (Real a, Vector v a, Fractional b) => CStack v a -> b
-mean c = realToFrac (sum c) / fromIntegral (curSize c)
-
--- | Compute the product of the elements on the stack.
---
--- For reasons of efficiency, commutativity of the combining function is
--- assumed. That is, the order of elements of the stack must not matter.
+-- | Compute the product of the elements on the stack. O(m).
 product :: (Num a, Vector v a) => CStack v a -> a
 product (CStack v i m)
   | m == n = V.product v