module MathObj.Matrix where
import qualified Algebra.Module as Module
import qualified Algebra.Vector as Vector
import qualified Algebra.Ring as Ring
import qualified Algebra.Additive as Additive
import Algebra.Module((*>))
import Algebra.Ring((*), fromInteger, scalarProduct)
import Algebra.Additive((+), (), zero, subtract)
import Data.Array (Array, listArray, elems, bounds, (!), ixmap, range)
import qualified Data.List as List
import Control.Monad (liftM2)
import Control.Exception (assert)
import NumericPrelude.List (outerProduct)
import NumericPrelude(Integer)
import PreludeBase hiding (zipWith)
data
T a = Cons (Array (Integer, Integer) a) deriving (Eq,Ord,Read)
twist :: (Integer,Integer) -> (Integer,Integer)
twist (x,y) = (y,x)
transpose :: T a -> T a
transpose (Cons m) =
let (lower,upper) = bounds m
in Cons (ixmap (twist lower, twist upper) twist m)
rows :: T a -> [[a]]
rows (Cons m) =
let ((lr,lc), (ur,uc)) = bounds m
in outerProduct (curry(m!)) (range (lr,ur)) (range (lc,uc))
columns :: T a -> [[a]]
columns (Cons m) =
let ((lr,lc), (ur,uc)) = bounds m
in outerProduct (curry(m!)) (range (lc,uc)) (range (lr,ur))
fromList :: Integer -> Integer -> [a] -> T a
fromList m n xs = Cons (listArray ((1,1),(m,n)) xs)
instance (Ring.C a, Show a) => Show (T a) where
show m = List.unlines $ map (\r -> "(" ++ r ++ ")")
$ map (unwords . map show) $ rows m
dimension :: T a -> (Integer,Integer)
dimension (Cons m) = uncurry subtract (bounds m) + (1,1)
numRows :: T a -> Integer
numRows = fst . dimension
numColumns :: T a -> Integer
numColumns = snd . dimension
instance (Additive.C a) => Additive.C (T a) where
(+) = zipWith (+)
() = zipWith ()
zero = zeroMatrix 1 1
zipWith :: (a -> b -> c) -> T a -> T b -> T c
zipWith op mM@(Cons m) nM@(Cons n) =
let d = dimension mM
em = elems m
en = elems n
in assert (d == dimension nM) $
uncurry fromList d (List.zipWith op em en)
zeroMatrix :: (Additive.C a) => Integer -> Integer -> T a
zeroMatrix m n = fromList m n $
List.repeat zero
instance (Ring.C a) => Ring.C (T a) where
mM * nM = assert (numRows mM == numColumns nM) $
fromList (numColumns mM) (numRows nM)
(liftM2 scalarProduct (rows mM) (columns nM))
fromInteger n = fromList 1 1 [fromInteger n]
instance Functor T where
fmap f (Cons m) = Cons (fmap f m)
instance Vector.C T where
zero = zero
(<+>) = (+)
(*>) = Vector.functorScale
instance Module.C a b => Module.C a (T b) where
x *> m = fmap (x*>) m
preimage :: (Ring.C a) => T a -> T a -> Maybe (T a)
preimage a y = assert
(numRows a == numRows y &&
numColumns y == 1)
Nothing