diff --git a/EigensystemNum.hs b/EigensystemNum.hs
--- a/EigensystemNum.hs
+++ b/EigensystemNum.hs
@@ -10,19 +10,20 @@
 matSqr x = mult x x
 
 powerIter :: (Fractional a, Ord a) => [[a]] -> [([[a]],[[a]])]
-powerIter x = tail (iterate
+powerIter x =
+  tail $ iterate
     (\(_,z)->let s=normalize (matSqr z) in (s,(mult x s)))
     ([],x)
-  )
 
 normalize :: (Fractional a, Ord a) => [[a]] -> [[a]]
-normalize x = map (map (/(matnorm1 x))) x
+normalize x = map (map (/ matnorm1 x)) x
 
 getGrowth :: (Fractional a, Ord a) => ([[a]],[[a]]) -> a
-getGrowth (x,y) = uncurry (/) (maximumBy
+getGrowth (x,y) =
+  uncurry (/) $
+  maximumBy
     (\(_,xc) (_,xa) -> compare (abs xc) (abs xa))
     (concat (zipWith zip y x))
-  )
 
 specRadApprox :: (Fractional a, Ord a) => [[a]] -> [a]
 specRadApprox = map getGrowth . powerIter
@@ -31,7 +32,8 @@
 eigenValuesApprox = map diagonals . iterate similar_to
 
 limit :: (Num a, Ord a) => a -> [a] -> a
-limit tol (x0:x1:xs) = if abs (x1-x0) < tol * abs x0
-                       then x0
-		       else limit tol (x1:xs)
+limit tol (x0:x1:xs) =
+    if abs (x1-x0) < tol * abs x0
+    then x0
+    else limit tol (x1:xs)
 limit _ _ = error "Only infinite sequences are allowed"
diff --git a/Makefile b/Makefile
--- a/Makefile
+++ b/Makefile
@@ -3,3 +3,14 @@
 
 %.html:	%.lhs
 	ln -s $< $@
+
+
+
+run-test:	update-test
+	runhaskell Setup.lhs configure --user --enable-tests
+	runhaskell Setup.lhs build
+	runhaskell Setup.lhs haddock
+	./dist/build/numeric-quest-test/numeric-quest-test
+
+update-test:
+	doctest-extract-0.1 -o test/ --executable-main=Test/Main.hs RowEchelon
diff --git a/RowEchelon.hs b/RowEchelon.hs
new file mode 100644
--- /dev/null
+++ b/RowEchelon.hs
@@ -0,0 +1,412 @@
+{- |
+<https://en.wikipedia.org/wiki/Row_echelon_form>
+
+Duplicate of @htam@ package.
+-}
+module RowEchelon (
+   Matrix(Matrix), matrixRows, matrixWidth, matrixHeight,
+   shortBesidesTall, (|||),
+   matrixValid, matrixFromRows,
+   identity, matrixProduct,
+   nullspace, reducedRowEchelon, rowReduction, layoutEchelonBlocks,
+   scatter,
+   Zipper0(Zipper0), Zipper1(Zipper1), altOuter,
+   ) where
+
+import Orthogonals (matrix_matrix)
+
+import qualified Data.Foldable as Fold
+import qualified Data.NonEmpty as NonEmpty
+import qualified Data.List.HT as ListHT
+import qualified Data.List as List
+import Data.NonEmpty ((!:))
+import Data.Maybe (fromMaybe)
+
+
+{- $setup
+>>> import RowEchelon as Matrix
+>>> import qualified Test.QuickCheck as QC
+>>> import Test.QuickCheck ((===))
+>>> import Control.Monad (replicateM)
+>>> import qualified Data.Foldable as Fold
+>>> import qualified Data.NonEmpty.Class as NonEmptyC
+>>> import qualified Data.NonEmpty as NonEmpty
+>>> import qualified Data.List.Match as Match
+>>> import qualified Data.List.HT as ListHT
+>>> import Data.NonEmpty ((!:))
+>>> import Data.Ratio ((%))
+>>>
+>>> genElementUniform, genElementNearZero :: QC.Gen Integer
+>>> genElementUniform = QC.choose (-10,10)
+>>> genElementNearZero = fmap (flip rem 11) $ QC.arbitrary
+>>>
+>>> genMatrixForSize :: Int -> Int -> QC.Gen (Matrix Integer)
+>>> genMatrixForSize m n = do
+>>>    fmap (Matrix n) $ replicateM m $ replicateM n genElementNearZero
+>>>
+>>> genMatrix :: QC.Gen (Matrix Integer)
+>>> genMatrix = do
+>>>    m <- QC.choose (0,10)
+>>>    n <- QC.choose (0,10)
+>>>    genMatrixForSize m n
+>>>
+>>> shrinkMatrix :: (Eq a) => Matrix a -> [Matrix a]
+>>> shrinkMatrix matrix@(Matrix width rows) =
+>>>    filter (matrix/=) $
+>>>    map (Matrix width . snd) (ListHT.removeEach rows)
+>>>    ++
+>>>    [Matrix (width-1) $ map (drop 1) rows]
+>>>
+>>> forMatrix :: (QC.Testable test) => (Matrix Rational -> test) -> QC.Property
+>>> forMatrix prop =
+>>>    QC.forAllShrink genMatrix shrinkMatrix (prop . rationalMatrix)
+>>>
+>>> rationalMatrix :: Matrix Integer -> Matrix Rational
+>>> rationalMatrix = fmap (%1)
+-}
+
+
+data Zipper0 a = Zipper0 [a] [a]
+   deriving (Eq, Show)
+
+data Zipper1 a = Zipper1 [a] a [a]
+   deriving (Eq, Show)
+
+instance Functor Zipper0 where
+   fmap f (Zipper0 xs ys) = Zipper0 (fmap f xs) (fmap f ys)
+
+instance Functor Zipper1 where
+   fmap f (Zipper1 xs y zs) = Zipper1 (fmap f xs) (f y) (fmap f zs)
+
+zipper0FromList :: [a] -> Zipper0 a
+zipper0FromList = Zipper0 []
+
+
+type List0 = []
+type List1 = NonEmpty.T List0
+
+-- data Matrix a = Matrix {matrixWidth :: Int, matrixRows :: [[a]]}
+data Matrix a = Matrix Int [[a]]
+   deriving (Eq, Show)
+
+instance Functor Matrix where
+   fmap f (Matrix width rows) = Matrix width $ map (map f) rows
+
+instance Foldable Matrix where
+   foldMap f (Matrix _width rows) = foldMap (foldMap f) rows
+
+matrixValid :: Matrix a -> Bool
+matrixValid (Matrix width rows) = all ((width==) . length) rows
+
+matrixWidth :: Matrix a -> Int
+matrixWidth (Matrix width _) = width
+
+matrixHeight :: Matrix a -> Int
+matrixHeight = length . matrixRows
+
+matrixRows :: Matrix a -> [[a]]
+matrixRows (Matrix _ rows) = rows
+
+matrixFromRows :: [[a]] -> Matrix a
+matrixFromRows rows =
+   Matrix (ListHT.switchL 0 (const . length) rows) rows
+
+infixr 3 |||
+
+-- | requires that both matrices have the same height
+(|||) :: Matrix a -> Matrix a -> Matrix a
+Matrix widthA a ||| Matrix widthB b =
+   Matrix (widthA+widthB) (zipWith (++) a b)
+
+{- |
+The expression @shortBesidesTall a b@ means:
+
+> /A  B_upper\
+> |          |
+> \0  B_lower/
+
+Matrix @a@ must be at most as tall as @b@.
+-}
+shortBesidesTall :: (Num a) => Matrix a -> Matrix a -> Matrix a
+shortBesidesTall a b =
+   Matrix (matrixWidth a + matrixWidth b) $
+   zipWith (++)
+      (matrixRows a ++ repeat (replicate (matrixWidth a) 0))
+      (matrixRows b)
+
+matrixProduct :: (Num a) => Matrix a -> Matrix a -> Matrix a
+matrixProduct (Matrix _ x) (Matrix width y) =
+   Matrix width $ matrix_matrix x $ List.transpose y
+
+{- |
+>>> nullspace (Matrix 2 []) :: Matrix Rational
+Matrix 2 [[1 % 1,0 % 1],[0 % 1,1 % 1]]
+
+
+prop> forMatrix $ matrixValid . nullspace
+
+
+prop> :{
+   forMatrix $ \matrix ->
+      matrixWidth matrix == matrixHeight (nullspace matrix)
+:}
+
+prop> :{
+   forMatrix $ \matrix ->
+      matrixWidth (nullspace matrix) <= matrixWidth matrix
+:}
+
+max 0 (width matrix - height matrix) <= width nullspace
+
+prop> :{
+   forMatrix $ \matrix ->
+      matrixWidth matrix <= matrixWidth (nullspace matrix) + matrixHeight matrix
+:}
+
+prop> :{
+   forMatrix $ \matrix ->
+      Fold.all (0==) $ matrixProduct matrix (nullspace matrix)
+:}
+-}
+nullspace :: (RealFrac a) => Matrix a -> Matrix a
+nullspace matrix =
+   let echelon = reducedRowEchelon matrix in
+   let flat = layoutEchelonBlocks echelon in
+   let nullDim = matrixWidth flat in
+   Matrix nullDim $
+   scatter
+      (fmap matrixWidth $ altOuter echelon)
+      (matrixRows $ fmap negate flat) $
+   matrixRows $ identity nullDim
+
+identity :: (Num a) => Int -> Matrix a
+identity n = Matrix n $ take n $ map (take n) $ iterate (0:) $ 1 : repeat 0
+
+
+-- cf. event-list
+data AlternatingList a b = AlternatingList a [(b,a)]
+   deriving (Eq)
+
+instance (Show a, Show b) => Show (AlternatingList a b) where
+   showsPrec p xs =
+      showParen (p>=5) $
+      flip
+         (altFoldr
+            (\a -> showsPrec 5 a . showString " ./ ")
+            (\b -> showsPrec 5 b . showString " /. "))
+         xs
+      .
+      showString "[]"
+
+altSingleton :: a -> AlternatingList a b
+altSingleton a = AlternatingList a []
+
+infixr 5 /. , ./
+
+(./) = altConsOuter
+(/.) = altConsInner
+
+altConsOuter, (./) :: a -> List0 (b,a) -> AlternatingList a b
+altConsOuter = AlternatingList
+
+altConsInner, (/.) :: b -> AlternatingList a b -> List0 (b,a)
+altConsInner b ~(AlternatingList a bas) = (b,a) : bas
+
+altOuter :: AlternatingList a b -> List1 a
+altOuter (AlternatingList a bas) = a !: map snd bas
+
+altFoldr :: (a -> c -> d) -> (b -> d -> c) -> c -> AlternatingList a b -> d
+altFoldr f g c (AlternatingList a0 bas) =
+   f a0 $ foldr (\(b,a) -> g b . f a) c bas
+
+
+{- |
+>>> scatter (3 !: 1 : 4 : 1 : []) "012" ['a'..'k']
+"abc0d1efgh2ijk"
+-}
+scatter :: List1 Int -> [a] -> [a] -> [a]
+scatter (NonEmpty.Cons k_ ks_) =
+   let go _k [] [] dst = dst
+       go k0 (k1:ks) (a:as) dst =
+         case splitAt k0 dst of
+            (prefix,suffix) -> prefix ++ a : go k1 ks as suffix
+       go _ _ _ _ = error "scatter: inconsistent lengths"
+   in go k_ ks_
+
+
+{- |
+> mapM_ print $ layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[1,3,0,5],[0,0,1,7::Rational]]
+
+>>> layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[0::Rational]]
+Matrix 1 []
+
+>>> layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[1,3,0,5::Rational]]
+Matrix 3 [[3 % 1,0 % 1,5 % 1]]
+
+>>> layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[0,1,3,5::Rational]]
+Matrix 3 [[0 % 1,3 % 1,5 % 1]]
+
+>>> layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[1,3,0,5],[0,0,1,7::Rational]]
+Matrix 2 [[3 % 1,5 % 1],[0 % 1,7 % 1]]
+
+prop> :{
+   forMatrix $ \matrix ->
+      (layoutEchelonBlocks $ reducedRowEchelon $
+         identity (matrixHeight matrix) ||| matrix)
+      ===
+      matrix
+:}
+
+Construct a matrix that is already in Echelon form
+and check that it is preserved by Echelon decomposition.
+
+prop> :{
+   QC.forAll (fmap (flip mod 10) . NonEmpty.mapTail (take 9) <$> QC.arbitrary) $
+      \blockWidths ->
+   QC.forAll (traverse (uncurry genMatrixForSize) $
+               NonEmptyC.zip (0!:[1..]) blockWidths) $
+      \blocksInt ->
+
+      let blocks :: NonEmpty.T [] (Matrix Rational)
+          blocks = fmap (fmap (%1)) blocksInt
+          unitBesidesBlock block =
+            Matrix (matrixWidth block + 1) $
+            zipWith (:)
+               (drop 1 $ Match.replicate (matrixRows block) 0 ++ [1])
+               (matrixRows block)
+          echelonWithInterleavedEye =
+            Fold.foldr1 shortBesidesTall $
+            NonEmpty.mapTail (map unitBesidesBlock) blocks
+      in
+
+      layoutEchelonBlocks (reducedRowEchelon echelonWithInterleavedEye)
+      ===
+      Fold.foldr1 shortBesidesTall blocks
+:}
+-}
+layoutEchelonBlocks :: (Num a) => AlternatingList (Matrix a) a -> Matrix a
+layoutEchelonBlocks = Fold.foldr1 shortBesidesTall . altOuter
+
+
+{-
+/I B\           /0\
+|   | * P * x = | |
+\0 0/           \0/
+
+B is a rectangular upper block triangular matrix with non-uniform block size.
+
+The solution is
+
+      /-B\
+P*X = |  |
+      \ I/
+
+such that
+
+        /-B\
+(I B) * |  | = I*(-B) + B*I = 0
+        \ I/
+
+This is the nullspace:
+
+      /-B\
+P^T * |  |
+      \ I/
+
+The null-space vectors have the form
+/* * *\
+|1 0 0|
+|* * *|
+|* * *|
+|* * *|
+|0 1 0|
+|* * *|
+|0 0 1|
+|* * *|
+\* * */
+-}
+{- |
+>>> reducedRowEchelon $ matrixFromRows [[2,1,3::Rational]]
+Matrix 0 [] ./ 2 % 1 /. Matrix 2 [[1 % 2,3 % 2]] ./ []
+
+>>> reducedRowEchelon $ matrixFromRows [[2],[1],[3::Rational]]
+Matrix 0 [] ./ 3 % 1 /. Matrix 0 [[]] ./ []
+
+>>> reducedRowEchelon $ matrixFromRows [[1,0,2],[0,1,3::Rational]]
+Matrix 0 [] ./ 1 % 1 /. Matrix 0 [[]] ./ 1 % 1 /. Matrix 1 [[2 % 1],[3 % 1]] ./ []
+
+>>> reducedRowEchelon $ matrixFromRows [[1,3,0,5],[0,0,1,7::Rational]]
+Matrix 0 [] ./ 1 % 1 /. Matrix 1 [[3 % 1]] ./ 1 % 1 /. Matrix 1 [[5 % 1],[7 % 1]] ./ []
+
+prop> :{
+   forMatrix $ \matrix ->
+      Fold.all matrixValid $ altOuter $ reducedRowEchelon matrix
+:}
+-}
+reducedRowEchelon :: (RealFrac a) => Matrix a -> AlternatingList (Matrix a) a
+reducedRowEchelon matrix@(Matrix _ []) = altSingleton matrix
+reducedRowEchelon (Matrix _ rows) =
+   let go matrix@(Zipper0 upper _lower) =
+         case rowReduction matrix of
+            (n, Nothing) -> altSingleton (Matrix n $ reverse upper)
+            (n, Just (hd, block, remaining)) ->
+               Matrix n block ./ hd /. go remaining
+   in go $ zipper0FromList rows
+
+{- |
+>>> rowReduction $ Zipper0 [[1,2,3,4,5],[6,7,8,9,10::Rational]] [[0,2,0,0,0],[3,0,0,0,0]]
+(0,Just (3 % 1,[[],[]],Zipper0 [[0 % 1,0 % 1,0 % 1,0 % 1],[2 % 1,3 % 1,4 % 1,5 % 1],[7 % 1,8 % 1,9 % 1,10 % 1]] [[2 % 1,0 % 1,0 % 1,0 % 1]]))
+-}
+rowReduction ::
+   (RealFrac a) => Zipper0 [a] -> (Int, Maybe (a, [[a]], Zipper0 [a]))
+rowReduction (Zipper0 upper []) =
+   (ListHT.switchL 0 (\row _rows -> length row) upper, Nothing)
+rowReduction (Zipper0 upper (focus:lower)) =
+   case shiftUntilNonZeroColumn 0 $ Zipper1 upper focus lower of
+      (n, Nothing) -> (n, Nothing)
+      (n, Just (Zipper1 upperNE (NonEmpty.Cons maxRowHead maxRow) lowerNE)) ->
+         (n, Just $
+            let maxRowNormalized = map (/maxRowHead) maxRow in
+            (maxRowHead,
+             map (take n) $ reverse upper,
+             Zipper0
+               (maxRowNormalized :
+                map (cancelRow maxRowNormalized) upperNE)
+               (map (cancelRow maxRowNormalized) lowerNE)))
+
+cancelRow :: (Num a) => List0 a -> List1 a -> List0 a
+cancelRow xs (NonEmpty.Cons y0 ys) = zipWith (-) ys $ map (y0*) xs
+
+
+shiftUntilNonZeroColumn ::
+   (Real a) => Int -> Zipper1 [a] -> (Int, Maybe (Zipper1 (List1 a)))
+shiftUntilNonZeroColumn n matrix =
+   case checkHeads matrix of
+      Nothing -> (n, Nothing)
+      Just matrixNE@(Zipper1 upperNE foucsNE lowerNE) ->
+         case maxHead (foucsNE!:lowerNE) of
+            (maxRow, lowerWithoutMax) ->
+               if NonEmpty.head maxRow == 0
+               then (shiftUntilNonZeroColumn $! (n+1)) $
+                    fmap NonEmpty.tail matrixNE
+               else (n, Just (Zipper1 upperNE maxRow lowerWithoutMax))
+
+maxHead :: (Real a) => List1 (List1 a) -> (List1 a, [List1 a])
+maxHead = NonEmpty.maximumKey (abs . NonEmpty.head . fst) . NonEmpty.removeEach
+
+checkHeads :: Zipper1 [a] -> Maybe (Zipper1 (List1 a))
+checkHeads (Zipper1 upper focus lower) =
+   case NonEmpty.fetch focus of
+      Nothing ->
+         if all null (upper++lower) then Nothing else error "row lengths differ"
+      Just xs -> Just $
+         Zipper1
+            (fromMaybe (error "upper row lengths differ") $ allNonEmpty upper)
+            xs
+            (fromMaybe (error "lower row lengths differ") $ allNonEmpty lower)
+
+allNonEmpty :: [List0 a] -> Maybe [List1 a]
+allNonEmpty xs =
+   case ListHT.partitionMaybe NonEmpty.fetch xs of
+      (ne,[]) -> Just ne
+      _ -> Nothing
diff --git a/numeric-quest.cabal b/numeric-quest.cabal
--- a/numeric-quest.cabal
+++ b/numeric-quest.cabal
@@ -1,5 +1,5 @@
 Name:           numeric-quest
-Version:        0.2.0.3
+Version:        0.2.1
 License:        GPL
 License-File:   LICENSE
 Author:         Jan Skibinski
@@ -25,7 +25,7 @@
   README
 
 Source-Repository this
-  Tag:         0.2.0.3
+  Tag:         0.2.1
   Type:        darcs
   Location:    http://code.haskell.org/~thielema/numeric-quest/
 
@@ -35,6 +35,8 @@
 
 Library
   Build-Depends:
+    non-empty >=0.3.2 && <0.4,
+    utility-ht >=0.0.2 && <0.1,
     prelude-compat >=0.0.0.1 && <0.1,
     array >=0.1 && <0.6,
     base >=3 && <5
@@ -50,4 +52,22 @@
      Orthogonals
      QuantumVector
      Roots
+     RowEchelon
      Tensor
+
+Test-Suite numeric-quest-test
+  Type: exitcode-stdio-1.0
+  Default-Language: Haskell98
+  GHC-Options:    -Wall
+  Hs-Source-Dirs: test
+  Other-Modules:  Test.RowEchelon
+  Main-Is: Test/Main.hs
+
+  Build-Depends:
+    numeric-quest,
+    doctest-exitcode-stdio >=0.0 && <0.1,
+    doctest-lib >=0.1 && <0.2,
+    QuickCheck >=2 && <3,
+    non-empty,
+    utility-ht,
+    base
diff --git a/test/Test/Main.hs b/test/Test/Main.hs
new file mode 100644
--- /dev/null
+++ b/test/Test/Main.hs
@@ -0,0 +1,10 @@
+-- Do not edit! Automatically created with doctest-extract.
+module Main where
+
+import qualified Test.RowEchelon
+
+import qualified Test.DocTest.Driver as DocTest
+
+main :: IO ()
+main = DocTest.run $ do
+    Test.RowEchelon.test
diff --git a/test/Test/RowEchelon.hs b/test/Test/RowEchelon.hs
new file mode 100644
--- /dev/null
+++ b/test/Test/RowEchelon.hs
@@ -0,0 +1,212 @@
+-- Do not edit! Automatically created with doctest-extract from RowEchelon.hs
+{-# LINE 26 "RowEchelon.hs" #-}
+
+module Test.RowEchelon where
+
+import Test.DocTest.Base
+import qualified Test.DocTest.Driver as DocTest
+
+{-# LINE 27 "RowEchelon.hs" #-}
+import     RowEchelon as Matrix
+import     qualified Test.QuickCheck as QC
+import     Test.QuickCheck ((===))
+import     Control.Monad (replicateM)
+import     qualified Data.Foldable as Fold
+import     qualified Data.NonEmpty.Class as NonEmptyC
+import     qualified Data.NonEmpty as NonEmpty
+import     qualified Data.List.Match as Match
+import     qualified Data.List.HT as ListHT
+import     Data.NonEmpty ((!:))
+import     Data.Ratio ((%))
+
+genElementUniform,     genElementNearZero :: QC.Gen Integer
+genElementUniform     = QC.choose (-10,10)
+genElementNearZero     = fmap (flip rem 11) $ QC.arbitrary
+
+genMatrixForSize     :: Int -> Int -> QC.Gen (Matrix Integer)
+genMatrixForSize     m n = do
+       fmap (Matrix n) $ replicateM m $ replicateM n genElementNearZero
+
+genMatrix     :: QC.Gen (Matrix Integer)
+genMatrix     = do
+       m <- QC.choose (0,10)
+       n <- QC.choose (0,10)
+       genMatrixForSize m n
+
+shrinkMatrix     :: (Eq a) => Matrix a -> [Matrix a]
+shrinkMatrix     matrix@(Matrix width rows) =
+       filter (matrix/=) $
+       map (Matrix width . snd) (ListHT.removeEach rows)
+       ++
+       [Matrix (width-1) $ map (drop 1) rows]
+
+forMatrix     :: (QC.Testable test) => (Matrix Rational -> test) -> QC.Property
+forMatrix     prop =
+       QC.forAllShrink genMatrix shrinkMatrix (prop . rationalMatrix)
+
+rationalMatrix     :: Matrix Integer -> Matrix Rational
+rationalMatrix     = fmap (%1)
+
+test :: DocTest.T ()
+test = do
+ DocTest.printPrefix "RowEchelon:142: "
+{-# LINE 142 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 142 "RowEchelon.hs" #-}
+    nullspace (Matrix 2 []) :: Matrix Rational
+  )
+  [ExpectedLine [LineChunk "Matrix 2 [[1 % 1,0 % 1],[0 % 1,1 % 1]]"]]
+ DocTest.printPrefix "RowEchelon:146: "
+{-# LINE 146 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 146 "RowEchelon.hs" #-}
+      forMatrix $ matrixValid . nullspace
+  )
+ DocTest.printPrefix "RowEchelon:149: "
+{-# LINE 149 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 149 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      matrixWidth matrix == matrixHeight (nullspace matrix)
+  )
+ DocTest.printPrefix "RowEchelon:154: "
+{-# LINE 154 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 154 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      matrixWidth (nullspace matrix) <= matrixWidth matrix
+  )
+ DocTest.printPrefix "RowEchelon:161: "
+{-# LINE 161 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 161 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      matrixWidth matrix <= matrixWidth (nullspace matrix) + matrixHeight matrix
+  )
+ DocTest.printPrefix "RowEchelon:166: "
+{-# LINE 166 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 166 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      Fold.all (0==) $ matrixProduct matrix (nullspace matrix)
+  )
+ DocTest.printPrefix "RowEchelon:224: "
+{-# LINE 224 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 224 "RowEchelon.hs" #-}
+    scatter (3 !: 1 : 4 : 1 : []) "012" ['a'..'k']
+  )
+  [ExpectedLine [LineChunk "\"abc0d1efgh2ijk\""]]
+ DocTest.printPrefix "RowEchelon:240: "
+{-# LINE 240 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 240 "RowEchelon.hs" #-}
+    layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[0::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 1 []"]]
+ DocTest.printPrefix "RowEchelon:243: "
+{-# LINE 243 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 243 "RowEchelon.hs" #-}
+    layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[1,3,0,5::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 3 [[3 % 1,0 % 1,5 % 1]]"]]
+ DocTest.printPrefix "RowEchelon:246: "
+{-# LINE 246 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 246 "RowEchelon.hs" #-}
+    layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[0,1,3,5::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 3 [[0 % 1,3 % 1,5 % 1]]"]]
+ DocTest.printPrefix "RowEchelon:249: "
+{-# LINE 249 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 249 "RowEchelon.hs" #-}
+    layoutEchelonBlocks $ reducedRowEchelon $ matrixFromRows [[1,3,0,5],[0,0,1,7::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 2 [[3 % 1,5 % 1],[0 % 1,7 % 1]]"]]
+ DocTest.printPrefix "RowEchelon:252: "
+{-# LINE 252 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 252 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      (layoutEchelonBlocks $ reducedRowEchelon $
+         identity (matrixHeight matrix) ||| matrix)
+      ===
+      matrix
+  )
+ DocTest.printPrefix "RowEchelon:263: "
+{-# LINE 263 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 263 "RowEchelon.hs" #-}
+        
+   QC.forAll (fmap (flip mod 10) . NonEmpty.mapTail (take 9) <$> QC.arbitrary) $
+      \blockWidths ->
+   QC.forAll (traverse (uncurry genMatrixForSize) $
+               NonEmptyC.zip (0!:[1..]) blockWidths) $
+      \blocksInt ->
+
+      let blocks :: NonEmpty.T [] (Matrix Rational)
+          blocks = fmap (fmap (%1)) blocksInt
+          unitBesidesBlock block =
+            Matrix (matrixWidth block + 1) $
+            zipWith (:)
+               (drop 1 $ Match.replicate (matrixRows block) 0 ++ [1])
+               (matrixRows block)
+          echelonWithInterleavedEye =
+            Fold.foldr1 shortBesidesTall $
+            NonEmpty.mapTail (map unitBesidesBlock) blocks
+      in
+
+      layoutEchelonBlocks (reducedRowEchelon echelonWithInterleavedEye)
+      ===
+      Fold.foldr1 shortBesidesTall blocks
+  )
+ DocTest.printPrefix "RowEchelon:329: "
+{-# LINE 329 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 329 "RowEchelon.hs" #-}
+    reducedRowEchelon $ matrixFromRows [[2,1,3::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 0 [] ./ 2 % 1 /. Matrix 2 [[1 % 2,3 % 2]] ./ []"]]
+ DocTest.printPrefix "RowEchelon:332: "
+{-# LINE 332 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 332 "RowEchelon.hs" #-}
+    reducedRowEchelon $ matrixFromRows [[2],[1],[3::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 0 [] ./ 3 % 1 /. Matrix 0 [[]] ./ []"]]
+ DocTest.printPrefix "RowEchelon:335: "
+{-# LINE 335 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 335 "RowEchelon.hs" #-}
+    reducedRowEchelon $ matrixFromRows [[1,0,2],[0,1,3::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 0 [] ./ 1 % 1 /. Matrix 0 [[]] ./ 1 % 1 /. Matrix 1 [[2 % 1],[3 % 1]] ./ []"]]
+ DocTest.printPrefix "RowEchelon:338: "
+{-# LINE 338 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 338 "RowEchelon.hs" #-}
+    reducedRowEchelon $ matrixFromRows [[1,3,0,5],[0,0,1,7::Rational]]
+  )
+  [ExpectedLine [LineChunk "Matrix 0 [] ./ 1 % 1 /. Matrix 1 [[3 % 1]] ./ 1 % 1 /. Matrix 1 [[5 % 1],[7 % 1]] ./ []"]]
+ DocTest.printPrefix "RowEchelon:341: "
+{-# LINE 341 "RowEchelon.hs" #-}
+ DocTest.property(
+{-# LINE 341 "RowEchelon.hs" #-}
+        
+   forMatrix $ \matrix ->
+      Fold.all matrixValid $ altOuter $ reducedRowEchelon matrix
+  )
+ DocTest.printPrefix "RowEchelon:357: "
+{-# LINE 357 "RowEchelon.hs" #-}
+ DocTest.example(
+{-# LINE 357 "RowEchelon.hs" #-}
+    rowReduction $ Zipper0 [[1,2,3,4,5],[6,7,8,9,10::Rational]] [[0,2,0,0,0],[3,0,0,0,0]]
+  )
+  [ExpectedLine [LineChunk "(0,Just (3 % 1,[[],[]],Zipper0 [[0 % 1,0 % 1,0 % 1,0 % 1],[2 % 1,3 % 1,4 % 1,5 % 1],[7 % 1,8 % 1,9 % 1,10 % 1]] [[2 % 1,0 % 1,0 % 1,0 % 1]]))"]]