{-# LANGUAGE RankNTypes #-}
{-# LANGUAGE ScopedTypeVariables #-}


module Graphics.HaGL.Numerical (
    Mat, Vec, RowVec,
    m0, (%|), (%-),
    horConcat, vertConcat,
    toList, 
    fromMapping,
    fromList,
    eltAt, eltsAt, matCol,
    Graphics.HaGL.Numerical.length,
    distance,
    dot,
    cross,
    normalize,
    transpose,
    matMult,
    outerProduct,
    determinant,
    inverse
) where

import Data.Array
import Data.List (intercalate, length)
import Data.Proxy
import GHC.TypeLits

import qualified Data.Matrix as Mat (fromList, toList, inverse, detLU)

import Graphics.HaGL.Util.Types (FinList, flToList)
import Graphics.HaGL.Util (warn)

-- | A matrix with @p@ rows, @q@ columns, and element type @t@
data Mat (p :: Nat) (q :: Nat) (t :: *) where
    Mat :: Array (Int, Int) t -> Mat p q t
    
-- | A column vector with @n@ elements and element type @t@
type Vec n t = Mat n 1 t

type RowVec n t = Mat 1 n t

m0 :: Mat 1 0 t
m0 :: forall t. Mat 1 0 t
m0 = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((Int
0, Int
0),(-Int
1, -Int
1)) []

infixr 9 %|
infixr 8 %-

(%|) :: t -> RowVec q t -> RowVec (q + 1) t
t
x %| :: forall t (q :: Nat). t -> RowVec q t -> RowVec (q + 1) t
%| Mat Array (Int, Int) t
xs = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
0, Int
q forall a. Num a => a -> a -> a
+ Int
1)) (t
x forall a. a -> [a] -> [a]
: forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs)  where
    ((Int
0, Int
0), (Int
_, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs

(%-) :: RowVec q t -> Mat p q t -> Mat (p + 1) q t
Mat Array (Int, Int) t
xs %- :: forall (q :: Nat) t (p :: Nat).
RowVec q t -> Mat p q t -> Mat (p + 1) q t
%- Mat Array (Int, Int) t
ys = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
p forall a. Num a => a -> a -> a
+ Int
1, Int
q)) (forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs forall a. [a] -> [a] -> [a]
++ forall i e. Array i e -> [e]
elems Array (Int, Int) t
ys) where
    ((Int
0, Int
0), (Int
p, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
ys

horConcat :: Mat p q1 t -> Mat p q2 t -> Mat p (q1 + q2) t
horConcat :: forall (p :: Nat) (q1 :: Nat) t (q2 :: Nat).
Mat p q1 t -> Mat p q2 t -> Mat p (q1 + q2) t
horConcat Mat p q1 t
m1 Mat p q2 t
m2 = forall (p :: Nat) (q :: Nat) t. Mat p q t -> Mat q p t
transpose (forall (p1 :: Nat) (q :: Nat) t (p2 :: Nat).
Mat p1 q t -> Mat p2 q t -> Mat (p1 + p2) q t
vertConcat (forall (p :: Nat) (q :: Nat) t. Mat p q t -> Mat q p t
transpose Mat p q1 t
m1) (forall (p :: Nat) (q :: Nat) t. Mat p q t -> Mat q p t
transpose Mat p q2 t
m2))

vertConcat :: Mat p1 q t -> Mat p2 q t -> Mat (p1 + p2) q t
vertConcat :: forall (p1 :: Nat) (q :: Nat) t (p2 :: Nat).
Mat p1 q t -> Mat p2 q t -> Mat (p1 + p2) q t
vertConcat (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
p1 forall a. Num a => a -> a -> a
+ Int
p2 forall a. Num a => a -> a -> a
+ Int
1, Int
q)) (forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs forall a. [a] -> [a] -> [a]
++ forall i e. Array i e -> [e]
elems Array (Int, Int) t
ys) where
    ((Int
0, Int
0), (Int
p1, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs
    ((Int
0, Int
0), (Int
p2, Int
_)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
ys

instance Eq t => Eq (Mat p q t) where
    (Mat Array (Int, Int) t
xs) == :: Mat p q t -> Mat p q t -> Bool
== (Mat Array (Int, Int) t
ys) = Array (Int, Int) t
xs forall a. Eq a => a -> a -> Bool
== Array (Int, Int) t
ys

instance Functor (Mat p q) where
    fmap :: forall a b. (a -> b) -> Mat p q a -> Mat p q b
fmap a -> b
f (Mat Array (Int, Int) a
xs) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat (forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap a -> b
f Array (Int, Int) a
xs)

instance Foldable (Mat p q) where
    foldr :: forall a b. (a -> b -> b) -> b -> Mat p q a -> b
foldr a -> b -> b
f b
z (Mat Array (Int, Int) a
xs) = forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr a -> b -> b
f b
z Array (Int, Int) a
xs

instance (KnownNat p, KnownNat q) => Applicative (Mat p q) where
    pure :: forall a. a -> Mat p q a
pure a
x = let p :: Int
p = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy p)
                 q :: Int
q = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy q) 
             in forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0),(Int
pforall a. Num a => a -> a -> a
-Int
1, Int
qforall a. Num a => a -> a -> a
-Int
1)) (forall a. Int -> a -> [a]
replicate (Int
pforall a. Num a => a -> a -> a
*Int
q) a
x)
    Mat Array (Int, Int) (a -> b)
fs <*> :: forall a b. Mat p q (a -> b) -> Mat p q a -> Mat p q b
<*> Mat Array (Int, Int) a
xs = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray (forall i e. Array i e -> (i, i)
bounds Array (Int, Int) (a -> b)
fs) 
        [(Array (Int, Int) (a -> b)
fs forall i e. Ix i => Array i e -> i -> e
! (Int, Int)
ind) (Array (Int, Int) a
xs forall i e. Ix i => Array i e -> i -> e
! (Int, Int)
ind) | (Int, Int)
ind <- forall a. Ix a => (a, a) -> [a]
range (forall i e. Array i e -> (i, i)
bounds Array (Int, Int) (a -> b)
fs)]

instance (KnownNat p, KnownNat q, Num t) => Num (Mat p q t) where
    Mat p q t
m1 + :: Mat p q t -> Mat p q t -> Mat p q t
+ Mat p q t
m2 = forall a. Num a => a -> a -> a
(+) forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mat p q t
m1 forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Mat p q t
m2
    Mat p q t
m1 - :: Mat p q t -> Mat p q t -> Mat p q t
- Mat p q t
m2 = (-) forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mat p q t
m1 forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Mat p q t
m2
    Mat p q t
m1 * :: Mat p q t -> Mat p q t -> Mat p q t
* Mat p q t
m2 = forall a. Num a => a -> a -> a
(*) forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mat p q t
m1 forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Mat p q t
m2
    negate :: Mat p q t -> Mat p q t
negate = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Num a => a -> a
negate
    abs :: Mat p q t -> Mat p q t
abs = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Num a => a -> a
abs
    signum :: Mat p q t -> Mat p q t
signum = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Num a => a -> a
signum
    fromInteger :: Integer -> Mat p q t
fromInteger = forall (f :: * -> *) a. Applicative f => a -> f a
pure forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Num a => Integer -> a
fromInteger

instance (KnownNat p, KnownNat q, Fractional t) => Fractional (Mat p q t) where
    Mat p q t
m1 / :: Mat p q t -> Mat p q t -> Mat p q t
/ Mat p q t
m2 = forall a. Fractional a => a -> a -> a
(/) forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mat p q t
m1 forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Mat p q t
m2
    fromRational :: Rational -> Mat p q t
fromRational = forall (f :: * -> *) a. Applicative f => a -> f a
pure forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Fractional a => Rational -> a
fromRational

-- TODO: Integral support?

instance Show t => Show (Mat p q t) where
    show :: Mat p q t -> String
show (Mat Array (Int, Int) t
xs) = forall a. [a] -> [[a]] -> [a]
intercalate String
"\n" [Int -> String
showRow Int
i | Int
i <- [Int
0..Int
p]] forall a. [a] -> [a] -> [a]
++ String
"\n" where
        showRow :: Int -> String
showRow Int
i = forall a. [a] -> [[a]] -> [a]
intercalate String
" " [Int -> Int -> String
showElt Int
i Int
j | Int
j <- [Int
0..Int
q]]
        showElt :: Int -> Int -> String
showElt Int
i Int
j = forall a. Show a => a -> String
show forall a b. (a -> b) -> a -> b
$ Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
j)
        ((Int
0, Int
0), (Int
p, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs


toList :: Mat p q t -> [t]
toList :: forall (p :: Nat) (q :: Nat) a. Mat p q a -> [a]
toList (Mat Array (Int, Int) t
xs) = forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs

-- | Construct a matrix from a mapping that maps indices @(i, j)@ to
-- the element at row @i@ and column @j@
fromMapping :: forall p q t. (KnownNat p, KnownNat q) => ((Int, Int) -> t) -> Mat p q t
fromMapping :: forall (p :: Nat) (q :: Nat) t.
(KnownNat p, KnownNat q) =>
((Int, Int) -> t) -> Mat p q t
fromMapping (Int, Int) -> t
f =
    let p :: Int
p = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy p)
        q :: Int
q = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy q)  
    in forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((Int
0, Int
0), (Int
pforall a. Num a => a -> a -> a
-Int
1, Int
qforall a. Num a => a -> a -> a
-Int
1))
        [((Int, Int)
ind, (Int, Int) -> t
f (Int, Int)
ind) | (Int, Int)
ind <- forall a. Ix a => (a, a) -> [a]
range ((Int
0, Int
0), (Int
pforall a. Num a => a -> a -> a
-Int
1, Int
qforall a. Num a => a -> a -> a
-Int
1))]

-- | Construct a matrix from a list of the matrix elements in row-major order
fromList :: forall p q t. (KnownNat p, KnownNat q) => [t] -> Mat p q t
fromList :: forall (p :: Nat) (q :: Nat) t.
(KnownNat p, KnownNat q) =>
[t] -> Mat p q t
fromList [t]
xs =
    let p :: Int
p = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy p)
        q :: Int
q = forall a. Num a => Integer -> a
fromInteger forall a b. (a -> b) -> a -> b
$ forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy q)  
    in forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
pforall a. Num a => a -> a -> a
-Int
1, Int
qforall a. Num a => a -> a -> a
-Int
1)) [t]
xs

eltAt :: Mat p q t -> (Int, Int) -> t
eltAt :: forall (p :: Nat) (q :: Nat) t. Mat p q t -> (Int, Int) -> t
eltAt (Mat Array (Int, Int) t
xs) (Int
i, Int
j) = Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
j)

eltsAt :: Mat p q t -> FinList m (Int, Int) -> Vec m t
eltsAt :: forall (p :: Nat) (q :: Nat) t (m :: Nat).
Mat p q t -> FinList m (Int, Int) -> Vec m t
eltsAt (Mat Array (Int, Int) t
xs) FinList m (Int, Int)
fl = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
mforall a. Num a => a -> a -> a
-Int
1, Int
0)) [t]
vs where
    vs :: [t]
vs = [Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int, Int)
ind | (Int, Int)
ind <- [(Int, Int)]
inds]
    inds :: [(Int, Int)]
inds = forall (n :: Nat) t. FinList n t -> [t]
flToList FinList m (Int, Int)
fl
    m :: Int
m = forall (t :: * -> *) a. Foldable t => t a -> Int
Data.List.length [(Int, Int)]
inds

matCol :: Mat p q t -> Int -> Vec p t
matCol :: forall (p :: Nat) (q :: Nat) t. Mat p q t -> Int -> Vec p t
matCol (Mat Array (Int, Int) t
xs) Int
j = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
p, Int
0)) [t]
vs where
    vs :: [t]
vs = [Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
j) | Int
i <- [Int
0..Int
p]]
    ((Int
0, Int
0), (Int
p, Int
_)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs


length :: Floating t => Vec n t -> t
length :: forall t (n :: Nat). Floating t => Vec n t -> t
length (Mat Array (Int, Int) t
xs) = forall a. Floating a => a -> a
sqrt forall a b. (a -> b) -> a -> b
$ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum
    [(Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
0)) forall a. Floating a => a -> a -> a
** t
2 | (Int
i, Int
0) <- forall a. Ix a => (a, a) -> [a]
range forall a b. (a -> b) -> a -> b
$ forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs]

distance :: Floating t => Vec n t -> Vec n t -> t
distance :: forall t (n :: Nat). Floating t => Vec n t -> Vec n t -> t
distance (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = forall a. Floating a => a -> a
sqrt forall a b. (a -> b) -> a -> b
$ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum
    [((Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
0)) forall a. Num a => a -> a -> a
- (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
0))) forall a. Floating a => a -> a -> a
** t
2 | (Int
i, Int
0) <- forall a. Ix a => (a, a) -> [a]
range forall a b. (a -> b) -> a -> b
$ forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs]

dot :: Num t => Vec n t -> Vec n t -> t
dot :: forall t (n :: Nat). Num t => Vec n t -> Vec n t -> t
dot (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = 
    forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [(Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
0)) | (Int
i, Int
0) <- forall a. Ix a => (a, a) -> [a]
range forall a b. (a -> b) -> a -> b
$ forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs]

cross :: Num t => Vec 3 t -> Vec 3 t -> Vec 3 t
cross :: forall t. Num t => Vec 3 t -> Vec 3 t -> Vec 3 t
cross (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray (forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs) [t]
zs where
    zs :: [t]
zs = [(Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
1, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
2, Int
0)) forall a. Num a => a -> a -> a
- (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
2, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
1, Int
0)),
          (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
0, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
2, Int
0)) forall a. Num a => a -> a -> a
- (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
2, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
0, Int
0)),
          (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
0, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
1, Int
0)) forall a. Num a => a -> a -> a
- (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
1, Int
0)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
0, Int
0))]

normalize :: Floating t => Vec n t -> Vec n t
normalize :: forall t (n :: Nat). Floating t => Vec n t -> Vec n t
normalize v :: Vec n t
v@(Mat Array (Int, Int) t
xs) = 
    forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray (forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs) (forall a b. (a -> b) -> [a] -> [b]
map (forall a. Fractional a => a -> a -> a
/ forall t (n :: Nat). Floating t => Vec n t -> t
Graphics.HaGL.Numerical.length Vec n t
v) (forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs))


transpose :: Mat p q t -> Mat q p t
transpose :: forall (p :: Nat) (q :: Nat) t. Mat p q t -> Mat q p t
transpose (Mat Array (Int, Int) t
xs) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((Int
0, Int
0), (Int
q, Int
p)) [((Int, Int), t)]
xst where
    xst :: [((Int, Int), t)]
xst =  [((Int
j, Int
i), Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
j)) | (Int
i, Int
j) <- forall a. Ix a => (a, a) -> [a]
range ((Int
0, Int
0), (Int
p, Int
q))]
    ((Int
0, Int
0), (Int
p, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs

matMult :: Num t => Mat p q t -> Mat q r t -> Mat p r t
matMult :: forall t (p :: Nat) (q :: Nat) (r :: Nat).
Num t =>
Mat p q t -> Mat q r t -> Mat p r t
matMult (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((Int
0, Int
0), (Int
p, Int
r)) [((Int, Int), t)]
zs where
    zs :: [((Int, Int), t)]
zs = [((Int
i, Int
j), Int -> Int -> t
matElt Int
i Int
j) | (Int
i, Int
j) <- forall a. Ix a => (a, a) -> [a]
range ((Int
0, Int
0), (Int
p, Int
r))]
    matElt :: Int -> Int -> t
matElt Int
i Int
j = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [(Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
k)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
k, Int
j)) | Int
k <- [Int
0..Int
q]]
    ((Int
0, Int
0), (Int
p, Int
q)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs
    ((Int
0, Int
0), (Int
_, Int
r)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
ys

outerProduct :: Num t => Vec p t -> Vec q t -> Mat p q t
outerProduct :: forall t (p :: Nat) (q :: Nat).
Num t =>
Vec p t -> Vec q t -> Mat p q t
outerProduct (Mat Array (Int, Int) t
xs) (Mat Array (Int, Int) t
ys) = forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((Int
0, Int
0), (Int
p, Int
q)) [((Int, Int), t)]
zs where
    zs :: [((Int, Int), t)]
zs = [((Int
i, Int
j), (Array (Int, Int) t
xs forall i e. Ix i => Array i e -> i -> e
! (Int
i, Int
1)) forall a. Num a => a -> a -> a
* (Array (Int, Int) t
ys forall i e. Ix i => Array i e -> i -> e
! (Int
j, Int
1))) | (Int
i, Int
j) <- forall a. Ix a => (a, a) -> [a]
range ((Int
0, Int
0), (Int
p, Int
q))]
    ((Int
0, Int
0), (Int
p, Int
1)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs
    ((Int
0, Int
0), (Int
q, Int
1)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
ys

determinant :: (Fractional t, Ord t) => Mat p p t -> t
determinant :: forall t (p :: Nat). (Fractional t, Ord t) => Mat p p t -> t
determinant (Mat Array (Int, Int) t
xs) = forall a. (Ord a, Fractional a) => Matrix a -> a
Mat.detLU forall a b. (a -> b) -> a -> b
$ forall a. Int -> Int -> [a] -> Matrix a
Mat.fromList (Int
p forall a. Num a => a -> a -> a
+ Int
1) (Int
p forall a. Num a => a -> a -> a
+ Int
1) (forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs) where
    ((Int
0, Int
0), (Int
p, Int
_)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs

inverse :: (Fractional t, Eq t) => Mat p p t -> Mat p p t
inverse :: forall t (p :: Nat). (Fractional t, Eq t) => Mat p p t -> Mat p p t
inverse (Mat Array (Int, Int) t
xs) =
    case forall a.
(Fractional a, Eq a) =>
Matrix a -> Either String (Matrix a)
Mat.inverse forall a b. (a -> b) -> a -> b
$ forall a. Int -> Int -> [a] -> Matrix a
Mat.fromList (Int
p forall a. Num a => a -> a -> a
+ Int
1) (Int
p forall a. Num a => a -> a -> a
+ Int
1) (forall i e. Array i e -> [e]
elems Array (Int, Int) t
xs) of
        Right Matrix t
m -> forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat forall a b. (a -> b) -> a -> b
$ forall i e. Ix i => (i, i) -> [e] -> Array i e
listArray ((Int
0, Int
0), (Int
p, Int
p)) forall a b. (a -> b) -> a -> b
$ forall a. Matrix a -> [a]
Mat.toList Matrix t
m
        Left String
msg -> forall a. String -> a -> a
warn String
msg (forall t (p :: Nat) (q :: Nat). Array (Int, Int) t -> Mat p q t
Mat Array (Int, Int) t
xs)
    where ((Int
0, Int
0), (Int
p, Int
_)) = forall i e. Array i e -> (i, i)
bounds Array (Int, Int) t
xs