{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE BangPatterns #-}  -- for Debug.Trace

-- |
-- Copyright   : (c) Johannes Kropp
-- License     : BSD 3-Clause
-- Maintainer  : Johannes Kropp <jodak932@gmail.com>

module Math.Nuha.Base where

import qualified Prelude as P
import Prelude hiding (map, replicate, (!!))
import Control.Monad (zipWithM_)
-- import qualified Debug.Trace as D
import Foreign.Storable (Storable, sizeOf)
import Data.Vector.Unboxed (Vector, Unbox)
import qualified Data.Vector.Unboxed as V
import qualified Data.Vector.Unboxed.Mutable as VM

import Math.Nuha.Types
import Math.Nuha.Internal


-- Naming conventions for indices:
-- idx : index for hValues :: Int
-- idcs : list of indices for hValues :: [Int]
-- mIdx : multiindex for a holor :: [Int]
-- mIdcs : list of multiindices for a holor :: [[Int]]


{- Some usefull links:
https://cheatsheets.quantecon.org/
http://hyperpolyglot.org/numerical-analysis
https://docs.julialang.org/en/v1/manual/arrays/
https://github.com/Magalame/fastest-matrices (benchmarks)
-}

{- Why unboxed immutable Vectors? :
https://stackoverflow.com/questions/34692809/lists-boxed-vectors-and-unboxed-vectors-for-heavy-scientific-computations
https://github.com/lehins/massiv
https://stackoverflow.com/ questions/40176678/differences-between-storable-and-unboxed-vectors
https://www.schoolofhaskell.com/user/commercial/content/vector
-}


-------------------
-- ** Constructions

{- | Basic function for creating a holor

>>> holor [2,2] [1,2,3,4]
  1.0  2.0
  3.0  4.0
>>> holor [2,2] [1,2,3,4] :: Holor Int
  1  2
  3  4
>>> holor [2,2,4] [1 .. 16]
[0,:,:] =
  1.0  2.0  3.0  4.0
  5.0  6.0  7.0  8.0
[1,:,:] =
  9.0  10.0  11.0  12.0
  13.0  14.0  15.0  16.0

-}
holor :: Unbox a
  => [Int] -- ^ shape
  -> [a] -- ^ values
  -> Holor a
{-# INLINE holor #-}
holor :: [Int] -> [a] -> Holor a
holor [Int]
shape [a]
values
    | [a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
values Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"holor : a holor can't be empty"
    | [Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
shape Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
2 = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"holor : a holor must have at least a dimension of 2"
    | (Int -> Int -> Int) -> Int -> [Int] -> Int
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr Int -> Int -> Int
forall a. Num a => a -> a -> a
(*) Int
1 [Int]
shape Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= [a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
values = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"holor : length of values "
        [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Int -> [Char]
forall a. Show a => a -> [Char]
show ([a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
values) [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
" doesn't match to shape " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Int] -> [Char]
forall a. Show a => a -> [Char]
show [Int]
shape
    | Bool
otherwise = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int]
shape ([Int] -> [Int]
fromShapeToStrides [Int]
shape) ([a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a]
values)

{- | Creates a column-vector (2d-holor with shape [m,1])

>>> vector [1,2,3,4]
 1.0
 2.0
 3.0
 4.0
-}
vector :: Unbox a => [a] -> Holor a
{-# INLINE vector #-}
vector :: [a] -> Holor a
vector [a]
values = [Int] -> [a] -> Holor a
forall a. Unbox a => [Int] -> [a] -> Holor a
holor [[a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
values, Int
1] [a]
values

{- | Creates a matrix (2d-holor)

>>> matrix [[1,2],[3,4]]
  1.0  2.0
  3.0  4.0
-}
matrix :: Unbox a => [[a]] -> Holor a
{-# INLINE matrix #-}
matrix :: [[a]] -> Holor a
matrix [[a]]
valuesList
    | Bool
lengths1dAreEqual = [Int] -> [a] -> Holor a
forall a. Unbox a => [Int] -> [a] -> Holor a
holor [Int]
shape [a]
values
    | Bool
otherwise = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"matrix : shape mismatch in nested list"
    where
        lengths1d :: [Int]
lengths1d = [[a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
v1 | [a]
v1 <- [[a]]
valuesList]
        lengths1dAreEqual :: Bool
lengths1dAreEqual = (Int -> Bool) -> [Int] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== [Int] -> Int
forall a. [a] -> a
head [Int]
lengths1d) ([Int] -> [Int]
forall a. [a] -> [a]
tail [Int]
lengths1d)
        values :: [a]
values = [a
entry | [a]
row <- [[a]]
valuesList , a
entry <- [a]
row]
        shape :: [Int]
shape = [[[a]] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [[a]]
valuesList, [Int] -> Int
forall a. [a] -> a
head [Int]
lengths1d]

-- | Creates a 2x1 holor from a tuple
vector2 :: Unbox a => T2 a -> Holor a
{-# INLINE vector2 #-}
vector2 :: T2 a -> Holor a
vector2 (a
t1,a
t2) = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
2,Int
1] [Int
1,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t1,a
t2]

-- | Creates a 3x1 holor from a tuple
vector3 :: Unbox a => T3 a -> Holor a
{-# INLINE vector3 #-}
vector3 :: T3 a -> Holor a
vector3 (a
t1,a
t2,a
t3) = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
3,Int
1] [Int
1,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t1,a
t2,a
t3]

-- | Creates a 4x1 holor from a tuple
vector4 :: Unbox a => T4 a -> Holor a
{-# INLINE vector4 #-}
vector4 :: T4 a -> Holor a
vector4 (a
t1,a
t2,a
t3,a
t4) = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
4,Int
1] [Int
1,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t1,a
t2,a
t3,a
t4]

-- | Creates a 2x2 holor from nested tuples
matrix22 :: Unbox a => T22 a -> Holor a
{-# INLINE matrix22 #-}
matrix22 :: T22 a -> Holor a
matrix22 ((a
t11,a
t12),(a
t21,a
t22)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
2,Int
2] [Int
2,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t21,a
t22]

-- | Creates a 3x2 holor from nested tuples
matrix32 :: Unbox a => T32 a -> Holor a
{-# INLINE matrix32 #-}
matrix32 :: T32 a -> Holor a
matrix32 ((a
t11,a
t12),(a
t21,a
t22),(a
t31,a
t32)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
3,Int
2] [Int
2,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t21,a
t22,a
t31,a
t32]

-- | Creates a 4x2 holor from nested tuples
matrix42 :: Unbox a => T42 a -> Holor a
{-# INLINE matrix42 #-}
matrix42 :: T42 a -> Holor a
matrix42 ((a
t11,a
t12),(a
t21,a
t22),(a
t31,a
t32),(a
t41,a
t42)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
4,Int
2] [Int
2,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t21,a
t22,a
t31,a
t32,a
t41,a
t42]

-- | Creates a 2x3 holor from nested tuples
matrix23 :: Unbox a => T23 a -> Holor a
{-# INLINE matrix23 #-}
matrix23 :: T23 a -> Holor a
matrix23 ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
2,Int
3] [Int
3,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23]

-- | Creates a 3x3 holor from nested tuples
matrix33 :: Unbox a => T33 a -> Holor a
{-# INLINE matrix33 #-}
matrix33 :: T33 a -> Holor a
matrix33 ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23),(a
t31,a
t32,a
t33)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
3,Int
3] [Int
3,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23,a
t31,a
t32,a
t33]

-- | Creates a 4x3 holor from nested tuples
matrix43 :: Unbox a => T43 a -> Holor a
{-# INLINE matrix43 #-}
matrix43 :: T43 a -> Holor a
matrix43 ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23),(a
t31,a
t32,a
t33),(a
t41,a
t42,a
t43)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
4,Int
3] [Int
3,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23,a
t31,a
t32,a
t33,a
t41,a
t42,a
t43]

-- | Creates a 2x4 holor from nested tuples
matrix24 :: Unbox a => T24 a -> Holor a
{-# INLINE matrix24 #-}
matrix24 :: T24 a -> Holor a
matrix24 ((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
2,Int
4] [Int
4,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24]

-- | Creates a 3x4 holor from nested tuples
matrix34 :: Unbox a => T34 a -> Holor a
{-# INLINE matrix34 #-}
matrix34 :: T34 a -> Holor a
matrix34 ((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24),(a
t31,a
t32,a
t33,a
t34)) =
    [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int
3,Int
4] [Int
4,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24,a
t31,a
t32,a
t33,a
t34]

-- | Creates a 4x4 holor from nested tuples
matrix44 :: Unbox a => T44 a -> Holor a
{-# INLINE matrix44 #-}
matrix44 :: T44 a -> Holor a
matrix44 ((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24),(a
t31,a
t32,a
t33,a
t34),(a
t41,a
t42,a
t43,a
t44)) = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor
    [Int
4,Int
4] [Int
4,Int
1] (Vector a -> Holor a) -> Vector a -> Holor a
forall a b. (a -> b) -> a -> b
$ [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24,a
t31,a
t32,a
t33,a
t34,a
t41,a
t42,a
t43,a
t44]

{- | Creates a holor for a given shape with all entries the same

>>> replicate [1,3] 0
  0.0  0.0  0.0
-}
replicate :: Unbox a => [Int] -> a -> Holor a
{-# INLINE replicate #-}
replicate :: [Int] -> a -> Holor a
replicate [Int]
shape a
value = [Int] -> [a] -> Holor a
forall a. Unbox a => [Int] -> [a] -> Holor a
holor [Int]
shape (Int -> a -> [a]
forall a. Int -> a -> [a]
P.replicate ([Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Int]
shape) a
value)


----------------
-- ** Properties

-- | Shape of a holor
shape :: Unbox a => Holor a -> [Int]
{-# INLINE shape #-}
shape :: Holor a -> [Int]
shape = Holor a -> [Int]
forall a. Holor a -> [Int]
hShape

-- | Number of holor entrys
numElems :: Unbox a => Holor a -> Int
{-# INLINE numElems #-}
numElems :: Holor a -> Int
numElems Holor a
hlr = [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr)

-- | Size of all holor elements in bytes (not equal to consumed memory of the holor)
sizeOfElems :: (Unbox a, Storable a) => Holor a -> Int
{-# INLINE sizeOfElems #-}
sizeOfElems :: Holor a -> Int
sizeOfElems Holor a
hlr = (Holor a -> Int
forall a. Unbox a => Holor a -> Int
numElems Holor a
hlr) Int -> Int -> Int
forall a. Num a => a -> a -> a
* (Holor a -> Int
forall a. (Unbox a, Storable a) => Holor a -> Int
sizeOfElem Holor a
hlr)

-- | Size of a single holor element in bytes
sizeOfElem :: forall a . (Unbox a, Storable a) => Holor a -> Int
{-# INLINE sizeOfElem #-}
sizeOfElem :: Holor a -> Int
sizeOfElem Holor a
hlr = a -> Int
forall a. Storable a => a -> Int
sizeOf (a
forall a. HasCallStack => a
undefined::a) -- requires "ScopedTypeVariables" and "forall a"

{- | Dimension of a holor (length of shape)

>>> dim $ holor [2,2,2] [0..7]
3
-}
dim :: Unbox a => Holor a -> Int
{-# INLINE dim #-}
dim :: Holor a -> Int
dim = [Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length ([Int] -> Int) -> (Holor a -> [Int]) -> Holor a -> Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Holor a -> [Int]
forall a. Holor a -> [Int]
hShape


---------------
-- ** Accessors

infixl 9 !
infixl 9 !!
infixl 9 |!
infixl 9 |!!
infixl 9 ||!
infixl 9 ||!!
infixl 9 |||!
infixl 9 |||!!

{- | Indexing a holor with shape (n,1) (column-vector) or (1,n) (row-vector)

>>> v = vector [5,6,7]
>>> v!1
6.0
-}
(!) :: Unbox a => Holor a -> Int -> a
{-# INLINE (!) #-}
(!) Holor a
vec Int
idx
    | Bool -> Bool
not (Bool -> Bool) -> Bool -> Bool
forall a b. (a -> b) -> a -> b
$ [Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
shape Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
2 Bool -> Bool -> Bool
&& ([Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Int]
shape Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
len)
        = [Char] -> a
forall a. HasCallStack => [Char] -> a
error ([Char] -> a) -> [Char] -> a
forall a b. (a -> b) -> a -> b
$ [Char]
"(!) : indexing with single Int only possible for column or row vectors"
    | Bool -> Bool
not (Bool -> Bool) -> Bool -> Bool
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Bool
isValidIdx Int
len Int
idx
        = [Char] -> a
forall a. HasCallStack => [Char] -> a
error ([Char] -> a) -> [Char] -> a
forall a b. (a -> b) -> a -> b
$ [Char]
"(!) : index out of range"
    | Bool
otherwise = Vector a
values Vector a -> Int -> a
forall a. Unbox a => Vector a -> Int -> a
V.! Int
idx
    where
        shape :: [Int]
shape = Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
vec
        values :: Vector a
values = Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
vec
        len :: Int
len = Vector a -> Int
forall a. Unbox a => Vector a -> Int
V.length Vector a
values

-- | Unsafe version of (!) without checking for row or column vector and not checking bounds
(!!) :: Unbox a => Holor a -> Int -> a
{-# INLINE (!!) #-}
!! :: Holor a -> Int -> a
(!!) Holor a
vec Int
idx = (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
vec) Vector a -> Int -> a
forall a. Unbox a => Vector a -> Int -> a
V.! Int
idx


{- | Indexing a single entry of a holor

>>> m = matrix [[2,3,4],[5,6,3],[9,0,2]]
>>> m|![0,2]
4.0
-}
(|!) :: Unbox a => Holor a -> [Int] -> a
{-# INLINE (|!) #-}
|! :: Holor a -> [Int] -> a
(|!) Holor a
hlr [Int]
mIdx
    | [Int] -> [Int] -> Bool
isValidMIdx (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr) [Int]
mIdx = Holor a
hlrHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|!![Int]
mIdx
    | Bool
otherwise = [Char] -> a
forall a. HasCallStack => [Char] -> a
error ([Char] -> a) -> [Char] -> a
forall a b. (a -> b) -> a -> b
$ [Char]
"(|!) : multiindex " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Int] -> [Char]
forall a. Show a => a -> [Char]
show [Int]
mIdx [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
" out of bounds"

-- | Unsafe version of (|!) without checking bounds
(|!!) :: Unbox a => Holor a -> [Int] -> a
{-# INLINE (|!!) #-}
|!! :: Holor a -> [Int] -> a
(|!!) Holor a
hlr [Int]
mIdx = Vector a -> Int -> a
forall a. Unbox a => Vector a -> Int -> a
V.unsafeIndex (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr) Int
vectorindex where
    vectorindex :: Int
vectorindex = [Int] -> [Int] -> Int
fromMultiIndexToIndex (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) [Int]
mIdx

{- | Extract a subholor by indices in the first dimension. Applied to a vector the result is a vector with the elements at the given element-indices. Applied to a matrix the result is a matrix with the rows at the given row-indices.

>>> h = holor [4,4] [1..16]
>>> h
  1.0  2.0  3.0  4.0
  5.0  6.0  7.0  8.0
  9.0  10.0  11.0  12.0
  13.0  14.0  15.0  16.0
>>> h||![1,3]
  5.0  6.0  7.0  8.0
  13.0  14.0  15.0  16.0
-}
(||!) :: Unbox a => Holor a -> [Int] -> Holor a
{-# INLINE (||!) #-}
||! :: Holor a -> [Int] -> Holor a
(||!) Holor a
hlr [Int]
idcs = case (Int -> Bool) -> [Int] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
any (\Int
i -> Int
iInt -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>(Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1) Bool -> Bool -> Bool
|| Int
iInt -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<Int
0) [Int]
idcs Bool -> Bool -> Bool
|| [Int]
idcs[Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
==[] of
    Bool
True -> [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"(||!) : invalid indices " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Int] -> [Char]
forall a. Show a => a -> [Char]
show [Int]
idcs
    Bool
False -> Holor a
hlrHolor a -> [Int] -> Holor a
forall a. Unbox a => Holor a -> [Int] -> Holor a
||!![Int]
idcs
    where
        m :: Int
m = [Int] -> Int
forall a. [a] -> a
head ([Int] -> Int) -> [Int] -> Int
forall a b. (a -> b) -> a -> b
$ Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr

-- | Unsafe version of (||!) without checking indices
(||!!) :: Unbox a => Holor a -> [Int] -> Holor a
{-# INLINE (||!!) #-}
||!! :: Holor a -> [Int] -> Holor a
(||!!) Holor a
hlr [Int]
idcs = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int]
shapeOut [Int]
stridesOut Vector a
valuesOut where
    valuesIn :: Vector a
valuesIn = Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr
    shapeInTail :: [Int]
shapeInTail = [Int] -> [Int]
forall a. [a] -> [a]
tail ([Int] -> [Int]) -> [Int] -> [Int]
forall a b. (a -> b) -> a -> b
$ Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr
    sliceLength :: Int
sliceLength = [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Int]
shapeInTail
    shapeOut :: [Int]
shapeOut = [Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
idcs Int -> [Int] -> [Int]
forall a. a -> [a] -> [a]
: [Int]
shapeInTail
    stridesOut :: [Int]
stridesOut = [Int] -> [Int]
fromShapeToStrides [Int]
shapeOut
    valuesOut :: Vector a
valuesOut = case [Int]
idcs of
        [] -> [Char] -> Vector a
forall a. HasCallStack => [Char] -> a
error [Char]
"(||!!) : invalid empty indices"
        [Int
idx] -> Int -> Int -> Vector a -> Vector a
forall a. Unbox a => Int -> Int -> Vector a -> Vector a
V.unsafeSlice (Int
idxInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
sliceLength) Int
sliceLength Vector a
valuesIn
        [Int]
_ -> [Vector a] -> Vector a
forall a. Unbox a => [Vector a] -> Vector a
V.concat [Int -> Int -> Vector a -> Vector a
forall a. Unbox a => Int -> Int -> Vector a -> Vector a
V.unsafeSlice (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
sliceLength) Int
sliceLength Vector a
valuesIn | Int
i <- [Int]
idcs]


{- | Multidimensional indexing. For each dimension a list of indices define the axes to extract. This means multiindexing is done by a nested list of indices which has the same length as the shape. The dimension of the extracted subholor remains the same as that of the original holor.

>>> h = holor [4,4] [1..16]
>>> h
  1.0  2.0  3.0  4.0
  5.0  6.0  7.0  8.0
  9.0  10.0  11.0  12.0
  13.0  14.0  15.0  16.0
>>> h|||![[1,2,3],[0,1,1,2,3]]
  5.0  6.0  6.0  7.0  8.0
  9.0  10.0  10.0  11.0  12.0
  13.0  14.0  14.0  15.0  16.0
-}
(|||!) :: Unbox a => Holor a -> [[Int]] -> Holor a
{-# INLINE (|||!) #-}
|||! :: Holor a -> [[Int]] -> Holor a
(|||!) Holor a
hlr [[Int]]
mIdcs
    | [Int] -> [[Int]] -> Bool
isValidMIdcs (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr) [[Int]]
mIdcs = Holor a
hlrHolor a -> [[Int]] -> Holor a
forall a. Unbox a => Holor a -> [[Int]] -> Holor a
|||!![[Int]]
mIdcs
    | Bool
otherwise = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"(|||!) : multiindices " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [[Int]] -> [Char]
forall a. Show a => a -> [Char]
show [[Int]]
mIdcs [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
" invalid"

(|||!!) :: Unbox a => Holor a -> [[Int]] -> Holor a
{-# INLINE (|||!!) #-}
-- | Unsafe version of (|||!) without checking multiindices
|||!! :: Holor a -> [[Int]] -> Holor a
(|||!!) Holor a
hlr [[Int]]
mIdcs = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int]
shp ([Int] -> [Int]
fromShapeToStrides [Int]
shp) ([a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [a]
values) where
    shp :: [Int]
shp = ([Int] -> Int) -> [[Int]] -> [Int]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap [Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [[Int]]
mIdcs
    values :: [a]
values = [ Holor a
hlrHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|!![Int]
indices | [Int]
indices <- [[Int]] -> [[Int]]
cartProd [[Int]]
mIdcs]



----------------
-- ** Operations

-- | Map a function over all entries of a holor
map :: (Unbox a, Unbox b) => (a -> b) -> Holor a -> Holor b
{-# INLINE map #-}
map :: (a -> b) -> Holor a -> Holor b
map a -> b
f Holor a
hlr = [Int] -> [Int] -> Vector b -> Holor b
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr) (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) ((a -> b) -> Vector a -> Vector b
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
V.map a -> b
f (Vector a -> Vector b) -> Vector a -> Vector b
forall a b. (a -> b) -> a -> b
$ Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)

{- | Filter the values of a holor by a condition. The remaining elements are returned in a list.

>>> h = holor [3,3,3] [1 .. 27] :: Holor Int
>>> filter (\e -> mod e 3 == 0) h
[3,6,9,12,15,18,21,24,27]
-}
filter :: (Unbox a) => (a -> Bool) -> Holor a -> [a]
{-# INLINE filter #-}
filter :: (a -> Bool) -> Holor a -> [a]
filter a -> Bool
cond Holor a
hlr = Vector a -> [a]
forall a. Unbox a => Vector a -> [a]
V.toList (Vector a -> [a]) -> Vector a -> [a]
forall a b. (a -> b) -> a -> b
$ (a -> Bool) -> Vector a -> Vector a
forall a. Unbox a => (a -> Bool) -> Vector a -> Vector a
V.filter a -> Bool
cond (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)

-- TODO: could the filter function be used for this?
-- TODO: may be working also with V.findIndices!
{- | Select all multiindices of a holor whose corresponding elements fulfill a condition.

>>> m = matrix [[-1,1],[-2,3]]
>>> selectBy (>0) m
  [[0,1],[1,1]]
-}
selectBy :: Unbox a => (a -> Bool) -> Holor a -> [[Int]]
{-# INLINE selectBy #-}
selectBy :: (a -> Bool) -> Holor a -> [[Int]]
selectBy a -> Bool
cond Holor a
hlr = Int -> [[Int]] -> [[Int]]
fillList Int
0 [] where
    fillList :: Int -> [[Int]] -> [[Int]]
fillList Int
idx [[Int]]
list
        | Int
idx Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Vector a -> Int
forall a. Unbox a => Vector a -> Int
V.length (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)
            = case () of
              ()
_ | a -> Bool
cond a
val -> ([Int] -> Int -> [Int]
fromIndexToMultiIndex (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) Int
idx) [Int] -> [[Int]] -> [[Int]]
forall a. a -> [a] -> [a]
: (Int -> [[Int]] -> [[Int]]
fillList (Int
idxInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) [[Int]]
list)
                | Bool
otherwise -> Int -> [[Int]] -> [[Int]]
fillList (Int
idxInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) [[Int]]
list
        | Bool
otherwise = [[Int]]
list
        where val :: a
val = (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)Vector a -> Int -> a
forall a. Unbox a => Vector a -> Int -> a
V.!Int
idx

-- | Count elements of a holor that fulfill a condition
countBy :: Unbox a => (a -> Bool) -> Holor a -> Int
{-# INLINE countBy #-}
countBy :: (a -> Bool) -> Holor a -> Int
countBy a -> Bool
cond Holor a
hlr = Vector Int -> Int
forall a. (Unbox a, Num a) => Vector a -> a
V.sum (Vector Int -> Int) -> Vector Int -> Int
forall a b. (a -> b) -> a -> b
$ (Bool -> Int) -> Vector Bool -> Vector Int
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
V.map Bool -> Int
forall a. Enum a => a -> Int
fromEnum ((a -> Bool) -> Vector a -> Vector Bool
forall a b. (Unbox a, Unbox b) => (a -> b) -> Vector a -> Vector b
V.map a -> Bool
cond (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr))

{- | Accumulate the values of a holor by a function.

>>> v = vector [2,3,4,1] :: Holor Int
>>> accumulate (+) v
  [2,5,9,10]
>>> accumulate max v
  [2,3,4,4]
-}
accumulate :: Unbox a => (a -> a -> a) -> Holor a -> [a]
{-# INLINE accumulate #-}
accumulate :: (a -> a -> a) -> Holor a -> [a]
accumulate a -> a -> a
f Holor a
hlr = Vector a -> [a]
forall a. Unbox a => Vector a -> [a]
V.toList (Vector a -> [a]) -> Vector a -> [a]
forall a b. (a -> b) -> a -> b
$ (a -> a -> a) -> Vector a -> Vector a
forall a. Unbox a => (a -> a -> a) -> Vector a -> Vector a
V.scanl1 a -> a -> a
f (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)

-- | Flatten a holor to a column vector
flatten :: Unbox a => Holor a -> Holor a
{-# INLINE flatten #-}
flatten :: Holor a -> Holor a
flatten Holor a
hlr = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [[Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product ([Int] -> Int) -> [Int] -> Int
forall a b. (a -> b) -> a -> b
$ Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr, Int
1] [Int
1,Int
1] (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)

-- | Transpose a holor
transpose :: Unbox a => Holor a -> Holor a
{-# INLINE transpose #-}
transpose :: Holor a -> Holor a
transpose Holor a
hlr = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int]
shapeOut [Int]
stridesOut Vector a
valuesOut where
    shapeOut :: [Int]
shapeOut = [Int] -> [Int]
forall a. [a] -> [a]
reverse (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr)
    stridesOut :: [Int]
stridesOut = [Int] -> [Int]
fromShapeToStrides [Int]
shapeOut
    valuesOut :: Vector a
valuesOut = [a] -> Vector a
forall a. Unbox a => [a] -> Vector a
V.fromList [ Holor a
hlrHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|!!([Int] -> [Int]
forall a. [a] -> [a]
reverse [Int]
idcs) | [Int]
idcs <- [Int] -> [[Int]]
fromShapeToMultiIndices [Int]
shapeOut]

{-
-- TODO:
-- Idea: Use V.unsafeSlice. Should be easier than flipCols since values are in row-major order.
-- | flip rows in reversed order for a given dimension
flipRows :: Int -> Holor a -> Holor a
{-# INLINE flipRows #-}
flipRows dim hlr = undefined
    -- = Holor (shape hlr) (hStrides hlr) values where
    -- values = undefined

-- TODO:
-- Idea: first transpose holor, then use flipRows, then transpose again?
-- | flip columns in reversed order for a given dimension
flipCols :: Int -> Holor a -> Holor a
{-# INLINE flipCols #-}
flipCols dim hlr = undefined
-}

{-
-- TODO:
tile = undefined

-- TODO:
concat = undefined
-}

{- | Reshape a holor without changing the order of the underlying values

>>> h = holor [2,2,2] [1 .. 8]
>>> reshape [2,4] h
  1.0  2.0  3.0  4.0
  5.0  6.0  7.0  8.0
-}
reshape :: Unbox a => [Int] -> Holor a -> Holor a
{-# INLINE reshape #-}
reshape :: [Int] -> Holor a -> Holor a
reshape [Int]
shp Holor a
hlr
    | [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr) Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Int]
shp = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor [Int]
shp ([Int] -> [Int]
fromShapeToStrides [Int]
shp) (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)
    | Bool
otherwise = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$
        [Char]
"reshape : cannot reshape from " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Int] -> [Char]
forall a. Show a => a -> [Char]
show (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr) [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
" to " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Int] -> [Char]
forall a. Show a => a -> [Char]
show [Int]
shp


-- | Get diagonal of a quadratical matrix (2d-holor) as a column vector
diagonal :: Unbox a => Holor a -> Holor a
{-# INLINE diagonal #-}
diagonal :: Holor a -> Holor a
diagonal Holor a
hlr
    | Holor a -> Bool
forall a. Unbox a => Holor a -> Bool
isSquare Holor a
hlr = [a] -> Holor a
forall a. Unbox a => [a] -> Holor a
vector [Holor a
hlr Holor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|!! [Int
i,Int
i] | Int
i <- [Int
0 .. [Int] -> Int
forall a. [a] -> a
head (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr)Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]]
    | Bool
otherwise = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"diagonal : Holor must be in quadratical matrix form"



-------------------
-- ** Conversions

-- | Returns the holor values in a list
toList :: Unbox a => Holor a -> [a]
{-# INLINE toList #-}
toList :: Holor a -> [a]
toList Holor a
hlr = Vector a -> [a]
forall a. Unbox a => Vector a -> [a]
V.toList (Vector a -> [a]) -> Vector a -> [a]
forall a b. (a -> b) -> a -> b
$ Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr

-- | Returns the values of an 2d-holor in a list of lists
toList2 :: Unbox a => Holor a -> [[a]]
{-# INLINE toList2 #-}
toList2 :: Holor a -> [[a]]
toList2 Holor a
hlr
    | Bool -> Bool
not ([Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
shp Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
2) = [Char] -> [[a]]
forall a. HasCallStack => [Char] -> a
error ([Char] -> [[a]]) -> [Char] -> [[a]]
forall a b. (a -> b) -> a -> b
$ [Char]
"toList2 : not a 2d-holor"
    | Bool
otherwise = [ [ Holor a
hlrHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|!![Int
i,Int
j] | Int
j <- [Int
0 .. Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]] | Int
i <- [Int
0 .. Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]]
    where
        shp :: [Int]
shp = Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr
        [Int
m,Int
n] = [Int]
shp

-- | Converts a holor with a single element to the element itself
toScalar :: Unbox a => Holor a -> a
{-# INLINE toScalar #-}
toScalar :: Holor a -> a
toScalar Holor a
hlr
    | Vector a -> Int
forall a. Unbox a => Vector a -> Int
V.length (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr) Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 = Vector a -> a
forall a. Unbox a => Vector a -> a
V.head (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)
    | Bool
otherwise = [Char] -> a
forall a. HasCallStack => [Char] -> a
error ([Char] -> a) -> [Char] -> a
forall a b. (a -> b) -> a -> b
$ [Char]
"toScalar : not a holor with a single element"


-- | Converts a holor to a 2-tuple if possible
toT2 :: Unbox a => Holor a -> T2 a
{-# INLINE toT2 #-}
toT2 :: Holor a -> T2 a
toT2 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
2,Int
1] Bool -> Bool -> Bool
|| Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
1,Int
2] = (a
t1,a
t2)
    | Bool
otherwise = [Char] -> T2 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T2 a) -> [Char] -> T2 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT2 : shape mismatch"
    where
        [a
t1,a
t2] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to a 3-tuple if possible
toT3 :: Unbox a => Holor a -> T3 a
{-# INLINE toT3 #-}
toT3 :: Holor a -> T3 a
toT3 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
3,Int
1] Bool -> Bool -> Bool
|| Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
1,Int
3] = (a
t1,a
t2,a
t3)
    | Bool
otherwise = [Char] -> T3 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T3 a) -> [Char] -> T3 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT3 : shape mismatch"
    where
        [a
t1,a
t2,a
t3] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to a 4-tuple if possible
toT4 :: Unbox a => Holor a -> T4 a
{-# INLINE toT4 #-}
toT4 :: Holor a -> T4 a
toT4 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
4,Int
1] Bool -> Bool -> Bool
|| Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
1,Int
4] = (a
t1,a
t2,a
t3,a
t4)
    | Bool
otherwise = [Char] -> T4 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T4 a) -> [Char] -> T4 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT4 : shape mismatch"
    where
        [a
t1,a
t2,a
t3,a
t4] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 2,2-tuple if possible
toT22 :: Unbox a => Holor a -> T22 a
{-# INLINE toT22 #-}
toT22 :: Holor a -> T22 a
toT22 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
2,Int
2] = ((a
t11,a
t12),(a
t21,a
t22))
    | Bool
otherwise = [Char] -> T22 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T22 a) -> [Char] -> T22 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT22 : shape mismatch"
    where
        [a
t11,a
t12,a
t21,a
t22] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 3,2-tuple if possible
toT32 :: Unbox a => Holor a -> T32 a
{-# INLINE toT32 #-}
toT32 :: Holor a -> T32 a
toT32 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
3,Int
2] = ((a
t11,a
t12),(a
t21,a
t22),(a
t31,a
t32))
    | Bool
otherwise = [Char] -> T32 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T32 a) -> [Char] -> T32 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT32 : shape mismatch"
    where
        [a
t11,a
t12,a
t21,a
t22,a
t31,a
t32] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 4,2-tuple if possible
toT42 :: Unbox a => Holor a -> T42 a
{-# INLINE toT42 #-}
toT42 :: Holor a -> T42 a
toT42 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
4,Int
2] = ((a
t11,a
t12),(a
t21,a
t22),(a
t31,a
t32),(a
t41,a
t42))
    | Bool
otherwise = [Char] -> T42 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T42 a) -> [Char] -> T42 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT42 : shape mismatch"
    where
        [a
t11,a
t12,a
t21,a
t22,a
t31,a
t32,a
t41,a
t42] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 2,3-tuple if possible
toT23 :: Unbox a => Holor a -> T23 a
{-# INLINE toT23 #-}
toT23 :: Holor a -> T23 a
toT23 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
2,Int
3] = ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23))
    | Bool
otherwise = [Char] -> T23 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T23 a) -> [Char] -> T23 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT23 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 3,3-tuple if possible
toT33 :: Unbox a => Holor a -> T33 a
{-# INLINE toT33 #-}
toT33 :: Holor a -> T33 a
toT33 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
3,Int
3] = ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23),(a
t31,a
t32,a
t33))
    | Bool
otherwise = [Char] -> T33 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T33 a) -> [Char] -> T33 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT33 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23,a
t31,a
t32,a
t33] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 4,3-tuple if possible
toT43 :: Unbox a => Holor a -> T43 a
{-# INLINE toT43 #-}
toT43 :: Holor a -> T43 a
toT43 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
4,Int
3] = ((a
t11,a
t12,a
t13),(a
t21,a
t22,a
t23),(a
t31,a
t32,a
t33),(a
t41,a
t42,a
t43))
    | Bool
otherwise = [Char] -> T43 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T43 a) -> [Char] -> T43 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT43 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t21,a
t22,a
t23,a
t31,a
t32,a
t33,a
t41,a
t42,a
t43] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 2,4-tuple if possible
toT24 :: Unbox a => Holor a -> T24 a
{-# INLINE toT24 #-}
toT24 :: Holor a -> T24 a
toT24 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
2,Int
4] = ((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24))
    | Bool
otherwise = [Char] -> T24 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T24 a) -> [Char] -> T24 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT24 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 3,4-tuple if possible
toT34 :: Unbox a => Holor a -> T34 a
{-# INLINE toT34 #-}
toT34 :: Holor a -> T34 a
toT34 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
3,Int
4] = ((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24),(a
t31,a
t32,a
t33,a
t34))
    | Bool
otherwise = [Char] -> T34 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T34 a) -> [Char] -> T34 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT34 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24,a
t31,a
t32,a
t33,a
t34] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr

-- | Converts a holor to nested 4,4-tuple if possible
toT44 :: Unbox a => Holor a -> T44 a
{-# INLINE toT44 #-}
toT44 :: Holor a -> T44 a
toT44 Holor a
hlr
    | Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
==[Int
4,Int
4] =((a
t11,a
t12,a
t13,a
t14),(a
t21,a
t22,a
t23,a
t24),(a
t31,a
t32,a
t33,a
t34),(a
t41,a
t42,a
t43,a
t44))
    | Bool
otherwise = [Char] -> T44 a
forall a. HasCallStack => [Char] -> a
error ([Char] -> T44 a) -> [Char] -> T44 a
forall a b. (a -> b) -> a -> b
$ [Char]
"toT44 : shape mismatch"
    where
        [a
t11,a
t12,a
t13,a
t14,a
t21,a
t22,a
t23,a
t24,a
t31,a
t32,a
t33,a
t34,a
t41,a
t42,a
t43,a
t44] = Holor a -> [a]
forall a. Unbox a => Holor a -> [a]
toList Holor a
hlr


-------------------
-- ** Modifications

{- | Modify a single entry of a holor

>>> v = vector [1,2,3]
>>> setElem [1,0] 22 v
  1.0  22.0  3.0
-}
setElem :: Unbox a => [Int] -> a -> Holor a -> Holor a
{-# INLINE setElem #-}
setElem :: [Int] -> a -> Holor a -> Holor a
setElem [Int]
mIdx a
value Holor a
hlr
    | Bool -> Bool
not ([Int] -> [Int] -> Bool
isValidMIdx (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr) [Int]
mIdx) = [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"setElem : multiindex invalid"
    | Bool
otherwise = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr) (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) Vector a
valuesOut
    where
    valuesOut :: Vector a
valuesOut = (forall s. ST s (MVector s a)) -> Vector a
forall a. Unbox a => (forall s. ST s (MVector s a)) -> Vector a
V.create ((forall s. ST s (MVector s a)) -> Vector a)
-> (forall s. ST s (MVector s a)) -> Vector a
forall a b. (a -> b) -> a -> b
$ do
        -- first yields mutable copy, then modifies and freezes it:
        MVector s a
vec <- Vector a -> ST s (MVector (PrimState (ST s)) a)
forall a (m :: * -> *).
(Unbox a, PrimMonad m) =>
Vector a -> m (MVector (PrimState m) a)
V.thaw (Vector a -> ST s (MVector (PrimState (ST s)) a))
-> Vector a -> ST s (MVector (PrimState (ST s)) a)
forall a b. (a -> b) -> a -> b
$ Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr
        let idx :: Int
idx = [Int] -> [Int] -> Int
fromMultiIndexToIndex (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) [Int]
mIdx
        MVector (PrimState (ST s)) a -> Int -> a -> ST s ()
forall (m :: * -> *) a.
(PrimMonad m, Unbox a) =>
MVector (PrimState m) a -> Int -> a -> m ()
VM.unsafeWrite MVector s a
MVector (PrimState (ST s)) a
vec Int
idx a
value
        MVector s a -> ST s (MVector s a)
forall (m :: * -> *) a. Monad m => a -> m a
return MVector s a
vec


{- | Modify many elements of a holor at once

>>> h = zeros [3,3] :: Holor Int
>>> setElems [[0,0],[1,1],[2,2]] [4,5,6] h
  4  0  0
  0  5  0
  0  0  6
-}
setElems :: Unbox a => [[Int]] -> [a] -> Holor a -> Holor a
{-# INLINE setElems #-}
-- TODO: Is there an faster alternative with V.(//) or V.unsafeUpd ?
setElems :: [[Int]] -> [a] -> Holor a -> Holor a
setElems [[Int]]
mIdxList [a]
values Holor a
hlr
    | [[Int]] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [[Int]]
mIdxList Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= [a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
values =
        [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"setElems : length of list with multiindices unequal to length of values"
    | Bool -> Bool
not ([Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and (([Int] -> Bool) -> [[Int]] -> [Bool]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap ([Int] -> [Int] -> Bool
isValidMIdx (Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr)) [[Int]]
mIdxList)) =
        [Char] -> Holor a
forall a. HasCallStack => [Char] -> a
error ([Char] -> Holor a) -> [Char] -> Holor a
forall a b. (a -> b) -> a -> b
$ [Char]
"setElems : list of multiindices invalid"
    | Bool
otherwise = [Int] -> [Int] -> Vector a -> Holor a
forall a. [Int] -> [Int] -> Vector a -> Holor a
Holor (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr) (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr) Vector a
valuesOut
    where
    valuesOut :: Vector a
valuesOut = (forall s. ST s (MVector s a)) -> Vector a
forall a. Unbox a => (forall s. ST s (MVector s a)) -> Vector a
V.create ((forall s. ST s (MVector s a)) -> Vector a)
-> (forall s. ST s (MVector s a)) -> Vector a
forall a b. (a -> b) -> a -> b
$ do
        -- first yields mutable copy, then modifies and freezes it:
        MVector s a
vec <- Vector a -> ST s (MVector (PrimState (ST s)) a)
forall a (m :: * -> *).
(Unbox a, PrimMonad m) =>
Vector a -> m (MVector (PrimState m) a)
V.thaw (Vector a -> ST s (MVector (PrimState (ST s)) a))
-> Vector a -> ST s (MVector (PrimState (ST s)) a)
forall a b. (a -> b) -> a -> b
$ Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr
        let idcs :: [Int]
idcs = ([Int] -> Int) -> [[Int]] -> [Int]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap ([Int] -> [Int] -> Int
fromMultiIndexToIndex (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr)) [[Int]]
mIdxList
        (Int -> a -> ST s ()) -> [Int] -> [a] -> ST s ()
forall (m :: * -> *) a b c.
Applicative m =>
(a -> b -> m c) -> [a] -> [b] -> m ()
zipWithM_ (MVector (PrimState (ST s)) a -> Int -> a -> ST s ()
forall (m :: * -> *) a.
(PrimMonad m, Unbox a) =>
MVector (PrimState m) a -> Int -> a -> m ()
VM.unsafeWrite MVector s a
MVector (PrimState (ST s)) a
vec) [Int]
idcs [a]
values
        MVector s a -> ST s (MVector s a)
forall (m :: * -> *) a. Monad m => a -> m a
return MVector s a
vec


-------------------
-- ** Checks

-- | Tests if a holor is instanced correctly, i.e. shape, strides and values fit together
isValidHolor :: Unbox a => Holor a -> Bool
{-# INLINE isValidHolor #-}
isValidHolor :: Holor a -> Bool
isValidHolor Holor a
holor = Vector a -> Int
forall a. Unbox a => Vector a -> Int
V.length (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
holor) Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
holor)
    Bool -> Bool -> Bool
&& [Int] -> [Int]
fromShapeToStrides (Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
holor) [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
holor

-- | Tests if a holor has shape [n,1] (column-vector) or [1,n] (row-vector)
isVector :: Unbox a => Holor a -> Bool
{-# INLINE isVector #-}
isVector :: Holor a -> Bool
isVector Holor a
hlr =
    ([Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
shp Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
2) Bool -> Bool -> Bool
&& ([Int] -> Int
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Int]
shp Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Vector a -> Int
forall a. Unbox a => Vector a -> Int
V.length (Holor a -> Vector a
forall a. Holor a -> Vector a
hValues Holor a
hlr)) Bool -> Bool -> Bool
&& (Holor a -> [Int]
forall a. Holor a -> [Int]
hStrides Holor a
hlr [Int] -> [Int] -> Bool
forall a. Eq a => a -> a -> Bool
== [Int
1,Int
1])
    where
        shp :: [Int]
shp = Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
hlr

-- | Tests if a holor has the shape of a square matrix
isSquare :: Unbox a => Holor a -> Bool
{-# INLINE isSquare #-}
isSquare :: Holor a -> Bool
isSquare Holor a
hlr = ([Int] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Int]
shp Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
2) Bool -> Bool -> Bool
&& ([Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and ([Bool] -> Bool) -> [Bool] -> Bool
forall a b. (a -> b) -> a -> b
$ (Int -> Bool) -> [Int] -> [Bool]
forall a b. (a -> b) -> [a] -> [b]
P.map (Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== [Int] -> Int
forall a. [a] -> a
head [Int]
shp) ([Int] -> [Int]
forall a. [a] -> [a]
tail [Int]
shp)) where
    shp :: [Int]
shp = Holor a -> [Int]
forall a. Unbox a => Holor a -> [Int]
shape Holor a
hlr

-- | Tests if a matrix A with shape [m,n] is of upper triangular form. If __m<n__, false is returned. If __m>=n__, all elements on the first diagonal must be unequal to zero and all elements below the first diagonal must be zero.
isUpperTri :: (Unbox a, Eq a, Num a) => Holor a -> Bool
{-# INLINE isUpperTri #-}
isUpperTri :: Holor a -> Bool
isUpperTri Holor a
_A = case Holor a -> [Int]
forall a. Holor a -> [Int]
hShape Holor a
_A of
    [Int
m,Int
n] -> case Int
mInt -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>=Int
n of
        Bool
False -> Bool
False
        Bool
True -> (a -> Bool) -> [a] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (a -> a -> Bool
forall a. Eq a => a -> a -> Bool
==a
0) [a]
belowDiagonalEntries Bool -> Bool -> Bool
&& (a -> Bool) -> [a] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/=a
0) [a]
diagonalEntries
        where
            diagonalEntries :: [a]
diagonalEntries = [ Holor a
_AHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|![Int
i,Int
i] | Int
i <- [Int
0..Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]]
            belowDiagonalEntries :: [a]
belowDiagonalEntries = [ Holor a
_AHolor a -> [Int] -> a
forall a. Unbox a => Holor a -> [Int] -> a
|![Int
i,Int
j] | Int
j<-[Int
0..Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1], Int
i<-[Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1..Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1]]
    [Int]
_ -> Bool
False