@@ -0,0 +1,471 @@
1+{-# LANGUAGE AllowAmbiguousTypes #-}
2+{-# LANGUAGE DataKinds #-}
3+{-# LANGUAGE FlexibleContexts #-}
4+{-# LANGUAGE FlexibleInstances #-}
5+{-# LANGUAGE GADTs #-}
6+{-# LANGUAGE InstanceSigs #-}
7+{-# LANGUAGE MultiParamTypeClasses #-}
8+{-# LANGUAGE NoStarIsType #-}
9+{-# LANGUAGE ScopedTypeVariables #-}
10+{-# LANGUAGE StandaloneDeriving #-}
11+{-# LANGUAGE TypeApplications #-}
12+{-# LANGUAGE TypeFamilies #-}
13+{-# LANGUAGE TypeOperators #-}
14+{-# LANGUAGE UndecidableInstances #-}
15+
16+module Matrix.Type
17+ ( -- | LAoP (Linear Algebra of Programming) Inductive Matrix definition.
18+ --
19+ -- LAoP generalises relations and functions treating them as
20+ -- Boolean matrices and in turn consider these as arrows.
21+ -- This library offers many of the combinators mentioned in the work of
22+ -- Macedo (2012) and Oliveira (2012).
23+ --
24+ -- This definition makes use of the fact that 'Void' is
25+ -- isomorphic to 0 and '()' to 1 and captures matrix
26+ -- dimensions as stacks of 'Either's.
27+ --
28+ -- There exists two type families that make it easier to write
29+ -- matrix dimensions: 'FromNat' and 'Count'. This approach
30+ -- leads to a very straightforward implementation
31+ -- of LAoP combinators.
32+
33+ -- * Type safe matrix representation
34+ Matrix (..),
35+
36+ -- * Auxiliary type families
37+ FromNat,
38+ Count,
39+ Normalize,
40+
41+ -- * Matrix construction and conversion
42+ FromLists,
43+ fromLists,
44+ toLists,
45+ toList,
46+ matrixBuilder,
47+ row,
48+ col,
49+ zeros,
50+ ones,
51+ bang,
52+ constant,
53+
54+ -- * Misc
55+ -- ** Get dimensions
56+ columns,
57+ rows,
58+
59+ -- ** Matrix Transposition
60+ tr,
61+
62+ -- ** Matrix "abiding"
63+ abideJS,
64+ abideSJ,
65+
66+ -- * Biproduct approach
67+ -- ** Split
68+ split,
69+ -- *** Projections
70+ p1,
71+ p2,
72+ -- ** Junc
73+ junc,
74+ -- *** Injections
75+ i1,
76+ i2,
77+ -- ** Bifunctors
78+ (-|-),
79+ (><),
80+
81+ -- ** Applicative matrix combinators
82+
83+ -- | Note that given the restrictions imposed it is not possible to
84+ -- implement the standard type classes present in standard Haskell.
85+ -- *** Matrix pairing projections
86+ KhatriP1,
87+ kp1,
88+ KhatriP2,
89+ kp2,
90+
91+ -- *** Matrix pairing
92+ khatri,
93+
94+ -- * Matrix composition and lifting
95+
96+ -- ** Arrow matrix combinators
97+
98+ -- | Note that given the restrictions imposed it is not possible to
99+ -- implement the standard type classes present in standard Haskell.
100+ identity,
101+ comp,
102+ fromF,
103+
104+ -- * Matrix printing
105+ pretty,
106+ prettyPrint
107+ )
108+ where
109+
110+import Data.Bool
111+import Data.Kind
112+import Data.List
113+import Data.Proxy
114+import Data.Void
115+import GHC.TypeLits
116+
117+-- | Matrix data type
118+data Matrix e cols rows where
119+ Empty :: Matrix e Void Void
120+ One :: e -> Matrix e () ()
121+ Junc :: Matrix e a rows -> Matrix e b rows -> Matrix e (Either a b) rows
122+ Split :: Matrix e cols a -> Matrix e cols b -> Matrix e cols (Either a b)
123+
124+deriving instance (Show e) => Show (Matrix e cols rows)
125+
126+-- | Type family that computes the cardinality of a given type dimension
127+type family Count (d :: Type) where
128+ Count Void = 0
129+ Count () = 1
130+ Count (Either a b) = (+) (Count a) (Count b)
131+ Count (a, b) = (*) (Count a) (Count b)
132+ Count (a -> b) = (^) (Count b) (Count a)
133+
134+-- | Type family that computes of a given type dimension from a fiven natural
135+type family FromNat (n :: Nat) where
136+ FromNat 0 = Void
137+ FromNat 1 = ()
138+ FromNat n = Either () (FromNat (n - 1))
139+
140+-- | Type family that normalizes the representation of a given data
141+-- structure
142+type family Normalize (d :: Type) where
143+ Normalize d = FromNat (Count d)
144+
145+instance Eq e => Eq (Matrix e cols rows) where
146+ Empty == Empty = True
147+ (One a) == (One b) = a == b
148+ (Junc a b) == (Junc c d) = a == c && b == d
149+ (Split a b) == (Split c d) = a == c && b == d
150+
151+instance Num e => Num (Matrix e cols rows) where
152+
153+ Empty + Empty = Empty
154+ (One a) + (One b) = One (a + b)
155+ (Junc a b) + (Junc c d) = Junc (a + c) (b + d)
156+ (Split a b) + (Split c d) = Split (a + c) (b + d)
157+
158+ Empty - Empty = Empty
159+ (One a) - (One b) = One (a - b)
160+ (Junc a b) - (Junc c d) = Junc (a - c) (b - d)
161+ (Split a b) - (Split c d) = Split (a - c) (b - d)
162+
163+ Empty * Empty = Empty
164+ (One a) * (One b) = One (a * b)
165+ (Junc a b) * (Junc c d) = Junc (a * c) (b * d)
166+ (Split a b) * (Split c d) = Split (a * c) (b * d)
167+
168+ abs Empty = Empty
169+ abs (One a) = One (abs a)
170+ abs (Junc a b) = Junc (abs a) (abs b)
171+ abs (Split a b) = Split (abs a) (abs b)
172+
173+ signum Empty = Empty
174+ signum (One a) = One (signum a)
175+ signum (Junc a b) = Junc (signum a) (signum b)
176+ signum (Split a b) = Split (signum a) (signum b)
177+
178+-- Primitives
179+
180+empty :: Matrix e Void Void
181+empty = Empty
182+
183+one :: e -> Matrix e () ()
184+one = One
185+
186+junc :: Matrix e a rows -> Matrix e b rows -> Matrix e (Either a b) rows
187+junc = Junc
188+
189+split :: Matrix e cols a -> Matrix e cols b -> Matrix e cols (Either a b)
190+split = Split
191+
192+-- Construction
193+
194+class FromLists e cols rows where
195+ fromLists :: [[e]] -> Matrix e cols rows
196+
197+instance FromLists e Void Void where
198+ fromLists [] = Empty
199+ fromLists _ = error "Wrong dimensions"
200+
201+instance {-# OVERLAPPING #-} FromLists e () () where
202+ fromLists [[e]] = One e
203+ fromLists _ = error "Wrong dimensions"
204+
205+instance (FromLists e cols ()) => FromLists e (Either () cols) () where
206+ fromLists [(h : t)] = Junc (One h) (fromLists [t])
207+ fromLists _ = error "Wrong dimensions"
208+
209+instance (FromLists e () rows) => FromLists e () (Either () rows) where
210+ fromLists ([h] : t) = Split (One h) (fromLists t)
211+ fromLists _ = error "Wrong dimensions"
212+
213+instance {-# OVERLAPPABLE #-} (FromLists e cols c, FromLists e cols d) => FromLists e cols (Either c d) where
214+ fromLists (h : t) =
215+ let lh = length h
216+ condition = all (== lh) (map length t)
217+ in if lh > 0 && condition
218+ then Split (fromLists [h]) (fromLists t)
219+ else error "Not all rows have the same length"
220+
221+matrixBuilder ::
222+ forall e cols rows.
223+ ( FromLists e cols rows,
224+ KnownNat (Count cols),
225+ KnownNat (Count rows)
226+ ) =>
227+ ((Int, Int) -> e) ->
228+ Matrix e cols rows
229+matrixBuilder f =
230+ let c = fromInteger $ natVal (Proxy :: Proxy (Count cols))
231+ r = fromInteger $ natVal (Proxy :: Proxy (Count rows))
232+ positions = [(a, b) | a <- [0 .. (r - 1)], b <- [0 .. (c - 1)]]
233+ in fromLists . map (map f) . groupBy (\(x, _) (w, _) -> x == w) $ positions
234+
235+col :: (FromLists e () rows) => [e] -> Matrix e () rows
236+col = fromLists . map (: [])
237+
238+row :: (FromLists e cols ()) => [e] -> Matrix e cols ()
239+row = fromLists . (: [])
240+
241+fromF ::
242+ forall a b cols rows e.
243+ ( Bounded a,
244+ Bounded b,
245+ Enum a,
246+ Enum b,
247+ Eq b,
248+ Num e,
249+ Ord e,
250+ KnownNat (Count cols),
251+ KnownNat (Count rows),
252+ FromLists e rows cols
253+ ) =>
254+ (a -> b) ->
255+ Matrix e cols rows
256+fromF f =
257+ let minA = minBound @a
258+ maxA = maxBound @a
259+ minB = minBound @b
260+ maxB = maxBound @b
261+ ccols = fromInteger $ natVal (Proxy :: Proxy (Count cols))
262+ rrows = fromInteger $ natVal (Proxy :: Proxy (Count rows))
263+ elementsA = take ccols [minA .. maxA]
264+ elementsB = take rrows [minB .. maxB]
265+ combinations = (,) <$> elementsA <*> elementsB
266+ combAp = map snd . sort . map (\(a, b) -> if f a == b
267+ then ((fromEnum a, fromEnum b), 1)
268+ else ((fromEnum a, fromEnum b), 0)) $ combinations
269+ mList = buildList combAp rrows
270+ in tr $ fromLists mList
271+ where
272+ buildList [] _ = []
273+ buildList l r = take r l : buildList (drop r l) r
274+
275+-- Conversion
276+
277+toLists :: Matrix e cols rows -> [[e]]
278+toLists Empty = []
279+toLists (One e) = [[e]]
280+toLists (Split l r) = toLists l ++ toLists r
281+toLists (Junc l r) = zipWith (++) (toLists l) (toLists r)
282+
283+toList :: Matrix e cols rows -> [e]
284+toList = concat . toLists
285+
286+-- Zeros Matrix
287+
288+zeros :: (Num e, FromLists e cols rows, KnownNat (Count cols), KnownNat (Count rows)) => Matrix e cols rows
289+zeros = matrixBuilder (const 0)
290+
291+-- Ones Matrix
292+
293+ones :: (Num e, FromLists e cols rows, KnownNat (Count cols), KnownNat (Count rows)) => Matrix e cols rows
294+ones = matrixBuilder (const 1)
295+
296+-- Const Matrix
297+
298+constant :: (Num e, FromLists e cols rows, KnownNat (Count cols), KnownNat (Count rows)) => e -> Matrix e cols rows
299+constant e = matrixBuilder (const e)
300+
301+-- Bang Matrix
302+
303+bang :: forall e cols. (Num e, Enum e, FromLists e cols (), KnownNat (Count cols)) => Matrix e cols ()
304+bang =
305+ let c = fromInteger $ natVal (Proxy :: Proxy (Count cols))
306+ in fromLists [take c [1, 1 ..]]
307+
308+-- Identity Matrix
309+
310+identity :: (Num e, FromLists e cols cols, KnownNat (Count cols)) => Matrix e cols cols
311+identity = matrixBuilder (bool 0 1 . uncurry (==))
312+
313+-- Matrix composition (MMM)
314+
315+comp :: (Num e) => Matrix e cr rows -> Matrix e cols cr -> Matrix e cols rows
316+comp Empty Empty = Empty
317+comp (One a) (One b) = One (a * b)
318+comp (Junc a b) (Split c d) = comp a c + comp b d
319+comp (Split a b) c = Split (comp a c) (comp b c)
320+comp c (Junc a b) = Junc (comp c a) (comp c b)
321+
322+-- Projections
323+
324+p1 :: forall e m n. (Num e, KnownNat (Count n), KnownNat (Count m), FromLists e n m, FromLists e m m) => Matrix e (Either m n) m
325+p1 =
326+ let iden = identity :: Matrix e m m
327+ zero = zeros :: Matrix e n m
328+ in junc iden zero
329+
330+p2 :: forall e m n. (Num e, KnownNat (Count n), KnownNat (Count m), FromLists e m n, FromLists e n n) => Matrix e (Either m n) n
331+p2 =
332+ let iden = identity :: Matrix e n n
333+ zero = zeros :: Matrix e m n
334+ in junc zero iden
335+
336+-- Injections
337+
338+i1 :: (Num e, KnownNat (Count n), KnownNat (Count m), FromLists e n m, FromLists e m m) => Matrix e m (Either m n)
339+i1 = tr p1
340+
341+i2 :: (Num e, KnownNat (Count n), KnownNat (Count m), FromLists e m n, FromLists e n n) => Matrix e n (Either m n)
342+i2 = tr p2
343+
344+-- Dimensions
345+
346+rows :: forall e cols rows. (KnownNat (Count rows)) => Matrix e cols rows -> Int
347+rows _ = fromInteger $ natVal (Proxy :: Proxy (Count rows))
348+
349+columns :: forall e cols rows. (KnownNat (Count cols)) => Matrix e cols rows -> Int
350+columns _ = fromInteger $ natVal (Proxy :: Proxy (Count cols))
351+
352+-- Coproduct Bifunctor
353+
354+infixl 5 -|-
355+
356+(-|-) ::
357+ forall e n k m j.
358+ ( Num e,
359+ KnownNat (Count j),
360+ KnownNat (Count k),
361+ FromLists e k k,
362+ FromLists e j k,
363+ FromLists e k j,
364+ FromLists e j j
365+ ) =>
366+ Matrix e n k ->
367+ Matrix e m j ->
368+ Matrix e (Either n m) (Either k j)
369+(-|-) a b = Junc (comp i1 a) (comp i2 b)
370+
371+-- Khatri Rao Product and projections
372+
373+class KhatriP1 e m k where
374+ kp1 :: Matrix e (Normalize (m, k)) m
375+
376+instance KhatriP1 e Void k where
377+ kp1 = Empty
378+
379+instance
380+ ( Num e,
381+ KnownNat (Count k),
382+ FromLists e (FromNat (Count m * Count k)) m,
383+ KnownNat (Count m),
384+ KnownNat (Count (Normalize (m, k)))
385+ ) =>
386+ KhatriP1 e m k
387+ where
388+ kp1 = matrixBuilder f
389+ where
390+ offset = fromInteger (natVal (Proxy :: Proxy (Count k)))
391+ f (x, y)
392+ | y >= (x * offset) && y <= (x * offset + offset - 1) = 1
393+ | otherwise = 0
394+
395+class KhatriP2 e m k where
396+ kp2 :: Matrix e (Normalize (m, k)) k
397+
398+instance KhatriP2 e m Void where
399+ kp2 = Empty
400+
401+instance
402+ ( Num e,
403+ KnownNat (Count k),
404+ FromLists e (FromNat (Count m * Count k)) k,
405+ KnownNat (Count m),
406+ KnownNat (Count (Normalize (m, k)))
407+ ) =>
408+ KhatriP2 e m k
409+ where
410+ kp2 = matrixBuilder f
411+ where
412+ offset = fromInteger (natVal (Proxy :: Proxy (Count k)))
413+ f (x, y)
414+ | x == y || mod (y - x) offset == 0 = 1
415+ | otherwise = 0
416+
417+khatri :: forall e cols a b. (Num e, KhatriP1 e a b, KhatriP2 e a b) => Matrix e cols a -> Matrix e cols b -> Matrix e cols (Normalize (a, b))
418+khatri a b =
419+ let kp1' = kp1 @e @a @b
420+ kp2' = kp2 @e @a @b
421+ in comp (tr kp1') a * comp (tr kp2') b
422+
423+-- Product Bifunctor (Kronecker)
424+
425+infixl 4 ><
426+
427+(><) :: forall e m p n q. (Num e, KhatriP1 e m n, KhatriP2 e m n, KhatriP1 e p q, KhatriP2 e p q) => Matrix e m p -> Matrix e n q -> Matrix e (Normalize (m, n)) (Normalize (p, q))
428+(><) a b =
429+ let kp1' = kp1 @e @m @n
430+ kp2' = kp2 @e @m @n
431+ in khatri (comp a kp1') (comp b kp2')
432+
433+-- Matrix abide Junc Split
434+
435+abideJS :: Matrix e cols rows -> Matrix e cols rows
436+abideJS (Junc (Split a c) (Split b d)) = Split (Junc (abideJS a) (abideJS b)) (Junc (abideJS c) (abideJS d))
437+abideJS Empty = Empty
438+abideJS (One e) = (One e)
439+abideJS (Junc a b) = Junc (abideJS a) (abideJS b)
440+abideJS (Split a b) = Split (abideJS a) (abideJS b)
441+
442+-- Matrix abide Split Junc
443+
444+abideSJ :: Matrix e cols rows -> Matrix e cols rows
445+abideSJ (Split (Junc a b) (Junc c d)) = Junc (Split (abideSJ a) (abideSJ c)) (Split (abideSJ b) (abideSJ d))
446+abideSJ Empty = Empty
447+abideSJ (One e) = (One e)
448+abideSJ (Junc a b) = Junc (abideSJ a) (abideSJ b)
449+abideSJ (Split a b) = Split (abideSJ a) (abideSJ b)
450+
451+-- Matrix transposition
452+
453+tr :: Matrix e cols rows -> Matrix e rows cols
454+tr Empty = Empty
455+tr (One e) = One e
456+tr (Junc a b) = Split (tr a) (tr b)
457+tr (Split a b) = Junc (tr a) (tr b)
458+
459+-- Pretty print
460+
461+prettyAux :: Show e => [[e]] -> String
462+prettyAux [] = ""
463+prettyAux [[e]] = show e
464+prettyAux [h] = concat $ intersperse " " (map show h)
465+prettyAux (h : t) = concat (intersperse " " (map show h)) ++ "\n" ++ prettyAux t
466+
467+pretty :: Show e => Matrix e cols rows -> String
468+pretty = prettyAux . toLists
469+
470+prettyPrint :: Show e => Matrix e cols rows -> IO ()
471+prettyPrint = putStrLn . pretty