HNumeric (empty) → 0.2.0.0

raw patch · 9 files changed

+760/−0 lines, 9 filesdep +HNumericdep +basedep +normaldistributionsetup-changed

Dependencies added: HNumeric, base, normaldistribution

Files

• ChangeLog.md +3/−0
• HNumeric.cabal +67/−0
• LICENSE +30/−0
• README.md +193/−0
• Setup.hs +2/−0
• app/Main.hs +12/−0
• src/HNum/Stats.hs +48/−0
• src/HNum/Vector.hs +403/−0
• test/Spec.hs +2/−0

@@ -0,0 +1,3 @@+# Changelog for HNumeric++## Unreleased changes

@@ -0,0 +1,67 @@+-- This file has been generated from package.yaml by hpack version 0.28.2.+--+-- see: https://github.com/sol/hpack+--+-- hash: 0a4c04fe65db3e146424c0f3450ab7549b8db8ec8c21cb1d08bc26575227481f++name:           HNumeric+version:        0.2.0.0+synopsis:       Haskell Numeric Library with pure functionality, R & MATLAB Syntax.+description:    Please see the README on GitHub at <https://github.com/Axect/HNumeric#readme>+category:       HNum, library, Numeric, LinearAlgebra, Statistics, bsd3+homepage:       https://github.com/Axect/HNumeric#readme+bug-reports:    https://github.com/Axect/HNumeric/issues+author:         Tae Geun Kim+maintainer:     edeftg@gmail.com+copyright:      2018 Tae Geun Kim+license:        BSD3+license-file:   LICENSE+build-type:     Simple+cabal-version:  >= 1.10+extra-source-files:+    ChangeLog.md+    README.md++source-repository head+  type: git+  location: https://github.com/Axect/HNumeric++library+  exposed-modules:+      HNum.Stats+      HNum.Vector+  other-modules:+      Paths_HNumeric+  hs-source-dirs:+      src+  build-depends:+      base >=4.7 && <5+    , normaldistribution+  default-language: Haskell2010++executable HNumeric-exe+  main-is: Main.hs+  other-modules:+      Paths_HNumeric+  hs-source-dirs:+      app+  ghc-options: -threaded -rtsopts -with-rtsopts=-N+  build-depends:+      HNumeric+    , base >=4.7 && <5+    , normaldistribution+  default-language: Haskell2010++test-suite HNumeric-test+  type: exitcode-stdio-1.0+  main-is: Spec.hs+  other-modules:+      Paths_HNumeric+  hs-source-dirs:+      test+  ghc-options: -threaded -rtsopts -with-rtsopts=-N+  build-depends:+      HNumeric+    , base >=4.7 && <5+    , normaldistribution+  default-language: Haskell2010

@@ -0,0 +1,30 @@+Copyright Tae Geun Kim (c) 2018++All rights reserved.++Redistribution and use in source and binary forms, with or without+modification, are permitted provided that the following conditions are met:++    * Redistributions of source code must retain the above copyright+      notice, this list of conditions and the following disclaimer.++    * Redistributions in binary form must reproduce the above+      copyright notice, this list of conditions and the following+      disclaimer in the documentation and/or other materials provided+      with the distribution.++    * Neither the name of Author name here nor the names of other+      contributors may be used to endorse or promote products derived+      from this software without specific prior written permission.++THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

@@ -0,0 +1,193 @@+# HNumeric++![Travis](https://travis-ci.org/Axect/HNumeric.svg?branch=master)++## Packages++* HNum.Vector : Contain vector, matrix, linear algebra+* HNum.Stats : Contain statistical functions++## Installation++### 1. Native Use++You can use this package just change `app/Main.hs`+Then, just type next command++```bash+git clone https://github.com/Axect/HNumeric+```++### 2. Cabal Install++First, install prerequisite.++```sh+cabal install normaldistribution+```+  +Second, install tarball & unpack.++```sh+export VERSION=0.2.0.0++wget -O HNumeric-${VERSION}.tar.gz https://github.com/Axect/HNumeric/blob/master/dist/HNumeric-${VERSION}.tar.gz\?raw\=true++tar -xvf HNumeric-${VERSION}.tar.gz+```++Finally, build & install++```sh+cabal build+cabal install+```++Finish!++### 3. Import to Stack project++If you use this package to your own project, then you should change `stack.yaml` and `package.yaml`++#### 1) Change `stack.yaml`++```yaml+# In stack.yaml+extra-deps:+  - git: https://github.com/Axect/HNumeric.git+    commit: [Latest Commit]+  - normaldistribution-1.1.0,3+```++* Replace [Latest Commit] with latest commit in [HNumeric Commit](https://github.com/Axect/HNumeric/commits/master)++#### 2) Change `package.yaml`++```yaml+# In package.yaml+dependecies:+- base+- HNumeric+- normaldistribution+```++Then enjoy!++## Usage++### Import Module++* HNum.Vector+* HNum.Stats++### Basic Vector Usage++```haskell+-- HNumeric-0.2.0.0 Documentation++let a = vector [1,2,3] -- Vector declaration+let b = vector [4,5,6]++-- Print Vector+print a++-- You can (+1) by fmap (Vector is functor)+(+1) <$> a ++-- Or MATLAB-like operator (.+, .-, .*, ./, .^)+a .+ 1++-- You can make list from vector+toList a -- [1, 2, 3]++-- You can make vector from list+fromList [1,2,3] -- Vector [1,2,3]++-- You can add (subtract, multiply, divide) vectors+a + b -- Vector [5,7,9]++-- Also dot product is here.+a .*. b -- 1*4 + 2*5 + 3*6 = 32++-- Declare Matrix (Syntactic Sugar)+let c = matrix [[1,2],[3,4]]++-- or Declare using R Syntax+let d = Matrix {val = Vector [5,6,7,8], row = 2, col = 2, byRow = True}++-- Determinant+det c++-- Inverse+inv c++-- Transpose+transpose c++-- Matrix ops with Constant (+, -, *, /, ^)+c .+ 1 -- Matrix [[2,3],[4,5]]++-- Matrix ops with Matrix (+, -)+c + c -- Matrix [[2,4],[6,8]]++-- Matrix Multiplication+c %*% d++-- Matrix - Inverse Multiplication+c %/% d++-- Vector Concatenate+hcat a b -- Vector [1,2,3,4,5,6]+vcat a b -- Matrix [[1,2,3],[4,5,6]]++-- Matrix Concatenate+hcat c d -- Matrix [[1,2,5,6],[3,4,7,8]]+vcat c d -- Matrix [[1,2],[3,4],[5,6],[7,8]]++-- Insert Vector to Matrix+Vector [1, 2] .: c -- Matrix [[1,2],[1,2],[3,4]]+```++### Basic Stats Usage++```haskell+-- Sample Vector (import Vector)+v = Vector [1..10]+w = Vector [10, 9 .. 1]++-- Mean+mean v++-- Var+var v++-- Std+std v++-- Cov Matrix+cov v w++-- Only Cov+cov' v w++-- Linear Fit+(intercept, slope) = lm v w -- (11.0, -1.0) -- (Intercept, Slope)++-- Linear Fit function+lineFit (intercept, slope) (Vector [1 .. 20])++-- RSS+rss v w++-- RSE+rse v w+```++### TODO++* ~~Effective Matrix Multiplication~~+* Write Vector to CSV+* ~~Haddock~~+* DataFrame using Map+* ~~Fix Matrix Implementation~~+* ~~Numeric Class Define~~

@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain

@@ -0,0 +1,12 @@+module Main where++import           HNum.Vector+import           HNum.Stats+import           Data.Random.Normal++main :: IO ()+main = do+  let a = Vector [1, 3, 4]+  print a+  print $ map ($ a) [mean, var, std]+

@@ -0,0 +1,48 @@+module HNum.Stats where++import           HNum.Vector++type Coeff a = (a, a)++-- | Expectation Value+mean :: Fractional a => Vector a -> a+mean v = sum v / fromIntegral (length v)++-- | Covariance (Single-Valued)+cov' :: Floating a => Vector a -> Vector a -> a+cov' x y+  | length x <= 1 || length y <= 1 = error "Samples are not enough"+  | length x /= length y = error "Length is not same"+  | otherwise = ((x .- mean x) .*. (y .- mean y)) / fromIntegral (length x - 1)++-- | Variance+var :: Floating a => Vector a -> a+var v = cov' v v++-- | Standard Deviation+std :: Floating a => Vector a -> a+std = sqrt . var++-- | Covariance Matrix+cov :: Floating a => Vector a -> Vector a -> Matrix a+cov x y = matrix [[var x, cov' x y], [cov' y x, var y]]++-- | Least Square Method - (Intercept, Slope)+lm :: Floating a => Vector a -> Vector a -> Coeff a+lm x y = (my - b1 * mx, b1)+ where+  mx = mean x+  my = mean y+  b1 = (x .- mx) .*. (y .- my) / ((x .- mx) .*. (x .- mx))++-- | Line Fitting with (Intercept, Slope) & Range of x+lineFit :: Floating a => Coeff a -> Vector a -> Vector a+lineFit (n, m) x = x .* m .+ n++-- | Residual Sum of Squares+rss :: Floating a => Vector a -> Vector a -> a+rss x y = sum ((y - lineFit (lm x y) x) .^ 2)++-- | Relative Standard Error+rse :: Floating a => Vector a -> Vector a -> a+rse x y = sqrt (1 / fromIntegral (length x - 2) * rss x y)

@@ -0,0 +1,403 @@+{-|+Module      : HNumeric.Vector+Description : Haskell Vector & Matrix & Linear Algebra Library to do machine learning+CopyRight   : (c) Tae Geun Kim, 2018+License     : GPL-3+Maintainer  : edeftg@gmail.com+Stability   : Experimental+-}+module HNum.Vector where++import           Data.Functor                   ( )+import           Control.Applicative            ( )++---------------------------------------------------+-- Vector+--------------------------------------------------+-- Type Section+newtype Vector a = Vector [a] deriving (Show, Eq)++vector :: [a] -> Vector a+vector = Vector++-- Instance Section+instance Functor Vector where+  fmap f (Vector x) = Vector (fmap f x)++instance Applicative Vector where+  pure a                  = Vector []+  Vector fs <*> Vector xs = Vector (zipWith ($) fs xs)++instance (Num a) => Num (Vector a) where+  negate v      = negate <$> v+  (+) v1 v2     = (+) <$> v1 <*> v2+  (*) v1 v2     = (*) <$> v1 <*> v2+  fromInteger n = fromInteger <$> Vector [n]+  signum v      = signum <$> v+  abs v         = abs <$> v++instance (Fractional a) => Fractional (Vector a) where+  recip v        = recip <$> v+  (/) v1 v2      = (*) <$> v1 <*> recip v2+  fromRational n = fromRational <$> Vector [n]++instance (Floating a) => Floating (Vector a) where+  pi = Vector [pi]+  exp v = exp <$> v+  log v = log <$> v+  sqrt v = sqrt <$> v+  sin v = sin <$> v+  cos v = cos <$> v+  tan v = tan <$> v+  asin v = asin <$> v+  acos v = acos <$> v+  atan v = atan <$> v+  sinh v = sinh <$> v+  cosh v = cosh <$> v+  tanh v = tanh <$> v+  asinh v = asinh <$> v+  acosh v = acosh <$> v+  atanh v = atanh <$> v++instance Foldable Vector where+  foldr _ z (Vector []) = z+  foldr f z (Vector xs) = foldr f z xs++  foldl _ z (Vector []) = z+  foldl f z (Vector xs) = foldl f z xs++-- Class for Vector with List+class List m where+  toList :: m a -> [a]+  fromList :: [a] -> m a++instance List Vector where+  toList (Vector xs) = xs+  fromList = Vector++---------------------------------------------------+-- Matrix+---------------------------------------------------++-- |Definition of Matrix+data Matrix a = Matrix {val :: Vector a, row :: Int, col :: Int, byRow :: Bool} deriving (Eq)++-- |matrix is syntactic sugar to create Matrix+matrix :: [[a]] -> Matrix a+matrix = formMat++-- |Matrices is necessary class for Matrix+class Matrices m where+  matForm :: m a -> [[a]]+  formMat :: [[a]] -> m a++instance Matrices Matrix where+  matForm (Matrix (Vector v) r c b)+    | r*c /= length v = error "Matrix Dimension mismatch!"+    | b = ctake c v+    | otherwise = dtake c v+      where ctake :: Int -> [a] -> [[a]]+            ctake _ [] = []+            ctake n m = take n m : ctake n (drop n m)+            dtake :: Int -> [a] -> [[a]]+            dtake _ [] = []+            dtake n m = [ptake n m r | r <- [0..(length m `div` n - 1)]]+            ptake n v r = [v !! x | x <- idx v, x `mod` (length v `div` n) == r]+            idx v = take (length v) [0..]+  formMat [] = Matrix (Vector []) 0 0 True+  formMat xs = Matrix (Vector (concat xs)) (length xs) (length (head xs)) True++instance Show a => Show (Matrix a) where+  show m = "Matrix " ++ show (matForm m)++instance Functor Matrix where+  fmap f mat = mat { val = fmap f (val mat) }++instance Applicative Matrix where+  pure a = matrix []+  mf <*> mx = mx { val = val mf <*> val mx }++instance Num a => Num (Matrix a) where+  negate m = negate <$> m+  (+) m n = (+) <$> m <*> n+  (*) m n = (*) <$> m <*> n+  fromInteger a = fromInteger <$> matrix [[a]]+  signum m = signum <$> m+  abs m = abs <$> m++instance Fractional a => Fractional (Matrix a) where+  recip m = recip <$> m+  (/) m n = (*) <$> m <*> recip n+  fromRational n = fromRational <$> matrix [[n]]++instance (Floating a) => Floating (Matrix a) where+  pi = matrix [[pi]]+  exp v = exp <$> v+  log v = log <$> v+  sqrt v = sqrt <$> v+  sin v = sin <$> v+  cos v = cos <$> v+  tan v = tan <$> v+  asin v = asin <$> v+  acos v = acos <$> v+  atan v = atan <$> v+  sinh v = sinh <$> v+  cosh v = cosh <$> v+  tanh v = tanh <$> v+  asinh v = asinh <$> v+  acosh v = acosh <$> v+  atanh v = atanh <$> v++instance Foldable Matrix where+  foldr _ z (Matrix (Vector []) _ _ _) = z+  foldr f z (Matrix (Vector xs) _ _ _) = foldr f z xs++  foldl _ z (Matrix (Vector []) _ _ _) = z+  foldl f z (Matrix (Vector xs) _ _ _) = foldl f z xs++---------------------------------------------------+-- Operation+---------------------------------------------------+{-|+   (.<ops>) is an operation Vector(or Matrix) with Constant.+   Dot means position of Vector.+   Example: a .* 2  = twice whole elements of a+            a .*. b = Dot product+-}+class Functor f => VecOps f where+  (.+) :: Num a => f a -> a -> f a+  (.-) :: Num a => f a -> a -> f a+  (.*) :: Num a => f a -> a -> f a+  (./) :: Fractional a => f a -> a -> f a+  (.^) :: Floating a => f a -> a -> f a+  (.*.) :: Num a => f a -> f a -> a+  norm :: Floating a => f a -> a++{-+   MatOps is just additional operations for Matrices.+-}+class Functor f => MatOps f where+  (%*%) :: Num a => f a -> f a -> f a+  (%/%) :: Fractional a => f a -> f a -> f a+  det :: Fractional a => f a -> a+  inv :: Fractional a => f a -> f a+  transpose :: f a -> f a++instance VecOps Vector where+  v .+ n = (+ n) <$> v+  v .- n = (+ negate n) <$> v+  v .* n = (* n) <$> v+  v ./ n = (/ n) <$> v+  v .^ n = (** n) <$> v+  v .*. w = sum $ v * w+  norm v = sqrt $ v .*. v++instance VecOps Matrix where+  v .+ n = (+ n) <$> v+  v .- n = (+ negate n) <$> v+  v .* n = (* n) <$> v+  v ./ n = (/ n) <$> v+  v .^ n = (** n) <$> v+  v .*. w = sum $ v * w+  norm v = sqrt $ v .*. v++instance MatOps Matrix where+  m %*% n | col m /= row n = error "Can't Multiply - Dimension mismatch!"+          | otherwise      = matrix $ matForm m %-*-% matForm n+  m %/% n = m %*% inv n+  det m | col m /= row m = error "Can't calculate determinant of non-square matrix"+        | otherwise = detMat (matForm m)+  inv m | col m /= row m = error "Can't calculate inverse of non-square matrix"+        | otherwise = (matrix . invMat . matForm) m+  transpose m = m {row = col m, col = row m, byRow = not (byRow m)}++-- |Block Partitioning+bp :: Int -> Matrix a -> Matrix a+bp n m = matrix $ bpMat n (matForm m)++---------------------------------------------------+-- Concatenate+---------------------------------------------------+class Functor f => Concatable f where+  hcat :: f a -> f a -> f a+  vcat :: f a -> f a -> Matrix a++instance Concatable Vector where+  hcat v w = fromList (toList v ++ toList w)+  vcat v w = matrix (toList v : [toList w])++instance Concatable Matrix where+  hcat m n | row m == row n = matrix (zipWith (++) mf nf)+           | otherwise = error "Can't concatenate matrices horizontally which have different row"+           where mf = matForm m+                 nf = matForm n+  vcat m n | col m == col n = m {val = hcat (val m) (val n), row = row m + row n}+           | otherwise = error "Can't concatenate matrices vertically which have different col"++-- |(.:) inserts vector to head of matrix.+(.:) :: Vector a -> Matrix a -> Matrix a+v .: m | length v == col m = matrix (toList v : matForm m)+       | otherwise         = error "Can't insert length(Vector) /= col(Matrix)"++---------------------------------------------------+-- Backend Functions (Do not Understand)+---------------------------------------------------+-- Transpose+transposeMat :: [[a]] -> [[a]]+transposeMat m = map (\l -> map (!! l) m) [0 .. (length (head m) - 1)]++-- indexMat+indexMat :: [[a]] -> [[(Int, Int)]]+indexMat m@(xs : xss) = do+  i <- [0 .. (length m - 1)]+  [zip (replicate (length xs) i) [0 .. (length xs - 1)]]++-- dropAtMat+dropAtMat :: Int -> Int -> [[a]] -> [[a]]+dropAtMat i j mat = map (dropAt j) $ dropAt i mat++-- postSplitAt+postSplitAt (x, y) = x ++ tail y++-- dropAt+dropAt :: Int -> [a] -> [a]+dropAt i = postSplitAt . splitAt i++-- dropAtMat' : drop nth array+dropAtMat' :: Int -> [[a]] -> [[a]]+dropAtMat' n mat | n /= (length mat - 1) = dropAt n mat+                 | otherwise             = take n mat++-- Block Partitioning+bpMat :: Int -> [[a]] -> [[a]]+bpMat _ [] = []+bpMat n m | n == 1    = (map (take sl) . take sl) m+          | n == 2    = (map (drop sl) . take sl) m+          | n == 3    = (map (take sl) . drop sl) m+          | n == 4    = (map (drop sl) . drop sl) m+          | otherwise = error "Please input 1 ~ 4"+ where+  l  = length m+  sl = (floor . sqrt . fromIntegral) l++-- Matrix + Matrix+(%-+-%) :: Num a => [[a]] -> [[a]] -> [[a]]+m    %-+-% []   = m+[]   %-+-% m    = m+[[]] %-+-% m    = m+m    %-+-% [[]] = m+m    %-+-% n    = zipWith (zipWith (+)) m n++negMap :: Num a => [[a]] -> [[a]]+negMap = map (map negate)++-- Matrix - Matrix+(%---%) :: Num a => [[a]] -> [[a]] -> [[a]]+m    %---% []   = m+[]   %---% m    = map (map negate) m+[[]] %---% m    = map (map negate) m+m    %---% [[]] = m+m    %---% n    = zipWith (zipWith (-)) m n++-- Matrix Multiplication+(%-*-%) :: Num a => [[a]] -> [[a]] -> [[a]]+_     %-*-% []    = []+[]    %-*-% _     = []+_     %-*-% [[]]  = [[]]+[[] ] %-*-% _     = [[]]+[[x]] %-*-% [[y]] = [[x * y]]+m     %-*-% n     = zipWith (++) a11 a12 ++ zipWith (++) a21 a22+ where+  (m11, n11) = (bpMat 1 m, bpMat 1 n)+  (m12, n12) = (bpMat 2 m, bpMat 2 n)+  (m21, n21) = (bpMat 3 m, bpMat 3 n)+  (m22, n22) = (bpMat 4 m, bpMat 4 n)+  a11        = (m11 %-*-% n11) %-+-% (m12 %-*-% n21)+  a12        = (m11 %-*-% n12) %-+-% (m12 %-*-% n22)+  a21        = (m21 %-*-% n11) %-+-% (m22 %-*-% n21)+  a22        = (m21 %-*-% n12) %-+-% (m22 %-*-% n22)++zerosVec :: Int -> [Int]+zerosVec n = take n [0, 0 ..]++eyeMat :: Int -> [[Int]]+eyeMat n = [ basisVec x n | x <- [0 .. (n - 1)] ]++-- Position -> Length +basisVec :: Int -> Int -> [Int]+basisVec n m = zerosVec n ++ [1] ++ zerosVec (m - n - 1)++permMat :: Int -> Int -> [[a]] -> [[Int]]+permMat i j m+  | i < j+  = take i idx+    ++ [idx !! j]+    ++ take (j - i - 1) (drop (i + 1) idx)+    ++ [idx !! i]+    ++ drop (j + 1) idx+  | otherwise+  = permMat j i m+  where idx = eyeMat (length m)++whichMax :: Ord a => [a] -> Int+whichMax v = whichMax' v 0 m+ where+  m = maximum v+  whichMax' :: Ord a => [a] -> Int -> a -> Int+  whichMax' (x : xs) n m' = if x == m' then n else whichMax' xs (n + 1) m'++colMat :: [[a]] -> Int -> [a]+colMat m n = map (!! n) m++colMaxIdx :: Ord a => [[a]] -> Int -> Int+colMaxIdx m n = whichMax $ colMat m n++-- | Another Block Partitioning+bpMat' :: Int -> [[a]] -> [[a]]+bpMat' _ []  = []+bpMat' _ [x] = [x]+bpMat' n m | n == 1 = (map (take l) . take l) m+           | n == 2 = (map (drop 1) . take l) m+           | n == 3 = (map (take l) . drop 1) m+           | n == 4 = (map (drop 1) . drop 1) m+           | n == 0 = (map (drop 1 . take l) . drop 1 . take l) m+  where l = length m - 1++-- | Determinant for Double List - Order ~ 4^n+detMat :: Fractional a => [[a]] -> a+detMat [[x]] = x+detMat m+  | l == 2    = detMat m11 * detMat m22 - detMat m12 * detMat m21+  | otherwise = (detMat m11 * detMat m22 - detMat m12 * detMat m21) / detMat m00+ where+  l   = length m+  m11 = bpMat' 1 m+  m12 = bpMat' 2 m+  m21 = bpMat' 3 m+  m22 = bpMat' 4 m+  m00 = bpMat' 0 m++-- | Inverse for Double List - Order ~ n * 2^n+invMat :: Fractional a => [[a]] -> [[a]]+invMat []    = []+invMat [[] ] = [[]]+invMat [[x]] = [[x]]+invMat m+  | length m == 2+  = map (map (/ detMat m))+    $  zipWith (++) m22          (negMap m12)+    ++ zipWith (++) (negMap m21) m11+  | otherwise+  = zipWith (++) a11 a12 ++ zipWith (++) a21 a22+ where+  m11 = bpMat 1 m+  m12 = bpMat 2 m+  m21 = bpMat 3 m+  m22 = bpMat 4 m+  a00 = invMat m11+  s   = m22 %---% (m21 %-*-% a00 %-*-% m12)+  s00 = invMat s+  a11 = a00 %-+-% (a00 %-*-% m12 %-*-% s00 %-*-% m21 %-*-% a00)+  a12 = negMap a00 %-*-% m12 %-*-% s00+  a21 = negMap s00 %-*-% m21 %-*-% a00+  a22 = s00

@@ -0,0 +1,2 @@+main :: IO ()+main = putStrLn "Test suite not yet implemented"