{-# OPTIONS_GHC -fplugin=GHC.TypeLits.KnownNat.Solver -fplugin=GHC.TypeLits.Normalise -fconstraint-solver-iterations=10 #-}
-- | Vectors and Matrices with statically typed dimensions based on storable vectors and using HMatrix where possible.
module Goal.Core.Vector.Storable
    ( -- * Vector
      module Data.Vector.Storable.Sized
      -- ** Construction
    , doubleton
    , range
      -- ** Deconstruction
    , concat
    , breakEvery
    , toPair
    -- ** Computation
    , average
    , zipFold
    -- * Matrix
    , Matrix
    , nRows
    , nColumns
    -- ** Construction
    , fromRows
    , fromColumns
    , matrixIdentity
    , outerProduct
    , sumOuterProduct
    , averageOuterProduct
    , weightedAverageOuterProduct
    , diagonalMatrix
    , fromLowerTriangular
    -- ** Deconstruction
    , toRows
    , toColumns
    , lowerTriangular
    , takeDiagonal
    -- ** Manipulation
    , columnVector
    , rowVector
    , combineTriangles
    , diagonalConcat
    , horizontalConcat
    , verticalConcat
    -- ** Computation
    , trace
    , withMatrix
    -- *** BLAS
    , scale
    , add
    , dotProduct
    , dotMap
    , matrixVectorMultiply
    , matrixMatrixMultiply
    , matrixMap
    , eigens
    , isSemiPositiveDefinite
    , determinant
    , inverseLogDeterminant
    , inverse
    , pseudoInverse
    , matrixRoot
    , transpose
    -- *** Least Squares
    , linearLeastSquares
    , meanSquaredError
    , rSquared
    , l2Norm
    , unsafeCholesky
    -- *** Convolutions
    , crossCorrelate2d
    , convolve2d
    , kernelOuterProduct
    , kernelTranspose
    -- ** Miscellaneous
    , prettyPrintMatrix
    ) where


--- Imports ---


-- Goal --

import Goal.Core.Util hiding (average,breakEvery,range)

-- Unqualified --

import Data.Proxy
import Data.Complex
import Foreign.Storable
import Data.Vector.Storable.Sized
import Numeric.LinearAlgebra (Field,Numeric)
import GHC.TypeNats
import Prelude hiding (concat,foldr1,concatMap,replicate,(++),length,map,sum,zip,and)

-- Qualified --

import qualified Prelude
import qualified Data.Vector.Storable as S
import qualified Goal.Core.Vector.Generic as G
import qualified Data.Vector.Generic.Sized.Internal as G
import qualified Numeric.LinearAlgebra as H
import qualified Data.List as L


--- Generic ---


-- | Matrices with static dimensions (storable).
type Matrix = G.Matrix S.Vector

-- | A fold over pairs of elements of 'Vector's of equal length.
zipFold :: (KnownNat n, Storable x, Storable y)
        => (z -> x -> y -> z)
        -> z
        -> Vector n x
        -> Vector n y
        -> z
{-# INLINE zipFold #-}
zipFold :: (z -> x -> y -> z) -> z -> Vector n x -> Vector n y -> z
zipFold z -> x -> y -> z
f z
z0 Vector n x
xs Vector n y
ys =
    let n :: Int
n = Vector n x -> Int
forall (n :: Nat) a. KnownNat n => Vector n a -> Int
length Vector n x
xs
        foldfun :: z -> Int -> z
foldfun z
z Int
i = z -> x -> y -> z
f z
z (Vector n x -> Int -> x
forall (n :: Nat) a. Storable a => Vector n a -> Int -> a
unsafeIndex Vector n x
xs Int
i) (Vector n y -> Int -> y
forall (n :: Nat) a. Storable a => Vector n a -> Int -> a
unsafeIndex Vector n y
ys Int
i)
     in (z -> Int -> z) -> z -> [Int] -> z
forall (t :: Type -> Type) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
L.foldl' z -> Int -> z
foldfun z
z0 [Int
0..Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]

-- | Concatenates a vector of vectors into a single vector.
concat :: (KnownNat n, Storable x) => Vector m (Vector n x) -> Vector (m*n) x
{-# INLINE concat #-}
concat :: Vector m (Vector n x) -> Vector (m * n) x
concat = Vector m (Vector n x) -> Vector (m * n) x
forall (n :: Nat) (v :: Type -> Type) x (m :: Nat).
(KnownNat n, Vector v x, Vector v (Vector v n x)) =>
Vector v m (Vector v n x) -> Vector v (m * n) x
G.concat

-- | Collect two values into a length 2 'Vector'.
doubleton :: Storable x => x -> x -> Vector 2 x
{-# INLINE doubleton #-}
doubleton :: x -> x -> Vector 2 x
doubleton = x -> x -> Vector 2 x
forall (v :: Type -> Type) x. Vector v x => x -> x -> Vector v 2 x
G.doubleton

-- | The number of rows in the 'Matrix'.
nRows :: forall m n a . KnownNat m => Matrix m n a -> Int
{-# INLINE nRows #-}
nRows :: Matrix m n a -> Int
nRows = Matrix m n a -> Int
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
KnownNat m =>
Matrix v m n a -> Int
G.nRows

-- | The columns of columns in the 'Matrix'.
nColumns :: forall m n a . KnownNat n => Matrix m n a -> Int
{-# INLINE nColumns #-}
nColumns :: Matrix m n a -> Int
nColumns = Matrix m n a -> Int
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
KnownNat n =>
Matrix v m n a -> Int
G.nColumns

-- | Convert a 'Matrix' into a 'Vector' of 'Vector's of rows.
toRows :: (KnownNat m, KnownNat n, Storable x) => Matrix m n x -> Vector m (Vector n x)
{-# INLINE toRows #-}
toRows :: Matrix m n x -> Vector m (Vector n x)
toRows = Matrix m n x -> Vector m (Vector n x)
forall (v :: Type -> Type) a (n :: Nat) (m :: Nat).
(Vector v a, Vector v (Vector v n a), KnownNat n, KnownNat m) =>
Matrix v m n a -> Vector v m (Vector v n a)
G.toRows

-- | Turn a 'Vector' into a single column 'Matrix'.
columnVector :: Vector n a -> Matrix n 1 a
{-# INLINE columnVector #-}
columnVector :: Vector n a -> Matrix n 1 a
columnVector = Vector n a -> Matrix n 1 a
forall (v :: Type -> Type) (n :: Nat) a.
Vector v n a -> Matrix v n 1 a
G.columnVector

-- | Turn a 'Vector' into a single row 'Matrix'.
rowVector :: Vector n a -> Matrix 1 n a
{-# INLINE rowVector #-}
rowVector :: Vector n a -> Matrix 1 n a
rowVector = Vector n a -> Matrix 1 n a
forall (v :: Type -> Type) (n :: Nat) a.
Vector v n a -> Matrix v 1 n a
G.rowVector

-- | Create a 'Matrix' from a 'Vector' of 'Vector's which represent the rows.
fromRows :: (KnownNat n, Storable x) => Vector m (Vector n x) -> Matrix m n x
{-# INLINE fromRows #-}
fromRows :: Vector m (Vector n x) -> Matrix m n x
fromRows = Vector m (Vector n x) -> Matrix m n x
forall (v :: Type -> Type) x (n :: Nat) (m :: Nat).
(Vector v x, Vector v (Vector v n x), KnownNat n) =>
Vector v m (Vector v n x) -> Matrix v m n x
G.fromRows

-- | Uniform partition of an interval into a 'Vector'.
range :: (KnownNat n, Fractional x, Storable x) => x -> x -> Vector n x
{-# INLINE range #-}
range :: x -> x -> Vector n x
range = x -> x -> Vector n x
forall (v :: Type -> Type) (n :: Nat) x.
(Vector v x, KnownNat n, Fractional x) =>
x -> x -> Vector v n x
G.range

-- | Reshapes a length 2 'Vector' into a pair of values.
toPair :: Storable x => Vector 2 x -> (x,x)
{-# INLINE toPair #-}
toPair :: Vector 2 x -> (x, x)
toPair = Vector 2 x -> (x, x)
forall (v :: Type -> Type) a. Vector v a => Vector v 2 a -> (a, a)
G.toPair


--- HMatrix ---


-- | Converts a pure, Storable-based 'Matrix' into an HMatrix matrix.
toHMatrix
    :: forall m n x . (KnownNat n, KnownNat m, H.Element x, Storable x)
    => Matrix m n x
    -> H.Matrix x
{-# INLINE toHMatrix #-}
toHMatrix :: Matrix m n x -> Matrix x
toHMatrix (G.Matrix Vector Vector (m * n) x
mtx) =
    let n :: Int
n = Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n)
        m :: Int
m = Proxy m -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy m
forall k (t :: k). Proxy t
Proxy :: Proxy m)
     in if Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0
           then [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromRows ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Int -> Vector x -> [Vector x]
forall a. Int -> a -> [a]
Prelude.replicate Int
m Vector x
forall a. Storable a => Vector a
S.empty
           else Int -> Vector x -> Matrix x
forall t. Storable t => Int -> Vector t -> Matrix t
H.reshape Int
n (Vector x -> Matrix x) -> Vector x -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector Vector (m * n) x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector Vector (m * n) x
mtx

-- | Converts an HMatrix matrix into a pure, Storable-based 'Matrix'.
fromHMatrix :: Numeric x => H.Matrix x -> Matrix m n x
{-# INLINE fromHMatrix #-}
fromHMatrix :: Matrix x -> Matrix m n x
fromHMatrix = Vector Vector (m * n) x -> Matrix m n x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (m * n) x -> Matrix m n x)
-> (Matrix x -> Vector Vector (m * n) x)
-> Matrix x
-> Matrix m n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector x -> Vector Vector (m * n) x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector Vector (m * n) x)
-> (Matrix x -> Vector x) -> Matrix x -> Vector Vector (m * n) x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.flatten

-- | Convert a 'Matrix' into a 'Vector' of 'Vector's of columns.
toColumns :: (KnownNat m, KnownNat n, Numeric x) => Matrix m n x -> Vector n (Vector m x)
{-# INLINE toColumns #-}
toColumns :: Matrix m n x -> Vector n (Vector m x)
toColumns = Matrix n m x -> Vector n (Vector m x)
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Storable x) =>
Matrix m n x -> Vector m (Vector n x)
toRows (Matrix n m x -> Vector n (Vector m x))
-> (Matrix m n x -> Matrix n m x)
-> Matrix m n x
-> Vector n (Vector m x)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix m n x -> Matrix n m x
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Numeric x) =>
Matrix m n x -> Matrix n m x
transpose

-- | Create a 'Matrix' from a 'Vector' of 'Vector's which represent the columns.
fromColumns :: (KnownNat m, KnownNat n, Numeric x) => Vector n (Vector m x) -> Matrix m n x
{-# INLINE fromColumns #-}
fromColumns :: Vector n (Vector m x) -> Matrix m n x
fromColumns = Matrix n m x -> Matrix m n x
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Numeric x) =>
Matrix m n x -> Matrix n m x
transpose (Matrix n m x -> Matrix m n x)
-> (Vector n (Vector m x) -> Matrix n m x)
-> Vector n (Vector m x)
-> Matrix m n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector n (Vector m x) -> Matrix n m x
forall (n :: Nat) x (m :: Nat).
(KnownNat n, Storable x) =>
Vector m (Vector n x) -> Matrix m n x
fromRows

-- | Breaks a 'Vector' into a Vector of Vectors.
breakEvery :: forall n k a . (KnownNat n, KnownNat k, Storable a) => Vector (n*k) a -> Vector n (Vector k a)
{-# INLINE breakEvery #-}
breakEvery :: Vector (n * k) a -> Vector n (Vector k a)
breakEvery Vector (n * k) a
v0 =
    let k :: Int
k = Proxy k -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy k
forall k (t :: k). Proxy t
Proxy :: Proxy k)
        v :: Vector a
v = Vector (n * k) a -> Vector a
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector (n * k) a
v0
     in (Finite n -> Vector k a) -> Vector n (Vector k a)
forall (n :: Nat) a.
(KnownNat n, Storable a) =>
(Finite n -> a) -> Vector n a
generate (\Finite n
i -> Vector a -> Vector k a
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector a -> Vector k a) -> Vector a -> Vector k a
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Vector a -> Vector a
forall a. Storable a => Int -> Int -> Vector a -> Vector a
S.unsafeSlice (Finite n -> Int
forall (n :: Nat). KnownNat n => Finite n -> Int
finiteInt Finite n
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k) Int
k Vector a
v)


--- BLAS ---


-- | The sum of two 'Vector's.
add :: Numeric x => Vector n x -> Vector n x -> Vector n x
{-# INLINE add #-}
add :: Vector n x -> Vector n x -> Vector n x
add (G.Vector Vector x
v1) (G.Vector Vector x
v2) = Vector x -> Vector n x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector x -> Vector x
forall c. Additive c => c -> c -> c
H.add Vector x
v1 Vector x
v2)

-- | Scalar multiplication of a 'Vector'.
scale :: Numeric x => x -> Vector n x -> Vector n x
{-# INLINE scale #-}
scale :: x -> Vector n x -> Vector n x
scale x
x (G.Vector Vector x
v) = Vector x -> Vector n x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (x -> Vector x -> Vector x
forall t (c :: Type -> Type). Linear t c => t -> c t -> c t
H.scale x
x Vector x
v)

-- | Apply a 'Vector' operation to a 'Matrix'.
withMatrix :: (Vector (n*m) x -> Vector (n*m) x) -> Matrix n m x -> Matrix n m x
{-# INLINE withMatrix #-}
withMatrix :: (Vector (n * m) x -> Vector (n * m) x)
-> Matrix n m x -> Matrix n m x
withMatrix Vector (n * m) x -> Vector (n * m) x
f (G.Matrix Vector (n * m) x
v) = Vector (n * m) x -> Matrix n m x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector (n * m) x -> Matrix n m x)
-> Vector (n * m) x -> Matrix n m x
forall a b. (a -> b) -> a -> b
$ Vector (n * m) x -> Vector (n * m) x
f Vector (n * m) x
v

-- | Returns the lower triangular part of a square matrix.
lowerTriangular :: forall n x . (Storable x, H.Element x, KnownNat n) => Matrix n n x -> Vector (Triangular n) x
{-# INLINE lowerTriangular #-}
lowerTriangular :: Matrix n n x -> Vector (Triangular n) x
lowerTriangular Matrix n n x
mtx =
    let hmtx :: Matrix x
hmtx = Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n n x
mtx
        rws :: [Vector x]
rws = Matrix x -> [Vector x]
forall t. Element t => Matrix t -> [Vector t]
H.toRows Matrix x
hmtx
        rws' :: [Vector x]
rws' = (Int -> Vector x -> Vector x) -> [Int] -> [Vector x] -> [Vector x]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
Prelude.zipWith Int -> Vector x -> Vector x
forall a. Storable a => Int -> Vector a -> Vector a
S.take [Int
1..] [Vector x]
rws
     in Vector x -> Vector (Triangular n) x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector (Triangular n) x)
-> Vector x -> Vector (Triangular n) x
forall a b. (a -> b) -> a -> b
$ [Vector x] -> Vector x
forall a. Storable a => [Vector a] -> Vector a
S.concat [Vector x]
rws'
--    let n = natValInt (Proxy :: Proxy n)
--        idxs = G.Vector . S.fromList
--            $ Prelude.concat [ from2Index n <$> Prelude.zip (repeat k) [0..k] | k <- [0..n-1] ]
--     in backpermute xs idxs

-- | Constructs a `Matrix` from a lower triangular part.
fromLowerTriangular :: forall n x . (Storable x, KnownNat n) => Vector (Triangular n) x -> Matrix n n x
{-# INLINE fromLowerTriangular #-}
fromLowerTriangular :: Vector (Triangular n) x -> Matrix n n x
fromLowerTriangular Vector (Triangular n) x
xs =
    let n :: Int
n = Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n)
        idxs :: Vector (n * n) Int
idxs = (Finite (n * n) -> Int) -> Vector (n * n) Int
forall (n :: Nat) a.
(KnownNat n, Storable a) =>
(Finite n -> a) -> Vector n a
generate ((Int, Int) -> Int
toTriangularIndex ((Int, Int) -> Int)
-> (Finite (n * n) -> (Int, Int)) -> Finite (n * n) -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Int -> (Int, Int)
to2Index Int
n (Int -> (Int, Int))
-> (Finite (n * n) -> Int) -> Finite (n * n) -> (Int, Int)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Finite (n * n) -> Int
forall (n :: Nat). KnownNat n => Finite n -> Int
finiteInt)
     in Vector Vector (n * n) x -> Matrix n n x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (n * n) x -> Matrix n n x)
-> Vector Vector (n * n) x -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ Vector (Triangular n) x
-> Vector (n * n) Int -> Vector Vector (n * n) x
forall a (m :: Nat) (n :: Nat).
Storable a =>
Vector m a -> Vector n Int -> Vector n a
backpermute Vector (Triangular n) x
xs Vector (n * n) Int
idxs

-- | Build a matrix with the given diagonal, lower triangular part given by the
-- first matrix, and upper triangular part given by the second matrix.
combineTriangles
    :: (KnownNat k, Storable x)
    => Vector k x -- ^ Diagonal
    -> Matrix k k x -- ^ Lower triangular source
    -> Matrix k k x -- ^ Upper triangular source
    -> Matrix k k x
{-# INLINE combineTriangles #-}
combineTriangles :: Vector k x -> Matrix k k x -> Matrix k k x -> Matrix k k x
combineTriangles (G.Vector Vector x
diag) Matrix k k x
crs1 Matrix k k x
crs2 =
    Vector k (Vector k x) -> Matrix k k x
forall (n :: Nat) x (m :: Nat).
(KnownNat n, Storable x) =>
Vector m (Vector n x) -> Matrix m n x
fromRows (Vector k (Vector k x) -> Matrix k k x)
-> Vector k (Vector k x) -> Matrix k k x
forall a b. (a -> b) -> a -> b
$ (Finite k -> Vector k x) -> Vector k (Vector k x)
forall (n :: Nat) a.
(KnownNat n, Storable a) =>
(Finite n -> a) -> Vector n a
generate (Vector k (Vector k x)
-> Vector k (Vector k x) -> Finite k -> Vector k x
generator (Matrix k k x -> Vector k (Vector k x)
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Storable x) =>
Matrix m n x -> Vector m (Vector n x)
toRows Matrix k k x
crs1) (Matrix k k x -> Vector k (Vector k x)
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Storable x) =>
Matrix m n x -> Vector m (Vector n x)
toRows Matrix k k x
crs2))
        where
            generator :: Vector k (Vector k x)
-> Vector k (Vector k x) -> Finite k -> Vector k x
generator Vector k (Vector k x)
rws1 Vector k (Vector k x)
rws2 Finite k
fnt =
                let (G.Vector Vector x
rw1) = Vector k (Vector k x) -> Finite k -> Vector k x
forall (n :: Nat) a. Storable a => Vector n a -> Finite n -> a
index Vector k (Vector k x)
rws1 Finite k
fnt
                    (G.Vector Vector x
rw2) = Vector k (Vector k x) -> Finite k -> Vector k x
forall (n :: Nat) a. Storable a => Vector n a -> Finite n -> a
index Vector k (Vector k x)
rws2 Finite k
fnt
                    i :: Int
i = Finite k -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Finite k
fnt
                 in Vector x -> Vector k x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector k x) -> Vector x -> Vector k x
forall a b. (a -> b) -> a -> b
$ Int -> Vector x -> Vector x
forall a. Storable a => Int -> Vector a -> Vector a
S.take Int
i Vector x
rw1 Vector x -> Vector x -> Vector x
forall a. Storable a => Vector a -> Vector a -> Vector a
S.++ x -> Vector x -> Vector x
forall a. Storable a => a -> Vector a -> Vector a
S.cons (Vector x
diag Vector x -> Int -> x
forall a. Storable a => Vector a -> Int -> a
S.! Int
i) (Int -> Vector x -> Vector x
forall a. Storable a => Int -> Vector a -> Vector a
S.drop (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) Vector x
rw2)

-- | The average of a 'Vector' of elements.
average :: (Numeric x, Fractional x) => Vector n x -> x
{-# INLINE average #-}
average :: Vector n x -> x
average (G.Vector Vector x
v) = Vector x -> x
forall (c :: Type -> Type) e. Container c e => c e -> e
H.sumElements Vector x
v x -> x -> x
forall a. Fractional a => a -> a -> a
/ Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Vector x -> Int
forall a. Storable a => Vector a -> Int
S.length Vector x
v)

-- | The dot product of two numerical 'Vector's.
dotProduct :: Numeric x => Vector n x -> Vector n x -> x
{-# INLINE dotProduct #-}
dotProduct :: Vector n x -> Vector n x -> x
dotProduct Vector n x
v1 Vector n x
v2 = Vector x -> Vector x -> x
forall t. Numeric t => Vector t -> Vector t -> t
H.dot (Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v1) (Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v2)

-- | The determinant of a 'Matrix'.
diagonalMatrix :: forall n x . (KnownNat n, Field x) => Vector n x -> Matrix n n x
{-# INLINE diagonalMatrix #-}
diagonalMatrix :: Vector n x -> Matrix n n x
diagonalMatrix Vector n x
v =
    let n :: Int
n = Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n)
     in Matrix x -> Matrix n n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix n n x) -> Matrix x -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ x -> Vector x -> Int -> Int -> Matrix x
forall t. Storable t => t -> Vector t -> Int -> Int -> Matrix t
H.diagRect x
0 (Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v) Int
n Int
n

-- | The determinant of a 'Matrix'.
takeDiagonal :: (KnownNat n, Field x) => Matrix n n x -> Vector n x
{-# INLINE takeDiagonal #-}
takeDiagonal :: Matrix n n x -> Vector n x
takeDiagonal = Vector x -> Vector n x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector n x)
-> (Matrix n n x -> Vector x) -> Matrix n n x -> Vector n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.takeDiag (Matrix x -> Vector x)
-> (Matrix n n x -> Matrix x) -> Matrix n n x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix

-- | The determinant of a 'Matrix'.
trace :: (KnownNat n, Field x) => Matrix n n x -> x
{-# INLINE trace #-}
trace :: Matrix n n x -> x
trace = Vector x -> x
forall a. (Storable a, Num a) => Vector a -> a
S.sum (Vector x -> x) -> (Matrix n n x -> Vector x) -> Matrix n n x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.takeDiag (Matrix x -> Vector x)
-> (Matrix n n x -> Matrix x) -> Matrix n n x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix

-- | Returns the eigenvalues and eigenvectors 'Matrix'.
eigens :: (KnownNat n, Field x) => Matrix n n x -> (Vector n (Complex Double), Vector n (Vector n (Complex Double)))
{-# INLINE eigens #-}
eigens :: Matrix n n x
-> (Vector n (Complex Double),
    Vector n (Vector n (Complex Double)))
eigens Matrix n n x
mtx =
    let (Vector (Complex Double)
exs,Matrix (Complex Double)
evs) = Matrix x -> (Vector (Complex Double), Matrix (Complex Double))
forall t.
Field t =>
Matrix t -> (Vector (Complex Double), Matrix (Complex Double))
H.eig (Matrix x -> (Vector (Complex Double), Matrix (Complex Double)))
-> Matrix x -> (Vector (Complex Double), Matrix (Complex Double))
forall a b. (a -> b) -> a -> b
$ Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n n x
mtx
     in (Vector (Complex Double) -> Vector n (Complex Double)
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector Vector (Complex Double)
exs, Vector (Vector n (Complex Double))
-> Vector n (Vector n (Complex Double))
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector (Vector n (Complex Double))
 -> Vector n (Vector n (Complex Double)))
-> ([Vector n (Complex Double)]
    -> Vector (Vector n (Complex Double)))
-> [Vector n (Complex Double)]
-> Vector n (Vector n (Complex Double))
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Vector n (Complex Double)] -> Vector (Vector n (Complex Double))
forall a. Storable a => [a] -> Vector a
S.fromList ([Vector n (Complex Double)]
 -> Vector n (Vector n (Complex Double)))
-> [Vector n (Complex Double)]
-> Vector n (Vector n (Complex Double))
forall a b. (a -> b) -> a -> b
$ Vector (Complex Double) -> Vector n (Complex Double)
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector (Complex Double) -> Vector n (Complex Double))
-> [Vector (Complex Double)] -> [Vector n (Complex Double)]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> Matrix (Complex Double) -> [Vector (Complex Double)]
forall t. Element t => Matrix t -> [Vector t]
H.toColumns Matrix (Complex Double)
evs)

-- | Test if the matrix is semi-positive definite.
isSemiPositiveDefinite :: (KnownNat n, Field x) => Matrix n n x -> Bool
{-# INLINE isSemiPositiveDefinite #-}
isSemiPositiveDefinite :: Matrix n n x -> Bool
isSemiPositiveDefinite =
    Vector n Bool -> Bool
forall (n :: Nat). Vector n Bool -> Bool
and (Vector n Bool -> Bool)
-> (Matrix n n x -> Vector n Bool) -> Matrix n n x -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Complex Double -> Bool)
-> Vector n (Complex Double) -> Vector n Bool
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map ((Double
0 Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<=) (Double -> Bool)
-> (Complex Double -> Double) -> Complex Double -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Complex Double -> Double
forall a. Complex a -> a
realPart) (Vector n (Complex Double) -> Vector n Bool)
-> (Matrix n n x -> Vector n (Complex Double))
-> Matrix n n x
-> Vector n Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Vector n (Complex Double), Vector n (Vector n (Complex Double)))
-> Vector n (Complex Double)
forall a b. (a, b) -> a
fst ((Vector n (Complex Double), Vector n (Vector n (Complex Double)))
 -> Vector n (Complex Double))
-> (Matrix n n x
    -> (Vector n (Complex Double),
        Vector n (Vector n (Complex Double))))
-> Matrix n n x
-> Vector n (Complex Double)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix n n x
-> (Vector n (Complex Double),
    Vector n (Vector n (Complex Double)))
forall (n :: Nat) x.
(KnownNat n, Field x) =>
Matrix n n x
-> (Vector n (Complex Double),
    Vector n (Vector n (Complex Double)))
eigens

-- | Returns the inverse, the logarithm of the absolute value of the
-- determinant, and the sign of the determinant of a given matrix.
inverseLogDeterminant :: (KnownNat n, Field x) => Matrix n n x -> (Matrix n n x, x, x)
{-# INLINE inverseLogDeterminant #-}
inverseLogDeterminant :: Matrix n n x -> (Matrix n n x, x, x)
inverseLogDeterminant Matrix n n x
mtx =
    let (Matrix x
imtx,(x
ldet,x
sgn)) = Matrix x -> (Matrix x, (x, x))
forall t. Field t => Matrix t -> (Matrix t, (t, t))
H.invlndet (Matrix x -> (Matrix x, (x, x))) -> Matrix x -> (Matrix x, (x, x))
forall a b. (a -> b) -> a -> b
$ Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n n x
mtx
     in (Matrix x -> Matrix n n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix Matrix x
imtx, x
ldet, x
sgn)

-- | The determinant of a 'Matrix'.
determinant :: (KnownNat n, Field x) => Matrix n n x -> x
{-# INLINE determinant #-}
determinant :: Matrix n n x -> x
determinant = Matrix x -> x
forall t. Field t => Matrix t -> t
H.det (Matrix x -> x) -> (Matrix n n x -> Matrix x) -> Matrix n n x -> x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix

-- | Transpose a 'Matrix'.
transpose
    :: forall m n x . (KnownNat m, KnownNat n, Numeric x)
    => Matrix m n x
    -> Matrix n m x
{-# INLINE transpose #-}
transpose :: Matrix m n x -> Matrix n m x
transpose (G.Matrix Vector Vector (m * n) x
mtx) =
    Vector Vector (n * m) x -> Matrix n m x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (n * m) x -> Matrix n m x)
-> Vector Vector (n * m) x -> Matrix n m x
forall a b. (a -> b) -> a -> b
$ (Vector x -> Vector x)
-> Vector Vector (m * n) x -> Vector Vector (m * n) x
forall a b (n :: Nat).
(Vector a -> Vector b) -> Vector n a -> Vector n b
withVectorUnsafe (Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.flatten (Matrix x -> Vector x)
-> (Vector x -> Matrix x) -> Vector x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Matrix x
forall m mt. Transposable m mt => m -> mt
H.tr (Matrix x -> Matrix x)
-> (Vector x -> Matrix x) -> Vector x -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Vector x -> Matrix x
forall t. Storable t => Int -> Vector t -> Matrix t
H.reshape (Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n))) Vector Vector (m * n) x
mtx

-- | Diagonally concatenate two matrices, padding the gaps with zeroes.
diagonalConcat
    :: (KnownNat n, KnownNat m, KnownNat o, KnownNat p, Numeric x)
    => Matrix n m x -> Matrix o p x -> Matrix (n+o) (m+p) x
{-# INLINE diagonalConcat #-}
diagonalConcat :: Matrix n m x -> Matrix o p x -> Matrix (n + o) (m + p) x
diagonalConcat Matrix n m x
mtx10 Matrix o p x
mtx20 =
    let mtx1 :: Matrix x
mtx1 = Matrix n m x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n m x
mtx10
        mtx2 :: Matrix x
mtx2 = Matrix o p x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix o p x
mtx20
     in Matrix x -> Matrix (n + o) (m + p) x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix (n + o) (m + p) x)
-> Matrix x -> Matrix (n + o) (m + p) x
forall a b. (a -> b) -> a -> b
$ [Matrix x] -> Matrix x
forall t. (Element t, Num t) => [Matrix t] -> Matrix t
H.diagBlock [Matrix x
mtx1,Matrix x
mtx2]

-- | Diagonally concatenate two matrices, padding the gaps with zeroes.
horizontalConcat
    :: (KnownNat n, KnownNat m, KnownNat o, Numeric x)
    => Matrix n m x -> Matrix n o x -> Matrix n (m+o) x
{-# INLINE horizontalConcat #-}
horizontalConcat :: Matrix n m x -> Matrix n o x -> Matrix n (m + o) x
horizontalConcat Matrix n m x
mtx10 Matrix n o x
mtx20 =
    let mtx1 :: Matrix x
mtx1 = Matrix n m x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n m x
mtx10
        mtx2 :: Matrix x
mtx2 = Matrix n o x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n o x
mtx20
     in Matrix x -> Matrix n (m + o) x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix n (m + o) x) -> Matrix x -> Matrix n (m + o) x
forall a b. (a -> b) -> a -> b
$ Matrix x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Element t => Matrix t -> Matrix t -> Matrix t
H.||| Matrix x
mtx2

-- | Diagonally concatenate two matrices, padding the gaps with zeroes.
verticalConcat
    :: (KnownNat n, KnownNat m, KnownNat o, Numeric x)
    => Matrix n o x -> Matrix m o x -> Matrix (n+m) o x
{-# INLINE verticalConcat #-}
verticalConcat :: Matrix n o x -> Matrix m o x -> Matrix (n + m) o x
verticalConcat Matrix n o x
mtx10 Matrix m o x
mtx20 =
    let mtx1 :: Matrix x
mtx1 = Matrix n o x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n o x
mtx10
        mtx2 :: Matrix x
mtx2 = Matrix m o x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix m o x
mtx20
     in Matrix x -> Matrix (n + m) o x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix (n + m) o x) -> Matrix x -> Matrix (n + m) o x
forall a b. (a -> b) -> a -> b
$ Matrix x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Element t => Matrix t -> Matrix t -> Matrix t
H.=== Matrix x
mtx2

-- | Invert a 'Matrix'.
inverse :: forall n x . (KnownNat n, Field x) => Matrix n n x -> Matrix n n x
{-# INLINE inverse #-}
inverse :: Matrix n n x -> Matrix n n x
inverse (G.Matrix Vector Vector (n * n) x
mtx) =
    Vector Vector (n * n) x -> Matrix n n x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (n * n) x -> Matrix n n x)
-> Vector Vector (n * n) x -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ (Vector x -> Vector x)
-> Vector Vector (n * n) x -> Vector Vector (n * n) x
forall a b (n :: Nat).
(Vector a -> Vector b) -> Vector n a -> Vector n b
withVectorUnsafe (Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.flatten (Matrix x -> Vector x)
-> (Vector x -> Matrix x) -> Vector x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Matrix x
forall t. Field t => Matrix t -> Matrix t
H.inv (Matrix x -> Matrix x)
-> (Vector x -> Matrix x) -> Vector x -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Vector x -> Matrix x
forall t. Storable t => Int -> Vector t -> Matrix t
H.reshape (Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n))) Vector Vector (n * n) x
mtx

-- | Pseudo-Invert a 'Matrix'.
pseudoInverse :: forall n x . (KnownNat n, Field x) => Matrix n n x -> Matrix n n x
{-# INLINE pseudoInverse #-}
pseudoInverse :: Matrix n n x -> Matrix n n x
pseudoInverse (G.Matrix Vector Vector (n * n) x
mtx) =
    Vector Vector (n * n) x -> Matrix n n x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (n * n) x -> Matrix n n x)
-> Vector Vector (n * n) x -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ (Vector x -> Vector x)
-> Vector Vector (n * n) x -> Vector Vector (n * n) x
forall a b (n :: Nat).
(Vector a -> Vector b) -> Vector n a -> Vector n b
withVectorUnsafe (Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.flatten (Matrix x -> Vector x)
-> (Vector x -> Matrix x) -> Vector x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Matrix x
forall t. Field t => Matrix t -> Matrix t
H.pinv (Matrix x -> Matrix x)
-> (Vector x -> Matrix x) -> Vector x -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Vector x -> Matrix x
forall t. Storable t => Int -> Vector t -> Matrix t
H.reshape (Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n))) Vector Vector (n * n) x
mtx

-- | Square root of a 'Matrix'.
matrixRoot :: forall n x . (KnownNat n, Field x) => Matrix n n x -> Matrix n n x
{-# INLINE matrixRoot #-}
matrixRoot :: Matrix n n x -> Matrix n n x
matrixRoot (G.Matrix Vector Vector (n * n) x
mtx) =
    Vector Vector (n * n) x -> Matrix n n x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector Vector (n * n) x -> Matrix n n x)
-> Vector Vector (n * n) x -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ (Vector x -> Vector x)
-> Vector Vector (n * n) x -> Vector Vector (n * n) x
forall a b (n :: Nat).
(Vector a -> Vector b) -> Vector n a -> Vector n b
withVectorUnsafe (Matrix x -> Vector x
forall t. Element t => Matrix t -> Vector t
H.flatten (Matrix x -> Vector x)
-> (Vector x -> Matrix x) -> Vector x -> Vector x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Matrix x
forall t. Field t => Matrix t -> Matrix t
H.sqrtm (Matrix x -> Matrix x)
-> (Vector x -> Matrix x) -> Vector x -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Vector x -> Matrix x
forall t. Storable t => Int -> Vector t -> Matrix t
H.reshape (Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n))) Vector Vector (n * n) x
mtx

-- | The outer product of two 'Vector's.
outerProduct :: (KnownNat m, KnownNat n, Numeric x) => Vector m x -> Vector n x -> Matrix m n x
{-# INLINE outerProduct #-}
outerProduct :: Vector m x -> Vector n x -> Matrix m n x
outerProduct Vector m x
v1 Vector n x
v2 =
    Matrix x -> Matrix m n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix m n x) -> Matrix x -> Matrix m n x
forall a b. (a -> b) -> a -> b
$ Vector x -> Vector x -> Matrix x
forall t. Product t => Vector t -> Vector t -> Matrix t
H.outer (Vector m x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector m x
v1) (Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v2)

-- | The summed outer product of two lists of 'Vector's.
sumOuterProduct :: (KnownNat m, KnownNat n, Fractional x, Numeric x) => [(Vector m x,Vector n x)] -> Matrix m n x
{-# INLINE sumOuterProduct #-}
sumOuterProduct :: [(Vector m x, Vector n x)] -> Matrix m n x
sumOuterProduct [(Vector m x, Vector n x)]
v12s =
    let ([Vector m x]
v1s,[Vector n x]
v2s) = [(Vector m x, Vector n x)] -> ([Vector m x], [Vector n x])
forall a b. [(a, b)] -> ([a], [b])
L.unzip [(Vector m x, Vector n x)]
v12s
        mtx1 :: Matrix x
mtx1 = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromColumns ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector m x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector m x -> Vector x) -> [Vector m x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector m x]
v1s
        mtx2 :: Matrix x
mtx2 = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromRows ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector n x -> Vector x) -> [Vector n x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector n x]
v2s
     in Matrix x -> Matrix m n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
H.<> Matrix x
mtx2)

-- | The average outer product of two lists of 'Vector's.
averageOuterProduct :: (KnownNat m, KnownNat n, Fractional x, Numeric x) => [(Vector m x,Vector n x)] -> Matrix m n x
{-# INLINE averageOuterProduct #-}
averageOuterProduct :: [(Vector m x, Vector n x)] -> Matrix m n x
averageOuterProduct [(Vector m x, Vector n x)]
v12s =
    let ([Vector m x]
v1s,[Vector n x]
v2s) = [(Vector m x, Vector n x)] -> ([Vector m x], [Vector n x])
forall a b. [(a, b)] -> ([a], [b])
L.unzip [(Vector m x, Vector n x)]
v12s
        mtx1 :: Matrix x
mtx1 = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromColumns ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector m x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector m x -> Vector x) -> [Vector m x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector m x]
v1s
        (Int
_,Int
n) = Matrix x -> IndexOf Matrix
forall (c :: Type -> Type) t. Container c t => c t -> IndexOf c
H.size Matrix x
mtx1
        mtx2 :: Matrix x
mtx2 = x -> Matrix x -> Matrix x
forall t (c :: Type -> Type). Linear t c => t -> c t -> c t
H.scale (x
1x -> x -> x
forall a. Fractional a => a -> a -> a
/Int -> x
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n) (Matrix x -> Matrix x)
-> ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromRows ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector n x -> Vector x) -> [Vector n x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector n x]
v2s
     in Matrix x -> Matrix m n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
H.<> Matrix x
mtx2)

-- | The average outer product of two lists of 'Vector's.
weightedAverageOuterProduct
    :: ( KnownNat m, KnownNat n, Fractional x, Numeric x )
    => [(x,Vector m x,Vector n x)]
    -> Matrix m n x
{-# INLINE weightedAverageOuterProduct #-}
weightedAverageOuterProduct :: [(x, Vector m x, Vector n x)] -> Matrix m n x
weightedAverageOuterProduct [(x, Vector m x, Vector n x)]
wv12s =
    let ([x]
ws,[Vector m x]
v1s,[Vector n x]
v2s) = [(x, Vector m x, Vector n x)] -> ([x], [Vector m x], [Vector n x])
forall a b c. [(a, b, c)] -> ([a], [b], [c])
L.unzip3 [(x, Vector m x, Vector n x)]
wv12s
        v1s' :: [Vector x]
v1s' = (x -> Vector x -> Vector x) -> [x] -> [Vector x] -> [Vector x]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
L.zipWith x -> Vector x -> Vector x
forall t (c :: Type -> Type). Linear t c => t -> c t -> c t
H.scale [x]
ws ([Vector x] -> [Vector x]) -> [Vector x] -> [Vector x]
forall a b. (a -> b) -> a -> b
$ Vector m x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector m x -> Vector x) -> [Vector m x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector m x]
v1s
        mtx1 :: Matrix x
mtx1 = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromColumns [Vector x]
v1s'
        mtx2 :: Matrix x
mtx2 = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromRows ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector n x -> Vector x) -> [Vector n x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector n x]
v2s
     in Matrix x -> Matrix m n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
H.<> Matrix x
mtx2)

-- | The identity 'Matrix'.
matrixIdentity :: forall n x . (KnownNat n, Numeric x, Num x) => Matrix n n x
{-# INLINE matrixIdentity #-}
matrixIdentity :: Matrix n n x
matrixIdentity =
    Matrix x -> Matrix n n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix n n x)
-> (Int -> Matrix x) -> Int -> Matrix n n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Matrix x
forall a. (Num a, Element a) => Int -> Matrix a
H.ident (Int -> Matrix n n x) -> Int -> Matrix n n x
forall a b. (a -> b) -> a -> b
$ Proxy n -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt (Proxy n
forall k (t :: k). Proxy t
Proxy :: Proxy n)

-- | The dot products of one vector with a list of vectors.
dotMap :: (KnownNat n, Numeric x) => Vector n x -> [Vector n x] -> [x]
{-# INLINE dotMap #-}
dotMap :: Vector n x -> [Vector n x] -> [x]
dotMap Vector n x
v [Vector n x]
vs =
    let mtx' :: Matrix x
mtx' = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromRows ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector n x -> Vector x) -> [Vector n x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector n x]
vs
     in Vector x -> [x]
forall a. Storable a => Vector a -> [a]
H.toList (Vector x -> [x]) -> Vector x -> [x]
forall a b. (a -> b) -> a -> b
$ Matrix x
mtx' Matrix x -> Vector x -> Vector x
forall t. Numeric t => Matrix t -> Vector t -> Vector t
H.#> Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v
--     in if S.null w
--           then replicate 0
--           else fmap G.Vector . H.toColumns $ toHMatrix mtx H.<> mtx'

-- | Map a linear transformation over a list of 'Vector's.
matrixMap :: (KnownNat m, KnownNat n, Numeric x)
                     => Matrix m n x -> [Vector n x] -> [Vector m x]
{-# INLINE matrixMap #-}
matrixMap :: Matrix m n x -> [Vector n x] -> [Vector m x]
matrixMap Matrix m n x
mtx [Vector n x]
vs =
    let mtx' :: Matrix x
mtx' = [Vector x] -> Matrix x
forall t. Element t => [Vector t] -> Matrix t
H.fromColumns ([Vector x] -> Matrix x) -> [Vector x] -> Matrix x
forall a b. (a -> b) -> a -> b
$ Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector n x -> Vector x) -> [Vector n x] -> [Vector x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector n x]
vs
     in (Vector x -> Vector m x) -> [Vector x] -> [Vector m x]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
fmap Vector x -> Vector m x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector ([Vector x] -> [Vector m x])
-> (Matrix x -> [Vector x]) -> Matrix x -> [Vector m x]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> [Vector x]
forall t. Element t => Matrix t -> [Vector t]
H.toColumns (Matrix x -> [Vector m x]) -> Matrix x -> [Vector m x]
forall a b. (a -> b) -> a -> b
$ Matrix m n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix m n x
mtx Matrix x -> Matrix x -> Matrix x
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
H.<> Matrix x
mtx'
--     in if S.null w
--           then replicate 0
--           else fmap G.Vector . H.toColumns $ toHMatrix mtx H.<> mtx'


-- | Apply a linear transformation to a 'Vector'.
matrixVectorMultiply :: (KnownNat m, KnownNat n, Numeric x)
                     => Matrix m n x -> Vector n x -> Vector m x
{-# INLINE matrixVectorMultiply #-}
matrixVectorMultiply :: Matrix m n x -> Vector n x -> Vector m x
matrixVectorMultiply Matrix m n x
mtx Vector n x
v =
    Vector x -> Vector m x
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector x -> Vector m x) -> Vector x -> Vector m x
forall a b. (a -> b) -> a -> b
$ Matrix m n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix m n x
mtx Matrix x -> Vector x -> Vector x
forall t. Numeric t => Matrix t -> Vector t -> Vector t
H.#> Vector n x -> Vector x
forall (n :: Nat) a. Vector n a -> Vector a
fromSized Vector n x
v
--    let w = toHMatrix mtx H.#> fromSized v
--     in if S.null w
--           then replicate 0
--           else G.Vector w

-- | Multiply a 'Matrix' with a second 'Matrix'.
matrixMatrixMultiply
    :: (KnownNat m, KnownNat n, KnownNat o, Numeric x)
    => Matrix m n x
    -> Matrix n o x
    -> Matrix m o x
{-# INLINE matrixMatrixMultiply #-}
matrixMatrixMultiply :: Matrix m n x -> Matrix n o x -> Matrix m o x
matrixMatrixMultiply Matrix m n x
mtx1 Matrix n o x
mtx2 = Matrix x -> Matrix m o x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix m o x) -> Matrix x -> Matrix m o x
forall a b. (a -> b) -> a -> b
$ Matrix m n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix m n x
mtx1 Matrix x -> Matrix x -> Matrix x
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
H.<> Matrix n o x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix Matrix n o x
mtx2

-- | Pretty print the values of a 'Matrix' (for extremely simple values of pretty).
prettyPrintMatrix :: (KnownNat m, KnownNat n, Numeric a, Show a) => Matrix m n a -> IO ()
prettyPrintMatrix :: Matrix m n a -> IO ()
prettyPrintMatrix = Matrix a -> IO ()
forall a. Show a => a -> IO ()
print (Matrix a -> IO ())
-> (Matrix m n a -> Matrix a) -> Matrix m n a -> IO ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix m n a -> Matrix a
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix

-- | The Mean Squared difference between two vectors.
meanSquaredError
    :: KnownNat k
    => Vector k Double
    -> Vector k Double
    -> Double
{-# INLINE meanSquaredError #-}
meanSquaredError :: Vector k Double -> Vector k Double -> Double
meanSquaredError Vector k Double
ys Vector k Double
yhts = Vector k Double -> Double
forall x (n :: Nat). (Numeric x, Fractional x) => Vector n x -> x
average (Vector k Double -> Double) -> Vector k Double -> Double
forall a b. (a -> b) -> a -> b
$ (Double -> Double) -> Vector k Double -> Vector k Double
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map Double -> Double
forall x. Floating x => x -> x
square (Vector k Double
ys Vector k Double -> Vector k Double -> Vector k Double
forall a. Num a => a -> a -> a
- Vector k Double
yhts)

-- | L2 length of a vector.
l2Norm
    :: KnownNat k
    => Vector k Double
    -> Double
{-# INLINE l2Norm #-}
l2Norm :: Vector k Double -> Double
l2Norm (G.Vector Vector Double
xs) = Vector Double -> Double
forall a. Normed a => a -> Double
H.norm_2 Vector Double
xs

-- | Computes the coefficient of determintation for the given outputs and model
-- predictions.
rSquared
    :: KnownNat k
    => Vector k Double -- ^ Dependent variable observations
    -> Vector k Double -- ^ Predicted Values
    -> Double -- ^ R-squared
{-# INLINE rSquared #-}
rSquared :: Vector k Double -> Vector k Double -> Double
rSquared Vector k Double
ys Vector k Double
yhts =
    let ybr :: Double
ybr = Vector k Double -> Double
forall x (n :: Nat). (Numeric x, Fractional x) => Vector n x -> x
average Vector k Double
ys
        ssres :: Double
ssres = Vector k Double -> Double
forall a (n :: Nat). (Storable a, Num a) => Vector n a -> a
sum (Vector k Double -> Double) -> Vector k Double -> Double
forall a b. (a -> b) -> a -> b
$ (Double -> Double) -> Vector k Double -> Vector k Double
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map Double -> Double
forall x. Floating x => x -> x
square (Vector k Double
ys Vector k Double -> Vector k Double -> Vector k Double
forall a. Num a => a -> a -> a
- Vector k Double
yhts)
        sstot :: Double
sstot = Vector k Double -> Double
forall a (n :: Nat). (Storable a, Num a) => Vector n a -> a
sum (Vector k Double -> Double) -> Vector k Double -> Double
forall a b. (a -> b) -> a -> b
$ (Double -> Double) -> Vector k Double -> Vector k Double
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map (Double -> Double
forall x. Floating x => x -> x
square (Double -> Double) -> (Double -> Double) -> Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> Double -> Double
forall a. Num a => a -> a -> a
subtract Double
ybr) Vector k Double
ys
     in Double
1 Double -> Double -> Double
forall a. Num a => a -> a -> a
- (Double
ssresDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
sstot)

-- | Solves the linear least squares problem.
linearLeastSquares
    :: KnownNat l
    => [Vector l Double] -- ^ Independent variable observations
    -> [Double] -- ^ Dependent variable observations
    -> Vector l Double -- ^ Parameter estimates
{-# INLINE linearLeastSquares #-}
linearLeastSquares :: [Vector l Double] -> [Double] -> Vector l Double
linearLeastSquares [Vector l Double]
as [Double]
xs =
    Vector Double -> Vector l Double
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector Double -> Vector l Double)
-> Vector Double -> Vector l Double
forall a b. (a -> b) -> a -> b
$ [Vector Double] -> Matrix Double
forall t. Element t => [Vector t] -> Matrix t
H.fromRows (Vector l Double -> Vector Double
forall (n :: Nat) a. Vector n a -> Vector a
fromSized (Vector l Double -> Vector Double)
-> [Vector l Double] -> [Vector Double]
forall (f :: Type -> Type) a b. Functor f => (a -> b) -> f a -> f b
<$> [Vector l Double]
as) Matrix Double -> Vector Double -> Vector Double
forall (c :: Type -> Type) t.
(LSDiv c, Field t) =>
Matrix t -> c t -> c t
H.<\> [Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
S.fromList [Double]
xs


unsafeCholesky
    :: (KnownNat n, Field x, Storable x)
    => Matrix n n x
    -> Matrix n n x
unsafeCholesky :: Matrix n n x -> Matrix n n x
unsafeCholesky =
    Matrix n n x -> Matrix n n x
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Numeric x) =>
Matrix m n x -> Matrix n m x
transpose (Matrix n n x -> Matrix n n x)
-> (Matrix n n x -> Matrix n n x) -> Matrix n n x -> Matrix n n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Matrix n n x
forall x (m :: Nat) (n :: Nat).
Numeric x =>
Matrix x -> Matrix m n x
fromHMatrix (Matrix x -> Matrix n n x)
-> (Matrix n n x -> Matrix x) -> Matrix n n x -> Matrix n n x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Herm x -> Matrix x
forall t. Field t => Herm t -> Matrix t
H.chol (Herm x -> Matrix x)
-> (Matrix n n x -> Herm x) -> Matrix n n x -> Matrix x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix x -> Herm x
forall t. Matrix t -> Herm t
H.trustSym (Matrix x -> Herm x)
-> (Matrix n n x -> Matrix x) -> Matrix n n x -> Herm x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix n n x -> Matrix x
forall (m :: Nat) (n :: Nat) x.
(KnownNat n, KnownNat m, Element x, Storable x) =>
Matrix m n x -> Matrix x
toHMatrix


--- Convolutions ---


-- | 2d cross-correlation of a kernel over a matrix of values.
crossCorrelate2d
    :: forall nk rdkr rdkc mr mc md x
    . ( KnownNat rdkr, KnownNat rdkc, KnownNat md, KnownNat mr, KnownNat mc
      , KnownNat nk, Numeric x, Storable x )
      => Proxy rdkr -- ^ Number of Kernel rows
      -> Proxy rdkc -- ^ Number of Kernel columns
      -> Proxy mr -- ^ Number of Matrix/Image rows
      -> Proxy mc -- ^ Number of Kernel/Image columns
      -> Matrix nk (md*(2*rdkr+1)*(2*rdkc+1)) x -- ^ Kernels (nk is their number)
      -> Matrix md (mr*mc) x -- ^ Image (md is the depth)
      -> Matrix nk (mr*mc) x -- ^ Cross-correlated image
{-# INLINE crossCorrelate2d #-}
crossCorrelate2d :: Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md (mr * mc) x
-> Matrix nk (mr * mc) x
crossCorrelate2d Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krns (G.Matrix Vector Vector (md * (mr * mc)) x
v) =
    let pmd :: Proxy md
pmd = Proxy md
forall k (t :: k). Proxy t
Proxy :: Proxy md
        mtx :: Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
mtx = Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
forall (rdkr :: Nat) (rdkc :: Nat) (md :: Nat) (mr :: Nat)
       (mc :: Nat) x.
(KnownNat rdkr, KnownNat rdkc, KnownNat mc, KnownNat md,
 KnownNat mr, Num x, Storable x) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
im2col Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
pmd Proxy mr
pmr Proxy mc
pmc Vector Vector (md * (mr * mc)) x
Vector ((md * mr) * mc) x
v
     in Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
-> Matrix nk (mr * mc) x
forall (m :: Nat) (n :: Nat) (o :: Nat) x.
(KnownNat m, KnownNat n, KnownNat o, Numeric x) =>
Matrix m n x -> Matrix n o x -> Matrix m o x
matrixMatrixMultiply Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krns Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
mtx

-- | The transpose of a convolutional kernel.
kernelTranspose
    :: (KnownNat nk, KnownNat md, KnownNat rdkr, KnownNat rdkc, Numeric x, Storable x)
    => Proxy nk
    -> Proxy md
    -> Proxy rdkr
    -> Proxy rdkc
    -> Matrix nk (md*(2*rdkr+1)*(2*rdkc+1)) x -- ^ Kernels (nk is their number)
    -> Matrix md (nk*(2*rdkr+1)*(2*rdkc+1)) x -- ^ Kernels (nk is their number)
{-# INLINE kernelTranspose #-}
kernelTranspose :: Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
kernelTranspose Proxy nk
pnk Proxy md
pmd Proxy rdkr
prdkr Proxy rdkc
prdkc (G.Matrix Vector Vector (nk * ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1))) x
kv) = Vector Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector
   Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
 -> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x)
-> (Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
    -> Vector
         Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x)
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector Vector (nk * ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1))) x
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
-> Vector
     Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall a (m :: Nat) (n :: Nat).
Storable a =>
Vector m a -> Vector n Int -> Vector n a
backpermute Vector Vector (nk * ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1))) x
kv (Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
 -> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x)
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall a b. (a -> b) -> a -> b
$ Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
forall (nk :: Nat) (md :: Nat) (rdkr :: Nat) (rdkc :: Nat).
(KnownNat nk, KnownNat md, KnownNat rdkr, KnownNat rdkc) =>
Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
kernelTransposeIndices Proxy nk
pnk Proxy md
pmd Proxy rdkr
prdkr Proxy rdkc
prdkc

-- | 2d convolution of a kernel over a matrix of values. This is the adjoint of crossCorrelate2d.
convolve2d
    :: forall nk rdkr rdkc md mr mc x
    . ( KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc
      , KnownNat md, KnownNat nk, Numeric x, Storable x )
      => Proxy rdkr -- ^ Number of Kernel rows
      -> Proxy rdkc -- ^ Number of Kernel columns
      -> Proxy mr -- ^ Number of Matrix/Image rows
      -> Proxy mc -- ^ Number of Kernel/Image columns
      -> Matrix nk (md*(2*rdkr+1)*(2*rdkc+1)) x -- ^ Kernels (nk is their number)
      -> Matrix nk (mr*mc) x -- ^ Dual image (nk is its depth)
      -> Matrix md (mr*mc) x -- ^ Convolved image
{-# INLINE convolve2d #-}
convolve2d :: Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix nk (mr * mc) x
-> Matrix md (mr * mc) x
convolve2d Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krn Matrix nk (mr * mc) x
mtxs =
    let pnk :: Proxy nk
pnk = Proxy nk
forall k (t :: k). Proxy t
Proxy :: Proxy nk
        pmd :: Proxy md
pmd = Proxy md
forall k (t :: k). Proxy t
Proxy :: Proxy md
        krn' :: Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krn' = Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall (nk :: Nat) (md :: Nat) (rdkr :: Nat) (rdkc :: Nat) x.
(KnownNat nk, KnownNat md, KnownNat rdkr, KnownNat rdkc, Numeric x,
 Storable x) =>
Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
kernelTranspose Proxy nk
pnk Proxy md
pmd Proxy rdkr
prdkr Proxy rdkc
prdkc Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krn
     in Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix nk (mr * mc) x
-> Matrix md (mr * mc) x
forall (nk :: Nat) (rdkr :: Nat) (rdkc :: Nat) (mr :: Nat)
       (mc :: Nat) (md :: Nat) x.
(KnownNat rdkr, KnownNat rdkc, KnownNat md, KnownNat mr,
 KnownNat mc, KnownNat nk, Numeric x, Storable x) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix md (mr * mc) x
-> Matrix nk (mr * mc) x
crossCorrelate2d Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Matrix md ((nk * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
krn' Matrix nk (mr * mc) x
mtxs

-- | The outer product of an image and a dual image to produce a convolutional kernel.
kernelOuterProduct
    :: forall nk rdkr rdkc md mr mc x
    . ( KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc
      , KnownNat md, KnownNat nk, Numeric x, Storable x )
      => Proxy rdkr -- ^ Number of Kernel rows
      -> Proxy rdkc -- ^ Number of Kernel columns
      -> Proxy mr -- ^ Number of Matrix/Image rows
      -> Proxy mc -- ^ Number of Kernel/Image columns
      -> Matrix nk (mr*mc) x -- ^ Dual image (nk is its depth)
      -> Matrix md (mr*mc) x -- ^ Image (md is the depth)
      -> Matrix nk (md*(2*rdkr+1)*(2*rdkc+1)) x -- ^ Kernels
{-# INLINE kernelOuterProduct #-}
kernelOuterProduct :: Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Matrix nk (mr * mc) x
-> Matrix md (mr * mc) x
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
kernelOuterProduct Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Matrix nk (mr * mc) x
omtx (G.Matrix Vector Vector (md * (mr * mc)) x
v) =
    let pmd :: Proxy md
pmd = Proxy md
forall k (t :: k). Proxy t
Proxy :: Proxy md
        imtx :: Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
imtx = Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
forall (rdkr :: Nat) (rdkc :: Nat) (md :: Nat) (mr :: Nat)
       (mc :: Nat) x.
(KnownNat rdkr, KnownNat rdkc, KnownNat mc, KnownNat md,
 KnownNat mr, Num x, Storable x) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
im2col Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
pmd Proxy mr
pmr Proxy mc
pmc Vector Vector (md * (mr * mc)) x
Vector ((md * mr) * mc) x
v
     in Matrix nk (mr * mc) x
-> Matrix (mr * mc) ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall (m :: Nat) (n :: Nat) (o :: Nat) x.
(KnownNat m, KnownNat n, KnownNat o, Numeric x) =>
Matrix m n x -> Matrix n o x -> Matrix m o x
matrixMatrixMultiply Matrix nk (mr * mc) x
omtx (Matrix (mr * mc) ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
 -> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x)
-> Matrix (mr * mc) ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
-> Matrix nk ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall a b. (a -> b) -> a -> b
$ Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
-> Matrix (mr * mc) ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) x
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Numeric x) =>
Matrix m n x -> Matrix n m x
transpose Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
imtx


--- Internal ---


toTriangularIndex :: (Int,Int) -> Int
toTriangularIndex :: (Int, Int) -> Int
toTriangularIndex (Int
i,Int
j)
    | Int
i Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
j = Int -> Int
triangularNumber Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
j
    | Bool
otherwise = (Int, Int) -> Int
toTriangularIndex (Int
j,Int
i)

to2Index :: Int -> Int -> (Int,Int)
to2Index :: Int -> Int -> (Int, Int)
to2Index Int
nj Int
ij = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
ij Int
nj

to3Index :: Int -> Int -> Int -> (Int,Int,Int)
{-# INLINE to3Index #-}
to3Index :: Int -> Int -> Int -> (Int, Int, Int)
to3Index Int
nj Int
nk Int
ijk =
    let nj' :: Int
nj' = Int
njInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk
        (Int
i,Int
jk) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
ijk Int
nj'
        (Int
j,Int
k) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
jk Int
nk
     in (Int
i,Int
j,Int
k)

from3Index :: Int -> Int -> (Int,Int,Int) -> Int
{-# INLINE from3Index #-}
from3Index :: Int -> Int -> (Int, Int, Int) -> Int
from3Index Int
nj Int
nk (Int
i,Int
j,Int
k) =
    let nj' :: Int
nj' = Int
njInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk
     in Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nj' Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
k

to4Index :: Int -> Int -> Int -> Int -> (Int,Int,Int,Int)
{-# INLINE to4Index #-}
to4Index :: Int -> Int -> Int -> Int -> (Int, Int, Int, Int)
to4Index Int
nj Int
nk Int
nl Int
ijkl =
    let nk' :: Int
nk' = Int
nlInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk
        nj' :: Int
nj' = Int
njInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk'
        (Int
i,Int
jkl) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
ijkl Int
nj'
        (Int
j,Int
kl) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
jkl Int
nk'
        (Int
k,Int
l) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
kl Int
nl
     in (Int
i,Int
j,Int
k,Int
l)

from4Index :: Int -> Int -> Int -> (Int,Int,Int,Int) -> Int
{-# INLINE from4Index #-}
from4Index :: Int -> Int -> Int -> (Int, Int, Int, Int) -> Int
from4Index Int
nj Int
nk Int
nl (Int
i,Int
j,Int
k,Int
l) =
    let nk' :: Int
nk' = Int
nlInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk
        nj' :: Int
nj' = Int
njInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk'
     in Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nj' Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nk' Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
kInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
nl Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
l

kernelTransposeIndices
    :: (KnownNat nk, KnownNat md, KnownNat rdkr, KnownNat rdkc)
    => Proxy nk
    -> Proxy md
    -> Proxy rdkr
    -> Proxy rdkc
    -> Vector (nk*md*(2*rdkr+1)*(2*rdkc+1)) Int
{-# INLINE kernelTransposeIndices #-}
kernelTransposeIndices :: Proxy nk
-> Proxy md
-> Proxy rdkr
-> Proxy rdkc
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
kernelTransposeIndices Proxy nk
pnk Proxy md
pmd Proxy rdkr
prdkr Proxy rdkc
prdkc =
    let nkrn :: Int
nkrn = Proxy nk -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy nk
pnk
        md :: Int
md = Proxy md -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy md
pmd
        rdkr :: Int
rdkr = Proxy rdkr -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkr
prdkr
        rdkc :: Int
rdkc = Proxy rdkc -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkc
prdkc
        dmkr :: Int
dmkr = Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkrInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
        dmkc :: Int
dmkc = Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkcInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
        nl :: Int
nl = Int
dmkc
        nk :: Int
nk = Int
dmkr
        nj :: Int
nj = Int
nkrn
        nl' :: Int
nl' = Int
dmkc
        nk' :: Int
nk' = Int
dmkr
        nj' :: Int
nj' = Int
md
        reIndex :: Int -> Int
reIndex Int
idx =
            let (Int
i,Int
j,Int
k,Int
l) = Int -> Int -> Int -> Int -> (Int, Int, Int, Int)
to4Index Int
nj Int
nk Int
nl Int
idx
             in Int -> Int -> Int -> (Int, Int, Int, Int) -> Int
from4Index Int
nj' Int
nk' Int
nl' (Int
j,Int
i,Int
nkInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
k,Int
nlInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
l)
     in (Finite (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) -> Int)
-> Vector (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) Int
forall (n :: Nat) a.
(KnownNat n, Storable a) =>
(Finite n -> a) -> Vector n a
generate (Int -> Int
reIndex (Int -> Int)
-> (Finite (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1))
    -> Int)
-> Finite (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1))
-> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Finite (((nk * md) * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral)

im2colIndices
    :: forall rdkr rdkc mr mc md
     . (KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc, KnownNat md)
    => Proxy rdkr
    -> Proxy rdkc
    -> Proxy md
    -> Proxy mr
    -> Proxy mc
    -> Vector (((2*rdkr+1)*(2*rdkc+1)*md)*(mr*mc)) Int
{-# INLINE im2colIndices #-}
im2colIndices :: Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
im2colIndices Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
_ Proxy mr
pmr Proxy mc
pmc =
    let rdkr :: Int
rdkr = Proxy rdkr -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkr
prdkr
        rdkc :: Int
rdkc = Proxy rdkc -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkc
prdkc
        nj :: Int
nj = (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkr Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
1)
        nk :: Int
nk = (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkc Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
1)
        reWindow :: Int -> Vector (mr * mc) Int
reWindow Int
idx =
            let (Int
i,Int
j,Int
k) = Int -> Int -> Int -> (Int, Int, Int)
to3Index Int
nj Int
nk Int
idx
             in Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Int
-> Int
-> Int
-> Vector (mr * mc) Int
forall (rdkr :: Nat) (rdkc :: Nat) (mr :: Nat) (mc :: Nat).
(KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Int
-> Int
-> Int
-> Vector (mr * mc) Int
windowIndices Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Int
i Int
j Int
k
          in ((Int -> Vector (mr * mc) Int)
-> Vector ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) Int
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
forall a b (m :: Nat) (n :: Nat).
(Storable a, Storable b) =>
(a -> Vector m b) -> Vector n a -> Vector (n * m) b
concatMap Int -> Vector (mr * mc) Int
reWindow :: Vector ((2*rdkr+1)*(2*rdkc+1)*md) Int -> Vector (((2*rdkr+1)*(2*rdkc+1)*md)*(mr*mc)) Int) (Vector ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) Int
 -> Vector
      (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int)
-> Vector ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) Int
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
forall a b. (a -> b) -> a -> b
$ (Finite ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) -> Int)
-> Vector ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) Int
forall (n :: Nat) a.
(KnownNat n, Storable a) =>
(Finite n -> a) -> Vector n a
generate Finite ((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) -> Int
forall (n :: Nat). KnownNat n => Finite n -> Int
finiteInt

im2col
    :: forall rdkr rdkc md mr mc x
    . (KnownNat rdkr, KnownNat rdkc, KnownNat mc, KnownNat md, KnownNat mr, Num x, Storable x)
    => Proxy rdkr
    -> Proxy rdkc
    -> Proxy md
    -> Proxy mr
    -> Proxy mc
    -> Vector (md*mr*mc) x
    -> Matrix (md*(2*rdkr+1)*(2*rdkc+1)) (mr*mc) x
{-# INLINE im2col #-}
im2col :: Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
im2col Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
pmd Proxy mr
pmr Proxy mc
pmc Vector ((md * mr) * mc) x
mtx =
    let idxs :: Vector
  (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
idxs = Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
forall (rdkr :: Nat) (rdkc :: Nat) (mr :: Nat) (mc :: Nat)
       (md :: Nat).
(KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc,
 KnownNat md) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
im2colIndices Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
pmd Proxy mr
pmr Proxy mc
pmc
        mtx' :: Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
mtx' = Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
forall (rdkr :: Nat) (rdkc :: Nat) (mr :: Nat) (mc :: Nat)
       (md :: Nat) x.
(KnownNat rdkr, KnownNat rdkc, KnownNat md, KnownNat mr,
 KnownNat mc, Num x, Storable x) =>
Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
padMatrix Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy md
pmd Proxy mr
pmr Proxy mc
pmc Vector ((md * mr) * mc) x
mtx
     in Vector
  Vector (((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) * (mr * mc)) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector
   Vector (((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) * (mr * mc)) x
 -> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x)
-> Vector
     Vector (((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) * (mr * mc)) x
-> Matrix ((md * ((2 * rdkr) + 1)) * ((2 * rdkc) + 1)) (mr * mc) x
forall a b. (a -> b) -> a -> b
$ Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
-> Vector
     (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) x
forall a (m :: Nat) (n :: Nat).
Storable a =>
Vector m a -> Vector n Int -> Vector n a
backpermute Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
mtx' Vector
  (((((2 * rdkr) + 1) * ((2 * rdkc) + 1)) * md) * (mr * mc)) Int
idxs

windowIndices
    :: forall rdkr rdkc mr mc . (KnownNat rdkr, KnownNat rdkc, KnownNat mr, KnownNat mc)
    => Proxy rdkr
    -> Proxy rdkc
    -> Proxy mr
    -> Proxy mc
    -> Int
    -> Int
    -> Int
    -> Vector (mr*mc) Int
{-# INLINE windowIndices #-}
windowIndices :: Proxy rdkr
-> Proxy rdkc
-> Proxy mr
-> Proxy mc
-> Int
-> Int
-> Int
-> Vector (mr * mc) Int
windowIndices Proxy rdkr
prdkr Proxy rdkc
prdkc Proxy mr
pmr Proxy mc
pmc Int
kd Int
kr Int
kc =
    let rdkr :: Int
rdkr = Proxy rdkr -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkr
prdkr
        rdkc :: Int
rdkc = Proxy rdkc -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy rdkc
prdkc
        mr :: Int
mr = Proxy mr -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy mr
pmr
        mc :: Int
mc = Proxy mc -> Int
forall (n :: Nat). KnownNat n => Proxy n -> Int
natValInt Proxy mc
pmc
        mrc :: Int
mrc = Int
mrInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
mc
        nj' :: Int
nj' = Int
mr Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkr
        nk' :: Int
nk' = Int
mc Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
rdkc
        reIndex :: Int -> Int
reIndex Int
idx =
            let (Int
j,Int
k) = Int -> Int -> (Int, Int)
forall a. Integral a => a -> a -> (a, a)
divMod Int
idx Int
mc
             in Int -> Int -> (Int, Int, Int) -> Int
from3Index Int
nj' Int
nk' (Int
kd,Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kr,Int
kInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kc)
     in Vector Int -> Vector (mr * mc) Int
forall (v :: Type -> Type) (n :: Nat) a. v a -> Vector v n a
G.Vector (Vector Int -> Vector (mr * mc) Int)
-> Vector Int -> Vector (mr * mc) Int
forall a b. (a -> b) -> a -> b
$ Int -> (Int -> Int) -> Vector Int
forall a. Storable a => Int -> (Int -> a) -> Vector a
S.generate Int
mrc Int -> Int
reIndex

padMatrix
    :: forall rdkr rdkc mr mc md x
    . (KnownNat rdkr, KnownNat rdkc, KnownNat md, KnownNat mr, KnownNat mc, Num x, Storable x)
    => Proxy rdkr
    -> Proxy rdkc
    -> Proxy md
    -> Proxy mr
    -> Proxy mc
    -> Vector (md*mr*mc) x
    -> Vector (md*(mr + 2*rdkr)*(mc + 2*rdkc)) x
{-# INLINE padMatrix #-}
padMatrix :: Proxy rdkr
-> Proxy rdkc
-> Proxy md
-> Proxy mr
-> Proxy mc
-> Vector ((md * mr) * mc) x
-> Vector ((md * (mr + (2 * rdkr))) * (mc + (2 * rdkc))) x
padMatrix Proxy rdkr
_ Proxy rdkc
_ Proxy md
_ Proxy mr
_ Proxy mc
_ Vector ((md * mr) * mc) x
v =
    let mtxs :: Vector md (Matrix mr mc x)
        mtxs :: Vector md (Matrix mr mc x)
mtxs = (Vector (mr * mc) x -> Matrix mr mc x)
-> Vector md (Vector (mr * mc) x) -> Vector md (Matrix mr mc x)
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map Vector (mr * mc) x -> Matrix mr mc x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Vector v (m * n) a -> Matrix v m n a
G.Matrix (Vector md (Vector (mr * mc) x) -> Vector md (Matrix mr mc x))
-> Vector md (Vector (mr * mc) x) -> Vector md (Matrix mr mc x)
forall a b. (a -> b) -> a -> b
$ Vector (md * (mr * mc)) x -> Vector md (Vector (mr * mc) x)
forall (n :: Nat) (k :: Nat) a.
(KnownNat n, KnownNat k, Storable a) =>
Vector (n * k) a -> Vector n (Vector k a)
breakEvery Vector (md * (mr * mc)) x
Vector ((md * mr) * mc) x
v
        pdrs :: Vector rdkr (Vector mc x)
        pdrs :: Vector rdkr (Vector mc x)
pdrs = Vector mc x -> Vector rdkr (Vector mc x)
forall (n :: Nat) a. (KnownNat n, Storable a) => a -> Vector n a
replicate (Vector mc x -> Vector rdkr (Vector mc x))
-> Vector mc x -> Vector rdkr (Vector mc x)
forall a b. (a -> b) -> a -> b
$ x -> Vector mc x
forall (n :: Nat) a. (KnownNat n, Storable a) => a -> Vector n a
replicate x
0
        mtxs' :: Vector md (Matrix ((rdkr + mr) + rdkr) mc x)
mtxs' = (Matrix mr mc x -> Matrix ((rdkr + mr) + rdkr) mc x)
-> Vector md (Matrix mr mc x)
-> Vector md (Matrix ((rdkr + mr) + rdkr) mc x)
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map (\Matrix mr mc x
mtx -> Vector ((rdkr + mr) + rdkr) (Vector mc x)
-> Matrix ((rdkr + mr) + rdkr) mc x
forall (n :: Nat) x (m :: Nat).
(KnownNat n, Storable x) =>
Vector m (Vector n x) -> Matrix m n x
fromRows (Vector ((rdkr + mr) + rdkr) (Vector mc x)
 -> Matrix ((rdkr + mr) + rdkr) mc x)
-> Vector ((rdkr + mr) + rdkr) (Vector mc x)
-> Matrix ((rdkr + mr) + rdkr) mc x
forall a b. (a -> b) -> a -> b
$ Vector rdkr (Vector mc x)
pdrs Vector rdkr (Vector mc x)
-> Vector mr (Vector mc x) -> Vector (rdkr + mr) (Vector mc x)
forall (n :: Nat) (m :: Nat) a.
Storable a =>
Vector n a -> Vector m a -> Vector (n + m) a
++ Matrix mr mc x -> Vector mr (Vector mc x)
forall (m :: Nat) (n :: Nat) x.
(KnownNat m, KnownNat n, Storable x) =>
Matrix m n x -> Vector m (Vector n x)
toRows Matrix mr mc x
mtx Vector (rdkr + mr) (Vector mc x)
-> Vector rdkr (Vector mc x)
-> Vector ((rdkr + mr) + rdkr) (Vector mc x)
forall (n :: Nat) (m :: Nat) a.
Storable a =>
Vector n a -> Vector m a -> Vector (n + m) a
++ Vector rdkr (Vector mc x)
pdrs) Vector md (Matrix mr mc x)
mtxs
        pdcs :: Vector rdkc (Vector (mr + 2*rdkr) x)
        pdcs :: Vector rdkc (Vector (mr + (2 * rdkr)) x)
pdcs = Vector ((rdkr + mr) + rdkr) x
-> Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
forall (n :: Nat) a. (KnownNat n, Storable a) => a -> Vector n a
replicate (Vector ((rdkr + mr) + rdkr) x
 -> Vector rdkc (Vector ((rdkr + mr) + rdkr) x))
-> Vector ((rdkr + mr) + rdkr) x
-> Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
forall a b. (a -> b) -> a -> b
$ x -> Vector ((rdkr + mr) + rdkr) x
forall (n :: Nat) a. (KnownNat n, Storable a) => a -> Vector n a
replicate x
0
     in (Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x
 -> Vector (((rdkr + mr) + rdkr) * ((rdkc + mc) + rdkc)) x)
-> Vector
     md (Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
-> Vector (md * (((rdkr + mr) + rdkr) * ((rdkc + mc) + rdkc))) x
forall a b (m :: Nat) (n :: Nat).
(Storable a, Storable b) =>
(a -> Vector m b) -> Vector n a -> Vector (n * m) b
concatMap Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x
-> Vector (((rdkr + mr) + rdkr) * ((rdkc + mc) + rdkc)) x
forall (v :: Type -> Type) (m :: Nat) (n :: Nat) a.
Matrix v m n a -> Vector v (m * n) a
G.toVector (Vector
   md (Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
 -> Vector (md * (((rdkr + mr) + rdkr) * ((rdkc + mc) + rdkc))) x)
-> Vector
     md (Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
-> Vector (md * (((rdkr + mr) + rdkr) * ((rdkc + mc) + rdkc))) x
forall a b. (a -> b) -> a -> b
$ (Matrix ((rdkr + mr) + rdkr) mc x
 -> Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
-> Vector md (Matrix ((rdkr + mr) + rdkr) mc x)
-> Vector
     md (Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
forall a b (n :: Nat).
(Storable a, Storable b) =>
(a -> b) -> Vector n a -> Vector n b
map (\Matrix ((rdkr + mr) + rdkr) mc x
mtx' -> Vector Vector ((rdkc + mc) + rdkc) (Vector ((rdkr + mr) + rdkr) x)
-> Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x
forall (v :: Type -> Type) x (n :: Nat) (m :: Nat).
(Vector v x, Vector v Int, Vector v (Vector v n x),
 Vector v (Vector v m x), KnownNat n, KnownNat m) =>
Vector v n (Vector v m x) -> Matrix v m n x
G.fromColumns (Vector Vector ((rdkc + mc) + rdkc) (Vector ((rdkr + mr) + rdkr) x)
 -> Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x)
-> Vector
     Vector ((rdkc + mc) + rdkc) (Vector ((rdkr + mr) + rdkr) x)
-> Matrix Vector ((rdkr + mr) + rdkr) ((rdkc + mc) + rdkc) x
forall a b. (a -> b) -> a -> b
$ Vector rdkc (Vector (mr + (2 * rdkr)) x)
Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
pdcs Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
-> Vector mc (Vector ((rdkr + mr) + rdkr) x)
-> Vector (rdkc + mc) (Vector ((rdkr + mr) + rdkr) x)
forall (n :: Nat) (m :: Nat) a.
Storable a =>
Vector n a -> Vector m a -> Vector (n + m) a
++ Matrix ((rdkr + mr) + rdkr) mc x
-> Vector mc (Vector ((rdkr + mr) + rdkr) x)
forall (v :: Type -> Type) a (m :: Nat) (n :: Nat).
(Vector v a, Vector v (Vector v m a), KnownNat m, KnownNat n,
 Vector v Int) =>
Matrix v m n a -> Vector v n (Vector v m a)
G.toColumns Matrix ((rdkr + mr) + rdkr) mc x
mtx' Vector (rdkc + mc) (Vector ((rdkr + mr) + rdkr) x)
-> Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
-> Vector
     Vector ((rdkc + mc) + rdkc) (Vector ((rdkr + mr) + rdkr) x)
forall (n :: Nat) (m :: Nat) a.
Storable a =>
Vector n a -> Vector m a -> Vector (n + m) a
++ Vector rdkc (Vector (mr + (2 * rdkr)) x)
Vector rdkc (Vector ((rdkr + mr) + rdkr) x)
pdcs) Vector md (Matrix ((rdkr + mr) + rdkr) mc x)
mtxs'