Skip to content
This repository has been archived by the owner on Nov 18, 2023. It is now read-only.

Add a matrix type using an Unboxed vector underneath. #56

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 10 additions & 4 deletions Data/Matrix.hs
Original file line number Diff line number Diff line change
Expand Up @@ -621,7 +621,7 @@ ref :: (Fractional a, Eq a) => Matrix a -> Either String (Matrix a)
ref mtx
| nrows mtx == 1
= clearedLeft
| otherwise = do
| otherwise = do
(tl, tr, bl, br) <- (splitBlocks 1 1 <$> clearedLeft)
br' <- ref br
return ((tl <|> tr) <-> (bl <|> br'))
Expand All @@ -633,7 +633,7 @@ ref mtx
Nothing -> Left "Attempt to invert a non-invertible matrix"
Just x -> return x
normalizedFirstRow = (\sigAtTop' -> scaleRow (1 / getElem 1 1 sigAtTop') 1 sigAtTop') <$> sigAtTop
clearedLeft = do
clearedLeft = do
comb <- mapM combinator [2..nrows mtx]
firstRow <- normalizedFirstRow
return $ (foldr (.) id comb) firstRow
Expand Down Expand Up @@ -857,9 +857,9 @@ Four methods are provided for matrix multiplication.

* 'multStrassen':
Matrix multiplication following the Strassen's algorithm.
Complexity grows slower but also some work is added
Complexity grows slower but memory requirement is @O(n^2)@, and some work is added
partitioning the matrix. Also, it only works on square
matrices of order @2^n@, so if this condition is not
matrices of order @2^k@, so if this condition is not
met, it is zero-padded until this is accomplished.
Therefore, its use is not recommended.

Expand Down Expand Up @@ -956,6 +956,9 @@ first f = go

-- | Strassen's algorithm over square matrices of order @2^n@.
strassen :: Num a => Matrix a -> Matrix a -> Matrix a
{-# SPECIALIZE strassen :: Matrix Double -> Matrix Double -> Matrix Double #-}
{-# SPECIALIZE strassen :: Matrix Int -> Matrix Int -> Matrix Int #-}
{-# SPECIALIZE strassen :: Matrix Rational -> Matrix Rational -> Matrix Rational #-}
-- Trivial 1x1 multiplication.
strassen a@(M 1 1 _ _ _ _) b@(M 1 1 _ _ _ _) = M 1 1 0 0 1 $ V.singleton $ (a ! (1,1)) * (b ! (1,1))
-- General case guesses that the input matrices are square matrices
Expand Down Expand Up @@ -983,6 +986,9 @@ strassen a b = joinBlocks (c11,c12,c21,c22)

-- | Strassen's matrix multiplication.
multStrassen :: Num a => Matrix a -> Matrix a -> Matrix a
{-# SPECIALIZE multStrassen :: Matrix Double -> Matrix Double -> Matrix Double #-}
{-# SPECIALIZE multStrassen :: Matrix Int -> Matrix Int -> Matrix Int #-}
{-# SPECIALIZE multStrassen :: Matrix Rational -> Matrix Rational -> Matrix Rational #-}
multStrassen a1@(M n m _ _ _ _) a2@(M n' m' _ _ _ _)
| m /= n' = error $ "Multiplication of " ++ sizeStr n m ++ " and "
++ sizeStr n' m' ++ " matrices."
Expand Down
Loading