packages feed

qlinear (empty) → 0.1.0.0

raw patch · 21 files changed

+1243/−0 lines, 21 filesdep +basedep +haskell-src-extsdep +haskell-src-meta

Dependencies added: base, haskell-src-exts, haskell-src-meta, hspec, linear, parsec, split, template-haskell

Files

+ LICENSE view
@@ -0,0 +1,30 @@+Copyright Author name here (c) 2020++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.
+ README.md view
@@ -0,0 +1,37 @@+# QLinear++QLinear is type safe library for linear algebra based on `"macro-constructors"`++## Constructors:+ * matrix+    > ```[matrix| a b ; c d |]``` builds matrix __2x2__  `[ [a, b], [c, d] ]`  +    +    Example: +    ```haskell+    [matrix| 1 2; 3 4 |] ~+~ [matrix| 2 3; 4 5 |] == [matrix| 3 5; 7 9 |] +    ```   +    > Also you can't, for example, add two matrix with different size. +    ```haskell+   [matrix| 1 2; 3 4 |] ~+~ [matrix| 1 2 3; 4 5 6; 7 8 9 |] -- will not be compiled+    ```++  * vector+    > ```[vector| a b c d |]``` builds matrix __4x1__ `[ [a], [b], [c], [d] ]`  + * operator+    > ```[operator| (x, y) => (y, x) |]``` builds matrix __2x2__ of operator `[ [0, 1], [1, 0] ]` that swaps coodrinates  +    > ```[operator| (x, y) => (2 * x, y) |]``` builds matrix __2x2__ of operator `[ [2, 0], [0, 1] ]` that doubles `x` coordinate  +   +   Example: +    ```haskell+    [operator| (x, y) => (3 * y, x / 2) |] ~*~ [vector| 2 8 |] == [vector| 24 1 |]+    ```+## Syntax:+   * matrix: `val11 val12 .. val1n; val21 val22 .. val2n; ..; valm1 valm2 .. valmn `+   * vector: `val1 val2 .. valn`+   * operator: `(var1, var2, .., varn) =>  (exp1, exp2, .., expn)`  +     where +     * `val` is `number literal`, `variable` or `any Haskell expression` between `(` and `)`+     * `var` is Haskell variable   +     * `exp` is any Haskell expression+  ---+Also there are basic operations as determinant, transposition, etc.
+ qlinear.cabal view
@@ -0,0 +1,82 @@+cabal-version: 3.0+name:           qlinear+version:        0.1.0.0+synopsis:       Typesafe library for linear algebra+description:    Please see the README on GitHub at <https://github.com/JuniorGarbageCollector/QLinear>+category:       Math+homepage:       https://github.com/JuniorGarbageCollector/QLinear#readme+bug-reports:    https://github.com/JuniorGarbageCollector/QLinear/issues+author:         JuniorGarbageCollector+maintainer:     GooseDB@yandex.ru+license:        BSD-3-Clause+license-file:   LICENSE+build-type:     Simple+extra-source-files:+    README.md++source-repository head+  type: git+  location: https://github.com/JuniorGarbageCollector/QLinear++library+  exposed-modules:+      QLinear+      QLinear.Constructor.Matrix+      QLinear.Constructor.Operator+      QLinear.Identity+      QLinear.Index+      QLinear.Integration.Linear.From+      QLinear.Integration.Linear.To+      QLinear.Matrix+      QLinear.Operations+  other-modules:+      Internal.Determinant+      Internal.Matrix+      Internal.Quasi.Matrix.Parser+      Internal.Quasi.Matrix.Quasi+      Internal.Quasi.Operator.Parser+      Internal.Quasi.Operator.Quasi+      Internal.Quasi.Parser+      Internal.Quasi.Quasi+  hs-source-dirs:+      src+  build-depends:+      base >=4.7 && <5+    , haskell-src-exts >=1.23.1 && <1.24+    , haskell-src-meta >=0.8.5 && <0.9+    , linear >=1.21.1 && <1.22+    , parsec >=3.1.13 && <3.2+    , split >=0.2.3 && <0.3+    , template-haskell >= 2.14.0 && < 2.15+  default-language: Haskell2010++test-suite UnitTests+  type:                exitcode-stdio-1.0+  main-is:             Main.hs+  hs-source-dirs:      test, src+  build-depends:       base,+                       hspec,+                       split,+                       template-haskell,+                       parsec,+                       haskell-src-exts,+                       haskell-src-meta,+                       linear+  other-modules:       QLinear+                       QLinear.Constructor.Matrix+                       QLinear.Constructor.Operator+                       QLinear.Identity+                       QLinear.Index+                       QLinear.Integration.Linear.From+                       QLinear.Integration.Linear.To+                       QLinear.Matrix+                       QLinear.Operations+                       Internal.Determinant+                       Internal.Matrix+                       Internal.Quasi.Matrix.Parser+                       Internal.Quasi.Matrix.Quasi+                       Internal.Quasi.Operator.Parser+                       Internal.Quasi.Operator.Quasi+                       Internal.Quasi.Parser+                       Internal.Quasi.Quasi+  default-language:    Haskell2010
+ src/Internal/Determinant.hs view
@@ -0,0 +1,98 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}++module Internal.Determinant where++import qualified Data.List as List+import Data.List.Split+import Data.Proxy+import qualified GHC.Natural as Natural+import GHC.TypeNats+import Internal.Matrix+import QLinear.Index++-- | Determinant of matrix+--+-- >>> det [matrix| 1 0; 0 1|]+-- 1+-- >>> det [matrix| 1 3; 4 2|]+-- -10+det :: Num a => Matrix n n a -> a+det (Matrix (1, _) [[a]]) = a+det (Matrix (2, _) [[a, b], [c, d]]) = a * d - b * c+det (Matrix (3, _) [[a, b, c], [d, e, f], [g, h, i]]) =+  a * (e * i - f * h) - d * (b * i - c * h) + g * (b * f - c * e)+det (Matrix (n, _) matrix) = sum $ map calcElem indices+  where+    calcElem index@(i, j) = ((!! (j - 1)) $ matrix !! (i - 1)) * unsafeAlgebraicComplement matrix n index+    indices = zip (repeat 1) [1 .. n]++-- | Typesafe algebraic complement+--+-- To use it you have to know i and j at compile time+--+-- >>> algebraicComplement [matrix| 1 2; 3 4 |] (Index @1 @1)+-- 4+-- >>> algebraicComplement [matrix| 1 2 3; 4 5 6; 7 8 9 |] (Index @1 @1)+-- -3+algebraicComplement ::+  forall n a i j.+  (KnownNat i, KnownNat j, KnownNat n, Num a, i <= n, j <= n) =>+  Matrix n n a ->+  Index i j ->+  a+algebraicComplement (Matrix (n, _) matrix) _ = unsafeAlgebraicComplement matrix n (i, j)+  where+    i = (Natural.naturalToInt $ natVal $ Proxy @i)+    j = (Natural.naturalToInt $ natVal $ Proxy @j)++-- | Algebraic complement.+--+-- Use it if you don't know indices at compile time+--+-- >>> algebraicComplement' [matrix| 1 2; 3 4 |] (1, 1)+-- Just 4+--+-- >>> algebraicComplement' [matrix| 1 2; 3 4 |] (34, 43)+-- Nothing+--+-- >>> algebraicComplement' [matrix| 1 2 3; 4 5 6; 7 8 9 |] (1, 1)+-- Just (-3)+algebraicComplement' :: Num a => Matrix n n a -> (Int, Int) -> Maybe a+algebraicComplement' (Matrix (n, _) matrix) ij@(i, j)+  | i <= n && j <= n = Just $ unsafeAlgebraicComplement matrix n ij+  | otherwise = Nothing++-- | Adjugate matrix+--+-- >>> adjugate [matrix| 1 2; 3 4|]+-- [4,-2]+-- [-3,1]+adjugate :: Num a => Matrix n n a -> Matrix n n a+adjugate (Matrix size@(n, _) matrix) = Matrix size $ chunksOf n $ adj+  where+    adj = map (unsafeAlgebraicComplement matrix n) [(i, j) | j <- [1 .. n], i <- [1 .. n]]++unsafeAlgebraicComplement :: forall a. Num a => [[a]] -> Int -> (Int, Int) -> a+unsafeAlgebraicComplement matrix n (i, j) = k * det (Matrix (n - 1, n - 1) minor)+  where+    minor = cross matrix (i - 1, j - 1)+    k = fromIntegral $ (-1) ^ (j + i)++cross :: [[a]] -> (Int, Int) -> [[a]]+cross matrix (i, j) = withoutColumn+  where+    withoutLine = cut i matrix+    withoutColumn = map (cut j) withoutLine++cut' :: Int -> [a] -> (a, [a])+cut' 0 (x : xs) = (x, xs)+cut' n (x : xs) = (x :) <$> cut' (n - 1) xs++cut :: Int -> [a] -> [a]+cut 0 (x : xs) = xs+cut n (x : xs) = x : cut (n - 1) xs
+ src/Internal/Matrix.hs view
@@ -0,0 +1,39 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE GADTs #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE RankNTypes #-}++module Internal.Matrix where++import qualified Data.List as List+import Data.List.Split+import Data.Proxy+import Data.Tuple+import qualified GHC.Natural as Natural+import GHC.TypeNats+import Prelude hiding (length)++data Matrix (m :: Nat) (n :: Nat) a where+  Matrix :: forall m n a. (Int, Int) -> ![[a]] -> Matrix m n a++instance Show a => Show (Matrix m n a) where+  show (Matrix _ matrix) = List.intercalate "\n" $ map show matrix++instance Functor (Matrix m n) where+  fmap f (Matrix size a) = Matrix size $ map (map f) a++instance Applicative (Matrix m n) where+  pure = Matrix (1, 1) . pure . pure+  Matrix size fs <*> (Matrix _ as) =+    Matrix size $ map (uncurry (<*>)) $ zip fs as++instance (Eq a) => Eq (Matrix m n a) where+  Matrix _ a == Matrix _ b = a == b++type Vector n a = Matrix n 1 a++size :: Integral b => Matrix m n a -> (b, b)+size (Matrix (m, n) _) = (fromIntegral m, fromIntegral n)++value :: Matrix m n a -> [[a]]+value (Matrix _ value) = value
+ src/Internal/Quasi/Matrix/Parser.hs view
@@ -0,0 +1,37 @@+module Internal.Quasi.Matrix.Parser where++import Internal.Quasi.Parser+import Language.Haskell.TH.Syntax++matrix :: Parser [[Exp]]+matrix = (line `sepBy` char ';') <* eof++vector :: Parser [[Exp]]+vector = map pure <$> (line <* eof)++line :: Parser [Exp]+line = spaces >> unit `endBy1` spaces++unit :: Parser Exp+unit = (var <|> num <|> inBrackets) >>= expr++num :: Parser String+num = do+  neg <- (char' '-') <|> pure []+  beforeDot <- ((many1 outer <> many inner) <|> char' '0')+  afterDot <- char' '.' <> many1 inner <|> mempty+  pure $ neg <> beforeDot <> afterDot+  where+    outer = oneOf ['1' .. '9']+    inner = char '0' <|> outer++inBrackets :: Parser String+inBrackets = nested '(' ')'++nested :: Char -> Char -> Parser String+nested open close = char' open <> scan 1+  where+    scan 0 = pure mempty+    scan n = many (noneOf [open, close]) <> (char' open <> inc n <|> char' close <> dec n)+    inc = scan . (+ 1)+    dec = scan . (\n -> n - 1)
+ src/Internal/Quasi/Matrix/Quasi.hs view
@@ -0,0 +1,81 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE TypeApplications #-}++module Internal.Quasi.Matrix.Quasi (matrix, vector) where++import Internal.Matrix+import qualified Internal.Quasi.Matrix.Parser as Parser+import qualified Internal.Quasi.Parser as Parser+import Internal.Quasi.Quasi+import Language.Haskell.TH.Quote+import Language.Haskell.TH.Syntax++-- | Macro constructor for 'QLinear.Matrix.Matrix'+--+-- >>> [matrix| 1 2; 3 4 |]+-- [1,2]+-- [3,4]+-- >>> :t [matrix| 1 2; 3 4|]+-- [matrix| 1 2; 3 4|] :: Num a => Matrix 2 2 a+matrix :: QuasiQuoter+matrix =+  QuasiQuoter+    { quoteExp = expr Parser.matrix,+      quotePat = notDefined "Pattern",+      quoteType = notDefined "Type",+      quoteDec = notDefined "Declaration"+    }+  where+    notDefined = isNotDefinedAs "matrix"++-- | Macro constructor for 'QLinear.Matrix.Vector'.+--+-- >>> [vector| 1 2 3 4 |]+-- [1]+-- [2]+-- [3]+-- [4]+-- >>> :t [vector| 1 2 3 4 |]+-- [vector| 1 2 3 4 |] :: Num a => Vector 4 a+vector :: QuasiQuoter+vector =+  QuasiQuoter+    { quoteExp = vectorExpr,+      quotePat = notDefined "Pattern",+      quoteType = notDefined "Type",+      quoteDec = notDefined "Declaration"+    }+  where+    notDefined = isNotDefinedAs "vector"++vectorExpr :: String -> Q Exp+vectorExpr source = f <$> expr Parser.vector source+  where+    f = AppE (VarE 'toVector)++expr :: Parser.Parser [[Exp]] -> String -> Q Exp+expr parser source = do+  let (matrix, (m, n)) = unwrap $ parse source parser+  let sizeType = LitT . NumTyLit+  let constructor = foldl AppTypeE (ConE 'Matrix) [sizeType m, sizeType n, WildCardT]+  let size = TupE $ map (LitE . IntegerL) [m, n]+  let value = ListE $ map ListE $ matrix+  pure $ foldl AppE constructor [size, value]++parse :: String -> Parser.Parser [[a]] -> Either [String] ([[a]], (Integer, Integer))+parse source parser = do+  matrix <- Parser.parse parser "QLinear" source+  size <- checkSize matrix+  pure (matrix, size)++checkSize :: [[a]] -> Either [String] (Integer, Integer)+checkSize [] = Left ["Matrix cannot be empty"]+checkSize matrix =+  let lines@(l : ls) = map length matrix+   in if all (== l) ls+        then Right (fromIntegral $ length matrix, fromIntegral l)+        else Left ["All lines must be the same length"]++toVector :: Matrix n 1 a -> Vector n a+toVector = id
+ src/Internal/Quasi/Operator/Parser.hs view
@@ -0,0 +1,24 @@+module Internal.Quasi.Operator.Parser where++import Internal.Quasi.Parser+import Language.Haskell.TH.Syntax++definition :: Parser ([Pat], [Exp])+definition = do+  params <- spaces *> parameters <* spaces+  string "=>"+  lams <- spaces *> lambdas <* spaces+  eof+  pure (params, lams)++lambdas :: Parser [Exp]+lambdas = do+  result <- many1 anyChar >>= expr+  case result of+    TupE elems -> pure elems+    otherwise -> parserFail $ show otherwise++parameters :: Parser [Pat]+parameters = char '(' *> inner <* char ')'+  where+    inner = map (VarP . mkName) <$> var `sepBy` (spaces >> char ',' >> spaces)
+ src/Internal/Quasi/Operator/Quasi.hs view
@@ -0,0 +1,77 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE MonoLocalBinds #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE TypeApplications #-}++module Internal.Quasi.Operator.Quasi where++import Data.List.Split (chunksOf)+import Data.Proxy+import qualified GHC.Natural as Natural+import GHC.TypeNats+import Internal.Matrix+import qualified Internal.Quasi.Operator.Parser as Parser+import qualified Internal.Quasi.Parser as Parser+import Internal.Quasi.Quasi+import Language.Haskell.TH.Quote+import Language.Haskell.TH.Syntax+import QLinear.Identity++{- | Macro constructor for operator++>>> [operator| (x, y) => (y, x) |]+[0,1]+[1,0]+>>> [operator| (x, y) => (2 * x, y + x) |] ~*~ [vector| 3 4 |]+[6]+[7]++Do note,constructor __doesn't prove__ linearity.+It just builds matrix of given operator.++-}+operator :: QuasiQuoter+operator =+  QuasiQuoter+    { quoteExp = expr,+      quotePat = notDefined "Pattern",+      quoteType = notDefined "Type",+      quoteDec = notDefined "Declaration"+    }+  where+    notDefined = isNotDefinedAs "operator"++expr :: String -> Q Exp+expr source = do+  let (params, lams, n) = unwrap $ parse source+  let sizeType = LitT . NumTyLit+  let size = TupE $ map (LitE . IntegerL) [n, 1]+  let func = VarE 'matrixOfOperator+  let constructor = foldl AppTypeE (ConE 'Matrix) [sizeType n, sizeType 1, WildCardT]+  let value = ListE $ map (ListE . pure . LamE [ListP params]) lams+  pure $ AppE func $ foldl AppE constructor [size, value]+  where++parse :: String -> Either [String] ([Pat], [Exp], Integer)+parse source = do+  (params, lams) <- Parser.parse Parser.definition "QLinear" source+  size <- checkSize (params, lams)+  pure (params, lams, size)++checkSize :: ([Pat], [Exp]) -> Either [String] Integer+checkSize ([], _) = Left ["Parameters of operator cannot be empty"]+checkSize (_, []) = Left ["Body of operator cannot be empty"]+checkSize (names, exprs) =+  let namesLength = length names+      exprsLength = length exprs+   in if namesLength == exprsLength+        then Right $ fromIntegral namesLength+        else Left ["Number of arguments and number of lambdas must be equal"]++matrixOfOperator :: forall n a b. (KnownNat n, HasIdentity a) => Matrix n 1 ([a] -> b) -> Matrix n n b+matrixOfOperator (Matrix _ fs) = Matrix (n, n) $ chunksOf n [f line | f <- concat fs, line <- identity]+  where+    (Matrix _ identity) = e :: Matrix n n a+    n = Natural.naturalToInt $ natVal (Proxy @n)
+ src/Internal/Quasi/Parser.hs view
@@ -0,0 +1,40 @@+module Internal.Quasi.Parser (module Parsec, module Char, Parser, parse, var, expr, satisfyOneOf, char') where++import Control.Monad+import Data.Char as Char+import Language.Haskell.Exts.Parser (ParseResult (..), parseExp)+import Language.Haskell.Meta.Syntax.Translate (toExp)+import Language.Haskell.TH.Syntax+import Text.Parsec as Parsec hiding (parse)+import qualified Text.Parsec as P+import Text.Parsec.Char as Parsec+import Text.Parsec.Combinator as Parsec+import Text.Parsec.Error as Parsec++type Parser a = Parsec String () a++parse :: Parser a -> SourceName -> String -> Either [String] a+parse parser name source = case P.parse parser name source of+  Right a -> Right a+  Left err -> Left $ map messageString $ errorMessages $ err++satisfyOneOf :: [Char -> Bool] -> Parser Char+satisfyOneOf ps = satisfy (or . ap ps . pure)++expr :: String -> Parser Exp+expr source =+  case parseExp source of+    ParseOk exp -> pure $ toExp exp+    ParseFailed _ e -> parserFail e++char' :: Char -> Parser String+char' = fmap pure . char++anyChar' :: Parser String+anyChar' = fmap pure $ anyChar++var :: Parser String+var = ((many1 $ satisfyOneOf outer) <> (many $ satisfyOneOf inner))+  where+    outer = [isAlpha, (== '_')]+    inner = isDigit : (== '\'') : outer
+ src/Internal/Quasi/Quasi.hs view
@@ -0,0 +1,10 @@+module Internal.Quasi.Quasi where++import qualified Data.List as List++isNotDefinedAs :: String -> String -> a+isNotDefinedAs name as = error $ "You cannot use " <> name <> " quasi as " <> as++unwrap :: Either [String] a -> a+unwrap (Left a) = error $ List.intercalate ", " a+unwrap (Right a ) = a
+ src/QLinear.hs view
@@ -0,0 +1,11 @@+module QLinear+  ( module E,+  )+where++import QLinear.Constructor.Matrix as E+import QLinear.Constructor.Operator as E+import QLinear.Identity as E+import QLinear.Index as E+import QLinear.Matrix as E+import QLinear.Operations as E
+ src/QLinear/Constructor/Matrix.hs view
@@ -0,0 +1,3 @@+module QLinear.Constructor.Matrix (matrix, vector) where++import Internal.Quasi.Matrix.Quasi
+ src/QLinear/Constructor/Operator.hs view
@@ -0,0 +1,3 @@+module QLinear.Constructor.Operator (operator) where++import Internal.Quasi.Operator.Quasi
+ src/QLinear/Identity.hs view
@@ -0,0 +1,48 @@+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MonoLocalBinds #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE UndecidableInstances #-}++module QLinear.Identity (e, Identity, HasIdentity (..)) where++import Data.Proxy+import qualified GHC.Natural as Natural+import GHC.TypeNats+import Internal.Matrix++type Identity n a = Matrix n n a++class HasIdentity a where+  zero :: a+  one :: a++instance (Num a) => HasIdentity a where+  zero = 0+  one = 1++-- | Polymirphic identity matrix+--+-- Identity matrix can udjust to other matrix with known size. If size is unknown, just set it yourself+--+-- >>> e :: Identity 4 Int+-- [1,0,0,0]+-- [0,1,0,0]+-- [0,0,1,0]+-- [0,0,0,1]+-- >>> e ~+~ [matrix| 1 2; 3 4 |]+-- [2,2]+-- [3,5]+e :: forall n a. (KnownNat n, HasIdentity a) => Identity n a+e = Matrix (n, n) $ finiteIdentityList (n, n) one zero+  where+    n = Natural.naturalToInt $ natVal (Proxy @n)++infiniteIdentityList :: a -> a -> [[a]]+infiniteIdentityList o z = stream (o : repeat z)+  where+    stream seed = seed : stream (z : seed)++finiteIdentityList :: (Int, Int) -> a -> a -> [[a]]+finiteIdentityList (m, n) o z = map (take n) $ take m $ infiniteIdentityList o z
+ src/QLinear/Index.hs view
@@ -0,0 +1,25 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeApplications #-}++module QLinear.Index where++import Data.Proxy+import qualified GHC.Natural as Natural+import GHC.TypeNats++-- | Typesafe index. To construct it use TypeApplications+--+-- >>> Index @1 @3+-- Index 1 3+-- >>> :t Index @1 @3+-- Index @1 @3 :: Index 1 3+data Index (i :: Nat) (j :: Nat) = Index++instance forall i j. (KnownNat i, KnownNat j) => Show (Index i j) where+  show _ = "Index " <> show i <> " " <> show j+    where+      i = natVal (Proxy @i)+      j = natVal (Proxy @j)
+ src/QLinear/Integration/Linear/From.hs view
@@ -0,0 +1,155 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE QuasiQuotes #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}++module QLinear.Integration.Linear.From where++import Data.List+import Internal.Matrix+import Internal.Quasi.Matrix.Quasi+import Linear (V1 (..), V2 (..), V3 (..), V4 (..))++class FromLinear m where+  type Q m :: *+  fromLinear :: m -> Q m++{- GENERATED CODE -}++instance FromLinear (V1 (V1 a)) where+  type Q (V1 (V1 a)) = Matrix 1 1 a+  fromLinear (V1 (V1 a1)) =+    [matrix| +     a1+    |]++instance FromLinear (V1 (V2 a)) where+  type Q (V1 (V2 a)) = Matrix 1 2 a+  fromLinear (V1 (V2 a2 a3)) =+    [matrix| +     a2  a3+    |]++instance FromLinear (V1 (V3 a)) where+  type Q (V1 (V3 a)) = Matrix 1 3 a+  fromLinear (V1 (V3 a3 a4 a5)) =+    [matrix| +     a3  a4  a5+    |]++instance FromLinear (V1 (V4 a)) where+  type Q (V1 (V4 a)) = Matrix 1 4 a+  fromLinear (V1 (V4 a4 a5 a6 a7)) =+    [matrix| +     a4  a5  a6  a7+    |]++instance FromLinear (V2 (V1 a)) where+  type Q (V2 (V1 a)) = Matrix 2 1 a+  fromLinear (V2 (V1 a1) (V1 a2)) =+    [matrix| +     a1;+     a2+    |]++instance FromLinear (V2 (V2 a)) where+  type Q (V2 (V2 a)) = Matrix 2 2 a+  fromLinear (V2 (V2 a2 a3) (V2 a4 a5)) =+    [matrix| +     a2  a3;+     a4  a5+    |]++instance FromLinear (V2 (V3 a)) where+  type Q (V2 (V3 a)) = Matrix 2 3 a+  fromLinear (V2 (V3 a3 a4 a5) (V3 a6 a7 a8)) =+    [matrix| +     a3  a4  a5;+     a6  a7  a8+    |]++instance FromLinear (V2 (V4 a)) where+  type Q (V2 (V4 a)) = Matrix 2 4 a+  fromLinear (V2 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11)) =+    [matrix| +     a4  a5  a6  a7;+     a8  a9  a10  a11+    |]++instance FromLinear (V3 (V1 a)) where+  type Q (V3 (V1 a)) = Matrix 3 1 a+  fromLinear (V3 (V1 a1) (V1 a2) (V1 a3)) =+    [matrix| +     a1;+     a2;+     a3+    |]++instance FromLinear (V3 (V2 a)) where+  type Q (V3 (V2 a)) = Matrix 3 2 a+  fromLinear (V3 (V2 a2 a3) (V2 a4 a5) (V2 a6 a7)) =+    [matrix| +     a2  a3;+     a4  a5;+     a6  a7+    |]++instance FromLinear (V3 (V3 a)) where+  type Q (V3 (V3 a)) = Matrix 3 3 a+  fromLinear (V3 (V3 a3 a4 a5) (V3 a6 a7 a8) (V3 a9 a10 a11)) =+    [matrix| +     a3  a4  a5;+     a6  a7  a8;+     a9  a10  a11+    |]++instance FromLinear (V3 (V4 a)) where+  type Q (V3 (V4 a)) = Matrix 3 4 a+  fromLinear (V3 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11) (V4 a12 a13 a14 a15)) =+    [matrix| +     a4  a5  a6  a7;+     a8  a9  a10  a11;+     a12  a13  a14  a15+    |]++instance FromLinear (V4 (V1 a)) where+  type Q (V4 (V1 a)) = Matrix 4 1 a+  fromLinear (V4 (V1 a1) (V1 a2) (V1 a3) (V1 a4)) =+    [matrix| +     a1;+     a2;+     a3;+     a4+    |]++instance FromLinear (V4 (V2 a)) where+  type Q (V4 (V2 a)) = Matrix 4 2 a+  fromLinear (V4 (V2 a2 a3) (V2 a4 a5) (V2 a6 a7) (V2 a8 a9)) =+    [matrix| +     a2  a3;+     a4  a5;+     a6  a7;+     a8  a9+    |]++instance FromLinear (V4 (V3 a)) where+  type Q (V4 (V3 a)) = Matrix 4 3 a+  fromLinear (V4 (V3 a3 a4 a5) (V3 a6 a7 a8) (V3 a9 a10 a11) (V3 a12 a13 a14)) =+    [matrix| +     a3  a4  a5;+     a6  a7  a8;+     a9  a10  a11;+     a12  a13  a14+    |]++instance FromLinear (V4 (V4 a)) where+  type Q (V4 (V4 a)) = Matrix 4 4 a+  fromLinear (V4 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11) (V4 a12 a13 a14 a15) (V4 a16 a17 a18 a19)) =+    [matrix| +     a4  a5  a6  a7;+     a8  a9  a10  a11;+     a12  a13  a14  a15;+     a16  a17  a18  a19+    |]
+ src/QLinear/Integration/Linear/To.hs view
@@ -0,0 +1,87 @@+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE QuasiQuotes #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE TypeApplications #-}+{-# LANGUAGE TypeFamilies #-}++module QLinear.Integration.Linear.To where++import Data.List+import Internal.Matrix+import Linear (V1 (..), V2 (..), V3 (..), V4 (..))++class ToLinear m where+  type L m :: *+  toLinear :: m -> L m++{- GENERATED CODE -}++instance ToLinear (Matrix 1 1 a) where+  type L (Matrix 1 1 a) = V1 (V1 a)+  toLinear (Matrix _ [[a1]]) = (V1 (V1 a1))++instance ToLinear (Matrix 1 2 a) where+  type L (Matrix 1 2 a) = V1 (V2 a)+  toLinear (Matrix _ [[a2, a3]]) = (V1 (V2 a2 a3))++instance ToLinear (Matrix 1 3 a) where+  type L (Matrix 1 3 a) = V1 (V3 a)+  toLinear (Matrix _ [[a3, a4, a5]]) = (V1 (V3 a3 a4 a5))++instance ToLinear (Matrix 1 4 a) where+  type L (Matrix 1 4 a) = V1 (V4 a)+  toLinear (Matrix _ [[a4, a5, a6, a7]]) = (V1 (V4 a4 a5 a6 a7))++instance ToLinear (Matrix 2 1 a) where+  type L (Matrix 2 1 a) = V2 (V1 a)+  toLinear (Matrix _ [[a1], [a2]]) = (V2 (V1 a1) (V1 a2))++instance ToLinear (Matrix 2 2 a) where+  type L (Matrix 2 2 a) = V2 (V2 a)+  toLinear (Matrix _ [[a2, a3], [a4, a5]]) = (V2 (V2 a2 a3) (V2 a4 a5))++instance ToLinear (Matrix 2 3 a) where+  type L (Matrix 2 3 a) = V2 (V3 a)+  toLinear (Matrix _ [[a3, a4, a5], [a6, a7, a8]]) = (V2 (V3 a3 a4 a5) (V3 a6 a7 a8))++instance ToLinear (Matrix 2 4 a) where+  type L (Matrix 2 4 a) = V2 (V4 a)+  toLinear (Matrix _ [[a4, a5, a6, a7], [a8, a9, a10, a11]]) = (V2 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11))++instance ToLinear (Matrix 3 1 a) where+  type L (Matrix 3 1 a) = V3 (V1 a)+  toLinear (Matrix _ [[a1], [a2], [a3]]) = (V3 (V1 a1) (V1 a2) (V1 a3))++instance ToLinear (Matrix 3 2 a) where+  type L (Matrix 3 2 a) = V3 (V2 a)+  toLinear (Matrix _ [[a2, a3], [a4, a5], [a6, a7]]) = (V3 (V2 a2 a3) (V2 a4 a5) (V2 a6 a7))++instance ToLinear (Matrix 3 3 a) where+  type L (Matrix 3 3 a) = V3 (V3 a)+  toLinear (Matrix _ [[a3, a4, a5], [a6, a7, a8], [a9, a10, a11]]) =+    (V3 (V3 a3 a4 a5) (V3 a6 a7 a8) (V3 a9 a10 a11))++instance ToLinear (Matrix 3 4 a) where+  type L (Matrix 3 4 a) = V3 (V4 a)+  toLinear (Matrix _ [[a4, a5, a6, a7], [a8, a9, a10, a11], [a12, a13, a14, a15]]) =+    (V3 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11) (V4 a12 a13 a14 a15))++instance ToLinear (Matrix 4 1 a) where+  type L (Matrix 4 1 a) = V4 (V1 a)+  toLinear (Matrix _ [[a1], [a2], [a3], [a4]]) = (V4 (V1 a1) (V1 a2) (V1 a3) (V1 a4))++instance ToLinear (Matrix 4 2 a) where+  type L (Matrix 4 2 a) = V4 (V2 a)+  toLinear (Matrix _ [[a2, a3], [a4, a5], [a6, a7], [a8, a9]]) =+    (V4 (V2 a2 a3) (V2 a4 a5) (V2 a6 a7) (V2 a8 a9))++instance ToLinear (Matrix 4 3 a) where+  type L (Matrix 4 3 a) = V4 (V3 a)+  toLinear (Matrix _ [[a3, a4, a5], [a6, a7, a8], [a9, a10, a11], [a12, a13, a14]]) =+    (V4 (V3 a3 a4 a5) (V3 a6 a7 a8) (V3 a9 a10 a11) (V3 a12 a13 a14))++instance ToLinear (Matrix 4 4 a) where+  type L (Matrix 4 4 a) = V4 (V4 a)+  toLinear (Matrix _ [[a4, a5, a6, a7], [a8, a9, a10, a11], [a12, a13, a14, a15], [a16, a17, a18, a19]]) =+    (V4 (V4 a4 a5 a6 a7) (V4 a8 a9 a10 a11) (V4 a12 a13 a14 a15) (V4 a16 a17 a18 a19))
+ src/QLinear/Matrix.hs view
@@ -0,0 +1,3 @@+module QLinear.Matrix (Matrix, Vector, size, value) where++import Internal.Matrix
+ src/QLinear/Operations.hs view
@@ -0,0 +1,149 @@+{-# LANGUAGE RankNTypes #-}++module QLinear.Operations+  ( length,+    mulMatricesWith,+    neg,+    transpose,+    zipMatricesWith,+    det,+    algebraicComplement,+    algebraicComplement',+    adjugate,+    inverted,+    (*~),+    (~*~),+    (~+),+    (+~),+    (~+~),+    (~-~),+  )+where++import qualified Data.List as List+import Data.List.Split (chunksOf)+import Data.Tuple+import Internal.Determinant (adjugate, algebraicComplement, algebraicComplement', det)+import Internal.Matrix+import Internal.Quasi.Matrix.Quasi+import Prelude hiding (length)++-- | Adds two matrices+--+-- >>> [matrix| 1 2 |] ~+~ [matrix| 2 3 |]+-- [3,5]+(~+~) :: Num a => Matrix m n a -> Matrix m n a -> Matrix m n a+(~+~) = zipMatricesWith (+)++-- | Multuplies all elements of matrix __m__ by __k__+--+-- >>> 5 *~ [matrix| 1 2 3; 4 5 6 |]+-- [5,10,15]+-- [20,25,30]+(*~) ::+  Num a =>+  -- | k+  a ->+  -- | m+  Matrix m n a ->+  Matrix m n a+(*~) n = fmap (n *)++-- | Adds __a__ to all elements of matrix __m__+--+-- >>> [matrix| 1 2 3 |] ~+ 8+-- [9,10,11]+(~+) ::+  Num a =>+  Matrix m n a ->+  a ->+  Matrix m n a+(~+) m n = (+ n) <$> m++-- | Flipped __~+__ :)+(+~) :: Num a => a -> Matrix m n a -> Matrix m n a+(+~) = flip (~+)++-- | Substracts second matrix from first one+--+-- >>> [matrix| 1 2 3 |] ~-~ [matrix| 3 2 1 |]+-- [-2,0,2]+(~-~) :: Num a => Matrix m n a -> Matrix m n a -> Matrix m n a+(~-~) = zipMatricesWith (-)++-- | Multiplies two matrix+--+-- >>> [matrix| 1 2; 3 4 |] ~*~ [matrix| 1; 2 |]+-- [5]+-- [11]+(~*~) :: Num a => Matrix m n a -> Matrix n k a -> Matrix m k a+(~*~) = mulMatricesWith (*) (+)++-- | Generalized matrices multiplication+mulMatricesWith ::+  -- | operation "__*__"+  (a -> b -> c) ->+  -- | operation "__+__"+  (c -> c -> c) ->+  Matrix m n a ->+  Matrix n k b ->+  Matrix m k c+mulMatricesWith mul add (Matrix (m, _) left) (Matrix (_, k) right) =+  Matrix (m, k) $+    chunksOf k [foldl1 add $ zipWith mul line column | line <- left, column <- List.transpose right]++-- | Generalized matrices addition+zipMatricesWith ::+  -- | operation "__+__"+  (a -> b -> c) ->+  Matrix m n a ->+  Matrix m n b ->+  Matrix m n c+zipMatricesWith op (Matrix size l) (Matrix _ r) =+  Matrix size $ zipWith (zipWith op) l r++-- | Transposes matrix+--+-- >>> transpose [matrix| 1 2 3; 4 5 6 |]+-- [1,4]+-- [2,5]+-- [3,6]+transpose :: Matrix m n a -> Matrix n m a+transpose (Matrix size matrix) = Matrix (swap size) (List.transpose matrix)++-- | Nagates all elements of matrix+--+-- >>> neg [matrix| 1 2 3 |]+-- [-1,-2,-3]+neg :: Num a => Matrix m n a -> Matrix m n a+neg = ((-1) *~)++-- | Length of vector+--+-- >>> length [vector| 3 4 |]+-- 5.0+-- >>> length [vector| 1 1 |]+-- 1.4142135623730951+length :: (Real a, Floating b) => Vector n a -> b+length (Matrix _ matrix) = sqrt $ sum $ squares+  where+    toFloating = realToFrac :: (Real a, Floating b) => a -> b+    squares = map ((** 2) . toFloating) $ concat matrix++-- | Inverted matrix+--+-- >>> inverted [matrix| 1 2; 3 4|]+-- Just [-2.0,1.0]+--      [1.5,-0.5]+-- >>> inverted [matrix| 1 4; 1 4|]+-- Nothing+inverted :: forall a b n. (Fractional b, Eq a, Real a) => Matrix n n a -> Maybe (Matrix n n b)+inverted (Matrix size@(1, 1) [[a]]) = if a /= 0 then Just (Matrix size [[1.0 / toFloating a]]) else Nothing+  where+    toFloating = realToFrac :: (Real a, Fractional b) => a -> b+inverted matrix = if determinant /= 0 then Just $ ((invertedDet *) . toFloating) <$> adj else Nothing+  where+    toFloating = realToFrac :: (Real a, Fractional b) => a -> b+    determinant = det matrix+    invertedDet = 1.0 / toFloating determinant+    adj = adjugate matrix
+ test/Main.hs view
@@ -0,0 +1,204 @@+{-# LANGUAGE BlockArguments #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE IncoherentInstances #-}+{-# LANGUAGE QuasiQuotes #-}+{-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE TypeApplications #-}++module Main where++import Data.Maybe+import Data.Typeable+import Internal.Determinant+import Internal.Matrix+import Linear (V2 (..), V3 (..), V4 (..))+import Linear.Matrix (det22, det33, det44)+import QLinear.Constructor.Matrix+import QLinear.Constructor.Operator+import QLinear.Identity+import QLinear.Index+import QLinear.Operations+import Test.Hspec++newtype TestMatrix m n a = TestMatrix (Matrix m n a)++instance {-# OVERLAPPING #-} Eq (TestMatrix m n Double) where+  TestMatrix (Matrix ls lv) == TestMatrix (Matrix rs rv) = ls == rs && (all id $ zipWith (\l r -> abs (l - r) < threshold) (concat lv) (concat rv))+    where+      threshold = 0.0000001++instance Eq a => Eq (TestMatrix m n a) where+  TestMatrix (Matrix ls lv) == TestMatrix (Matrix rs rv) = ls == rs && lv == rv++instance Show a => Show (TestMatrix m n a) where+  show (TestMatrix (Matrix s l)) = show s <> " " <> show l++matrixEq :: (Show a, Eq a) => Matrix m n a -> Matrix m n a -> Expectation+matrixEq l r = TestMatrix l `shouldBe` TestMatrix r++matrixEqDouble :: Matrix m n Double -> Matrix m n Double -> Expectation+matrixEqDouble l r = TestMatrix l `shouldBe` TestMatrix r++main :: IO ()+main = hspec do+  describe "QuasiQuoter" do+    it "build matrix by quasi qouter" do+      {- will not be compiled -}+      -- [matrix| 1; 2 3; 4 5 6 |] `matrixEq` undefined+      [matrix| 1 |] `matrixEq` (Matrix (1, 1) [[1]] :: Matrix 1 1 Int)+      [matrix| 1 2 |] `matrixEq` (Matrix (1, 2) [[1, 2]] :: Matrix 1 2 Int)+      [matrix| 1 2 3 |] `matrixEq` (Matrix (1, 3) [[1, 2, 3]] :: Matrix 1 3 Int)+      [matrix| 1 2 3 4|] `matrixEq` (Matrix (1, 4) [[1, 2, 3, 4]] :: Matrix 1 4 Int)+      [matrix| 1; 2; 3 |] `matrixEq` (Matrix (3, 1) [[1], [2], [3]] :: Matrix 3 1 Int)+      [matrix| (1 + 2) -4 5.0 |] `matrixEq` (Matrix (1, 3) [[3, -4, 5]] :: Matrix 1 3 Double)+      [matrix| ("foo") ("bar") |] `matrixEq` (Matrix (1, 2) [["foo", "bar"]] :: Matrix 1 2 String)+    it "build vector by quasi qouter" do+      [vector| 1 2 3 |] `matrixEq` (Matrix (3, 1) [[1], [2], [3]] :: Vector 3 Int)+    it "build operator by quasi quoter" do+      [operator| (x, y) => (x, y) |] `matrixEq` (Matrix (2, 2) [[1, 0], [0, 1]] :: Matrix 2 2 Int)+      [operator| (x, y) => (2 * x, y) |] `matrixEq` (Matrix (2, 2) [[2, 0], [0, 1]] :: Matrix 2 2 Int)+      [operator| (x, y) => (y, x) |] `matrixEq` (Matrix (2, 2) [[0, 1], [1, 0]] :: Matrix 2 2 Int)+      [operator| (x, y) => (x + 2 * y, y) |] `matrixEq` (Matrix (2, 2) [[1, 2], [0, 1]] :: Matrix 2 2 Int)++  describe "Addition" do+    it "two matrices" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~+~ [matrix| 1 2 3 |] `matrixEq` undefined+      [matrix| 1 2 |] ~+~ [matrix| 3 4 |] `matrixEq` [matrix| 4 6 |]+      [matrix| 1 2; 3 4 |] ~+~ [matrix| 3 4; 5 6 |] `matrixEq` [matrix| 4 6; 8 10 |]+    it "matrix and identity matrix" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~+~ e `matrixEq` undefined+      [matrix| 1 2; 3 4 |] ~+~ e `matrixEq` [matrix| 2 2; 3 5 |]+      e ~+~ [matrix| 1 2; 3 4 |] `matrixEq` [matrix| 2 2; 3 5 |]+    it "two identity matrices" do+      {- will not be compiled -}+      -- (e :: Matrix 3 3 Int) ~+~ (e :: Matrix 4 4 Int) `matrixEq` undefined+      (e :: Matrix 3 3 Int) ~+~ e `matrixEq` [matrix| 2 0 0; 0 2 0; 0 0 2 |]+    it "two vectors" do+      {- will not be compiled -}+      -- [vector| 1 2 |] ~+~ [vector| 1 2 3 |] `matrixEq` undefined+      [vector| 1 2 |] ~+~ [vector| 3 4 |] `matrixEq` [vector| 4 6 |]++  describe "Substaction" do+    it "two matrices" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~-~ [matrix| 1 2 3 |] `matrixEq` undefined+      [matrix| 1 2 |] ~-~ [matrix| 3 4 |] `matrixEq` [matrix| -2 -2 |]+      [matrix| 1 2; 3 4 |] ~-~ [matrix| 3 4; 5 6 |] `matrixEq` [matrix| -2 -2; -2 -2 |]+    it "matrix and identity matrix" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~-~ e `matrixEq` undefined+      [matrix| 1 2; 3 4 |] ~-~ e `matrixEq` [matrix| 0 2; 3 3 |]+      e ~-~ [matrix| 1 2; 3 4 |] `matrixEq` [matrix| 0 -2; -3 -3 |]+    it "two identity matrices" do+      {- will not be compiled -}+      -- (e :: Matrix 3 3 Int) ~-~ (e :: Matrix 4 4 Int) `matrixEq` undefined+      (e :: Matrix 3 3 Int) ~-~ e `matrixEq` [matrix| 0 0 0; 0 0 0; 0 0 0 |]+    it "two vectors" do+      {- will not be compiled -}+      -- [vector| 1 2 |] ~-~ [vector| 1 2 3 |] `matrixEq` undefined+      [vector| 1 2 |] ~-~ [vector| 3 4 |] `matrixEq` [vector| -2 -2 |]++  describe "Multiplication" do+    it "two matrices" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~*~ [matrix| 1 2 |] `matrixEq` undefined+      [matrix| 2 |] ~*~ [matrix| 3 |] `matrixEq` [matrix| 6 |]+      [matrix| 1 2 |] ~*~ [matrix| 1; 2 |] `matrixEq` [matrix| 5 |]+      [matrix| 1 2 |] ~*~ [matrix| 1 2; 3 4 |] `matrixEq` [matrix| 7 10 |]+      [matrix| 1 2; 3 4 |] ~*~ [matrix| 2 3; 4 5 |] `matrixEq` [matrix| 10 13; 22 29 |]+    it "matrix and identity matrix" do+      {- will not be compiled -}+      -- [matrix| 1 2 |] ~*~ (e :: Matrix 3 3 Int) `matrixEq` undefined+      [matrix| 1 2 |] ~*~ e `matrixEq` [matrix| 1 2 |]+      [matrix| 1 2; 3 4 |] ~*~ e `matrixEq` [matrix| 1 2; 3 4 |]+      e ~*~ [matrix| 1 2; 3 4 |] `matrixEq` [matrix| 1 2; 3 4 |]+    it "two identity matrices" do+      {- will not be compiled -}+      -- (e :: Matrix 3 3 Int) ~*~ (e :: Matrix 4 4 Int) `matrixEq` undefined+      (e :: Matrix 3 3 Int) ~*~ e `matrixEq` e++  describe "Determinant" do+    it "not identity matrix" do+      let [a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p] = [1, 3, 5, 7, 8, 2, 4, 6, 9, 10, 11, 13, 15, 17, 12, 14]+      {- will not be compiled -}+      --  det [matrix| 1 2 |] `shouldBe` undefined+      det [matrix| a |] `shouldBe` a+      det [matrix| a b; c d |] `shouldBe` det22 (V2 (V2 a b) (V2 c d))+      det [matrix| a b c; d e f; g h i |] `shouldBe` det33 (V3 (V3 a b c) (V3 d e f) (V3 g h i))+      det+        [matrix| +        a b c d; +        e f g h; +        i j k l; +        m n o p |]+        `shouldBe` det44 (V4 (V4 a b c d) (V4 e f g h) (V4 i j k l) (V4 m n o p))+      det+        [matrix| +        -1  2  3  6   5; +        10  7 14  9  18; +        11 22 13 26  15; +        30 17 34 19  38; +        21 42 23 46 -25 |]+        `shouldBe` (-125200)++  describe "Identity matrix" do+    it "multiplication" do+      let m = [matrix| 1 2 3 4; 5 6 7 8 |]+      let i = e :: Matrix 3 3 Int+      (m ~*~ e) `matrixEq` (e ~*~ m)+      (i ~*~ e) `matrixEq` e+      (e ~*~ i) `matrixEq` e+    it "determinant" do+      det (e :: Matrix 1 1 Int) `shouldBe` 1+      det (e :: Matrix 2 2 Int) `shouldBe` 1+      det (e :: Matrix 3 3 Int) `shouldBe` 1+      det (e :: Matrix 4 4 Int) `shouldBe` 1+      det (e :: Matrix 5 5 Int) `shouldBe` 1++  describe "Algebraic complement" do+    it "typesafe" do+      {- will not be compiled-}+      -- algebraicComplement [matrix| 1 2; 3 4 |] (Index @1 @3) `shouldBe` 4+      algebraicComplement [matrix| 1 2; 3 4 |] (Index @1 @1) `shouldBe` 4+      algebraicComplement [matrix| 1 2; 3 4 |] (Index @1 @2) `shouldBe` (-3)+    it "not typesafe" do+      algebraicComplement' [matrix| 1 2; 3 4 |] (1, 3) `shouldBe` Nothing+      algebraicComplement' [matrix| 1 2; 3 4 |] (1, 1) `shouldBe` Just 4+      algebraicComplement' [matrix| 1 2; 3 4 |] (1, 2) `shouldBe` Just (-3)++  describe "Inverted matrix" do+    let a = [matrix| 1 -2 3; 0 4 -1; 7 8 9 |] :: Matrix 3 3 Double+    let b = [matrix| 1 8 3; 50 4 -1; 5 0 1 |] :: Matrix 3 3 Double+    let Just invA = inverted a+    let Just invB = inverted b+    it "(A^(-1))^(-1) = A" do+      let Just res = flip matrixEqDouble a <$> inverted invA+      res+    it "A^(-1) * A = E" do+      let Just inv = inverted a+      a ~*~ invA `matrixEqDouble` e+      invA ~*~ a `matrixEqDouble` e+    it "(AB)^(-1) = B^(-1)A^(-1)" do+      let Just invAB = inverted $ a ~*~ b+      invAB `matrixEqDouble` (invB ~*~ invA)+    it "E(^1) = E" do+      let Just invE1 = inverted (e :: Matrix 1 1 Double)+      let Just invE2 = inverted (e :: Matrix 2 2 Double)+      let Just invE3 = inverted (e :: Matrix 3 3 Double)+      let Just invE4 = inverted (e :: Matrix 4 4 Double)+      let Just invE5 = inverted (e :: Matrix 5 5 Double)+      invE1 `matrixEqDouble` e+      invE2 `matrixEqDouble` e+      invE3 `matrixEqDouble` e+      invE4 `matrixEqDouble` e+      invE5 `matrixEqDouble` e+    it "(A^T)^(-1) = (A^(-1))^T" do+      let Just invAT = inverted $ transpose a+      invAT `matrixEq` (transpose invA)+    it "(kA)^(-1) = k^(-1)A^(-1)" do+      let k = 5.5 :: Double+      let Just invkA = inverted $ k *~ a+      invkA `matrixEqDouble` ((1 / k) *~ invA)