diff --git a/CHANGELOG.md b/CHANGELOG.md
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,3 +1,46 @@
+# 0.4.0
+
+* Made `Construct` a super class of `Mutable`
+* Reimplement a safe version of `makeLoadArray`, that is parallelizable.
+* Switch from `EltRepr r ix` to much simpler `R r`
+* Remove `Construct` instance for `M` representation.
+* `unsafeLinearSet` - length argument now accepts `Sz1` instead of an `Int`
+* Renamed:
+  * `forPrimM_` -> `forPrimM`
+  * `iforPrimM_` -> `iforPrimM`
+  * `iforLinearPrimM_` -> `iforLinearPrimM`
+* Introduced new functions that do not mutate the original array: `forPrimM_`,
+  `iforPrimM_` and `iforLinearPrimM_`
+* Addition of `readM`, `writeM`, `modifyM`, `swapM`, `modifyM_`, `swapM_`
+* Add an orphan instance of `MonadThrow` for `ST` monad for older versions of
+  `exceptions`. See [ekmett/exceptions#72](https://github.com/ekmett/exceptions/pull/72)
+* Deprecation of `read'`, `write'` `modify'` and `swap'`
+* Make `modify` accept a monadic action, rather than a pure function. Also not it returns
+  the old element.
+* Make `swap` return the swapped elements.
+* Addition of `unsafeLinearSwap` and `unsafeSwap`
+* Expose `unsafeLinearModify` and `unsafeModify`
+* Expose `Data.Massiv.Core.List`
+* Expose `indexWith`, so macro `INDEX_CHECK` from `massiv.h` could be used outside massiv.
+* Addition of `liftSz`
+* Fixed `expand*` functions by making them accept `Sz1` instead of an `Int`
+* Addition of `expandWithinM`
+* Bunch of minor fixes to `Show` instances
+* Extracted test-suite into it's own package.
+* Stop accepting computation strategy for all functions that can be performed sequentially only:
+  * `iterateN`
+  * `iiterateN`
+  * `unfoldrS_`
+  * `iunfoldrS_`
+  * `unfoldlS_`
+  * `iunfoldlS_`
+  * `makeArrayA`
+  * `makeArrayAR`
+  * `generateArrayLinearS`
+  * `generateArrayS`
+* Redefined most of the numeric operators with `Numeric` and `NumericFloat`. Will be
+  required for SIMD operations.
+
 # 0.3.6
 
 * Addition of `unsafeArrayLinearCopy`, `unsafeLinearCopy`, `unsafeLinearShrink`, `unsafeLinearGrow`
diff --git a/massiv.cabal b/massiv.cabal
--- a/massiv.cabal
+++ b/massiv.cabal
@@ -1,5 +1,5 @@
 name:                massiv
-version:             0.3.6.0
+version:             0.4.0.0
 synopsis:            Massiv (Массив) is an Array Library.
 description:         Multi-dimensional Arrays with fusion, stencils and parallel computation.
 homepage:            https://github.com/lehins/massiv
@@ -42,6 +42,8 @@
                      , Data.Massiv.Array.Unsafe
                      , Data.Massiv.Core
                      , Data.Massiv.Core.Index
+                     , Data.Massiv.Core.List
+                     , Data.Massiv.Core.Operations
 
   other-modules:       Data.Massiv.Array.Delayed.Interleaved
                      , Data.Massiv.Array.Delayed.Pull
@@ -69,7 +71,6 @@
                      , Data.Massiv.Core.Index.Stride
                      , Data.Massiv.Core.Index.Tuple
                      , Data.Massiv.Core.Iterator
-                     , Data.Massiv.Core.List
   build-depends:       base >= 4.9 && < 5
                      , bytestring
                      , data-default-class
@@ -80,9 +81,6 @@
                      , unliftio-core
                      , vector
 
-  if impl(ghc < 8.4)
-    build-depends: ghc-prim
-
   include-dirs: include
   install-includes: massiv.h
 
@@ -94,52 +92,6 @@
                       -Wincomplete-record-updates
                       -Wincomplete-uni-patterns
                       -Wredundant-constraints
-
-test-suite tests
-  type:               exitcode-stdio-1.0
-  hs-source-dirs:     tests
-  main-is:            Main.hs
-  other-modules:      Data.Massiv.ArraySpec
-                    , Data.Massiv.Array.DelayedSpec
-                    , Data.Massiv.Array.Delayed.InterleavedSpec
-                    , Data.Massiv.Array.Delayed.WindowedSpec
-                    , Data.Massiv.Array.Delayed.PushSpec
-                    , Data.Massiv.Array.ManifestSpec
-                    , Data.Massiv.Array.Manifest.VectorSpec
-                    , Data.Massiv.Array.MutableSpec
-                    , Data.Massiv.Array.Numeric.IntegralSpec
-                    , Data.Massiv.Array.Ops.ConstructSpec
-                    , Data.Massiv.Array.Ops.FoldSpec
-                    , Data.Massiv.Array.Ops.SliceSpec
-                    , Data.Massiv.Array.Ops.SortSpec
-                    , Data.Massiv.Array.Ops.TransformSpec
-                    , Data.Massiv.Array.Ops.MapSpec
-                    , Data.Massiv.Array.StencilSpec
-                    , Data.Massiv.CoreArbitrary
-                    , Data.Massiv.Core.IndexSpec
-                    , Data.Massiv.Core.SchedulerSpec
-                    , Spec
-  build-depends:      base
-                    , async
-                    , bytestring
-                    , containers
-                    , data-default
-                    , deepseq
-                    , massiv
-                    , hspec
-                    , scheduler
-                    , QuickCheck
-                    , unliftio
-                    , vector
-
-  default-language:   Haskell2010
-  ghc-options:       -Wall
-                     -Wincomplete-record-updates
-                     -Wincomplete-uni-patterns
-                     -Wredundant-constraints
-                     -fno-warn-orphans
-                     -threaded
-                     -with-rtsopts=-N2
 
 test-suite doctests
   type:             exitcode-stdio-1.0
diff --git a/src/Data/Massiv/Array.hs b/src/Data/Massiv/Array.hs
--- a/src/Data/Massiv/Array.hs
+++ b/src/Data/Massiv/Array.hs
@@ -89,7 +89,7 @@
   , convert
   , convertAs
   , convertProxy
-  , fromRaggedArray
+  , fromRaggedArrayM
   , fromRaggedArray'
   -- * Size
   , size
@@ -106,7 +106,6 @@
   , borderIndex
   , evaluateM
   , evaluate'
-  , evaluateAt
   -- * Mapping
   , module Data.Massiv.Array.Ops.Map
   -- * Folding
diff --git a/src/Data/Massiv/Array/Delayed/Interleaved.hs b/src/Data/Massiv/Array/Delayed/Interleaved.hs
--- a/src/Data/Massiv/Array/Delayed/Interleaved.hs
+++ b/src/Data/Massiv/Array/Delayed/Interleaved.hs
@@ -24,11 +24,9 @@
 import Data.Massiv.Core.List (L, showArrayList, showsArrayPrec)
 
 
--- | Delayed array that will be loaded in an interleaved fasion during parallel
+-- | Delayed array that will be loaded in an interleaved fashion during parallel
 -- computation.
 data DI = DI
-
-type instance EltRepr DI ix = DI
 
 newtype instance Array DI ix e = DIArray
   { diArray :: Array D ix e
diff --git a/src/Data/Massiv/Array/Delayed/Pull.hs b/src/Data/Massiv/Array/Delayed/Pull.hs
--- a/src/Data/Massiv/Array/Delayed/Pull.hs
+++ b/src/Data/Massiv/Array/Delayed/Pull.hs
@@ -21,13 +21,12 @@
   , delay
   , eq
   , ord
-  , liftArray
-  , liftArray2
   ) where
 
 import qualified Data.Foldable as F
 import Data.Massiv.Array.Ops.Fold.Internal as A
 import Data.Massiv.Core.Common
+import Data.Massiv.Core.Operations
 import Data.Massiv.Core.List (L, showArrayList, showsArrayPrec)
 import GHC.Base (build)
 import Prelude hiding (zipWith)
@@ -40,7 +39,6 @@
 data instance Array D ix e = DArray { dComp :: !Comp
                                     , dSize :: !(Sz ix)
                                     , dIndex :: ix -> e }
-type instance EltRepr D ix = D
 
 instance (Ragged L ix e, Show e) => Show (Array D ix e) where
   showsPrec = showsArrayPrec id
@@ -70,6 +68,8 @@
 instance Index ix => Source D ix e where
   unsafeIndex = INDEX_CHECK("(Source D ix e).unsafeIndex", size, dIndex)
   {-# INLINE unsafeIndex #-}
+  -- unsafeLinearSlice ix sz arr = unsafeExtract ix sz (unsafeResize sz arr)
+  -- {-# INLINE unsafeLinearSlice #-}
 
 
 instance ( Index ix
@@ -114,7 +114,7 @@
   {-# INLINE pure #-}
   (<*>) (DArray c1 (SafeSz sz1) uIndex1) (DArray c2 (SafeSz sz2) uIndex2) =
     DArray (c1 <> c2) (SafeSz (liftIndex2 min sz1 sz2)) $ \ !ix ->
-      (uIndex1 ix) (uIndex2 ix)
+      uIndex1 ix (uIndex2 ix)
   {-# INLINE (<*>) #-}
 
 
@@ -152,21 +152,21 @@
 instance Index ix => StrideLoad D ix e
 
 instance (Index ix, Num e) => Num (Array D ix e) where
-  (+)         = liftArray2 (+)
+  (+)         = unsafeLiftArray2 (+)
   {-# INLINE (+) #-}
-  (-)         = liftArray2 (-)
+  (-)         = unsafeLiftArray2 (-)
   {-# INLINE (-) #-}
-  (*)         = liftArray2 (*)
+  (*)         = unsafeLiftArray2 (*)
   {-# INLINE (*) #-}
-  abs         = liftArray abs
+  abs         = unsafeLiftArray abs
   {-# INLINE abs #-}
-  signum      = liftArray signum
+  signum      = unsafeLiftArray signum
   {-# INLINE signum #-}
   fromInteger = singleton . fromInteger
   {-# INLINE fromInteger #-}
 
 instance (Index ix, Fractional e) => Fractional (Array D ix e) where
-  (/)          = liftArray2 (/)
+  (/)          = unsafeLiftArray2 (/)
   {-# INLINE (/) #-}
   fromRational = singleton . fromRational
   {-# INLINE fromRational #-}
@@ -175,33 +175,77 @@
 instance (Index ix, Floating e) => Floating (Array D ix e) where
   pi    = singleton pi
   {-# INLINE pi #-}
-  exp   = liftArray exp
+  exp   = unsafeLiftArray exp
   {-# INLINE exp #-}
-  log   = liftArray log
+  log   = unsafeLiftArray log
   {-# INLINE log #-}
-  sin   = liftArray sin
+  sin   = unsafeLiftArray sin
   {-# INLINE sin #-}
-  cos   = liftArray cos
+  cos   = unsafeLiftArray cos
   {-# INLINE cos #-}
-  asin  = liftArray asin
+  asin  = unsafeLiftArray asin
   {-# INLINE asin #-}
-  atan  = liftArray atan
+  atan  = unsafeLiftArray atan
   {-# INLINE atan #-}
-  acos  = liftArray acos
+  acos  = unsafeLiftArray acos
   {-# INLINE acos #-}
-  sinh  = liftArray sinh
+  sinh  = unsafeLiftArray sinh
   {-# INLINE sinh #-}
-  cosh  = liftArray cosh
+  cosh  = unsafeLiftArray cosh
   {-# INLINE cosh #-}
-  asinh = liftArray asinh
+  asinh = unsafeLiftArray asinh
   {-# INLINE asinh #-}
-  atanh = liftArray atanh
+  atanh = unsafeLiftArray atanh
   {-# INLINE atanh #-}
-  acosh = liftArray acosh
+  acosh = unsafeLiftArray acosh
   {-# INLINE acosh #-}
 
 
+instance Num e => Numeric D e where
+  -- plusScalar arr e = unsafeLiftArray (+ e) arr
+  -- {-# INLINE plusScalar #-}
+  -- minusScalar arr e = unsafeLiftArray (subtract e) arr
+  -- {-# INLINE minusScalar #-}
+  -- multiplyScalar arr e = unsafeLiftArray (* e) arr
+  -- {-# INLINE multiplyScalar #-}
+  -- absPointwise = unsafeLiftArray abs
+  -- {-# INLINE absPointwise #-}
+  -- additionPointwise = unsafeLiftArray2 (+)
+  -- {-# INLINE additionPointwise #-}
+  -- subtractionPointwise = unsafeLiftArray2 (-)
+  -- {-# INLINE subtractionPointwise #-}
+  -- multiplicationPointwise = unsafeLiftArray2 (*)
+  -- {-# INLINE multiplicationPointwise #-}
+  -- powerPointwise arr pow = unsafeLiftArray (^ pow) arr
+  -- {-# INLINE powerPointwise #-}
+  -- powerSumArray arr = sumArray . powerPointwise arr
+  -- {-# INLINE powerSumArray #-}
+  -- unsafeDotProduct a1 a2 = sumArray (multiplicationPointwise a1 a2)
+  -- {-# INLINE unsafeDotProduct #-}
+  unsafeLiftArray f arr = arr {dIndex = f . dIndex arr}
+  {-# INLINE unsafeLiftArray #-}
+  unsafeLiftArray2 f a1 a2 =
+    DArray (dComp a1 <> dComp a2) (SafeSz (liftIndex2 min (unSz (dSize a1)) (unSz (dSize a2)))) $ \i ->
+      f (dIndex a1 i) (dIndex a2 i)
+  {-# INLINE unsafeLiftArray2 #-}
 
+
+instance Floating e => NumericFloat D e where
+  -- recipPointwise = liftDArray recip
+  -- {-# INLINE recipPointwise #-}
+  -- sqrtPointwise = liftDArray sqrt
+  -- {-# INLINE sqrtPointwise #-}
+  -- floorPointwise = liftDArray floor
+  -- {-# INLINE floorPointwise #-}
+  -- ceilingPointwise = liftDArray ceiling
+  -- {-# INLINE ceilingPointwise #-}
+  -- divisionPointwise = liftDArray2 (/)
+  -- {-# INLINE divisionPointwise #-}
+  -- divideScalar arr e = liftDArray (/ e) arr
+  -- {-# INLINE divideScalar #-}
+
+
+
 -- | /O(1)/ Conversion from a source array to `D` representation.
 delay :: Source r ix e => Array r ix e -> Array D ix e
 delay arr = DArray (getComp arr) (size arr) (unsafeIndex arr)
@@ -234,28 +278,29 @@
        f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))
 {-# INLINE ord #-}
 
--- | The usual map.
-liftArray :: Source r ix b => (b -> e) -> Array r ix b -> Array D ix e
-liftArray f !arr = DArray (getComp arr) (size arr) (f . unsafeIndex arr)
-{-# INLINE liftArray #-}
 
--- | Similar to `Data.Massiv.Array.zipWith`, except dimensions of both arrays either have to be the
--- same, or at least one of the two array must be a singleton array, in which case it will behave as
--- a `Data.Massiv.Array.map`.
---
--- @since 0.1.4
-liftArray2
-  :: (Source r1 ix a, Source r2 ix b)
-  => (a -> b -> e) -> Array r1 ix a -> Array r2 ix b -> Array D ix e
-liftArray2 f !arr1 !arr2
-  | sz1 == oneSz = liftArray (f (unsafeIndex arr1 zeroIndex)) arr2
-  | sz2 == oneSz = liftArray (`f` unsafeIndex arr2 zeroIndex) arr1
-  | sz1 == sz2 =
-    DArray (getComp arr1 <> getComp arr2) sz1 (\ !ix -> f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))
-  | otherwise = throw $ SizeMismatchException (size arr1) (size arr2)
-  where
-    sz1 = size arr1
-    sz2 = size arr2
-{-# INLINE liftArray2 #-}
+-- -- | The usual map.
+-- liftArray :: Source r ix b => (b -> e) -> Array r ix b -> Array D ix e
+-- liftArray f !arr = DArray (getComp arr) (size arr) (f . unsafeIndex arr)
+-- {-# INLINE liftArray #-}
+
+-- -- | Similar to `Data.Massiv.Array.zipWith`, except dimensions of both arrays either have to be the
+-- -- same, or at least one of the two array must be a singleton array, in which case it will behave as
+-- -- a `Data.Massiv.Array.map`.
+-- --
+-- -- @since 0.1.4
+-- liftArray2
+--   :: (Source r1 ix a, Source r2 ix b)
+--   => (a -> b -> e) -> Array r1 ix a -> Array r2 ix b -> Array D ix e
+-- liftArray2 f !arr1 !arr2
+--   | sz1 == oneSz = liftArray (f (unsafeIndex arr1 zeroIndex)) arr2
+--   | sz2 == oneSz = liftArray (`f` unsafeIndex arr2 zeroIndex) arr1
+--   | sz1 == sz2 =
+--     DArray (getComp arr1 <> getComp arr2) sz1 (\ !ix -> f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))
+--   | otherwise = throw $ SizeMismatchException (size arr1) (size arr2)
+--   where
+--     sz1 = size arr1
+--     sz2 = size arr2
+-- {-# INLINE liftArray2 #-}
 
 
diff --git a/src/Data/Massiv/Array/Delayed/Push.hs b/src/Data/Massiv/Array/Delayed/Push.hs
--- a/src/Data/Massiv/Array/Delayed/Push.hs
+++ b/src/Data/Massiv/Array/Delayed/Push.hs
@@ -26,11 +26,12 @@
   , fromStrideLoad
   ) where
 
+import Control.Monad
 import Data.Massiv.Core.Common
 import Data.Massiv.Core.Index.Internal (Sz(SafeSz))
-import qualified Data.Semigroup as Semigroup
 import Prelude hiding (map, zipWith)
-import Control.Applicative
+import Control.Scheduler as S (traverse_)
+import Data.Foldable as F
 
 #include "massiv.h"
 
@@ -49,8 +50,6 @@
               -> m ()
   }
 
-type instance EltRepr DL ix = DL
-
 instance Index ix => Construct DL ix e where
   setComp c arr = arr {dlComp = c}
   {-# INLINE setComp #-}
@@ -64,27 +63,63 @@
   {-# INLINE unsafeResize #-}
 
 instance Semigroup (Array DL Ix1 e) where
-  (<>) (DLArray c1 sz1 def1 load1) (DLArray c2 sz2 def2 load2) =
-    DLArray
-      {dlComp = c1 <> c2, dlSize = SafeSz (k + unSz sz2), dlDefault = def1 <|> def2, dlLoad = load}
-    where
-      !k = unSz sz1
-      load :: Monad m => Scheduler m () -> Int -> (Int -> e -> m ()) -> m ()
-      load scheduler startAt dlWrite = do
-        load1 scheduler startAt dlWrite
-        load2 scheduler (startAt + k) dlWrite
-      {-# INLINE load #-}
+  (<>) = mappendDL
   {-# INLINE (<>) #-}
 
-
 instance Monoid (Array DL Ix1 e) where
-  mempty = makeArray Seq zeroSz (const (throwImpossible Uninitialized))
+  mempty =
+    DLArray
+      {dlComp = mempty, dlSize = Sz zeroIndex, dlDefault = Nothing, dlLoad = \_ _ _ -> pure ()}
   {-# INLINE mempty #-}
-
-  mappend = (Semigroup.<>)
+  mappend = mappendDL
   {-# INLINE mappend #-}
+  mconcat [] = mempty
+  mconcat [x] = x
+  mconcat [x, y] = x <> y
+  mconcat xs = mconcatDL xs
+  {-# INLINE mconcat #-}
 
--- | Describe how an array should be loaded into memory
+mconcatDL :: forall e . [Array DL Ix1 e] -> Array DL Ix1 e
+mconcatDL !arrs =
+  DLArray {dlComp = foldMap getComp arrs, dlSize = SafeSz k, dlDefault = Nothing, dlLoad = load}
+  where
+    !k = F.foldl' (+) 0 (unSz . size <$> arrs)
+    load :: Monad m => Scheduler m () -> Int -> (Int -> e -> m ()) -> m ()
+    load scheduler startAt dlWrite =
+      let loadArr !startAtCur DLArray {dlSize = SafeSz kCur, dlDefault, dlLoad} = do
+            let !endAtCur = startAtCur + kCur
+            scheduleWork_ scheduler $ do
+              S.traverse_
+                (\def -> loopM_ startAtCur (< endAtCur) (+ 1) (`dlWrite` def))
+                dlDefault
+              dlLoad scheduler startAtCur dlWrite
+            pure endAtCur
+          {-# INLINE loadArr #-}
+       in foldM_ loadArr startAt arrs
+    {-# INLINE load #-}
+{-# INLINE mconcatDL #-}
+
+
+mappendDL :: forall e . Array DL Ix1 e -> Array DL Ix1 e -> Array DL Ix1 e
+mappendDL (DLArray c1 sz1 mDef1 load1) (DLArray c2 sz2 mDef2 load2) =
+  DLArray {dlComp = c1 <> c2, dlSize = SafeSz (k1 + k2), dlDefault = Nothing, dlLoad = load}
+  where
+    !k1 = unSz sz1
+    !k2 = unSz sz2
+    load :: Monad m => Scheduler m () -> Int -> (Int -> e -> m ()) -> m ()
+    load scheduler startAt dlWrite = do
+      scheduleWork_ scheduler $ do
+        S.traverse_ (\def1 -> loopM_ startAt (< k1) (+ 1) (`dlWrite` def1)) mDef1
+        load1 scheduler startAt dlWrite
+      scheduleWork_ scheduler $ do
+        let startAt2 = startAt + k1
+        S.traverse_ (\def2 -> loopM_ startAt2 (< startAt2 + k2) (+ 1) (`dlWrite` def2)) mDef2
+        load2 scheduler startAt2 dlWrite
+    {-# INLINE load #-}
+{-# INLINE mappendDL #-}
+
+-- | Describe how an array should be loaded into memory sequentially. For parallelizable
+-- version see `makeLoadArray`.
 --
 -- @since 0.3.1
 makeLoadArrayS ::
@@ -92,7 +127,7 @@
      Sz ix
   -- ^ Size of the resulting array
   -> e
-  -- ^ Default value to use for all cells that have possibly been ommitted by the writing function
+  -- ^ Default value to use for all cells that might have been ommitted by the writing function
   -> (forall m. Monad m => (ix -> e -> m Bool) -> m ())
   -- ^ Writing function that described which elements to write into the target array.
   -> Array DL ix e
@@ -105,28 +140,60 @@
      in writer safeWrite
 {-# INLINE makeLoadArrayS #-}
 
--- | Specify how an array can be loaded/computed through creation of a `DL` array.
+-- | Specify how an array should be loaded into memory. Unlike `makeLoadArrayS`, loading
+-- function accepts a scheduler, thus can be parallelized. If you need an unsafe version
+-- of this function see `unsafeMakeLoadArray`.
 --
--- @since 0.3.0
+-- @since 0.4.0
 makeLoadArray ::
-     Comp
+     Index ix
+  => Comp
+  -- ^ Computation strategy to use. Directly affects the scheduler that gets created for
+  -- the loading function.
   -> Sz ix
-  -> (forall m. Monad m => Scheduler m () -> Int -> (Int -> e -> m ()) -> m ())
+  -- ^ Size of the resulting array
+  -> e
+  -- ^ Default value to use for all cells that might have been ommitted by the writing function
+  -> (forall m. Monad m =>
+                  Scheduler m () -> (ix -> e -> m Bool) -> m ())
+  -- ^ Writing function that described which elements to write into the target array. It
+  -- accepts a scheduler, that can be used for parallelization, as well as a safe element
+  -- writing function.
   -> Array DL ix e
-makeLoadArray comp sz = DLArray comp sz Nothing
+makeLoadArray comp sz defVal writer =
+  DLArray comp sz (Just defVal) $ \scheduler !startAt uWrite ->
+    let safeWrite !ix !e
+          | isSafeIndex sz ix = uWrite (startAt + toLinearIndex sz ix) e >> pure True
+          | otherwise = pure False
+        {-# INLINE safeWrite #-}
+     in writer scheduler safeWrite
 {-# INLINE makeLoadArray #-}
-{-# DEPRECATED makeLoadArray "In favor of equivalent `unsafeMakeLoadArray` and safe `makeLoadArrayS`" #-}
 
 -- | Specify how an array can be loaded/computed through creation of a `DL` array. Unlike
--- `makeLoadArrayS` this function is unsafe since there is no guarantee that all elements will be
--- initialized and in case of parallel scheduler there is a possibility of non-determinism.
+-- `makeLoadArrayS` or `makeLoadArray` this function is unsafe, since there is no
+-- guarantee that all elements will be initialized and the supplied element writing
+-- function does not perform any bounds checking.
 --
 -- @since 0.3.1
 unsafeMakeLoadArray ::
      Comp
+  -- ^ Computation strategy to use. Directly affects the scheduler that gets created for
+  -- the loading function.
   -> Sz ix
+  -- ^ Size of the array
   -> Maybe e
-  -> (forall m. Monad m => Scheduler m () -> Int -> (Int -> e -> m ()) -> m ())
+  -- ^ An element to use for initialization of the mutable array that will be created in
+  -- the future
+  -> (forall m. Monad m =>
+                  Scheduler m () -> Int -> (Int -> e -> m ()) -> m ())
+  -- ^ This function accepts:
+  --
+  -- * A scheduler that can be used for parallelization of loading
+  --
+  -- * Linear index at which this load array will start (an offset that should be added to
+  --   the linear writng function)
+  --
+  -- * Linear element writing function
   -> Array DL ix e
 unsafeMakeLoadArray = DLArray
 {-# INLINE unsafeMakeLoadArray #-}
diff --git a/src/Data/Massiv/Array/Delayed/Windowed.hs b/src/Data/Massiv/Array/Delayed/Windowed.hs
--- a/src/Data/Massiv/Array/Delayed/Windowed.hs
+++ b/src/Data/Massiv/Array/Delayed/Windowed.hs
@@ -42,8 +42,6 @@
 -- | Delayed Windowed Array representation.
 data DW = DW
 
-type instance EltRepr DW ix = D
-
 data Window ix e = Window { windowStart     :: !ix
                           -- ^ Index of where window will start at.
                           , windowSize      :: !(Sz ix)
diff --git a/src/Data/Massiv/Array/Manifest/Boxed.hs b/src/Data/Massiv/Array/Manifest/Boxed.hs
--- a/src/Data/Massiv/Array/Manifest/Boxed.hs
+++ b/src/Data/Massiv/Array/Manifest/Boxed.hs
@@ -78,8 +78,6 @@
 -- spine strict, but elements are strict to Weak Head Normal Form (WHNF) only.
 data B = B deriving Show
 
-type instance EltRepr B ix = M
-
 data instance Array B ix e = BArray { bComp :: !Comp
                                     , bSize :: !(Sz ix)
                                     , bData :: {-# UNPACK #-} !(A.Array e)
@@ -185,6 +183,7 @@
   {-# INLINE unsafeLinearWrite #-}
 
 instance Index ix => Load B ix e where
+  type R B = M
   size = bSize
   {-# INLINE size #-}
   getComp = bComp
@@ -246,8 +245,6 @@
 -- instance is required.
 data N = N deriving Show
 
-type instance EltRepr N ix = M
-
 newtype instance Array N ix e = NArray { bArray :: Array B ix e }
 
 instance (Ragged L ix e, Show e, NFData e) => Show (Array N ix e) where
@@ -354,6 +351,7 @@
   {-# INLINE unsafeLinearWrite #-}
 
 instance (Index ix, NFData e) => Load N ix e where
+  type R N = M
   size = bSize . bArray
   {-# INLINE size #-}
   getComp = bComp . bArray
diff --git a/src/Data/Massiv/Array/Manifest/Internal.hs b/src/Data/Massiv/Array/Manifest/Internal.hs
--- a/src/Data/Massiv/Array/Manifest/Internal.hs
+++ b/src/Data/Massiv/Array/Manifest/Internal.hs
@@ -35,7 +35,6 @@
   , gcastArr
   , fromRaggedArrayM
   , fromRaggedArray'
-  , fromRaggedArray
   , sizeofArray
   , sizeofMutableArray
   , iterateUntil
@@ -53,8 +52,6 @@
 import Data.Massiv.Core.List
 import Data.Maybe (fromMaybe)
 import Data.Typeable
-import qualified Data.Vector as V
-import qualified Data.Vector.Mutable as MV
 import GHC.Base hiding (ord)
 import System.IO.Unsafe (unsafePerformIO)
 
@@ -81,7 +78,6 @@
 data instance Array M ix e = MArray { mComp :: !Comp
                                     , mSize :: !(Sz ix)
                                     , mLinearIndex :: Int -> e }
-type instance EltRepr M ix = M
 
 instance (Ragged L ix e, Show e) => Show (Array M ix e) where
   showsPrec = showsArrayPrec id
@@ -96,20 +92,7 @@
   compare = ord compare
   {-# INLINE compare #-}
 
-instance Index ix => Construct M ix e where
-  setComp c arr = arr {mComp = c}
-  {-# INLINE setComp #-}
-  makeArrayLinear !comp !sz f =
-    unsafePerformIO $ do
-      let !k = totalElem sz
-      mv <- MV.unsafeNew k
-      withScheduler_ comp $ \scheduler ->
-        splitLinearlyWithM_ scheduler k (pure . f) (MV.unsafeWrite mv)
-      v <- V.unsafeFreeze mv
-      pure $ MArray comp sz (V.unsafeIndex v)
-  {-# INLINE makeArrayLinear #-}
 
-
 -- | /O(1)/ - Conversion of `Manifest` arrays to `M` representation.
 toManifest :: Manifest r ix e => Array r ix e -> Array M ix e
 toManifest !arr = MArray (getComp arr) (size arr) (unsafeLinearIndexM arr)
@@ -307,25 +290,10 @@
 {-# INLINE convertProxy #-}
 
 
--- | Convert a ragged array into a usual rectangular shaped one.
-fromRaggedArray :: (Mutable r ix e, Ragged r' ix e, Load r' ix e) =>
-                   Array r' ix e -> Either ShapeException (Array r ix e)
-fromRaggedArray arr =
-  unsafePerformIO $ do
-    let !sz = edgeSize arr
-        !comp = getComp arr
-    mArr <- unsafeNew sz
-    try $ do
-      withScheduler_ comp $ \scheduler ->
-        loadRagged (scheduleWork scheduler) (unsafeLinearWrite mArr) 0 (totalElem sz) sz arr
-      unsafeFreeze comp mArr
-{-# INLINE fromRaggedArray #-}
-{-# DEPRECATED fromRaggedArray "In favor of a more general `fromRaggedArrayM`" #-}
-
 -- | Convert a ragged array into a common array with rectangular shape. Throws `ShapeException`
 -- whenever supplied ragged array does not have a rectangular shape.
 --
--- @since 0.3.0
+-- @since 0.4.0
 fromRaggedArrayM ::
      forall r ix e r' m . (Mutable r ix e, Ragged r' ix e, Load r' ix e, MonadThrow m)
   => Array r' ix e
@@ -341,8 +309,10 @@
 {-# INLINE fromRaggedArrayM #-}
 
 
--- | Same as `fromRaggedArray`, but will throw an error if its shape is not
+-- | Same as `fromRaggedArrayM`, but will throw a pure exception if its shape is not
 -- rectangular.
+--
+-- @since 0.1.1
 fromRaggedArray' ::
      forall r ix e r'. (Mutable r ix e, Load r' ix e, Ragged r' ix e)
   => Array r' ix e
diff --git a/src/Data/Massiv/Array/Manifest/List.hs b/src/Data/Massiv/Array/Manifest/List.hs
--- a/src/Data/Massiv/Array/Manifest/List.hs
+++ b/src/Data/Massiv/Array/Manifest/List.hs
@@ -19,7 +19,6 @@
     fromList
   , fromListsM
   , fromLists'
-  , fromLists
   , toList
   , toLists
   , toLists2
@@ -77,6 +76,8 @@
 --   , [ [4,5] ]
 --   ]
 -- )
+-- >>> fromListsM Seq [[[1,2,3]],[[4,5]]] :: Maybe (Array B Ix3 Int)
+-- Nothing
 -- >>> fromListsM Seq [[[1,2,3]],[[4,5]]] :: IO (Array B Ix3 Int)
 -- *** Exception: DimTooShortException: expected (Sz1 3), got (Sz1 2)
 --
@@ -86,19 +87,12 @@
 fromListsM comp = fromRaggedArrayM . setComp comp . throughNested
 {-# INLINE fromListsM #-}
 
--- | Similar to `fromListsM`, but less general.
-fromLists :: (Nested LN ix e, Ragged L ix e, Mutable r ix e)
-         => Comp -> [ListItem ix e] -> Maybe (Array r ix e)
-fromLists comp = fromRaggedArrayM . setComp comp . throughNested
-{-# INLINE fromLists #-}
-{-# DEPRECATED fromLists "In favor of a more general `fromListsM`" #-}
-
 -- TODO: Figure out QuickCheck properties. Best guess idea so far IMHO is to add it as dependency
 -- and move Arbitrary instances int the library
 --
 -- prop> fromLists' Seq xs == fromList xs
 --
--- | Same as `fromLists`, but will throw an error on irregular shaped lists.
+-- | Same as `fromListsM`, but will throw a pure error on irregular shaped lists.
 --
 -- __Note__: This function is the same as if you would turn on @{-\# LANGUAGE OverloadedLists #-}@
 -- extension. For that reason you can also use `GHC.Exts.fromList`.
diff --git a/src/Data/Massiv/Array/Manifest/Primitive.hs b/src/Data/Massiv/Array/Manifest/Primitive.hs
--- a/src/Data/Massiv/Array/Manifest/Primitive.hs
+++ b/src/Data/Massiv/Array/Manifest/Primitive.hs
@@ -52,9 +52,6 @@
 import Data.Primitive.Types
 import GHC.Base (Int(..))
 import GHC.Exts as GHC
-#if !MIN_VERSION_base(4,11,0)
-import GHC.Prim as GHC
-#endif
 import Prelude hiding (mapM)
 import System.IO.Unsafe (unsafePerformIO)
 
@@ -63,8 +60,6 @@
 -- | Representation for `Prim`itive elements
 data P = P deriving Show
 
-type instance EltRepr P ix = M
-
 data instance Array P ix e = PArray { pComp :: !Comp
                                     , pSize :: !(Sz ix)
                                     , pData :: {-# UNPACK #-} !ByteArray
@@ -194,7 +189,7 @@
                 Sz . elemsMBA _mpa, writeByteArray) ma
   {-# INLINE unsafeLinearWrite #-}
 
-  unsafeLinearSet (MPArray _ ma) = setByteArray ma
+  unsafeLinearSet (MPArray _ ma) offset (SafeSz sz) = setByteArray ma offset sz
   {-# INLINE unsafeLinearSet #-}
 
   unsafeLinearCopy (MPArray _ maFrom) iFrom (MPArray _ maTo) iTo (Sz k) =
@@ -213,11 +208,12 @@
   {-# INLINE unsafeLinearShrink #-}
 
   unsafeLinearGrow (MPArray _ ma) sz =
-    MPArray sz <$> resizeMutableByteArray ma (totalElem sz * sizeOf (undefined :: e))
+    MPArray sz <$> resizeMutableByteArrayCompat ma (totalElem sz * sizeOf (undefined :: e))
   {-# INLINE unsafeLinearGrow #-}
 
 
 instance (Prim e, Index ix) => Load P ix e where
+  type R P = M
   size = pSize
   {-# INLINE size #-}
   getComp = pComp
@@ -271,13 +267,12 @@
     arr = PArray comp sz ba
 {-# INLINE fromByteArrayM #-}
 
--- | See `fromByteArrayM`.
+-- | /O(1)/ - Construct a flat Array from `ByteArray`
 --
--- @since 0.2.1
-fromByteArray :: (Index ix, Prim e) => Comp -> Sz ix -> ByteArray -> Maybe (Array P ix e)
-fromByteArray = fromByteArrayM
+-- @since 0.4.0
+fromByteArray :: forall e . Prim e => Comp -> ByteArray -> Array P Ix1 e
+fromByteArray comp ba = PArray comp (SafeSz (elemsBA (Proxy :: Proxy e) ba)) ba
 {-# INLINE fromByteArray #-}
-{-# DEPRECATED fromByteArray "In favor of more general `fromByteArrayM`" #-}
 
 
 -- | /O(1)/ - Extract the internal `MutableByteArray`.
@@ -288,8 +283,8 @@
 {-# INLINE toMutableByteArray #-}
 
 
--- | /O(1)/ - Construct a primitive mutable array from the `MutableByteArray`. Will return `Nothing`
--- if number of elements doesn't match.
+-- | /O(1)/ - Construct a primitive mutable array from the `MutableByteArray`. Will throw
+-- `SizeElementsMismatchException` if number of elements doesn't match.
 --
 -- @since 0.3.0
 fromMutableByteArrayM ::
@@ -300,13 +295,12 @@
     marr = MPArray sz mba
 {-# INLINE fromMutableByteArrayM #-}
 
--- | See `fromMutableByteArray`.
+-- | /O(1)/ - Construct a flat Array from `MutableByteArray`
 --
--- @since 0.2.1
-fromMutableByteArray :: (Index ix, Prim e) => Sz ix -> MutableByteArray s -> Maybe (MArray s P ix e)
-fromMutableByteArray = fromMutableByteArrayM
+-- @since 0.4.0
+fromMutableByteArray :: forall e s . Prim e => MutableByteArray s -> MArray s P Ix1 e
+fromMutableByteArray mba = MPArray (SafeSz (elemsMBA (Proxy :: Proxy e) mba)) mba
 {-# INLINE fromMutableByteArray #-}
-{-# DEPRECATED fromMutableByteArray "In favor of more general `fromMutableByteArrayM`" #-}
 
 
 -- | Atomically read an `Int` element from the array
@@ -483,22 +477,24 @@
 {-# INLINE unsafeAtomicXorIntArray #-}
 
 
-
-#if !MIN_VERSION_primitive(0,6,4)
-resizeMutableByteArray ::
-     PrimMonad m => MutableByteArray (PrimState m) -> Int -> m (MutableByteArray (PrimState m))
-resizeMutableByteArray (MutableByteArray arr#) (I# n#) =
-  primitive
-    (\s# ->
-       case resizeMutableByteArray# arr# n# s# of
-         (# s'#, arr'# #) -> (# s'#, MutableByteArray arr'# #))
-{-# INLINE resizeMutableByteArray #-}
-#endif
-
 shrinkMutableByteArray :: forall m. (PrimMonad m)
   => MutableByteArray (PrimState m)
   -> Int -- ^ new size
   -> m ()
-{-# INLINE shrinkMutableByteArray #-}
 shrinkMutableByteArray (MutableByteArray arr#) (I# n#)
   = primitive_ (shrinkMutableByteArray# arr# n#)
+{-# INLINE shrinkMutableByteArray #-}
+
+
+resizeMutableByteArrayCompat ::
+  PrimMonad m => MutableByteArray (PrimState m) -> Int -> m (MutableByteArray (PrimState m))
+#if MIN_VERSION_primitive(0,6,4)
+resizeMutableByteArrayCompat = resizeMutableByteArray
+#else
+resizeMutableByteArrayCompat (MutableByteArray arr#) (I# n#) =
+  primitive
+    (\s# ->
+       case resizeMutableByteArray# arr# n# s# of
+         (# s'#, arr'# #) -> (# s'#, MutableByteArray arr'# #))
+#endif
+{-# INLINE resizeMutableByteArrayCompat #-}
diff --git a/src/Data/Massiv/Array/Manifest/Storable.hs b/src/Data/Massiv/Array/Manifest/Storable.hs
--- a/src/Data/Massiv/Array/Manifest/Storable.hs
+++ b/src/Data/Massiv/Array/Manifest/Storable.hs
@@ -59,8 +59,6 @@
 -- | Representation for `Storable` elements
 data S = S deriving Show
 
-type instance EltRepr S ix = M
-
 data instance Array S ix e = SArray { sComp :: !Comp
                                     , sSize :: !(Sz ix)
                                     , sData :: !(VS.Vector e)
@@ -202,6 +200,7 @@
 
 
 instance (Index ix, VS.Storable e) => Load S ix e where
+  type R S = M
   size = sSize
   {-# INLINE size #-}
   getComp = sComp
diff --git a/src/Data/Massiv/Array/Manifest/Unboxed.hs b/src/Data/Massiv/Array/Manifest/Unboxed.hs
--- a/src/Data/Massiv/Array/Manifest/Unboxed.hs
+++ b/src/Data/Massiv/Array/Manifest/Unboxed.hs
@@ -42,8 +42,6 @@
 -- | Representation for `Unbox`ed elements
 data U = U deriving Show
 
-type instance EltRepr U ix = M
-
 data instance Array U ix e = UArray { uComp :: !Comp
                                     , uSize :: !(Sz ix)
                                     , uData :: !(VU.Vector e)
@@ -90,6 +88,7 @@
   {-# INLINE unsafeExtract #-}
 
 instance (VU.Unbox e, Index ix) => Load U ix e where
+  type R U = M
   size = uSize
   {-# INLINE size #-}
   getComp = uComp
diff --git a/src/Data/Massiv/Array/Manifest/Vector.hs b/src/Data/Massiv/Array/Manifest/Vector.hs
--- a/src/Data/Massiv/Array/Manifest/Vector.hs
+++ b/src/Data/Massiv/Array/Manifest/Vector.hs
@@ -15,7 +15,6 @@
 module Data.Massiv.Array.Manifest.Vector
   ( fromVectorM
   , fromVector'
-  , fromVector
   , castFromVector
   , toVector
   , castToVector
@@ -95,7 +94,6 @@
      , Typeable v
      , VG.Vector v a
      , Mutable (ARepr v) ix a
-     , Construct r ix a
      , Mutable r ix a
      )
   => Comp
@@ -115,23 +113,13 @@
 --
 -- @since 0.3.0
 fromVector' ::
-     (Typeable v, VG.Vector v a, Mutable (ARepr v) ix a, Construct r ix a, Mutable r ix a)
+     (Typeable v, VG.Vector v a, Mutable (ARepr v) ix a, Mutable r ix a)
   => Comp
   -> Sz ix -- ^ Resulting size of the array
   -> v a -- ^ Source Vector
   -> Array r ix a
 fromVector' comp sz = either throw id . fromVectorM comp sz
 {-# INLINE fromVector' #-}
-
-fromVector ::
-     (Typeable v, VG.Vector v a, Mutable (ARepr v) ix a, Construct r ix a, Mutable r ix a)
-  => Comp
-  -> Sz ix -- ^ Resulting size of the array
-  -> v a -- ^ Source Vector
-  -> Array r ix a
-fromVector comp sz = either throw id . fromVectorM comp sz
-{-# INLINE fromVector #-}
-{-# DEPRECATED fromVector "In favor of safer `fromVectorM`" #-}
 
 -- | /O(1)/ - conversion from `Mutable` array to a corresponding vector. Will
 -- return `Nothing` only if source array representation was not one of `B`, `N`,
diff --git a/src/Data/Massiv/Array/Mutable.hs b/src/Data/Massiv/Array/Mutable.hs
--- a/src/Data/Massiv/Array/Mutable.hs
+++ b/src/Data/Massiv/Array/Mutable.hs
@@ -17,12 +17,18 @@
     msize
     -- ** Element-wise mutation
   , read
+  , readM
   , read'
   , write
+  , writeM
   , write'
   , modify
+  , modifyM
+  , modifyM_
   , modify'
   , swap
+  , swapM
+  , swapM_
   , swap'
   -- ** Operations on @MArray@
   -- *** Immutable conversion
@@ -61,8 +67,11 @@
   , unfoldlPrimM
   , iunfoldlPrimM
   -- *** Mapping
+  , forPrimM
   , forPrimM_
+  , iforPrimM
   , iforPrimM_
+  , iforLinearPrimM
   , iforLinearPrimM_
   -- *** Modify
   , withMArray
@@ -82,7 +91,7 @@
 
 -- TODO: add fromListM, et al.
 
-import Control.Monad (when, unless)
+import Control.Monad (void, when, unless, (>=>))
 import Control.Monad.ST
 import Control.Scheduler
 import Data.Massiv.Core.Common
@@ -110,7 +119,7 @@
 --   [ [ 0, 0, 0, 0, 0, 0 ]
 --   , [ 0, 0, 0, 0, 0, 0 ]
 --   ]
--- >>> new @B @_ @Int (Sz2 2 6) >>= (`read'` 1)
+-- >>> new @B @_ @Int (Sz2 2 6) >>= (`readM` 1)
 -- *** Exception: Uninitialized
 --
 -- @since 0.1.0
@@ -130,8 +139,8 @@
 -- >>> :set -XTypeApplications
 -- >>> arr <- fromListsM @U @Ix2 @Double Par [[12,21],[13,31]]
 -- >>> marr <- thaw arr
--- >>> modify marr (+ 10) (1 :. 0)
--- True
+-- >>> modify marr (pure . (+ 10)) (1 :. 0)
+-- Just 13.0
 -- >>> freeze Par marr
 -- Array U Par (Sz (2 :. 2))
 --   [ [ 12.0, 21.0 ]
@@ -164,7 +173,7 @@
 -- >>> :set -XOverloadedLists
 -- >>> thawS @P @Ix1 @Double [1..10]
 -- >>> marr <- thawS @P @Ix1 @Double [1..10]
--- >>> write' marr 5 100
+-- >>> writeM marr 5 100
 -- >>> freezeS marr
 -- Array P Seq (Sz1 10)
 --   [ 1.0, 2.0, 3.0, 4.0, 5.0, 100.0, 7.0, 8.0, 9.0, 10.0 ]
@@ -236,11 +245,7 @@
 
 newMaybeInitialized ::
      (Load r' ix e, Mutable r ix e, PrimMonad m) => Array r' ix e -> m (MArray (PrimState m) r ix e)
-newMaybeInitialized !arr = do
-  let !sz = size arr
-  marr <- unsafeNew sz
-  mapM_ (unsafeLinearSet marr 0 (totalElem sz)) $ defaultElement arr
-  pure marr
+newMaybeInitialized !arr = initializeNew (defaultElement arr) (size arr)
 {-# INLINE newMaybeInitialized #-}
 
 
@@ -356,7 +361,7 @@
 --
 -- >>> :set -XTypeApplications
 -- >>> import Data.Massiv.Array
--- >>> createArray_ @P @_ @Int Seq (Sz1 2) (\ s marr -> scheduleWork s (write' marr 0 10) >> scheduleWork s (write' marr 1 11))
+-- >>> createArray_ @P @_ @Int Seq (Sz1 2) (\ s marr -> scheduleWork s (writeM marr 0 10) >> scheduleWork s (writeM marr 1 11))
 -- Array P Seq (Sz1 2)
 --   [ 10, 11 ]
 --
@@ -464,11 +469,11 @@
 {-# INLINE createArrayST #-}
 
 
--- | Sequentially generate a pure array. Much like `makeArray` creates a pure array this function
--- will use `Mutable` interface to generate a pure `Array` in the end, except that computation
--- strategy is ignored. Element producing function no longer has to be pure but is a stateful
--- action, since it is restricted to `PrimMonad` and allows for sharing the state between
--- computation of each element, which could be arbitrary effects if that monad is `IO`.
+-- | Sequentially generate a pure array. Much like `makeArray` creates a pure array this
+-- function will use `Mutable` interface to generate a pure `Array` in the end, except that
+-- computation strategy is set to `Seq`. Element producing function no longer has to be pure
+-- but is a stateful action, becuase it is restricted to `PrimMonad` thus allows for sharing
+-- the state between computation of each element.
 --
 -- @since 0.2.6
 --
@@ -477,7 +482,7 @@
 -- >>> import Data.Massiv.Array
 -- >>> import Data.IORef
 -- >>> ref <- newIORef (0 :: Int)
--- >>> generateArray Seq (Sz1 6) (\ i -> modifyIORef' ref (+i) >> print i >> pure i) :: IO (Array U Ix1 Int)
+-- >>> generateArrayS (Sz1 6) (\ i -> modifyIORef' ref (+i) >> print i >> pure i) :: IO (Array U Ix1 Int)
 -- 0
 -- 1
 -- 2
@@ -491,26 +496,24 @@
 --
 generateArrayS ::
      forall r ix e m. (Mutable r ix e, PrimMonad m)
-  => Comp -- ^ Computation strategy (ingored during generation)
-  -> Sz ix -- ^ Resulting size of the array
+  => Sz ix -- ^ Resulting size of the array
   -> (ix -> m e) -- ^ Element producing generator
   -> m (Array r ix e)
-generateArrayS comp sz gen = generateArrayLinearS comp sz (gen . fromLinearIndex sz)
+generateArrayS sz gen = generateArrayLinearS sz (gen . fromLinearIndex sz)
 {-# INLINE generateArrayS #-}
 
--- | Same as `generateArray` but with action takes row-major linear index.
+-- | Same as `generateArray` but with action that accepts row-major linear index.
 --
 -- @since 0.3.0
 generateArrayLinearS ::
      forall r ix e m. (Mutable r ix e, PrimMonad m)
-  => Comp -- ^ Computation strategy (ingored during generation)
-  -> Sz ix -- ^ Resulting size of the array
+  => Sz ix -- ^ Resulting size of the array
   -> (Int -> m e) -- ^ Element producing generator
   -> m (Array r ix e)
-generateArrayLinearS comp sz gen = do
+generateArrayLinearS sz gen = do
   marr <- unsafeNew sz
   loopM_ 0 (< totalElem (msize marr)) (+ 1) $ \i -> gen i >>= unsafeLinearWrite marr i
-  unsafeFreeze comp marr
+  unsafeFreeze Seq marr
 {-# INLINE generateArrayLinearS #-}
 
 
@@ -640,7 +643,7 @@
      in iterLinearM sz' 0 (totalElem sz') 1 (<) acc0 $ \i ix acc -> do
           (e, acc') <- gen acc ix
           unsafeLinearWrite marr i e
-          return acc'
+          pure $! acc'
 {-# INLINE iunfoldrPrimM #-}
 
 -- | Just like `iunfoldrPrimM`, but do the unfolding with index aware function.
@@ -659,7 +662,7 @@
      in loopM 0 (< totalElem sz') (+1) acc0 $ \i acc -> do
           (e, acc') <- gen acc
           unsafeLinearWrite marr i e
-          return acc'
+          pure $! acc'
 {-# INLINE unfoldrPrimM #-}
 
 -- | Sequentially unfold an array from the left.
@@ -727,7 +730,7 @@
      in iterLinearM sz' (totalElem sz' - 1) 0 (negate 1) (>=) acc0 $ \i ix acc -> do
           (acc', e) <- gen acc ix
           unsafeLinearWrite marr i e
-          return acc'
+          pure $! acc'
 {-# INLINE iunfoldlPrimM #-}
 
 -- | Just like `iunfoldlPrimM`, but do the unfolding with index aware function.
@@ -746,39 +749,66 @@
      in loopDeepM 0 (< totalElem sz') (+1) acc0 $ \i acc -> do
           (acc', e) <- gen acc
           unsafeLinearWrite marr i e
-          return acc'
+          pure $! acc'
 {-# INLINE unfoldlPrimM #-}
 
---TODO: in 0.4.0:
---  * forPrimM_  does not modify each element, but simply iterates over
---  * forPrimM - does what forPrimM_ does now.
--- | Sequentially loop over a mutable array while modifying each element with an action.
+-- | Sequentially loop over a mutable array while reading each element and applying an
+-- action to it. There is no mutation to the array, unless the action itself modifies it.
 --
--- @since 0.3.0
-forPrimM_ :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (e -> m e) -> m ()
+-- @since 0.4.0
+forPrimM_ :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (e -> m ()) -> m ()
 forPrimM_ marr f =
-  loopM_ 0 (< totalElem (msize marr)) (+1) (unsafeLinearModify marr (const f))
+  loopM_ 0 (< totalElem (msize marr)) (+1) (unsafeLinearRead marr >=> f)
 {-# INLINE forPrimM_ #-}
 
+-- | Sequentially loop over a mutable array while modifying each element with an action.
+--
+-- @since 0.4.0
+forPrimM :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (e -> m e) -> m ()
+forPrimM marr f =
+  loopM_ 0 (< totalElem (msize marr)) (+1) (unsafeLinearModify marr f)
+{-# INLINE forPrimM #-}
 
--- | Sequentially loop over a mutable array while modifying each element with an index aware action.
+
+-- | Sequentially loop over a mutable array while reading each element and applying an
+-- index aware action to it. There is no mutation to the array, unless the
+-- action itself modifies it.
 --
--- @since 0.3.0
+-- @since 0.4.0
 iforPrimM_ ::
-     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (ix -> e -> m e) -> m ()
+     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (ix -> e -> m ()) -> m ()
 iforPrimM_ marr f = iforLinearPrimM_ marr (f . fromLinearIndex (msize marr))
 {-# INLINE iforPrimM_ #-}
 
-
 -- | Sequentially loop over a mutable array while modifying each element with an index aware action.
 --
--- @since 0.3.0
+-- @since 0.4.0
+iforPrimM ::
+     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (ix -> e -> m e) -> m ()
+iforPrimM marr f = iforLinearPrimM marr (f . fromLinearIndex (msize marr))
+{-# INLINE iforPrimM #-}
+
+
+-- | Sequentially loop over a mutable array while reading each element and applying a
+-- linear index aware action to it. There is no mutation to the array, unless the action
+-- itself modifies it.
+--
+-- @since 0.4.0
 iforLinearPrimM_ ::
-     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (Int -> e -> m e) -> m ()
-iforLinearPrimM_ marr f = loopM_ 0 (< totalElem (msize marr)) (+ 1) (unsafeLinearModify marr f)
+     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (Int -> e -> m ()) -> m ()
+iforLinearPrimM_ marr f =
+  loopM_ 0 (< totalElem (msize marr)) (+ 1) (\i -> unsafeLinearRead marr i >>= f i)
 {-# INLINE iforLinearPrimM_ #-}
 
--- FIXME: supply Scheduler instead of numWorkers and schedule action
+-- | Sequentially loop over a mutable array while modifying each element with an index aware action.
+--
+-- @since 0.4.0
+iforLinearPrimM ::
+     (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (Int -> e -> m e) -> m ()
+iforLinearPrimM marr f =
+  loopM_ 0 (< totalElem (msize marr)) (+ 1) (\i -> unsafeLinearModify marr (f i) i)
+{-# INLINE iforLinearPrimM #-}
+
 -- | Create a copy of a pure array, mutate it in place and return its frozen version. The big
 -- difference between `withMArrayS` is that it's not only gonna respect the computation strategy
 -- supplied to it while making a copy, but it will also pass extra argumens to the action that
@@ -795,12 +825,11 @@
 withMArray ::
      (Mutable r ix e, MonadUnliftIO m)
   => Array r ix e
-  -> (Int -> (m () -> m ()) -> MArray RealWorld r ix e -> m a)
+  -> (Scheduler m () -> MArray RealWorld r ix e -> m a)
   -> m (Array r ix e)
 withMArray arr action = do
   marr <- thaw arr
-  withScheduler_ (getComp arr) $ \scheduler ->
-    action (numWorkers scheduler) (scheduleWork scheduler) marr
+  withScheduler_ (getComp arr) (`action` marr)
   liftIO $ unsafeFreeze (getComp arr) marr
 {-# INLINE withMArray #-}
 
@@ -849,6 +878,18 @@
 
 -- | /O(1)/ - Same as `read`, but throws `IndexOutOfBoundsException` on an invalid index.
 --
+-- @since 0.4.0
+readM :: (Mutable r ix e, PrimMonad m, MonadThrow m) =>
+        MArray (PrimState m) r ix e -> ix -> m e
+readM marr ix =
+  read marr ix >>= \case
+    Just e -> pure e
+    Nothing -> throwM $ IndexOutOfBoundsException (msize marr) ix
+{-# INLINE readM #-}
+
+
+-- | /O(1)/ - Same as `read`, but throws `IndexOutOfBoundsException` on an invalid index.
+--
 -- @since 0.1.0
 read' :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> m e
 read' marr ix =
@@ -856,10 +897,13 @@
     Just e -> pure e
     Nothing -> throw $ IndexOutOfBoundsException (msize marr) ix
 {-# INLINE read' #-}
+{-# DEPRECATED read' "In favor of more general `readM`" #-}
 
 
 -- | /O(1)/ - Write an element into the cell of a mutable array. Returns `False` when index is out
 -- of bounds.
+--
+-- @since 0.1.0
 write :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> e -> m Bool
 write marr ix e =
   if isSafeIndex (msize marr) ix
@@ -867,7 +911,16 @@
   else pure False
 {-# INLINE write #-}
 
+-- | /O(1)/ - Same as `write`, but throws `IndexOutOfBoundsException` on an invalid index.
+--
+-- @since 0.4.0
+writeM ::
+     (Mutable r ix e, PrimMonad m, MonadThrow m) => MArray (PrimState m) r ix e -> ix -> e -> m ()
+writeM marr ix e =
+  write marr ix e >>= (`unless` throwM (IndexOutOfBoundsException (msize marr) ix))
+{-# INLINE writeM #-}
 
+
 -- | /O(1)/ - Same as `write`, but lives in IO and throws `IndexOutOfBoundsException` on invalid
 -- index.
 --
@@ -876,58 +929,130 @@
      (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> e -> m ()
 write' marr ix e = write marr ix e >>= (`unless` throw (IndexOutOfBoundsException (msize marr) ix))
 {-# INLINE write' #-}
+{-# DEPRECATED write' "In favor of more general `writeM`" #-}
 
--- TODO: switch to `... -> m (Maybe e)`
--- | /O(1)/ - Modify an element in the cell of a mutable array with a supplied function. Returns
--- `False` when index is out of bounds.
+-- | /O(1)/ - Modify an element in the cell of a mutable array with a supplied
+-- action. Returns the previous value, if index was not out of bounds.
 --
 -- @since 0.1.0
-modify :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> (e -> e) -> ix -> m Bool
+modify ::
+     (Mutable r ix e, PrimMonad m)
+  => MArray (PrimState m) r ix e -- ^ Array to mutate.
+  -> (e -> m e) -- ^ Monadic action that modifies the element
+  -> ix -- ^ Index at which to perform modification.
+  -> m (Maybe e)
 modify marr f ix =
   if isSafeIndex (msize marr) ix
-  then do
-    val <- unsafeRead marr ix
-    unsafeWrite marr ix $ f val
-    return True
-  else return False
+    then Just <$> unsafeModify marr f ix
+    else return Nothing
 {-# INLINE modify #-}
 
+-- | /O(1)/ - Modify an element in the cell of a mutable array with a supplied
+-- action. Throws an `IndexOutOfBoundsException` exception for invalid index and returns
+-- the previous value otherwise.
+--
+-- @since 0.4.0
+modifyM ::
+     (Mutable r ix e, PrimMonad m, MonadThrow m)
+  => MArray (PrimState m) r ix e -- ^ Array to mutate.
+  -> (e -> m e) -- ^ Monadic action that modifies the element
+  -> ix -- ^ Index at which to perform modification.
+  -> m e
+modifyM marr f ix
+  | isSafeIndex (msize marr) ix = unsafeModify marr f ix
+  | otherwise = throwM (IndexOutOfBoundsException (msize marr) ix)
+{-# INLINE modifyM #-}
 
+-- | /O(1)/ - Same as `modifyM`, but discard the returned element
+--
+-- ====__Examples__
+--
+-- >>> :set -XTypeApplications
+-- >>> import Control.Monad.ST
+-- >>> import Data.Massiv.Array
+-- >>> runST $ new @P @Ix1 @Int (Sz1 3) >>= (\ma -> modifyM_ ma (pure . (+10)) 1 >> freezeS ma)
+-- Array P Seq (Sz1 3)
+--   [ 0, 10, 0 ]
+--
+-- @since 0.4.0
+modifyM_ ::
+     (Mutable r ix e, PrimMonad m, MonadThrow m)
+  => MArray (PrimState m) r ix e -- ^ Array to mutate.
+  -> (e -> m e) -- ^ Monadic action that modifies the element
+  -> ix -- ^ Index at which to perform modification.
+  -> m ()
+modifyM_ marr f ix = void $ modifyM marr f ix
+{-# INLINE modifyM_ #-}
+
+
 -- | /O(1)/ - Same as `modify`, but throws an error if index is out of bounds.
 --
 -- @since 0.1.0
 modify' :: (Mutable r ix e, PrimMonad m) =>
         MArray (PrimState m) r ix e -> (e -> e) -> ix -> m ()
 modify' marr f ix =
-  modify marr f ix >>= (`unless` throw (IndexOutOfBoundsException (msize marr) ix))
+  modify marr (pure . f) ix >>= \case
+    Just _ -> pure ()
+    Nothing -> throw (IndexOutOfBoundsException (msize marr) ix)
 {-# INLINE modify' #-}
+{-# DEPRECATED modify' "In favor of more general `modifyM`" #-}
 
 
--- | /O(1)/ - Swap two elements in a mutable array by supplying their indices. Returns `False` when
--- either one of the indices is out of bounds.
+-- | /O(1)/ - Same as `swapM`, but instead of thropwing an exception returns `Nothing` when
+-- either one of the indices is out of bounds and `Just` elements under those indices
+-- otherwise.
 --
 -- @since 0.1.0
-swap :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> ix -> m Bool
-swap marr ix1 ix2 = do
+swap :: (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> ix -> m (Maybe (e, e))
+swap marr ix1 ix2 =
   let sz = msize marr
-  if isSafeIndex sz ix1 && isSafeIndex sz ix2
-  then do
-    val <- unsafeRead marr ix1
-    unsafeRead marr ix2 >>= unsafeWrite marr ix1
-    unsafeWrite marr ix2 val
-    return True
-  else return False
+   in if isSafeIndex sz ix1 && isSafeIndex sz ix2
+        then Just <$> unsafeSwap marr ix1 ix2
+        else pure Nothing
 {-# INLINE swap #-}
 
+-- | /O(1)/ - Swap two elements in a mutable array under the supplied indices. Throws an
+-- `IndexOutOfBoundsException` when either one of the indices is out of bounds and
+-- elements under those indices otherwise.
+--
+-- @since 0.4.0
+swapM ::
+     (Mutable r ix e, PrimMonad m, MonadThrow m)
+  => MArray (PrimState m) r ix e
+  -> ix -- ^ Index for the first element, which will be returned as the first element in the
+        -- tuple.
+  -> ix -- ^ Index for the second element, which will be returned as the second element in
+        -- the tuple.
+  -> m (e, e)
+swapM marr ix1 ix2
+  | not (isSafeIndex sz ix1) = throwM $ IndexOutOfBoundsException (msize marr) ix1
+  | not (isSafeIndex sz ix2) = throwM $ IndexOutOfBoundsException (msize marr) ix2
+  | otherwise = unsafeSwap marr ix1 ix2
+  where
+    !sz = msize marr
+{-# INLINE swapM #-}
 
+
+-- | /O(1)/ - Same as `swapM`, but discard the returned elements
+--
+-- @since 0.4.0
+swapM_ ::
+     (Mutable r ix e, PrimMonad m, MonadThrow m) => MArray (PrimState m) r ix e -> ix -> ix -> m ()
+swapM_ marr ix1 ix2 = void $ swapM marr ix1 ix2
+{-# INLINE swapM_ #-}
+
+
 -- | /O(1)/ - Same as `swap`, but throws an `IndexOutOfBoundsException` on invalid indices.
 --
 -- @since 0.1.0
 swap' ::
      (Mutable r ix e, PrimMonad m) => MArray (PrimState m) r ix e -> ix -> ix -> m ()
 swap' marr ix1 ix2 =
-  swap marr ix1 ix2 >>=
-    (`unless` if isSafeIndex (msize marr) ix1
-                then throw $ IndexOutOfBoundsException (msize marr) ix2
-                else throw $ IndexOutOfBoundsException (msize marr) ix1)
+  swap marr ix1 ix2 >>= \case
+    Just _ -> pure ()
+    Nothing ->
+      if isSafeIndex (msize marr) ix1
+        then throw $ IndexOutOfBoundsException (msize marr) ix2
+        else throw $ IndexOutOfBoundsException (msize marr) ix1
 {-# INLINE swap' #-}
+{-# DEPRECATED swap' "In favor of more general `swapM`" #-}
diff --git a/src/Data/Massiv/Array/Mutable/Atomic.hs b/src/Data/Massiv/Array/Mutable/Atomic.hs
--- a/src/Data/Massiv/Array/Mutable/Atomic.hs
+++ b/src/Data/Massiv/Array/Mutable/Atomic.hs
@@ -71,7 +71,8 @@
 {-# INLINE casIntArray #-}
 
 
--- | Atomically modify an `Int` element of the array. Returns the old value.
+-- | Atomically modify an `Int` element of the array. Returns the old value, unless the
+-- supplied index was out of bounds.
 --
 -- @since 0.3.0
 atomicModifyIntArray ::
diff --git a/src/Data/Massiv/Array/Numeric.hs b/src/Data/Massiv/Array/Numeric.hs
--- a/src/Data/Massiv/Array/Numeric.hs
+++ b/src/Data/Massiv/Array/Numeric.hs
@@ -13,9 +13,15 @@
 --
 module Data.Massiv.Array.Numeric
   ( -- * Num
-    (.+)
+    (.+.)
+  , (.+)
+  , (+.)
+  , (.-.)
   , (.-)
+  , (-.)
+  , (.*.)
   , (.*)
+  , (*.)
   , (.^)
   , (|*|)
   , multiplyTransposed
@@ -32,6 +38,7 @@
   , quotRemA
   , divModA
   -- * Fractional
+  , (./.)
   , (./)
   , (.^^)
   , recipA
@@ -72,13 +79,14 @@
 import Data.Massiv.Array.Ops.Transform as A
 import Data.Massiv.Core
 import Data.Massiv.Core.Common
+import Data.Massiv.Core.Operations
 import Data.Massiv.Core.Index.Internal (Sz(SafeSz))
 import Prelude as P
 
 
 infixr 8  .^, .^^
-infixl 7  .*, ./, `quotA`, `remA`, `divA`, `modA`
-infixl 6  .+, .-
+infixl 7  .*., .*, *., ./., ./, `quotA`, `remA`, `divA`, `modA`
+infixl 6  .+., .+, +., .-., .-, -.
 
 liftArray2Matching
   :: (Source r1 ix a, Source r2 ix b)
@@ -95,38 +103,105 @@
     sz2 = size arr2
 {-# INLINE liftArray2Matching #-}
 
+liftArray2M ::
+     (Load r ix e, Numeric r e, MonadThrow m)
+  => (e -> e -> e)
+  -> Array r ix e
+  -> Array r ix e
+  -> m (Array r ix e)
+liftArray2M f a1 a2
+  | size a1 == size a2 = pure $ unsafeLiftArray2 f a1 a2
+  | otherwise = throwM $ SizeMismatchException (size a1) (size a2)
+{-# INLINE liftArray2M #-}
 
-(.+)
-  :: (Source r1 ix e, Source r2 ix e, Num e)
-  => Array r1 ix e -> Array r2 ix e -> Array D ix e
-(.+) = liftArray2Matching (+)
+
+liftNumericArray2M ::
+     (Load r ix e, MonadThrow m)
+  => (Array r ix e -> Array r ix e -> Array r ix e)
+  -> Array r ix e
+  -> Array r ix e
+  -> m (Array r ix e)
+liftNumericArray2M f a1 a2
+  | size a1 == size a2 = pure $ f a1 a2
+  | otherwise = throwM $ SizeMismatchException (size a1) (size a2)
+{-# INLINE liftNumericArray2M #-}
+
+
+-- | Add two arrays together pointwise. Throws `SizeMismatchException` if arrays sizes do
+-- not match.
+--
+-- @since 0.4.0
+(.+.) ::
+     (Load r ix e, Numeric r e, MonadThrow m) => Array r ix e -> Array r ix e -> m (Array r ix e)
+(.+.) = liftNumericArray2M additionPointwise
+{-# INLINE (.+.) #-}
+
+-- | Add a scalar to each element of the array. Array is on the left.
+--
+-- @since 0.1.0
+(.+) :: (Index ix, Numeric r e) => Array r ix e -> e -> Array r ix e
+(.+) = plusScalar
 {-# INLINE (.+) #-}
 
-(.-)
-  :: (Source r1 ix e, Source r2 ix e, Num e)
-  => Array r1 ix e -> Array r2 ix e -> Array D ix e
-(.-) = liftArray2Matching (-)
+-- | Add a scalar to each element of the array. Array is on the right.
+--
+-- @since 0.4.0
+(+.) :: (Index ix, Numeric r e) => e -> Array r ix e -> Array r ix e
+(+.) = flip plusScalar
+{-# INLINE (+.) #-}
+
+-- | Subtract two arrays pointwise. Throws `SizeMismatchException` if arrays sizes do not
+-- match.
+--
+-- @since 0.4.0
+(.-.) ::
+     (Load r ix e, Numeric r e, MonadThrow m) => Array r ix e -> Array r ix e -> m (Array r ix e)
+(.-.) = liftNumericArray2M subtractionPointwise
+{-# INLINE (.-.) #-}
+
+
+-- | Subtract a scalar from each element of the array. Array is on the left.
+--
+-- @since 0.1.0
+(.-) :: (Index ix, Numeric r e) => Array r ix e -> e -> Array r ix e
+(.-) = minusScalar
 {-# INLINE (.-) #-}
 
-(.*)
-  :: (Source r1 ix e, Source r2 ix e, Num e)
-  => Array r1 ix e -> Array r2 ix e -> Array D ix e
-(.*) = liftArray2Matching (*)
+-- | Subtract a scalar from each element of the array. Array is on the right.
+--
+-- @since 0.4.0
+(-.) :: (Index ix, Numeric r e) => e -> Array r ix e -> Array r ix e
+(-.) = flip minusScalar
+{-# INLINE (-.) #-}
+
+
+-- | Multiply two arrays together pointwise.
+--
+-- @since 0.4.0
+(.*.) ::
+     (Load r ix e, Numeric r e, MonadThrow m) => Array r ix e -> Array r ix e -> m (Array r ix e)
+(.*.) = liftNumericArray2M multiplicationPointwise
+{-# INLINE (.*.) #-}
+
+(.*) :: (Index ix, Numeric r e) => Array r ix e -> e -> Array r ix e
+(.*) = multiplyScalar
 {-# INLINE (.*) #-}
 
-(.^)
-  :: (Source r ix e, Num e, Integral b)
-  => Array r ix e -> b -> Array D ix e
-(.^) arr n = liftArray (^ n) arr
+(*.) :: (Index ix, Numeric r e) => e -> Array r ix e -> Array r ix e
+(*.) = flip multiplyScalar
+{-# INLINE (*.) #-}
+
+(.^) :: (Index ix, Numeric r e) => Array r ix e -> Int -> Array r ix e
+(.^) = powerPointwise
 {-# INLINE (.^) #-}
 
 -- | Perform matrix multiplication. Inner dimensions must agree, otherwise `SizeMismatchException`.
 (|*|) ::
-     (Mutable r Ix2 e, Source r' Ix2 e, OuterSlice r Ix2 e, Source (EltRepr r Ix2) Ix1 e, Num e)
+     (Mutable r Ix2 e, Source r' Ix2 e, OuterSlice r Ix2 e, Source (R r) Ix1 e, Num e, MonadThrow m)
   => Array r Ix2 e
   -> Array r' Ix2 e
-  -> Array r Ix2 e
-(|*|) a1 = compute . multArrs a1
+  -> m (Array r Ix2 e)
+(|*|) a1 a2 = compute <$> multArrs a1 a2
 {-# INLINE [1] (|*|) #-}
 
 {-# RULES
@@ -137,24 +212,26 @@
 multiplyTransposedFused ::
      ( Mutable r Ix2 e
      , OuterSlice r Ix2 e
-     , Source (EltRepr r Ix2) Ix1 e
+     , Source (R r) Ix1 e
      , Num e
+     , MonadThrow m
      )
   => Array r Ix2 e
   -> Array r Ix2 e
-  -> Array r Ix2 e
-multiplyTransposedFused arr1 arr2 = compute (multiplyTransposed arr1 arr2)
+  -> m (Array r Ix2 e)
+multiplyTransposedFused arr1 arr2 = compute <$> multiplyTransposed arr1 arr2
 {-# INLINE multiplyTransposedFused #-}
 
 
-multArrs :: forall r r' e.
+multArrs :: forall r r' e m.
             ( Mutable r Ix2 e
             , Source r' Ix2 e
             , OuterSlice r Ix2 e
-            , Source (EltRepr r Ix2) Ix1 e
+            , Source (R r) Ix1 e
             , Num e
+            , MonadThrow m
             )
-         => Array r Ix2 e -> Array r' Ix2 e -> Array D Ix2 e
+         => Array r Ix2 e -> Array r' Ix2 e -> m (Array D Ix2 e)
 multArrs arr1 arr2 = multiplyTransposed arr1 arr2'
   where
     arr2' :: Array r Ix2 e
@@ -166,15 +243,17 @@
 multiplyTransposed ::
      ( Manifest r Ix2 e
      , OuterSlice r Ix2 e
-     , Source (EltRepr r Ix2) Ix1 e
+     , Source (R r) Ix1 e
      , Num e
+     , MonadThrow m
      )
   => Array r Ix2 e
   -> Array r Ix2 e
-  -> Array D Ix2 e
+  -> m (Array D Ix2 e)
 multiplyTransposed arr1 arr2
-  | n1 /= m2 = throw $ SizeMismatchException (size arr1) (size arr2)
+  | n1 /= m2 = throwM $ SizeMismatchException (size arr1) (size arr2)
   | otherwise =
+    pure $
     DArray (getComp arr1 <> getComp arr2) (SafeSz (m1 :. n2)) $ \(i :. j) ->
       A.foldlS (+) 0 (A.zipWith (*) (unsafeOuterSlice arr1 i) (unsafeOuterSlice arr2 j))
   where
@@ -184,7 +263,7 @@
 
 -- | Create an indentity matrix.
 --
--- ====___Example__
+-- ==== __Example__
 --
 -- >>> import Data.Massiv.Array
 -- >>> identityMatrix 5
@@ -197,51 +276,47 @@
 --   ]
 --
 -- @since 0.3.6
-identityMatrix :: Int -> Array DL Ix2 Int
-identityMatrix n = makeLoadArrayS (Sz2 n n) 0 $ \ w -> loopM_ 0 (< n) (+1) $ \ i -> w (i :. i) 1
+identityMatrix :: Sz1 -> Array DL Ix2 Int
+identityMatrix (Sz n) = makeLoadArrayS (Sz2 n n) 0 $ \ w -> loopM_ 0 (< n) (+1) $ \ i -> w (i :. i) 1
 {-# INLINE identityMatrix #-}
 
 
-negateA
-  :: (Source r ix e, Num e)
-  => Array r ix e -> Array D ix e
-negateA = liftArray negate
+negateA :: (Index ix, Numeric r e) => Array r ix e -> Array r ix e
+negateA = unsafeLiftArray negate
 {-# INLINE negateA #-}
 
-absA
-  :: (Source r ix e, Num e)
-  => Array r ix e -> Array D ix e
-absA = liftArray abs
+absA :: (Index ix, Numeric r e) => Array r ix e -> Array r ix e
+absA = absPointwise
 {-# INLINE absA #-}
 
-signumA
-  :: (Source r ix e, Num e)
-  => Array r ix e -> Array D ix e
-signumA = liftArray signum
+signumA :: (Index ix, Numeric r e) => Array r ix e -> Array r ix e
+signumA = unsafeLiftArray signum
 {-# INLINE signumA #-}
 
-fromIntegerA
-  :: (Index ix, Num e)
-  => Integer -> Array D ix e
+fromIntegerA :: (Index ix, Num e) => Integer -> Array D ix e
 fromIntegerA = singleton . fromInteger
 {-# INLINE fromIntegerA #-}
 
-(./)
-  :: (Source r1 ix e, Source r2 ix e, Fractional e)
-  => Array r1 ix e -> Array r2 ix e -> Array D ix e
-(./) = liftArray2Matching (/)
+(./.) ::
+     (Load r ix e, NumericFloat r e, MonadThrow m)
+  => Array r ix e
+  -> Array r ix e
+  -> m (Array r ix e)
+(./.) = liftNumericArray2M divisionPointwise
+{-# INLINE (./.) #-}
+
+(./) ::(Index ix,  NumericFloat r e) => Array r ix e -> e -> Array r ix e
+(./) = divideScalar
 {-# INLINE (./) #-}
 
 (.^^)
-  :: (Source r ix e, Fractional e, Integral b)
-  => Array r ix e -> b -> Array D ix e
-(.^^) arr n = liftArray (^^ n) arr
+  :: (Index ix, Numeric r e, Fractional e, Integral b)
+  => Array r ix e -> b -> Array r ix e
+(.^^) arr n = unsafeLiftArray (^^ n) arr
 {-# INLINE (.^^) #-}
 
-recipA
-  :: (Source r ix e, Fractional e)
-  => Array r ix e -> Array D ix e
-recipA = liftArray recip
+recipA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+recipA = recipPointwise
 {-# INLINE recipA #-}
 
 
@@ -257,22 +332,16 @@
 piA = singleton pi
 {-# INLINE piA #-}
 
-expA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-expA = liftArray exp
+expA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+expA = unsafeLiftArray exp
 {-# INLINE expA #-}
 
-sqrtA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-sqrtA = liftArray sqrt
+sqrtA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+sqrtA = unsafeLiftArray sqrt
 {-# INLINE sqrtA #-}
 
-logA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-logA = liftArray log
+logA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+logA = unsafeLiftArray log
 {-# INLINE logA #-}
 
 logBaseA
@@ -289,102 +358,78 @@
 
 
 
-sinA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-sinA = liftArray sin
+sinA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+sinA = unsafeLiftArray sin
 {-# INLINE sinA #-}
 
-cosA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-cosA = liftArray cos
+cosA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+cosA = unsafeLiftArray cos
 {-# INLINE cosA #-}
 
-tanA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-tanA = liftArray cos
+tanA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+tanA = unsafeLiftArray cos
 {-# INLINE tanA #-}
 
-asinA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-asinA = liftArray asin
+asinA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+asinA = unsafeLiftArray asin
 {-# INLINE asinA #-}
 
-atanA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-atanA = liftArray atan
+atanA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+atanA = unsafeLiftArray atan
 {-# INLINE atanA #-}
 
-acosA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-acosA = liftArray acos
+acosA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+acosA = unsafeLiftArray acos
 {-# INLINE acosA #-}
 
-sinhA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-sinhA = liftArray sinh
+sinhA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+sinhA = unsafeLiftArray sinh
 {-# INLINE sinhA #-}
 
-tanhA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-tanhA = liftArray cos
+tanhA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+tanhA = unsafeLiftArray cos
 {-# INLINE tanhA #-}
 
-coshA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-coshA = liftArray cosh
+coshA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+coshA = unsafeLiftArray cosh
 {-# INLINE coshA #-}
 
-asinhA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-asinhA = liftArray asinh
+asinhA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+asinhA = unsafeLiftArray asinh
 {-# INLINE asinhA #-}
 
-acoshA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-acoshA = liftArray acosh
+acoshA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+acoshA = unsafeLiftArray acosh
 {-# INLINE acoshA #-}
 
-atanhA
-  :: (Source r ix e, Floating e)
-  => Array r ix e -> Array D ix e
-atanhA = liftArray atanh
+atanhA :: (Index ix, NumericFloat r e) => Array r ix e -> Array r ix e
+atanhA = unsafeLiftArray atanh
 {-# INLINE atanhA #-}
 
 
 quotA
   :: (Source r1 ix e, Source r2 ix e, Integral e)
   => Array r1 ix e -> Array r2 ix e -> Array D ix e
-quotA = liftArray2Matching (quot)
+quotA = liftArray2Matching quot
 {-# INLINE quotA #-}
 
 
 remA
   :: (Source r1 ix e, Source r2 ix e, Integral e)
   => Array r1 ix e -> Array r2 ix e -> Array D ix e
-remA = liftArray2Matching (rem)
+remA = liftArray2Matching rem
 {-# INLINE remA #-}
 
 divA
   :: (Source r1 ix e, Source r2 ix e, Integral e)
   => Array r1 ix e -> Array r2 ix e -> Array D ix e
-divA = liftArray2Matching (div)
+divA = liftArray2Matching div
 {-# INLINE divA #-}
 
 modA
   :: (Source r1 ix e, Source r2 ix e, Integral e)
   => Array r1 ix e -> Array r2 ix e -> Array D ix e
-modA = liftArray2Matching (mod)
+modA = liftArray2Matching mod
 {-# INLINE modA #-}
 
 
@@ -405,35 +450,30 @@
 
 
 truncateA
-  :: (Source r ix a, RealFrac a, Integral b)
-  => Array r ix a -> Array D ix b
-truncateA = liftArray truncate
+  :: (Index ix, Numeric r e, RealFrac a, Integral e)
+  => Array r ix a -> Array r ix e
+truncateA = unsafeLiftArray truncate
 {-# INLINE truncateA #-}
 
 
-roundA
-  :: (Source r ix a, RealFrac a, Integral b)
-  => Array r ix a -> Array D ix b
-roundA = liftArray round
+roundA :: (Index ix, Numeric r e, RealFrac a, Integral e) => Array r ix a -> Array r ix e
+roundA = unsafeLiftArray round
 {-# INLINE roundA #-}
 
 
-ceilingA
-  :: (Source r ix a, RealFrac a, Integral b)
-  => Array r ix a -> Array D ix b
-ceilingA = liftArray ceiling
+ceilingA :: (Index ix, Numeric r e, RealFrac a, Integral e) => Array r ix a -> Array r ix e
+ceilingA = unsafeLiftArray ceiling
 {-# INLINE ceilingA #-}
 
 
-floorA
-  :: (Source r ix a, RealFrac a, Integral b)
-  => Array r ix a -> Array D ix b
-floorA = liftArray floor
+floorA :: (Index ix, Numeric r e, RealFrac a, Integral e) => Array r ix a -> Array r ix e
+floorA = unsafeLiftArray floor
 {-# INLINE floorA #-}
 
-atan2A
-  :: (Source r ix e, RealFloat e)
-  => Array r ix e -> Array r ix e -> Array D ix e
-atan2A = liftArray2Matching atan2
+atan2A ::
+     (Load r ix e, Numeric r e, RealFloat e, MonadThrow m)
+  => Array r ix e
+  -> Array r ix e
+  -> m (Array r ix e)
+atan2A = liftArray2M atan2
 {-# INLINE atan2A #-}
-
diff --git a/src/Data/Massiv/Array/Numeric/Integral.hs b/src/Data/Massiv/Array/Numeric/Integral.hs
--- a/src/Data/Massiv/Array/Numeric/Integral.hs
+++ b/src/Data/Massiv/Array/Numeric/Integral.hs
@@ -217,7 +217,7 @@
   -> Int -- ^ @n@ - Scaling factor, i.e. number of sample points per cell.
   -> Array D ix e
 fromFunction comp f a d (Sz sz) n =
-  (f scale) <$> rangeInclusive comp zeroIndex (liftIndex (n *) sz)
+  f scale <$> rangeInclusive comp zeroIndex (liftIndex (n *) sz)
   where
     nFrac = fromIntegral n
     scale i = a + d * fromIntegral i / nFrac
@@ -244,7 +244,7 @@
   :: (Index ix, Fractional a) =>
      Comp -> ((Int -> a) -> ix -> e) -> a -> a -> Sz ix -> Int -> Array D ix e
 fromFunctionMidpoint comp f a d (Sz sz) n =
-  (\ix -> f scale ix) <$> rangeInclusive comp zeroIndex (liftIndex (\i -> n * i - 1) sz)
+  f scale <$> rangeInclusive comp zeroIndex (liftIndex (\i -> n * i - 1) sz)
   where
     nFrac = fromIntegral n
     dx2 = d / nFrac / 2
diff --git a/src/Data/Massiv/Array/Ops/Construct.hs b/src/Data/Massiv/Array/Ops/Construct.hs
--- a/src/Data/Massiv/Array/Ops/Construct.hs
+++ b/src/Data/Massiv/Array/Ops/Construct.hs
@@ -47,7 +47,6 @@
   , range
   , rangeStepM
   , rangeStep'
-  , rangeStep
   , rangeInclusive
   , rangeStepInclusiveM
   , rangeStepInclusive'
@@ -57,6 +56,7 @@
   , enumFromStepN
     -- ** Expansion
   , expandWithin
+  , expandWithinM
   , expandWithin'
   , expandOuter
   , expandInner
@@ -123,8 +123,7 @@
 runSTA !sz (STA m) = runST (unsafeNew sz >>= m)
 {-# INLINE runSTA  #-}
 
--- | Similar to `makeArray`, but construct the array sequentially using an `Applicative` interface
--- disregarding the supplied `Comp`.
+-- | Similar to `makeArray`, but construct the array sequentially using an `Applicative` interface.
 --
 -- /Note/ - using `Data.Massiv.Array.Mutable.generateArray` or
 -- `Data.Massiv.Array.Mutable.generateArrayS` will always be faster, althought not always possible.
@@ -134,11 +133,10 @@
 --
 makeArrayA ::
      forall r ix e f. (Mutable r ix e, Applicative f)
-  => Comp
-  -> Sz ix
+  => Sz ix
   -> (ix -> f e)
   -> f (Array r ix e)
-makeArrayA !comp !sz f =
+makeArrayA !sz f =
   let n = totalElem sz
       go !i
         | i < n =
@@ -146,7 +144,7 @@
             (\e (STA st) -> STA (\ma -> unsafeLinearWrite ma i e >> st ma))
             (f (fromLinearIndex sz i))
             (go (i + 1))
-        | otherwise = pure (STA (unsafeFreeze comp))
+        | otherwise = pure (STA (unsafeFreeze Seq))
    in runSTA sz <$> go 0
 {-# INLINE makeArrayA  #-}
 
@@ -158,7 +156,6 @@
 makeArrayAR ::
      forall r ix e f. (Mutable r ix e, Applicative f)
   => r
-  -> Comp
   -> Sz ix
   -> (ix -> f e)
   -> f (Array r ix e)
@@ -172,40 +169,47 @@
 -- ==== __Example__
 --
 -- >>> import Data.Massiv.Array
--- >>> iterateN Seq (Sz2 2 10) succ (10 :: Int)
+-- >>> iterateN (Sz2 2 10) succ (10 :: Int)
 -- Array DL Seq (Sz (2 :. 10))
 --   [ [ 11, 12, 13, 14, 15, 16, 17, 18, 19, 20 ]
 --   , [ 21, 22, 23, 24, 25, 26, 27, 28, 29, 30 ]
 --   ]
 --
 -- @since 0.3.0
-iterateN :: forall ix e . Index ix => Comp -> Sz ix -> (e -> e) -> e -> Array DL ix e
-iterateN comp sz f = unfoldrS_ comp sz $ \a -> let !a' = f a in (a', a')
+iterateN :: forall ix e . Index ix => Sz ix -> (e -> e) -> e -> Array DL ix e
+iterateN sz f = unfoldrS_ sz $ \a -> let !a' = f a in (a', a')
 {-# INLINE iterateN #-}
 
 -- | Same as `iterateN`, but with index aware function.
 --
 -- @since 0.3.0
-iiterateN :: forall ix e . Index ix => Comp -> Sz ix -> (e -> ix -> e) -> e -> Array DL ix e
-iiterateN comp sz f = iunfoldrS_ comp sz $ \a ix -> let !a' = f a ix in (a', a')
+iiterateN :: forall ix e . Index ix => Sz ix -> (e -> ix -> e) -> e -> Array DL ix e
+iiterateN sz f = iunfoldrS_ sz $ \a ix -> let !a' = f a ix in (a', a')
 {-# INLINE iiterateN #-}
 
 
--- |
+-- | Right unfold of a delayed load array. For the inverse direction use `unfoldlS_`.
 --
+-- ==== __Examples__
+--
+-- >>> import Data.Massiv.Array
+-- >>> unfoldrS_ (Sz1 10) (\xs -> (head xs, tail xs)) ([10 ..] :: [Int])
+-- Array DL Seq (Sz1 10)
+--   [ 10, 11, 12, 13, 14, 15, 16, 17, 18, 19 ]
+--
 -- @since 0.3.0
-unfoldrS_ :: forall ix e a . Construct DL ix e => Comp -> Sz ix -> (a -> (e, a)) -> a -> Array DL ix e
-unfoldrS_ comp sz f = iunfoldrS_ comp sz (\a _ -> f a)
+unfoldrS_ :: forall ix e a . Construct DL ix e => Sz ix -> (a -> (e, a)) -> a -> Array DL ix e
+unfoldrS_ sz f = iunfoldrS_ sz (\a _ -> f a)
 {-# INLINE unfoldrS_ #-}
 
--- |
+-- | Right unfold of a delayed load array with index aware function
 --
 -- @since 0.3.0
 iunfoldrS_
-  :: Construct DL ix e => Comp -> Sz ix -> (a -> ix -> (e, a)) -> a -> Array DL ix e
-iunfoldrS_ comp sz f acc0 =
+  :: Construct DL ix e => Sz ix -> (a -> ix -> (e, a)) -> a -> Array DL ix e
+iunfoldrS_ sz f acc0 =
   DLArray
-    { dlComp = comp
+    { dlComp = Seq
     , dlSize = sz
     , dlDefault = Nothing
     , dlLoad =
@@ -219,23 +223,23 @@
 {-# INLINE iunfoldrS_ #-}
 
 
--- | Unfold sequentially from the end. There is no way to save the accumulator after unfolding is
--- done, since resulting array is delayed, but it's possible to use
+-- | Unfold sequentially from the end. There is no way to save the accumulator after
+-- unfolding is done, since resulting array is delayed, but it's possible to use
 -- `Data.Massiv.Array.Mutable.unfoldlPrimM` to achive such effect.
 --
 -- @since 0.3.0
-unfoldlS_ :: Construct DL ix e => Comp -> Sz ix -> (a -> (a, e)) -> a -> Array DL ix e
-unfoldlS_ comp sz f = iunfoldlS_ comp sz (const f)
+unfoldlS_ :: Construct DL ix e => Sz ix -> (a -> (a, e)) -> a -> Array DL ix e
+unfoldlS_ sz f = iunfoldlS_ sz (const f)
 {-# INLINE unfoldlS_ #-}
 
 -- | Unfold sequentially from the right with an index aware function.
 --
 -- @since 0.3.0
 iunfoldlS_
-  :: Construct DL ix e => Comp -> Sz ix -> (ix -> a -> (a, e)) -> a -> Array DL ix e
-iunfoldlS_ comp sz f acc0 =
+  :: Construct DL ix e => Sz ix -> (ix -> a -> (a, e)) -> a -> Array DL ix e
+iunfoldlS_ sz f acc0 =
   DLArray
-    { dlComp = comp
+    { dlComp = Seq
     , dlSize = sz
     , dlDefault = Nothing
     , dlLoad =
@@ -443,17 +447,8 @@
 range comp !from !to = rangeSize comp from (Sz (liftIndex2 (-) to from))
 {-# INLINE range #-}
 
-
 -- | Same as `range`, but with a custom step.
 --
--- @since 0.1.0
-rangeStep :: Index ix => Comp -> ix -> ix -> ix -> Maybe (Array D ix ix)
-rangeStep = rangeStepM
-{-# INLINE rangeStep #-}
-{-# DEPRECATED rangeStep "In favor of more general `rangeStepM`" #-}
-
--- | Same as `range`, but with a custom step.
---
 -- ==== __Examples__
 --
 -- >>> import Data.Massiv.Array
@@ -629,17 +624,17 @@
 expandWithin ::
      forall ix e r n a. (IsIndexDimension ix n, Manifest r (Lower ix) a)
   => Dimension n
-  -> Int
-  -> (a -> Int -> e)
+  -> Sz1
+  -> (a -> Ix1 -> e)
   -> Array r (Lower ix) a
   -> Array D ix e
-expandWithin dim k f arr =
+expandWithin dim (Sz k) f arr =
   makeArray (getComp arr) sz $ \ix ->
     let (i, ixl) = pullOutDimension ix dim
      in f (unsafeIndex arr ixl) i
   where
     szl = unSz (size arr)
-    sz = Sz (insertDimension szl dim k)
+    sz = SafeSz (insertDimension szl dim k)
 {-# INLINE expandWithin #-}
 
 -- | Similar to `expandWithin`, except that dimension is specified at a value level, which means it
@@ -649,26 +644,39 @@
 expandWithin'
   :: (Index ix, Manifest r (Lower ix) a)
   => Dim
-  -> Int
-  -> (a -> Int -> b)
+  -> Sz1
+  -> (a -> Ix1 -> b)
   -> Array r (Lower ix) a
   -> Array D ix b
-expandWithin' dim k f arr =
-  makeArray (getComp arr) sz $ \ix ->
-    let (i, ixl) = pullOutDim' ix dim
-     in f (unsafeIndex arr ixl) i
-  where
-    szl = unSz (size arr)
-    sz = Sz (insertDim' szl dim k)
+expandWithin' dim k f arr = either throw id $ expandWithinM dim k f arr
 {-# INLINE expandWithin' #-}
 
+-- | Similar to `expandWithin`, except that dimension is specified at a value level, which means it
+-- will throw an exception on an invalid dimension.
+--
+-- @since 0.4.0
+expandWithinM
+  :: (Index ix, Manifest r (Lower ix) a, MonadThrow m)
+  => Dim
+  -> Sz1
+  -> (a -> Ix1 -> b)
+  -> Array r (Lower ix) a
+  -> m (Array D ix b)
+expandWithinM dim k f arr = do
+  sz <- insertSzM (size arr) dim k
+  pure $
+    makeArray (getComp arr) sz $ \ix ->
+      let (i, ixl) = pullOutDim' ix dim -- dim has been checked above
+       in f (unsafeIndex arr ixl) i
+{-# INLINE expandWithinM #-}
+
 -- | Similar to `expandWithin`, except it uses the outermost dimension.
 --
 -- @since 0.2.6
 expandOuter
   :: (Index ix, Manifest r (Lower ix) a)
-  => Int
-  -> (a -> Int -> b)
+  => Sz1
+  -> (a -> Ix1 -> b)
   -> Array r (Lower ix) a
   -> Array D ix b
 expandOuter k f arr =
@@ -677,7 +685,7 @@
      in f (unsafeIndex arr ixl) i
   where
     szl = size arr
-    sz = consSz (Sz k) szl
+    sz = consSz k szl
 {-# INLINE expandOuter #-}
 
 -- | Similar to `expandWithin`, except it uses the innermost dimension.
@@ -685,8 +693,8 @@
 -- @since 0.2.6
 expandInner
   :: (Index ix, Manifest r (Lower ix) a)
-  => Int
-  -> (a -> Int -> b)
+  => Sz1
+  -> (a -> Ix1 -> b)
   -> Array r (Lower ix) a
   -> Array D ix b
 expandInner k f arr =
@@ -695,5 +703,5 @@
      in f (unsafeIndex arr ixl) i
   where
     szl = size arr
-    sz = snocSz szl (Sz k)
+    sz = snocSz szl k
 {-# INLINE expandInner #-}
diff --git a/src/Data/Massiv/Array/Ops/Fold.hs b/src/Data/Massiv/Array/Ops/Fold.hs
--- a/src/Data/Massiv/Array/Ops/Fold.hs
+++ b/src/Data/Massiv/Array/Ops/Fold.hs
@@ -23,10 +23,8 @@
   , foldSemi
   , minimumM
   , minimum'
-  , minimum
   , maximumM
   , maximum'
-  , maximum
   , sum
   , product
   , and
@@ -89,15 +87,14 @@
   , ifoldrP
   , ifoldlIO
   , ifoldrIO
+  -- , splitReduce
   ) where
 
 import Data.Massiv.Array.Delayed.Pull
 import Data.Massiv.Array.Ops.Fold.Internal
 import Data.Massiv.Core
 import Data.Massiv.Core.Common
-import Data.Massiv.Core.Index.Internal (Sz(..))
-import Prelude hiding (all, and, any, foldl, foldr, map, maximum, minimum, or,
-                product, sum)
+import Prelude hiding (all, and, any, foldl, foldr, map, maximum, minimum, or, product, sum)
 
 -- | /O(n)/ - Monoidal fold over an array with an index aware function. Also known as reduce.
 --
@@ -299,14 +296,6 @@
 
 -- | /O(n)/ - Compute maximum of all elements.
 --
--- @since 0.1.0
-maximum :: (Source r ix e, Ord e) => Array r ix e -> e
-maximum = maximum'
-{-# INLINE maximum #-}
-{-# DEPRECATED maximum "In favor of a safer `maximumM` or an equivalent `maximum'`" #-}
-
--- | /O(n)/ - Compute maximum of all elements.
---
 -- @since 0.3.0
 maximum' :: (Source r ix e, Ord e) => Array r ix e -> e
 maximum' = either throw id . maximumM
@@ -331,13 +320,16 @@
 minimum' = either throw id . minimumM
 {-# INLINE minimum' #-}
 
--- | /O(n)/ - Compute minimum of all elements.
---
--- @since 0.1.0
-minimum :: (Source r ix e, Ord e) => Array r ix e -> e
-minimum = minimum'
-{-# INLINE minimum #-}
-{-# DEPRECATED minimum "In favor of a safer `minimumM` or an equivalent `minimum'`" #-}
+
+-- -- | /O(n)/ - Compute sum of all elements.
+-- --
+-- -- @since 0.1.0
+-- sum' ::
+--      forall r ix e. (Source r ix e, Numeric r e)
+--   => Array r ix e
+--   -> IO e
+-- sum' = splitReduce (\_ -> pure . sumArray) (\x y -> pure (x + y)) 0
+-- {-# INLINE sum' #-}
 
 -- | /O(n)/ - Compute sum of all elements.
 --
diff --git a/src/Data/Massiv/Array/Ops/Fold/Internal.hs b/src/Data/Massiv/Array/Ops/Fold/Internal.hs
--- a/src/Data/Massiv/Array/Ops/Fold/Internal.hs
+++ b/src/Data/Massiv/Array/Ops/Fold/Internal.hs
@@ -41,10 +41,11 @@
   , ifoldrP
   , ifoldlIO
   , ifoldrIO
+  -- , splitReduce
   ) where
 
-import Control.Scheduler
 import Control.Monad (void, when)
+import Control.Scheduler
 import qualified Data.Foldable as F
 import Data.Functor.Identity (runIdentity)
 import Data.Massiv.Core.Common
@@ -326,6 +327,32 @@
             f acc ix (unsafeLinearIndex arr i)
   F.foldlM g tAcc results
 {-# INLINE ifoldlIO #-}
+
+-- -- | Split an array into linear row-major vector chunks and apply an action to each of
+-- -- them. Number of chunks will depend on the computation strategy. Results of each action
+-- -- will be combined with a folding function.
+-- --
+-- -- @since 0.4.1
+-- splitReduce ::
+--      (MonadUnliftIO m, Source r ix e)
+--   => (Scheduler m a -> Array r Ix1 e -> m a)
+--   -> (b -> a -> m b) -- ^ Folding action that is applied to the results of a parallel fold
+--   -> b -- ^ Accumulator for chunks folding
+--   -> Array r ix e
+--   -> m b
+-- splitReduce f g !tAcc !arr = do
+--   let !sz = size arr
+--       !totalLength = totalElem sz
+--   results <-
+--     withScheduler (getComp arr) $ \scheduler ->
+--       splitLinearly (numWorkers scheduler) totalLength $ \chunkLength slackStart -> do
+--         loopM_ 0 (< slackStart) (+ chunkLength) $ \ !start ->
+--           scheduleWork scheduler $ f scheduler $ unsafeLinearSlice start (SafeSz chunkLength) arr
+--         when (slackStart < totalLength) $
+--           scheduleWork scheduler $
+--           f scheduler $ unsafeLinearSlice slackStart (SafeSz (totalLength - slackStart)) arr
+--   F.foldlM g tAcc results
+-- {-# INLINE splitReduce #-}
 
 
 
diff --git a/src/Data/Massiv/Array/Ops/Map.hs b/src/Data/Massiv/Array/Ops/Map.hs
--- a/src/Data/Massiv/Array/Ops/Map.hs
+++ b/src/Data/Massiv/Array/Ops/Map.hs
@@ -201,10 +201,10 @@
   -> Array r2 ix e2
   -> f (Array r ix e)
 izipWithA f arr1 arr2 =
+  setComp (getComp arr1 <> getComp arr2) <$>
   makeArrayA
-    (getComp arr1 <> getComp arr2)
-    (SafeSz (liftIndex2 min (coerce (size arr1)) (coerce (size arr2)))) $ \ !ix ->
-    f ix (unsafeIndex arr1 ix) (unsafeIndex arr2 ix)
+    (SafeSz (liftIndex2 min (coerce (size arr1)) (coerce (size arr2))))
+    (\ !ix -> f ix (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))
 {-# INLINE izipWithA #-}
 
 -- | Same as `zipWithA`, but for three arrays.
@@ -231,8 +231,8 @@
   -> Array r3 ix e3
   -> f (Array r ix e)
 izipWith3A f arr1 arr2 arr3 =
-  makeArrayA (getComp arr1 <> getComp arr2 <> getComp arr3) sz $ \ !ix ->
-    f ix (unsafeIndex arr1 ix) (unsafeIndex arr2 ix) (unsafeIndex arr3 ix)
+  setComp (getComp arr1 <> getComp arr2 <> getComp arr3) <$>
+  makeArrayA sz (\ !ix -> f ix (unsafeIndex arr1 ix) (unsafeIndex arr2 ix) (unsafeIndex arr3 ix))
   where
     sz =
       SafeSz $
@@ -240,6 +240,27 @@
 {-# INLINE izipWith3A #-}
 
 
+
+-- | Similar to `Data.Massiv.Array.zipWith`, except dimensions of both arrays either have to be the
+-- same, or at least one of the two array must be a singleton array, in which case it will behave as
+-- a `Data.Massiv.Array.map`.
+--
+-- @since 0.1.4
+liftArray2
+  :: (Source r1 ix a, Source r2 ix b)
+  => (a -> b -> e) -> Array r1 ix a -> Array r2 ix b -> Array D ix e
+liftArray2 f !arr1 !arr2
+  | sz1 == oneSz = map (f (unsafeIndex arr1 zeroIndex)) arr2
+  | sz2 == oneSz = map (`f` unsafeIndex arr2 zeroIndex) arr1
+  | sz1 == sz2 =
+    DArray (getComp arr1 <> getComp arr2) sz1 (\ !ix -> f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))
+  | otherwise = throw $ SizeMismatchException (size arr1) (size arr2)
+  where
+    sz1 = size arr1
+    sz2 = size arr2
+{-# INLINE liftArray2 #-}
+
+
 --------------------------------------------------------------------------------
 -- traverse --------------------------------------------------------------------
 --------------------------------------------------------------------------------
@@ -255,7 +276,7 @@
   => (a -> f e)
   -> Array r' ix a
   -> f (Array r ix e)
-traverseA f arr = makeArrayA (getComp arr) (size arr) (f . unsafeIndex arr)
+traverseA f arr = setComp (getComp arr) <$> makeArrayA (size arr) (f . unsafeIndex arr)
 {-# INLINE traverseA #-}
 
 -- | Traverse sequentially over a source array, while discarding the result.
@@ -293,7 +314,8 @@
   => (ix -> a -> f e)
   -> Array r' ix a
   -> f (Array r ix e)
-itraverseA f arr = makeArrayA (getComp arr) (size arr) $ \ !ix -> f ix (unsafeIndex arr ix)
+itraverseA f arr =
+  setComp (getComp arr) <$> makeArrayA (size arr) (\ !ix -> f ix (unsafeIndex arr ix))
 {-# INLINE itraverseA #-}
 
 
@@ -343,11 +365,11 @@
 -- @since 0.3.0
 --
 traversePrim ::
-     (Source r' ix a, Mutable r ix b, PrimMonad m)
+     forall r ix b r' a m . (Source r' ix a, Mutable r ix b, PrimMonad m)
   => (a -> m b)
   -> Array r' ix a
   -> m (Array r ix b)
-traversePrim f arr = generateArrayS (getComp arr) (size arr) (f . unsafeIndex arr)
+traversePrim f = itraversePrim (const f)
 {-# INLINE traversePrim #-}
 
 -- | Same as `traversePrim`, but traverse with index aware action.
@@ -355,11 +377,17 @@
 -- @since 0.3.0
 --
 itraversePrim ::
-     (Source r' ix a, Mutable r ix b, PrimMonad m)
+     forall r ix b r' a m . (Source r' ix a, Mutable r ix b, PrimMonad m)
   => (ix -> a -> m b)
   -> Array r' ix a
   -> m (Array r ix b)
-itraversePrim f arr = generateArrayS (getComp arr) (size arr) (\ !ix -> f ix (unsafeIndex arr ix))
+itraversePrim f arr =
+  setComp (getComp arr) <$>
+  generateArrayLinearS
+    (size arr)
+    (\ !i ->
+       let ix = fromLinearIndex (size arr) i
+        in f ix (unsafeLinearIndex arr i))
 {-# INLINE itraversePrim #-}
 
 
diff --git a/src/Data/Massiv/Array/Ops/Slice.hs b/src/Data/Massiv/Array/Ops/Slice.hs
--- a/src/Data/Massiv/Array/Ops/Slice.hs
+++ b/src/Data/Massiv/Array/Ops/Slice.hs
@@ -101,7 +101,7 @@
 -- >>> arr !?> 2 ??> 0 ??> -1 :: Maybe Ix3T
 -- Nothing
 -- >>> arr !?> 2 ??> -10 ?? 1
--- *** Exception: IndexOutOfBoundsException: -10 not safe for (Sz1 2)
+-- *** Exception: IndexOutOfBoundsException: -10 is not safe for (Sz1 2)
 --
 -- @since 0.1.0
 (??>) :: (MonadThrow m, OuterSlice r ix e) => m (Array r ix e) -> Int -> m (Elt r ix e)
diff --git a/src/Data/Massiv/Array/Ops/Sort.hs b/src/Data/Massiv/Array/Ops/Sort.hs
--- a/src/Data/Massiv/Array/Ops/Sort.hs
+++ b/src/Data/Massiv/Array/Ops/Sort.hs
@@ -62,20 +62,9 @@
 -- @since 0.3.2
 quicksort ::
      (Mutable r Ix1 e, Ord e) => Array r Ix1 e -> Array r Ix1 e
-quicksort arr = unsafePerformIO $ withMArray' arr quicksortM_
+quicksort arr = unsafePerformIO $ withMArray arr quicksortM_
 {-# INLINE quicksort #-}
 
--- TODO: switch to the regular withMArray, once it is fixed.
-withMArray' ::
-     (Mutable r ix e, MonadUnliftIO m)
-  => Array r ix e
-  -> (Scheduler m () -> MArray RealWorld r ix e -> m a)
-  -> m (Array r ix e)
-withMArray' arr action = do
-  marr <- thaw arr
-  withScheduler_ (getComp arr) $ \scheduler -> action scheduler marr
-  liftIO $ unsafeFreeze (getComp arr) marr
-{-# INLINE withMArray' #-}
 
 
 -- | Mutable version of `quicksort`
diff --git a/src/Data/Massiv/Array/Ops/Transform.hs b/src/Data/Massiv/Array/Ops/Transform.hs
--- a/src/Data/Massiv/Array/Ops/Transform.hs
+++ b/src/Data/Massiv/Array/Ops/Transform.hs
@@ -20,18 +20,14 @@
   -- ** Backpermute
   , backpermuteM
   , backpermute'
-  , backpermute
   -- ** Resize
   , resizeM
   , resize'
-  , resize
   , flatten
   -- ** Extract
   , extractM
-  , extract
   , extract'
   , extractFromToM
-  , extractFromTo
   , extractFromTo'
   , deleteRowsM
   , deleteColumnsM
@@ -42,12 +38,10 @@
   , snoc
   , unsnocM
   , appendM
-  , append
   , append'
   , concatM
   , concat'
   , splitAtM
-  , splitAt
   , splitAt'
   , splitExtractM
   -- ** Upsample/Downsample
@@ -60,9 +54,6 @@
   , transform'
   , transform2M
   , transform2'
-  -- ** Traverse (deprecated)
-  , traverse
-  , traverse2
   ) where
 
 import Control.Scheduler (traverse_)
@@ -86,7 +77,7 @@
          => ix -- ^ Starting index
          -> Sz ix -- ^ Size of the resulting array
          -> Array r ix e -- ^ Source array
-         -> m (Array (EltRepr r ix) ix e)
+         -> m (Array (R r) ix e)
 extractM !sIx !newSz !arr
   | isSafeIndex sz1 sIx && isSafeIndex eIx1 sIx && isSafeIndex sz1 eIx =
     pure $ unsafeExtract sIx newSz arr
@@ -105,29 +96,11 @@
         => ix -- ^ Starting index
         -> Sz ix -- ^ Size of the resulting array
         -> Array r ix e -- ^ Source array
-        -> Array (EltRepr r ix) ix e
+        -> Array (R r) ix e
 extract' sIx newSz = either throw id . extractM sIx newSz
 {-# INLINE extract' #-}
 
--- | Extract a sub-array from within a larger source array. Array that is being extracted must be
--- fully encapsulated in a source array, otherwise `Nothing` is returned,
-extract :: Extract r ix e
-        => ix -- ^ Starting index
-        -> Sz ix -- ^ Size of the resulting array
-        -> Array r ix e -- ^ Source array
-        -> Maybe (Array (EltRepr r ix) ix e)
-extract !sIx !newSz !arr
-  | isSafeIndex sz1 sIx && isSafeIndex eIx1 sIx && isSafeIndex sz1 eIx =
-    Just $ unsafeExtract sIx newSz arr
-  | otherwise = Nothing
-  where
-    sz1 = Sz (liftIndex (+1) (unSz (size arr)))
-    eIx1 = Sz (liftIndex (+1) eIx)
-    eIx = liftIndex2 (+) sIx $ unSz newSz
-{-# INLINE extract #-}
-{-# DEPRECATED extract "In favor of a more general `extractM`" #-}
 
-
 -- | Similar to `extractM`, except it takes starting and ending index. Result array will not include
 -- the ending index.
 --
@@ -136,21 +109,10 @@
                   ix -- ^ Starting index
                -> ix -- ^ Index up to which elements should be extracted.
                -> Array r ix e -- ^ Source array.
-               -> m (Array (EltRepr r ix) ix e)
+               -> m (Array (R r) ix e)
 extractFromToM sIx eIx = extractM sIx (Sz (liftIndex2 (-) eIx sIx))
 {-# INLINE extractFromToM #-}
 
--- | Similar to `extract`, except it takes starting and ending index. Result array will not include
--- the ending index.
-extractFromTo :: Extract r ix e =>
-                 ix -- ^ Starting index
-              -> ix -- ^ Index up to which elmenets should be extracted.
-              -> Array r ix e -- ^ Source array.
-              -> Maybe (Array (EltRepr r ix) ix e)
-extractFromTo sIx eIx = extract sIx $ Sz (liftIndex2 (-) eIx sIx)
-{-# INLINE extractFromTo #-}
-{-# DEPRECATED extractFromTo "In favor of a more general `extractFromToM`" #-}
-
 -- | Same as `extractFromTo`, but throws an error on invalid indices.
 --
 -- @since 0.2.4
@@ -158,23 +120,13 @@
                  ix -- ^ Starting index
               -> ix -- ^ Index up to which elmenets should be extracted.
               -> Array r ix e -- ^ Source array.
-              -> Array (EltRepr r ix) ix e
+              -> Array (R r) ix e
 extractFromTo' sIx eIx = extract' sIx $ Sz (liftIndex2 (-) eIx sIx)
 {-# INLINE extractFromTo' #-}
 
 
 -- | /O(1)/ - Changes the shape of an array. Returns `Nothing` if total
 -- number of elements does not match the source array.
-resize ::
-     (Index ix', Load r ix e, Resize r ix) => Sz ix' -> Array r ix e -> Maybe (Array r ix' e)
-resize !sz !arr
-  | totalElem sz == totalElem (size arr) = Just $ unsafeResize sz arr
-  | otherwise = Nothing
-{-# INLINE resize #-}
-{-# DEPRECATED resize "In favor of a more general `resizeM`" #-}
-
--- | /O(1)/ - Changes the shape of an array. Returns `Nothing` if total
--- number of elements does not match the source array.
 --
 -- @since 0.3.0
 resizeM ::
@@ -367,7 +319,7 @@
 backpermuteM sz ixF !arr = generateArray (getComp arr) sz (evaluateM arr . ixF)
 {-# INLINE backpermuteM #-}
 
--- | Similar to `backpermuteM`, with few notable differences:
+-- | Similar to `backpermuteM`, with a few notable differences:
 --
 -- * Creates a delayed array, instead of manifest, therefore it can be fused
 -- * Respects computation strategy, so it can be parallelized
@@ -382,19 +334,7 @@
 backpermute' sz ixF !arr = makeArray (getComp arr) sz (evaluate' arr . ixF)
 {-# INLINE backpermute' #-}
 
--- | See `backpermute'`.
---
--- @since 0.1.0
-backpermute :: (Source r' ix' e, Index ix) =>
-               Sz ix -- ^ Size of the result array
-            -> (ix -> ix') -- ^ A function that maps indices of the new array into the source one.
-            -> Array r' ix' e -- ^ Source array.
-            -> Array D ix e
-backpermute = backpermute'
-{-# INLINE backpermute #-}
-{-# DEPRECATED backpermute "In favor of a safe `backpermuteM` or an equivalent `backpermute'`" #-}
 
-
 -- | /O(1)/ - Add an element to the vector from the left side
 --
 -- @since 0.3.0
@@ -516,14 +456,7 @@
       }
 {-# INLINE appendM #-}
 
--- | Append two arrays together along a specified dimension.
-append :: (Source r1 ix e, Source r2 ix e) =>
-          Dim -> Array r1 ix e -> Array r2 ix e -> Maybe (Array DL ix e)
-append = appendM
-{-# INLINE append #-}
-{-# DEPRECATED append "In favor of a more general `appendM`" #-}
 
-
 -- | Same as `appendM`, but will throw an exception in pure code on mismatched sizes.
 --
 -- @since 0.3.0
@@ -585,11 +518,11 @@
 --
 -- @since 0.3.0
 splitAtM ::
-     (MonadThrow m, Extract r ix e, r' ~ EltRepr r ix)
+     (MonadThrow m, Extract r ix e)
   => Dim -- ^ Dimension along which to split
   -> Int -- ^ Index along the dimension to split at
   -> Array r ix e -- ^ Source array
-  -> m (Array r' ix e, Array r' ix e)
+  -> m (Array (R r) ix e, Array (R r) ix e)
 splitAtM dim i arr = do
   let Sz sz = size arr
   eIx <- setDimM sz dim i
@@ -599,23 +532,12 @@
   return (arr1, arr2)
 {-# INLINE splitAtM #-}
 
--- | /O(1)/ - Split an array at an index along a specified dimension.
-splitAt ::
-     (Extract r ix e, r' ~ EltRepr r ix)
-  => Dim -- ^ Dimension along which to split
-  -> Int -- ^ Index along the dimension to split at
-  -> Array r ix e -- ^ Source array
-  -> Maybe (Array r' ix e, Array r' ix e)
-splitAt = splitAtM
-{-# INLINE splitAt #-}
-{-# DEPRECATED splitAt "In favor of a more general `splitAtM`" #-}
-
 -- | Same as `splitAt`, but will throw an error instead of returning `Nothing` on wrong dimension
 -- and index out of bounds.
 --
 -- @since 0.1.0
-splitAt' :: (Extract r ix e, r' ~ EltRepr r ix) =>
-           Dim -> Int -> Array r ix e -> (Array r' ix e, Array r' ix e)
+splitAt' :: Extract r ix e =>
+            Dim -> Int -> Array r ix e -> (Array (R r) ix e, Array (R r) ix e)
 splitAt' dim i arr = either throw id $ splitAtM dim i arr
 {-# INLINE splitAt' #-}
 
@@ -624,12 +546,12 @@
 --
 -- @since 0.3.5
 splitExtractM ::
-     (MonadThrow m, Extract r ix e, Source r' ix e, r' ~ EltRepr r ix)
+     (MonadThrow m, Extract r ix e, Source (R r) ix e)
   => Dim -- ^ Dimension along which to do the extraction
   -> Ix1 -- ^ Start index along the dimension that needs to be extracted
   -> Sz Ix1 -- ^ Size of the extracted array along the dimension that it will be extracted
   -> Array r ix e
-  -> m (Array r' ix e, Array r' ix e, Array r' ix e)
+  -> m (Array (R r) ix e, Array (R r) ix e, Array (R r) ix e)
 splitExtractM dim startIx1 (Sz extractSzIx1) arr = do
   let Sz szIx = size arr
   midStartIx <- setDimM zeroIndex dim startIx1
@@ -667,7 +589,7 @@
 --
 -- @since 0.3.5
 deleteRegionM ::
-     (MonadThrow m, Extract r ix e, Source (EltRepr r ix) ix e)
+     (MonadThrow m, Extract r ix e, Source (R r) ix e)
   => Dim -- ^ Along which axis should the removal happen
   -> Ix1 -- ^ At which index to start dropping slices
   -> Sz Ix1 -- ^ Number of slices to drop
@@ -699,7 +621,7 @@
 --
 -- @since 0.3.5
 deleteRowsM ::
-     (MonadThrow m, Extract r ix e, Source (EltRepr r ix) ix e, Index (Lower ix))
+     (MonadThrow m, Extract r ix e, Source (R r) ix e, Index (Lower ix))
   => Ix1
   -> Sz Ix1
   -> Array r ix e
@@ -728,7 +650,7 @@
 --
 -- @since 0.3.5
 deleteColumnsM ::
-     (MonadThrow m, Extract r ix e, Source (EltRepr r ix) ix e)
+     (MonadThrow m, Extract r ix e, Source (R r) ix e)
   => Ix1
   -> Sz Ix1
   -> Array r ix e
@@ -794,31 +716,6 @@
     !newsz = SafeSz (timesStride $ unSz sz)
 {-# INLINE upsample #-}
 
-
--- | Create an array by traversing a source array.
-traverse
-  :: (Source r1 ix1 e1, Index ix)
-  => Sz ix -- ^ Size of the result array
-  -> ((ix1 -> e1) -> ix -> e) -- ^ Function that will receive a source array safe index function and
-                              -- an index for an element it should return a value of.
-  -> Array r1 ix1 e1 -- ^ Source array
-  -> Array D ix e
-traverse sz f arr1 = makeArray (getComp arr1) sz (f (evaluate' arr1))
-{-# INLINE traverse #-}
-{-# DEPRECATED traverse "In favor of more general `transform'`" #-}
-
--- | Create an array by traversing two source arrays.
-traverse2
-  :: (Source r1 ix1 e1, Source r2 ix2 e2, Index ix)
-  => Sz ix
-  -> ((ix1 -> e1) -> (ix2 -> e2) -> ix -> e)
-  -> Array r1 ix1 e1
-  -> Array r2 ix2 e2
-  -> Array D ix e
-traverse2 sz f arr1 arr2 =
-  makeArray (getComp arr1 <> getComp arr2) sz (f (evaluate' arr1) (evaluate' arr2))
-{-# INLINE traverse2 #-}
-{-# DEPRECATED traverse2 "In favor of more general `transform2'`" #-}
 
 -- | General array transformation, that forces computation and produces a manifest array.
 --
diff --git a/src/Data/Massiv/Array/Unsafe.hs b/src/Data/Massiv/Array/Unsafe.hs
--- a/src/Data/Massiv/Array/Unsafe.hs
+++ b/src/Data/Massiv/Array/Unsafe.hs
@@ -26,24 +26,32 @@
   , unsafeExtract
   , unsafeTransform
   , unsafeTransform2
-    -- ** Deprecated
-  , unsafeTraverse
-  , unsafeTraverse2
     -- * Slicing
   , unsafeSlice
   , unsafeOuterSlice
   , unsafeInnerSlice
+  -- , unsafeLinearSlice
     -- * Mutable interface
   , unsafeThaw
   , unsafeFreeze
   , unsafeNew
+    -- ** Read
   , unsafeRead
   , unsafeLinearRead
+    -- ** Write
   , unsafeWrite
   , unsafeLinearWrite
+    -- ** Modify
+  , unsafeModify
+  , unsafeLinearModify
+    -- ** Swap
+  , unsafeSwap
+  , unsafeLinearSwap
+    -- ** Range modification
   , unsafeLinearSet
   , unsafeLinearCopy
   , unsafeArrayLinearCopy
+    -- ** Resizing
   , unsafeLinearShrink
   , unsafeLinearGrow
     -- * Pointer access
@@ -83,31 +91,6 @@
 unsafeBackpermute !sz ixF !arr =
   makeArray (getComp arr) sz $ \ !ix -> unsafeIndex arr (ixF ix)
 {-# INLINE unsafeBackpermute #-}
-
-
-unsafeTraverse
-  :: (Source r ix' e', Index ix)
-  => Sz ix
-  -> ((ix' -> e') -> ix -> e)
-  -> Array r ix' e'
-  -> Array D ix e
-unsafeTraverse sz f arr1 =
-  makeArray (getComp arr1) sz (f (unsafeIndex arr1))
-{-# INLINE unsafeTraverse #-}
-{-# DEPRECATED unsafeTraverse "In favor of more general `unsafeTransform'`" #-}
-
-
-unsafeTraverse2
-  :: (Source r1 ix1 e1, Source r2 ix2 e2, Index ix)
-  => Sz ix
-  -> ((ix1 -> e1) -> (ix2 -> e2) -> ix -> e)
-  -> Array r1 ix1 e1
-  -> Array r2 ix2 e2
-  -> Array D ix e
-unsafeTraverse2 sz f arr1 arr2 =
-  makeArray (getComp arr1 <> getComp arr2) sz (f (unsafeIndex arr1) (unsafeIndex arr2))
-{-# INLINE unsafeTraverse2 #-}
-{-# DEPRECATED unsafeTraverse2 "In favor of more general `unsafeTransform2'`" #-}
 
 -- | Same `Data.Array.transform'`, except no bounds checking is performed, thus making it faster,
 -- but unsafe.
diff --git a/src/Data/Massiv/Core.hs b/src/Data/Massiv/Core.hs
--- a/src/Data/Massiv/Core.hs
+++ b/src/Data/Massiv/Core.hs
@@ -9,9 +9,8 @@
 module Data.Massiv.Core
   ( Array(List, unList)
   , Elt
-  , EltRepr
   , Construct
-  , Load(loadArrayM, defaultElement)
+  , Load(R, loadArrayM, defaultElement)
   , Source
   , Resize
   , Extract
diff --git a/src/Data/Massiv/Core/Common.hs b/src/Data/Massiv/Core/Common.hs
--- a/src/Data/Massiv/Core/Common.hs
+++ b/src/Data/Massiv/Core/Common.hs
@@ -16,7 +16,6 @@
 module Data.Massiv.Core.Common
   ( Array
   , Elt
-  , EltRepr
   , Construct(..)
   , Source(..)
   , Load(..)
@@ -36,7 +35,10 @@
   , WorkerStates
   , unsafeRead
   , unsafeWrite
+  , unsafeModify
   , unsafeLinearModify
+  , unsafeSwap
+  , unsafeLinearSwap
   , unsafeDefaultLinearShrink
   , Ragged(..)
   , Nested(..)
@@ -51,7 +53,6 @@
   , (!?)
   , index
   , indexM
-  , indexWith
   , (!)
   , index'
   , (??)
@@ -59,7 +60,6 @@
   , borderIndex
   , evaluateM
   , evaluate'
-  , evaluateAt
   , module Data.Massiv.Core.Index
   -- * Common Operations
   , imapM_
@@ -71,6 +71,7 @@
   , SizeException(..)
   , ShapeException(..)
   , module Data.Massiv.Core.Exception
+  , Proxy(..)
   -- * Stateful Monads
   , MonadUnliftIO
   , MonadIO(liftIO)
@@ -96,14 +97,12 @@
 -- | The array family. Representations @r@ describes how data is arranged or computed. All arrays
 -- have a common property that each index @ix@ always maps to the same unique element, even if that
 -- element does not exist in memory and has to be computed upon lookup. Data is always arranged in a
--- nested fasion, depth of which is controlled by @`Rank` ix@.
+-- nested fashion, depth of which is controlled by @`Rank` ix@.
 data family Array r ix e :: *
 
-type family EltRepr r ix :: *
-
 type family Elt r ix e :: * where
   Elt r Ix1 e = e
-  Elt r ix  e = Array (EltRepr r ix) (Lower ix) e
+  Elt r ix  e = Array (R r) (Lower ix) e
 
 type family NestedStruct r ix e :: *
 
@@ -183,7 +182,7 @@
 class Load r ix e => Extract r ix e where
   -- | /O(1)/ - Extract a portion of an array. Staring index and new size are
   -- not validated.
-  unsafeExtract :: ix -> Sz ix -> Array r ix e -> Array (EltRepr r ix) ix e
+  unsafeExtract :: ix -> Sz ix -> Array r ix e -> Array (R r) ix e
 
 
 -- | Arrays that can be used as source to practically any manipulation function.
@@ -192,20 +191,31 @@
 
   -- | Lookup element in the array. No bounds check is performed and access of
   -- arbitrary memory is possible when invalid index is supplied.
+  --
+  -- @since 0.1.0
   unsafeIndex :: Array r ix e -> ix -> e
   unsafeIndex =
     INDEX_CHECK("(Source r ix e).unsafeIndex",
                 size, \ !arr -> unsafeLinearIndex arr . toLinearIndex (size arr))
   {-# INLINE unsafeIndex #-}
 
-  -- | Lookup element in the array using flat index in a row-major fasion. No
+  -- | Lookup element in the array using flat index in a row-major fashion. No
   -- bounds check is performed
+  --
+  -- @since 0.1.0
   unsafeLinearIndex :: Array r ix e -> Int -> e
   unsafeLinearIndex !arr = unsafeIndex arr . fromLinearIndex (size arr)
   {-# INLINE unsafeLinearIndex #-}
 
+  -- -- | Source arrays also give us ability to look at their linear slices
+  -- --
+  -- -- @since 0.4.0
+  -- unsafeLinearSlice :: Ix1 -> Sz1 -> Array r ix e -> Array r Ix1 e
+
 -- | Any array that can be computed and loaded into memory
 class (Typeable r, Index ix) => Load r ix e where
+  type family R r :: *
+  type instance R r = r
 
   -- | Get computation strategy of this array
   --
@@ -280,7 +290,7 @@
   unsafeLinearIndexM :: Array r ix e -> Int -> e
 
 
-class Manifest r ix e => Mutable r ix e where
+class (Construct r ix e, Manifest r ix e) => Mutable r ix e where
   data MArray s r ix e :: *
 
   -- | Get the size of a mutable array.
@@ -321,12 +331,11 @@
   -- | Create new mutable array while initializing all elements to some default value.
   --
   -- @since 0.3.0
-  initializeNew :: PrimMonad m =>
-                   Maybe e -> Sz ix -> m (MArray (PrimState m) r ix e)
+  initializeNew :: PrimMonad m => Maybe e -> Sz ix -> m (MArray (PrimState m) r ix e)
   initializeNew mdef sz = do
     marr <- unsafeNew sz
     case mdef of
-      Just val -> unsafeLinearSet marr 0 (totalElem sz) val
+      Just val -> unsafeLinearSet marr 0 (SafeSz (totalElem sz)) val
       Nothing  -> initialize marr
     return marr
   {-# INLINE initializeNew #-}
@@ -335,16 +344,16 @@
   --
   -- @since 0.3.0
   unsafeLinearSet :: PrimMonad m =>
-                     MArray (PrimState m) r ix e -> Ix1 -> Int -> e -> m ()
+                     MArray (PrimState m) r ix e -> Ix1 -> Sz1 -> e -> m ()
   unsafeLinearSet marr offset len e =
-    loopM_ offset (< (offset + len)) (+1) (\i -> unsafeLinearWrite marr i e)
+    loopM_ offset (< (offset + unSz len)) (+1) (\i -> unsafeLinearWrite marr i e)
   {-# INLINE unsafeLinearSet #-}
 
   -- | Copy part of one mutable array into another
   --
   -- @since 0.3.6
-  unsafeLinearCopy :: PrimMonad m =>
-                      MArray (PrimState m) r ix e -- ^ Source mutable array
+  unsafeLinearCopy :: (Mutable r ix' e, PrimMonad m) =>
+                      MArray (PrimState m) r ix' e -- ^ Source mutable array
                    -> Ix1 -- ^ Starting index at source array
                    -> MArray (PrimState m) r ix e -- ^ Target mutable array
                    -> Ix1 -- ^ Starting index at target array
@@ -359,8 +368,8 @@
   -- | Copy a part of a pure array into a mutable array
   --
   -- @since 0.3.6
-  unsafeArrayLinearCopy :: PrimMonad m =>
-                           Array r ix e -- ^ Source pure array
+  unsafeArrayLinearCopy :: (Mutable r ix' e, PrimMonad m) =>
+                           Array r ix' e -- ^ Source pure array
                         -> Ix1 -- ^ Starting index at source array
                         -> MArray (PrimState m) r ix e -- ^ Target mutable array
                         -> Ix1 -- ^ Starting index at target array
@@ -372,19 +381,23 @@
       unsafeLinearWrite marrTo (i + delta) (unsafeLinearIndex arrFrom i)
   {-# INLINE unsafeArrayLinearCopy #-}
 
-  -- | Linearly reduce the size of an array. Total number of elements should be smaller.
+  -- | Linearly reduce the size of an array. Total number of elements should be smaller or
+  -- equal. There is no guarantee that the original array is left unchanged, so it should
+  -- no longer be used.
   --
   -- @since 0.3.6
   unsafeLinearShrink :: PrimMonad m =>
-                     MArray (PrimState m) r ix e -> Sz ix -> m (MArray (PrimState m) r ix e)
+                        MArray (PrimState m) r ix e -> Sz ix -> m (MArray (PrimState m) r ix e)
   unsafeLinearShrink = unsafeDefaultLinearShrink
   {-# INLINE unsafeLinearShrink #-}
 
-  -- | Linearly increase the size of an array. Total number of elements should be larger.
+  -- | Linearly increase the size of an array. Total number of elements should be larger
+  -- or equal. There is no guarantee that the original array is left unchanged, so it
+  -- should no longer be used.
   --
   -- @since 0.3.6
   unsafeLinearGrow :: PrimMonad m =>
-                     MArray (PrimState m) r ix e -> Sz ix -> m (MArray (PrimState m) r ix e)
+                      MArray (PrimState m) r ix e -> Sz ix -> m (MArray (PrimState m) r ix e)
   unsafeLinearGrow marr sz = do
     marr' <- unsafeNew sz
     unsafeLinearCopy marr 0 marr' 0 $ SafeSz (totalElem (msize marr))
@@ -405,28 +418,68 @@
 
 
 -- | Read an array element
+--
+-- @since 0.1.0
 unsafeRead :: (Mutable r ix e, PrimMonad m) =>
                MArray (PrimState m) r ix e -> ix -> m e
 unsafeRead !marr !ix = unsafeLinearRead marr (toLinearIndex (msize marr) ix)
 {-# INLINE unsafeRead #-}
 
 -- | Write an element into array
+--
+-- @since 0.1.0
 unsafeWrite :: (Mutable r ix e, PrimMonad m) =>
                MArray (PrimState m) r ix e -> ix -> e -> m ()
 unsafeWrite !marr !ix = unsafeLinearWrite marr (toLinearIndex (msize marr) ix)
 {-# INLINE unsafeWrite #-}
 
 
--- | Modify an element in the array with an index aware action.
+-- | Modify an element in the array with a monadic action. Returns the previous value.
+--
+-- @since 0.4.0
 unsafeLinearModify :: (Mutable r ix e, PrimMonad m) =>
-                      MArray (PrimState m) r ix e -> (Int -> e -> m e) -> Int -> m ()
+                      MArray (PrimState m) r ix e -> (e -> m e) -> Int -> m e
 unsafeLinearModify !marr f !i = do
   v <- unsafeLinearRead marr i
-  v' <- f i v
+  v' <- f v
   unsafeLinearWrite marr i v'
+  pure v
 {-# INLINE unsafeLinearModify #-}
 
+-- | Modify an element in the array with a monadic action. Returns the previous value.
+--
+-- @since 0.4.0
+unsafeModify :: (Mutable r ix e, PrimMonad m) =>
+                MArray (PrimState m) r ix e -> (e -> m e) -> ix -> m e
+unsafeModify marr f ix = unsafeLinearModify marr f (toLinearIndex (msize marr) ix)
+{-# INLINE unsafeModify #-}
 
+-- | Swap two elements in a mutable array under the supplied indices. Returns the previous
+-- values.
+--
+-- @since 0.4.0
+unsafeSwap :: (Mutable r ix e, PrimMonad m) =>
+                    MArray (PrimState m) r ix e -> ix -> ix -> m (e, e)
+unsafeSwap !marr !ix1 !ix2 = unsafeLinearSwap marr (toLinearIndex sz ix1) (toLinearIndex sz ix2)
+  where sz = msize marr
+{-# INLINE unsafeSwap #-}
+
+
+-- | Swap two elements in a mutable array under the supplied linear indices. Returns the
+-- previous values.
+--
+-- @since 0.4.0
+unsafeLinearSwap :: (Mutable r ix e, PrimMonad m) =>
+                    MArray (PrimState m) r ix e -> Int -> Int -> m (e, e)
+unsafeLinearSwap !marr !i1 !i2 = do
+  val1 <- unsafeLinearRead marr i1
+  val2 <- unsafeLinearRead marr i2
+  unsafeLinearWrite marr i1 val2
+  unsafeLinearWrite marr i2 val1
+  return (val1, val2)
+{-# INLINE unsafeLinearSwap #-}
+
+
 class Nested r ix e where
   fromNested :: NestedStruct r ix e -> Array r ix e
 
@@ -511,8 +564,10 @@
 
 -- | Infix version of `index'`.
 --
+-- ==== __Examples__
+--
 -- >>> import Data.Massiv.Array as A
--- >>> a = computeAs U $ iterateN Seq (Sz (2 :. 3)) succ (0 :: Int)
+-- >>> a = computeAs U $ iterateN (Sz (2 :. 3)) succ (0 :: Int)
 -- >>> a
 -- Array U Seq (Sz (2 :. 3))
 --   [ [ 1, 2, 3 ]
@@ -521,7 +576,7 @@
 -- >>> a ! 0 :. 2
 -- 3
 -- >>> a ! 0 :. 3
--- *** Exception: IndexOutOfBoundsException: (0 :. 3) not safe for (Sz (2 :. 3))
+-- *** Exception: IndexOutOfBoundsException: (0 :. 3) is not safe for (Sz (2 :. 3))
 --
 -- @since 0.1.0
 (!) :: Manifest r ix e => Array r ix e -> ix -> e
@@ -529,7 +584,7 @@
 {-# INLINE (!) #-}
 
 
--- | Infix version of `index`.
+-- | Infix version of `indexM`.
 --
 -- ==== __Examples__
 --
@@ -544,7 +599,7 @@
 -- >>> a !? 0 :. 2
 -- 3
 -- >>> a !? 0 :. 3
--- *** Exception: IndexOutOfBoundsException: (0 :. 3) not safe for (Sz (2 :. 3))
+-- *** Exception: IndexOutOfBoundsException: (0 :. 3) is not safe for (Sz (2 :. 3))
 -- >>> a !? 0 :. 3 :: Maybe Int
 -- Nothing
 --
@@ -624,6 +679,8 @@
 -- | /O(1)/ - Lookup an element in the array. Use a border resolution technique
 -- when index is out of bounds.
 --
+-- ==== __Examples__
+--
 -- >>> import Data.Massiv.Array as A
 -- >>> :set -XOverloadedLists
 -- >>> xs = [0..100] :: Array U Ix1 Int
@@ -647,7 +704,7 @@
 -- >>> index' xs 50
 -- 50
 -- >>> index' xs 150
--- *** Exception: IndexOutOfBoundsException: 150 not safe for (Sz1 101)
+-- *** Exception: IndexOutOfBoundsException: 150 is not safe for (Sz1 101)
 --
 -- @since 0.1.0
 index' :: Manifest r ix e => Array r ix e -> ix -> e
@@ -666,7 +723,7 @@
 -- >>> evaluateM (range Seq (Ix2 10 20) (100 :. 210)) 50 :: Either SomeException Ix2
 -- Right (60 :. 70)
 -- >>> evaluateM (range Seq (Ix2 10 20) (100 :. 210)) 150 :: Either SomeException Ix2
--- Left (IndexOutOfBoundsException: (150 :. 150) not safe for (Sz (90 :. 190)))
+-- Left (IndexOutOfBoundsException: (150 :. 150) is not safe for (Sz (90 :. 190)))
 --
 -- @since 0.3.0
 evaluateM :: (Source r ix e, MonadThrow m) => Array r ix e -> ix -> m e
@@ -686,7 +743,7 @@
 -- >>> evaluate' (range Seq (Ix2 10 20) (100 :. 210)) 50
 -- 60 :. 70
 -- >>> evaluate' (range Seq (Ix2 10 20) (100 :. 210)) 150
--- *** Exception: IndexOutOfBoundsException: (150 :. 150) not safe for (Sz (90 :. 190))
+-- *** Exception: IndexOutOfBoundsException: (150 :. 150) is not safe for (Sz (90 :. 190))
 --
 -- @since 0.3.0
 evaluate' :: Source r ix e => Array r ix e -> ix -> e
@@ -697,39 +754,6 @@
     (unsafeIndex arr)
     ix
 {-# INLINE evaluate' #-}
-
--- | See `evaluate'`.
---
--- @since 0.1.0
-evaluateAt :: Source r ix e => Array r ix e -> ix -> e
-evaluateAt = evaluate'
-{-# INLINE evaluateAt #-}
-{-# DEPRECATED evaluateAt "In favor of a safe `evaluateM` or an equivalent `evaluate'`" #-}
-
-
--- | This is only used together with the @unsafe-checks@ cabal flag
-indexWith ::
-     Index ix
-  => String -- ^ Source file name, eg. __FILE__
-  -> Int -- ^ Line number in th source file, eg. __LINE__
-  -> String
-  -> (arr -> Sz ix) -- ^ Get size of the array
-  -> (arr -> ix -> e) -- ^ Indexing function
-  -> arr -- ^ Array
-  -> ix -- ^ Index
-  -> e
-indexWith fileName lineNo funName getSize' f arr ix
-  | isSafeIndex (getSize' arr) ix = f arr ix
-  | otherwise = errorIx ("<" ++ fileName ++ ":" ++ show lineNo ++ "> " ++ funName) (getSize' arr) ix
-{-# NOINLINE indexWith #-}
-
--- | Helper function for throwing out of bounds error. Used by `indexWith`
-errorIx :: (Show ix, Show ix') => String -> ix -> ix' -> a
-errorIx fName sz ix =
-  error $
-  fName ++
-  ": Index out of bounds: (" ++ show ix ++ ") for Array of size: (" ++ show sz ++ ")"
-{-# NOINLINE errorIx #-}
 
 
 -- | Map a monadic index aware function over an array sequentially, while discarding the result.
diff --git a/src/Data/Massiv/Core/Exception.hs b/src/Data/Massiv/Core/Exception.hs
--- a/src/Data/Massiv/Core/Exception.hs
+++ b/src/Data/Massiv/Core/Exception.hs
@@ -1,4 +1,6 @@
+{-# LANGUAGE CPP #-}
 {-# LANGUAGE GADTs #-}
+{-# OPTIONS_GHC -fno-warn-orphans #-}
 module Data.Massiv.Core.Exception
   ( ImpossibleException(..)
   , throwImpossible
@@ -10,6 +12,8 @@
 import Control.Monad
 import Control.Monad.Catch
 import Data.Massiv.Core.Index.Internal
+import Control.Monad.ST (ST)
+import Control.Monad.ST.Unsafe (unsafeIOToST)
 
 newtype ImpossibleException =
   ImpossibleException SomeException
@@ -42,3 +46,9 @@
 guardNumberOfElements sz sz' =
   unless (totalElem sz == totalElem sz') $ throwM $ SizeElementsMismatchException sz sz'
 {-# INLINE guardNumberOfElements #-}
+
+#if !MIN_VERSION_exceptions(0, 10, 3)
+-- | Orphan instance in "massiv"
+instance MonadThrow (ST s) where
+  throwM = unsafeIOToST . throwIO
+#endif
diff --git a/src/Data/Massiv/Core/Index.hs b/src/Data/Massiv/Core/Index.hs
--- a/src/Data/Massiv/Core/Index.hs
+++ b/src/Data/Massiv/Core/Index.hs
@@ -31,6 +31,7 @@
   , unSz
   , zeroSz
   , oneSz
+  , liftSz
   , consSz
   , unconsSz
   , snocSz
@@ -38,19 +39,23 @@
   , setSzM
   , insertSzM
   , pullOutSzM
+  -- ** Dimension
   , Dim(..)
   , Dimension(Dim1, Dim2, Dim3, Dim4, Dim5, DimN)
   , IsIndexDimension
-  , Lower
-  , Index(..)
+  -- ** Stride
   , Stride(Stride)
   , unStride
   , toLinearIndexStride
   , strideStart
   , strideSize
   , oneStride
+  -- ** Border
   , Border(..)
   , handleBorderIndex
+  -- ** Index functions
+  , Lower
+  , Index(..)
   , zeroIndex
   , oneIndex
   , isNonEmpty
@@ -59,16 +64,11 @@
   , lastDim
   , initDim
   , getDim'
-  , getDim
   , setDim'
-  , setDim
   , dropDimM
   , dropDim'
-  , dropDim
   , pullOutDim'
-  , pullOutDim
   , insertDim'
-  , insertDim
   , fromDimension
   , getDimension
   , setDimension
@@ -86,6 +86,7 @@
   , SizeException(..)
   , ShapeException(..)
   , guardNumberOfElements
+  , indexWith
   ) where
 
 import Control.DeepSeq
@@ -216,7 +217,7 @@
 {-# INLINE [1] isNonEmpty #-}
 -- TODO: benchmark against (also adjust `isEmpty` with fastest):
 -- - foldlIndex (*) 1 (unSz sz) /= 0
--- - foldlIndex ((&&) . (==0)) True (unSz sz)
+-- - foldlIndex (\a x -> a && x /= 0) True (unSz sz)
 -- - totalElem sz == 0
 
 -- | Get the outmost dimension of the index.
@@ -280,14 +281,6 @@
 setDim' ix dim = either throw id . setDimM ix dim
 {-# INLINE [1] setDim' #-}
 
--- | See `setDimM`
---
--- @since 0.2.4
-setDim :: Index ix => ix -> Dim -> Int -> Maybe ix
-setDim = setDimM
-{-# INLINE [1] setDim #-}
-{-# DEPRECATED setDim "In favor of more general `setDimM`" #-}
-
 -- | Change the value from a specific dimension within the index. Throws `IndexException`. See
 -- `getDimM` for a safer version and `getDimension` for a type safe version.
 --
@@ -296,21 +289,13 @@
 -- >>> getDim' (2 :> 3 :> 4 :. 5) 3
 -- 3
 -- >>> getDim' (2 :> 3 :> 4 :. 5) 0
--- *** Exception: IndexDimensionException: (Dim 0) for 3 :> 4 :. 5
+-- *** Exception: IndexDimensionException: (Dim 0) for (2 :> 3 :> 4 :. 5)
 --
 -- @since 0.2.4
 getDim' :: Index ix => ix -> Dim -> Int
 getDim' ix = either throw id . getDimM ix
 {-# INLINE [1] getDim' #-}
 
--- | See `getDimM`
---
--- @since 0.2.4
-getDim :: Index ix => ix -> Dim -> Maybe Int
-getDim = getDimM
-{-# INLINE [1] getDim #-}
-{-# DEPRECATED getDim "In favor of more general `getDimM`" #-}
-
 -- | Remove a dimension from the index.
 --
 -- ==== __Examples__
@@ -325,14 +310,6 @@
 dropDimM ix = fmap snd . pullOutDimM ix
 {-# INLINE [1] dropDimM #-}
 
--- | See `dropDimM`
---
--- @since 0.1.0
-dropDim :: Index ix => ix -> Dim -> Maybe (Lower ix)
-dropDim = dropDimM
-{-# INLINE [1] dropDim #-}
-{-# DEPRECATED dropDim "In favor of more general `dropDimM`" #-}
-
 -- | Remove a dimension from the index.
 --
 -- ==== __Examples__
@@ -340,7 +317,7 @@
 -- >>> dropDim' (2 :> 3 :> 4 :. 5) 3
 -- 2 :> 4 :. 5
 -- >>> dropDim' (2 :> 3 :> 4 :. 5) 6
--- *** Exception: IndexDimensionException: (Dim 6) for 3 :> 4 :. 5
+-- *** Exception: IndexDimensionException: (Dim 6) for (2 :> 3 :> 4 :. 5)
 --
 -- @since 0.2.4
 dropDim' :: Index ix => ix -> Dim -> Lower ix
@@ -360,14 +337,6 @@
 pullOutDim' ix = either throw id . pullOutDimM ix
 {-# INLINE [1] pullOutDim' #-}
 
--- | See `pullOutDimM`
---
--- @since 0.2.4
-pullOutDim :: Index ix => ix -> Dim -> Maybe (Int, Lower ix)
-pullOutDim = pullOutDimM
-{-# INLINE [1] pullOutDim #-}
-{-# DEPRECATED pullOutDim "In favor of more general `pullOutDimM`" #-}
-
 -- | Raise the dimension of the index by inserting one in the specified dimension. Throws
 -- `IndexException`. See `insertDimM` for a safer version and `insertDimension` for a type safe
 -- version.
@@ -377,21 +346,13 @@
 -- >>> insertDim' (2 :> 3 :> 4 :. 5) 3 10 :: Ix5
 -- 2 :> 3 :> 10 :> 4 :. 5
 -- >>> insertDim' (2 :> 3 :> 4 :. 5) 11 10 :: Ix5
--- *** Exception: IndexDimensionException: (Dim 11) for 4 :. 5
+-- *** Exception: IndexDimensionException: (Dim 11) for (2 :> 3 :> 4 :. 5)
 --
 -- @since 0.2.4
 insertDim' :: Index ix => Lower ix -> Dim -> Int -> ix
 insertDim' ix dim = either throw id . insertDimM ix dim
 {-# INLINE [1] insertDim' #-}
 
--- | See `insertDimM`
---
--- @since 0.2.4
-insertDim :: Index ix => Lower ix -> Dim -> Int -> Maybe ix
-insertDim = insertDimM
-{-# INLINE [1] insertDim #-}
-{-# DEPRECATED insertDim "In favor of more general `insertDimM`" #-}
-
 -- | Get the value level `Dim` from the type level equivalent.
 --
 -- ==== __Examples__
@@ -416,7 +377,7 @@
 --
 -- @since 0.2.4
 setDimension :: IsIndexDimension ix n => ix -> Dimension n -> Int -> ix
-setDimension ix d = setDim' ix (fromDimension d)
+setDimension ix = setDim' ix . fromDimension
 {-# INLINE [1] setDimension #-}
 
 -- | Type safe way to extract value of index at a particular dimension.
@@ -428,7 +389,7 @@
 --
 -- @since 0.2.4
 getDimension :: IsIndexDimension ix n => ix -> Dimension n -> Int
-getDimension ix d = getDim' ix (fromDimension d)
+getDimension ix = getDim' ix . fromDimension
 {-# INLINE [1] getDimension #-}
 
 
@@ -442,7 +403,7 @@
 --
 -- @since 0.2.4
 dropDimension :: IsIndexDimension ix n => ix -> Dimension n -> Lower ix
-dropDimension ix d = dropDim' ix (fromDimension d)
+dropDimension ix = dropDim' ix . fromDimension
 {-# INLINE [1] dropDimension #-}
 
 -- | Type safe way of pulling out a particular dimension, thus lowering index
@@ -455,7 +416,7 @@
 --
 -- @since 0.2.4
 pullOutDimension :: IsIndexDimension ix n => ix -> Dimension n -> (Int, Lower ix)
-pullOutDimension ix d = pullOutDim' ix (fromDimension d)
+pullOutDimension ix = pullOutDim' ix . fromDimension
 {-# INLINE [1] pullOutDimension #-}
 
 -- | Type safe way of inserting a particular dimension, thus raising index dimensionality.
@@ -475,7 +436,7 @@
 --
 -- @since 0.2.5
 insertDimension :: IsIndexDimension ix n => Lower ix -> Dimension n -> Int -> ix
-insertDimension ix d = insertDim' ix (fromDimension d)
+insertDimension ix = insertDim' ix . fromDimension
 {-# INLINE [1] insertDimension #-}
 
 -- | Row-major iterator for the index. Same as `iterM`, but pure.
@@ -516,10 +477,10 @@
 -- @since 0.1.0
 iterLinearM :: (Index ix, Monad m)
             => Sz ix -- ^ Size
-            -> Int -- ^ Linear start
-            -> Int -- ^ Linear end
-            -> Int -- ^ Increment
-            -> (Int -> Int -> Bool) -- ^ Continuation condition (continue if True)
+            -> Int -- ^ Linear start (must be non-negative)
+            -> Int -- ^ Linear end (must be less than or equal to @`totalElem` sz@)
+            -> Int -- ^ Increment (must not be zero)
+            -> (Int -> Int -> Bool) -- ^ Continuation condition (continue if @True@)
             -> a -- ^ Accumulator
             -> (Int -> ix -> a -> m a)
             -> m a
@@ -540,12 +501,42 @@
 -- @since 0.1.0
 iterLinearM_ :: (Index ix, Monad m) =>
                 Sz ix -- ^ Size
-             -> Int -- ^ Start
+             -> Int -- ^ Start (must be non-negative)
              -> Int -- ^ End
-             -> Int -- ^ Increment
-             -> (Int -> Int -> Bool) -- ^ Continuation condition
+             -> Int -- ^ Increment (must not be zero)
+             -> (Int -> Int -> Bool) -- ^ Continuation condition (continue if @True@)
              -> (Int -> ix -> m ()) -- ^ Monadic action that takes index in both forms
              -> m ()
 iterLinearM_ sz !k0 !k1 !inc cond f =
   loopM_ k0 (`cond` k1) (+ inc) $ \ !i -> f i (fromLinearIndex sz i)
 {-# INLINE iterLinearM_ #-}
+
+
+-- | This is used by @INDEX_CHECK@ macro and thus used whenever the @unsafe-checks@ cabal
+-- flag is on.
+--
+-- @since 0.4.0
+indexWith ::
+     Index ix
+  => String -- ^ Source file name, eg. __FILE__
+  -> Int -- ^ Line number in th source file, eg. __LINE__
+  -> String
+  -> (arr -> Sz ix) -- ^ Get size of the array
+  -> (arr -> ix -> e) -- ^ Indexing function
+  -> arr -- ^ Array
+  -> ix -- ^ Index
+  -> e
+indexWith fileName lineNo funName getSize f arr ix
+  | isSafeIndex sz ix = f arr ix
+  | otherwise = errorIx ("<" ++ fileName ++ ":" ++ show lineNo ++ "> " ++ funName) sz ix
+  where
+    sz = getSize arr
+
+-- | Helper function for throwing out of bounds error. Used by `indexWith`
+errorIx :: (Show ix, Show ix') => String -> ix -> ix' -> a
+errorIx fName sz ix =
+  error $
+  fName ++
+  ": Index out of bounds: (" ++ show ix ++ ") for Array of size: (" ++ show sz ++ ")"
+{-# NOINLINE errorIx #-}
+
diff --git a/src/Data/Massiv/Core/Index/Internal.hs b/src/Data/Massiv/Core/Index/Internal.hs
--- a/src/Data/Massiv/Core/Index/Internal.hs
+++ b/src/Data/Massiv/Core/Index/Internal.hs
@@ -32,6 +32,7 @@
   , unSz
   , zeroSz
   , oneSz
+  , liftSz
   , consSz
   , unconsSz
   , snocSz
@@ -55,6 +56,7 @@
   , IndexException(..)
   , SizeException(..)
   , ShapeException(..)
+  , showsPrecWrapped
   ) where
 
 import Control.DeepSeq
@@ -102,16 +104,13 @@
 
 
 instance Index ix => Show (Sz ix) where
-  showsPrec n sz@(SafeSz usz) s =
-    if n == 0
-      then str ++ s
-      else '(' : str ++ ')' : s
+  showsPrec n sz@(SafeSz usz) = showsPrecWrapped n (str ++)
     where
       str =
         "Sz" ++
         case unDim (dimensions sz) of
           1 -> "1 " ++ show usz
-          _ -> " (" ++ show usz ++ ")"
+          _ -> " (" ++ shows usz ")"
 
 instance (Num ix, Index ix) => Num (Sz ix) where
   (+) x y = Sz (coerce x + coerce y)
@@ -170,6 +169,20 @@
 {-# INLINE oneSz #-}
 
 
+-- | Same as `liftIndex`, but for `Sz`
+--
+-- ==== __Example__
+--
+-- >>> import Data.Massiv.Core.Index
+-- >>> liftSz succ (Sz2 2 3)
+-- Sz (3 :. 4)
+--
+-- @since 0.4.0
+liftSz :: Index ix => (Int -> Int) -> Sz ix -> Sz ix
+liftSz f (SafeSz ix) = Sz (liftIndex f ix)
+{-# INLINE liftSz #-}
+
+
 -- | Same as `consDim`, but for `Sz`
 --
 -- ==== __Example__
@@ -205,7 +218,7 @@
 -- >>> setSzM (Sz2 2 3) 2 (Sz1 1) :: IO Sz2
 -- Sz (1 :. 3)
 -- >>> setSzM (Sz2 2 3) 3 (Sz1 1) :: IO Sz2
--- *** Exception: IndexDimensionException: (Dim 3) for 2 :. 3
+-- *** Exception: IndexDimensionException: (Dim 3) for (2 :. 3)
 --
 -- @since 0.3.0
 setSzM :: (MonadThrow m, Index ix) => Sz ix -> Dim -> Sz Int -> m (Sz ix)
@@ -220,7 +233,7 @@
 -- >>> insertSzM (Sz2 2 3) 3 (Sz1 1) :: IO Sz3
 -- Sz (1 :> 2 :. 3)
 -- >>> insertSzM (Sz2 2 3) 4 (Sz1 1) :: IO Sz3
--- *** Exception: IndexDimensionException: (Dim 4) for 2 :. 3
+-- *** Exception: IndexDimensionException: (Dim 4) for (2 :. 3)
 --
 -- @since 0.3.0
 insertSzM :: (MonadThrow m, Index ix) => Sz (Lower ix) -> Dim -> Sz Int -> m (Sz ix)
@@ -259,7 +272,7 @@
 -- >>> pullOutSzM (Sz3 1 2 3) 3
 -- (Sz1 1,Sz (2 :. 3))
 -- >>> pullOutSzM (Sz3 1 2 3) 0
--- *** Exception: IndexDimensionException: (Dim 0) for 1 :> 2 :. 3
+-- *** Exception: IndexDimensionException: (Dim 0) for (1 :> 2 :. 3)
 --
 -- @since 0.3.0
 pullOutSzM :: (MonadThrow m, Index ix) => Sz ix -> Dim -> m (Sz Ix1, Sz (Lower ix))
@@ -336,6 +349,8 @@
       , Ord (Lower ix)
       , Show (Lower ix)
       , NFData (Lower ix)
+      , 1 <= Dimensions ix
+      , KnownNat (Dimensions ix)
       ) =>
       Index ix
   where
@@ -566,6 +581,7 @@
 -- @since 0.1.0
 pattern Ix1 :: Int -> Ix1
 pattern Ix1 i = i
+{-# COMPLETE Ix1 #-}
 
 type instance Lower Int = Ix0
 
@@ -637,13 +653,12 @@
   IndexOutOfBoundsException :: Index ix => !(Sz ix) -> !ix -> IndexException
 
 instance Show IndexException where
-  show (IndexZeroException ix) = "IndexZeroException: " ++ show ix
+  show (IndexZeroException ix) = "IndexZeroException: " ++ showsPrec 1 ix ""
   show (IndexDimensionException ix dim) =
-    "IndexDimensionException: " ++ show dim ++ " for " ++ show ix
+    "IndexDimensionException: " ++ showsPrec 1 dim " for " ++ showsPrec 1 ix ""
   show (IndexOutOfBoundsException sz ix) =
-    "IndexOutOfBoundsException: " ++ showsPrec 1 ix " not safe for (" ++ show sz ++ ")"
-  showsPrec 0 arr s = show arr ++ s
-  showsPrec _ arr s = '(' : show arr ++ ")" ++ s
+    "IndexOutOfBoundsException: " ++ showsPrec 1 ix " is not safe for " ++ showsPrec 1 sz ""
+  showsPrec n exc = showsPrecWrapped n (show exc ++)
 
 instance Eq IndexException where
   e1 == e2 =
@@ -695,8 +710,7 @@
     show sz' ++ ") is to small for " ++ show ix ++ " (" ++ show sz ++ ")"
   show (SizeEmptyException sz) =
     "SizeEmptyException: (" ++ show sz ++ ") corresponds to an empty array"
-  showsPrec 0 arr s = show arr ++ s
-  showsPrec _ arr s = '(' : show arr ++ ")" ++ s
+  showsPrec n exc = showsPrecWrapped n (show exc ++)
 
 -- | Exception that can happen upon conversion of a ragged type array into the rectangular kind. Which
 -- means conversion from lists is susceptible to this exception.
@@ -708,11 +722,16 @@
   deriving Eq
 
 instance Show ShapeException where
-  show (DimTooShortException sz sz') =
-    "DimTooShortException: expected (" ++ show sz ++ "), got (" ++ show sz' ++ ")"
-  show DimTooLongException =
-    "DimTooLongException"
-  showsPrec 0 arr s = show arr ++ s
-  showsPrec _ arr s = '(' : show arr ++ ")" ++ s
+  showsPrec _ DimTooLongException = ("DimTooLongException" ++)
+  showsPrec n (DimTooShortException sz sz') =
+    showsPrecWrapped
+      n
+      (("DimTooShortException: expected (" ++) . shows sz . ("), got (" ++) . shows sz' . (")" ++))
 
 instance Exception ShapeException
+
+
+showsPrecWrapped :: Int -> ShowS -> ShowS
+showsPrecWrapped n inner
+  | n < 1 = inner
+  | otherwise = ('(':) . inner . (")" ++)
diff --git a/src/Data/Massiv/Core/Index/Ix.hs b/src/Data/Massiv/Core/Index/Ix.hs
--- a/src/Data/Massiv/Core/Index/Ix.hs
+++ b/src/Data/Massiv/Core/Index/Ix.hs
@@ -69,6 +69,7 @@
 -- @since 0.1.0
 pattern Ix2 :: Int -> Int -> Ix2
 pattern Ix2 i2 i1 = i2 :. i1
+{-# COMPLETE Ix2 #-}
 
 -- | 2-dimensional size type synonym.
 --
@@ -93,6 +94,7 @@
 -- @since 0.1.0
 pattern Ix3 :: Int -> Int -> Int -> Ix3
 pattern Ix3 i3 i2 i1 = i3 :> i2 :. i1
+{-# COMPLETE Ix3 #-}
 
 -- | 3-dimensional size type synonym.
 --
@@ -116,6 +118,7 @@
 -- @since 0.1.0
 pattern Ix4 :: Int -> Int -> Int -> Int -> Ix4
 pattern Ix4 i4 i3 i2 i1 = i4 :> i3 :> i2 :. i1
+{-# COMPLETE Ix4 #-}
 
 -- | 4-dimensional size type synonym.
 --
@@ -139,6 +142,7 @@
 -- @since 0.1.0
 pattern Ix5 :: Int -> Int -> Int -> Int -> Int -> Ix5
 pattern Ix5 i5 i4 i3 i2 i1 = i5 :> i4 :> i3 :> i2 :. i1
+{-# COMPLETE Ix5 #-}
 
 -- | 5-dimensional size type synonym.
 --
@@ -172,20 +176,10 @@
 
 
 instance Show Ix2 where
-  showsPrec n (i :. j) =
-    if n < 1
-      then inner
-      else ("(" ++) . inner . (")" ++)
-    where
-      inner = shows i . (" :. " ++) . shows j
+  showsPrec n (i :. j) = showsPrecWrapped n (shows i . (" :. " ++) . shows j)
 
 instance Show (Ix (n - 1)) => Show (IxN n) where
-  showsPrec n (i :> ix) =
-    if n < 1
-      then inner
-      else ("(" ++) . inner . (")" ++)
-    where
-      inner = shows i . (" :> " ++) . shows ix
+  showsPrec n (i :> ix) = showsPrecWrapped n (shows i . (" :> " ++) . shows ix)
 
 
 instance Num Ix2 where
@@ -221,11 +215,14 @@
   {-# INLINE [1] fromInteger #-}
 
 
-instance {-# OVERLAPPABLE #-} (4 <= n,
-          KnownNat n,
-          Index (Ix (n - 1)),
-          IxN (n - 1) ~ Ix (n - 1)
-          ) => Num (IxN n) where
+instance {-# OVERLAPPABLE #-} ( 1 <= n
+                              , 4 <= n
+                              , KnownNat n
+                              , KnownNat (n - 1)
+                              , Index (Ix (n - 1))
+                              , IxN (n - 1) ~ Ix (n - 1)
+                              ) =>
+                              Num (IxN n) where
   (+) = liftIndex2 (+)
   {-# INLINE [1] (+) #-}
   (-) = liftIndex2 (-)
@@ -255,11 +252,14 @@
   maxBound = pureIndex maxBound
   {-# INLINE maxBound #-}
 
-instance {-# OVERLAPPABLE #-} (4 <= n,
-          KnownNat n,
-          Index (Ix (n - 1)),
-          IxN (n - 1) ~ Ix (n - 1)
-          ) => Bounded (IxN n) where
+instance {-# OVERLAPPABLE #-} ( 1 <= n
+                              , 4 <= n
+                              , KnownNat n
+                              , KnownNat (n - 1)
+                              , Index (Ix (n - 1))
+                              , IxN (n - 1) ~ Ix (n - 1)
+                              ) =>
+                              Bounded (IxN n) where
   minBound = pureIndex minBound
   {-# INLINE minBound #-}
   maxBound = pureIndex maxBound
@@ -385,11 +385,14 @@
     repairIndex (SafeSz n) i rBelow rOver :> repairIndex (SafeSz szL) ixL rBelow rOver
   {-# INLINE [1] repairIndex #-}
 
-instance {-# OVERLAPPABLE #-} (4 <= n,
-          KnownNat n,
-          Index (Ix (n - 1)),
-          IxN (n - 1) ~ Ix (n - 1)
-          ) => Index (IxN n) where
+instance {-# OVERLAPPABLE #-} ( 1 <= n
+                              , 4 <= n
+                              , KnownNat n
+                              , KnownNat (n - 1)
+                              , Index (Ix (n - 1))
+                              , IxN (n - 1) ~ Ix (n - 1)
+                              ) =>
+                              Index (IxN n) where
   type Dimensions (IxN n) = n
   dimensions _ = fromInteger $ natVal (Proxy :: Proxy n)
   {-# INLINE [1] dimensions #-}
@@ -401,20 +404,27 @@
   {-# INLINE [1] unconsDim #-}
   snocDim (i :> ixl) i1 = i :> snocDim ixl i1
   {-# INLINE [1] snocDim #-}
-  unsnocDim (i :> ixl) = case unsnocDim ixl of
-                          (ix, i1) -> (i :> ix, i1)
+  unsnocDim (i :> ixl) =
+    case unsnocDim ixl of
+      (ix, i1) -> (i :> ix, i1)
   {-# INLINE [1] unsnocDim #-}
-  getDimM ix@(i :> ixl) d | d == dimensions (Just ix) = pure i
-                         | otherwise = getDimM ixl d
+  getDimM ix@(i :> ixl) d
+    | d == dimensions (Proxy :: Proxy (IxN n)) = pure i
+    | otherwise = maybe (throwM $ IndexDimensionException ix d) pure (getDimM ixl d)
   {-# INLINE [1] getDimM #-}
-  setDimM ix@(i :> ixl) d di | d == dimensions (Just ix) = pure (di :> ixl)
-                            | otherwise = (i :>) <$> setDimM ixl d di
+  setDimM ix@(i :> ixl) d di
+    | d == dimensions (Proxy :: Proxy (IxN n)) = pure (di :> ixl)
+    | otherwise = maybe (throwM $ IndexDimensionException ix d) (pure . (i :>)) (setDimM ixl d di)
   {-# INLINE [1] setDimM #-}
-  pullOutDimM ix@(i :> ixl) d | d == dimensions (Just ix) = pure (i, ixl)
-                             | otherwise = fmap (i :>) <$> pullOutDimM ixl d
+  pullOutDimM ix@(i :> ixl) d
+    | d == dimensions (Proxy :: Proxy (IxN n)) = pure (i, ixl)
+    | otherwise =
+      maybe (throwM $ IndexDimensionException ix d) (pure . fmap (i :>)) (pullOutDimM ixl d)
   {-# INLINE [1] pullOutDimM #-}
-  insertDimM ix@(i :> ixl) d di | d == dimensions (Just ix) + 1 = pure (di :> ix)
-                               | otherwise = (i :>) <$> insertDimM ixl d di
+  insertDimM ix@(i :> ixl) d di
+    | d == dimensions (Proxy :: Proxy (IxN n)) = pure (di :> ix)
+    | otherwise =
+      maybe (throwM $ IndexDimensionException ix d) (pure . (i :>)) (insertDimM ixl d di)
   {-# INLINE [1] insertDimM #-}
   pureIndex i = i :> (pureIndex i :: Ix (n - 1))
   {-# INLINE [1] pureIndex #-}
diff --git a/src/Data/Massiv/Core/Index/Stride.hs b/src/Data/Massiv/Core/Index/Stride.hs
--- a/src/Data/Massiv/Core/Index/Stride.hs
+++ b/src/Data/Massiv/Core/Index/Stride.hs
@@ -51,12 +51,14 @@
 -- * In case of two dimensions, if what you want is to keep all rows divisible by 5, but keep every
 --   column intact then you'd use @Stride (5 :. 1)@.
 --
-
+-- @since 0.2.1
 newtype Stride ix = SafeStride ix deriving (Eq, Ord, NFData)
 
 
 -- | A safe bidirectional pattern synonym for `Stride` construction that will make sure stride
 -- elements are always positive.
+--
+-- @since 0.2.1
 pattern Stride :: Index ix => ix -> Stride ix
 pattern Stride ix <- SafeStride ix where
         Stride ix = SafeStride (liftIndex (max 1) ix)
@@ -64,15 +66,19 @@
 
 
 instance Index ix => Show (Stride ix) where
-  show (SafeStride ix) = "Stride (" ++ show ix ++ ")"
+  showsPrec n (SafeStride ix) = showsPrecWrapped n (("Stride " ++) . showsPrec 1 ix)
 
 
 -- | Just a helper function for unwrapping `Stride`.
+--
+-- @since 0.2.1
 unStride :: Stride ix -> ix
 unStride (SafeStride ix) = ix
 {-# INLINE unStride #-}
 
 -- | Adjust starting index according to the stride
+--
+-- @since 0.2.1
 strideStart :: Index ix => Stride ix -> ix -> ix
 strideStart (SafeStride stride) ix =
   liftIndex2
@@ -82,14 +88,19 @@
 {-# INLINE strideStart #-}
 
 -- | Adjust size according to the stride.
+--
+-- @since 0.2.1
 strideSize :: Index ix => Stride ix -> Sz ix -> Sz ix
 strideSize (SafeStride stride) (SafeSz sz) =
   SafeSz (liftIndex (+ 1) $ liftIndex2 div (liftIndex (subtract 1) sz) stride)
 {-# INLINE strideSize #-}
 
--- | Compute an index with stride using the original size and index
-toLinearIndexStride :: Index ix =>
-  Stride ix -- ^ Stride
+-- | Compute linear index with stride using the original size and index
+--
+-- @since 0.2.1
+toLinearIndexStride ::
+     Index ix
+  => Stride ix -- ^ Stride
   -> Sz ix -- ^ Size
   -> ix -- ^ Index
   -> Int
@@ -98,6 +109,8 @@
 
 
 -- | A default stride of @1@, where all elements are kept
+--
+-- @since 0.2.1
 oneStride :: Index ix => Stride ix
 oneStride = SafeStride (pureIndex 1)
 {-# INLINE oneStride #-}
diff --git a/src/Data/Massiv/Core/List.hs b/src/Data/Massiv/Core/List.hs
--- a/src/Data/Massiv/Core/List.hs
+++ b/src/Data/Massiv/Core/List.hs
@@ -33,15 +33,12 @@
 import Data.Foldable (foldr')
 import qualified Data.List as L
 import Data.Massiv.Core.Common
-import Data.Proxy
 import Data.Typeable
 import GHC.Exts
 import System.IO.Unsafe (unsafePerformIO)
 
 data LN
 
-type instance EltRepr LN ix = LN
-
 type family ListItem ix e :: * where
   ListItem Ix1 e = e
   ListItem ix  e = [ListItem (Lower ix) e]
@@ -77,7 +74,6 @@
 
 
 data L = L
-type instance EltRepr L ix = L
 
 type instance NestedStruct L ix e = Array LN ix e
 
@@ -144,7 +140,20 @@
   {-# INLINE loadArrayM #-}
 
 
+instance (Index ix, Load L ix e, Ragged L ix e) => Load LN ix e where
+  size = edgeSize . LArray Seq
+  {-# INLINE size #-}
+  getComp _ = Seq
+  {-# INLINE getComp #-}
+  loadArrayM scheduler arr uWrite =
+    loadRagged (scheduleWork scheduler) uWrite 0 (totalElem sz) sz arrL
+    where
+      !arrL = LArray Seq arr
+      !sz = size arrL
+  {-# INLINE loadArrayM #-}
 
+
+
 outerLength :: Array L ix e -> Sz Int
 outerLength = SafeSz . length . unList . lData
 
@@ -272,6 +281,9 @@
 {-# INLINE unsafeGenerateN #-}
 
 
+-- | Construct an array backed by linked lists from any source array
+--
+-- @since 0.4.0
 toListArray :: (Construct L ix e, Source r ix e)
             => Array r ix e
             -> Array L ix e
@@ -315,7 +327,9 @@
         else (('(':), ("\n)" ++))
     lnarr = toNested arr
 
-
+-- | Helper function for declaring `Show` instances for arrays
+--
+-- @since 0.4.0
 showsArrayPrec ::
      forall r r' ix ix' e. (Ragged L ix' e, Load r ix e, Source r' ix' e, Show e)
   => (Array r ix e -> Array r' ix' e) -- ^ Modifier
@@ -328,6 +342,9 @@
     larr = makeArray (getComp arr') (size arr') (evaluate' arr') :: Array L ix' e
 
 
+-- | Helper function for declaring `Show` instances for arrays
+--
+-- @since 0.4.0
 showArrayList
   :: Show arr => [arr] -> String -> String
 showArrayList arrs = ('[':) . go arrs . (']':)
diff --git a/src/Data/Massiv/Core/Operations.hs b/src/Data/Massiv/Core/Operations.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Massiv/Core/Operations.hs
@@ -0,0 +1,130 @@
+{-# LANGUAGE DefaultSignatures #-}
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE MultiParamTypeClasses #-}
+{-# LANGUAGE ScopedTypeVariables #-}
+{-# LANGUAGE TypeFamilies #-}
+-- |
+-- Module      : Data.Massiv.Core.Operations
+-- Copyright   : (c) Alexey Kuleshevich 2018-2019
+-- License     : BSD3
+-- Maintainer  : Alexey Kuleshevich <lehins@yandex.ru>
+-- Stability   : experimental
+-- Portability : non-portable
+module Data.Massiv.Core.Operations
+  ( Numeric(..)
+  , NumericFloat(..)
+  ) where
+
+import Data.Massiv.Core.Common
+-- import Data.Massiv.Array.Ops.Fold.Internal
+
+
+class Num e => Numeric r e where
+
+  {-# MINIMAL unsafeLiftArray, unsafeLiftArray2 #-}
+
+  -- sumArray :: Array r Ix1 e -> e
+  -- default sumArray :: Source r Ix1 e => Array r Ix1 e -> e
+  -- sumArray = foldlS (+) 0
+  -- {-# INLINE sumArray #-}
+
+  -- productArray :: Array r Ix1 e -> e
+  -- default productArray :: Source r Ix1 e => Array r Ix1 e -> e
+  -- productArray = foldlS (*) 1
+  -- {-# INLINE productArray #-}
+
+  -- -- | Raise each element in the array to some non-negative power and sum the results
+  -- powerSumArray :: Array r Ix1 e -> Int -> e
+
+  -- unsafeDotProduct :: Array r Ix1 e -> Array r Ix1 e -> e
+
+  plusScalar :: Index ix => Array r ix e -> e -> Array r ix e
+  plusScalar arr e = unsafeLiftArray (+ e) arr
+  {-# INLINE plusScalar #-}
+
+  minusScalar :: Index ix => Array r ix e -> e -> Array r ix e
+  minusScalar arr e = unsafeLiftArray (subtract e) arr
+  {-# INLINE minusScalar #-}
+
+  multiplyScalar :: Index ix => Array r ix e -> e -> Array r ix e
+  multiplyScalar arr e = unsafeLiftArray (* e) arr
+  {-# INLINE multiplyScalar #-}
+
+  absPointwise :: Index ix => Array r ix e -> Array r ix e
+  absPointwise = unsafeLiftArray abs
+  {-# INLINE absPointwise #-}
+
+  additionPointwise :: Index ix => Array r ix e -> Array r ix e -> Array r ix e
+  additionPointwise = unsafeLiftArray2 (+)
+  {-# INLINE additionPointwise #-}
+
+  subtractionPointwise :: Index ix => Array r ix e -> Array r ix e -> Array r ix e
+  subtractionPointwise = unsafeLiftArray2 (-)
+  {-# INLINE subtractionPointwise #-}
+
+  multiplicationPointwise :: Index ix => Array r ix e -> Array r ix e -> Array r ix e
+  multiplicationPointwise = unsafeLiftArray2 (*)
+  {-# INLINE multiplicationPointwise #-}
+
+  -- | Raise each element of the array to the power
+  powerPointwise :: Index ix => Array r ix e -> Int -> Array r ix e
+  powerPointwise arr pow = unsafeLiftArray (^ pow) arr
+  {-# INLINE powerPointwise #-}
+
+
+  unsafeLiftArray :: Index ix => (a -> e) -> Array r ix a -> Array r ix e
+
+  unsafeLiftArray2 :: Index ix => (a -> b -> e) -> Array r ix a -> Array r ix b -> Array r ix e
+
+
+
+class (Numeric r e, Floating e) => NumericFloat r e where
+
+  divideScalar :: Index ix => Array r ix e -> e -> Array r ix e
+  divideScalar arr e = unsafeLiftArray (/ e) arr
+  {-# INLINE divideScalar #-}
+
+  divisionPointwise :: Index ix => Array r ix e -> Array r ix e -> Array r ix e
+  divisionPointwise = unsafeLiftArray2 (/)
+  {-# INLINE divisionPointwise #-}
+
+  recipPointwise :: Index ix => Array r ix e -> Array r ix e
+  recipPointwise = unsafeLiftArray recip
+  {-# INLINE recipPointwise #-}
+
+  sqrtPointwise :: Index ix => Array r ix e -> Array r ix e
+  sqrtPointwise = unsafeLiftArray sqrt
+  {-# INLINE sqrtPointwise #-}
+
+  -- floorPointwise :: (Index ix, Integral a) => Array r ix e -> Array r ix a
+  -- floorPointwise = unsafeLiftArray floor
+  -- {-# INLINE floorPointwise #-}
+
+  -- ceilingPointwise :: (Index ix, Integral a) => Array r ix e -> Array r ix a
+  -- ceilingPointwise = unsafeLiftArray ceiling
+  -- {-# INLINE ceilingPointwise #-}
+
+
+-- class Equality r e where
+
+--   unsafeEq :: Index ix => Array r ix e -> Array r ix e -> Bool
+
+--   unsafeEqPointwise :: Index ix => Array r ix e -> Array r ix e -> Array r ix Bool
+
+
+-- class Relation r e where
+
+--   unsafePointwiseLT :: Array r ix e -> Array r ix e -> Array r ix Bool
+--   unsafePointwiseLTE :: Array r ix e -> Array r ix e -> Array r ix Bool
+
+--   unsafePointwiseGT :: Array r ix e -> Array r ix e -> Array r ix Bool
+--   unsafePointwiseGTE :: Array r ix e -> Array r ix e -> Array r ix Bool
+
+--   unsafePointwiseMin :: Array r ix e -> Array r ix e -> Array r ix e
+--   unsafePointwiseMax :: Array r ix e -> Array r ix e -> Array r ix e
+
+--   unsafeMinimum :: Array r ix e -> e
+
+--   unsafeMaximum :: Array r ix e -> e
+
+
diff --git a/tests/Data/Massiv/Array/Delayed/InterleavedSpec.hs b/tests/Data/Massiv/Array/Delayed/InterleavedSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Delayed/InterleavedSpec.hs
+++ /dev/null
@@ -1,40 +0,0 @@
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.Delayed.InterleavedSpec
-  ( spec
-  ) where
-
-import Data.Massiv.CoreArbitrary as A
-
-
-prop_EqDelayed ::
-     (Ragged L ix Int, Load D ix Int, Load DI ix Int)
-  => Array D ix Int
-  -> Property
-prop_EqDelayed arr = computeAs P arr === computeAs P (toInterleaved arr)
-
-
-prop_Resize ::
-     (Ragged L ix Int, Load D ix Int, Load DI ix Int)
-  => Array DI ix Int
-  -> Property
-prop_Resize arr =
-  computeAs P (resize' k arr) === computeAs P (resize' k arrD)
-  where
-    arrD = fromInterleaved arr
-    k = Sz (totalElem (size arr))
-
-spec :: Spec
-spec =
-  describe "Interleaved same as Delayed" $ do
-    it "EqDelayed Ix1" $ property $ prop_EqDelayed @Ix1
-    it "EqDelayed Ix2" $ property $ prop_EqDelayed @Ix2
-    it "EqDelayed Ix3" $ property $ prop_EqDelayed @Ix3
-    it "EqDelayed Ix4" $ property $ prop_EqDelayed @Ix4
-    it "EqDelayed Ix5" $ property $ prop_EqDelayed @Ix5
-    it "Resize Ix1" $ property $ prop_Resize @Ix1
-    it "Resize Ix2" $ property $ prop_Resize @Ix2
-    it "Resize Ix3" $ property $ prop_Resize @Ix3
-    it "Resize Ix4" $ property $ prop_Resize @Ix4
-    it "Resize Ix5" $ property $ prop_Resize @Ix5
diff --git a/tests/Data/Massiv/Array/Delayed/PushSpec.hs b/tests/Data/Massiv/Array/Delayed/PushSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Delayed/PushSpec.hs
+++ /dev/null
@@ -1,31 +0,0 @@
-{-# LANGUAGE DataKinds #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-module Data.Massiv.Array.Delayed.PushSpec (spec) where
-
-import Data.Massiv.Array.Delayed
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary as A
-import Data.Typeable
-
-
--- prop_upsampleDifferentDefault ::
---      Proxy ix -> Comp -> SzIx ix -> Int -> Maybe Int -> Property
--- prop_upsampleDifferentDefault _ comp (SzIx sz ix) v mDef =
---   computeAs P (unsafeMakeLoadArray comp sz mDef $ \ put -> put ix v)
-
-
-spec :: Spec
-spec = pure ()
-  -- describe "upsampleDifferentDefault" $ do
-  --   it "Ix1" $ property $ prop_upsampleDifferentDefault (Proxy :: Proxy Ix1)
-  --   it "Ix2" $ property $ prop_upsampleDifferentDefault (Proxy :: Proxy Ix2)
-  --   it "Ix3" $ property $ prop_upsampleDifferentDefault (Proxy :: Proxy Ix3)
-  --   it "Ix4" $ property $ prop_upsampleDifferentDefault (Proxy :: Proxy Ix4)
-  --   it "Ix5" $ property $ prop_upsampleDifferentDefault (Proxy :: Proxy Ix5)
-
-
-identityDL :: Int -> Array DL Ix2 Int
-identityDL n = makeLoadArrayS (Sz2 n n) 0 $ \ writeCell -> do
-  let f i = writeCell (i :. i) 1
-  A.mapM_ f (0 ... n - 1)
diff --git a/tests/Data/Massiv/Array/Delayed/WindowedSpec.hs b/tests/Data/Massiv/Array/Delayed/WindowedSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Delayed/WindowedSpec.hs
+++ /dev/null
@@ -1,75 +0,0 @@
-{-# LANGUAGE DataKinds #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE NamedFieldPuns #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE UndecidableInstances #-}
-module Data.Massiv.Array.Delayed.WindowedSpec (spec) where
-
-import Data.Massiv.Array.Delayed
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary as A
-import Data.Typeable
-
-
-data ArrDW ix e = ArrDW (Array D ix e) (Array DW ix e)
-
-instance (Show ix, Index ix, Show (Array D ix e), Show (Array DW ix e)) => Show (ArrDW ix e) where
-  show (ArrDW d dw) =
-    "Delayed:\n" ++
-    show d ++
-    "\nCorresponding Windowed:\n" ++
-    --show dw ++
-    windowInfo
-    where
-      windowInfo =
-        maybe
-          "\n No Window"
-          (\Window {windowStart, windowSize} ->
-             "\n With Window starting index (" ++
-             show windowStart ++ ") and size (" ++ show windowSize ++ ")") $
-        getWindow dw
-
-instance (Arbitrary ix, CoArbitrary ix, Index ix, Arbitrary e, Typeable e) =>
-         Arbitrary (ArrDW ix e) where
-  arbitrary = do
-    ArrTiny (arr :: Array D ix e) <- arbitrary
-    let sz = size arr
-    ArrDW arr <$>
-      if totalElem sz == 0
-        then return (makeArray (getComp arr) sz (unsafeIndex arr))
-        else do
-          wix <- flip (liftIndex2 mod) (unSz sz) <$> arbitrary
-          wsz <- liftIndex (+1) . flip (liftIndex2 mod) (liftIndex2 (-) (unSz sz) wix) <$> arbitrary
-          return $ makeWindowedArray arr wix (Sz wsz) (unsafeIndex arr)
-
-
-prop_EqDelayed ::
-     (Ragged L ix Int, Load DW ix Int) => Proxy ix -> ArrDW ix Int -> Property
-prop_EqDelayed _ (ArrDW arrD arrDW) =
-  computeAs P arrD === computeAs P arrDW
-
-prop_EqDelayedStride ::
-     (Ragged L ix Int, StrideLoad DW ix Int) => Proxy ix -> Stride ix -> ArrDW ix Int -> Property
-prop_EqDelayedStride _ stride (ArrDW arrD arrDW) =
-  computeWithStrideAs P stride arrD === computeWithStrideAs P stride arrDW
-
-
-spec :: Spec
-spec = do
-  describe "Equivalency with Delayed" $ do
-    it "Ix1" $ property $ prop_EqDelayed (Proxy :: Proxy Ix1)
-    it "Ix2" $ property $ prop_EqDelayed (Proxy :: Proxy Ix2)
-    it "Ix3" $ property $ prop_EqDelayed (Proxy :: Proxy Ix3)
-    it "Ix4" $ property $ prop_EqDelayed (Proxy :: Proxy Ix4)
-    it "Ix5" $ property $ prop_EqDelayed (Proxy :: Proxy Ix5)
-    it "Ix2T" $ property $ prop_EqDelayed (Proxy :: Proxy Ix2T)
-  describe "Equivalency with Stride With Delayed" $ do
-    it "Ix1" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix1)
-    it "Ix2" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix2)
-    it "Ix3" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix3)
-    it "Ix4" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix4)
-    it "Ix5" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix5)
-    it "Ix2T" $ property $ prop_EqDelayedStride (Proxy :: Proxy Ix2T)
diff --git a/tests/Data/Massiv/Array/DelayedSpec.hs b/tests/Data/Massiv/Array/DelayedSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/DelayedSpec.hs
+++ /dev/null
@@ -1,66 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-module Data.Massiv.Array.DelayedSpec (spec) where
-
-import Data.Massiv.Array
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary as A
-import Data.Proxy
-
-
-downsampleArr :: Source r ix e => Stride ix -> Array r ix e -> Array D ix e
-downsampleArr stride arr =
-  unsafeBackpermute (strideSize stride (size arr)) (liftIndex2 (*) (unStride stride)) arr
-
-prop_computeWithStrideEqDownsample ::
-     Ragged L ix Int
-  => proxy ix
-  -> Stride ix
-  -> Array D ix Int
-  -> Property
-prop_computeWithStrideEqDownsample _ stride arr =
-  computeWithStride stride arr === computeAs U (downsampleArr stride arr)
-
-
-prop_computeWithStrideInterleavedEqDownsample ::
-     Ragged L ix Int
-  => proxy ix
-  -> Stride ix
-  -> Array D ix Int
-  -> Property
-prop_computeWithStrideInterleavedEqDownsample _ stride arr =
-  computeWithStride stride (toInterleaved arr) === computeAs U (downsampleArr stride arr)
-
-prop_computeWithStrideWindowedEqDownsample ::
-     (Ragged L ix Int, StrideLoad DW ix Int)
-  => proxy ix
-  -> Stride ix
-  -> ArrIx D ix Int
-  -> Property
-prop_computeWithStrideWindowedEqDownsample _ stride (ArrIx arr _) =
-  computeWithStride stride (insertWindow arr (Window zeroIndex (size arr) (unsafeIndex arr) Nothing)) ===
-  -- Below triggers a bug in ghc-8.0 which results in a deadlock.
-  -- computeWithStride stride (makeWindowedArray arr zeroIndex (size arr) (unsafeIndex arr)) ===
-  computeAs U (downsampleArr stride arr)
-
-
-delayedSpec ::
-     (Arbitrary ix, CoArbitrary ix, StrideLoad DW ix Int, Ragged L ix Int)
-  => String
-  -> proxy ix
-  -> Spec
-delayedSpec dimName proxy =
-  describe dimName $ do
-    it "computeWithStrideEqDownsample" $ property $ prop_computeWithStrideEqDownsample proxy
-    it "computeWithStrideInterleavedEqDownsample" $
-      property $ prop_computeWithStrideInterleavedEqDownsample proxy
-    it "computeWithStrideWindowedEqDownsample" $
-      property $ prop_computeWithStrideWindowedEqDownsample proxy
-
-spec :: Spec
-spec = do
-  delayedSpec "Ix1" (Proxy :: Proxy Ix1)
-  delayedSpec "Ix2" (Proxy :: Proxy Ix2)
-  delayedSpec "Ix3" (Proxy :: Proxy Ix3)
-  delayedSpec "Ix4" (Proxy :: Proxy Ix4)
diff --git a/tests/Data/Massiv/Array/Manifest/VectorSpec.hs b/tests/Data/Massiv/Array/Manifest/VectorSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Manifest/VectorSpec.hs
+++ /dev/null
@@ -1,79 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE TypeFamilies #-}
-module Data.Massiv.Array.Manifest.VectorSpec (spec) where
-
-import Data.Massiv.Array.Manifest.Vector
-import Data.Massiv.CoreArbitrary
-import Data.Proxy
-import Data.Typeable
-import qualified Data.Vector as VB
-import qualified Data.Vector.Generic as VG
-import qualified Data.Vector.Primitive as VP
-import qualified Data.Vector.Storable as VS
-import qualified Data.Vector.Unboxed as VU
-
-prop_castToFromVector
-  :: ( VG.Vector (VRepr r) Int
-     , Mutable r ix Int
-     , Typeable (VRepr r)
-     , ARepr (VRepr r) ~ r
-     , Eq (Array r ix Int)
-     , Show (Array r ix Int)
-     )
-  => proxy ix -> r -> Arr r ix Int -> Property
-prop_castToFromVector _ _ (Arr arr) =
-  Just arr === (castToVector arr >>= castFromVector (getComp arr) (size arr))
-
-
-prop_toFromVector ::
-     forall r ix v.
-     ( Construct r ix Int
-     , Mutable r ix Int
-     , Mutable (ARepr v) ix Int
-     , VRepr (ARepr v) ~ v
-     , Eq (Array r ix Int)
-     , VG.Vector v Int
-     , Show (Array r ix Int)
-     , Typeable v
-     )
-  => Proxy v
-  -> Proxy ix
-  -> r
-  -> Arr r ix Int
-  -> Property
-prop_toFromVector _ _ _ (Arr arr) =
-  arr === fromVector' (getComp arr) (size arr) (toVector arr :: v Int)
-
-
-toFromVectorSpec :: Spec
-toFromVectorSpec = do
-  it_prop "Unboxed" U
-  it_prop "Primitive" P
-  it_prop "Storable" S
-  it_prop "BoxedStrict" B
-  where
-    it_prop name r =
-      describe name $ do
-        describe "CastToFrom" $ do
-          it "Ix1" $ property $ prop_castToFromVector (Proxy :: Proxy Ix1) r
-          it "Ix2" $ property $ prop_castToFromVector (Proxy :: Proxy Ix2) r
-          it "Ix3" $ property $ prop_castToFromVector (Proxy :: Proxy Ix3) r
-        describe "Through Boxed Vector" $ do
-          it "Ix1" $ property $ prop_toFromVector (Proxy :: Proxy VB.Vector) (Proxy :: Proxy Ix1) r
-          it "Ix2" $ property $ prop_toFromVector (Proxy :: Proxy VB.Vector) (Proxy :: Proxy Ix2) r
-        describe "Through Unboxed Vector" $ do
-          it "Ix1" $ property $ prop_toFromVector (Proxy :: Proxy VU.Vector) (Proxy :: Proxy Ix1) r
-          it "Ix2" $ property $ prop_toFromVector (Proxy :: Proxy VU.Vector) (Proxy :: Proxy Ix2) r
-        describe "Through Primitive Vector" $ do
-          it "Ix1" $ property $ prop_toFromVector (Proxy :: Proxy VP.Vector) (Proxy :: Proxy Ix1) r
-          it "Ix2" $ property $ prop_toFromVector (Proxy :: Proxy VP.Vector) (Proxy :: Proxy Ix2) r
-        describe "Through Storable Vector" $ do
-          it "Ix1" $ property $ prop_toFromVector (Proxy :: Proxy VS.Vector) (Proxy :: Proxy Ix1) r
-          it "Ix2" $ property $ prop_toFromVector (Proxy :: Proxy VS.Vector) (Proxy :: Proxy Ix2) r
-
-
-spec :: Spec
-spec = describe "toFromVector" toFromVectorSpec
diff --git a/tests/Data/Massiv/Array/ManifestSpec.hs b/tests/Data/Massiv/Array/ManifestSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/ManifestSpec.hs
+++ /dev/null
@@ -1,42 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.ManifestSpec (spec) where
-
-import Data.ByteString as S
-import Data.ByteString.Builder as S
-import Data.ByteString.Lazy as SL
-import Data.Massiv.CoreArbitrary as A
-import Data.Word (Word8)
-
-
--- ByteString
-prop_toFromByteString :: Manifest r Ix1 Word8 => Array r Ix1 Word8 -> Property
-prop_toFromByteString arr = toManifest arr === fromByteString (getComp arr) (toByteString arr)
-
-prop_castToFromByteString :: Array S Ix1 Word8 -> Property
-prop_castToFromByteString arr = arr === castFromByteString (getComp arr) (castToByteString arr)
-
-
-prop_fromToByteString :: Comp -> [Word8] -> Property
-prop_fromToByteString comp ls = bs === toByteString (fromByteString comp bs)
-  where bs = S.pack ls
-
-prop_toBuilder :: Array P Ix1 Word8 -> Property
-prop_toBuilder arr = bs === SL.toStrict (S.toLazyByteString (toBuilder S.word8 arr))
-  where bs = toByteString arr
-
-conversionSpec :: Spec
-conversionSpec =
-  describe "ByteString" $ do
-    it "castTo/TromByteString" $ property prop_castToFromByteString
-    it "to/from ByteString P" $ property (prop_toFromByteString @P)
-    it "to/from ByteString S" $ property (prop_toFromByteString @S)
-    it "from/to ByteString" $ property prop_fromToByteString
-    it "toBuilder" $ property prop_toBuilder
-
-
-spec :: Spec
-spec = describe "Conversion" conversionSpec
diff --git a/tests/Data/Massiv/Array/MutableSpec.hs b/tests/Data/Massiv/Array/MutableSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/MutableSpec.hs
+++ /dev/null
@@ -1,217 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.MutableSpec (spec) where
-
-import Control.Concurrent.Async
-import Control.Monad (when)
-import Control.Monad.ST
-import Data.Functor.Identity
-import Data.List as L
-import Data.Massiv.Array.Mutable.Atomic
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary as A
-import Data.Proxy
-
-
-prop_MapMapM :: (Show (Array r ix Int), Eq (Array r ix Int), Mutable r ix Int) =>
-                r -> Proxy ix -> Fun Int Int -> ArrTiny D ix Int -> Property
-prop_MapMapM r _ f (ArrTiny arr) =
-  computeAs r (A.map (apply f) arr) === runIdentity (A.mapMR r (return . apply f) arr)
-
-prop_iMapiMapM :: (Show (Array r ix Int), Eq (Array r ix Int), Mutable r ix Int) =>
-                r -> Proxy ix -> Fun (ix, Int) Int -> ArrTiny D ix Int -> Property
-prop_iMapiMapM r _ f (ArrTiny arr) =
-  computeAs r (A.imap (curry (apply f)) arr) ===
-  runIdentity (A.imapMR r (\ix e -> return $ apply f (ix, e)) arr)
-
-
-prop_generateMakeST ::
-     (Show (Array r ix Int), Eq (Array r ix Int), Mutable r ix Int)
-  => r
-  -> Proxy ix
-  -> Arr r ix Int
-  -> Property
-prop_generateMakeST _ _ (Arr arr) =
-  arr === runST (generateArrayS (getComp arr) (size arr) (return . evaluate' arr))
-
-prop_generateMakeIO :: (Show (Array r ix Int), Eq (Array r ix Int), Mutable r ix Int) =>
-                             r -> Proxy ix -> Arr r ix Int -> Property
-prop_generateMakeIO _ _ (Arr arr) = monadicIO $ do
-  arr' <- run $ generateArray (getComp arr) (size arr) (evaluateM arr)
-  return (arr === arr')
-
-prop_shrinkIO ::
-     ( Mutable r ix Int
-     , Resize r ix
-     , Source r Ix1 Int
-     )
-  => r
-  -> Proxy ix
-  -> ArrIx r ix Int
-  -> Property
-prop_shrinkIO _ _ (ArrIx arr ix) =
-  monadicIO $
-  run $ do
-    marr <- thaw arr
-    sarr <- unsafeFreeze (getComp arr) =<< unsafeLinearShrink marr (Sz ix)
-    pure (A.foldlS (.&&.) (property True) $ A.zipWith (===) (flatten arr) (flatten sarr))
-
-prop_growShrinkIO ::
-     ( Show (Array r ix Int)
-     , Eq (Array r ix Int)
-     , Mutable r ix Int
-     , Extract r ix Int
-     , Num ix
-     , Load (EltRepr r ix) ix Int
-     )
-  => r
-  -> Proxy ix
-  -> Arr r ix Int
-  -> NonNegative Int
-  -> Property
-prop_growShrinkIO _ _ (Arr arr) (NonNegative delta) =
-  monadicIO $
-  run $ do
-    marr <- thaw arr
-    let sz = size arr
-    k <- getDimM (unSz sz) (dimensions sz)
-    -- increase the outer most dimension, just so the structure doesn't change
-    newSz <- Sz <$> setDimM (unSz sz) (dimensions sz) (k + delta)
-    gMarr <- unsafeLinearGrow marr newSz
-    -- Make sure we can write into the newly allocated area
-    when (delta > 0) $ write' gMarr (unSz newSz - 1) delta
-    garr <- compute . extract' 0 sz <$> unsafeFreeze (getComp arr) gMarr
-    sarr <- freeze (getComp arr) =<< unsafeLinearShrink gMarr sz
-    pure (garr === arr .&&. sarr === arr)
-
-prop_atomicModifyIntArrayMany :: ArrIx P Ix2 Int -> Array B Ix1 Int -> Property
-prop_atomicModifyIntArrayMany (ArrIx arr ix) barr =
-  monadicIO $ do
-    xs <-
-      run $ do
-        marr <- thaw arr
-        mbarr' <- mapConcurrently (atomicModifyIntArray marr ix . const) barr
-        x <- A.read' marr ix
-        let xs = maybe (error "atomicModifyIntArray read'") toList (Prelude.sequenceA mbarr')
-        pure (x : xs)
-    return (L.sort (index' arr ix : toList barr) === L.sort xs)
-
-
-
-prop_atomicReadIntArrayMany :: Array P Ix2 Int -> Array B Ix1 Ix2 -> Property
-prop_atomicReadIntArrayMany arr bix =
-  monadicIO $
-  run $ do
-    marr <- thaw arr
-    as :: Array N Ix1 (Maybe Int) <- forM bix (A.read marr)
-    as' <- forM bix (atomicReadIntArray marr)
-    pure (as === as')
-
-
-prop_atomicWriteIntArrayMany :: Array P Ix2 Int -> Array B Ix1 Ix2 -> Fun Ix2 Int -> Property
-prop_atomicWriteIntArrayMany arr bix f =
-  monadicIO $
-  run $ do
-    marr <- thaw arr
-    marr' <- unsafeThaw arr
-    bs :: Array N Ix1 Bool <- forM bix (\ix -> write marr ix (apply f ix))
-    bs' <- forM bix (\ix -> atomicWriteIntArray marr' ix (apply f ix))
-    arrRes <- unsafeFreeze (getComp arr) marr
-    arrRes' <- unsafeFreeze (getComp arr) marr'
-    pure (bs === bs' .&&. arrRes === arrRes')
-
-
-
-prop_unfoldrList :: Sz1 -> Fun Word (Int, Word) -> Word -> Property
-prop_unfoldrList sz1 f i =
-  conjoin $
-  L.zipWith
-    (===)
-    (A.toList (runST $ unfoldrPrimM_ @P Seq sz1 (pure . apply f) i))
-    (L.unfoldr (Just . apply f) i)
-
-prop_unfoldrReverseUnfoldl :: Sz1 -> Fun Word (Int, Word) -> Word -> Property
-prop_unfoldrReverseUnfoldl sz1 f i =
-  runST (unfoldrPrimM_ @P Seq sz1 (pure . apply f) i) ===
-  rev (runST (unfoldlPrimM_ @P Seq sz1 (pure . swapTuple . apply f) i))
-    where swapTuple (x, y) = (y, x)
-          rev a = computeAs P $ backpermute' sz1 (\ix1 -> unSz sz1 - ix1 - 1) a
-
-
-
-mutableSpec ::
-     ( Show r
-     , Show (Array r Ix1 Int)
-     , Show (Array r Ix2 Int)
-     , Show (Array r Ix3 Int)
-     , Eq (Array r Ix1 Int)
-     , Eq (Array r Ix2 Int)
-     , Eq (Array r Ix3 Int)
-     , Mutable r Ix1 Int
-     , Mutable r Ix2 Int
-     , Mutable r Ix3 Int
-     , Construct r Ix1 Int
-     , Construct r Ix2 Int
-     , Construct r Ix3 Int
-     , Extract r Ix1 Int
-     , Extract r Ix2 Int
-     , Extract r Ix3 Int
-     , Resize r Ix1
-     , Resize r Ix2
-     , Resize r Ix3
-     , Load (EltRepr r Ix1) Ix1 Int
-     , Load (EltRepr r Ix2) Ix2 Int
-     , Load (EltRepr r Ix3) Ix3 Int
-     )
-  => r
-  -> SpecWith ()
-mutableSpec r =
-  describe (show r) $ do
-    describe "map == mapM" $ do
-      it "Ix1" $ property $ prop_MapMapM r (Proxy :: Proxy Ix1)
-      it "Ix2" $ property $ prop_MapMapM r (Proxy :: Proxy Ix2)
-      it "Ix3" $ property $ prop_MapMapM r (Proxy :: Proxy Ix3)
-    describe "imap == imapM" $ do
-      it "Ix1" $ property $ prop_iMapiMapM r (Proxy :: Proxy Ix1)
-      it "Ix2T" $ property $ prop_iMapiMapM r (Proxy :: Proxy Ix2)
-      it "Ix3T" $ property $ prop_iMapiMapM r (Proxy :: Proxy Ix3)
-    describe "makeArray == generateArrayST" $ do
-      it "Ix1" $ property $ prop_generateMakeST r (Proxy :: Proxy Ix1)
-      it "Ix2" $ property $ prop_generateMakeST r (Proxy :: Proxy Ix2)
-      it "Ix3" $ property $ prop_generateMakeST r (Proxy :: Proxy Ix3)
-    describe "makeArray == generateArrayIO" $ do
-      it "Ix1" $ property $ prop_generateMakeIO r (Proxy :: Proxy Ix1)
-      it "Ix2" $ property $ prop_generateMakeIO r (Proxy :: Proxy Ix2)
-      it "Ix3" $ property $ prop_generateMakeIO r (Proxy :: Proxy Ix3)
-    describe "shrink" $ do
-      it "Ix1" $ property $ prop_shrinkIO r (Proxy :: Proxy Ix1)
-      it "Ix2" $ property $ prop_shrinkIO r (Proxy :: Proxy Ix2)
-      it "Ix3" $ property $ prop_shrinkIO r (Proxy :: Proxy Ix3)
-    describe "grow+shrink" $ do
-      it "Ix1" $ property $ prop_growShrinkIO r (Proxy :: Proxy Ix1)
-      it "Ix2" $ property $ prop_growShrinkIO r (Proxy :: Proxy Ix2)
-      it "Ix3" $ property $ prop_growShrinkIO r (Proxy :: Proxy Ix3)
-
-generateSpec :: Spec
-generateSpec = do
-  mutableSpec P
-  mutableSpec S
-  mutableSpec U
-  mutableSpec B
-  mutableSpec N
-
-
-spec :: Spec
-spec = do
-  describe "GenerateM" generateSpec
-  describe "AtomicIntArraySpec" $ do
-    it "atomicReadIntArrayMany" $ property prop_atomicReadIntArrayMany
-    it "atomicWriteIntArrayMany" $ property prop_atomicWriteIntArrayMany
-    it "atomicModifyIntArrayMany" $ property prop_atomicModifyIntArrayMany
-  describe "Unfolding" $ do
-    it "unfoldrList" $ property prop_unfoldrList
-    it "unfoldrReverseUnfoldl" $ property prop_unfoldrReverseUnfoldl
diff --git a/tests/Data/Massiv/Array/Numeric/IntegralSpec.hs b/tests/Data/Massiv/Array/Numeric/IntegralSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Numeric/IntegralSpec.hs
+++ /dev/null
@@ -1,37 +0,0 @@
-module Data.Massiv.Array.Numeric.IntegralSpec
-  ( spec
-  ) where
-
-import Data.Massiv.Array as A
-import Data.Massiv.Array.Numeric.Integral
-import Test.Hspec
-
-gaussian :: Float -> Float
-gaussian x = exp (x ^ (2 :: Int))
-
-spec :: Spec
-spec = do
-  let (a, b) = (0, 2)
-      integrator rule = rule Seq N (\ scale -> gaussian . scale) a b (Sz1 1)
-  describe "Integral Approximation" $ do
-    it "Midpoint Rule" $ do
-      integrator midpointRule 4 ! 0 `shouldBe` 14.485613
-      integrator midpointRule 8 ! 0 `shouldBe` 15.905677
-      integrator midpointRule 16 ! 0 `shouldBe` 16.311854
-      integrator midpointRule 32 ! 0 `shouldBe` 16.417171
-      integrator midpointRule 64 ! 0 `shouldBe` 16.443748
-      integrator midpointRule 128 ! 0 `shouldBe` 16.450407
-    it "Trapezoid Rule" $ do
-      integrator trapezoidRule 4 ! 0 `shouldBe` 20.644558
-      integrator trapezoidRule 8 ! 0 `shouldBe` 17.565086
-      integrator trapezoidRule 16 ! 0 `shouldBe` 16.735381
-      integrator trapezoidRule 32 ! 0 `shouldBe` 16.523618
-      integrator trapezoidRule 64 ! 0 `shouldBe` 16.470394
-      integrator trapezoidRule 128 ! 0 `shouldBe` 16.457073
-    it "Simspon's Rule" $ do
-      integrator simpsonsRule 4 ! 0 `shouldBe` 17.353626
-      integrator simpsonsRule 8 ! 0 `shouldBe` 16.538595
-      integrator simpsonsRule 16 ! 0 `shouldBe` 16.458815
-      integrator simpsonsRule 32 ! 0 `shouldBe` 16.453030
-      integrator simpsonsRule 64 ! 0 `shouldBe` 16.452653
-      integrator simpsonsRule 128 ! 0 `shouldBe` 16.452629
diff --git a/tests/Data/Massiv/Array/Ops/ConstructSpec.hs b/tests/Data/Massiv/Array/Ops/ConstructSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/ConstructSpec.hs
+++ /dev/null
@@ -1,130 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE GADTs #-}
-{-# LANGUAGE LambdaCase #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-module Data.Massiv.Array.Ops.ConstructSpec (spec) where
-
-import Data.List as L
-import Data.Massiv.CoreArbitrary as A
-import Data.Proxy
-import qualified GHC.Exts as GHC (IsList(..))
-import Prelude as P
-
-prop_rangeEqRangeStep1 :: Int -> Int -> Property
-prop_rangeEqRangeStep1 from to = range Seq from to === rangeStep' Par from 1 to
-
-prop_rangeEqEnumFromN :: Int -> Int -> Property
-prop_rangeEqEnumFromN from to = range Seq from to === enumFromN Par from (Sz (to - from))
-
-prop_rangeStepEqEnumFromStepN :: Int -> NonZero Int -> Int -> Property
-prop_rangeStepEqEnumFromStepN from (NonZero step) sz =
-  rangeStep' Seq from step (from + step * sz) === enumFromStepN Par from step (Sz sz)
-
-
-prop_rangeStepExc :: Int -> Int -> Property
-prop_rangeStepExc from to =
-  assertException
-    (\case
-       IndexZeroException _ -> True
-       _ -> False)
-    (computeAs U (rangeStep' Seq from 0 to))
-
-prop_toFromListIsList ::
-     (Show (Array U ix Int), GHC.IsList (Array U ix Int), Index ix)
-  => Proxy ix
-  -> Arr U ix Int
-  -> Property
-prop_toFromListIsList _ (Arr arr) = arr === GHC.fromList (GHC.toList arr)
-
-
-prop_toFromList ::
-  forall ix . (Show (Array U ix Int), Nested LN ix Int, Ragged L ix Int)
-  => Proxy ix
-  -> Arr U ix Int
-  -> Property
-prop_toFromList _ (Arr arr) = comp === comp' .&&. arr === arr'
-  where comp = getComp arr
-        arr' = fromLists' comp (toLists arr)
-        comp' = getComp arr'
-
-
-prop_excFromToListIx2 :: Comp -> [[Int]] -> Property
-prop_excFromToListIx2 comp ls2 =
-  if P.null lsL || P.all (head lsL ==) lsL
-     then label "Expected Success" $ resultLs === ls2
-     else label "Expected Failure" $ assertSomeException resultLs
-  where
-    lsL = P.map P.length ls2
-    resultLs = toLists (fromLists' comp ls2 :: Array U Ix2 Int)
-
-
-prop_excFromToListIx3 :: Comp -> [[[Int]]] -> Property
-prop_excFromToListIx3 comp ls3
-  | P.null (P.concat (P.concat ls3)) =
-    classify True "Expected Success" $ counterexample (show arr) $ totalElem (size arr) === 0
-  | P.all (head lsL ==) lsL && P.all (P.all (head (head lsLL) ==)) lsLL =
-    classify True "Expected Success" $ counterexample (show arr) $ resultLs === ls3
-  | otherwise = classify True "Expected Failure" $ assertSomeException resultLs
-  where
-    arr = fromLists' comp ls3 :: Array U Ix3 Int
-    resultLs = toLists arr
-    lsL = P.map P.length ls3
-    lsLL = P.map (P.map P.length) ls3
-
-
-specIx1 :: Spec
-specIx1 = do
-  it "toFromList" $ property (prop_toFromList (Proxy :: Proxy Ix1))
-  it "toFromListIsList" $ property (prop_toFromListIsList (Proxy :: Proxy Ix1))
-  it "rangeEqRangeStep1" $ property prop_rangeEqRangeStep1
-  it "rangeEqEnumFromN" $ property prop_rangeEqEnumFromN
-  it "rangeStepEqEnumFromStepN" $ property prop_rangeStepEqEnumFromStepN
-  it "rangeStepExc" $ property prop_rangeStepExc
-
-specIx2 :: Spec
-specIx2 = do
-  it "toFromList" $ property (prop_toFromList (Proxy :: Proxy Ix2))
-  it "toFromListIsList" $ property (prop_toFromListIsList (Proxy :: Proxy Ix2))
-  it "excFromToListIx2" $ property prop_excFromToListIx2
-
-specIx3 :: Spec
-specIx3 = do
-  it "toFromList" $ property (prop_toFromList (Proxy :: Proxy Ix3))
-  it "toFromListIsList" $ property (prop_toFromListIsList (Proxy :: Proxy Ix3))
-  it "excFromToListIx3" $ property prop_excFromToListIx3
-
-mkIntermediate :: Int -> Array U Ix1 Int
-mkIntermediate t = A.fromList Seq [t + 50, t + 75]
-
-initArr :: Array N Ix1 (Array U Ix1 Int)
-initArr = makeArray Seq (Sz1 3) mkIntermediate
-
-initArr2 :: Array N Ix2 (Array U Ix1 Int)
-initArr2 = makeArray Seq (Sz 2) (\ (x :. y) -> mkIntermediate (x+y))
-
-prop_unfoldrList :: Sz1 -> Fun Word (Int, Word) -> Word -> Property
-prop_unfoldrList sz1 f i =
-  conjoin $
-  L.zipWith
-    (===)
-    (A.toList (computeAs P $ unfoldrS_ Seq sz1 (apply f) i))
-    (L.unfoldr (Just . apply f) i)
-
-specExpand :: Spec
-specExpand = do
-  it "expandOuter" $ compute (expandOuter 2 A.index' initArr :: Array D Ix2 Int) `shouldBe`
-    resize' (Sz2 2 3) (fromList Seq [50, 51, 52, 75, 76, 77] :: Array U Ix1 Int)
-  it "expandInner" $ compute (expandInner 2 A.index' initArr :: Array D Ix2 Int) `shouldBe`
-    resize' (Sz2 3 2) (fromList Seq [50, 75, 51, 76, 52, 77] :: Array U Ix1 Int)
-  it "expandwithin" $ compute (expandWithin Dim1 2 A.index' initArr2 :: Array D Ix3 Int) `shouldBe`
-    resize' (Sz 2) (fromList Seq [50, 75, 51, 76, 51, 76, 52, 77] :: Array U Ix1 Int)
-  it "expandwithin'" $ compute (expandWithin' 1 2 A.index' initArr2 :: Array D Ix3 Int) `shouldBe`
-    resize' (Sz 2) (fromList Seq [50, 75, 51, 76, 51, 76, 52, 77] :: Array U Ix1 Int)
-
-spec :: Spec
-spec = do
-  describe "Ix1" specIx1
-  describe "Ix2" specIx2
-  describe "Ix3" specIx3
-  describe "Expand" specExpand
-  describe "Unfolding" $ it "unfoldrS_" $ property prop_unfoldrList
diff --git a/tests/Data/Massiv/Array/Ops/FoldSpec.hs b/tests/Data/Massiv/Array/Ops/FoldSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/FoldSpec.hs
+++ /dev/null
@@ -1,57 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-module Data.Massiv.Array.Ops.FoldSpec (spec) where
-
-import qualified Data.Foldable as F
-import Data.Massiv.CoreArbitrary as A
-import Data.Semigroup
-import Prelude hiding (map, product, sum)
-
-
-
-prop_SumSEqSumP :: Index ix => proxy ix -> Array D ix Int -> Bool
-prop_SumSEqSumP _ arr = sum arr == sum (setComp Par arr)
-
-
-prop_ProdSEqProdP :: Index ix => proxy ix -> Array D ix Int -> Bool
-prop_ProdSEqProdP _ arr = product arr == product (setComp Par arr)
-
-prop_NestedFoldP :: Array D Ix1 (Array D Ix1 Int) -> Bool
-prop_NestedFoldP arr = sum (setComp Par (map sum $ setComp Par arr)) == sum (map sum arr)
-
-
-specFold ::
-     (Arbitrary ix, CoArbitrary ix, Index ix, Show (Array D ix Int))
-  => proxy ix
-  -> String
-  -> Spec
-specFold proxy dimStr =
-  describe dimStr $ do
-    it "sumS Eq sumP" $ property $ prop_SumSEqSumP proxy
-    it "prodS Eq prodP" $ property $ prop_ProdSEqProdP proxy
-
-foldOpsProp :: (Source P ix Int) => proxy ix -> Fun Int Bool -> ArrTiny1 P ix Int -> Property
-foldOpsProp _ f (ArrTiny1 arr) =
-  (A.maximum' arr === getMax (foldMono Max arr)) .&&.
-  (A.minimum' arr === getMin (foldSemi Min maxBound arr)) .&&.
-  (A.sum arr === F.sum ls) .&&.
-  (A.product (A.map ((+ 0.1) . (fromIntegral :: Int -> Double)) arr) ===
-   getProduct (foldMono (Product . (+ 0.1) . fromIntegral) arr)) .&&.
-  (A.all (apply f) arr === F.all (apply f) ls) .&&.
-  (A.any (apply f) arr === F.any (apply f) ls) .&&.
-  (A.or (A.map (apply f) arr) === F.or (fmap (apply f) ls)) .&&.
-  (A.and (A.map (apply f) arr) === F.and (fmap (apply f) ls))
-  where
-    ls = toList arr
-
-spec :: Spec
-spec = do
-  specFold (Nothing :: Maybe Ix1) "Ix1"
-  specFold (Nothing :: Maybe Ix2) "Ix2"
-  it "Nested Parallel Fold" $ property prop_NestedFoldP
-  describe "Foldable Props" $ do
-    it "Ix1" $ property $ foldOpsProp (Nothing :: Maybe Ix1)
-    it "Ix2" $ property $ foldOpsProp (Nothing :: Maybe Ix2)
-    it "Ix3" $ property $ foldOpsProp (Nothing :: Maybe Ix3)
diff --git a/tests/Data/Massiv/Array/Ops/MapSpec.hs b/tests/Data/Massiv/Array/Ops/MapSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/MapSpec.hs
+++ /dev/null
@@ -1,125 +0,0 @@
-{-# LANGUAGE AllowAmbiguousTypes #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE GADTs #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.Ops.MapSpec (spec) where
-
-import Data.IORef
-import Control.Monad.ST
-import Data.Foldable as F
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary as A
-import Prelude as P
-import Control.Scheduler.Internal
-
-prop_zipUnzip ::
-     (Index ix, Show (Array D ix Int))
-  => Array D ix Int
-  -> Array D ix Int
-  -> Property
-prop_zipUnzip arr1 arr2 =
-  (extract' zeroIndex sz arr1, extract' zeroIndex sz arr2) === A.unzip (A.zip arr1 arr2)
-  where sz = Sz (liftIndex2 min (unSz (size arr1)) (unSz (size arr2)))
-
-prop_zipFlip ::
-     (Index ix, Show (Array D ix (Int, Int)))
-  => Array D ix Int
-  -> Array D ix Int
-  -> Property
-prop_zipFlip arr1 arr2 =
-  A.zip arr1 arr2 ===
-  A.map (\(e2, e1) -> (e1, e2)) (A.zip arr2 arr1)
-
-prop_zipUnzip3 ::
-     (Index ix, Show (Array D ix Int))
-  => Array D ix Int
-  -> Array D ix Int
-  -> Array D ix Int
-  -> Property
-prop_zipUnzip3 arr1 arr2 arr3 =
-  (extract' zeroIndex sz arr1, extract' zeroIndex sz arr2, extract' zeroIndex sz arr3) ===
-  A.unzip3 (A.zip3 arr1 arr2 arr3)
-  where
-    sz =
-      Sz (liftIndex2 min (liftIndex2 min (unSz (size arr1)) (unSz (size arr2))) (unSz (size arr3)))
-
-prop_zipFlip3 ::
-     (Index ix, Show (Array D ix (Int, Int, Int)))
-  => Array D ix Int
-  -> Array D ix Int
-  -> Array D ix Int
-  -> Property
-prop_zipFlip3 arr1 arr2 arr3 =
-  A.zip3 arr1 arr2 arr3 === A.map (\(e3, e2, e1) -> (e1, e2, e3)) (A.zip3 arr3 arr2 arr1)
-
-
-
-prop_itraverseA ::
-     (Index ix, Show (Array U ix Int)) => Array D ix Int -> Fun (ix, Int) Int -> Property
-prop_itraverseA arr fun =
-  alt_imapM (\ix -> Just . applyFun2Compat fun ix) arr ===
-  itraverseAR U (\ix -> Just . applyFun2Compat fun ix) arr
-
-
-mapSpec ::
-     forall ix.
-     ( Arbitrary ix
-     , CoArbitrary ix
-     , Index ix
-     , Function ix
-     , Show (Array U ix Int)
-     , Show (Array D ix Int)
-     , Show (Array D ix (Int, Int))
-     , Show (Array D ix (Int, Int, Int))
-     )
-  => Spec
-mapSpec = do
-  describe "Zipping" $ do
-    it "zipUnzip" $ property $ prop_zipUnzip @ix
-    it "zipFlip" $ property $ prop_zipFlip @ix
-    it "zipUnzip3" $ property $ prop_zipUnzip3 @ix
-    it "zipFlip3" $ property $ prop_zipFlip3 @ix
-  describe "Traversing" $ do
-    it "itraverseA" $ property $ prop_itraverseA @ix
-  describe "StatefulMapping" $ do
-    it "mapWS" $ property $ prop_MapWS @ix
-
-spec :: Spec
-spec = do
-  describe "Ix1" $ mapSpec @Ix1
-  describe "Ix2" $ mapSpec @Ix2
-  describe "Ix3" $ mapSpec @Ix3
-  describe "Ix4" $ mapSpec @Ix4
-
-
-
-alt_imapM
-  :: (Applicative f, Mutable r2 t1 b, Source r1 t1 t2) =>
-     (t1 -> t2 -> f b) -> Array r1 t1 t2 -> f (Array r2 t1 b)
-alt_imapM f arr = fmap loadList $ P.traverse (uncurry f) $ foldrS (:) [] (zipWithIndex arr)
-  where
-    loadList xs =
-      runST $ do
-        marr <- unsafeNew (size arr)
-        _ <- F.foldlM (\i e -> unsafeLinearWrite marr i e >> return (i + 1)) 0 xs
-        unsafeFreeze (getComp arr) marr
-    {-# INLINE loadList #-}
-
-zipWithIndex :: forall r ix e . Source r ix e => Array r ix e -> Array D ix (ix, e)
-zipWithIndex arr = A.zip (range Seq zeroIndex (unSz (size arr))) arr
-{-# INLINE zipWithIndex #-}
-
-
-prop_MapWS :: (Show (Array U ix Int), Index ix) => Array U ix Int -> Property
-prop_MapWS arr =
-  monadicIO $
-  run $ do
-    states <- initWorkerStates (getComp arr) (\_ -> newIORef 0)
-    arr' <-
-      forWS states arr $ \e ref -> do
-        acc <- readIORef ref
-        writeIORef ref (acc + e)
-        pure e
-    accsArr <- A.mapM @P readIORef (evalArray Seq (_workerStatesArray states))
-    pure (A.sum arr' === A.sum accsArr .&&. arr === arr')
diff --git a/tests/Data/Massiv/Array/Ops/SliceSpec.hs b/tests/Data/Massiv/Array/Ops/SliceSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/SliceSpec.hs
+++ /dev/null
@@ -1,240 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE TypeFamilies #-}
-module Data.Massiv.Array.Ops.SliceSpec (spec) where
-
-import Control.Applicative ((<|>))
-import Control.Exception
-import Data.Massiv.Array.Unsafe
-import Data.Massiv.CoreArbitrary
------------
--- Size --
------------
-
-prop_ExtractEqualsExtractFromTo ::
-     ( Eq (Array (EltRepr r ix) ix e)
-     , Show (Array (EltRepr r ix) ix e)
-     , Extract r ix e
-     )
-  => proxy (r, ix, e)
-  -> SzIx ix
-  -> Array r ix e
-  -> Property
-prop_ExtractEqualsExtractFromTo _ (SzIx (Sz eIx) sIx) arr =
-  (extractFromToM sIx eIx arr <|> Nothing) === extractM sIx (Sz (liftIndex2 (-) eIx sIx)) arr
-
-
-specSizeN ::
-     ( Eq (Array (EltRepr r ix) ix e)
-     , Show (Array (EltRepr r ix) ix e)
-     , Arbitrary (Array r ix e)
-     , Show (Array r ix e)
-     , Arbitrary ix
-     , Extract r ix e
-     )
-  => proxy (r, ix, e)
-  -> Spec
-specSizeN proxy =
-  describe "extract" $
-    it "ExtractEqualsExtractFromTo" $ property $ prop_ExtractEqualsExtractFromTo proxy
-
-
------------
--- Slice --
------------
-
-
-prop_SliceRight ::
-     (Slice r ix e, OuterSlice r ix e, Eq (Elt r ix e), Show (Elt r ix e))
-  => proxy (r, ix, e)
-  -> Int
-  -> Array r ix e
-  -> Property
-prop_SliceRight _ i arr =
-  either (Left . displayException) Right (arr !?> i) ===
-  either (Left . displayException) Right (arr <!?> (dimensions (size arr), i))
-
-
-prop_SliceLeft ::
-     (Slice r ix e, InnerSlice r ix e, Eq (Elt r ix e), Show (Elt r ix e))
-  => proxy (r, ix, e)
-  -> Int
-  -> Array r ix e
-  -> Property
-prop_SliceLeft _ i arr =
-  either (Left . displayException) Right (arr <!? i) ===
-  either (Left . displayException) Right (arr <!?> (1, i))
-
-
-prop_SliceIndexDim2D :: ArrIx D Ix2 Int -> Property
-prop_SliceIndexDim2D (ArrIx arr ix@(i :. j)) =
-  val === evaluate' (arr <! j) i .&&.
-  val === evaluate' (arr !> i) j
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim2RankD :: ArrIx D Ix2 Int -> Property
-prop_SliceIndexDim2RankD (ArrIx arr ix@(i :. j)) =
-  val === evaluate' (arr <!> (2, i)) j .&&.
-  val === evaluate' (arr <!> (1, j)) i
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim3D :: ArrIx D Ix3 Int -> Property
-prop_SliceIndexDim3D (ArrIx arr ix@(i :> j :. k)) =
-  val === evaluate' (arr <! k <! j) i .&&.
-  val === evaluate' (arr !> i !> j) k .&&.
-  val === evaluate' (arr <! k !> i) j .&&.
-  val === evaluate' (arr !> i <! k) j
-  where
-    val = unsafeIndex arr ix
-
-prop_SliceIndexDim3RankD :: ArrIx D Ix3 Int -> Property
-prop_SliceIndexDim3RankD (ArrIx arr ix@(i :> j :. k)) =
-  val === evaluate' (arr <!> (3, i) <!> (2, j)) k .&&.
-  val === evaluate' (arr <!> (3, i) <!> (1, k)) j .&&.
-  val === evaluate' (arr <!> (2, j) <!> (2, i)) k .&&.
-  val === evaluate' (arr <!> (2, j) <!> (1, k)) i .&&.
-  val === evaluate' (arr <!> (1, k) <!> (2, i)) j .&&.
-  val === evaluate' (arr <!> (1, k) <!> (1, j)) i
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim2M :: ArrIx M Ix2 Int -> Property
-prop_SliceIndexDim2M (ArrIx arr ix@(i :. j)) =
-  val === (arr !> i ! j) .&&.
-  val === (arr <! j ! i)
-  where
-    val = unsafeIndex arr ix
-
-prop_SliceIndexDim2RankM :: ArrIx M Ix2 Int -> Property
-prop_SliceIndexDim2RankM (ArrIx arr ix@(i :. j)) =
-  val === (arr <!> (2, i) ! j) .&&.
-  val === (arr <!> (1, j) ! i)
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim3M :: ArrIx M Ix3 Int -> Property
-prop_SliceIndexDim3M (ArrIx arr ix@(i :> j :. k)) =
-  val === (arr <! k <! j ! i) .&&.
-  val === (arr !> i !> j ! k) .&&.
-  val === (arr <! k !> i ! j) .&&.
-  val === (arr !> i <! k ! j)
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim3RankM :: ArrIx M Ix3 Int -> Property
-prop_SliceIndexDim3RankM (ArrIx arr ix@(i :> j :. k)) =
-  val === (arr <!> (3, i) <!> (2, j) ! k) .&&.
-  val === (arr <!> (3, i) <!> (1, k) ! j) .&&.
-  val === (arr <!> (2, j) <!> (2, i) ! k) .&&.
-  val === (arr <!> (2, j) <!> (1, k) ! i) .&&.
-  val === (arr <!> (1, k) <!> (2, i) ! j) .&&.
-  val === (arr <!> (1, k) <!> (1, j) ! i)
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim4D :: ArrIx D Ix4 Int -> Property
-prop_SliceIndexDim4D (ArrIx arr ix@(i1 :> i2 :> i3 :. i4)) =
-  val === evaluate' (arr !> i1 !> i2 !> i3) i4 .&&.
-  val === evaluate' (arr !> i1 !> i2 <! i4) i3 .&&.
-  val === evaluate' (arr !> i1 <! i4 <! i3) i2 .&&.
-  val === evaluate' (arr !> i1 <! i4 !> i2) i3 .&&.
-  val === evaluate' (arr <! i4 !> i1 !> i2) i3 .&&.
-  val === evaluate' (arr <! i4 !> i1 <! i3) i2 .&&.
-  val === evaluate' (arr <! i4 <! i3 <! i2) i1 .&&.
-  val === evaluate' (arr <! i4 <! i3 !> i1) i2
-  where
-    val = unsafeIndex arr ix
-
-prop_SliceIndexDim4RankD :: ArrIx D Ix4 Int -> Property
-prop_SliceIndexDim4RankD (ArrIx arr ix@(i1 :> i2 :> i3 :. i4)) =
-  val === unsafeIndex (arr <!> (4, i1) <!> (3, i2) <!> (2, i3)) i4 .&&.
-  val === unsafeIndex (arr <!> (4, i1) <!> (2, i3) <! i4) i2 .&&.
-  val === unsafeIndex (arr <!> (3, i2) <!> (3, i1)) (i3 :. i4) .&&.
-  val === unsafeIndex (arr <!> (2, i3) <!> (2, i2)) (i1 :. i4) .&&.
-  val === unsafeIndex (arr <!> (2, i3) <!> (1, i4) !> i1) i2 .&&.
-  val === unsafeIndex (arr <!> (1, i4) !> i1 !> i2) i3
-  where
-    val = evaluate' arr ix
-
-
-prop_SliceIndexDim4RankM :: ArrIx M Ix4 Int -> Property
-prop_SliceIndexDim4RankM (ArrIx arr ix@(i1 :> i2 :> i3 :. i4)) =
-  val === (arr <!> (4, i1) <!> (3, i2) <!> (2, i3) ! i4) .&&.
-  val === (arr <!> (4, i1) <!> (2, i3) <! i4 ! i2) .&&.
-  val === (arr <!> (3, i2) <!> (3, i1) ! (i3 :. i4)) .&&.
-  val === (arr <!> (2, i3) <!> (2, i2) ! (i1 :. i4)) .&&.
-  val === (arr <!> (2, i3) <!> (1, i4) !> i1 ! i2) .&&.
-  val === (arr <!> (1, i4) !> i1 !> i2 ! i3)
-  where
-    val = unsafeIndex arr ix
-
-
-prop_SliceIndexDim4M :: ArrIx M Ix4 Int -> Property
-prop_SliceIndexDim4M (ArrIx arr ix@(i1 :> i2 :> i3 :. i4)) =
-  val === (arr !> i1 !> i2 !> i3 ! i4) .&&.
-  val === (arr !> i1 !> i2 <! i4 ! i3) .&&.
-  val === (arr !> i1 <! i4 <! i3 ! i2) .&&.
-  val === (arr !> i1 <! i4 !> i2 ! i3) .&&.
-  val === (arr <! i4 !> i1 !> i2 ! i3) .&&.
-  val === (arr <! i4 !> i1 <! i3 ! i2) .&&.
-  val === (arr <! i4 <! i3 <! i2 ! i1) .&&.
-  val === (arr <! i4 <! i3 !> i1 ! i2)
-  where
-    val = unsafeIndex arr ix
-
-
-
-specSliceN :: ( Arbitrary (Array r ix e)
-              , Show (Array r ix e)
-              , Slice r ix e
-              , OuterSlice r ix e
-              , InnerSlice r ix e
-              , Eq (Elt r ix e)
-              , Show (Elt r ix e)
-              )
-           => proxy (r, ix, e) -> Spec
-specSliceN proxy =
-  describe "Slice" $ do
-    it "SliceRight" $ property $ prop_SliceRight proxy
-    it "SliceLeft" $ property $ prop_SliceLeft proxy
-
-
-
-spec :: Spec
-spec = do
-  describe "Ix1" $
-    specSizeN (Nothing :: Maybe (D, Ix1, Int))
-  describe "Ix2" $ do
-    specSizeN (Nothing :: Maybe (D, Ix2, Int))
-    specSliceN (Nothing :: Maybe (D, Ix2, Int))
-    describe "SliceIndex" $ do
-      it "Delayed" $ property prop_SliceIndexDim2D
-      it "Rank - Delayed" $ property prop_SliceIndexDim2RankD
-      it "Manifest" $ property prop_SliceIndexDim2M
-      it "Rank - Manifest" $ property prop_SliceIndexDim2RankM
-  describe "Ix3" $ do
-    specSizeN (Nothing :: Maybe (D, Ix3, Int))
-    specSliceN (Nothing :: Maybe (D, Ix3, Int))
-    describe "SliceIndex" $ do
-      it "Delayed" $ property prop_SliceIndexDim3D
-      it "Rank - Delayed" $ property prop_SliceIndexDim3RankD
-      it "Manifest" $ property prop_SliceIndexDim3M
-      it "Rank - Manifest" $ property prop_SliceIndexDim3RankM
-  describe "Ix4" $ do
-    specSizeN (Nothing :: Maybe (D, Ix4, Int))
-    specSliceN (Nothing :: Maybe (D, Ix4, Int))
-    describe "SliceIndex" $ do
-      it "Delayed" $ property prop_SliceIndexDim4D
-      it "Rank - Delayed" $ property prop_SliceIndexDim4RankD
-      it "Manifest" $ property prop_SliceIndexDim4M
-      it "Rank - Manifest" $ property prop_SliceIndexDim4RankM
diff --git a/tests/Data/Massiv/Array/Ops/SortSpec.hs b/tests/Data/Massiv/Array/Ops/SortSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/SortSpec.hs
+++ /dev/null
@@ -1,17 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.Ops.SortSpec (spec) where
-
-import Data.List as L
-import Data.Massiv.CoreArbitrary as A
-
-
-prop_IsSorted :: (b -> b) -> ([Int] -> b) -> (b -> [Int]) -> [Int] -> Property
-prop_IsSorted sortWith from to xs =
-  to (sortWith (from xs)) === sort xs
-
-spec :: Spec
-spec =
-  describe "QuickSort" $ do
-    it "Seq" $ property $ prop_IsSorted (quicksort @P) (A.fromList Seq) A.toList
-    it "Par" $ property $ prop_IsSorted (quicksort @P) (A.fromList (ParN 4)) A.toList
diff --git a/tests/Data/Massiv/Array/Ops/TransformSpec.hs b/tests/Data/Massiv/Array/Ops/TransformSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/Ops/TransformSpec.hs
+++ /dev/null
@@ -1,157 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE LambdaCase #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE TypeApplications #-}
-module Data.Massiv.Array.Ops.TransformSpec (spec) where
-
-import Data.Massiv.CoreArbitrary as A
-import Data.Sequence as S
-import Prelude as P
-import Data.Foldable as F (foldl', toList)
-import Data.Maybe
-
-prop_transposeOuterInner :: Arr D Ix2 Int -> Property
-prop_transposeOuterInner (Arr arr) = transposeOuter arr === transpose arr
-
-prop_upsampleDownsample ::
-     (Show (Array P ix Int), Index ix) => ArrTiny P ix Int -> Stride ix -> Int -> Property
-prop_upsampleDownsample (ArrTiny arr) stride fill =
-  arr === compute (downsample stride (computeAs P (upsample fill stride arr)))
-
-prop_ExtractAppend
-  :: (Show (Array P ix Int), Index ix)
-  => DimIx ix -> ArrIx P ix Int -> Property
-prop_ExtractAppend (DimIx dim) (ArrIx arr ix) =
-  arr === compute (uncurry (append' dim) $ A.splitAt' dim (getDim' ix dim) arr)
-
-prop_SplitExtract
-  :: (Show (Array P ix Int), Show (Array M ix Int), Index ix)
-  => DimIx ix -> ArrIx P ix Int -> Positive Int -> Property
-prop_SplitExtract (DimIx dim) (ArrIx arr ix) (Positive n) =
-  (computeAs P <$> splitAt' dim i arr) === (left, computeAs P (append' dim center right)) .&&.
-  (computeAs P splitLeft, splitRight) === (computeAs P (append' dim left center), right)
-  where i = getDim' ix dim
-        k = getDim' (unSz (size arr)) dim
-        n' = n `mod` (k - i)
-        (left, center, right) = either throw id (splitExtractM dim i (Sz n') arr)
-        (splitLeft, splitRight) = splitAt' dim (i + n') arr
-
-prop_ConcatAppend
-  :: (Show (Array P ix Int), Index ix)
-  => DimIx ix -> Comp -> Sz ix -> NonEmptyList (Fun ix Int) -> Property
-prop_ConcatAppend (DimIx dim) comp sz (NonEmpty fns) =
-  foldl1 (\arr -> computeAs P . append' dim arr) arrs ===
-  computeAs P (concat' dim arrs)
-  where
-    arrs = P.map (makeArrayR P comp sz . apply) fns
-
-prop_AppendMappend
-  :: Array D Ix1 Int -> Array D Ix1 Int -> Property
-prop_AppendMappend arr1 arr2 =
-  computeAs P (append' 1 arr1 arr2) === computeAs P (toLoadArray arr1 <> toLoadArray arr2)
-
-prop_ConcatMconcat
-  :: [Array D Ix1 Int] -> Property
-prop_ConcatMconcat arrs =
-  computeAs P (concat' 1 (A.empty : arrs)) === computeAs P (mconcat (fmap toLoadArray arrs))
-
-prop_ExtractSizeMismatch ::
-     Index ix => ArrTiny P ix Int -> Positive Int -> Property
-prop_ExtractSizeMismatch (ArrTiny arr) (Positive n) =
-  assertExceptionIO (SizeElementsMismatchException sz sz' ==) $ resizeM sz' arr
-  where
-    sz = size arr
-    sz' = Sz (totalElem sz + n)
-
-spec :: Spec
-spec = do
-  it "transposeOuterInner" $ property prop_transposeOuterInner
-  describe "upsampleDownsample" $ do
-    it "Ix1" $ property (prop_upsampleDownsample @Ix1)
-    it "Ix2" $ property (prop_upsampleDownsample @Ix2)
-    it "Ix3" $ property (prop_upsampleDownsample @Ix3)
-    it "Ix4" $ property (prop_upsampleDownsample @Ix4)
-  describe "extractSizeMismatch" $ do
-    it "Ix1" $ property (prop_ExtractSizeMismatch @Ix1)
-    it "Ix2" $ property (prop_ExtractSizeMismatch @Ix2)
-    it "Ix3" $ property (prop_ExtractSizeMismatch @Ix3)
-    it "Ix4" $ property (prop_ExtractSizeMismatch @Ix4)
-  describe "ExtractAppend" $ do
-    it "Ix1" $ property (prop_ExtractAppend @Ix1)
-    it "Ix2" $ property (prop_ExtractAppend @Ix2)
-    it "Ix3" $ property (prop_ExtractAppend @Ix3)
-    it "Ix4" $ property (prop_ExtractAppend @Ix4)
-  describe "ExtractAppend" $ do
-    it "Ix1" $ property (prop_SplitExtract @Ix1)
-    it "Ix2" $ property (prop_SplitExtract @Ix2)
-    it "Ix3" $ property (prop_SplitExtract @Ix3)
-    it "Ix4" $ property (prop_SplitExtract @Ix4)
-  describe "ConcatAppend" $ do
-    it "Ix1" $ property (prop_ConcatAppend @Ix1)
-    it "Ix2" $ property (prop_ConcatAppend @Ix2)
-    it "Ix3" $ property (prop_ConcatAppend @Ix3)
-    it "Ix4" $ property (prop_ConcatAppend @Ix4)
-  describe "Monoid" $ do
-    it "Ix1" $ property prop_AppendMappend
-    it "Ix1" $ property prop_ConcatMconcat
-  describe "Sequence" $ do
-    it "ConsSnoc" $ property prop_ConsSnoc
-    it "UnconsUnsnoc" $ property prop_UnconsUnsnoc
-  describe "zoomWithGrid" $ do
-    it "Ix1" $ property (prop_zoomWithGridStrideCompute @Ix1)
-    it "Ix2" $ property (prop_zoomWithGridStrideCompute @Ix2)
-    it "Ix3" $ property (prop_zoomWithGridStrideCompute @Ix3)
-    it "Ix4" $ property (prop_zoomWithGridStrideCompute @Ix4)
-
-prop_zoomWithGridStrideCompute :: (Show (Array P ix Int), Index ix) => Array D ix Int -> Stride ix -> Int -> Property
-prop_zoomWithGridStrideCompute arr stride defVal =
-  (computeWithStrideAs P stride' arr' ===
-   A.replicate Seq (Sz (liftIndex (+ 1) $ unSz (size arr))) defVal) .&&.
-  (computeWithStrideAs P stride' (extract' (pureIndex 1) sz' arr') === compute arr)
-  where
-    arr' = computeAs P (zoomWithGrid defVal stride arr)
-    sz' = Sz (liftIndex (subtract 1) $ unSz (size arr'))
-    stride' = Stride (liftIndex (+ 1) $ unStride stride)
-
-
-prop_UnconsUnsnoc :: Array D Ix1 Int -> Bool -> Property
-prop_UnconsUnsnoc arr unconsFirst =
-  preJust $ do
-    (arr', u, s) <-
-      if unconsFirst
-        then do
-          (u, au) <- unconsM arr
-          (as, s) <- unsnocM au
-          pure (as, u, s)
-        else do
-          (as, s) <- unsnocM arr
-          (u, au) <- unconsM as
-          pure (au, u, s)
-    pure (computeAs U (A.snoc (A.cons u (toLoadArray (computeAs U arr'))) s) === compute arr)
-
-preJust :: Testable prop => Maybe prop -> Property
-preJust m = isJust m ==> fromJust m
-
-prop_ConsSnoc :: Array D Ix1 Int -> [SeqOp Int] -> Property
-prop_ConsSnoc arr ops =
-  A.toList (computeAs U (foldl' applyArraySeqOp (toLoadArray arr) ops)) ===
-  F.toList (foldl' applySequenceSeqOp (S.fromList (A.toList arr)) ops)
-
-data SeqOp e = Cons e | Snoc e deriving (Eq, Show)
-
-instance Arbitrary e => Arbitrary (SeqOp e) where
-  arbitrary = do
-    e <- arbitrary
-    elements [Cons e, Snoc e]
-
-applyArraySeqOp :: Array DL Ix1 e -> SeqOp e -> Array DL Ix1 e
-applyArraySeqOp arr = \case
-  Cons x -> A.cons x arr
-  Snoc x -> A.snoc arr x
-
-
-applySequenceSeqOp :: Seq a -> SeqOp a -> Seq a
-applySequenceSeqOp arr = \case
-  Cons x -> x <| arr
-  Snoc x -> arr |> x
diff --git a/tests/Data/Massiv/Array/StencilSpec.hs b/tests/Data/Massiv/Array/StencilSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Array/StencilSpec.hs
+++ /dev/null
@@ -1,196 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE OverloadedLists #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-module Data.Massiv.Array.StencilSpec (spec) where
-
-import Control.DeepSeq (deepseq)
-import Data.Default (Default(def))
-import Data.Massiv.CoreArbitrary as A
-import Data.Proxy
-
--- sum3x3Stencil :: Fractional a => Stencil Ix2 a a
--- sum3x3Stencil = makeConvolutionStencil (3 :. 3) (1 :. 1) $ \ get ->
---   get (-1 :. -1) 1 . get (-1 :. 0) 1 . get (-1 :. 1) 1 .
---   get ( 0 :. -1) 1 . get ( 0 :. 0) 1 . get ( 0 :. 1) 1 .
---   get ( 1 :. -1) 1 . get ( 1 :. 0) 1 . get ( 1 :. 1) 1
--- {-# INLINE sum3x3Stencil #-}
-
-
-singletonStencil :: (Index ix) => (Int -> Int) -> Stencil ix Int Int
-singletonStencil f =
-  makeStencil oneSz zeroIndex $ \ get -> fmap f (get zeroIndex)
-{-# INLINE singletonStencil #-}
-
-
-prop_MapSingletonStencil :: (Load DW ix Int, Manifest U ix Int) =>
-                            Proxy ix -> Fun Int Int -> Border Int -> ArrP U ix Int -> Bool
-prop_MapSingletonStencil _ f b (ArrP arr) =
-  computeAs U (mapStencil b (singletonStencil (apply f)) arr) == computeAs U (A.map (apply f) arr)
-
-prop_MapSingletonStencilWithStride :: (StrideLoad DW ix Int, Manifest U ix Int) =>
-                                      Proxy ix -> Fun Int Int -> Border Int -> ArrP U ix Int -> Bool
-prop_MapSingletonStencilWithStride _ f b (ArrP arr) =
-  computeWithStride oneStride (mapStencil b (singletonStencil (apply f)) arr) ==
-  computeAs U (A.map (apply f) arr)
-
--- Tests out of bounds stencil indexing
-prop_DangerousStencil ::
-     Index ix => Proxy ix -> NonZero Int -> DimIx ix -> SzIx ix -> Property
-prop_DangerousStencil _ (NonZero s) (DimIx r) (SzIx sz ix) =
-  ix' `deepseq` assertSomeException $ makeStencil sz ix $ \get -> get ix' :: Value Int
-  where
-    ix' = liftIndex (* signum s) (setDim' zeroIndex r (getDim' (unSz sz) r))
-
-
-stencilSpec :: Spec
-stencilSpec = do
-  describe "MapSingletonStencil" $ do
-    it "Ix1" $ property $ prop_MapSingletonStencil (Proxy :: Proxy Ix1)
-    it "Ix2" $ property $ prop_MapSingletonStencil (Proxy :: Proxy Ix2)
-    it "Ix3" $ property $ prop_MapSingletonStencil (Proxy :: Proxy Ix3)
-    it "Ix4" $ property $ prop_MapSingletonStencil (Proxy :: Proxy Ix4)
-    it "Ix2T" $ property $ prop_MapSingletonStencil (Proxy :: Proxy Ix2T)
-  describe "MapSingletonStencilWithStride" $ do
-    it "Ix1" $ property $ prop_MapSingletonStencilWithStride (Proxy :: Proxy Ix1)
-    it "Ix2" $ property $ prop_MapSingletonStencilWithStride (Proxy :: Proxy Ix2)
-    it "Ix3" $ property $ prop_MapSingletonStencilWithStride (Proxy :: Proxy Ix3)
-  describe "DangerousStencil" $ do
-    it "Ix1" $ property $ prop_DangerousStencil (Proxy :: Proxy Ix1)
-    it "Ix2" $ property $ prop_DangerousStencil (Proxy :: Proxy Ix2)
-    it "Ix3" $ property $ prop_DangerousStencil (Proxy :: Proxy Ix3)
-    it "Ix4" $ property $ prop_DangerousStencil (Proxy :: Proxy Ix4)
-
-
-stencilDirection :: (Default a, Unbox a, Manifest r Ix2 a) => Ix2 -> Array r Ix2 a -> Array U Ix2 a
-stencilDirection ix =
-  computeAs U . mapStencil (Fill def) (makeStencil (Sz 3) (1 :. 1) $ \f -> f ix)
-
-
-stencilCorners ::
-     (Default a, Unbox a, Manifest r Ix2 a) => Ix2 -> Ix2 -> Array r Ix2 a -> Array U Ix2 a
-stencilCorners ixC ix = computeAs U . mapStencil (Fill def) (makeStencil (Sz 3) ixC $ \f -> f ix)
-
-
-stencilConvolution :: Spec
-stencilConvolution = do
-  let xs3 :: Array U Ix1 Int
-      xs3 = [1, 2, 3]
-      xs3f f = f (-1) 1 . f 0 2 . f 1 3
-      xs4 :: Array U Ix1 Int
-      xs4 = [1, 2, 3, 4]
-      xs4f f = f (-2) 1 . f (-1) 2 . f 0 3 . f 1 4
-      ys :: Array U Ix1 Int
-      ys = [1, 2, 3, 4, 5]
-      ysConvXs3 = [4, 10, 16, 22, 22]
-      ysConvXs4 = [10, 20, 30, 34, 31]
-      ysCorrXs3 = [8, 14, 20, 26, 14]
-      ysCorrXs4 = [11, 20, 30, 40, 26]
-      ysConvXs4' = [4, 10, 20, 30, 34]
-      ysCorrXs4' = [20, 30, 40, 26, 14]
-      xs4f' f = f (-1) 1 . f 0 2 . f 1 3 . f 2 4
-      applyStencil s = computeAs U . mapStencil (Fill 0) s
-  describe "makeConvolutionStencilFromKernel" $ do
-    it "1x3" $ applyStencil (makeConvolutionStencilFromKernel xs3) ys `shouldBe` ysConvXs3
-    it "1x4" $ applyStencil (makeConvolutionStencilFromKernel xs4) ys `shouldBe` ysConvXs4
-  describe "makeCorrelationStencilFromKernel" $ do
-    it "1x3" $ applyStencil (makeCorrelationStencilFromKernel xs3) ys `shouldBe` ysCorrXs3
-    it "1x4" $ applyStencil (makeCorrelationStencilFromKernel xs4) ys `shouldBe` ysCorrXs4
-  describe "makeConvolutionStencil" $ do
-    it "1x3" $ applyStencil (makeConvolutionStencil (Sz1 3) 1 xs3f) ys `shouldBe` ysConvXs3
-    it "1x4" $ applyStencil (makeConvolutionStencil (Sz1 4) 2 xs4f) ys `shouldBe` ysConvXs4
-    it "1x4" $ applyStencil (makeConvolutionStencil (Sz1 4) 1 xs4f') ys `shouldBe` ysConvXs4'
-  describe "makeCorrelationStencil" $ do
-    it "1x3" $ applyStencil (makeCorrelationStencil (Sz1 3) 1 xs3f) ys `shouldBe` ysCorrXs3
-    it "1x4" $ applyStencil (makeCorrelationStencil (Sz1 4) 2 xs4f) ys `shouldBe` ysCorrXs4
-    it "1x4" $ applyStencil (makeCorrelationStencil (Sz1 4) 1 xs4f') ys `shouldBe` ysCorrXs4'
-  describe "makeConvolutionStencil == makeConvolutionStencilFromKernel" $ do
-    it "Sobel Horizontal" $
-      property $ \(arr :: Array U Ix2 Int) ->
-        applyStencil (makeConvolutionStencil (Sz 3) 1 sobelX) arr ===
-        applyStencil (makeConvolutionStencilFromKernel sobelKernelX) arr
-    it "1x3" $
-      property $ \(arr :: Array U Ix1 Int) ->
-        applyStencil (makeConvolutionStencil (Sz1 3) 1 xs3f) arr ===
-        applyStencil (makeConvolutionStencilFromKernel xs3) arr
-    it "1x4" $
-      property $ \(arr :: Array U Ix1 Int) ->
-        applyStencil (makeConvolutionStencil (Sz1 4) 2 xs4f) arr ===
-        applyStencil (makeConvolutionStencilFromKernel xs4) arr
-  describe "makeCorrelationStencil == makeCorrelationStencilFromKernel" $ do
-    it "Sobel Horizontal" $
-      property $ \(arr :: Array U Ix2 Int) ->
-        applyStencil (makeCorrelationStencil (Sz 3) 1 sobelX) arr ===
-        applyStencil (makeCorrelationStencilFromKernel sobelKernelX) arr
-    it "1x3" $
-      property $ \(arr :: Array U Ix1 Int) ->
-        applyStencil (makeCorrelationStencil (Sz1 3) 1 xs3f) arr ===
-        applyStencil (makeCorrelationStencilFromKernel xs3) arr
-    it "1x4" $
-      property $ \(arr :: Array U Ix1 Int) ->
-        applyStencil (makeCorrelationStencil (Sz1 4) 2 xs4f) arr ===
-        applyStencil (makeCorrelationStencilFromKernel xs4) arr
-  describe "makeConvolutionStencil == makeCorrelationStencil . rotate180" $ do
-    it "Sobel Horizontal" $
-      property $ \(arr :: Array U Ix2 Int) ->
-        applyStencil (makeConvolutionStencilFromKernel sobelKernelX) arr ===
-        applyStencil (makeCorrelationStencilFromKernel (rotate180 sobelKernelX)) arr
-    it "1x3" $
-      property $ \(arr :: Array U Ix1 Int) ->
-        applyStencil (makeConvolutionStencilFromKernel xs3) arr ===
-        applyStencil (makeCorrelationStencilFromKernel (rotate180 xs3)) arr
-    -- it "1x4" $
-    --   property $ \(arr :: Array U Ix1 Int) ->
-    --     applyStencil (makeConvolutionStencilFromKernel xs4) arr ===
-    --     applyStencil (makeCorrelationStencilFromKernel (rotate180 xs4)) arr
-
-spec :: Spec
-spec = do
-  describe "Stencil" $ do
-    stencilSpec
-    let arr = [[1, 2, 3], [4, 5, 6], [7, 8, 9]] :: Array U Ix2 Int
-    describe "Unit tests Ix2" $ do
-      it "Direction Left" $
-        stencilDirection (0 :. 1) arr `shouldBe` [[2, 3, 0], [5, 6, 0], [8, 9, 0]]
-      it "Direction Right" $
-        stencilDirection (0 :. -1) arr `shouldBe` [[0, 1, 2], [0, 4, 5], [0, 7, 8]]
-      it "Direction Down" $
-        stencilDirection (1 :. 0) arr `shouldBe` [[4, 5, 6], [7, 8, 9], [0, 0, 0]]
-      it "Direction Up" $
-        stencilDirection (-1 :. 0) arr `shouldBe` [[0, 0, 0], [1, 2, 3], [4, 5, 6]]
-      it "Direction Left/Top Corner" $
-        stencilCorners (0 :. 0) (2 :. 2) arr `shouldBe` [[9, 0, 0], [0, 0, 0], [0, 0, 0]]
-      it "Direction Right/Top Corner" $
-        stencilCorners (0 :. 2) (2 :. -2) arr `shouldBe` [[0, 0, 7], [0, 0, 0], [0, 0, 0]]
-      it "Direction Right/Bottom Corner" $
-        stencilCorners (2 :. 2) (-2 :. -2) arr `shouldBe` [[0, 0, 0], [0, 0, 0], [0, 0, 1]]
-      it "Direction Left/Bottom Corner" $
-        stencilCorners (2 :. 0) (-2 :. 2) arr `shouldBe` [[0, 0, 0], [0, 0, 0], [3, 0, 0]]
-    describe "mapStencil with stride" $ do
-      let kernel = [[-1, 0, 1], [0, 1, 0], [-1, 0, 1]] :: Array U Ix2 Int
-          stencil = makeConvolutionStencilFromKernel kernel
-          stride = Stride 2
-      it "map stencil with stride on small array" $
-        let strideArr = mapStencil (Fill 0) stencil arr
-         in computeWithStrideAs U stride strideArr `shouldBe` [[-4, 8], [2, 14]]
-      it "map stencil with stride on larger array" $
-        let largeArr = makeArrayR U Seq (Sz 5) (succ . toLinearIndex (Sz 5))
-            strideArr = mapStencil (Fill 0) stencil largeArr
-         in computeWithStrideAs U stride strideArr `shouldBe`
-            [[-6, 1, 14], [-13, 9, 43], [4, 21, 44]]
-  stencilConvolution
-
-sobelX :: Num e => (Ix2 -> e -> e -> e) -> e -> e
-sobelX f = f (-1 :. -1) (-1) . f (-1 :. 1) 1 .
-           f ( 0 :. -1) (-2) . f ( 0 :. 1) 2 .
-           f ( 1 :. -1) (-1) . f ( 1 :. 1) 1
-
-sobelKernelX :: Array U Ix2 Int
-sobelKernelX = [ [-1, 0, 1]
-               , [-2, 0, 2]
-               , [-1, 0, 1] ]
-
-rotate180 :: (Num ix, Index ix) => Array U ix Int -> Array U ix Int
-rotate180 = computeAs U . transform' (\sz -> (sz, sz)) (\(Sz sz) f ix -> f (sz - 1 - ix))
diff --git a/tests/Data/Massiv/ArraySpec.hs b/tests/Data/Massiv/ArraySpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/ArraySpec.hs
+++ /dev/null
@@ -1,141 +0,0 @@
-{-# LANGUAGE AllowAmbiguousTypes #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE MonoLocalBinds #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE TypeApplications #-}
--- Here are contained tests for all instances for all main classes
-module Data.Massiv.ArraySpec
-  ( spec
-  ) where
-
-import Data.Massiv.Core
-import Data.Massiv.CoreArbitrary as A
-
-
-prop_Construct_makeArray_Manifest ::
-     forall r ix. (Load D ix Int, Ragged L ix Int, Source r ix Int, Construct r ix Int)
-  => Comp
-  -> Sz ix
-  -> Fun Int Int
-  -> Property
-prop_Construct_makeArray_Manifest comp sz f =
-  makeArrayLinearR D comp sz (apply f) ===
-  delay (setComp Seq (makeArray comp sz (apply f . toLinearIndex sz) :: Array r ix Int))
-
-prop_Construct_makeArray_Delayed ::
-     forall r ix. (Load D ix Int, Ragged L ix Int, Load r ix Int, Construct r ix Int)
-  => Comp
-  -> Sz ix
-  -> Fun Int Int
-  -> Property
-prop_Construct_makeArray_Delayed comp sz f =
-  makeArrayLinearR P comp sz (apply f) ===
-  compute (setComp Seq (makeArrayLinear comp sz (apply f)) :: Array r ix Int)
-
-prop_Functor ::
-     forall r ix.
-     (Load D ix Int, Ragged L ix Int, Load r ix Int, Construct r ix Int, Functor (Array r ix))
-  => Comp
-  -> Sz ix
-  -> Fun Int Int
-  -> Fun Int Int
-  -> Property
-prop_Functor comp sz f g =
-  makeArrayLinearR P comp sz (apply g . apply f) ===
-  compute (fmap (apply g) (makeArrayLinear comp sz (apply f) :: Array r ix Int))
-
-prop_Extract ::
-     forall r ix.
-     ( Load D ix Int
-     , Ragged L ix Int
-     , Load (EltRepr r ix) ix Int
-     , Construct r ix Int
-     , Extract r ix Int
-     )
-  => Comp
-  -> Sz ix
-  -> Fun Int Int
-  -> ix
-  -> Sz ix
-  -> Property
-prop_Extract comp sz f start newSize =
-  (computeAs P <$> toStringException (extractM start newSize arrD))  ===
-  (compute <$> toStringException (extractM start newSize arr))
-  where
-    arrD = makeArrayLinearR D comp sz (apply f)
-    arr = makeArrayLinear comp sz (apply f) :: Array r ix Int
-
-prop_IxUnbox ::
-     forall ix.
-     ( Load D ix ix
-     , Ragged L ix ix
-     , Construct U ix ix
-     , Source U ix ix
-     )
-  => Comp
-  -> Sz ix
-  -> Fun Int ix
-  -> Property
-prop_IxUnbox comp sz f =
-  makeArrayLinearR D comp sz (apply f) ===
-  delay (makeArrayLinear comp sz (apply f) :: Array U ix ix)
-
-prop_computeWithStride ::
-     forall r ix. (Load D ix Int, Ragged L ix Int, StrideLoad r ix Int, Construct r ix Int)
-  => Comp
-  -> Sz ix
-  -> Fun Int Int
-  -> Stride ix
-  -> Property
-prop_computeWithStride comp sz f stride =
-  arr === computeWithStride stride arrL .&&.
-  arr === compute (fromStrideLoad stride arrL)
-  where
-    arrL = makeArrayLinear comp sz (apply f) :: Array r ix Int
-    arr = computeWithStrideAs P stride (makeArrayLinearR D comp sz (apply f))
-
-
-specCommon ::
-     forall ix.
-     (Arbitrary ix, Load D ix Int, StrideLoad DW ix Int, Ragged L ix Int, Ragged L ix ix, Unbox ix)
-  => Spec
-specCommon =
-  describe "Construct" $ do
-    it "Construct_makeArray B" $ property $ prop_Construct_makeArray_Manifest @B @ix
-    it "Construct_makeArray N" $ property $ prop_Construct_makeArray_Manifest @N @ix
-    it "Construct_makeArray S" $ property $ prop_Construct_makeArray_Manifest @S @ix
-    it "Construct_makeArray P" $ property $ prop_Construct_makeArray_Manifest @P @ix
-    it "Construct_makeArray U" $ property $ prop_Construct_makeArray_Manifest @U @ix
-    it "Construct_makeArray M" $ property $ prop_Construct_makeArray_Manifest @M @ix
-    it "Construct_makeArray_Delayed DI" $ property $ prop_Construct_makeArray_Delayed @DI @ix
-    it "Construct_makeArray_Delayed DL" $ property $ prop_Construct_makeArray_Delayed @DL @ix
-    it "Construct_makeArray_Delayed DW" $ property $ prop_Construct_makeArray_Delayed @DW @ix
-    it "Construct_makeArray_Delayed M" $ property $ prop_Construct_makeArray_Delayed @M @ix
-    it "Functor D" $ property $ prop_Functor @D @ix
-    it "Functor DI" $ property $ prop_Functor @DI @ix
-    it "Functor DL" $ property $ prop_Functor @DL @ix
-    it "Functor DW" $ property $ prop_Functor @DW @ix
-    it "Extract DI" $ property $ prop_Extract @DI @ix
-    it "Extract D" $ property $ prop_Extract @D @ix
-    it "Extract B" $ property $ prop_Extract @B @ix
-    it "Extract N" $ property $ prop_Extract @N @ix
-    it "Extract S" $ property $ prop_Extract @S @ix
-    it "Extract U" $ property $ prop_Extract @U @ix
-    it "Extract M" $ property $ prop_Extract @M @ix
-    it "computeWithStride DI" $ property $ prop_computeWithStride @DI @ix
-    it "computeWithStride DW" $ property $ prop_computeWithStride @DW @ix
-    it "computeWithStride B" $ property $ prop_computeWithStride @B @ix
-    it "computeWithStride N" $ property $ prop_computeWithStride @N @ix
-    it "computeWithStride S" $ property $ prop_computeWithStride @S @ix
-    it "computeWithStride U" $ property $ prop_computeWithStride @U @ix
-    it "computeWithStride M" $ property $ prop_computeWithStride @M @ix
-    it "IxUnbox" $ property $ prop_IxUnbox @ix
-
-
-spec :: Spec
-spec = do
-  specCommon @Ix1
-  specCommon @Ix2
-  specCommon @Ix3
-  specCommon @Ix4
-  specCommon @Ix5
diff --git a/tests/Data/Massiv/Core/IndexSpec.hs b/tests/Data/Massiv/Core/IndexSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Core/IndexSpec.hs
+++ /dev/null
@@ -1,415 +0,0 @@
-{-# LANGUAGE AllowAmbiguousTypes #-}
-{-# LANGUAGE ExplicitNamespaces #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE GADTs #-}
-{-# LANGUAGE PatternSynonyms #-}
-{-# LANGUAGE RankNTypes #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE TypeApplications #-}
-{-# LANGUAGE TypeOperators #-}
-module Data.Massiv.Core.IndexSpec (SzNE(..), SzIx(..), DimIx(..), spec) where
-
-import Control.Monad
-import Data.Functor.Identity
-import Data.Massiv.Core.Index
-import Data.Proxy
-import Test.Hspec
-import Test.QuickCheck
-import Test.QuickCheck.Function
-import Data.Typeable
-
--- | Size that will result in a non-empty array
-newtype SzNE ix = SzNE (Sz ix) deriving Show
-
--- | Dimension that is always within bounds of an index
-newtype DimIx ix = DimIx Dim deriving Show
-
--- | Non-empty size together with an index that is within bounds of that index.
-data SzIx ix = SzIx (Sz ix) ix deriving Show
-
-instance (Index ix, Arbitrary ix) => Arbitrary (Sz ix) where
-  arbitrary = do
-    sz <- Sz . liftIndex abs <$> arbitrary
-    if totalElem sz > 200000
-      then arbitrary
-      else return sz
-
-instance (Index ix, Arbitrary ix) => Arbitrary (SzNE ix) where
-  arbitrary = SzNE . Sz . liftIndex (+1) . unSz <$> arbitrary
-
-instance (Index ix, Arbitrary ix) => Arbitrary (Stride ix) where
-  arbitrary = do
-    Positive (Small x) <- arbitrary
-    Stride . liftIndex ((+1) . (`mod` min 6 x)) <$> arbitrary
-
-instance (Index ix, Arbitrary ix) => Arbitrary (SzIx ix) where
-  arbitrary = do
-    SzNE sz <- arbitrary
-    -- Make sure index is within bounds:
-    SzIx sz . flip (liftIndex2 mod) (unSz sz) <$> arbitrary
-
-
-instance Arbitrary e => Arbitrary (Border e) where
-  arbitrary =
-    oneof
-      [ Fill <$> arbitrary
-      , return Wrap
-      , return Edge
-      , return Reflect
-      , return Continue
-      ]
-
-
-instance Index ix => Arbitrary (DimIx ix) where
-  arbitrary = do
-    n <- arbitrary
-    return $ DimIx (1 + (Dim n `mod` dimensions (Proxy :: Proxy ix)))
-
-arbitraryIntIx :: Gen Int
-arbitraryIntIx = sized (\s -> resize (floor $ (sqrt :: Double -> Double) $ fromIntegral s) arbitrary)
-  -- Generators are quadratic in QuickCheck.
-
-
-instance Arbitrary Ix2 where
-  arbitrary = (:.) <$> arbitraryIntIx <*> arbitraryIntIx
-
-instance Arbitrary Ix3 where
-  arbitrary = (:>) <$> arbitraryIntIx <*> ((:.) <$> arbitraryIntIx <*> arbitraryIntIx)
-
-instance Arbitrary Ix4 where
-  arbitrary = (:>) <$> arbitraryIntIx <*> arbitrary
-
-instance Arbitrary Ix5 where
-  arbitrary = (:>) <$> arbitraryIntIx <*> arbitrary
-
-instance CoArbitrary Ix2 where
-  coarbitrary (i :. j) = coarbitrary i . coarbitrary j
-
-instance CoArbitrary Ix3 where
-  coarbitrary (i :> ix) = coarbitrary i . coarbitrary ix
-
-instance CoArbitrary Ix4 where
-  coarbitrary (i :> ix) = coarbitrary i . coarbitrary ix
-
-instance CoArbitrary Ix5 where
-  coarbitrary (i :> ix) = coarbitrary i . coarbitrary ix
-
-instance Function Ix2 where
-  function = functionMap fromIx2 toIx2
-
-instance Function Ix3 where
-  function = functionMap fromIx3 toIx3
-
-instance Function Ix4 where
-  function = functionMap fromIx4 toIx4
-
-instance Function Ix5 where
-  function = functionMap fromIx5 toIx5
-
-
-prop_IsSafeIx :: Index ix => proxy ix -> SzIx ix -> Bool
-prop_IsSafeIx _ (SzIx sz ix) = isSafeIndex sz ix
-
-prop_RepairSafeIx :: Index ix => proxy ix -> SzIx ix -> Bool
-prop_RepairSafeIx _ (SzIx sz ix) =
-  ix == repairIndex sz ix (error "Impossible") (error "Impossible")
-
-prop_UnconsCons :: Index ix => proxy ix -> ix -> Bool
-prop_UnconsCons _ ix = ix == uncurry consDim (unconsDim ix)
-
-prop_UnsnocSnoc :: Index ix => proxy ix -> ix -> Bool
-prop_UnsnocSnoc _ ix = ix == uncurry snocDim (unsnocDim ix)
-
-prop_ToFromLinearIndex :: Index ix => proxy ix -> SzIx ix -> Property
-prop_ToFromLinearIndex _ (SzIx sz ix) =
-  isSafeIndex sz ix ==> ix == fromLinearIndex sz (toLinearIndex sz ix)
-
-prop_FromToLinearIndex :: Index ix => proxy ix -> SzNE ix -> NonNegative Int -> Property
-prop_FromToLinearIndex _ (SzNE sz) (NonNegative i) =
-  totalElem sz >= i ==> i == toLinearIndex sz (fromLinearIndex sz i)
-
-prop_CountElements :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_CountElements _ thresh sz =
-  totalElem sz < thresh ==> totalElem sz ==
-  iter zeroIndex (unSz sz) (pureIndex 1) (<) 0 (const (+ 1))
-
-prop_IterMonotonic :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonic _ thresh sz =
-  totalElem sz < thresh ==> fst $
-  iter (liftIndex succ zeroIndex) (unSz sz) (pureIndex 1) (<) (True, zeroIndex) mono
-  where
-    mono curIx (prevMono, prevIx) =
-      let isMono = prevMono && prevIx < curIx
-       in isMono `seq` (isMono, curIx)
-
-
-prop_IterMonotonic' :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonic' _ thresh sz =
-  totalElem sz <
-  thresh ==>
-  if isM
-    then isM
-    else error (show a)
-  where
-    (isM, a, _) =
-      iter (liftIndex succ zeroIndex) (unSz sz) (pureIndex 1) (<) (True, [], zeroIndex) mono
-    mono curIx (prevMono, acc, prevIx) =
-      let nAcc = (prevIx, curIx, prevIx < curIx) : acc
-          isMono = prevMono && prevIx < curIx
-       in isMono `seq` (isMono, nAcc, curIx)
-
-
-prop_IterMonotonicBackwards' :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonicBackwards' _ thresh sz@(Sz szix) =
-  totalElem sz <
-  thresh ==>
-  if isM
-    then isM
-    else error (show a)
-  where
-    (isM, a, _) = iter (liftIndex pred szix) zeroIndex (pureIndex (-1)) (>=) (True, [], szix) mono
-    mono curIx (prevMono, acc, prevIx) =
-      let isMono = prevMono && prevIx > curIx
-          nAcc = (prevIx, curIx, prevIx > curIx) : acc
-       in isMono `seq` (isMono, nAcc, curIx)
-
-prop_IterMonotonicM :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonicM _ thresh sz =
-  totalElem sz < thresh ==> fst $
-  runIdentity $ iterM (liftIndex succ zeroIndex) (unSz sz) (pureIndex 1) (<) (True, zeroIndex) mono
-  where
-    mono curIx (prevMono, prevIx) =
-      let isMono = prevMono && prevIx < curIx
-       in return $ isMono `seq` (isMono, curIx)
-
-
-prop_IterMonotonicBackwards :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonicBackwards _ thresh sz@(Sz szix) =
-  totalElem sz < thresh ==> fst $
-  iter (liftIndex pred szix) zeroIndex (pureIndex (-1)) (>=) (True, szix) mono
-  where
-    mono curIx (prevMono, prevIx) =
-      let isMono = prevMono && prevIx > curIx
-       in isMono `seq` (isMono, curIx)
-
-prop_IterMonotonicBackwardsM :: Index ix => proxy ix -> Int -> Sz ix -> Property
-prop_IterMonotonicBackwardsM _ thresh sz@(Sz szix) =
-  totalElem sz < thresh ==> fst $
-  runIdentity $ iterM (liftIndex pred szix) zeroIndex (pureIndex (-1)) (>=) (True, szix) mono
-  where
-    mono curIx (prevMono, prevIx) =
-      let isMono = prevMono && prevIx > curIx
-       in return $ isMono `seq` (isMono, curIx)
-
-prop_LiftLift2 :: Index ix => proxy ix -> ix -> Int -> Bool
-prop_LiftLift2 _ ix delta = liftIndex2 (+) ix (liftIndex (+delta) zeroIndex) ==
-                            liftIndex (+delta) ix
-
-
-prop_BorderRepairSafe :: Index ix => proxy ix -> Border ix -> SzNE ix -> ix -> Property
-prop_BorderRepairSafe _ border@(Fill defIx) (SzNE sz) ix =
-  not (isSafeIndex sz ix) ==> handleBorderIndex border sz id ix == defIx
-prop_BorderRepairSafe _ border (SzNE sz) ix =
-  not (isSafeIndex sz ix) ==> isSafeIndex sz (handleBorderIndex border sz id ix)
-
-
-prop_GetDropInsert :: Index ix => proxy ix -> NonNegative Int -> ix -> Property
-prop_GetDropInsert _ (NonNegative d) ix =
-  expected === do
-    i <- getDimM ix dim
-    ixL <- dropDimM ix dim
-    insertDimM ixL dim i
-  where expected = if d >= 1 && dim <= dimensions (Just ix) then Just ix else Nothing
-        dim = Dim d
-
-prop_PullOutInsert :: Index ix => proxy ix -> NonNegative Int -> ix -> Property
-prop_PullOutInsert _ (NonNegative d) ix =
-  expected === do
-    (i, ixL) <- pullOutDimM ix dim
-    insertDimM ixL dim i
-  where expected = if d >= 1 && dim <= dimensions (Just ix) then Just ix else Nothing
-        dim = Dim d
-
-prop_UnconsGetDrop :: (Index (Lower ix), Index ix) => proxy ix -> ix -> Property
-prop_UnconsGetDrop _ ix =
-  Just (unconsDim ix) === do
-    i <- getDimM ix (dimensions (Just ix))
-    ixL <- dropDimM ix (dimensions (Just ix))
-    return (i, ixL)
-
-prop_UnsnocGetDrop :: (Index (Lower ix), Index ix) => proxy ix -> ix -> Property
-prop_UnsnocGetDrop _ ix =
-  Just (unsnocDim ix) === do
-    i <- getDimM ix 1
-    ixL <- dropDimM ix 1
-    return (ixL, i)
-
-prop_SetAll :: Index ix => proxy ix -> ix -> Int -> Bool
-prop_SetAll _ ix i =
-  foldM (\cix d -> setDimM cix d i) ix ([1 .. dimensions (Just ix)] :: [Dim]) ==
-  Just (pureIndex i)
-
-
-prop_SetGet :: Index ix => proxy ix -> ix -> DimIx ix -> Int -> Bool
-prop_SetGet _ ix (DimIx dim) n = Just n == (setDimM ix dim n >>= (`getDimM` dim))
-
-
-prop_BorderIx1 :: Positive Int -> Border Char -> Fun Ix1 Char -> SzNE Ix1 -> Ix1 -> Bool
-prop_BorderIx1 (Positive period) border getVal (SzNE sz) ix =
-  if isSafeIndex sz ix
-    then apply getVal ix == val
-    else case border of
-           Fill defVal -> defVal == val
-           Wrap ->
-             val ==
-             handleBorderIndex
-               border
-               sz
-               (apply getVal)
-               (liftIndex2 (+) (liftIndex (* period) (unSz sz)) ix)
-           Edge ->
-             if ix < 0
-               then val == apply getVal (liftIndex (max 0) ix)
-               else val ==
-                    apply getVal (liftIndex2 min (liftIndex (subtract 1) (unSz sz)) ix)
-           Reflect ->
-             val ==
-             handleBorderIndex
-               border
-               sz
-               (apply getVal)
-               (liftIndex2 (+) (liftIndex (* (2 * signum ix * period)) (unSz sz)) ix)
-           Continue ->
-             val ==
-             handleBorderIndex
-               Reflect
-               sz
-               (apply getVal)
-               (if ix < 0
-                  then ix - 1
-                  else ix + 1)
-  where
-    val = handleBorderIndex border sz (apply getVal) ix
-
-specDimN :: (Index ix, Arbitrary ix) => proxy ix -> Spec
-specDimN proxy = do
-  describe "Safety" $ do
-    it "isSafeIndex" $ property $ prop_IsSafeIx proxy
-    it "RepairSafeIx" $ property $ prop_RepairSafeIx proxy
-  describe "Lifting" $
-    it "Lift/Lift2" $ property $ prop_LiftLift2 proxy
-  describe "Linear" $ do
-    it "ToFromLinearIndex" $ property $ prop_ToFromLinearIndex proxy
-    it "FromToLinearIndex" $ property $ prop_FromToLinearIndex proxy
-  describe "Iterator" $ do
-    it "CountElements" $ property $ prop_CountElements proxy 2000000
-    it "Monotonic" $ property $ prop_IterMonotonic proxy 2000000
-    it "MonotonicBackwards" $ property $ prop_IterMonotonicBackwards proxy 2000000
-    it "MonotonicM" $ property $ prop_IterMonotonicM proxy 2000000
-    it "MonotonicBackwardsM" $ property $ prop_IterMonotonicBackwardsM proxy 2000000
-  describe "Border" $
-    it "BorderRepairSafe" $ property $ prop_BorderRepairSafe proxy
-  describe "SetGetDrop" $ do
-    it "SetAll" $ property $ prop_SetAll proxy
-    it "SetGet" $ property $ prop_SetGet proxy
-    it "GetDropInsert" $ property $ prop_GetDropInsert proxy
-    it "PullOutInsert" $ property $ prop_PullOutInsert proxy
-
-specDim2AndUp
-  :: (Index ix, Index (Lower ix), Arbitrary ix)
-  => proxy ix -> Spec
-specDim2AndUp proxy =
-  describe "Higher/Lower" $ do
-    it "UnconsCons" $ property $ prop_UnconsCons proxy
-    it "UnsnocSnoc" $ property $ prop_UnsnocSnoc proxy
-    it "UnconsGetDrop" $ property $ prop_UnconsGetDrop proxy
-    it "UnsnocGetDrop" $ property $ prop_UnsnocGetDrop proxy
-
-ixToList :: Index ix => ix -> [Int]
-ixToList = reverse . foldlIndex (flip (:)) []
-
-prop_BinaryNumIx ::
-  (Num ix, Index ix) => (forall n . Num n => n -> n -> n) -> ix -> ix -> Property
-prop_BinaryNumIx f ix1 ix2 = zipWith f (ixToList ix1) (ixToList ix2) === ixToList (f ix1 ix2)
-
-prop_UnaryNumIx ::
-  (Num ix, Index ix) => (forall n . Num n => n -> n) -> ix -> Property
-prop_UnaryNumIx f ix = map f (ixToList ix) === ixToList (f ix)
-
-prop_BinaryNumSz ::
-  (Num ix, Index ix) => (forall n . Num n => n -> n -> n) -> Sz ix -> Sz ix -> Property
-prop_BinaryNumSz f sz1 sz2 =
-  zipWith f' (ixToList (unSz sz1)) (ixToList (unSz sz2)) === ixToList (unSz (f sz1 sz2))
-  where
-    f' x y = max 0 (f x y)
-
-prop_UnaryNumSz ::
-  (Num ix, Index ix) => (forall n . Num n => n -> n) -> Sz ix -> Property
-prop_UnaryNumSz f sz = map f' (ixToList (unSz sz)) === ixToList (unSz (f sz))
-  where
-    f' = max 0 . f
-
-specClasses :: forall ix . (Typeable ix, Num ix, Index ix, Arbitrary ix) => Spec
-specClasses = do
-  describe ("Num (" ++ showsTypeRep (typeRep (Proxy :: Proxy ix)) ")") $ do
-    it "(+)" $ property $ prop_BinaryNumIx @ix (+)
-    it "(-)" $ property $ prop_BinaryNumIx @ix (-)
-    it "(*)" $ property $ prop_BinaryNumIx @ix (*)
-    it "negate" $ property $ prop_UnaryNumIx @ix negate
-    it "abs" $ property $ prop_UnaryNumIx @ix abs
-    it "signum" $ property $ prop_UnaryNumIx @ix signum
-    it "fromInteger" $ property $ \ (i :: Int) ->
-      (fromIntegral i :: ix) === liftIndex (const i) zeroIndex
-  describe ("Num (" ++ showsTypeRep (typeRep (Proxy :: Proxy (Sz ix))) ")") $ do
-    it "(+)" $ property $ prop_BinaryNumSz @ix (+)
-    it "(-)" $ property $ prop_BinaryNumSz @ix (-)
-    it "(*)" $ property $ prop_BinaryNumSz @ix (*)
-    it "negate" $ property $ prop_UnaryNumSz @ix negate
-    it "abs" $ property $ prop_UnaryNumSz @ix abs
-    it "signum" $ property $ prop_UnaryNumSz @ix signum
-    it "fromInteger" $ property $ \ (i :: Int) ->
-      (fromIntegral i :: ix) === liftIndex (const i) zeroIndex
-    it "fromIx" $ property $ \ (ix :: ix) -> unSz (Sz ix) == liftIndex (max 0) ix
-
-spec :: Spec
-spec = do
-  describe "Tuple based indices" $ do
-    describe "Ix1T" $ do
-      specDimN (Nothing :: Maybe Ix1T)
-      it "prop_BorderIx1" $ property prop_BorderIx1
-    describe "Ix2T" $ do
-      specDimN (Nothing :: Maybe Ix2T)
-      specDim2AndUp (Nothing :: Maybe Ix2T)
-    describe "Ix3T" $ do
-      specDimN (Nothing :: Maybe Ix3T)
-      specDim2AndUp (Nothing :: Maybe Ix3T)
-    describe "Ix4T" $ do
-      specDimN (Nothing :: Maybe Ix4T)
-      specDim2AndUp (Nothing :: Maybe Ix4T)
-    describe "Ix5T" $ do
-      specDimN (Nothing :: Maybe Ix5T)
-      specDim2AndUp (Nothing :: Maybe Ix5T)
-  describe "Specialized indices" $ do
-    describe "Ix2" $ do
-      -- These can be used to quickly debug monotonicity
-      it "Monotonic'" $
-        property $ prop_IterMonotonic' (Nothing :: Maybe Ix2) 20000
-      it "MonotonicBackwards'" $
-        property $ prop_IterMonotonicBackwards' (Nothing :: Maybe Ix2) 20000
-      specDimN (Nothing :: Maybe Ix2)
-      specDim2AndUp (Nothing :: Maybe Ix2)
-    describe "Ix3" $ do
-      specDimN (Nothing :: Maybe Ix3)
-      specDim2AndUp (Nothing :: Maybe Ix3)
-    describe "Ix4" $ do
-      specDimN (Nothing :: Maybe Ix4)
-      specDim2AndUp (Nothing :: Maybe Ix4)
-    describe "Ix5" $ do
-      specDimN (Nothing :: Maybe Ix5)
-      specDim2AndUp (Nothing :: Maybe Ix5)
-  specClasses @Ix1
-  specClasses @Ix2
-  specClasses @Ix3
-  specClasses @Ix4
-  specClasses @Ix5
diff --git a/tests/Data/Massiv/Core/SchedulerSpec.hs b/tests/Data/Massiv/Core/SchedulerSpec.hs
deleted file mode 100644
--- a/tests/Data/Massiv/Core/SchedulerSpec.hs
+++ /dev/null
@@ -1,47 +0,0 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-module Data.Massiv.Core.SchedulerSpec (spec) where
-
-import Control.Exception.Base (ArithException(DivideByZero))
-import Data.Massiv.CoreArbitrary as A
-import Prelude as P
-
-
--- | Ensure proper exception handling.
-prop_CatchDivideByZero :: ArrIx D Ix2 Int -> [Int] -> Property
-prop_CatchDivideByZero (ArrIx arr ix) caps =
-  assertException
-    (== DivideByZero)
-    (A.sum $
-     A.imap
-       (\ix' x ->
-          if ix == ix'
-            then x `div` 0
-            else x)
-       (setComp (ParOn caps) arr))
-
--- | Ensure proper exception handling in nested parallel computation
-prop_CatchNested :: ArrIx D Ix1 (ArrIxP D Ix1 Int) -> [Int] -> Property
-prop_CatchNested (ArrIx arr ix) caps =
-  assertException
-    (== DivideByZero)
-    (computeAs U $
-     A.map A.sum $
-     A.imap
-       (\ix' (ArrIxP iarr ixi) ->
-          if ix == ix'
-            then A.imap
-                   (\ixi' e ->
-                      if ixi == ixi'
-                        then e `div` 0
-                        else e)
-                   iarr
-            else iarr)
-       (setComp (ParOn caps) arr))
-
-
-spec :: Spec
-spec =
-  describe "Exceptions" $ do
-    it "CatchDivideByZero" $ property prop_CatchDivideByZero
-    it "CatchNested" $ property prop_CatchNested
diff --git a/tests/Data/Massiv/CoreArbitrary.hs b/tests/Data/Massiv/CoreArbitrary.hs
deleted file mode 100644
--- a/tests/Data/Massiv/CoreArbitrary.hs
+++ /dev/null
@@ -1,203 +0,0 @@
-{-# OPTIONS_GHC -Wno-redundant-constraints #-}
-{-# LANGUAGE CPP #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-{-# LANGUAGE GADTs               #-}
-{-# LANGUAGE ScopedTypeVariables #-}
-{-# LANGUAGE StandaloneDeriving #-}
-{-# LANGUAGE UndecidableInstances #-}
-module Data.Massiv.CoreArbitrary
-  ( Arr(..)
-  , ArrTiny(..)
-  , ArrTiny1(..)
-  , ArrIx(..)
-  , ArrP(..)
-  , ArrS(..)
-  , ArrIxP(..)
-  , SzIx(..)
-  , SzNE(..)
-  , DimIx(..)
-  , toIx
-  , assertException
-  , assertSomeException
-  , assertExceptionIO
-  , assertSomeExceptionIO
-  , toStringException
-  , Semigroup((<>))
-  , applyFun2Compat
-  , module X
-  ) where
-
-import Control.DeepSeq (NFData, deepseq)
-import UnliftIO.Exception (Exception, SomeException, catch, catchAny)
-import Data.Foldable as F
-import Data.Massiv.Array as X
-import Data.Massiv.Core.IndexSpec hiding (spec)
-import Data.Typeable
-import Test.Hspec as X
-import Test.QuickCheck as X hiding (resize)
-import Test.QuickCheck.Function as X
-import Test.QuickCheck.Monadic as X
-#if !MIN_VERSION_base(4,11,0)
-import Data.Semigroup
-#endif
-applyFun2Compat :: Fun (a, b) c -> (a -> b -> c)
-#if MIN_VERSION_QuickCheck(2,10,0)
-applyFun2Compat = applyFun2
-#else
-applyFun2Compat (Fun _ f) a b = f (a, b)
-instance Function Word where
-  function = functionMap fromIntegral fromInteger
-#endif
-
--- | Arbitrary non-empty array. Computation strategy can be either `Seq` or `Par`.
-newtype Arr r ix e = Arr {unArr :: Array r ix e}
-
--- | Arbitrary small and possibly empty array. Computation strategy can be either `Seq` or `Par`.
-newtype ArrTiny r ix e = ArrTiny {unArrTiny :: Array r ix e}
-
--- | Tiny but non-empty
-newtype ArrTiny1 r ix e = ArrTiny1 {unArrTiny1 :: Array r ix e}
-
-newtype ArrS r ix e = ArrS {unArrS :: Array r ix e}
-
-newtype ArrP r ix e = ArrP {unArrP :: Array r ix e}
-
--- | Arbitrary non-empty array with a valid index. Can be either `Seq` or `Par`
-data ArrIx r ix e = ArrIx (Array r ix e) ix
-
--- | Arbitrary non-empty array with a valid index and `Seq` computation strategy
-data ArrIxS r ix e = ArrIxS (Array r ix e) ix
-
--- | Arbitrary non-empty array with a valid index and `Par` computation strategy
-data ArrIxP r ix e = ArrIxP (Array r ix e) ix
-
-deriving instance (Show (Array r ix e)) => Show (Arr r ix e)
-deriving instance (Show (Array r ix e)) => Show (ArrTiny r ix e)
-deriving instance (Show (Array r ix e)) => Show (ArrTiny1 r ix e)
-deriving instance (Show (Array r ix e)) => Show (ArrS r ix e)
-deriving instance (Show (Array r ix e)) => Show (ArrP r ix e)
-deriving instance (Show (Array r ix e), Show ix) => Show (ArrIx r ix e)
-deriving instance (Show (Array r ix e), Show ix) => Show (ArrIxS r ix e)
-deriving instance (Show (Array r ix e), Show ix) => Show (ArrIxP r ix e)
-
-instance Arbitrary Comp where
-  arbitrary = oneof [pure Seq, fmap ParOn arbitrary]
-
-
--- | Arbitrary array
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (Array r ix e) where
-  arbitrary = do
-    sz <- arbitrary
-    func <- arbitrary
-    comp <- oneof [pure Seq, pure Par]
-    return $ makeArray comp sz func
-
-
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrTiny r ix e) where
-  arbitrary = do
-    sz <- unSz <$> arbitrary
-    func <- arbitrary
-    comp <- oneof [pure Seq, pure Par]
-    return $ ArrTiny $ makeArray comp (Sz (liftIndex (`mod` 10) sz)) func
-
--- | Arbitrary small and possibly empty array. Computation strategy can be either `Seq` or `Par`.
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrTiny1 r ix e) where
-  arbitrary = do
-    sz <- unSz <$> arbitrary
-    func <- arbitrary
-    comp <- oneof [pure Seq, pure Par]
-    return $ ArrTiny1 $ makeArray comp (Sz (liftIndex (succ . (`mod` 10)) sz)) func
-
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (Arr r ix e) where
-  arbitrary = do
-    SzNE sz <- arbitrary
-    func <- arbitrary
-    comp <- oneof [pure Seq, pure Par]
-    return $ Arr $ makeArray comp sz func
-
-
--- | Arbitrary non-empty array
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrS r ix e) where
-  arbitrary = do
-    SzNE sz <- arbitrary
-    ArrS . makeArray Seq sz <$> arbitrary
-
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrP r ix e) where
-  arbitrary = do
-    Arr arr <- arbitrary
-    return $ ArrP (setComp Par arr)
-
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrIx r ix e) where
-  arbitrary = do
-    SzIx sz ix <- arbitrary
-    func <- arbitrary
-    comp <- arbitrary
-    return $ ArrIx (makeArray comp sz func) ix
-
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrIxS r ix e) where
-  arbitrary = do
-    SzIx sz ix <- arbitrary
-    func <- arbitrary
-    return $ ArrIxS (makeArray Seq sz func) ix
-
-
--- | Arbitrary non-empty array with a valid index
-instance (CoArbitrary ix, Arbitrary ix, Typeable e, Construct r ix e, Arbitrary e) =>
-         Arbitrary (ArrIxP r ix e) where
-  arbitrary = do
-    ArrIx arrIx ix <- arbitrary
-    return $ ArrIxP (setComp Par arrIx) ix
-
-
-toIx :: (Dimensions ix' ~ Dimensions ix, Index ix', Index ix) => ix -> ix'
-toIx ix = F.foldl' setEachIndex zeroIndex [1.. dimensions (Sz ix)]
-  where setEachIndex ix' d = setDim' ix' d (getDim' ix d)
-
-
-assertException :: (NFData a, Exception exc) =>
-                   (exc -> Bool) -- ^ Return True if that is the exception that was expected
-                -> a -- ^ Value that should throw an exception, when fully evaluated
-                -> Property
-assertException isExc = assertExceptionIO isExc . pure
-
-
-assertSomeException :: NFData a => a -> Property
-assertSomeException = assertSomeExceptionIO . pure
-
-
-assertExceptionIO :: (NFData a, Exception exc) =>
-                     (exc -> Bool) -- ^ Return True if that is the exception that was expected
-                  -> IO a -- ^ IO Action that should throw an exception
-                  -> Property
-assertExceptionIO isExc action =
-  monadicIO $ do
-    assert =<<
-      run
-        (catch
-           (do res <- action
-               res `deepseq` return False)
-           (\exc -> displayException exc `deepseq` return (isExc exc)))
-
-assertSomeExceptionIO :: NFData a => IO a -> Property
-assertSomeExceptionIO action =
-  monadicIO $ do
-    assert =<<
-      run
-        (catchAny
-           (do res <- action
-               res `deepseq` return False)
-           (\exc -> displayException exc `deepseq` return True))
-
-
-toStringException :: Either SomeException a -> Either String a
-toStringException = either (Left . displayException) Right
diff --git a/tests/Main.hs b/tests/Main.hs
deleted file mode 100644
--- a/tests/Main.hs
+++ /dev/null
@@ -1,10 +0,0 @@
-module Main where
-
-import System.IO (BufferMode(LineBuffering), hSetBuffering, stdout)
-import Test.Hspec
-import Spec
-
-main :: IO ()
-main = do
-  hSetBuffering stdout LineBuffering
-  hspec spec
diff --git a/tests/Spec.hs b/tests/Spec.hs
deleted file mode 100644
--- a/tests/Spec.hs
+++ /dev/null
@@ -1,1 +0,0 @@
-{-# OPTIONS_GHC -F -pgmF hspec-discover -optF --no-main #-}