diff --git a/LICENSE b/LICENSE
new file mode 100644
--- /dev/null
+++ b/LICENSE
@@ -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.
diff --git a/README.md b/README.md
new file mode 100644
--- /dev/null
+++ b/README.md
@@ -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.
diff --git a/qlinear.cabal b/qlinear.cabal
new file mode 100644
--- /dev/null
+++ b/qlinear.cabal
@@ -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
diff --git a/src/Internal/Determinant.hs b/src/Internal/Determinant.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Determinant.hs
@@ -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
diff --git a/src/Internal/Matrix.hs b/src/Internal/Matrix.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Matrix.hs
@@ -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
diff --git a/src/Internal/Quasi/Matrix/Parser.hs b/src/Internal/Quasi/Matrix/Parser.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Matrix/Parser.hs
@@ -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)
diff --git a/src/Internal/Quasi/Matrix/Quasi.hs b/src/Internal/Quasi/Matrix/Quasi.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Matrix/Quasi.hs
@@ -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
diff --git a/src/Internal/Quasi/Operator/Parser.hs b/src/Internal/Quasi/Operator/Parser.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Operator/Parser.hs
@@ -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)
diff --git a/src/Internal/Quasi/Operator/Quasi.hs b/src/Internal/Quasi/Operator/Quasi.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Operator/Quasi.hs
@@ -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)
diff --git a/src/Internal/Quasi/Parser.hs b/src/Internal/Quasi/Parser.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Parser.hs
@@ -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
diff --git a/src/Internal/Quasi/Quasi.hs b/src/Internal/Quasi/Quasi.hs
new file mode 100644
--- /dev/null
+++ b/src/Internal/Quasi/Quasi.hs
@@ -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
diff --git a/src/QLinear.hs b/src/QLinear.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear.hs
@@ -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
diff --git a/src/QLinear/Constructor/Matrix.hs b/src/QLinear/Constructor/Matrix.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Constructor/Matrix.hs
@@ -0,0 +1,3 @@
+module QLinear.Constructor.Matrix (matrix, vector) where
+
+import Internal.Quasi.Matrix.Quasi
diff --git a/src/QLinear/Constructor/Operator.hs b/src/QLinear/Constructor/Operator.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Constructor/Operator.hs
@@ -0,0 +1,3 @@
+module QLinear.Constructor.Operator (operator) where
+
+import Internal.Quasi.Operator.Quasi
diff --git a/src/QLinear/Identity.hs b/src/QLinear/Identity.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Identity.hs
@@ -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
diff --git a/src/QLinear/Index.hs b/src/QLinear/Index.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Index.hs
@@ -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)
diff --git a/src/QLinear/Integration/Linear/From.hs b/src/QLinear/Integration/Linear/From.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Integration/Linear/From.hs
@@ -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
+    |]
diff --git a/src/QLinear/Integration/Linear/To.hs b/src/QLinear/Integration/Linear/To.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Integration/Linear/To.hs
@@ -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))
diff --git a/src/QLinear/Matrix.hs b/src/QLinear/Matrix.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Matrix.hs
@@ -0,0 +1,3 @@
+module QLinear.Matrix (Matrix, Vector, size, value) where
+
+import Internal.Matrix
diff --git a/src/QLinear/Operations.hs b/src/QLinear/Operations.hs
new file mode 100644
--- /dev/null
+++ b/src/QLinear/Operations.hs
@@ -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
diff --git a/test/Main.hs b/test/Main.hs
new file mode 100644
--- /dev/null
+++ b/test/Main.hs
@@ -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)