{-# LANGUAGE FlexibleContexts #-}

{-# OPTIONS_GHC -fno-warn-missing-signatures #-}

-----------------------------------------------------------------------------
-- |
-- Module      :  Internal.Chain
-- Copyright   :  (c) Vivian McPhail 2010
-- License     :  BSD3
--
-- Maintainer  :  Vivian McPhail <haskell.vivian.mcphail <at> gmail.com>
-- Stability   :  provisional
-- Portability :  portable
--
-- optimisation of association order for chains of matrix multiplication
--
-----------------------------------------------------------------------------

{-# LANGUAGE FlexibleContexts #-}

module Internal.Chain (
                      optimiseMult,
                     ) where

import Data.Maybe

import Internal.Matrix
import Internal.Numeric

import qualified Data.Array.IArray as A

-----------------------------------------------------------------------------
{- | 
     Provide optimal association order for a chain of matrix multiplications 
     and apply the multiplications.

     The algorithm is the well-known O(n\^3) dynamic programming algorithm
     that builds a pyramid of optimal associations.

> m1, m2, m3, m4 :: Matrix Double
> m1 = (10><15) [1..]
> m2 = (15><20) [1..]
> m3 = (20><5) [1..]
> m4 = (5><10) [1..]

> >>> optimiseMult [m1,m2,m3,m4]

will perform @((m1 `multiply` (m2 `multiply` m3)) `multiply` m4)@

The naive left-to-right multiplication would take @4500@ scalar multiplications
whereas the optimised version performs @2750@ scalar multiplications.  The complexity
in this case is 32 (= 4^3/2) * (2 comparisons, 3 scalar multiplications, 3 scalar additions,
5 lookups, 2 updates) + a constant (= three table allocations)
-}
optimiseMult :: Product t => [Matrix t] -> Matrix t
optimiseMult :: [Matrix t] -> Matrix t
optimiseMult = [Matrix t] -> Matrix t
forall a. Product a => [Matrix a] -> Matrix a
chain

-----------------------------------------------------------------------------

type Matrices a = A.Array Int (Matrix a)
type Sizes      = A.Array Int (Int,Int)
type Cost       = A.Array Int (A.Array Int (Maybe Int))
type Indexes    = A.Array Int (A.Array Int (Maybe ((Int,Int),(Int,Int))))

update :: A.Array Int (A.Array Int a) -> (Int,Int) -> a -> A.Array Int (A.Array Int a)
update :: Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int a)
a (Int
r,Int
c) a
e = Array Int (Array Int a)
a Array Int (Array Int a)
-> [(Int, Array Int a)] -> Array Int (Array Int a)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> [(i, e)] -> a i e
A.// [(Int
r,(Array Int (Array Int a)
a Array Int (Array Int a) -> Int -> Array Int a
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
r) Array Int a -> [(Int, a)] -> Array Int a
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> [(i, e)] -> a i e
A.// [(Int
c,a
e)])]

newWorkSpaceCost :: Int -> A.Array Int (A.Array Int (Maybe Int))
newWorkSpaceCost :: Int -> Array Int (Array Int (Maybe Int))
newWorkSpaceCost Int
n = (Int, Int)
-> [(Int, Array Int (Maybe Int))]
-> Array Int (Array Int (Maybe Int))
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [(i, e)] -> a i e
A.array (Int
1,Int
n) ([(Int, Array Int (Maybe Int))]
 -> Array Int (Array Int (Maybe Int)))
-> [(Int, Array Int (Maybe Int))]
-> Array Int (Array Int (Maybe Int))
forall a b. (a -> b) -> a -> b
$ (Int -> (Int, Array Int (Maybe Int)))
-> [Int] -> [(Int, Array Int (Maybe Int))]
forall a b. (a -> b) -> [a] -> [b]
map (\Int
i -> (Int
i, Int -> Array Int (Maybe Int)
forall (a :: * -> * -> *) a i.
(IArray a (Maybe a), Ix i, Num i) =>
i -> a i (Maybe a)
subArray Int
i)) [Int
1..Int
n]
   where subArray :: i -> a i (Maybe a)
subArray i
i = (i, i) -> [Maybe a] -> a i (Maybe a)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [e] -> a i e
A.listArray (i
1,i
i) (Maybe a -> [Maybe a]
forall a. a -> [a]
repeat Maybe a
forall a. Maybe a
Nothing)

newWorkSpaceIndexes :: Int -> A.Array Int (A.Array Int (Maybe ((Int,Int),(Int,Int))))
newWorkSpaceIndexes :: Int -> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
newWorkSpaceIndexes Int
n = (Int, Int)
-> [(Int, Array Int (Maybe ((Int, Int), (Int, Int))))]
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [(i, e)] -> a i e
A.array (Int
1,Int
n) ([(Int, Array Int (Maybe ((Int, Int), (Int, Int))))]
 -> Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> [(Int, Array Int (Maybe ((Int, Int), (Int, Int))))]
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
forall a b. (a -> b) -> a -> b
$ (Int -> (Int, Array Int (Maybe ((Int, Int), (Int, Int)))))
-> [Int] -> [(Int, Array Int (Maybe ((Int, Int), (Int, Int))))]
forall a b. (a -> b) -> [a] -> [b]
map (\Int
i -> (Int
i, Int -> Array Int (Maybe ((Int, Int), (Int, Int)))
forall (a :: * -> * -> *) a i.
(IArray a (Maybe a), Ix i, Num i) =>
i -> a i (Maybe a)
subArray Int
i)) [Int
1..Int
n]
   where subArray :: i -> a i (Maybe a)
subArray i
i = (i, i) -> [Maybe a] -> a i (Maybe a)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [e] -> a i e
A.listArray (i
1,i
i) (Maybe a -> [Maybe a]
forall a. a -> [a]
repeat Maybe a
forall a. Maybe a
Nothing)

matricesToSizes :: [Matrix a] -> Sizes
matricesToSizes :: [Matrix a] -> Sizes
matricesToSizes [Matrix a]
ms = (Int, Int) -> [(Int, Int)] -> Sizes
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [e] -> a i e
A.listArray (Int
1,[Matrix a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Matrix a]
ms) ([(Int, Int)] -> Sizes) -> [(Int, Int)] -> Sizes
forall a b. (a -> b) -> a -> b
$ (Matrix a -> (Int, Int)) -> [Matrix a] -> [(Int, Int)]
forall a b. (a -> b) -> [a] -> [b]
map (\Matrix a
m -> (Matrix a -> Int
forall t. Matrix t -> Int
rows Matrix a
m,Matrix a -> Int
forall t. Matrix t -> Int
cols Matrix a
m)) [Matrix a]
ms

chain :: Product a => [Matrix a] -> Matrix a
chain :: [Matrix a] -> Matrix a
chain []  = [Char] -> Matrix a
forall a. HasCallStack => [Char] -> a
error [Char]
"chain: zero matrices to multiply"
chain [Matrix a
m] = Matrix a
m
chain [Matrix a
ml,Matrix a
mr] = Matrix a
ml Matrix a -> Matrix a -> Matrix a
forall e. Product e => Matrix e -> Matrix e -> Matrix e
`multiply` Matrix a
mr
chain [Matrix a]
ms = let ln :: Int
ln = [Matrix a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Matrix a]
ms
               ma :: Array Int (Matrix a)
ma = (Int, Int) -> [Matrix a] -> Array Int (Matrix a)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
(i, i) -> [e] -> a i e
A.listArray (Int
1,Int
ln) [Matrix a]
ms
               mz :: Sizes
mz = [Matrix a] -> Sizes
forall a. [Matrix a] -> Sizes
matricesToSizes [Matrix a]
ms
               i :: Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
i = Sizes -> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
chain_cost Sizes
mz
           in (Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Array Int (Matrix a)
-> Matrix a
forall a.
Product a =>
(Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
chain_paren (Int
ln,Int
ln) Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
i Array Int (Matrix a)
ma

chain_cost :: Sizes -> Indexes
chain_cost :: Sizes -> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
chain_cost Sizes
mz = let (Int
_,Int
u) = Sizes -> (Int, Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> (i, i)
A.bounds Sizes
mz
                    cost :: Array Int (Array Int (Maybe Int))
cost = Int -> Array Int (Array Int (Maybe Int))
newWorkSpaceCost Int
u
                    ixes :: Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes = Int -> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
newWorkSpaceIndexes Int
u
                    (Sizes
_,Array Int (Array Int (Maybe Int))
_,Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
i) =  ((Sizes, Array Int (Array Int (Maybe Int)),
  Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
 -> (Int, Int)
 -> (Sizes, Array Int (Array Int (Maybe Int)),
     Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))))
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> [(Int, Int)]
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> (Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
chain_cost' (Sizes
mz,Array Int (Array Int (Maybe Int))
cost,Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes) (Int -> [(Int, Int)]
forall b. (Eq b, Num b, Enum b) => b -> [(b, b)]
order Int
u)
                in Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
i

chain_cost' :: (Sizes,Cost,Indexes) -> (Int,Int) -> (Sizes,Cost,Indexes)
chain_cost' :: (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> (Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
chain_cost' sci :: (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
sci@(Sizes
mz,Array Int (Array Int (Maybe Int))
cost,Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes) (Int
r,Int
c) 
    | Int
c Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1                     = let cost' :: Array Int (Array Int (Maybe Int))
cost' = Array Int (Array Int (Maybe Int))
-> (Int, Int) -> Maybe Int -> Array Int (Array Int (Maybe Int))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe Int))
cost (Int
r,Int
c) (Int -> Maybe Int
forall a. a -> Maybe a
Just Int
0)
                                       ixes' :: Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes' = Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> (Int, Int)
-> Maybe ((Int, Int), (Int, Int))
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes (Int
r,Int
c) (((Int, Int), (Int, Int)) -> Maybe ((Int, Int), (Int, Int))
forall a. a -> Maybe a
Just ((Int
r,Int
c),(Int
r,Int
c)))
                                       in (Sizes
mz,Array Int (Array Int (Maybe Int))
cost',Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes')
    | Bool
otherwise                  = (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> (Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
minimum_cost (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
sci (Int
r,Int
c)

minimum_cost :: (Sizes,Cost,Indexes) -> (Int,Int) -> (Sizes,Cost,Indexes)
minimum_cost :: (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> (Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
minimum_cost (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
sci (Int, Int)
fu = ((Sizes, Array Int (Array Int (Maybe Int)),
  Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
 -> ((Int, Int), (Int, Int))
 -> (Sizes, Array Int (Array Int (Maybe Int)),
     Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))))
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> [((Int, Int), (Int, Int))]
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl ((Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> ((Int, Int), (Int, Int))
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
smaller_cost (Int, Int)
fu) (Sizes, Array Int (Array Int (Maybe Int)),
 Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
sci ((Int, Int) -> [((Int, Int), (Int, Int))]
forall a. (Num a, Enum a) => (a, a) -> [((a, a), (a, a))]
fulcrum_order (Int, Int)
fu)

smaller_cost :: (Int,Int) -> (Sizes,Cost,Indexes) -> ((Int,Int),(Int,Int)) -> (Sizes,Cost,Indexes)
smaller_cost :: (Int, Int)
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
-> ((Int, Int), (Int, Int))
-> (Sizes, Array Int (Array Int (Maybe Int)),
    Array Int (Array Int (Maybe ((Int, Int), (Int, Int)))))
smaller_cost (Int
r,Int
c) (Sizes
mz,Array Int (Array Int (Maybe Int))
cost,Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes) ix :: ((Int, Int), (Int, Int))
ix@((Int
lr,Int
lc),(Int
rr,Int
rc)) =
    let op_cost :: Int
op_cost =   Maybe Int -> Int
forall a. HasCallStack => Maybe a -> a
fromJust ((Array Int (Array Int (Maybe Int))
cost Array Int (Array Int (Maybe Int)) -> Int -> Array Int (Maybe Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
lr) Array Int (Maybe Int) -> Int -> Maybe Int
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
lc)
               Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Maybe Int -> Int
forall a. HasCallStack => Maybe a -> a
fromJust ((Array Int (Array Int (Maybe Int))
cost Array Int (Array Int (Maybe Int)) -> Int -> Array Int (Maybe Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
rr) Array Int (Maybe Int) -> Int -> Maybe Int
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
rc)
               Int -> Int -> Int
forall a. Num a => a -> a -> a
+ (Int, Int) -> Int
forall a b. (a, b) -> a
fst (Sizes
mz Sizes -> Int -> (Int, Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! (Int
lrInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
lcInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1))
                 Int -> Int -> Int
forall a. Num a => a -> a -> a
* (Int, Int) -> Int
forall a b. (a, b) -> b
snd (Sizes
mz Sizes -> Int -> (Int, Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
lc)
                 Int -> Int -> Int
forall a. Num a => a -> a -> a
* (Int, Int) -> Int
forall a b. (a, b) -> b
snd (Sizes
mz Sizes -> Int -> (Int, Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
rr)
        cost' :: Maybe Int
cost' = (Array Int (Array Int (Maybe Int))
cost Array Int (Array Int (Maybe Int)) -> Int -> Array Int (Maybe Int)
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
r) Array Int (Maybe Int) -> Int -> Maybe Int
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
c
    in case Maybe Int
cost' of
               Maybe Int
Nothing -> let cost'' :: Array Int (Array Int (Maybe Int))
cost'' = Array Int (Array Int (Maybe Int))
-> (Int, Int) -> Maybe Int -> Array Int (Array Int (Maybe Int))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe Int))
cost (Int
r,Int
c) (Int -> Maybe Int
forall a. a -> Maybe a
Just Int
op_cost)
                              ixes'' :: Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes'' = Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> (Int, Int)
-> Maybe ((Int, Int), (Int, Int))
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes (Int
r,Int
c) (((Int, Int), (Int, Int)) -> Maybe ((Int, Int), (Int, Int))
forall a. a -> Maybe a
Just ((Int, Int), (Int, Int))
ix)
                          in (Sizes
mz,Array Int (Array Int (Maybe Int))
cost'',Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes'')
               Just Int
ct -> if Int
op_cost Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
ct then
                          let cost'' :: Array Int (Array Int (Maybe Int))
cost'' = Array Int (Array Int (Maybe Int))
-> (Int, Int) -> Maybe Int -> Array Int (Array Int (Maybe Int))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe Int))
cost (Int
r,Int
c) (Int -> Maybe Int
forall a. a -> Maybe a
Just Int
op_cost)
                              ixes'' :: Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes'' = Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> (Int, Int)
-> Maybe ((Int, Int), (Int, Int))
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
forall a.
Array Int (Array Int a)
-> (Int, Int) -> a -> Array Int (Array Int a)
update Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes (Int
r,Int
c) (((Int, Int), (Int, Int)) -> Maybe ((Int, Int), (Int, Int))
forall a. a -> Maybe a
Just ((Int, Int), (Int, Int))
ix)
                          in (Sizes
mz,Array Int (Array Int (Maybe Int))
cost'',Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes'')
                          else (Sizes
mz,Array Int (Array Int (Maybe Int))
cost,Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes)
                                                                         

fulcrum_order :: (a, a) -> [((a, a), (a, a))]
fulcrum_order (a
r,a
c) = let fs' :: [(a, a)]
fs' = [a] -> [a] -> [(a, a)]
forall a b. [a] -> [b] -> [(a, b)]
zip (a -> [a]
forall a. a -> [a]
repeat a
r) [a
1..(a
ca -> a -> a
forall a. Num a => a -> a -> a
-a
1)]
                      in ((a, a) -> ((a, a), (a, a))) -> [(a, a)] -> [((a, a), (a, a))]
forall a b. (a -> b) -> [a] -> [b]
map ((a, a) -> (a, a) -> ((a, a), (a, a))
forall b a. Num b => (b, b) -> (a, b) -> ((b, b), (a, b))
partner (a
r,a
c)) [(a, a)]
fs'

partner :: (b, b) -> (a, b) -> ((b, b), (a, b))
partner (b
r,b
c) (a
a,b
b) = ((b
rb -> b -> b
forall a. Num a => a -> a -> a
-b
b, b
cb -> b -> b
forall a. Num a => a -> a -> a
-b
b), (a
a,b
b))

order :: b -> [(b, b)]
order b
0 = []
order b
n = b -> [(b, b)]
order (b
nb -> b -> b
forall a. Num a => a -> a -> a
-b
1) [(b, b)] -> [(b, b)] -> [(b, b)]
forall a. [a] -> [a] -> [a]
++ [b] -> [b] -> [(b, b)]
forall a b. [a] -> [b] -> [(a, b)]
zip (b -> [b]
forall a. a -> [a]
repeat b
n) [b
1..b
n]

chain_paren :: Product a => (Int,Int) -> Indexes -> Matrices a -> Matrix a
chain_paren :: (Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
chain_paren (Int
r,Int
c) Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes Matrices a
ma = let ((Int
lr,Int
lc),(Int
rr,Int
rc)) = Maybe ((Int, Int), (Int, Int)) -> ((Int, Int), (Int, Int))
forall a. HasCallStack => Maybe a -> a
fromJust (Maybe ((Int, Int), (Int, Int)) -> ((Int, Int), (Int, Int)))
-> Maybe ((Int, Int), (Int, Int)) -> ((Int, Int), (Int, Int))
forall a b. (a -> b) -> a -> b
$ (Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Int -> Array Int (Maybe ((Int, Int), (Int, Int)))
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
r) Array Int (Maybe ((Int, Int), (Int, Int)))
-> Int -> Maybe ((Int, Int), (Int, Int))
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
c
                            in if Int
lr Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
rr Bool -> Bool -> Bool
&& Int
lc Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
rc then (Matrices a
ma Matrices a -> Int -> Matrix a
forall (a :: * -> * -> *) e i.
(IArray a e, Ix i) =>
a i e -> i -> e
A.! Int
lr)
                               else ((Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
forall a.
Product a =>
(Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
chain_paren (Int
lr,Int
lc) Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes Matrices a
ma) Matrix a -> Matrix a -> Matrix a
forall e. Product e => Matrix e -> Matrix e -> Matrix e
`multiply` ((Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
forall a.
Product a =>
(Int, Int)
-> Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
-> Matrices a
-> Matrix a
chain_paren (Int
rr,Int
rc) Array Int (Array Int (Maybe ((Int, Int), (Int, Int))))
ixes Matrices a
ma) 

--------------------------------------------------------------------------

{- TESTS

-- optimal association is ((m1*(m2*m3))*m4)
m1, m2, m3, m4 :: Matrix Double
m1 = (10><15) [1..]
m2 = (15><20) [1..]
m3 = (20><5) [1..]
m4 = (5><10) [1..]

-}