diff --git a/rings.cabal b/rings.cabal --- a/rings.cabal +++ b/rings.cabal @@ -1,7 +1,7 @@ name: rings -version: 0.0.3 +version: 0.0.3.1 synopsis: Ring-like objects. -description: Semirings, rings, division rings, modules, and algebras. +description: Semirings, rings, division rings, and modules. homepage: https://github.com/cmk/rings license: BSD3 license-file: LICENSE @@ -20,18 +20,16 @@ ghc-options: -Wall -optc-std=c99 exposed-modules: - Data.Algebra - , Data.Algebra.Quaternion - , Data.Semiring + Data.Semiring , Data.Semiring.Property , Data.Semifield - , Data.Semimodule - , Data.Semimodule.Vector - , Data.Semimodule.Matrix - , Data.Semimodule.Transform , Data.Semigroup.Additive , Data.Semigroup.Multiplicative , Data.Semigroup.Property + , Data.Semimodule + , Data.Semimodule.Free + , Data.Semimodule.Basis + , Data.Semimodule.Transform default-extensions: ScopedTypeVariables diff --git a/src/Data/Algebra.hs b/src/Data/Algebra.hs deleted file mode 100644 --- a/src/Data/Algebra.hs +++ /dev/null @@ -1,272 +0,0 @@ -{-# LANGUAGE CPP #-} -{-# LANGUAGE Safe #-} -{-# LANGUAGE PolyKinds #-} -{-# LANGUAGE ConstraintKinds #-} -{-# LANGUAGE DefaultSignatures #-} -{-# LANGUAGE DeriveFunctor #-} -{-# LANGUAGE DeriveGeneric #-} -{-# LANGUAGE FlexibleContexts #-} -{-# LANGUAGE FlexibleInstances #-} -{-# LANGUAGE NoImplicitPrelude #-} -{-# LANGUAGE RebindableSyntax #-} -{-# LANGUAGE TypeOperators #-} -{-# LANGUAGE TypeFamilies #-} - -module Data.Algebra ( - (><) - , (//) - , (.@.) - , unit - , norm - , conj - , triple - , reciprocal - , Algebra(..) - , Composition(..) - , Unital(..) - , Division(..) -) where - -import safe Data.Bool -import safe Data.Functor.Rep -import safe Data.Semifield -import safe Data.Semigroup.Additive as A -import safe Data.Semigroup.Multiplicative as M -import safe Data.Semimodule -import safe Data.Semiring hiding ((//)) -import safe Prelude hiding (Num(..), Fractional(..), sum, product) - --- | < https://en.wikipedia.org/wiki/Algebra_over_a_field#Generalization:_algebra_over_a_ring Algebra > over a semiring. --- --- Needn't be associative or unital. --- -class Semiring r => Algebra r a where - multiplyWith :: (a -> a -> r) -> a -> r - -infixl 7 >< - --- | Multiplication operator on a free algebra. --- --- In particular this is cross product on the 'I3' basis in /R^3/: --- --- >>> V3 1 0 0 >< V3 0 1 0 >< V3 0 1 0 :: V3 Int --- V3 (-1) 0 0 --- >>> V3 1 0 0 >< (V3 0 1 0 >< V3 0 1 0) :: V3 Int --- V3 0 0 0 --- --- /Caution/ in general (><) needn't be commutative, nor even associative. --- --- The cross product in particular satisfies the following properties: --- --- @ --- a '><' a = 'mempty' --- a '><' b = 'negate' ( b '><' a ) , --- a '><' ( b <> c ) = ( a '><' b ) <> ( a '><' c ) , --- ( r a ) '><' b = a '><' ( r b ) = r ( a '><' b ) . --- a '><' ( b '><' c ) <> b '><' ( c '><' a ) <> c '><' ( a '><' b ) = 'mempty' . --- @ --- --- See < https://en.wikipedia.org/wiki/Jacobi_identity Jacobi identity >. --- --- For associative algebras, use (*) instead for clarity: --- --- >>> (1 :+ 2) >< (3 :+ 4) :: Complex Int --- (-5) :+ 10 --- >>> (1 :+ 2) * (3 :+ 4) :: Complex Int --- (-5) :+ 10 --- >>> qi >< qj :: QuatM --- Quaternion 0.000000 (V3 0.000000 0.000000 1.000000) --- >>> qi * qj :: QuatM --- Quaternion 0.000000 (V3 0.000000 0.000000 1.000000) --- -(><) :: (Representable f, Algebra r (Rep f)) => f r -> f r -> f r -(><) x y = tabulate $ multiplyWith (\i1 i2 -> index x i1 * index y i2) - --- | Scalar triple product. --- --- @ --- 'triple' x y z = 'triple' z x y = 'triple' y z x --- 'triple' x y z = 'negate' '$' 'triple' x z y = 'negate' '$' 'triple' y x z --- 'triple' x x y = 'triple' x y y = 'triple' x y x = 'zero' --- ('triple' x y z) '*.' x = (x '><' y) '><' (x '><' z) --- @ --- --- >>> triple (V3 0 0 1) (V3 1 0 0) (V3 0 1 0) :: Double --- 1.0 --- -triple :: Free f => Foldable f => Algebra a (Rep f) => f a -> f a -> f a -> a -triple x y z = x .*. (y >< z) -{-# INLINE triple #-} - --- | < https://en.wikipedia.org/wiki/Composition_algebra Composition algebra > over a free semimodule. --- -class Algebra r a => Composition r a where - conjugateWith :: (a -> r) -> a -> r - - normWith :: (a -> r) -> r - - --- @ 'conj' a '><' 'conj' b = 'conj' (b >< a) @ ---prop_conj :: Representable f => Foldable f => Semiring b => Composition a (Rep f) => Rel a b -> f a -> f a -> b ---prop_conj (~~) p q = sum $ mzipWithRep (~~) (conj (p >< q)) (conj q >< conj p) - --- @ 'conj' a '><' 'conj' b = 'conj' (b >< a) @ -conj :: Representable f => Composition r (Rep f) => f r -> f r -conj = tabulate . conjugateWith . index - --- | Norm of a composition algebra. --- --- @ --- 'norm' x '*' 'norm' y = 'norm' (x >< y) --- 'norm' . 'norm'' $ x = 'norm' x '*' 'norm' x --- @ --- -norm :: (Representable f, Composition r (Rep f)) => f r -> r -norm x = normWith $ index x - ---norm' :: (Representable f, Composition r (Rep f)) => f r -> f r ---norm' x = x >< conj x - -class (Semiring r, Algebra r a) => Unital r a where - unitWith :: r -> a -> r - --- | Unit of a unital algebra. --- --- >>> unit :: Complex Int --- 1 :+ 0 --- >>> unit :: QuatD --- Quaternion 1.0 (V3 0.0 0.0 0.0) --- -unit :: Representable f => Unital r (Rep f) => f r -unit = tabulate $ unitWith one - --- | A (not necessarily associative) < https://en.wikipedia.org/wiki/Division_algebra division algebra >. --- -class (Semifield r, Unital r a) => Division r a where - --divideWith :: (a -> a -> r) -> a -> r - - reciprocalWith :: (a -> r) -> a -> r - - - - --- | @ 'reciprocal' x = (/ 'quadrance' x) '<$>' 'conj' x@ -reciprocal :: Representable f => Division a (Rep f) => f a -> f a -reciprocal = tabulate . reciprocalWith . index - --- reciprocal' x = (/ quadrance x) <$> conj x - - -infixl 7 // - --- | Division operator on a free division algebra. --- --- >>> (1 :+ 0) // (0 :+ 1) --- 0.0 :+ (-1.0) --- -(//) :: Representable f => Division r (Rep f) => f r -> f r -> f r -(//) x y = x >< reciprocal y - -infix 6 .@. --- | Bilinear form on a free composition algebra. --- --- >>> V2 1 2 .@. V2 1 2 --- 5.0 --- >>> V2 1 2 .@. V2 2 (-1) --- 0.0 --- >>> V3 1 1 1 .@. V3 1 1 (-2) --- 0.0 --- --- >>> (1 :+ 2) .@. (2 :+ (-1)) :: Double --- 0.0 --- --- >>> qi .@. qj :: Double --- 0.0 --- >>> qj .@. qk :: Double --- 0.0 --- >>> qk .@. qi :: Double --- 0.0 --- >>> qk .@. qk :: Double --- 1.0 --- -(.@.) :: Representable f => Composition a (Rep f) => Semigroup (f a) => Field a => f a -> f a -> a -x .@. y = prod / two where prod = norm (x <> y) - norm x - norm y - ---------------------------------------------------------------------- --- Instances ---------------------------------------------------------------------- - - ---instance (Semiring r, Unital r a) => Unital r (a -> r) where --- unitWith = unitWith one - ---instance (Semiring r, Division r a) => Division r (a -> r) where --- reciprocalWith = reciprocalWith - --- incoherent --- instance Unital () a where unitWith _ _ = () --- instance (Unital r a, Unital r b) => Unital (a -> r) b where unitWith f b a = unitWith (f a) b - -instance Semiring r => Algebra r () where - multiplyWith f = f () - -instance Semiring r => Unital r () where - unitWith r () = r - -instance (Algebra r a, Algebra r b) => Algebra r (a,b) where - multiplyWith f (a,b) = multiplyWith (\a1 a2 -> multiplyWith (\b1 b2 -> f (a1,b1) (a2,b2)) b) a - -instance (Algebra r a, Algebra r b, Algebra r c) => Algebra r (a,b,c) where - multiplyWith f (a,b,c) = multiplyWith (\a1 a2 -> multiplyWith (\b1 b2 -> multiplyWith (\c1 c2 -> f (a1,b1,c1) (a2,b2,c2)) c) b) a - -instance (Unital r a, Unital r b) => Unital r (a,b) where - unitWith r (a,b) = unitWith r a * unitWith r b - -instance (Unital r a, Unital r b, Unital r c) => Unital r (a,b,c) where - unitWith r (a,b,c) = unitWith r a * unitWith r b * unitWith r c - --- | Tensor algebra --- --- >>> multiplyWith (<>) [1..3 :: Int] --- [1,2,3,1,2,3,1,2,3,1,2,3] --- --- >>> multiplyWith (\f g -> fold (f ++ g)) [1..3] :: Int --- 24 --- -instance Semiring r => Algebra r [a] where - multiplyWith f = go [] where - go ls rrs@(r:rs) = f (reverse ls) rrs + go (r:ls) rs - go ls [] = f (reverse ls) [] - -instance Semiring r => Unital r [a] where - unitWith r [] = r - unitWith _ _ = zero - -type ComplexBasis = Bool - --- Complex basis ---instance Module r ComplexBasis => Algebra r ComplexBasis where -instance Ring r => Algebra r ComplexBasis where - multiplyWith f = f' where - fe = f False False - f True True - fi = f False True + f True False - f' False = fe - f' True = fi - ---instance Module r ComplexBasis => Composition r ComplexBasis where -instance Ring r => Composition r ComplexBasis where - conjugateWith f = f' where - afe = f False - nfi = negate (f True) - f' False = afe - f' True = nfi - - normWith f = flip multiplyWith zero $ \i1 i2 -> f i1 * conjugateWith f i2 - ---instance Module r ComplexBasis => Unital r ComplexBasis where -instance Ring r => Unital r ComplexBasis where - unitWith x False = x - unitWith _ _ = zero - -instance Field r => Division r ComplexBasis where - reciprocalWith f i = conjugateWith f i / normWith f diff --git a/src/Data/Algebra/Quaternion.hs b/src/Data/Algebra/Quaternion.hs deleted file mode 100644 --- a/src/Data/Algebra/Quaternion.hs +++ /dev/null @@ -1,247 +0,0 @@ -{-# LANGUAGE CPP #-} -{-# LANGUAGE Safe #-} -{-# LANGUAGE PolyKinds #-} -{-# LANGUAGE ConstraintKinds #-} -{-# LANGUAGE DefaultSignatures #-} -{-# LANGUAGE DeriveFunctor #-} -{-# LANGUAGE DeriveGeneric #-} -{-# LANGUAGE FlexibleContexts #-} -{-# LANGUAGE FlexibleInstances #-} ---{-# LANGUAGE RebindableSyntax #-} -{-# LANGUAGE TypeOperators #-} -{-# LANGUAGE TypeFamilies #-} - --- | See the /spatial-math/ package for usage. -module Data.Algebra.Quaternion where - -import safe Data.Algebra -import safe Data.Distributive -import safe Data.Fixed -import safe Data.Functor.Rep -import safe Data.Semifield -import safe Data.Semigroup.Foldable -import safe Data.Semimodule -import safe Data.Semimodule.Vector -import safe Data.Semiring -import safe GHC.Generics hiding (Rep) -import safe Prelude hiding (Num(..), Fractional(..), sum, product) - -{- need tolerances: -λ> prop_conj q12 (q3 :: QuatP) -False -λ> prop_conj q14 (q3 :: QuatP) -False - -prop_conj :: Ring a => (a -> a -> Bool) -> Quaternion a -> Quaternion a -> Bool -prop_conj (~~) p q = sum $ mzipWithRep (~~) (conj (p * q)) (conj q * conj p) - --- conj (p * q) = conj q * conj p --- conj q = (-0.5) * (q <> (i * q * i) <> (j * q * j) <> (k * q * k)) --- 2 * real q '==' q <> conj q --- 2 * imag q '==' q << conj q -conj :: Group a => Quaternion a -> Quaternion a -conj (Quaternion r v) = Quaternion r $ fmap negate v - --- TODO: add to Property module -prop_conj' :: Field a => Rel (Quaternion a) b -> Quaternion a -> b -prop_conj' (~~) q = (conj q) ~~ (conj' q) where - conj' q = ((one / negate two) *) <$> q <> (qi * q * qi) <> (qj * q * qj) <> (qk * q * qk) --} - - - - - -type QuatF = Quaternion Float -type QuatD = Quaternion Double -type QuatR = Quaternion Rational -type QuatM = Quaternion Micro -type QuatN = Quaternion Nano -type QuatP = Quaternion Pico - -data Quaternion a = Quaternion !a {-# UNPACK #-}! (V3 a) deriving (Eq, Ord, Show, Generic, Generic1) - --- | Obtain a 'Quaternion' from 4 base field elements. --- -quat :: a -> a -> a -> a -> Quaternion a -quat r x y z = Quaternion r (V3 x y z) - --- | Real or scalar part of a quaternion. --- -scal :: Quaternion a -> a -scal (Quaternion r _) = r - -vect :: Quaternion a -> V3 a -vect (Quaternion _ v) = v - --- | Use a quaternion to rotate a vector. --- --- >>> rotate qk . rotate qj $ V3 1 1 0 :: V3 Int --- V3 1 (-1) 0 --- -rotate :: Ring a => Quaternion a -> V3 a -> V3 a -rotate q v = v' where Quaternion _ v' = q * Quaternion zero v * conj q - --- | Scale a 'QuatD' to unit length. --- --- >>> normalize $ normalize $ quat 2.0 2.0 2.0 2.0 --- Quaternion 0.5 (V3 0.5 0.5 0.5) --- -normalize :: QuatD -> QuatD -normalize q = 1.0 / (sqrt $ norm q) *. q - -------------------------------------------------------------------------------- --- Standard quaternion basis elements -------------------------------------------------------------------------------- - --- | The real quaternion. --- --- Represents no rotation. --- --- 'qe' = 'unit' --- -qe :: Semiring a => Quaternion a -qe = idx Nothing - --- | The /i/ quaternion. --- --- Represents a \( \pi \) radian rotation about the /x/ axis. --- --- >>> rotate (qi :: QuatM) $ V3 1 0 0 --- V3 1.000000 0.000000 0.000000 --- >>> rotate (qi :: QuatM) $ V3 0 1 0 --- V3 0.000000 -1.000000 0.000000 --- >>> rotate (qi :: QuatM) $ V3 0 0 1 --- V3 0.000000 0.000000 -1.000000 --- --- >>> qi * qj --- Quaternion 0 (V3 0 0 1) --- -qi :: Semiring a => Quaternion a -qi = idx (Just I31) - --- | The /j/ quaternion. --- --- Represents a \( \pi \) radian rotation about the /y/ axis. --- --- >>> rotate (qj :: QuatM) $ V3 1 0 0 --- V3 -1.000000 0.000000 0.000000 --- >>> rotate (qj :: QuatM) $ V3 0 1 0 --- V3 0.000000 1.000000 0.000000 --- >>> rotate (qj :: QuatM) $ V3 0 0 1 --- V3 0.000000 0.000000 -1.000000 --- --- >>> qj * qk --- Quaternion 0 (V3 1 0 0) --- -qj :: Semiring a => Quaternion a -qj = idx (Just I32) - --- | The /k/ quaternion. --- --- Represents a \( \pi \) radian rotation about the /z/ axis. --- --- >>> rotate (qk :: QuatM) $ V3 1 0 0 --- V3 -1.000000 0.000000 0.000000 --- >>> rotate (qk :: QuatM) $ V3 0 1 0 --- V3 0.000000 -1.000000 0.000000 --- >>> rotate (qk :: QuatM) $ V3 0 0 1 --- V3 0.000000 0.000000 1.000000 --- --- >>> qk * qi --- Quaternion 0 (V3 0 1 0) --- >>> qi * qj * qk --- Quaternion (-1) (V3 0 0 0) --- -qk :: Semiring a => Quaternion a -qk = idx (Just I33) - -------------------------------------------------------------------------------- --- Instances -------------------------------------------------------------------------------- - -instance (Additive-Semigroup) a => Semigroup (Quaternion a) where - (<>) = mzipWithRep (+) - -instance (Additive-Monoid) a => Monoid (Quaternion a) where - mempty = pureRep zero - -instance (Additive-Group) a => Magma (Quaternion a) where - (<<) = mzipWithRep (-) - -instance (Additive-Group) a => Quasigroup (Quaternion a) - -instance (Additive-Group) a => Loop (Quaternion a) - -instance (Additive-Group) a => Group (Quaternion a) - -instance (Additive-Group) a => Magma (Additive (Quaternion a)) where - (<<) = mzipWithRep (<<) - -instance (Additive-Group) a => Quasigroup (Additive (Quaternion a)) - -instance (Additive-Group) a => Loop (Additive (Quaternion a)) - -instance (Additive-Group) a => Group (Additive (Quaternion a)) - -instance Semiring a => Semimodule a (Quaternion a) where - (*.) = multl - -instance (Additive-Semigroup) a => Semigroup (Additive (Quaternion a)) where - (<>) = mzipWithRep (<>) - -instance (Additive-Monoid) a => Monoid (Additive (Quaternion a)) where - mempty = pure mempty - -instance Ring a => Semigroup (Multiplicative (Quaternion a)) where - -- >>> qi * qj :: QuatM - -- Quaternion 0.000000 (V3 0.000000 0.000000 1.000000) - -- >>> qk * qi :: QuatM - -- Quaternion 0.000000 (V3 0.000000 1.000000 0.000000) - -- >>> qj * qk :: QuatM - -- Quaternion 0.000000 (V3 1.000000 0.000000 0.000000) - (<>) = mzipWithRep (><) - -instance Ring a => Monoid (Multiplicative (Quaternion a)) where - mempty = pure unit - -instance Ring a => Presemiring (Quaternion a) - -instance Ring a => Semiring (Quaternion a) - -instance Ring a => Ring (Quaternion a) - -instance Functor Quaternion where - fmap f (Quaternion r v) = Quaternion (f r) (fmap f v) - {-# INLINE fmap #-} - - a <$ _ = Quaternion a (V3 a a a) - {-# INLINE (<$) #-} - -instance Foldable Quaternion where - foldMap f (Quaternion e v) = f e <> foldMap f v - {-# INLINE foldMap #-} - foldr f z (Quaternion e v) = f e (foldr f z v) - {-# INLINE foldr #-} - null _ = False - length _ = 4 - -instance Foldable1 Quaternion where - foldMap1 f (Quaternion r v) = f r <> foldMap1 f v - {-# INLINE foldMap1 #-} - -instance Distributive Quaternion where - distribute f = Quaternion (fmap (\(Quaternion x _) -> x) f) $ V3 - (fmap (\(Quaternion _ (V3 y _ _)) -> y) f) - (fmap (\(Quaternion _ (V3 _ z _)) -> z) f) - (fmap (\(Quaternion _ (V3 _ _ w)) -> w) f) - {-# INLINE distribute #-} - -instance Representable Quaternion where - type Rep Quaternion = Maybe I3 - - tabulate f = Quaternion (f Nothing) (V3 (f $ Just I31) (f $ Just I32) (f $ Just I33)) - {-# INLINE tabulate #-} - - index (Quaternion r v) = maybe r (index v) - {-# INLINE index #-} diff --git a/src/Data/Semifield.hs b/src/Data/Semifield.hs --- a/src/Data/Semifield.hs +++ b/src/Data/Semifield.hs @@ -10,82 +10,88 @@ {-# LANGUAGE TypeOperators #-} module Data.Semifield ( - (/) - , (^^) + -- * Semifields + type SemifieldLaw, Semifield + , anan, pinf + , (/), (\\), (^^) , recip - , anan - , pinf + -- * Fields + , type FieldLaw, Field, Real , ninf - , type SemifieldLaw, Semifield - , type FieldLaw, Field ) where -import safe Data.Bool import safe Data.Complex import safe Data.Fixed -import safe Data.Foldable as Foldable (fold, foldl') -import safe Data.Int import safe Data.Semiring -import safe Data.Semigroup.Foldable as Foldable1 -import safe Data.Semigroup.Additive import safe Data.Semigroup.Multiplicative -import safe Data.Tuple -import safe Data.Word -import safe GHC.Real hiding (Fractional(..), (^^), (^), div) +import safe GHC.Real hiding (Real, Fractional(..), (^^), (^), div) import safe Numeric.Natural import safe Foreign.C.Types (CFloat(..),CDouble(..)) -import Prelude ( Eq(..), Ord(..), Show(..), Applicative(..), Functor(..), Monoid(..), Semigroup(..), (.), ($), flip, (<$>), Integer, fromInteger, Float, Double) -import qualified Prelude as P +import Prelude (Monoid(..) , Float, Double) -infixr 8 ^^ +------------------------------------------------------------------------------- +-- Semifields +------------------------------------------------------------------------------- --- @ 'one' '==' a '^^' 0 @ +type SemifieldLaw a = ((Additive-Monoid) a, (Multiplicative-Group) a) + +-- | A semifield, near-field, division ring, or associative division algebra. -- --- >>> 8 ^^ 0 :: Double --- 1.0 --- >>> 8 ^^ 0 :: Pico --- 1.000000000000 +-- Instances needn't have commutative multiplication or additive inverses, +-- however addition and multiplication must be associative as usual. -- -(^^) :: (Multiplicative-Group) a => a -> Integer -> a -a ^^ n = unMultiplicative $ greplicate n (Multiplicative a) - --- | Take the reciprocal of a multiplicative group element. +-- See also the wikipedia definitions of: -- --- >>> recip (3 :+ 4) :: Complex Rational --- 3 % 25 :+ (-4) % 25 --- >>> recip (3 :+ 4) :: Complex Double --- 0.12 :+ (-0.16) --- >>> recip (3 :+ 4) :: Complex Pico --- 0.120000000000 :+ -0.160000000000 +-- * < https://en.wikipedia.org/wiki/Semifield semifield > +-- * < https://en.wikipedia.org/wiki/Near-field_(mathematics) near-field > +-- * < https://en.wikipedia.org/wiki/Division_ring division ring > +-- * < https://en.wikipedia.org/wiki/Division_algebra division algebra >. -- -recip :: (Multiplicative-Group) a => a -> a -recip a = one / a -{-# INLINE recip #-} +class (Semiring a, SemifieldLaw a) => Semifield a +-- | The /NaN/ value of the semifield. +-- +-- @ 'anan' = 'zero' '/' 'zero' @ +-- anan :: Semifield a => a anan = zero / zero {-# INLINE anan #-} +-- | The positive infinity of the semifield. +-- +-- @ 'pinf' = 'one' '/' 'zero' @ +-- pinf :: Semifield a => a pinf = one / zero {-# INLINE pinf #-} -ninf :: Field a => a -ninf = negate one / zero -{-# INLINE ninf #-} +------------------------------------------------------------------------------- +-- Fields +------------------------------------------------------------------------------- --- Sometimes called a division ring -type SemifieldLaw a = ((Additive-Monoid) a, (Multiplicative-Group) a) +type FieldLaw a = ((Additive-Group) a, (Multiplicative-Group) a) --- | A semifield, near-field, division ring, or associative division algebra. +-- | A < https://en.wikipedia.org/wiki/Field_(mathematics) field >. -- --- Instances needn't have commutative multiplication or additive inverses. +class (Ring a, Semifield a, FieldLaw a) => Field a + +-- | A type modeling the real numbers. -- --- See also the wikipedia definitions of < https://en.wikipedia.org/wiki/Semifield semifield >, < https://en.wikipedia.org/wiki/Near-field_(mathematics) near-field >, < https://en.wikipedia.org/wiki/Division_ring division ring >, and < https://en.wikipedia.org/wiki/Division_algebra division algebra >. --- -class (Semiring a, SemifieldLaw a) => Semifield a +class Field a => Real a +-- | The 'one' '/' 'zero' value of the field. +-- +-- @ 'ninf' = 'negate' 'one' '/' 'zero' @ +-- +ninf :: Field a => a +ninf = negate one / zero +{-# INLINE ninf #-} + +------------------------------------------------------------------------------- +-- Instances +------------------------------------------------------------------------------- + instance Semifield () instance Semifield (Ratio Natural) instance Semifield Rational @@ -105,10 +111,7 @@ instance Field a => Semifield (Complex a) -type FieldLaw a = ((Additive-Group) a, (Multiplicative-Group) a) -class (Ring a, Semifield a, FieldLaw a) => Field a - instance Field () instance Field Rational @@ -127,9 +130,8 @@ instance Field a => Field (Complex a) -{- -class (Ord a, Field a) => Real a + instance Real Rational instance Real Uni @@ -144,5 +146,3 @@ instance Real Double instance Real CFloat instance Real CDouble --} - diff --git a/src/Data/Semigroup/Additive.hs b/src/Data/Semigroup/Additive.hs --- a/src/Data/Semigroup/Additive.hs +++ b/src/Data/Semigroup/Additive.hs @@ -20,31 +20,36 @@ import safe Data.Distributive import safe Data.Functor.Rep import safe Data.Fixed -import safe Data.Foldable hiding (sum) import safe Data.Group import safe Data.Int -import safe Data.List import safe Data.List.NonEmpty import safe Data.Ord import safe Data.Semigroup -import safe Data.Semigroup.Foldable import safe Data.Semigroup.Multiplicative -import safe Data.Tuple import safe Data.Word import safe Foreign.C.Types (CFloat(..),CDouble(..)) import safe GHC.Generics (Generic) import safe GHC.Real hiding (Fractional(..), div, (^^), (^), (%)) import safe Numeric.Natural -import safe Prelude ( Eq(..), Ord(..), Show, Ordering(..), Bounded(..), Applicative(..), Functor(..), Monoid(..), Semigroup(..), (.), ($), flip, (<$>), Integer, Float, Double) +import safe Prelude + ( Eq(..), Ord(..), Show, Applicative(..), Functor(..), Monoid(..), Semigroup(..) + , (.), ($), (<$>), Integer, Float, Double) import safe qualified Prelude as P -import qualified Data.Map as Map -import qualified Data.Set as Set -import qualified Data.IntMap as IntMap -import qualified Data.IntSet as IntSet +import safe qualified Data.Map as Map +import safe qualified Data.Set as Set +import safe qualified Data.IntMap as IntMap +import safe qualified Data.IntSet as IntSet +-- | A commutative 'Semigroup' under '+'. +newtype Additive a = Additive { unAdditive :: a } deriving (Eq, Generic, Ord, Show, Functor) + +zero :: (Additive-Monoid) a => a +zero = unAdditive mempty +{-# INLINE zero #-} + infixl 6 + -- >>> Dual [2] + Dual [3] :: Dual [Int] @@ -59,13 +64,24 @@ a - b = unAdditive (Additive a << Additive b) {-# INLINE (-) #-} -zero :: (Additive-Monoid) a => a -zero = unAdditive mempty -{-# INLINE zero #-} +negate :: (Additive-Group) a => a -> a +negate a = zero - a +{-# INLINE negate #-} --- | A commutative 'Semigroup' under '+'. -newtype Additive a = Additive { unAdditive :: a } deriving (Eq, Generic, Ord, Show, Functor) +-- | Absolute value of an element. +-- +-- @ 'abs' r = 'mul' r ('signum' r) @ +-- +-- https://en.wikipedia.org/wiki/Linearly_ordered_group +abs :: (Additive-Group) a => Ord a => a -> a +abs x = bool (negate x) x $ zero <= x +{-# INLINE abs #-} +------------------------------------------------------------------------------- +-- Instances +------------------------------------------------------------------------------- + + instance Applicative Additive where pure = Additive Additive f <*> Additive a = Additive (f a) @@ -84,7 +100,6 @@ - {- newtype Ordered a = Ordered { unOrdered :: a } deriving (Eq, Generic, Ord, Show, Functor) @@ -151,7 +166,6 @@ -} - --------------------------------------------------------------------- -- Num-based --------------------------------------------------------------------- @@ -505,22 +519,31 @@ -- instance (Meet-Monoid) (Down a) => Monoid (Meet (Down a)) where mempty = Down <$> mempty instance ((Additive-Semigroup) a, (Additive-Semigroup) b) => Semigroup (Additive (a, b)) where - Additive (x1, y1) <> Additive (x2, y2) = Additive (x1 + x2, y1 + y2) + (<>) = liftA2 $ \(x1,y1) (x2,y2) -> (x1+x2, y1+y2) +instance ((Additive-Monoid) a, (Additive-Monoid) b) => Monoid (Additive (a, b)) where + mempty = pure (zero, zero) + +instance ((Additive-Semigroup) a, (Additive-Semigroup) b, (Additive-Semigroup) c) => Semigroup (Additive (a, b, c)) where + (<>) = liftA2 $ \(x1,y1,z1) (x2,y2,z2) -> (x1+x2, y1+y2, z1+z2) + +instance ((Additive-Monoid) a, (Additive-Monoid) b, (Additive-Monoid) c) => Monoid (Additive (a, b, c)) where + mempty = pure (zero, zero, zero) + instance (Additive-Semigroup) a => Semigroup (Additive (Maybe a)) where Additive (Just x) <> Additive (Just y) = Additive . Just $ x + y Additive (x@Just{}) <> _ = Additive x Additive Nothing <> y = y +instance (Additive-Semigroup) a => Monoid (Additive (Maybe a)) where + mempty = Additive Nothing + instance ((Additive-Semigroup) a, (Additive-Semigroup) b) => Semigroup (Additive (Either a b)) where Additive (Right x) <> Additive (Right y) = Additive . Right $ x + y Additive(x@Right{}) <> _ = Additive x Additive (Left x) <> Additive (Left y) = Additive . Left $ x + y Additive (Left _) <> y = y - -instance (Additive-Semigroup) a => Monoid (Additive (Maybe a)) where - mempty = Additive Nothing instance Ord a => Semigroup (Additive (Set.Set a)) where (<>) = liftA2 Set.union diff --git a/src/Data/Semigroup/Multiplicative.hs b/src/Data/Semigroup/Multiplicative.hs --- a/src/Data/Semigroup/Multiplicative.hs +++ b/src/Data/Semigroup/Multiplicative.hs @@ -15,69 +15,96 @@ import safe Data.Ord import safe Control.Applicative import safe Data.Bool -import safe Data.Complex +import safe Data.Distributive +import safe Data.Functor.Rep import safe Data.Maybe import safe Data.Either import safe Data.Fixed -import safe Data.Foldable as Foldable (Foldable, foldr', foldl') import safe Data.Group import safe Data.Int -import safe Data.List -import safe Data.List.NonEmpty import safe Data.Semigroup -import safe Data.Semigroup.Foldable as Foldable1 -import safe Data.Tuple import safe Data.Word import safe Foreign.C.Types (CFloat(..),CDouble(..)) import safe GHC.Generics (Generic) import safe GHC.Real hiding (Fractional(..), div, (^^), (^)) import safe Numeric.Natural ---import safe Prelude ( Eq, Ord, Show, Applicative(..), Functor(..), Monoid(..), Semigroup(..), (.), ($), flip, (<$>), Integer, Float, Double) -import safe qualified Prelude as P -import safe Prelude ( Eq(..), Ord, Show, Ordering(..), Bounded(..), Applicative(..), Functor(..), Monoid(..), Semigroup(..), (.), ($), flip, (<$>), Integer, Float, Double) -import safe qualified Prelude as P +import safe Prelude + ( Eq(..), Ord, Show, Applicative(..), Functor(..), Monoid(..) + , Semigroup(..), (.), ($), flip, (<$>), Integer, Float, Double) -import qualified Data.Map as Map -import qualified Data.Set as Set -import qualified Data.IntMap as IntMap -import qualified Data.IntSet as IntSet -import qualified Data.Sequence as Seq +import safe qualified Prelude as P +import safe qualified Data.Map as Map +import safe qualified Data.Set as Set +import safe qualified Data.IntMap as IntMap +import safe qualified Data.IntSet as IntSet -import safe Data.Distributive -import safe Data.Functor.Rep - infixr 1 - -- | Hyphenation operator. type (g - f) a = f (g a) -infixl 7 * +-- | A (potentially non-commutative) 'Semigroup' under '+'. +newtype Multiplicative a = Multiplicative { unMultiplicative :: a } deriving (Eq, Generic, Ord, Show, Functor) +one :: (Multiplicative-Monoid) a => a +one = unMultiplicative mempty +{-# INLINE one #-} + +infixl 7 *, \\, / + -- >>> Dual [2] * Dual [3] :: Dual [Int] -- Dual {getDual = [5]} (*) :: (Multiplicative-Semigroup) a => a -> a -> a a * b = unMultiplicative (Multiplicative a <> Multiplicative b) {-# INLINE (*) #-} -infixl 7 / - (/) :: (Multiplicative-Group) a => a -> a -> a a / b = unMultiplicative (Multiplicative a << Multiplicative b) {-# INLINE (/) #-} +-- | Left division by a multiplicative group element. +-- +-- When '*' is commutative we must have: +-- +-- @ x '\\' y = y '/' x @ +-- +(\\) :: (Multiplicative-Group) a => a -> a -> a +(\\) x y = recip x * y +infixr 8 ^^ -one :: (Multiplicative-Monoid) a => a -one = unMultiplicative mempty -{-# INLINE one #-} +-- | Integral power of a multiplicative group element. +-- +-- @ 'one' '==' a '^^' 0 @ +-- +-- >>> 8 ^^ 0 :: Double +-- 1.0 +-- >>> 8 ^^ 0 :: Pico +-- 1.000000000000 +-- +(^^) :: (Multiplicative-Group) a => a -> Integer -> a +a ^^ n = unMultiplicative $ greplicate n (Multiplicative a) -div :: (Multiplicative-Group) a => a -> a -> a -a `div` b = unMultiplicative (Multiplicative a << Multiplicative b) -{-# INLINE div #-} +-- | Reciprocal of a multiplicative group element. +-- +-- @ +-- x '/' y = x '*' 'recip' y +-- x '\\' y = 'recip' x '*' y +-- @ +-- +-- >>> recip (3 :+ 4) :: Complex Rational +-- 3 % 25 :+ (-4) % 25 +-- >>> recip (3 :+ 4) :: Complex Double +-- 0.12 :+ (-0.16) +-- >>> recip (3 :+ 4) :: Complex Pico +-- 0.120000000000 :+ -0.160000000000 +-- +recip :: (Multiplicative-Group) a => a -> a +recip a = one / a +{-# INLINE recip #-} -newtype Multiplicative a = Multiplicative { unMultiplicative :: a } deriving (Eq, Generic, Ord, Show, Functor) instance Applicative Multiplicative where pure = Multiplicative @@ -96,7 +123,7 @@ {-# INLINE index #-} --------------------------------------------------------------------- --- Num-based instances +-- Instances --------------------------------------------------------------------- #define deriveMultiplicativeSemigroup(ty) \ @@ -327,7 +354,7 @@ instance (Multiplicative-Semigroup) a => Semigroup (Multiplicative (Maybe a)) where Multiplicative Nothing <> _ = Multiplicative Nothing - Multiplicative (x@Just{}) <> Multiplicative Nothing = Multiplicative Nothing + Multiplicative (Just{}) <> Multiplicative Nothing = Multiplicative Nothing Multiplicative (Just x) <> Multiplicative (Just y) = Multiplicative . Just $ x * y -- Mul a <> Mul b = Mul $ liftA2 (*) a b @@ -336,7 +363,7 @@ instance ((Multiplicative-Semigroup) a, (Multiplicative-Semigroup) b) => Semigroup (Multiplicative (Either a b)) where Multiplicative (Right x) <> Multiplicative (Right y) = Multiplicative . Right $ x * y - Multiplicative(x@Right{}) <> y = y + Multiplicative (Right{}) <> y = y Multiplicative (Left x) <> Multiplicative (Left y) = Multiplicative . Left $ x * y Multiplicative (x@Left{}) <> _ = Multiplicative x diff --git a/src/Data/Semigroup/Property.hs b/src/Data/Semigroup/Property.hs --- a/src/Data/Semigroup/Property.hs +++ b/src/Data/Semigroup/Property.hs @@ -1,31 +1,30 @@ {-# Language AllowAmbiguousTypes #-} {-# LANGUAGE Safe #-} -module Data.Semigroup.Property where -{- ( +module Data.Semigroup.Property ( -- * Required properties of semigroups associative_addition_on , associative_multiplication_on -- * Required properties of monoids , neutral_addition_on , neutral_multiplication_on - -- * Required properties of semigroup & monoid morphisms - , morphism_additive_on - , morphism_multiplicative_on - , morphism_additive_on' - , morphism_multiplicative_on' -- * Properties of commuative semigroups , commutative_addition_on , commutative_multiplication_on - -- * Properties of idempotent semigroups - , idempotent_addition_on - , idempotent_multiplication_on -- * Properties of cancellative semigroups , cancellative_addition_on , cancellative_multiplication_on + -- * Properties of idempotent semigroups + , idempotent_addition_on + , idempotent_multiplication_on + -- * Required properties of semigroup & monoid morphisms + , morphism_additive_on + , morphism_multiplicative_on + , morphism_additive_on' + , morphism_multiplicative_on' ) where --} + import safe Test.Logic (Rel) import safe Data.Semigroup.Additive import safe Data.Semigroup.Multiplicative @@ -34,7 +33,7 @@ import safe Prelude hiding (Num(..), sum) -{- + ------------------------------------------------------------------------------------ -- Required properties of semigroups @@ -45,7 +44,7 @@ -- This is a required property. -- associative_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> r -> r -> b -associative_addition_on (~~) = Prop.associative_on (~~) add +associative_addition_on (~~) = Prop.associative_on (~~) (+) -- | \( \forall a, b, c \in R: (a * b) * c \sim a * (b * c) \) -- @@ -54,7 +53,7 @@ -- This is a required property. -- associative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> r -> r -> b -associative_multiplication_on (~~) = Prop.associative_on (~~) mul +associative_multiplication_on (~~) = Prop.associative_on (~~) (*) ------------------------------------------------------------------------------------ -- Required properties of monoids @@ -64,44 +63,120 @@ -- A semigroup with a right-neutral additive identity must satisfy: -- -- @ --- 'neutral_addition' 'zero' ~~ const True +-- 'neutral_addition' 'zero' = const True -- @ -- -- Or, equivalently: -- -- @ --- 'zero' '+' r ~~ r +-- 'zero' '+' r = r -- @ -- -- This is a required property for additive monoids. -- neutral_addition_on :: (Additive-Monoid) r => Rel r b -> r -> b -neutral_addition_on (~~) = Prop.neutral_on (~~) add zero +neutral_addition_on (~~) = Prop.neutral_on (~~) (+) zero -- | \( \forall a \in R: (o * a) \sim a \) -- -- A semigroup with a right-neutral multiplicative identity must satisfy: -- -- @ --- 'neutral_multiplication' 'one' ~~ const True +-- 'neutral_multiplication' 'one' = const True -- @ -- -- Or, equivalently: -- -- @ --- 'one' '*' r ~~ r +-- 'one' '*' r = r -- @ -- -- This is a required propert for multiplicative monoids. -- neutral_multiplication_on :: (Multiplicative-Monoid) r => Rel r b -> r -> b -neutral_multiplication_on (~~) = Prop.neutral_on (~~) mul one +neutral_multiplication_on (~~) = Prop.neutral_on (~~) (*) one --} +------------------------------------------------------------------------------------ +-- Properties of commutative semigroups +-- | \( \forall a, b \in R: a + b \sim b + a \) +-- +-- This is a an /optional/ property for semigroups, and a /required/ property for semirings. +-- +commutative_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> r -> b +commutative_addition_on (~~) = Prop.commutative_on (~~) (+) + +-- | \( \forall a, b \in R: a * b \sim b * a \) +-- +-- This is a an /optional/ property for semigroups, and a /optional/ property for semirings. +-- It is a /required/ property for rings. +-- +commutative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> r -> b +commutative_multiplication_on (~~) = Prop.commutative_on (~~) (*) + ------------------------------------------------------------------------------------ +-- Properties of cancellative semigroups + +-- | \( \forall a, b, c \in R: b + a \sim c + a \Rightarrow b = c \) +-- +-- If /R/ is right-cancellative wrt addition then for all /a/ +-- the section /(a +)/ is injective. +-- +-- See < https://en.wikipedia.org/wiki/Cancellation_property > +-- +cancellative_addition_on :: (Additive-Semigroup) r => Rel r Bool -> r -> r -> r -> Bool +cancellative_addition_on (~~) a = Prop.injective_on (~~) (+ a) + +-- | \( \forall a, b, c \in R: b * a \sim c * a \Rightarrow b = c \) +-- +-- If /R/ is right-cancellative wrt multiplication then for all /a/ +-- the section /(a *)/ is injective. +-- +cancellative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r Bool -> r -> r -> r -> Bool +cancellative_multiplication_on (~~) a = Prop.injective_on (~~) (* a) + +------------------------------------------------------------------------------------ +-- Properties of idempotent semigroups + +-- | Idempotency property for additive semigroups. +-- +-- @ 'idempotent_addition' = 'absorbative_addition' 'one' @ +-- +-- See < https://en.wikipedia.org/wiki/Band_(mathematics) >. +-- +-- This is a required property for lattices. +-- +idempotent_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> b +idempotent_addition_on (~~) r = (r + r) ~~ r + +-- | Idempotency property for multplicative semigroups. +-- +-- @ 'idempotent_multiplication' = 'absorbative_multiplication' 'zero' @ +-- +-- See < https://en.wikipedia.org/wiki/Band_(mathematics) >. +-- +-- This is a an /optional/ property for semigroups, and a /optional/ property for semirings. +-- +-- This is a /required/ property for lattices. +-- +idempotent_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> b +idempotent_multiplication_on (~~) r = (r * r) ~~ r + +------------------------------------------------------------------------------------ -- Properties of semigroup morphisms +morphism_additive_on :: (Additive-Semigroup) r => (Additive-Semigroup) s => Rel s b -> (r -> s) -> r -> r -> b +morphism_additive_on (~~) f x y = (f $ x + y) ~~ (f x + f y) + +morphism_multiplicative_on :: (Multiplicative-Semigroup) r => (Multiplicative-Semigroup) s => Rel s b -> (r -> s) -> r -> r -> b +morphism_multiplicative_on (~~) f x y = (f $ x * y) ~~ (f x * f y) + +morphism_additive_on' :: (Additive-Monoid) r => (Additive-Monoid) s => Rel s b -> (r -> s) -> b +morphism_additive_on' (~~) f = (f zero) ~~ zero + +morphism_multiplicative_on' :: (Multiplicative-Monoid) r => (Multiplicative-Monoid) s => Rel s b -> (r -> s) -> b +morphism_multiplicative_on' (~~) f = (f one) ~~ one + {- morphism_additive_on :: (Additive-Semigroup) r => (Additive-Semigroup) s => Rel s b -> (r -> s) -> r -> r -> b morphism_additive_on (~~) f x y = (f $ x `add` y) ~~ (f x `add` f y) @@ -115,16 +190,7 @@ morphism_multiplicative_on' :: (Multiplicative-Monoid) r => (Multiplicative-Monoid) s => Rel s b -> (r -> s) -> b morphism_multiplicative_on' (~~) f = (f one) ~~ one ------------------------------------------------------------------------------------- --- Properties of commutative semigroups --- | \( \forall a, b \in R: a + b \sim b + a \) --- --- This is a an /optional/ property for semigroups, and a /required/ property for semirings. --- -commutative_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> r -> b -commutative_addition_on (~~) = Prop.commutative_on (~~) add - -- | \( \forall a, b \in R: a * b \sim b * a \) -- -- This is a an /optional/ property for semigroups, and a /optional/ property for semirings. @@ -134,7 +200,6 @@ commutative_multiplication_on (~~) = Prop.commutative_on (~~) mul -} ------------------------------------------------------------------------------------- --- Properties of idempotent dioids and predioids + diff --git a/src/Data/Semimodule.hs b/src/Data/Semimodule.hs --- a/src/Data/Semimodule.hs +++ b/src/Data/Semimodule.hs @@ -11,28 +11,43 @@ {-# LANGUAGE TypeOperators #-} {-# LANGUAGE TypeFamilies #-} -module Data.Semimodule where +module Data.Semimodule ( + -- * Types + type Free + , type Basis + , type Basis2 + , type Basis3 + , type FreeModule + , type FreeSemimodule + -- * Left modules + , type LeftModule + , LeftSemimodule(..) + , lscaleDef + , negateDef + , lerp + , (*.) + , (/.) + , (\.) + -- * Right modules + , type RightModule + , RightSemimodule(..) + , rscaleDef + , (.*) + , (./) + , (.\) + -- * Bimodules + , type Bimodule + , Bisemimodule(..) +) where -import safe Data.Bool import safe Data.Complex import safe Data.Semifield -import safe Data.Fixed -import safe Data.Functor.Compose import safe Data.Functor.Rep -import safe Data.Int import safe Data.Semiring -import safe Data.Semigroup.Foldable as Foldable1 -import safe Data.Tuple -import safe Data.Word import safe GHC.Real hiding (Fractional(..)) import safe Numeric.Natural -import safe Foreign.C.Types (CFloat(..),CDouble(..)) import safe Prelude hiding (Num(..), Fractional(..), sum, product) -import safe qualified Prelude as N -import safe Data.Semigroup.Additive as A -import safe Data.Semigroup.Multiplicative as M - import safe Prelude (fromInteger) @@ -40,74 +55,51 @@ type Basis b f = (Free f, Rep f ~ b) -{- - --- Semimodule over a semifield --- dioids -type DSpace r a = (Semifield r, Semimodule r a) - - --- | Free semimodule over a generating set. --- -type FreeSemimodule a f = (Free f, Semimodule a (f a)) - -type FreeModule a f = (Free f, Module a (f a)) - -type CommutativeGroup a = Module Integer a +type Basis2 b c f g = (Basis b f, Basis c g) --} +type Basis3 b c d f g h = (Basis b f, Basis c g, Basis d h) +type FreeModule a f = (Free f, Bimodule a a (f a)) ---instance (Unital (f a), Algebra (f a), Functor f) => Semifield (f a) where - --recip q = conj' q // norm' q --- recip q = ((recip . norm' $ q) ><) <$> conj' q +type FreeSemimodule a f = (Free f, Bisemimodule a a (f a)) -type Module r a = (Ring r, Group a, Semimodule r a) +------------------------------------------------------------------------------- +-- Left modules +------------------------------------------------------------------------------- -infixl 7 .*, *. +type LeftModule l a = (Ring l, (Additive-Group) a, LeftSemimodule l a) --- | < https://en.wikipedia.org/wiki/Semimodule Semimodule > over a commutative semiring. +-- | < https://en.wikipedia.org/wiki/Semimodule Left semimodule > over a commutative semiring. -- -- All instances must satisfy the following identities: -- --- @ r '*.' (x '<>' y) '==' r '*.' x '<>' r '*.' y @ --- --- @ (r '+' s) '*.' x '==' r '*.' x '<>' s '*.' x @ --- --- @ (r '*' s) '*.' x '==' r '*.' (s '*.' x) @ --- --- When the ring of coefficients /r/ is unital we must additionally have: --- --- @ 'one' '*.' x '==' x @ +-- @ +-- 'lscale' s (x '+' y) = 'lscale' s x '+' 'lscale' s y +-- 'lscale' (s1 '+' s2) x = 'lscale' s1 x '+' 'lscale' s2 x +-- 'lscale' (s1 '*' s2) = 'lscale' s1 . 'lscale' s2 +-- 'lscale' 'zero' = 'zero' +-- @ -- +-- When the ring of coefficients /s/ is unital we must additionally have: +-- @ +-- 'lscale' 'one' = 'id' +-- @ +-- -- See the properties module for a detailed specification of the laws. -- -class (Semiring r, Semigroup a) => Semimodule r a where +class (Semiring l, (Additive-Monoid) a) => LeftSemimodule l a where -- | Left-multiply by a scalar. -- - (*.) :: r -> a -> a - (*.) = flip (.*) - - -- | Right-multiply by a scalar. - -- - (.*) :: a -> r -> a - (.*) = flip (*.) - - - --- | Default definition of '(*.)' for a free module. --- -multl :: Semiring a => Functor f => a -> f a -> f a -multl a f = (a *) <$> f + lscale :: l -> a -> a --- | Default definition of '(.*)' for a free module. +-- | Default definition of 'lscale' for a free module. -- -multr :: Semiring a => Functor f => f a -> a -> f a -multr f a = (* a) <$> f +lscaleDef :: Semiring a => Functor f => a -> f a -> f a +lscaleDef a f = (a *) <$> f -- | Default definition of '<<' for a commutative group. -- -negateDef :: Semimodule Integer a => a -> a +negateDef :: LeftModule Integer a => a -> a negateDef a = (-1 :: Integer) *. a -- | Linearly interpolate between two vectors. @@ -118,107 +110,178 @@ -- >>> lerp r u v -- V3 (6 % 4) (12 % 4) (18 % 4) -- -lerp :: Module r a => r -> a -> a -> a -lerp r f g = r *. f <> (one - r) *. g +lerp :: LeftModule r a => r -> a -> a -> a +lerp r f g = r *. f + (one - r) *. g {-# INLINE lerp #-} -infix 6 .*. +infixr 7 *., \., /. --- | Dot product. --- --- >>> V3 1 2 3 .*. V3 1 2 3 --- 14 --- -(.*.) :: Free f => Foldable f => Semiring a => f a -> f a -> a -(.*.) x y = sum $ liftR2 (*) x y -{-# INLINE (.*.) #-} +(*.) :: LeftSemimodule l a => l -> a -> a +(*.) = lscale --- | Squared /l2/ norm of a vector. --- -quadrance :: Free f => Foldable f => Semiring a => f a -> a -quadrance f = f .*. f -{-# INLINE quadrance #-} +(/.) :: Semifield a => Functor f => a -> f a -> f a +a /. f = (a /) <$> f --- | Squared /l2/ norm of the difference between two vectors. --- -qd :: Free f => Foldable f => Module a (f a) => f a -> f a -> a -qd f g = quadrance $ f << g -{-# INLINE qd #-} +(\.) :: Semifield a => Functor f => a -> f a -> f a +a \. f = (a \\) <$> f --- | Dirac delta function. + +------------------------------------------------------------------------------- +-- Right modules +------------------------------------------------------------------------------- + +type RightModule r a = (Ring r, (Additive-Group) a, RightSemimodule r a) + +-- | < https://en.wikipedia.org/wiki/Semimodule Right semimodule > over a commutative semiring. -- -dirac :: Eq i => Semiring a => i -> i -> a -dirac i j = bool zero one (i == j) -{-# INLINE dirac #-} +-- The laws for right semimodules are analagous to those of left semimodules. +-- +-- See the properties module for a detailed specification. +-- +class (Semiring r, (Additive-Monoid) a) => RightSemimodule r a where --- | Create a unit vector at an index. + -- | Right-multiply by a scalar. + -- + rscale :: r -> a -> a + +-- | Default definition of 'rscale' for a free module. -- --- >>> idx I21 :: V2 Int --- V2 1 0 +rscaleDef :: Semiring a => Functor f => a -> f a -> f a +rscaleDef a f = (* a) <$> f + +infixl 7 .*, .\, ./ +(.*) :: RightSemimodule r a => a -> r -> a +(.*) = flip rscale + +(./) :: Semifield a => Functor f => f a -> a -> f a +(./) = flip (/.) + +(.\) :: Semifield a => Functor f => f a -> a -> f a +(.\) = flip (\.) + +------------------------------------------------------------------------------- +-- Bimodules +------------------------------------------------------------------------------- + +type Bimodule l r a = (LeftModule l a, RightModule r a, Bisemimodule l r a) + +-- | < https://en.wikipedia.org/wiki/Bimodule Bisemimodule > over a commutative semiring. -- --- >>> idx I42 :: V4 Int --- V4 0 1 0 0 +-- @ +-- 'lscale' l . 'rscale' r = 'rscale' r . 'lscale' l +-- @ -- -idx :: Free f => Semiring a => Rep f -> f a -idx i = tabulate $ dirac i -{-# INLINE idx #-} +class (LeftSemimodule l a, RightSemimodule r a) => Bisemimodule l r a where + discale :: l -> r -> a -> a + discale l r = lscale l . rscale r + ------------------------------------------------------------------------------- -- Instances ------------------------------------------------------------------------------- -instance Semiring r => Semimodule r () where - _ *. _ = () +instance Semiring l => LeftSemimodule l () where + lscale _ = const () -instance Semigroup a => Semimodule () a where - _ *. a = a +instance (Additive-Monoid) a => LeftSemimodule () a where + lscale _ = id -instance Monoid a => Semimodule Natural a where - (*.) = mreplicate +instance (Additive-Monoid) a => LeftSemimodule Natural a where + lscale l a = unAdditive $ mreplicate l (Additive a) -instance Group a => Semimodule Integer a where - (*.) = greplicate +instance ((Additive-Monoid) a, (Additive-Group) a) => LeftSemimodule Integer a where + lscale l a = unAdditive $ greplicate l (Additive a) -instance Semimodule r a => Semimodule r (e -> a) where - a *. f = (a *.) <$> f +instance LeftSemimodule l a => LeftSemimodule l (e -> a) where + lscale l = fmap (l *.) -instance (Semimodule r a, Semimodule r b) => Semimodule r (a, b) where - n *. (a, b) = (n *. a, n *. b) +instance (LeftSemimodule l a, LeftSemimodule l b) => LeftSemimodule l (a, b) where + lscale n (a, b) = (n *. a, n *. b) -instance (Semimodule r a, Semimodule r b, Semimodule r c) => Semimodule r (a, b, c) where - n *. (a, b, c) = (n *. a, n *. b, n *. c) +instance (LeftSemimodule l a, LeftSemimodule l b, LeftSemimodule l c) => LeftSemimodule l (a, b, c) where + lscale n (a, b, c) = (n *. a, n *. b, n *. c) -instance (Semiring a, Semimodule r a) => Semimodule r (Additive (Ratio a)) where - a *. (Additive (x :% y)) = Additive $ (a *. x) :% y +instance Semiring a => LeftSemimodule a (Ratio a) where + lscale l (x :% y) = (l * x) :% y -instance (Ring a, Semimodule r a) => Semimodule r (Additive (Complex a)) where - a *. (Additive (x :+ y)) = Additive $ (a *. x) :+ (a *. y) +instance Ring a => LeftSemimodule a (Complex a) where + lscale l (x :+ y) = (l * x) :+ (l * y) -#define deriveSemimodule(ty) \ -instance Semiring ty => Semimodule ty (Additive ty) where { \ - r *. (Additive a) = Additive $ r * a \ -; {-# INLINE (*.) #-} \ +instance Ring a => LeftSemimodule (Complex a) (Complex a) where + lscale = (*) + +{- +#define deriveLeftSemimodule(ty) \ +instance LeftSemimodule ty ty where { \ + lscale = (*) \ +; {-# INLINE lscale #-} \ } -deriveSemimodule(Bool) -deriveSemimodule(Int) -deriveSemimodule(Int8) -deriveSemimodule(Int16) -deriveSemimodule(Int32) -deriveSemimodule(Int64) -deriveSemimodule(Word) -deriveSemimodule(Word8) -deriveSemimodule(Word16) -deriveSemimodule(Word32) -deriveSemimodule(Word64) -deriveSemimodule(Uni) -deriveSemimodule(Deci) -deriveSemimodule(Centi) -deriveSemimodule(Milli) -deriveSemimodule(Micro) -deriveSemimodule(Nano) -deriveSemimodule(Pico) -deriveSemimodule(Float) -deriveSemimodule(Double) -deriveSemimodule(CFloat) -deriveSemimodule(CDouble) +deriveLeftSemimodule(Bool) +deriveLeftSemimodule(Int) +deriveLeftSemimodule(Int8) +deriveLeftSemimodule(Int16) +deriveLeftSemimodule(Int32) +deriveLeftSemimodule(Int64) +deriveLeftSemimodule(Word) +deriveLeftSemimodule(Word8) +deriveLeftSemimodule(Word16) +deriveLeftSemimodule(Word32) +deriveLeftSemimodule(Word64) +deriveLeftSemimodule(Uni) +deriveLeftSemimodule(Deci) +deriveLeftSemimodule(Centi) +deriveLeftSemimodule(Milli) +deriveLeftSemimodule(Micro) +deriveLeftSemimodule(Nano) +deriveLeftSemimodule(Pico) +deriveLeftSemimodule(Float) +deriveLeftSemimodule(Double) +deriveLeftSemimodule(CFloat) +deriveLeftSemimodule(CDouble) +-} + +instance Semiring r => RightSemimodule r () where + rscale _ = const () + +instance (Additive-Monoid) a => RightSemimodule () a where + rscale _ = id + +instance (Additive-Monoid) a => RightSemimodule Natural a where + rscale r a = unAdditive $ mreplicate r (Additive a) + +instance ((Additive-Monoid) a, (Additive-Group) a) => RightSemimodule Integer a where + rscale r a = unAdditive $ greplicate r (Additive a) + +instance RightSemimodule r a => RightSemimodule r (e -> a) where + rscale r = fmap (.* r) + +instance (RightSemimodule r a, RightSemimodule r b) => RightSemimodule r (a, b) where + rscale n (a, b) = (a .* n, b .* n) + +instance (RightSemimodule r a, RightSemimodule r b, RightSemimodule r c) => RightSemimodule r (a, b, c) where + rscale n (a, b, c) = (a .* n, b .* n, c .* n) + +instance Semiring a => RightSemimodule a (Ratio a) where + rscale r (x :% y) = (r * x) :% y + +instance Ring a => RightSemimodule a (Complex a) where + rscale r (x :+ y) = (r * x) :+ (r * y) + +instance Ring a => RightSemimodule (Complex a) (Complex a) where + rscale = (*) + +instance Semiring r => Bisemimodule r r () + +instance Bisemimodule r r a => Bisemimodule r r (e -> a) + +instance (Bisemimodule r r a, Bisemimodule r r b) => Bisemimodule r r (a, b) + +instance (Bisemimodule r r a, Bisemimodule r r b, Bisemimodule r r c) => Bisemimodule r r (a, b, c) + +instance Semiring a => Bisemimodule a a (Ratio a) + +instance Ring a => Bisemimodule a a (Complex a) + +instance Ring a => Bisemimodule (Complex a) (Complex a) (Complex a) + diff --git a/src/Data/Semimodule/Basis.hs b/src/Data/Semimodule/Basis.hs new file mode 100644 --- /dev/null +++ b/src/Data/Semimodule/Basis.hs @@ -0,0 +1,83 @@ +{-# LANGUAGE Safe #-} +{-# LANGUAGE RankNTypes #-} +{-# LANGUAGE TypeFamilies #-} +module Data.Semimodule.Basis ( + -- * Euclidean bases + E1(..), e1, fillE1 + , E2(..), e2, fillE2 + , E3(..), e3, fillE3 + , E4(..), e4, fillE4 + , E6(..) +) where + +import safe Data.Functor.Rep +import safe Data.Semimodule +import safe Prelude hiding (Num(..), Fractional(..), negate, sum, product) + + + +------------------------------------------------------------------------------- +-- Standard basis on one real dimension +------------------------------------------------------------------------------- + +data E1 = E1 deriving (Eq, Ord, Show) + +e1 :: a -> E1 -> a +e1 = const + +fillE1 :: Basis E1 f => a -> f a +fillE1 x = tabulate $ e1 x + +------------------------------------------------------------------------------- +-- Standard basis on two real dimensions +------------------------------------------------------------------------------- + +data E2 = E21 | E22 deriving (Eq, Ord, Show) + +e2 :: a -> a -> E2 -> a +e2 x _ E21 = x +e2 _ y E22 = y + +fillE2 :: Basis E2 f => a -> a -> f a +fillE2 x y = tabulate $ e2 x y + +------------------------------------------------------------------------------- +-- Standard basis on three real dimensions +------------------------------------------------------------------------------- + +data E3 = E31 | E32 | E33 deriving (Eq, Ord, Show) + +e3 :: a -> a -> a -> E3 -> a +e3 x _ _ E31 = x +e3 _ y _ E32 = y +e3 _ _ z E33 = z + +fillE3 :: Basis E3 f => a -> a -> a -> f a +fillE3 x y z = tabulate $ e3 x y z + +------------------------------------------------------------------------------- +-- Standard basis on four real dimensions +------------------------------------------------------------------------------- + +data E4 = E41 | E42 | E43 | E44 deriving (Eq, Ord, Show) + +e4 :: a -> a -> a -> a -> E4 -> a +e4 x _ _ _ E41 = x +e4 _ y _ _ E42 = y +e4 _ _ z _ E43 = z +e4 _ _ _ w E44 = w + +fillE4 :: Basis E4 f => a -> a -> a -> a -> f a +fillE4 x y z w = tabulate $ e4 x y z w + +------------------------------------------------------------------------------- +-- Standard basis on five real dimensions +------------------------------------------------------------------------------- + +data E5 = E51 | E52 | E53 | E54 | E55 deriving (Eq, Ord, Show) + +------------------------------------------------------------------------------- +-- Standard basis on six real dimensions +------------------------------------------------------------------------------- + +data E6 = E61 | E62 | E63 | E64 | E65 | E66 deriving (Eq, Ord, Show) diff --git a/src/Data/Semimodule/Free.hs b/src/Data/Semimodule/Free.hs new file mode 100644 --- /dev/null +++ b/src/Data/Semimodule/Free.hs @@ -0,0 +1,1389 @@ +{-# LANGUAGE CPP #-} +{-# LANGUAGE Safe #-} +{-# LANGUAGE PolyKinds #-} +{-# LANGUAGE ConstraintKinds #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE DeriveFunctor #-} +{-# LANGUAGE DeriveGeneric #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE NoImplicitPrelude #-} +{-# LANGUAGE RebindableSyntax #-} +{-# LANGUAGE TypeOperators #-} +{-# LANGUAGE TypeFamilies #-} +{-# LANGUAGE RankNTypes #-} + +module Data.Semimodule.Free ( + -- * Types + type Free + , type Basis + , type Basis2 + , type Basis3 + -- * Vector arithmetic + , (.*) + , (.#) + , (*.) + , (#.) + , dot + , lerp + , quadrance + , qd + , cross + , triple + -- * Vector accessors and constructors + , dirac + , idx + , elt + , lensRep + , grateRep + -- * Matrix arithmetic + , (.#.) + , trace + , transpose + , inv1 + , inv2 + , bdet2 + , det2 + , bdet3 + , det3 + , inv3 + , bdet4 + , det4 + , inv4 + -- * Matrix accessors and constructors + , tran + , elt2 + , row + , rows + , col + , cols + , diag + , outer + , identity + , diagonal + -- * Vector types + , V1(..) + , unV1 + , V2(..) + , V3(..) + , V4(..) + -- * Matrix types + , type M11 + , type M12 + , type M13 + , type M14 + , type M21 + , type M31 + , type M41 + , type M22 + , type M23 + , type M24 + , type M32 + , type M33 + , type M34 + , type M42 + , type M43 + , type M44 + , m11 + , m12 + , m13 + , m14 + , m21 + , m31 + , m41 + , m22 + , m23 + , m24 + , m32 + , m33 + , m34 + , m42 + , m43 + , m44 +) where + +import safe Control.Applicative +import safe Data.Bool +import safe Data.Distributive +import safe Data.Functor.Classes +import safe Data.Functor.Compose +import safe Data.Functor.Rep +import safe Data.Semifield +import safe Data.Semigroup.Foldable as Foldable1 +import safe Data.Semimodule +import safe Data.Semimodule.Basis +import safe Data.Semimodule.Transform +import safe Data.Semiring +import safe Prelude hiding (Num(..), Fractional(..), negate, sum, product) +import safe Prelude (fromInteger) + +------------------------------------------------------------------------------- +-- Vector Arithmetic +------------------------------------------------------------------------------- + + +infix 7 #., .# + +-- | Multiply a matrix on the left by a row vector. +-- +-- >>> V2 1 2 #. m23 3 4 5 6 7 8 +-- V3 15 18 21 +-- +-- >>> (V2 1 2 #. m23 3 4 5 6 7 8) #. m32 1 0 0 0 0 0 +-- V2 15 0 +-- +(#.) :: Semiring a => Foldable f => Basis2 b c f g => f a -> (f**g) a -> g a +x #. y = tabulate (\j -> x `dot` col j y) +{-# INLINE (#.) #-} + +-- | Multiply a matrix on the right by a column vector. +-- +-- @ ('.#') = 'app' . 'tran' @ +-- +-- >>> app (tran $ m23 1 2 3 4 5 6) (V3 7 8 9) :: V2 Int +-- V2 50 122 +-- >>> m23 1 2 3 4 5 6 .# V3 7 8 9 :: V2 Int +-- V2 50 122 +-- >>> m22 1 0 0 0 .# (m23 1 2 3 4 5 6 .# V3 7 8 9) +-- V2 50 0 +-- +(.#) :: Semiring a => Foldable g => Basis2 b c f g => (f**g) a -> g a -> f a +x .# y = tabulate (\i -> row i x `dot` y) +{-# INLINE (.#) #-} + + +infix 6 `dot` + +-- | Dot product. +-- +-- This is a variant of 'Data.Semiring.xmult' restricted to free functors. +-- +-- >>> V3 1 2 3 `dot` V3 1 2 3 +-- 14 +-- +dot :: Semiring a => Foldable f => Basis b f => f a -> f a -> a +dot x y = sum $ liftR2 (*) x y +{-# INLINE dot #-} + +-- | Squared /l2/ norm of a vector. +-- +quadrance :: Semiring a => Foldable f => Basis b f => f a -> a +quadrance x = x `dot` x +{-# INLINE quadrance #-} + +-- | Squared /l2/ norm of the difference between two vectors. +-- +qd :: FreeModule a f => Foldable f => f a -> f a -> a +qd x y = quadrance $ x - y +{-# INLINE qd #-} + +-- | Cross product +-- +-- @ +-- a `'cross'` a = 'zero' +-- a `'cross'` b = 'negate' ( b `'cross'` a ) , +-- a `'cross'` ( b '+' c ) = ( a `'cross'` b ) '+' ( a `'cross'` c ) , +-- ( r a ) `'cross'` b = a `'cross'` ( r b ) = r ( a `'cross'` b ) . +-- a `'cross'` ( b `'cross'` c ) '+' b `'cross'` ( c `'cross'` a ) '+' c `'cross'` ( a `'cross'` b ) = 'zero' . +-- @ +-- +-- See < https://en.wikipedia.org/wiki/Jacobi_identity Jacobi identity >. +-- +cross :: Ring a => V3 a -> V3 a -> V3 a +cross (V3 a b c) (V3 d e f) = V3 (b*f-c*e) (c*d-a*f) (a*e-b*d) +{-# INLINABLE cross #-} + +-- | Scalar triple product. +-- +-- @ +-- 'triple' x y z = 'triple' z x y = 'triple' y z x +-- 'triple' x y z = 'negate' '$' 'triple' x z y = 'negate' '$' 'triple' y x z +-- 'triple' x x y = 'triple' x y y = 'triple' x y x = 'zero' +-- ('triple' x y z) '*.' x = (x `'cross'` y) `'cross'` (x `'cross'` z) +-- @ +-- +-- >>> triple (V3 0 0 1) (V3 1 0 0) (V3 0 1 0) :: Double +-- 1.0 +-- +triple :: Ring a => V3 a -> V3 a -> V3 a -> a +triple x y z = dot x (cross y z) +{-# INLINE triple #-} + +------------------------------------------------------------------------------- +-- Vector Constructors & Acessors +------------------------------------------------------------------------------- + +-- | Dirac delta function. +-- +dirac :: Eq i => Semiring a => i -> i -> a +dirac i j = bool zero one (i == j) +{-# INLINE dirac #-} + +-- | Create a unit vector at an index. +-- +-- >>> idx E21 :: V2 Int +-- V2 1 0 +-- +-- >>> idx E42 :: V4 Int +-- V4 0 1 0 0 +-- +idx :: Semiring a => Basis b f => b -> f a +idx i = tabulate $ dirac i +{-# INLINE idx #-} + +-- | Retrieve an element of a vector. +-- +-- >>> elt E21 (V2 1 2) +-- 1 +-- +elt :: Basis b f => b -> f a -> a +elt = flip index +{-# INLINE elt #-} + +lensRep :: Basis b f => b -> forall g. Functor g => (a -> g a) -> f a -> (g**f) a +lensRep i f s = Compose $ setter s <$> f (getter s) + where getter = flip index i + setter s' b = tabulate $ \j -> bool (index s' j) b (i == j) +{-# INLINE lensRep #-} + +grateRep :: Basis b f => forall g. Functor g => (b -> g a1 -> a2) -> (g**f) a1 -> f a2 +grateRep iab s = tabulate $ \i -> iab i (fmap (`index` i) $ getCompose s) +{-# INLINE grateRep #-} + +------------------------------------------------------------------------------- +-- Matrix Arithmetic +------------------------------------------------------------------------------- + +infixr 7 .#. + +-- | Multiply two matrices. +-- +-- >>> m22 1 2 3 4 .#. m22 1 2 3 4 :: M22 Int +-- Compose (V2 (V2 7 10) (V2 15 22)) +-- +-- >>> m23 1 2 3 4 5 6 .#. m32 1 2 3 4 4 5 :: M22 Int +-- Compose (V2 (V2 19 25) (V2 43 58)) +-- +(.#.) :: Semiring a => Foldable g => Basis3 b c d f g h => (f**g) a -> (g**h) a -> (f**h) a +(.#.) x y = tabulate (\(i,j) -> row i x `dot` col j y) +{-# INLINE (.#.) #-} + +-- | Compute the trace of a matrix. +-- +-- >>> trace $ m22 1.0 2.0 3.0 4.0 +-- 5.0 +-- +trace :: Semiring a => Foldable f => Basis b f => (f**f) a -> a +trace = sum . diagonal +{-# INLINE trace #-} + +-- | 1x1 matrix inverse over a field. +-- +-- >>> inv1 $ m11 4.0 :: M11 Double +-- Compose (V1 (V1 0.25)) +-- +inv1 :: Field a => M11 a -> M11 a +inv1 = transpose . fmap recip + +-- | 2x2 matrix bdeterminant over a commutative semiring. +-- +-- >>> bdet2 $ m22 1 2 3 4 +-- (4,6) +-- +bdet2 :: Semiring a => Basis2 E2 E2 f g => (f**g) a -> (a, a) +bdet2 m = (elt2 E21 E21 m * elt2 E22 E22 m, elt2 E21 E22 m * elt2 E22 E21 m) +{-# INLINE bdet2 #-} + +-- | 2x2 matrix determinant over a commutative ring. +-- +-- @ +-- 'det2' = 'uncurry' ('-') . 'bdet2' +-- @ +-- +-- >>> det2 $ m22 1 2 3 4 :: Double +-- -2.0 +-- +det2 :: Ring a => Basis2 E2 E2 f g => (f**g) a -> a +det2 = uncurry (-) . bdet2 +{-# INLINE det2 #-} + +-- | 2x2 matrix inverse over a field. +-- +-- >>> inv2 $ m22 1 2 3 4 :: M22 Double +-- Compose (V2 (V2 (-2.0) 1.0) (V2 1.5 (-0.5))) +-- +inv2 :: Field a => M22 a -> M22 a +inv2 m = lscaleDef (recip $ det2 m) $ m22 d (-b) (-c) a where + a = elt2 E21 E21 m + b = elt2 E21 E22 m + c = elt2 E22 E21 m + d = elt2 E22 E22 m +{-# INLINE inv2 #-} + +-- | 3x3 matrix bdeterminant over a commutative semiring. +-- +-- >>> bdet3 (V3 (V3 1 2 3) (V3 4 5 6) (V3 7 8 9)) +-- (225, 225) +-- +bdet3 :: Semiring a => Basis2 E3 E3 f g => (f**g) a -> (a, a) +bdet3 m = (evens, odds) where + evens = a*e*i + g*b*f + d*h*c + odds = a*h*f + d*b*i + g*e*c + a = elt2 E31 E31 m + b = elt2 E31 E32 m + c = elt2 E31 E33 m + d = elt2 E32 E31 m + e = elt2 E32 E32 m + f = elt2 E32 E33 m + g = elt2 E33 E31 m + h = elt2 E33 E32 m + i = elt2 E33 E33 m +{-# INLINE bdet3 #-} + +-- | 3x3 double-precision matrix determinant. +-- +-- @ +-- 'det3' = 'uncurry' ('-') . 'bdet3' +-- @ +-- +-- Implementation uses a cofactor expansion to avoid loss of precision. +-- +-- >>> det3 $ m33 1 2 3 4 5 6 7 8 9 +-- 0 +-- +det3 :: Ring a => Basis2 E3 E3 f g => (f**g) a -> a +det3 m = a * (e*i-f*h) - d * (b*i-c*h) + g * (b*f-c*e) where + a = elt2 E31 E31 m + b = elt2 E31 E32 m + c = elt2 E31 E33 m + d = elt2 E32 E31 m + e = elt2 E32 E32 m + f = elt2 E32 E33 m + g = elt2 E33 E31 m + h = elt2 E33 E32 m + i = elt2 E33 E33 m +{-# INLINE det3 #-} + +-- | 3x3 matrix inverse. +-- +-- >>> inv3 $ m33 1 2 4 4 2 2 1 1 1 :: M33 Double +-- Compose (V3 (V3 0.0 0.5 (-1.0)) (V3 (-0.5) (-0.75) 3.5) (V3 0.5 0.25 (-1.5))) +-- +inv3 :: Field a => M33 a -> M33 a +inv3 m = lscaleDef (recip $ det3 m) $ m33 a' b' c' d' e' f' g' h' i' where + a = elt2 E31 E31 m + b = elt2 E31 E32 m + c = elt2 E31 E33 m + d = elt2 E32 E31 m + e = elt2 E32 E32 m + f = elt2 E32 E33 m + g = elt2 E33 E31 m + h = elt2 E33 E32 m + i = elt2 E33 E33 m + a' = cofactor (e,f,h,i) + b' = cofactor (c,b,i,h) + c' = cofactor (b,c,e,f) + d' = cofactor (f,d,i,g) + e' = cofactor (a,c,g,i) + f' = cofactor (c,a,f,d) + g' = cofactor (d,e,g,h) + h' = cofactor (b,a,h,g) + i' = cofactor (a,b,d,e) + cofactor (q,r,s,t) = det2 (m22 q r s t) +{-# INLINE inv3 #-} + +-- | 4x4 matrix bdeterminant over a commutative semiring. +-- +-- >>> bdet4 $ m44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 +-- (27728,27728) +-- +bdet4 :: Semiring a => Basis2 E4 E4 f g => (f**g) a -> (a, a) +bdet4 x = (evens, odds) where + evens = a * (f*k*p + g*l*n + h*j*o) + + b * (g*i*p + e*l*o + h*k*m) + + c * (e*j*p + f*l*m + h*i*n) + + d * (f*i*o + e*k*n + g*j*m) + odds = a * (g*j*p + f*l*o + h*k*n) + + b * (e*k*p + g*l*m + h*i*o) + + c * (f*i*p + e*l*n + h*j*m) + + d * (e*j*o + f*k*m + g*i*n) + a = elt2 E41 E41 x + b = elt2 E41 E42 x + c = elt2 E41 E43 x + d = elt2 E41 E44 x + e = elt2 E42 E41 x + f = elt2 E42 E42 x + g = elt2 E42 E43 x + h = elt2 E42 E44 x + i = elt2 E43 E41 x + j = elt2 E43 E42 x + k = elt2 E43 E43 x + l = elt2 E43 E44 x + m = elt2 E44 E41 x + n = elt2 E44 E42 x + o = elt2 E44 E43 x + p = elt2 E44 E44 x +{-# INLINE bdet4 #-} + +-- | 4x4 matrix determinant over a commutative ring. +-- +-- @ +-- 'det4' = 'uncurry' ('-') . 'bdet4' +-- @ +-- +-- This implementation uses a cofactor expansion to avoid loss of precision. +-- +-- >>> det4 $ m44 1 0 3 2 2 0 2 1 0 0 0 1 0 3 4 0 :: Rational +-- (-12) % 1 +-- +det4 :: Ring a => Basis2 E4 E4 f g => (f**g) a -> a +det4 x = s0 * c5 - s1 * c4 + s2 * c3 + s3 * c2 - s4 * c1 + s5 * c0 where + s0 = i00 * e11 - e10 * i01 + s1 = i00 * e12 - e10 * i02 + s2 = i00 * e13 - e10 * i03 + s3 = i01 * e12 - e11 * i02 + s4 = i01 * e13 - e11 * i03 + s5 = i02 * e13 - e12 * i03 + + c5 = e22 * e33 - e32 * e23 + c4 = e21 * e33 - e31 * e23 + c3 = e21 * e32 - e31 * e22 + c2 = e20 * e33 - e30 * e23 + c1 = e20 * e32 - e30 * e22 + c0 = e20 * e31 - e30 * e21 + + i00 = elt2 E41 E41 x + i01 = elt2 E41 E42 x + i02 = elt2 E41 E43 x + i03 = elt2 E41 E44 x + e10 = elt2 E42 E41 x + e11 = elt2 E42 E42 x + e12 = elt2 E42 E43 x + e13 = elt2 E42 E44 x + e20 = elt2 E43 E41 x + e21 = elt2 E43 E42 x + e22 = elt2 E43 E43 x + e23 = elt2 E43 E44 x + e30 = elt2 E44 E41 x + e31 = elt2 E44 E42 x + e32 = elt2 E44 E43 x + e33 = elt2 E44 E44 x +{-# INLINE det4 #-} + +-- | 4x4 matrix inverse. +-- +-- >>> row E41 . inv4 $ m44 1 0 3 2 2 0 2 1 0 0 0 1 0 3 4 0 :: V4 Rational +-- V4 (6 % (-12)) ((-9) % (-12)) ((-3) % (-12)) (0 % (-12)) +-- +inv4 :: Field a => M44 a -> M44 a +inv4 x = lscaleDef (recip det) $ x' where + i00 = elt2 E41 E41 x + i01 = elt2 E41 E42 x + i02 = elt2 E41 E43 x + i03 = elt2 E41 E44 x + e10 = elt2 E42 E41 x + e11 = elt2 E42 E42 x + e12 = elt2 E42 E43 x + e13 = elt2 E42 E44 x + e20 = elt2 E43 E41 x + e21 = elt2 E43 E42 x + e22 = elt2 E43 E43 x + e23 = elt2 E43 E44 x + e30 = elt2 E44 E41 x + e31 = elt2 E44 E42 x + e32 = elt2 E44 E43 x + e33 = elt2 E44 E44 x + + s0 = i00 * e11 - e10 * i01 + s1 = i00 * e12 - e10 * i02 + s2 = i00 * e13 - e10 * i03 + s3 = i01 * e12 - e11 * i02 + s4 = i01 * e13 - e11 * i03 + s5 = i02 * e13 - e12 * i03 + c5 = e22 * e33 - e32 * e23 + c4 = e21 * e33 - e31 * e23 + c3 = e21 * e32 - e31 * e22 + c2 = e20 * e33 - e30 * e23 + c1 = e20 * e32 - e30 * e22 + c0 = e20 * e31 - e30 * e21 + + det = s0 * c5 - s1 * c4 + s2 * c3 + s3 * c2 - s4 * c1 + s5 * c0 + + x' = m44 (e11 * c5 - e12 * c4 + e13 * c3) + (-i01 * c5 + i02 * c4 - i03 * c3) + (e31 * s5 - e32 * s4 + e33 * s3) + (-e21 * s5 + e22 * s4 - e23 * s3) + (-e10 * c5 + e12 * c2 - e13 * c1) + (i00 * c5 - i02 * c2 + i03 * c1) + (-e30 * s5 + e32 * s2 - e33 * s1) + (e20 * s5 - e22 * s2 + e23 * s1) + (e10 * c4 - e11 * c2 + e13 * c0) + (-i00 * c4 + i01 * c2 - i03 * c0) + (e30 * s4 - e31 * s2 + e33 * s0) + (-e20 * s4 + e21 * s2 - e23 * s0) + (-e10 * c3 + e11 * c1 - e12 * c0) + (i00 * c3 - i01 * c1 + i02 * c0) + (-e30 * s3 + e31 * s1 - e32 * s0) + (e20 * s3 - e21 * s1 + e22 * s0) +{-# INLINE inv4 #-} + +------------------------------------------------------------------------------- +-- Matrix constructors and accessors +------------------------------------------------------------------------------- + +-- | Lift a matrix into a linear transformation +-- +-- @ ('.#') = 'app' . 'tran' @ +-- +tran :: Semiring a => Basis2 b c f g => Foldable g => (f**g) a -> Tran a b c +tran m = Tran $ \f -> index $ m .# (tabulate f) + +-- | Retrieve an element of a matrix. +-- +-- >>> elt2 E21 E21 $ m22 1 2 3 4 +-- 1 +-- +elt2 :: Basis2 b c f g => b -> c -> (f**g) a -> a +elt2 i j = elt i . col j +{-# INLINE elt2 #-} + +-- | Retrieve a row of a matrix. +-- +-- >>> row E22 $ m23 1 2 3 4 5 6 +-- V3 4 5 6 +-- +row :: Basis b f => b -> (f**g) a -> g a +row i = elt i . getCompose +{-# INLINE row #-} + +-- | Retrieve a column of a matrix. +-- +-- >>> elt E22 . col E31 $ m23 1 2 3 4 5 6 +-- 4 +-- +col :: Basis2 b c f g => c -> (f**g) a -> f a +col j = elt j . distribute . getCompose +{-# INLINE col #-} + +-- | Obtain a diagonal matrix from a vector. +-- +-- >>> diag $ V2 2 3 +-- Compose (V2 (V2 2 0) (V2 0 3)) +-- +diag :: Semiring a => Basis b f => f a -> (f**f) a +diag f = Compose $ flip imapRep f $ \i x -> flip imapRep f (\j _ -> bool zero x $ i == j) +{-# INLINE diag #-} + +-- | Outer product of two vectors. +-- +-- >>> V2 1 1 `outer` V2 1 1 +-- Compose (V2 (V2 1 1) (V2 1 1)) +-- +outer :: Semiring a => Basis2 b c f g => f a -> g a -> (f**g) a +outer x y = Compose $ fmap (\z-> fmap (*z) y) x + +-- | Identity matrix. +-- +-- >>> identity :: M33 Int +-- Compose (V3 (V3 1 0 0) (V3 0 1 0) (V3 0 0 1)) +-- +identity :: Semiring a => Basis b f => (f**f) a +identity = diag $ pureRep one +{-# INLINE identity #-} + +-- | Obtain the diagonal of a matrix as a vector. +-- +-- >>> diagonal $ m22 1.0 2.0 3.0 4.0 +-- V2 1.0 4.0 +-- +diagonal :: Representable f => (f**f) a -> f a +diagonal = flip bindRep id . getCompose +{-# INLINE diagonal #-} + +-- | Construct a 1x1 matrix. +-- +-- >>> m11 1 :: M11 Int +-- Compose (V1 (V1 1)) +-- +m11 :: a -> M11 a +m11 a = Compose $ V1 (V1 a) +{-# INLINE m11 #-} + +-- | Construct a 1x2 matrix. +-- +-- >>> m12 1 2 :: M12 Int +-- Compose (V1 (V2 1 2)) +-- +m12 :: a -> a -> M12 a +m12 a b = Compose $ V1 (V2 a b) +{-# INLINE m12 #-} + +-- | Construct a 1x3 matrix. +-- +-- >>> m13 1 2 3 :: M13 Int +-- Compose (V1 (V3 1 2 3)) +-- +m13 :: a -> a -> a -> M13 a +m13 a b c = Compose $ V1 (V3 a b c) +{-# INLINE m13 #-} + +-- | Construct a 1x4 matrix. +-- +-- >>> m14 1 2 3 4 :: M14 Int +-- Compose (V1 (V4 1 2 3 4)) +-- +m14 :: a -> a -> a -> a -> M14 a +m14 a b c d = Compose $ V1 (V4 a b c d) +{-# INLINE m14 #-} + +-- | Construct a 2x1 matrix. +-- +-- >>> m21 1 2 :: M21 Int +-- Compose (V2 (V1 1) (V1 2)) +-- +m21 :: a -> a -> M21 a +m21 a b = Compose $ V2 (V1 a) (V1 b) +{-# INLINE m21 #-} + +-- | Construct a 3x1 matrix. +-- +-- >>> m31 1 2 3 :: M31 Int +-- Compose (V3 (V1 1) (V1 2) (V1 3)) +-- +m31 :: a -> a -> a -> M31 a +m31 a b c = Compose $ V3 (V1 a) (V1 b) (V1 c) +{-# INLINE m31 #-} + +-- | Construct a 4x1 matrix. +-- +-- >>> m41 1 2 3 4 :: M41 Int +-- Compose (V4 (V1 1) (V1 2) (V1 3) (V1 4)) +-- +m41 :: a -> a -> a -> a -> M41 a +m41 a b c d = Compose $ V4 (V1 a) (V1 b) (V1 c) (V1 d) +{-# INLINE m41 #-} + +-- | Construct a 2x2 matrix. +-- +-- Arguments are in row-major order. +-- +-- >>> m22 1 2 3 4 :: M22 Int +-- Compose (V2 (V2 1 2) (V2 3 4)) +-- +m22 :: a -> a -> a -> a -> M22 a +m22 a b c d = Compose $ V2 (V2 a b) (V2 c d) +{-# INLINE m22 #-} + +-- | Construct a 2x3 matrix. +-- +-- Arguments are in row-major order. +-- +m23 :: a -> a -> a -> a -> a -> a -> M23 a +m23 a b c d e f = Compose $ V2 (V3 a b c) (V3 d e f) +{-# INLINE m23 #-} + +-- | Construct a 2x4 matrix. +-- +-- Arguments are in row-major order. +-- +m24 :: a -> a -> a -> a -> a -> a -> a -> a -> M24 a +m24 a b c d e f g h = Compose $ V2 (V4 a b c d) (V4 e f g h) +{-# INLINE m24 #-} + +-- | Construct a 3x2 matrix. +-- +-- Arguments are in row-major order. +-- +m32 :: a -> a -> a -> a -> a -> a -> M32 a +m32 a b c d e f = Compose $ V3 (V2 a b) (V2 c d) (V2 e f) +{-# INLINE m32 #-} + +-- | Construct a 3x3 matrix. +-- +-- Arguments are in row-major order. +-- +m33 :: a -> a -> a -> a -> a -> a -> a -> a -> a -> M33 a +m33 a b c d e f g h i = Compose $ V3 (V3 a b c) (V3 d e f) (V3 g h i) +{-# INLINE m33 #-} + +-- | Construct a 3x4 matrix. +-- +-- Arguments are in row-major order. +-- +m34 :: a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> M34 a +m34 a b c d e f g h i j k l = Compose $ V3 (V4 a b c d) (V4 e f g h) (V4 i j k l) +{-# INLINE m34 #-} + +-- | Construct a 4x2 matrix. +-- +-- Arguments are in row-major order. +-- +m42 :: a -> a -> a -> a -> a -> a -> a -> a -> M42 a +m42 a b c d e f g h = Compose $ V4 (V2 a b) (V2 c d) (V2 e f) (V2 g h) +{-# INLINE m42 #-} + +-- | Construct a 4x3 matrix. +-- +-- Arguments are in row-major order. +-- +m43 :: a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> M43 a +m43 a b c d e f g h i j k l = Compose $ V4 (V3 a b c) (V3 d e f) (V3 g h i) (V3 j k l) +{-# INLINE m43 #-} + +-- | Construct a 4x4 matrix. +-- +-- Arguments are in row-major order. +-- +m44 :: a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> M44 a +m44 a b c d e f g h i j k l m n o p = Compose $ V4 (V4 a b c d) (V4 e f g h) (V4 i j k l) (V4 m n o p) +{-# INLINE m44 #-} + +------------------------------------------------------------------------------- +-- Matrix types +------------------------------------------------------------------------------- + +-- All matrices use row-major representation. + +-- | A 1x1 matrix. +type M11 = Compose V1 V1 + +-- | A 1x2 matrix. +type M12 = Compose V1 V2 + +-- | A 1x3 matrix. +type M13 = Compose V1 V3 + +-- | A 1x4 matrix. +type M14 = Compose V1 V4 + +-- | A 2x1 matrix. +type M21 = Compose V2 V1 + +-- | A 3x1 matrix. +type M31 = Compose V3 V1 + +-- | A 4x1 matrix. +type M41 = Compose V4 V1 + +-- | A 2x2 matrix. +type M22 = Compose V2 V2 + +-- | A 2x3 matrix. +type M23 = Compose V2 V3 + +-- | A 2x4 matrix. +type M24 = Compose V2 V4 + +-- | A 3x2 matrix. +type M32 = Compose V3 V2 + +-- | A 3x3 matrix. +type M33 = Compose V3 V3 + +-- | A 3x4 matrix. +type M34 = Compose V3 V4 + +-- | A 4x2 matrix. +type M42 = Compose V4 V2 + +-- | A 4x3 matrix. +type M43 = Compose V4 V3 + +-- | A 4x4 matrix. +type M44 = Compose V4 V4 + + + +------------------------------------------------------------------------------- +-- V1 +------------------------------------------------------------------------------- + +unV1 :: V1 a -> a +unV1 (V1 a) = a + +newtype V1 a = V1 a deriving (Eq,Ord,Show) + +instance Show1 V1 where + liftShowsPrec f _ d (V1 a) = showParen (d >= 10) $ showString "V1 " . f d a + +{- +instance Field a => Composition a V1 where + conj = id + + norm f = unV1 $ liftA2 (*) f f +-} + +instance Functor V1 where + fmap f (V1 a) = V1 (f a) + {-# INLINE fmap #-} + a <$ _ = V1 a + {-# INLINE (<$) #-} + +instance Applicative V1 where + pure = pureRep + liftA2 = liftR2 + +instance Foldable V1 where + foldMap f (V1 a) = f a + {-# INLINE foldMap #-} + null _ = False + length _ = one + +instance Foldable1 V1 where + foldMap1 f (V1 a) = f a + {-# INLINE foldMap1 #-} + +instance Distributive V1 where + distribute f = V1 $ fmap (\(V1 x) -> x) f + {-# INLINE distribute #-} + +instance Representable V1 where + type Rep V1 = E1 + tabulate f = V1 (f E1) + {-# INLINE tabulate #-} + + index (V1 x) E1 = x + {-# INLINE index #-} + +------------------------------------------------------------------------------- +-- V2 +------------------------------------------------------------------------------- + +data V2 a = V2 !a !a deriving (Eq,Ord,Show) + +instance Show1 V2 where + liftShowsPrec f _ d (V2 a b) = showsBinaryWith f f "V2" d a b + +instance Functor V2 where + fmap f (V2 a b) = V2 (f a) (f b) + {-# INLINE fmap #-} + a <$ _ = V2 a a + {-# INLINE (<$) #-} + +instance Applicative V2 where + pure = pureRep + liftA2 = liftR2 + +instance Foldable V2 where + foldMap f (V2 a b) = f a <> f b + {-# INLINE foldMap #-} + null _ = False + length _ = two + +instance Foldable1 V2 where + foldMap1 f (V2 a b) = f a <> f b + {-# INLINE foldMap1 #-} + +instance Distributive V2 where + distribute f = V2 (fmap (\(V2 x _) -> x) f) (fmap (\(V2 _ y) -> y) f) + {-# INLINE distribute #-} + +instance Representable V2 where + type Rep V2 = E2 + tabulate f = V2 (f E21) (f E22) + {-# INLINE tabulate #-} + + index (V2 x _) E21 = x + index (V2 _ y) E22 = y + {-# INLINE index #-} + +------------------------------------------------------------------------------- +-- V3 +------------------------------------------------------------------------------- + + +data V3 a = V3 !a !a !a deriving (Eq,Ord,Show) + +-- TODO add Prd1 and push instance downstream +instance Eq1 V3 where + liftEq k (V3 a b c) (V3 d e f) = k a d && k b e && k c f + +instance Show1 V3 where + liftShowsPrec f _ d (V3 a b c) = showParen (d > 10) $ + showString "V3 " . f 11 a . showChar ' ' . f 11 b . showChar ' ' . f 11 c + +instance Functor V3 where + fmap f (V3 a b c) = V3 (f a) (f b) (f c) + {-# INLINE fmap #-} + a <$ _ = V3 a a a + {-# INLINE (<$) #-} + +instance Applicative V3 where + pure = pureRep + liftA2 = liftR2 + +instance Foldable V3 where + foldMap f (V3 a b c) = f a <> f b <> f c + {-# INLINE foldMap #-} + null _ = False + --length _ = 3 + +instance Foldable1 V3 where + foldMap1 f (V3 a b c) = f a <> f b <> f c + {-# INLINE foldMap1 #-} + +instance Distributive V3 where + distribute f = V3 (fmap (\(V3 x _ _) -> x) f) (fmap (\(V3 _ y _) -> y) f) (fmap (\(V3 _ _ z) -> z) f) + {-# INLINE distribute #-} + +instance Representable V3 where + type Rep V3 = E3 + tabulate f = V3 (f E31) (f E32) (f E33) + {-# INLINE tabulate #-} + + index (V3 x _ _) E31 = x + index (V3 _ y _) E32 = y + index (V3 _ _ z) E33 = z + {-# INLINE index #-} + + + +------------------------------------------------------------------------------- +-- V4 +------------------------------------------------------------------------------- + +data V4 a = V4 !a !a !a !a deriving (Eq,Ord,Show) + +instance Show1 V4 where + liftShowsPrec f _ z (V4 a b c d) = showParen (z > 10) $ + showString "V4 " . f 11 a . showChar ' ' . f 11 b . showChar ' ' . f 11 c . showChar ' ' . f 11 d + +instance Functor V4 where + fmap f (V4 a b c d) = V4 (f a) (f b) (f c) (f d) + {-# INLINE fmap #-} + a <$ _ = V4 a a a a + {-# INLINE (<$) #-} + +instance Applicative V4 where + pure = pureRep + liftA2 = liftR2 + +instance Foldable V4 where + foldMap f (V4 a b c d) = f a <> f b <> f c <> f d + {-# INLINE foldMap #-} + null _ = False + length _ = two + two + +instance Foldable1 V4 where + foldMap1 f (V4 a b c d) = f a <> f b <> f c <> f d + {-# INLINE foldMap1 #-} + +instance Distributive V4 where + distribute f = V4 (fmap (\(V4 x _ _ _) -> x) f) (fmap (\(V4 _ y _ _) -> y) f) (fmap (\(V4 _ _ z _) -> z) f) (fmap (\(V4 _ _ _ w) -> w) f) + {-# INLINE distribute #-} + +instance Representable V4 where + type Rep V4 = E4 + tabulate f = V4 (f E41) (f E42) (f E43) (f E44) + {-# INLINE tabulate #-} + + index (V4 x _ _ _) E41 = x + index (V4 _ y _ _) E42 = y + index (V4 _ _ z _) E43 = z + index (V4 _ _ _ w) E44 = w + {-# INLINE index #-} + + +------------------------------------------------------------------------------- +-- Autogenerated instances +------------------------------------------------------------------------------- + + +#define deriveAdditiveSemigroup(ty) \ +instance (Additive-Semigroup) a => Semigroup (Additive (ty a)) where { \ + (<>) = liftA2 $ mzipWithRep (+) \ +; {-# INLINE (<>) #-} \ +} + +#define deriveAdditiveMonoid(ty) \ +instance (Additive-Monoid) a => Monoid (Additive (ty a)) where { \ + mempty = pure $ pureRep zero \ +; {-# INLINE mempty #-} \ +} + +#define deriveMultiplicativeSemigroup(ty) \ +instance (Multiplicative-Semigroup) a => Semigroup (Multiplicative (ty a)) where { \ + (<>) = liftA2 $ mzipWithRep (*) \ +; {-# INLINE (<>) #-} \ +} + +#define deriveMultiplicativeMonoid(ty) \ +instance (Multiplicative-Monoid) a => Monoid (Multiplicative (ty a)) where { \ + mempty = pure $ pureRep one \ +; {-# INLINE mempty #-} \ +} + +#define deriveMultiplicativeMatrixSemigroup(ty) \ +instance Semiring a => Semigroup (Multiplicative (ty a)) where { \ + (<>) = liftA2 $ (.#.) \ +; {-# INLINE (<>) #-} \ +} + +#define deriveMultiplicativeMatrixMonoid(ty) \ +instance Semiring a => Monoid (Multiplicative (ty a)) where { \ + mempty = pure identity \ +; {-# INLINE mempty #-} \ +} + +#define deriveAdditiveMagma(ty) \ +instance (Additive-Group) a => Magma (Additive (ty a)) where { \ + (<<) = liftA2 $ mzipWithRep (-) \ +; {-# INLINE (<<) #-} \ +} + +#define deriveAdditiveQuasigroup(ty) \ +instance (Additive-Group) a => Quasigroup (Additive (ty a)) \ + +#define deriveAdditiveLoop(ty) \ +instance (Additive-Group) a => Loop (Additive (ty a)) \ + +#define deriveAdditiveGroup(ty) \ +instance (Additive-Group) a => Group (Additive (ty a)) \ + +#define derivePresemiring(ty) \ +instance Semiring a => Presemiring (ty a) \ + +#define deriveSemiring(ty) \ +instance Semiring a => Semiring (ty a) \ + +#define deriveRing(ty) \ +instance Ring a => Ring (ty a) \ + +#define deriveFreeLeftSemimodule(ty) \ +instance Semiring a => LeftSemimodule a (ty a) where { \ + lscale = lscaleDef \ +; {-# INLINE lscale #-} \ +} + +#define deriveFreeRightSemimodule(ty) \ +instance Semiring a => RightSemimodule a (ty a) where { \ + rscale = rscaleDef \ +; {-# INLINE rscale #-} \ +} + +#define deriveFreeBisemimodule(ty) \ +instance Semiring a => Bisemimodule a a (ty a) \ + +#define deriveBisemimodule(tyl, tyr, ty) \ +instance Semiring a => Bisemimodule (tyl a) (tyr a) (ty a) \ + +#define deriveLeftSemimodule(tyl,ty) \ +instance Semiring a => LeftSemimodule (tyl a) (ty a) where { \ + lscale = (.#.) \ +; {-# INLINE lscale #-} \ +} + +#define deriveRightSemimodule(tyr,ty) \ +instance Semiring a => RightSemimodule (tyr a) (ty a) where { \ + rscale = flip (.#.) \ +; {-# INLINE rscale #-} \ +} + +#define deriveBisemimodule(tyl, tyr, ty) \ +instance Semiring a => Bisemimodule (tyl a) (tyr a) (ty a) \ + + + +-- V1 +deriveAdditiveSemigroup(V1) +deriveAdditiveMonoid(V1) + +deriveAdditiveMagma(V1) +deriveAdditiveQuasigroup(V1) +deriveAdditiveLoop(V1) +deriveAdditiveGroup(V1) + +deriveFreeLeftSemimodule(V1) +deriveFreeRightSemimodule(V1) +deriveFreeBisemimodule(V1) + + +-- V2 +deriveAdditiveSemigroup(V2) +deriveAdditiveMonoid(V2) + +deriveAdditiveMagma(V2) +deriveAdditiveQuasigroup(V2) +deriveAdditiveLoop(V2) +deriveAdditiveGroup(V2) + +deriveFreeLeftSemimodule(V2) +deriveFreeRightSemimodule(V2) +deriveFreeBisemimodule(V2) + + +-- V3 +deriveAdditiveSemigroup(V3) +deriveAdditiveMonoid(V3) + +deriveAdditiveMagma(V3) +deriveAdditiveQuasigroup(V3) +deriveAdditiveLoop(V3) +deriveAdditiveGroup(V3) + +deriveFreeLeftSemimodule(V3) +deriveFreeRightSemimodule(V3) +deriveFreeBisemimodule(V3) + +-- V4 +deriveAdditiveSemigroup(V4) +deriveAdditiveMonoid(V4) + +deriveAdditiveMagma(V4) +deriveAdditiveQuasigroup(V4) +deriveAdditiveLoop(V4) +deriveAdditiveGroup(V4) + +deriveFreeLeftSemimodule(V4) +deriveFreeRightSemimodule(V4) +deriveFreeBisemimodule(V4) + +-- M11 +deriveAdditiveSemigroup(M11) +deriveAdditiveMonoid(M11) + +deriveAdditiveMagma(M11) +deriveAdditiveQuasigroup(M11) +deriveAdditiveLoop(M11) +deriveAdditiveGroup(M11) + +deriveLeftSemimodule(M11, M11) +deriveRightSemimodule(M11, M11) +deriveBisemimodule(M11, M11, M11) + +deriveMultiplicativeMatrixSemigroup(M11) +deriveMultiplicativeMatrixMonoid(M11) + +derivePresemiring(M11) +deriveSemiring(M11) +deriveRing(M11) + +-- M21 +deriveAdditiveSemigroup(M21) +deriveAdditiveMonoid(M21) + +deriveAdditiveMagma(M21) +deriveAdditiveQuasigroup(M21) +deriveAdditiveLoop(M21) +deriveAdditiveGroup(M21) + +deriveLeftSemimodule(M22, M21) +deriveRightSemimodule(M11, M21) +deriveBisemimodule(M22, M11, M21) + + +-- M31 +deriveAdditiveSemigroup(M31) +deriveAdditiveMonoid(M31) + +deriveAdditiveMagma(M31) +deriveAdditiveQuasigroup(M31) +deriveAdditiveLoop(M31) +deriveAdditiveGroup(M31) + +deriveLeftSemimodule(M33, M31) +deriveRightSemimodule(M11, M31) +deriveBisemimodule(M33, M11, M31) + + +-- M41 +deriveAdditiveSemigroup(M41) +deriveAdditiveMonoid(M41) + +deriveAdditiveMagma(M41) +deriveAdditiveQuasigroup(M41) +deriveAdditiveLoop(M41) +deriveAdditiveGroup(M41) + +deriveLeftSemimodule(M44, M41) +deriveRightSemimodule(M11, M41) +deriveBisemimodule(M44, M11, M41) + + +-- M12 +deriveAdditiveSemigroup(M12) +deriveAdditiveMonoid(M12) + +deriveAdditiveMagma(M12) +deriveAdditiveQuasigroup(M12) +deriveAdditiveLoop(M12) +deriveAdditiveGroup(M12) + +deriveLeftSemimodule(M11, M12) +deriveRightSemimodule(M22, M12) +deriveBisemimodule(M11, M22, M12) + + +-- M22 +deriveAdditiveSemigroup(M22) +deriveAdditiveMonoid(M22) + +deriveAdditiveMagma(M22) +deriveAdditiveQuasigroup(M22) +deriveAdditiveLoop(M22) +deriveAdditiveGroup(M22) + +deriveLeftSemimodule(M22, M22) +deriveRightSemimodule(M22, M22) +deriveBisemimodule(M22, M22, M22) + +deriveMultiplicativeMatrixSemigroup(M22) +deriveMultiplicativeMatrixMonoid(M22) + +derivePresemiring(M22) +deriveSemiring(M22) +deriveRing(M22) + + +-- M32 +deriveAdditiveSemigroup(M32) +deriveAdditiveMonoid(M32) + +deriveAdditiveMagma(M32) +deriveAdditiveQuasigroup(M32) +deriveAdditiveLoop(M32) +deriveAdditiveGroup(M32) + +deriveLeftSemimodule(M33, M32) +deriveRightSemimodule(M22, M32) +deriveBisemimodule(M33, M22, M32) + + +-- M42 +deriveAdditiveSemigroup(M42) +deriveAdditiveMonoid(M42) + +deriveAdditiveMagma(M42) +deriveAdditiveQuasigroup(M42) +deriveAdditiveLoop(M42) +deriveAdditiveGroup(M42) + +deriveLeftSemimodule(M44, M42) +deriveRightSemimodule(M22, M42) +deriveBisemimodule(M44, M22, M42) + + +-- M13 +deriveAdditiveSemigroup(M13) +deriveAdditiveMonoid(M13) + +deriveAdditiveMagma(M13) +deriveAdditiveQuasigroup(M13) +deriveAdditiveLoop(M13) +deriveAdditiveGroup(M13) + +deriveLeftSemimodule(M11, M13) +deriveRightSemimodule(M33, M13) +deriveBisemimodule(M11, M33, M13) + + +-- M23 +deriveAdditiveSemigroup(M23) +deriveAdditiveMonoid(M23) + +deriveAdditiveMagma(M23) +deriveAdditiveQuasigroup(M23) +deriveAdditiveLoop(M23) +deriveAdditiveGroup(M23) + +deriveLeftSemimodule(M22, M23) +deriveRightSemimodule(M33, M23) +deriveBisemimodule(M22, M33, M23) + + +-- M33 +deriveAdditiveSemigroup(M33) +deriveAdditiveMonoid(M33) + +deriveAdditiveMagma(M33) +deriveAdditiveQuasigroup(M33) +deriveAdditiveLoop(M33) +deriveAdditiveGroup(M33) + +deriveLeftSemimodule(M33, M33) +deriveRightSemimodule(M33, M33) +deriveBisemimodule(M33, M33, M33) + +deriveMultiplicativeMatrixSemigroup(M33) +deriveMultiplicativeMatrixMonoid(M33) + +derivePresemiring(M33) +deriveSemiring(M33) +deriveRing(M33) + + +-- M43 +deriveAdditiveSemigroup(M43) +deriveAdditiveMonoid(M43) + +deriveAdditiveMagma(M43) +deriveAdditiveQuasigroup(M43) +deriveAdditiveLoop(M43) +deriveAdditiveGroup(M43) + +deriveLeftSemimodule(M44, M43) +deriveRightSemimodule(M33, M43) +deriveBisemimodule(M44, M33, M43) + + +-- M14 +deriveAdditiveSemigroup(M14) +deriveAdditiveMonoid(M14) + +deriveAdditiveMagma(M14) +deriveAdditiveQuasigroup(M14) +deriveAdditiveLoop(M14) +deriveAdditiveGroup(M14) + +deriveLeftSemimodule(M11, M14) +deriveRightSemimodule(M44, M14) +deriveBisemimodule(M11, M44, M14) + + +-- M24 +deriveAdditiveSemigroup(M24) +deriveAdditiveMonoid(M24) + +deriveAdditiveMagma(M24) +deriveAdditiveQuasigroup(M24) +deriveAdditiveLoop(M24) +deriveAdditiveGroup(M24) + +deriveLeftSemimodule(M22, M24) +deriveRightSemimodule(M44, M24) +deriveBisemimodule(M22, M44, M24) + + +-- M34 +deriveAdditiveSemigroup(M34) +deriveAdditiveMonoid(M34) + +deriveAdditiveMagma(M34) +deriveAdditiveQuasigroup(M34) +deriveAdditiveLoop(M34) +deriveAdditiveGroup(M34) + +deriveLeftSemimodule(M33, M34) +deriveRightSemimodule(M44, M34) +deriveBisemimodule(M33, M44, M34) + + +-- M44 +deriveAdditiveSemigroup(M44) +deriveAdditiveMonoid(M44) + +deriveAdditiveMagma(M44) +deriveAdditiveQuasigroup(M44) +deriveAdditiveLoop(M44) +deriveAdditiveGroup(M44) + +deriveLeftSemimodule(M44, M44) +deriveRightSemimodule(M44, M44) +deriveBisemimodule(M44, M44, M44) + +deriveMultiplicativeMatrixSemigroup(M44) +deriveMultiplicativeMatrixMonoid(M44) + +derivePresemiring(M44) +deriveSemiring(M44) +deriveRing(M44) diff --git a/src/Data/Semimodule/Matrix.hs b/src/Data/Semimodule/Matrix.hs deleted file mode 100644 --- a/src/Data/Semimodule/Matrix.hs +++ /dev/null @@ -1,552 +0,0 @@ -{-# LANGUAGE CPP #-} -{-# LANGUAGE Safe #-} -{-# LANGUAGE PolyKinds #-} -{-# LANGUAGE ConstraintKinds #-} -{-# LANGUAGE DefaultSignatures #-} -{-# LANGUAGE DeriveFunctor #-} -{-# LANGUAGE DeriveGeneric #-} -{-# LANGUAGE FlexibleContexts #-} -{-# LANGUAGE FlexibleInstances #-} -{-# LANGUAGE TypeOperators #-} -{-# LANGUAGE TypeFamilies #-} -{-# LANGUAGE RebindableSyntax #-} -{-# LANGUAGE RankNTypes #-} - -module Data.Semimodule.Matrix ( - type M22 - , type M23 - , type M24 - , type M32 - , type M33 - , type M34 - , type M42 - , type M43 - , type M44 - , lensRep - , grateRep - , tran - , row - , rows - , col - , cols - , (.#) - , (.*) - , (#.) - , (*.) - , (.#.) - , (.*.) - , outer - , scale - , dirac - , identity - , transpose - , trace - , diagonal - , bdet2 - , det2 - , inv2 - , bdet3 - , det3 - , inv3 - , bdet4 - , det4 - , inv4 - , m22 - , m23 - , m24 - , m32 - , m33 - , m34 - , m42 - , m43 - , m44 - ) where - -import safe Data.Bool -import safe Data.Distributive -import safe Data.Functor.Compose -import safe Data.Functor.Rep -import safe Data.Semifield -import safe Data.Semigroup.Additive -import safe Data.Semigroup.Multiplicative -import safe Data.Semimodule -import safe Data.Semimodule.Transform -import safe Data.Semimodule.Vector -import safe Data.Semiring -import safe Data.Tuple -import safe Prelude hiding (Num(..), Fractional(..), sum, negate) - - --- All matrices use row-major representation. - --- | A 2x2 matrix. -type M22 a = V2 (V2 a) - --- | A 2x3 matrix. -type M23 a = V2 (V3 a) - --- | A 2x4 matrix. -type M24 a = V2 (V4 a) - --- | A 3x2 matrix. -type M32 a = V3 (V2 a) - --- | A 3x3 matrix. -type M33 a = V3 (V3 a) - --- | A 3x4 matrix. -type M34 a = V3 (V4 a) - --- | A 4x2 matrix. -type M42 a = V4 (V2 a) - --- | A 4x3 matrix. -type M43 a = V4 (V3 a) - --- | A 4x4 matrix. -type M44 a = V4 (V4 a) - -lensRep :: Eq (Rep f) => Representable f => Rep f -> forall g. Functor g => (a -> g a) -> f a -> g (f a) -lensRep i f s = setter s <$> f (getter s) - where getter = flip index i - setter s' b = tabulate $ \j -> bool (index s' j) b (i == j) -{-# INLINE lensRep #-} - -grateRep :: Representable f => forall g. Functor g => (Rep f -> g a -> b) -> g (f a) -> f b -grateRep iab s = tabulate $ \i -> iab i (fmap (`index` i) s) -{-# INLINE grateRep #-} - --- @ ('.#') = 'app' . 'tran' @ -tran :: Semiring a => Basis b f => Basis c g => Foldable g => f (g a) -> Tran a b c -tran m = Tran $ \f -> index $ m .# (tabulate f) - --- | Retrieve a row of a row-major matrix or element of a row vector. --- --- >>> row I21 (V2 1 2) --- 1 --- -row :: Representable f => Rep f -> f a -> a -row = flip index -{-# INLINE row #-} - --- | Retrieve a column of a row-major matrix. --- --- >>> row I22 . col I31 $ V2 (V3 1 2 3) (V3 4 5 6) --- 4 --- -col :: Functor f => Representable g => Rep g -> f (g a) -> f a -col j = flip index j . distribute -{-# INLINE col #-} - --- | Outer product of two vectors. --- --- >>> V2 1 1 `outer` V2 1 1 --- V2 (V2 1 1) (V2 1 1) --- -outer :: Semiring a => Functor f => Functor g => f a -> g a -> f (g a) -outer x y = fmap (\z-> fmap (*z) y) x - -infixl 7 #. - --- | Multiply a matrix on the left by a row vector. --- --- >>> V2 1 2 #. m23 3 4 5 6 7 8 --- V3 15 18 21 --- --- >>> V2 1 2 #. m23 3 4 5 6 7 8 #. m32 1 0 0 0 0 0 --- V2 15 0 --- -(#.) :: (Semiring a, Free f, Foldable f, Free g) => f a -> f (g a) -> g a -x #. y = tabulate (\j -> x .*. col j y) -{-# INLINE (#.) #-} - -infixr 7 .#, .#. - --- | Multiply a matrix on the right by a column vector. --- --- @ ('.#') = 'app' . 'fromMatrix' @ --- --- >>> m23 1 2 3 4 5 6 .# V3 7 8 9 --- V2 50 122 --- --- >>> m22 1 0 0 0 .# m23 1 2 3 4 5 6 .# V3 7 8 9 --- V2 50 0 --- -(.#) :: (Semiring a, Free f, Free g, Foldable g) => f (g a) -> g a -> f a -x .# y = tabulate (\i -> row i x .*. y) -{-# INLINE (.#) #-} - --- | Multiply two matrices. --- --- >>> m22 1 2 3 4 .#. m22 1 2 3 4 :: M22 Int --- V2 (V2 7 10) (V2 15 22) --- --- >>> m23 1 2 3 4 5 6 .#. m32 1 2 3 4 4 5 :: M22 Int --- V2 (V2 19 25) (V2 43 58) --- -(.#.) :: (Semiring a, Free f, Free g, Free h, Foldable g) => f (g a) -> g (h a) -> f (h a) -(.#.) x y = getCompose $ tabulate (\(i,j) -> row i x .*. col j y) -{-# INLINE (.#.) #-} - --- | Obtain a diagonal matrix from a vector. --- --- >>> scale (V2 2 3) --- V2 (V2 2 0) (V2 0 3) --- -scale :: (Additive-Monoid) a => Free f => f a -> f (f a) -scale f = flip imapRep f $ \i x -> flip imapRep f (\j _ -> bool zero x $ i == j) -{-# INLINE scale #-} - --- | Identity matrix. --- --- >>> identity :: M44 Int --- V4 (V4 1 0 0 0) (V4 0 1 0 0) (V4 0 0 1 0) (V4 0 0 0 1) --- --- >>> identity :: V3 (V3 Int) --- V3 (V3 1 0 0) (V3 0 1 0) (V3 0 0 1) --- -identity :: Semiring a => Free f => f (f a) -identity = scale $ pureRep one -{-# INLINE identity #-} - --- | Compute the trace of a matrix. --- --- >>> trace (V2 (V2 a b) (V2 c d)) --- a <> d --- -trace :: Semiring a => Free f => Foldable f => f (f a) -> a -trace = sum . diagonal -{-# INLINE trace #-} - --- | Obtain the diagonal of a matrix as a vector. --- --- >>> diagonal (V2 (V2 a b) (V2 c d)) --- V2 a d --- -diagonal :: Representable f => f (f a) -> f a -diagonal = flip bindRep id -{-# INLINE diagonal #-} - -ij :: Representable f => Representable g => Rep f -> Rep g -> f (g a) -> a -ij i j = row i . col j - --- | 2x2 matrix bdeterminant over a commutative semiring. --- --- >>> bdet2 $ m22 1 2 3 4 --- (4,6) --- -bdet2 :: Semiring a => Basis I2 f => Basis I2 g => f (g a) -> (a, a) -bdet2 m = (ij I21 I21 m * ij I22 I22 m, ij I21 I22 m * ij I22 I21 m) -{-# INLINE bdet2 #-} - --- | 2x2 matrix determinant over a commutative ring. --- --- @ --- 'det2' '==' 'uncurry' ('-') . 'bdet2' --- @ --- --- >>> det2 $ m22 1 2 3 4 :: Double --- -2.0 --- -det2 :: Ring a => Basis I2 f => Basis I2 g => f (g a) -> a -det2 = uncurry (-) . bdet2 -{-# INLINE det2 #-} - --- | 2x2 matrix inverse over a field. --- --- >>> inv2 $ m22 1 2 3 4 :: M22 Double --- V2 (V2 (-2.0) 1.0) (V2 1.5 (-0.5)) --- -inv2 :: Field a => Basis I2 f => Basis I2 g => f (g a) -> g (f a) -inv2 m = multl (recip $ det2 m) <$> m22 d (-b) (-c) a where - a = ij I21 I21 m - b = ij I21 I22 m - c = ij I22 I21 m - d = ij I22 I22 m -{-# INLINE inv2 #-} - --- | 3x3 matrix bdeterminant over a commutative semiring. --- --- >>> bdet3 (V3 (V3 1 2 3) (V3 4 5 6) (V3 7 8 9)) --- (225, 225) --- -bdet3 :: Semiring a => Basis I3 f => Basis I3 g => f (g a) -> (a, a) -bdet3 m = (evens, odds) where - evens = a*e*i + g*b*f + d*h*c - odds = a*h*f + d*b*i + g*e*c - a = ij I31 I31 m - b = ij I31 I32 m - c = ij I31 I33 m - d = ij I32 I31 m - e = ij I32 I32 m - f = ij I32 I33 m - g = ij I33 I31 m - h = ij I33 I32 m - i = ij I33 I33 m -{-# INLINE bdet3 #-} - --- | 3x3 double-precision matrix determinant. --- --- @ --- 'det3' '==' 'uncurry' ('-') . 'bdet3' --- @ --- --- Implementation uses a cofactor expansion to avoid loss of precision. --- --- >>> det3 (V3 (V3 1 2 3) (V3 4 5 6) (V3 7 8 9)) --- 0 --- -det3 :: Ring a => Basis I3 f => Basis I3 g => f (g a) -> a -det3 m = a * (e*i-f*h) - d * (b*i-c*h) + g * (b*f-c*e) where - a = ij I31 I31 m - b = ij I31 I32 m - c = ij I31 I33 m - d = ij I32 I31 m - e = ij I32 I32 m - f = ij I32 I33 m - g = ij I33 I31 m - h = ij I33 I32 m - i = ij I33 I33 m -{-# INLINE det3 #-} - --- | 3x3 matrix inverse. --- --- >>> inv3 $ m33 1 2 4 4 2 2 1 1 1 :: M33 Double --- V3 (V3 0.0 0.5 (-1.0)) (V3 (-0.5) (-0.75) 3.5) (V3 0.5 0.25 (-1.5)) --- -inv3 :: forall a f g. Field a => Basis I3 f => Basis I3 g => f (g a) -> g (f a) -inv3 m = multl (recip $ det3 m) <$> m33 a' b' c' d' e' f' g' h' i' where - a = ij I31 I31 m - b = ij I31 I32 m - c = ij I31 I33 m - d = ij I32 I31 m - e = ij I32 I32 m - f = ij I32 I33 m - g = ij I33 I31 m - h = ij I33 I32 m - i = ij I33 I33 m - a' = cofactor (e,f,h,i) - b' = cofactor (c,b,i,h) - c' = cofactor (b,c,e,f) - d' = cofactor (f,d,i,g) - e' = cofactor (a,c,g,i) - f' = cofactor (c,a,f,d) - g' = cofactor (d,e,g,h) - h' = cofactor (b,a,h,g) - i' = cofactor (a,b,d,e) - cofactor (q,r,s,t) = det2 (m22 q r s t :: M22 a) -{-# INLINE inv3 #-} - --- | 4x4 matrix bdeterminant over a commutative semiring. --- --- >>> bdet4 (V4 (V4 1 2 3 4) (V4 5 6 7 8) (V4 9 10 11 12) (V4 13 14 15 16)) --- (27728,27728) --- -bdet4 :: Semiring a => Basis I4 f => Basis I4 g => f (g a) -> (a, a) -bdet4 x = (evens, odds) where - evens = a * (f*k*p + g*l*n + h*j*o) + - b * (g*i*p + e*l*o + h*k*m) + - c * (e*j*p + f*l*m + h*i*n) + - d * (f*i*o + e*k*n + g*j*m) - odds = a * (g*j*p + f*l*o + h*k*n) + - b * (e*k*p + g*l*m + h*i*o) + - c * (f*i*p + e*l*n + h*j*m) + - d * (e*j*o + f*k*m + g*i*n) - a = ij I41 I41 x - b = ij I41 I42 x - c = ij I41 I43 x - d = ij I41 I44 x - e = ij I42 I41 x - f = ij I42 I42 x - g = ij I42 I43 x - h = ij I42 I44 x - i = ij I43 I41 x - j = ij I43 I42 x - k = ij I43 I43 x - l = ij I43 I44 x - m = ij I44 I41 x - n = ij I44 I42 x - o = ij I44 I43 x - p = ij I44 I44 x -{-# INLINE bdet4 #-} - --- | 4x4 matrix determinant over a commutative ring. --- --- @ --- 'det4' '==' 'uncurry' ('-') . 'bdet4' --- @ --- --- This implementation uses a cofactor expansion to avoid loss of precision. --- --- >>> det4 (m44 1 0 3 2 2 0 2 1 0 0 0 1 0 3 4 0 :: M44 Rational) --- (-12) % 1 --- -det4 :: Ring a => Basis I4 f => Basis I4 g => f (g a) -> a -det4 x = s0 * c5 - s1 * c4 + s2 * c3 + s3 * c2 - s4 * c1 + s5 * c0 where - s0 = i00 * i11 - i10 * i01 - s1 = i00 * i12 - i10 * i02 - s2 = i00 * i13 - i10 * i03 - s3 = i01 * i12 - i11 * i02 - s4 = i01 * i13 - i11 * i03 - s5 = i02 * i13 - i12 * i03 - - c5 = i22 * i33 - i32 * i23 - c4 = i21 * i33 - i31 * i23 - c3 = i21 * i32 - i31 * i22 - c2 = i20 * i33 - i30 * i23 - c1 = i20 * i32 - i30 * i22 - c0 = i20 * i31 - i30 * i21 - - i00 = ij I41 I41 x - i01 = ij I41 I42 x - i02 = ij I41 I43 x - i03 = ij I41 I44 x - i10 = ij I42 I41 x - i11 = ij I42 I42 x - i12 = ij I42 I43 x - i13 = ij I42 I44 x - i20 = ij I43 I41 x - i21 = ij I43 I42 x - i22 = ij I43 I43 x - i23 = ij I43 I44 x - i30 = ij I44 I41 x - i31 = ij I44 I42 x - i32 = ij I44 I43 x - i33 = ij I44 I44 x -{-# INLINE det4 #-} - --- | 4x4 matrix inverse. --- --- >>> row I41 $ inv4 (m44 1 0 3 2 2 0 2 1 0 0 0 1 0 3 4 0 :: M44 Rational) --- V4 (6 % (-12)) ((-9) % (-12)) ((-3) % (-12)) (0 % (-12)) --- -inv4 :: forall a f g. Field a => Basis I4 f => Basis I4 g => f (g a) -> g (f a) -inv4 x = multl (recip det) <$> x' where - i00 = ij I41 I41 x - i01 = ij I41 I42 x - i02 = ij I41 I43 x - i03 = ij I41 I44 x - i10 = ij I42 I41 x - i11 = ij I42 I42 x - i12 = ij I42 I43 x - i13 = ij I42 I44 x - i20 = ij I43 I41 x - i21 = ij I43 I42 x - i22 = ij I43 I43 x - i23 = ij I43 I44 x - i30 = ij I44 I41 x - i31 = ij I44 I42 x - i32 = ij I44 I43 x - i33 = ij I44 I44 x - - s0 = i00 * i11 - i10 * i01 - s1 = i00 * i12 - i10 * i02 - s2 = i00 * i13 - i10 * i03 - s3 = i01 * i12 - i11 * i02 - s4 = i01 * i13 - i11 * i03 - s5 = i02 * i13 - i12 * i03 - c5 = i22 * i33 - i32 * i23 - c4 = i21 * i33 - i31 * i23 - c3 = i21 * i32 - i31 * i22 - c2 = i20 * i33 - i30 * i23 - c1 = i20 * i32 - i30 * i22 - c0 = i20 * i31 - i30 * i21 - - det = s0 * c5 - s1 * c4 + s2 * c3 + s3 * c2 - s4 * c1 + s5 * c0 - - x' = m44 (i11 * c5 - i12 * c4 + i13 * c3) - (-i01 * c5 + i02 * c4 - i03 * c3) - (i31 * s5 - i32 * s4 + i33 * s3) - (-i21 * s5 + i22 * s4 - i23 * s3) - (-i10 * c5 + i12 * c2 - i13 * c1) - (i00 * c5 - i02 * c2 + i03 * c1) - (-i30 * s5 + i32 * s2 - i33 * s1) - (i20 * s5 - i22 * s2 + i23 * s1) - (i10 * c4 - i11 * c2 + i13 * c0) - (-i00 * c4 + i01 * c2 - i03 * c0) - (i30 * s4 - i31 * s2 + i33 * s0) - (-i20 * s4 + i21 * s2 - i23 * s0) - (-i10 * c3 + i11 * c1 - i12 * c0) - (i00 * c3 - i01 * c1 + i02 * c0) - (-i30 * s3 + i31 * s1 - i32 * s0) - (i20 * s3 - i21 * s1 + i22 * s0) -{-# INLINE inv4 #-} - --- | Construct a 2x2 matrix. --- --- Arguments are in row-major order. --- --- >>> m22 1 2 3 4 :: M22 Int --- V2 (V2 1 2) (V2 3 4) --- --- @ 'm22' :: a -> a -> a -> a -> 'M22' a @ --- -m22 :: Basis I2 f => Basis I2 g => a -> a -> a -> a -> f (g a) -m22 a b c d = fillI2 (fillI2 a b) (fillI2 c d) -{-# INLINE m22 #-} - --- | Construct a 2x3 matrix. --- --- Arguments are in row-major order. --- --- @ 'm23' :: a -> a -> a -> a -> a -> a -> 'M23' a @ --- -m23 :: Basis I2 f => Basis I3 g => a -> a -> a -> a -> a -> a -> f (g a) -m23 a b c d e f = fillI2 (fillI3 a b c) (fillI3 d e f) -{-# INLINE m23 #-} - --- | Construct a 2x4 matrix. --- --- Arguments are in row-major order. --- -m24 :: Basis I2 f => Basis I4 g => a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m24 a b c d e f g h = fillI2 (fillI4 a b c d) (fillI4 e f g h) -{-# INLINE m24 #-} - --- | Construct a 3x2 matrix. --- --- Arguments are in row-major order. --- -m32 :: Basis I3 f => Basis I2 g => a -> a -> a -> a -> a -> a -> f (g a) -m32 a b c d e f = fillI3 (fillI2 a b) (fillI2 c d) (fillI2 e f) -{-# INLINE m32 #-} - --- | Construct a 3x3 matrix. --- --- Arguments are in row-major order. --- -m33 :: Basis I3 f => Basis I3 g => a -> a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m33 a b c d e f g h i = fillI3 (fillI3 a b c) (fillI3 d e f) (fillI3 g h i) -{-# INLINE m33 #-} - --- | Construct a 3x4 matrix. --- --- Arguments are in row-major order. --- -m34 :: Basis I3 f => Basis I4 g => a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m34 a b c d e f g h i j k l = fillI3 (fillI4 a b c d) (fillI4 e f g h) (fillI4 i j k l) -{-# INLINE m34 #-} - --- | Construct a 4x2 matrix. --- --- Arguments are in row-major order. --- -m42 :: Basis I4 f => Basis I2 g => a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m42 a b c d e f g h = fillI4 (fillI2 a b) (fillI2 c d) (fillI2 e f) (fillI2 g h) -{-# INLINE m42 #-} - --- | Construct a 4x3 matrix. --- --- Arguments are in row-major order. --- -m43 :: Basis I4 f => Basis I3 g => a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m43 a b c d e f g h i j k l = fillI4 (fillI3 a b c) (fillI3 d e f) (fillI3 g h i) (fillI3 j k l) -{-# INLINE m43 #-} - --- | Construct a 4x4 matrix. --- --- Arguments are in row-major order. --- -m44 :: Basis I4 f => Basis I4 g => a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> a -> f (g a) -m44 a b c d e f g h i j k l m n o p = fillI4 (fillI4 a b c d) (fillI4 e f g h) (fillI4 i j k l) (fillI4 m n o p) -{-# INLINE m44 #-} diff --git a/src/Data/Semimodule/Transform.hs b/src/Data/Semimodule/Transform.hs --- a/src/Data/Semimodule/Transform.hs +++ b/src/Data/Semimodule/Transform.hs @@ -1,6 +1,5 @@ {-# LANGUAGE CPP #-} {-# LANGUAGE Safe #-} -{-# LANGUAGE PolyKinds #-} {-# LANGUAGE ConstraintKinds #-} {-# LANGUAGE DefaultSignatures #-} {-# LANGUAGE DeriveFunctor #-} @@ -14,7 +13,42 @@ {-# LANGUAGE RankNTypes #-} -module Data.Semimodule.Transform where +module Data.Semimodule.Transform ( + -- * Types + type (**) + , type (++) + , type Dim + , type Endo + , Tran(..) + , app + , arr + , invmap + -- * Matrix combinators + , rows + , cols + , projl + , projr + , compl + , compr + , complr + , transpose + -- * Dimensional combinators + , braid + , sbraid + , first + , second + , left + , right + , (***) + , (+++) + , (&&&) + , (|||) + , ($$$) + , adivide + , adivide' + , aselect + , aselect' +) where import safe Control.Category (Category, (>>>)) import safe Data.Functor.Compose @@ -28,41 +62,28 @@ import safe qualified Control.Category as C import safe qualified Data.Bifunctor as B -{- -app' = app @I3 @V3 @I3 @V3 @Int - --- >>> app' foo $ V3 1 2 3 --- V3 2 1 3 --- >>> app' foo >>> app' foo $ V3 1 2 3 --- V3 1 2 3 --- >>> app' (foo >>> foo) $ V3 1 2 3 --- V3 1 2 3 --- -foo = Tran $ \f -> f . t - where - t I31 = I32 - t I32 = I31 - t I33 = I33 --} - +infixr 2 ** +infixr 1 ++ ---------------------------------------------------------------------- +type (f ** g) = Compose f g +type (f ++ g) = Product f g --- | A binary relation between two basis indices. +-- | A dimensional (binary) relation between two bases. -- --- @ 'Index' b c @ relations correspond to (compositions of) +-- @ 'Dim' b c @ relations correspond to (compositions of) -- permutation, projection, and embedding transformations. -- -- See also < https://en.wikipedia.org/wiki/Logical_matrix >. -- -type Index b c = forall a . Tran a b c +type Dim b c = forall a . Tran a b c --- | A general linear transformation between free semimodules indexed with bases /b/ and /c/. +-- | An endomorphism over a free semimodule. -- -newtype Tran a b c = Tran { runTran :: (c -> a) -> (b -> a) } deriving Functor +type Endo a b = Tran a b b -app :: Basis b f => Basis c g => Tran a b c -> g a -> f a -app t = tabulate . runTran t . index +-- | A morphism between free semimodules indexed with bases /b/ and /c/. +-- +newtype Tran a b c = Tran { runTran :: (c -> a) -> (b -> a) } deriving Functor instance Category (Tran a) where id = Tran id @@ -72,6 +93,20 @@ lmap f (Tran t) = Tran $ \ca -> t ca . f rmap = fmap +--------------------------------------------------------------------- + +-- | Apply a transformation to a vector. +-- +app :: Basis2 b c f g => Tran a b c -> g a -> f a +app t = tabulate . runTran t . index + +-- | Lift a function on basis indices into a transformation. +-- +-- @ 'arr' f = 'rmap' f 'C.id' @ +-- +arr :: (b -> c) -> Tran a b c +arr f = Tran (. f) + -- | /Tran a b c/ is an invariant functor on /a/. -- -- See also < http://comonad.com/reader/2008/rotten-bananas/ >. @@ -81,17 +116,13 @@ --------------------------------------------------------------------- --- | An endomorphism over a free semimodule. --- -type Endo a b = Tran a b b - -- | Obtain a matrix by stacking rows. -- -- >>> rows (V2 1 2) :: M22 Int -- V2 (V2 1 2) (V2 1 2) -- -rows :: Free f => Free g => g a -> f (g a) -rows = getCompose . app in1 +rows :: Basis2 b c f g => g a -> (f**g) a +rows = app $ arr snd {-# INLINE rows #-} -- | Obtain a matrix by stacking columns. @@ -99,42 +130,48 @@ -- >>> cols (V2 1 2) :: M22 Int -- V2 (V2 1 1) (V2 2 2) -- -cols :: Free f => Free g => f a -> f (g a) -cols = getCompose . app in2 +cols :: Basis2 b c f g => f a -> (f**g) a +cols = app $ arr fst {-# INLINE cols #-} -projl :: Free f => Free g => Product f g a -> f a -projl = app exl +-- | Project onto the left-hand component of a direct sum. +-- +projl :: Basis2 b c f g => (f++g) a -> f a +projl = app $ arr Left +{-# INLINE projl #-} -projr :: Free f => Free g => Product f g a -> g a -projr = app exr +-- | Project onto the right-hand component of a direct sum. +-- +projr :: Basis2 b c f g => (f++g) a -> g a +projr = app $ arr Right +{-# INLINE projr #-} -- | Left (post) composition with a linear transformation. -- -compl :: Basis b f1 => Basis c f2 => Free g => Index b c -> f2 (g a) -> f1 (g a) -compl f = getCompose . app (first f) . Compose +compl :: Basis3 b c d f1 f2 g => Dim b c -> (f2**g) a -> (f1**g) a +compl f = app (first f) -- | Right (pre) composition with a linear transformation. -- -compr :: Basis b g1 => Basis c g2 => Free f => Index b c -> f (g2 a) -> f (g1 a) -compr f = getCompose . app (second f) . Compose +compr :: Basis3 b c d f g1 g2 => Dim c d -> (f**g2) a -> (f**g1) a +compr f = app (second f) -- | Left and right composition with a linear transformation. -- --- @ 'complr f g' = 'compl f' . 'compr g' @ +-- @ 'complr' f g = 'compl' f '>>>' 'compr' g @ -- -- When /f . g = id/ this induces a similarity transformation: -- --- >>> perm1 = arr (+ I32) --- >>> perm2 = arr (+ I33) +-- >>> perm1 = arr (+ E32) +-- >>> perm2 = arr (+ E33) -- >>> m = m33 1 2 3 4 5 6 7 8 9 :: M33 Int --- >>> conjugate perm1 perm2 m :: M33 Int +-- >>> complr perm1 perm2 m :: M33 Int -- V3 (V3 5 6 4) (V3 8 9 7) (V3 2 3 1) -- -- See also < https://en.wikipedia.org/wiki/Matrix_similarity > & < https://en.wikipedia.org/wiki/Conjugacy_class >. -- -complr :: Basis b1 f1 => Basis c1 f2 => Basis b2 g1 => Basis c2 g2 => Index b1 c1 -> Index b2 c2 -> f2 (g2 a) -> f1 (g1 a) -complr f g = getCompose . app (f *** g) . Compose +complr :: Basis2 b1 c1 f1 f2 => Basis2 b2 c2 g1 g2 => Dim b1 c1 -> Dim b2 c2 -> (f2**g2) a -> (f1**g1) a +complr f g = app (f *** g) -- | Transpose a matrix. -- @@ -144,104 +181,99 @@ -- >>> transpose $ m23 1 2 3 4 5 6 :: M32 Int -- V3 (V2 1 4) (V2 2 5) (V2 3 6) -- -transpose :: Free f => Free g => f (g a) -> g (f a) -transpose = getCompose . app braid . Compose +transpose :: Basis2 b c f g => (f**g) a -> (g**f) a +transpose = app braid {-# INLINE transpose #-} --------------------------------------------------------------------- --- arr toI3 :: Dim3 e => Index e I3 - --- @ 'arr' f = 'rmap' f 'C.id' @ -arr :: (b -> c) -> Index b c -arr f = Tran (. f) - -in1 :: Index (a , b) b -in1 = arr snd -{-# INLINE in1 #-} - -in2 :: Index (a , b) a -in2 = arr fst -{-# INLINE in2 #-} - -exl :: Index a (a + b) -exl = arr Left -{-# INLINE exl #-} - -exr :: Index b (a + b) -exr = arr Right -{-# INLINE exr #-} - -braid :: Index (a , b) (b , a) +-- | Swap components of a tensor product. +-- +braid :: Dim (a , b) (b , a) braid = arr swap {-# INLINE braid #-} -ebraid :: Index (a + b) (b + a) -ebraid = arr eswap -{-# INLINE ebraid #-} +-- | Swap components of a direct sum. +-- +sbraid :: Dim (a + b) (b + a) +sbraid = arr eswap +{-# INLINE sbraid #-} -first :: Index b c -> Index (b , d) (c , d) +-- | Lift a transform into a transform on tensor products. +-- +first :: Dim b c -> Dim (b , d) (c , d) first (Tran caba) = Tran $ \cda -> cda . B.first (caba id) -second :: Index b c -> Index (d , b) (d , c) +-- | Lift a transform into a transform on tensor products. +-- +second :: Dim b c -> Dim (d , b) (d , c) second (Tran caba) = Tran $ \cda -> cda . B.second (caba id) -left :: Index b c -> Index (b + d) (c + d) +-- | Lift a transform into a transform on direct sums. +-- +left :: Dim b c -> Dim (b + d) (c + d) left (Tran caba) = Tran $ \cda -> cda . B.first (caba id) -right :: Index b c -> Index (d + b) (d + c) +-- | Lift a transform into a transform on direct sums. +-- +right :: Dim b c -> Dim (d + b) (d + c) right (Tran caba) = Tran $ \cda -> cda . B.second (caba id) infixr 3 *** -(***) :: Index a1 b1 -> Index a2 b2 -> Index (a1 , a2) (b1 , b2) +-- | Create a transform on a tensor product of semimodules. +-- +(***) :: Dim a1 b1 -> Dim a2 b2 -> Dim (a1 , a2) (b1 , b2) x *** y = first x >>> arr swap >>> first y >>> arr swap {-# INLINE (***) #-} infixr 2 +++ -(+++) :: Index a1 b1 -> Index a2 b2 -> Index (a1 + a2) (b1 + b2) +-- | Create a transform on a direct sum of semimodules. +-- +(+++) :: Dim a1 b1 -> Dim a2 b2 -> Dim (a1 + a2) (b1 + b2) x +++ y = left x >>> arr eswap >>> left y >>> arr eswap {-# INLINE (+++) #-} infixr 3 &&& -(&&&) :: Index a b1 -> Index a b2 -> Index a (b1 , b2) +(&&&) :: Dim a b1 -> Dim a b2 -> Dim a (b1 , b2) x &&& y = dimap fork id $ x *** y {-# INLINE (&&&) #-} infixr 2 ||| -(|||) :: Index a1 b -> Index a2 b -> Index (a1 + a2) b +(|||) :: Dim a1 b -> Dim a2 b -> Dim (a1 + a2) b x ||| y = dimap id join $ x +++ y {-# INLINE (|||) #-} infixr 0 $$$ -($$$) :: Index a (b -> c) -> Index a b -> Index a c +($$$) :: Dim a (b -> c) -> Dim a b -> Dim a c ($$$) f x = dimap fork apply (f *** x) {-# INLINE ($$$) #-} -adivide :: (a -> (a1 , a2)) -> Index a1 b -> Index a2 b -> Index a b +-- | +-- +-- @ 'adivide' 'fork' = 'C.id' @ +-- +adivide :: (a -> (a1 , a2)) -> Dim a1 b -> Dim a2 b -> Dim a b adivide f x y = dimap f fst $ x *** y {-# INLINE adivide #-} -adivide' :: Index a b -> Index a b -> Index a b -adivide' = adivide fork +adivide' :: Dim a1 b -> Dim a2 b -> Dim (a1 , a2) b +adivide' = adivide id {-# INLINE adivide' #-} -adivided :: Index a1 b -> Index a2 b -> Index (a1 , a2) b -adivided = adivide id -{-# INLINE adivided #-} - -aselect :: ((b1 + b2) -> b) -> Index a b1 -> Index a b2 -> Index a b +-- | +-- +-- @ 'aselect' 'join' = 'C.id' @ +-- +aselect :: ((b1 + b2) -> b) -> Dim a b1 -> Dim a b2 -> Dim a b aselect f x y = dimap Left f $ x +++ y {-# INLINE aselect #-} -aselect' :: Index a b -> Index a b -> Index a b -aselect' = aselect join +aselect' :: Dim a b1 -> Dim a b2 -> Dim a (b1 + b2) +aselect' = aselect id {-# INLINE aselect' #-} -aselected :: Index a b1 -> Index a b2 -> Index a (b1 + b2) -aselected = aselect id -{-# INLINE aselected #-} diff --git a/src/Data/Semimodule/Vector.hs b/src/Data/Semimodule/Vector.hs deleted file mode 100644 --- a/src/Data/Semimodule/Vector.hs +++ /dev/null @@ -1,461 +0,0 @@ -{-# LANGUAGE CPP #-} -{-# LANGUAGE Safe #-} -{-# LANGUAGE PolyKinds #-} -{-# LANGUAGE ConstraintKinds #-} -{-# LANGUAGE DefaultSignatures #-} -{-# LANGUAGE DeriveFunctor #-} -{-# LANGUAGE DeriveGeneric #-} -{-# LANGUAGE FlexibleContexts #-} -{-# LANGUAGE FlexibleInstances #-} -{-# LANGUAGE NoImplicitPrelude #-} -{-# LANGUAGE RebindableSyntax #-} -{-# LANGUAGE TypeOperators #-} -{-# LANGUAGE TypeFamilies #-} -{-# LANGUAGE RankNTypes #-} - -module Data.Semimodule.Vector ( - type Basis - , (*.) - , (.*) - , (.*.) - , (><) - , triple - , lerp - , quadrance - , qd - , dirac - , module Data.Semimodule.Vector -) where - -import safe Control.Applicative -import safe Data.Algebra -import safe Data.Bool -import safe Data.Distributive -import safe Data.Functor.Rep -import safe Data.Semifield -import safe Data.Semigroup.Foldable as Foldable1 -import safe Data.Semimodule -import safe Data.Semiring -import safe Prelude hiding (Num(..), Fractional(..), negate, sum, product) - -------------------------------------------------------------------------------- --- V2 -------------------------------------------------------------------------------- - -data V2 a = V2 !a !a deriving (Eq,Ord,Show) - --- | Vector addition. --- --- >>> V2 1 2 <> V2 3 4 --- V2 4 6 --- -instance (Additive-Semigroup) a => Semigroup (V2 a) where - (<>) = mzipWithRep (+) - --- | Matrix addition. --- --- >>> m23 1 2 3 4 5 6 <> m23 7 8 9 1 2 3 :: M23 Int --- V2 (V3 8 10 12) (V3 5 7 9) --- -instance (Additive-Semigroup) a => Semigroup (Additive (V2 a)) where - (<>) = liftA2 $ mzipWithRep (+) - -instance (Additive-Monoid) a => Monoid (V2 a) where - mempty = pureRep zero - -instance (Additive-Monoid) a => Monoid (Additive (V2 a)) where - mempty = pure $ pureRep zero - --- | Vector subtraction. --- --- >>> V2 1 2 << V2 3 4 --- V2 (-2) (-2) --- -instance (Additive-Group) a => Magma (V2 a) where - (<<) = mzipWithRep (-) - --- | Matrix subtraction. --- --- >>> m23 1 2 3 4 5 6 << m23 7 8 9 1 2 3 :: M23 Int --- V2 (V3 (-6) (-6) (-6)) (V3 3 3 3) --- -instance (Additive-Group) a => Magma (Additive (V2 a)) where - (<<) = liftA2 $ mzipWithRep (-) - -instance (Additive-Group) a => Quasigroup (V2 a) -instance (Additive-Group) a => Quasigroup (Additive (V2 a)) -instance (Additive-Group) a => Loop (V2 a) -instance (Additive-Group) a => Loop (Additive (V2 a)) -instance (Additive-Group) a => Group (V2 a) -instance (Additive-Group) a => Group (Additive (V2 a)) - -instance Semiring a => Semimodule a (V2 a) where - (*.) = multl - {-# INLINE (*.) #-} - -instance Functor V2 where - fmap f (V2 a b) = V2 (f a) (f b) - {-# INLINE fmap #-} - a <$ _ = V2 a a - {-# INLINE (<$) #-} - -instance Applicative V2 where - pure = pureRep - liftA2 = liftR2 - -instance Foldable V2 where - foldMap f (V2 a b) = f a <> f b - {-# INLINE foldMap #-} - null _ = False - length _ = two - -instance Foldable1 V2 where - foldMap1 f (V2 a b) = f a <> f b - {-# INLINE foldMap1 #-} - -instance Distributive V2 where - distribute f = V2 (fmap (\(V2 x _) -> x) f) (fmap (\(V2 _ y) -> y) f) - {-# INLINE distribute #-} - -instance Representable V2 where - type Rep V2 = I2 - tabulate f = V2 (f I21) (f I22) - {-# INLINE tabulate #-} - - index (V2 x _) I21 = x - index (V2 _ y) I22 = y - {-# INLINE index #-} - -------------------------------------------------------------------------------- --- Standard basis on two real dimensions -------------------------------------------------------------------------------- - -data I2 = I21 | I22 deriving (Eq, Ord, Show) - -i2 :: a -> a -> I2 -> a -i2 x _ I21 = x -i2 _ y I22 = y - -fillI2 :: Basis I2 f => a -> a -> f a -fillI2 x y = tabulate $ i2 x y - -instance Semigroup (Additive I2) where - Additive I21 <> x = x - x <> Additive I21 = x - - Additive I22 <> Additive I22 = Additive I21 - -instance Monoid (Additive I2) where - mempty = pure I21 - --- trivial diagonal algebra -instance Semiring r => Algebra r I2 where - multiplyWith f = f' where - fi = f I21 I21 - fj = f I22 I22 - - f' I21 = fi - f' I22 = fj - -instance Semiring r => Composition r I2 where - conjugateWith = id - - normWith f = flip multiplyWith I21 $ \ix1 ix2 -> - flip multiplyWith I22 $ \jx1 jx2 -> - f ix1 * f ix2 + f jx1 * f jx2 - -------------------------------------------------------------------------------- --- V3 -------------------------------------------------------------------------------- - - -data V3 a = V3 !a !a !a deriving (Eq,Ord,Show) - --- | Vector addition. --- --- >>> V3 1 2 3 <> V3 4 5 6 --- V3 5 7 9 --- -instance (Additive-Semigroup) a => Semigroup (V3 a) where - (<>) = mzipWithRep (+) - --- | Matrix addition. --- --- >>> V2 (V3 1 2 3) (V3 4 5 6) <> V2 (V3 7 8 9) (V3 1 2 3) --- V2 (V3 8 10 12) (V3 5 7 9) --- -instance (Additive-Semigroup) a => Semigroup (Additive (V3 a)) where - (<>) = liftA2 $ mzipWithRep (+) - -instance (Additive-Monoid) a => Monoid (V3 a) where - mempty = pureRep zero - -instance (Additive-Monoid) a => Monoid (Additive (V3 a)) where - mempty = pure $ pureRep zero - --- | Vector subtraction. --- --- >>> V3 1 2 3 << V3 4 5 6 --- V3 (-3) (-3) (-3) --- -instance (Additive-Group) a => Magma (V3 a) where - (<<) = mzipWithRep (-) - --- | Matrix subtraction. --- --- >>> V3 (V3 1 2 3) (V3 4 5 6) (V3 7 8 9) << V3 (V3 7 8 9) (V3 7 8 9) (V3 7 8 9) --- V3 (V3 (-6) (-6) (-6)) (V3 (-3) (-3) (-3)) (V3 0 0 0) --- -instance (Additive-Group) a => Magma (Additive (V3 a)) where - (<<) = liftA2 $ mzipWithRep (-) - -instance (Additive-Group) a => Quasigroup (V3 a) -instance (Additive-Group) a => Quasigroup (Additive (V3 a)) -instance (Additive-Group) a => Loop (V3 a) -instance (Additive-Group) a => Loop (Additive (V3 a)) -instance (Additive-Group) a => Group (V3 a) -instance (Additive-Group) a => Group (Additive (V3 a)) - -instance Semiring a => Semimodule a (V3 a) where - (*.) = multl - {-# INLINE (*.) #-} - -instance Functor V3 where - fmap f (V3 a b c) = V3 (f a) (f b) (f c) - {-# INLINE fmap #-} - a <$ _ = V3 a a a - {-# INLINE (<$) #-} - -instance Applicative V3 where - pure = pureRep - liftA2 = liftR2 - -instance Foldable V3 where - foldMap f (V3 a b c) = f a <> f b <> f c - {-# INLINE foldMap #-} - null _ = False - --length _ = 3 - -instance Foldable1 V3 where - foldMap1 f (V3 a b c) = f a <> f b <> f c - {-# INLINE foldMap1 #-} - -instance Distributive V3 where - distribute f = V3 (fmap (\(V3 x _ _) -> x) f) (fmap (\(V3 _ y _) -> y) f) (fmap (\(V3 _ _ z) -> z) f) - {-# INLINE distribute #-} - -instance Representable V3 where - type Rep V3 = I3 - tabulate f = V3 (f I31) (f I32) (f I33) - {-# INLINE tabulate #-} - - index (V3 x _ _) I31 = x - index (V3 _ y _) I32 = y - index (V3 _ _ z) I33 = z - {-# INLINE index #-} - -------------------------------------------------------------------------------- --- Standard basis on three real dimensions -------------------------------------------------------------------------------- - -data I3 = I31 | I32 | I33 deriving (Eq, Ord, Show) - -i3 :: a -> a -> a -> I3 -> a -i3 x _ _ I31 = x -i3 _ y _ I32 = y -i3 _ _ z I33 = z - -fillI3 :: Basis I3 f => a -> a -> a -> f a -fillI3 x y z = tabulate $ i3 x y z - -instance Semigroup (Additive I3) where - Additive I31 <> x = x - x <> Additive I31 = x - - Additive I32 <> Additive I33 = Additive I31 - Additive I33 <> Additive I32 = Additive I31 - - Additive I32 <> Additive I32 = Additive I33 - Additive I33 <> Additive I33 = Additive I32 - -instance Monoid (Additive I3) where - mempty = pure I31 - -instance Ring r => Algebra r I3 where - multiplyWith f = f' where - i31 = f I32 I33 - f I33 I32 - i32 = f I33 I31 - f I31 I33 - i33 = f I31 I32 - f I32 I31 - f' I31 = i31 - f' I32 = i32 - f' I33 = i33 - -instance Ring r => Composition r I3 where - conjugateWith = id - - normWith f = flip multiplyWith' I31 $ \ix1 ix2 -> - flip multiplyWith' I32 $ \jx1 jx2 -> - flip multiplyWith' I33 $ \kx1 kx2 -> - f ix1 * f ix2 + f jx1 * f jx2 + f kx1 * f kx2 - - where - multiplyWith' f1 = f1' where - i31 = f1 I31 I31 - i32 = f1 I32 I32 - i33 = f1 I33 I33 - f1' I31 = i31 - f1' I32 = i32 - f1' I33 = i33 - - -------------------------------------------------------------------------------- --- QuaternionBasis -------------------------------------------------------------------------------- - -type QuaternionBasis = Maybe I3 - -instance Ring r => Algebra r QuaternionBasis where - multiplyWith f = maybe fe f' where - e = Nothing - i = Just I31 - j = Just I32 - k = Just I33 - fe = f e e - (f i i + f j j + f k k) - fi = f e i + f i e + (f j k - f k j) - fj = f e j + f j e + (f k i - f i k) - fk = f e k + f k e + (f i j - f j i) - f' I31 = fi - f' I32 = fj - f' I33 = fk - -instance Ring r => Unital r QuaternionBasis where - unitWith x Nothing = x - unitWith _ _ = zero - -instance Ring r => Composition r QuaternionBasis where - conjugateWith f = maybe fe f' where - fe = f Nothing - f' I31 = negate . f $ Just I31 - f' I32 = negate . f $ Just I32 - f' I33 = negate . f $ Just I33 - - normWith f = flip multiplyWith zero $ \ix1 ix2 -> f ix1 * conjugateWith f ix2 - -instance Field r => Division r QuaternionBasis where - reciprocalWith f i = conjugateWith f i / normWith f -{- -reciprocal'' x = divq unit x - -divq (Quaternion r0 (V3 r1 r2 r3)) (Quaternion q0 (V3 q1 q2 q3)) = - (/denom) <$> Quaternion (r0*q0 + r1*q1 + r2*q2 + r3*q3) imag - where denom = q0*q0 + q1*q1 + q2*q2 + q3*q3 - imag = (V3 (r0*q1 + (negate r1*q0) + (negate r2*q3) + r3*q2) - (r0*q2 + r1*q3 + (negate r2*q0) + (negate r3*q1)) - (r0*q3 + (negate r1*q2) + r2*q1 + (negate r3*q0))) - --} - -------------------------------------------------------------------------------- --- V4 -------------------------------------------------------------------------------- - -data V4 a = V4 !a !a !a !a deriving (Eq,Ord,Show) - --- | Vector addition. --- --- >>> V4 1 2 3 4 <> V4 5 6 7 8 --- V4 6 8 10 12 --- -instance (Additive-Semigroup) a => Semigroup (V4 a) where - (<>) = mzipWithRep (+) - --- | Matrix addition. --- --- >>> m24 1 2 3 4 5 6 7 8 <> m24 1 2 3 4 5 6 7 8 :: M24 Int --- V2 (V4 2 4 6 8) (V4 10 12 14 16) --- -instance (Additive-Semigroup) a => Semigroup (Additive (V4 a)) where - (<>) = liftA2 $ mzipWithRep (+) - -instance (Additive-Monoid) a => Monoid (V4 a) where - mempty = pureRep zero - -instance (Additive-Monoid) a => Monoid (Additive (V4 a)) where - mempty = pure $ pureRep zero - --- | Vector subtraction. --- --- >>> V4 1 2 3 << V4 4 5 6 --- V4 (-3) (-3) (-3) --- -instance (Additive-Group) a => Magma (V4 a) where - (<<) = mzipWithRep (-) - --- | Matrix subtraction. --- --- >>> V4 (V4 1 2 3) (V4 4 5 6) (V4 7 8 9) << V4 (V4 7 8 9) (V4 7 8 9) (V4 7 8 9) --- V4 (V4 (-6) (-6) (-6)) (V4 (-3) (-3) (-3)) (V4 0 0 0) --- -instance (Additive-Group) a => Magma (Additive (V4 a)) where - (<<) = liftA2 $ mzipWithRep (-) - -instance (Additive-Group) a => Quasigroup (V4 a) -instance (Additive-Group) a => Quasigroup (Additive (V4 a)) -instance (Additive-Group) a => Loop (V4 a) -instance (Additive-Group) a => Loop (Additive (V4 a)) -instance (Additive-Group) a => Group (V4 a) -instance (Additive-Group) a => Group (Additive (V4 a)) - -instance Semiring a => Semimodule a (V4 a) where - (*.) = multl - {-# INLINE (*.) #-} - -instance Functor V4 where - fmap f (V4 a b c d) = V4 (f a) (f b) (f c) (f d) - {-# INLINE fmap #-} - a <$ _ = V4 a a a a - {-# INLINE (<$) #-} - -instance Applicative V4 where - pure = pureRep - liftA2 = liftR2 - -instance Foldable V4 where - foldMap f (V4 a b c d) = f a <> f b <> f c <> f d - {-# INLINE foldMap #-} - null _ = False - length _ = two + two - -instance Foldable1 V4 where - foldMap1 f (V4 a b c d) = f a <> f b <> f c <> f d - {-# INLINE foldMap1 #-} - -instance Distributive V4 where - distribute f = V4 (fmap (\(V4 x _ _ _) -> x) f) (fmap (\(V4 _ y _ _) -> y) f) (fmap (\(V4 _ _ z _) -> z) f) (fmap (\(V4 _ _ _ w) -> w) f) - {-# INLINE distribute #-} - -instance Representable V4 where - type Rep V4 = I4 - tabulate f = V4 (f I41) (f I42) (f I43) (f I44) - {-# INLINE tabulate #-} - - index (V4 x _ _ _) I41 = x - index (V4 _ y _ _) I42 = y - index (V4 _ _ z _) I43 = z - index (V4 _ _ _ w) I44 = w - {-# INLINE index #-} - -------------------------------------------------------------------------------- --- Standard basis on four real dimensions -------------------------------------------------------------------------------- - -data I4 = I41 | I42 | I43 | I44 deriving (Eq, Ord, Show) - -i4 :: a -> a -> a -> a -> I4 -> a -i4 x _ _ _ I41 = x -i4 _ y _ _ I42 = y -i4 _ _ z _ I43 = z -i4 _ _ _ w I44 = w - -fillI4 :: Basis I4 f => a -> a -> a -> a -> f a -fillI4 x y z w = tabulate $ i4 x y z w diff --git a/src/Data/Semiring.hs b/src/Data/Semiring.hs --- a/src/Data/Semiring.hs +++ b/src/Data/Semiring.hs @@ -11,23 +11,35 @@ {-# LANGUAGE MonoLocalBinds #-} module Data.Semiring ( + -- * Types type (-) - , zero, one, two, (+), (*), (-), (^) - , sum, sum1, sumWith, sumWith1 - , product, product1, productWith, productWith1 - , cross, cross1 - , eval, evalWith, eval1, evalWith1 - , negate, abs, signum + -- * Presemirings , type PresemiringLaw, Presemiring + , (+), (*) + , sum1, sumWith1 + , product1, productWith1 + , xmult1 + , eval1, evalWith1 + -- * Semirings , type SemiringLaw, Semiring + , zero, one, two + , (^) + , sum, sumWith + , product, productWith + , xmult + , eval, evalWith + -- * Rings , type RingLaw, Ring + , (-) + , negate, abs, signum + -- * Re-exports + , mreplicate , Additive(..) , Multiplicative(..) , Magma(..) - , Quasigroup(..) - , Loop(..) + , Quasigroup + , Loop , Group(..) - , mreplicate ) where import safe Control.Applicative @@ -68,7 +80,7 @@ -- /Distributivity/ -- -- @ --- (a '+' b) '*' c '==' (a '*' c) '+' (b '*' c) +-- (a '+' b) '*' c = (a '*' c) '+' (b '*' c) -- @ -- -- Note that addition and multiplication needn't be commutative. @@ -79,6 +91,71 @@ class PresemiringLaw a => Presemiring a +-- | Evaluate a non-empty presemiring sum. +-- +sum1 :: Presemiring a => Foldable1 f => f a -> a +sum1 = sumWith1 id + +-- | Evaluate a non-empty presemiring sum using a given presemiring. +-- +-- >>> evalWith1 Max $ (1 :| [2..5 :: Int]) :| [1 :| [2..5 :: Int]] +-- | Fold over a non-empty collection using the additive operation of an arbitrary semiring. +-- +-- >>> sumWith1 First $ (1 :| [2..5 :: Int]) * (1 :| [2..5 :: Int]) +-- First {getFirst = 1} +-- >>> sumWith1 First $ Nothing :| [Just (5 :: Int), Just 6, Nothing] +-- First {getFirst = Nothing} +-- >>> sumWith1 Just $ 1 :| [2..5 :: Int] +-- Just 15 +-- +sumWith1 :: Foldable1 t => Presemiring a => (b -> a) -> t b -> a +sumWith1 f = unAdditive . foldMap1 (Additive . f) +{-# INLINE sumWith1 #-} + +-- | Evaluate a non-empty presemiring product. +-- +product1 :: Presemiring a => Foldable1 f => f a -> a +product1 = productWith1 id + +-- | Evaluate a non-empty presemiring product using a given presemiring. +-- +-- As the collection is non-empty this does not require a distinct multiplicative unit: +-- +-- >>> productWith1 Just $ 1 :| [2..5 :: Int] +-- Just 120 +-- >>> productWith1 First $ 1 :| [2..(5 :: Int)] +-- First {getFirst = 15} +-- >>> productWith1 First $ Nothing :| [Just (5 :: Int), Just 6, Nothing] +-- First {getFirst = Just 11} +-- +productWith1 :: Foldable1 t => Presemiring a => (b -> a) -> t b -> a +productWith1 f = unMultiplicative . foldMap1 (Multiplicative . f) +{-# INLINE productWith1 #-} + +-- | Cross-multiply two non-empty collections. +-- +-- >>> xmult1 (Right 2 :| [Left "oops"]) (Right 2 :| [Right 3]) :: Either [Char] Int +-- Right 4 +-- +xmult1 :: Foldable1 f => Apply f => Presemiring a => f a -> f a -> a +xmult1 a b = sum1 $ liftF2 (*) a b +{-# INLINE xmult1 #-} + +-- | Evaluate a presemiring expression. +-- +eval1 :: Presemiring a => Functor f => Foldable1 f => Foldable1 g => f (g a) -> a +eval1 = sum1 . fmap product1 + +-- | Evaluate a presemiring expression using a given presemiring. +-- +-- >>> evalWith1 (Max . Down) $ (1 :| [2..5 :: Int]) :| [-5 :| [2..5 :: Int]] +-- Max {getMax = Down 9} +-- >>> evalWith1 Max $ (1 :| [2..5 :: Int]) :| [-5 :| [2..5 :: Int]] +-- Max {getMax = 15} +-- +evalWith1 :: Presemiring r => Functor f => Functor g => Foldable1 f => Foldable1 g => (a -> r) -> f (g a) -> r +evalWith1 f = sum1 . fmap product1 . (fmap . fmap) f + ------------------------------------------------------------------------------- -- Semiring ------------------------------------------------------------------------------- @@ -94,76 +171,68 @@ -- /Neutrality/ -- -- @ --- 'zero' '+' r '==' r --- 'one' '*' r '==' r +-- 'zero' '+' r = r +-- 'one' '*' r = r -- @ -- -- /Absorbtion/ -- -- @ --- 'zero' '*' a '==' 'zero' +-- 'zero' '*' a = 'zero' -- @ -- class (Presemiring a, SemiringLaw a) => Semiring a -two :: (Additive-Semigroup) a => (Multiplicative-Monoid) a => a +-- | +-- +-- @ +-- 'two' = 'one' '+' 'one' +-- @ +-- +two :: Semiring a => a two = one + one {-# INLINE two #-} - infixr 8 ^ --- @ 'one' == a '^' 0 @ +-- | Evaluate a natural-numbered power of a semiring element. -- +-- @ 'one' = a '^' 0 @ +-- -- >>> 8 ^ 0 :: Int -- 1 -- (^) :: Semiring a => a -> Natural -> a a ^ n = unMultiplicative $ mreplicate (P.fromIntegral n) (Multiplicative a) +-- | Evaluate a semiring sum. +-- -- >>> sum [1..5 :: Int] -- 15 +-- sum :: (Additive-Monoid) a => Presemiring a => Foldable f => f a -> a sum = sumWith id -sum1 :: Presemiring a => Foldable1 f => f a -> a -sum1 = sumWith1 id - +-- | Evaluate a semiring sum using a given semiring. +-- sumWith :: (Additive-Monoid) a => Presemiring a => Foldable t => (b -> a) -> t b -> a sumWith f = foldr' ((+) . f) zero {-# INLINE sumWith #-} --- >>> evalWith1 Max $ (1 :| [2..5 :: Int]) :| [1 :| [2..5 :: Int]] --- | Fold over a non-empty collection using the additive operation of an arbitrary semiring. --- --- >>> sumWith1 First $ (1 :| [2..5 :: Int]) * (1 :| [2..5 :: Int]) --- First {getFirst = 1} --- >>> sumWith1 First $ Nothing :| [Just (5 :: Int), Just 6, Nothing] --- First {getFirst = Nothing} --- >>> sumWith1 Just $ 1 :| [2..5 :: Int] --- Just 15 +-- | Evaluate a semiring product. -- -sumWith1 :: Foldable1 t => Presemiring a => (b -> a) -> t b -> a -sumWith1 f = unAdditive . foldMap1 (Additive . f) -{-# INLINE sumWith1 #-} - -- >>> product [1..5 :: Int] -- 120 +-- product :: (Multiplicative-Monoid) a => Presemiring a => Foldable f => f a -> a product = productWith id +-- | Evaluate a semiring product using a given semiring. -- --- | The product of at a list of semiring elements (of length at least one) -product1 :: Presemiring a => Foldable1 f => f a -> a -product1 = productWith1 id - --- | Fold over a collection using the multiplicative operation of an arbitrary semiring. --- -- @ --- 'product' f '==' 'Data.foldr'' ((*) . f) 'one' +-- 'product' f = 'Data.foldr'' (('*') . f) 'one' -- @ -- --- -- >>> productWith Just [1..5 :: Int] -- Just 120 -- @@ -171,43 +240,18 @@ productWith f = foldr' ((*) . f) one {-# INLINE productWith #-} - --- | Fold over a non-empty collection using the multiplicative operation of a semiring. --- --- As the collection is non-empty this does not require a distinct multiplicative unit: --- --- >>> productWith1 Just $ 1 :| [2..5 :: Int] --- Just 120 --- >>> productWith1 First $ 1 :| [2..(5 :: Int)] --- First {getFirst = 15} --- >>> productWith1 First $ Nothing :| [Just (5 :: Int), Just 6, Nothing] --- First {getFirst = Just 11} --- -productWith1 :: Foldable1 t => Presemiring a => (b -> a) -> t b -> a -productWith1 f = unMultiplicative . foldMap1 (Multiplicative . f) -{-# INLINE productWith1 #-} - -- | Cross-multiply two collections. -- --- >>> cross (V3 1 2 3) (V3 1 2 3) +-- >>> xmult (V3 1 2 3) (V3 1 2 3) -- 14 --- >>> cross [1,2,3 :: Int] [1,2,3] +-- >>> xmult [1,2,3 :: Int] [1,2,3] -- 36 --- >>> cross [1,2,3 :: Int] [] +-- >>> xmult [1,2,3 :: Int] [] -- 0 -- -cross :: Foldable f => Applicative f => Presemiring a => (Additive-Monoid) a => f a -> f a -> a -cross a b = sum $ liftA2 (*) a b -{-# INLINE cross #-} - --- | Cross-multiply two non-empty collections. --- --- >>> cross1 (Right 2 :| [Left "oops"]) (Right 2 :| [Right 3]) :: Either [Char] Int --- Right 4 --- -cross1 :: Foldable1 f => Apply f => Presemiring a => f a -> f a -> a -cross1 a b = sum1 $ liftF2 (*) a b -{-# INLINE cross1 #-} +xmult :: Foldable f => Applicative f => Presemiring a => (Additive-Monoid) a => f a -> f a -> a +xmult a b = sum $ liftA2 (*) a b +{-# INLINE xmult #-} -- | Evaluate a semiring expression. -- @@ -221,22 +265,14 @@ eval :: Semiring a => Functor f => Foldable f => Foldable g => f (g a) -> a eval = sum . fmap product +-- | Evaluate a semiring expression using a given semiring. +-- -- >>> evalWith Max [[1..4 :: Int], [0..2 :: Int]] -- Max {getMax = 24} +-- evalWith :: Semiring r => Functor f => Functor g => Foldable f => Foldable g => (a -> r) -> f (g a) -> r evalWith f = sum . fmap product . (fmap . fmap) f -eval1 :: Presemiring a => Functor f => Foldable1 f => Foldable1 g => f (g a) -> a -eval1 = sum1 . fmap product1 - --- >>> evalWith1 (Max . Down) $ (1 :| [2..5 :: Int]) :| [-5 :| [2..5 :: Int]] --- Max {getMax = Down 9} --- >>> evalWith1 Max $ (1 :| [2..5 :: Int]) :| [-5 :| [2..5 :: Int]] --- Max {getMax = 15} --- -evalWith1 :: Presemiring r => Functor f => Functor g => Foldable1 f => Foldable1 g => (a -> r) -> f (g a) -> r -evalWith1 f = sum1 . fmap product1 . (fmap . fmap) f - ------------------------------------------------------------------------------- -- Ring ------------------------------------------------------------------------------- @@ -249,9 +285,9 @@ -- -- The basic properties of a ring follow immediately from the axioms: -- --- @ r '*' 'zero' '==' 'zero' '==' 'zero' '*' r @ +-- @ r '*' 'zero' = 'zero' = 'zero' '*' r @ -- --- @ 'negate' 'one' '*' r '==' 'negate' r @ +-- @ 'negate' 'one' '*' r = 'negate' r @ -- -- Furthermore, the binomial formula holds for any commuting pair of elements (that is, any /a/ and /b/ such that /a * b = b * a/). -- @@ -270,177 +306,12 @@ -- class (Semiring a, RingLaw a) => Ring a where -negate :: (Additive-Group) a => a -> a -negate a = zero - a -{-# INLINE negate #-} - --- | Absolute value of an element. --- --- @ 'abs' r '==' 'mul' r ('signum' r) @ --- --- https://en.wikipedia.org/wiki/Linearly_ordered_group -abs :: (Additive-Group) a => Ord a => a -> a -abs x = bool (negate x) x $ zero <= x -{-# INLINE abs #-} - -- satisfies trichotomy law: -- Exactly one of the following is true: a is positive, -a is positive, or a = 0. -- This property follows from the fact that ordered rings are abelian, linearly ordered groups with respect to addition. -signum :: RingLaw a => Ord a => a -> a +signum :: Ring a => Ord a => a -> a signum x = bool (negate one) one $ zero <= x {-# INLINE signum #-} - -{- --- | Default implementation of 'fromBoolean' given a multiplicative unit. --- -fromBooleanDef :: Unital a => a -> Bool -> a -fromBooleanDef _ False = mempty -fromBooleanDef o True = o -{-# INLINE fromBooleanDef #-} - --- | Multiplicative unit. --- --- Note that 'one' needn't be distinct from 'mempty' for a semiring to be valid. --- -one :: Unital a => a -one = fromBoolean True -{-# INLINE one #-} - - -infixr 8 ^ - -(^) :: Unital a => a -> Natural -> a -(^) = flip sinnum' -{-# INLINE (^) #-} - --- | A generalization of 'Data.List.replicate' to an arbitrary 'Monoid'. --- --- Adapted from . --- -sinnum :: Monoid a => Natural -> a -> a -sinnum n a - | n == 0 = mempty - | otherwise = f a n - where - f x y - | even y = f (x <> x) (y `quot` 2) - | y == 1 = x - | otherwise = g (x <> x) ((y N.- 1) `quot` 2) x - g x y z - | even y = g (x <> x) (y `quot` 2) z - | y == 1 = x <> z - | otherwise = g (x <> x) ((y N.- 1) `quot` 2) (x <> z) -{-# INLINE sinnum #-} - -sinnum' :: Unital a => Natural -> a -> a -sinnum' n a = getProd $ sinnum n (Prod a) -{-# INLINE sinnum' #-} - -powers :: Unital a => Natural -> a -> a -powers n a = foldr' (<>) one . flip unfoldr n $ \m -> - if m == 0 then Nothing else Just (a^m,m N.- 1) -{-# INLINE powers #-} - -------------------------------------------------------------------------------- --- Pre-semirings -------------------------------------------------------------------------------- - -instance Semigroup a => Semiring (Either e a) where - (*) = liftA2 (<>) - {-# INLINE (*) #-} - - -instance Semiring Ordering where - LT * LT = LT - LT * GT = LT - _ * EQ = EQ - EQ * _ = EQ - GT * x = x - - fromBoolean = fromBooleanDef GT - - - - fromBoolean = const . fromBoolean - -instance Unital a => Semiring (Op a b) where - Op f * Op g = Op $ \x -> f x * g x - {-# INLINE (*) #-} - - fromBoolean = fromBooleanDef $ Op (const one) - -instance (Unital a, Unital b) => Semiring (a, b) where - (a, b) * (c, d) = (a*c, b*d) - {-# INLINE (*) #-} - - fromBoolean = liftA2 (,) fromBoolean fromBoolean - -instance (Unital a, Unital b, Unital c) => Semiring (a, b, c) where - (a, b, c) * (d, e, f) = (a*d, b*e, c*f) - {-# INLINE (*) #-} - - fromBoolean = liftA3 (,,) fromBoolean fromBoolean fromBoolean - - - - -{- ---------------------------------------------------------------------- --- Instances (contravariant) ---------------------------------------------------------------------- - --- Note that due to the underlying 'Monoid' instance this instance --- has 'All' semiring semantics rather than 'Any'. -instance Semiring (Predicate a) where - Predicate f * Predicate g = Predicate $ \x -> f x || g x - {-# INLINE (*) #-} - - --Note that the truth values are flipped here to create a - --valid semiring homomorphism. Users should precompose with 'not' - --where necessary. - fromBoolean False = Predicate $ const True - fromBoolean True = Predicate $ const False - - --- Note that due to the underlying 'Monoid' instance this instance --- has 'All' semiring semantics rather than 'Any'. -instance Semiring (Equivalence a) where - Equivalence f * Equivalence g = Equivalence $ \x y -> f x y || g x y - {-# INLINE (*) #-} - - --Note that the truth values are flipped here to create a - --valid semiring homomorphism. Users should precompose with 'not' - --where necessary. - fromBoolean False = Equivalence $ \_ _ -> True - fromBoolean True = Equivalence $ \_ _ -> False --} - ---------------------------------------------------------------------- --- Instances (containers) ---------------------------------------------------------------------- - -instance Ord a => Semiring (Set.Set a) where - (*) = Set.intersection - -instance Monoid a => Semiring (Seq.Seq a) where - (*) = liftA2 (<>) - {-# INLINE (*) #-} - - fromBoolean = fromBooleanDef $ Seq.singleton mempty - -instance (Ord k, Monoid k, Monoid a) => Semiring (Map.Map k a) where - xs * ys = foldMap (flip Map.map xs . (<>)) ys - {-# INLINE (*) #-} - - fromBoolean = fromBooleanDef $ Map.singleton mempty mempty - -instance Monoid a => Semiring (IntMap.IntMap a) where - xs * ys = foldMap (flip IntMap.map xs . (<>)) ys - {-# INLINE (*) #-} - - fromBoolean = fromBooleanDef $ IntMap.singleton 0 mempty - --} --------------------------------------------------------------------- -- Instances diff --git a/src/Data/Semiring/Property.hs b/src/Data/Semiring/Property.hs --- a/src/Data/Semiring/Property.hs +++ b/src/Data/Semiring/Property.hs @@ -17,8 +17,8 @@ , annihilative_multiplication_on , distributive_finite_on -- * Left-distributive presemirings and semirings - , distributive_cross_on - , distributive_cross1_on + , distributive_xmult_on + , distributive_xmult1_on -- * Commutative presemirings & semirings , commutative_multiplication_on -- * Cancellative presemirings & semirings @@ -32,10 +32,7 @@ import Data.Foldable (Foldable) import Data.Functor.Apply (Apply) import Data.Semigroup.Foldable (Foldable1) -import Data.Semigroup.Additive -import Data.Semigroup.Multiplicative import Data.Semigroup.Property -import qualified Test.Function as Prop import qualified Test.Operation as Prop import Prelude hiding (Num(..), sum) @@ -71,32 +68,7 @@ morphism_distribitive_on (==) f x y z ------------------------------------------------------------------------------------ --- Required properties of semigroups - --- | \( \forall a, b, c \in R: (a + b) + c \sim a + (b + c) \) --- --- All semigroups must right-associate addition. --- --- This is a required property. --- -associative_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> r -> r -> b -associative_addition_on (~~) = Prop.associative_on (~~) (+) - --- | \( \forall a, b, c \in R: (a * b) * c \sim a * (b * c) \) --- --- All semigroups must right-associate multiplication. --- --- This is a required property. --- -associative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> r -> r -> b -associative_multiplication_on (~~) = Prop.associative_on (~~) (*) - --- | \( \forall a, b \in R: a + b \sim b + a \) --- --- This is a an /optional/ property for semigroups, and a /required/ property for semirings. --- -commutative_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> r -> b -commutative_addition_on (~~) = Prop.commutative_on (~~) (+) +-- Required properties of semirings -- | \( \forall a, b, c \in R: (a + b) * c \sim (a * c) + (b * c) \) -- @@ -135,18 +107,6 @@ ------------------------------------------------------------------------------------ -- Required properties of semirings -morphism_additive_on :: (Additive-Semigroup) r => (Additive-Semigroup) s => Rel s b -> (r -> s) -> r -> r -> b -morphism_additive_on (~~) f x y = (f $ x + y) ~~ (f x + f y) - -morphism_multiplicative_on :: (Multiplicative-Semigroup) r => (Multiplicative-Semigroup) s => Rel s b -> (r -> s) -> r -> r -> b -morphism_multiplicative_on (~~) f x y = (f $ x * y) ~~ (f x * f y) - -morphism_additive_on' :: (Additive-Monoid) r => (Additive-Monoid) s => Rel s b -> (r -> s) -> b -morphism_additive_on' (~~) f = (f zero) ~~ zero - -morphism_multiplicative_on' :: (Multiplicative-Monoid) r => (Multiplicative-Monoid) s => Rel s b -> (r -> s) -> b -morphism_multiplicative_on' (~~) f = (f one) ~~ one - -- | Semiring morphisms are monoidal presemiring morphisms. -- -- This is a required property for semiring morphisms. @@ -157,45 +117,6 @@ morphism_additive_on' (==) f && morphism_multiplicative_on' (==) f - --- | \( \forall a \in R: (z + a) \sim a \) --- --- A semigroup with a right-neutral additive identity must satisfy: --- --- @ --- 'neutral_addition' 'zero' ~~ const True --- @ --- --- Or, equivalently: --- --- @ --- 'zero' '+' r ~~ r --- @ --- --- This is a required property for additive monoids. --- -neutral_addition_on :: (Additive-Monoid) r => Rel r b -> r -> b -neutral_addition_on (~~) = Prop.neutral_on (~~) (+) zero - --- | \( \forall a \in R: (o * a) \sim a \) --- --- A semigroup with a right-neutral multiplicative identity must satisfy: --- --- @ --- 'neutral_multiplication' 'one' ~~ const True --- @ --- --- Or, equivalently: --- --- @ --- 'one' '*' r ~~ r --- @ --- --- This is a required propert for multiplicative monoids. --- -neutral_multiplication_on :: (Multiplicative-Monoid) r => Rel r b -> r -> b -neutral_multiplication_on (~~) = Prop.neutral_on (~~) (*) one - -- | \( \forall a \in R: (z * a) \sim u \) -- -- A /R/ is semiring then its addititive one must be right-annihilative, i.e.: @@ -237,70 +158,14 @@ -- | \( \forall M,N \geq 0; a_1 \dots a_M, b_1 \dots b_N \in R: (\sum_{i=1}^M a_i) * (\sum_{j=1}^N b_j) \sim \sum_{i=1 j=1}^{i=M j=N} a_i * b_j \) -- --- If /R/ is also left-distributive then it supports cross-multiplication. +-- If /R/ is also left-distributive then it supports xmult-multiplication. -- -distributive_cross_on :: Semiring r => Applicative f => Foldable f => Rel r b -> f r -> f r -> b -distributive_cross_on (~~) as bs = (sum as * sum bs) ~~ (cross as bs) +distributive_xmult_on :: Semiring r => Applicative f => Foldable f => Rel r b -> f r -> f r -> b +distributive_xmult_on (~~) as bs = (sum as * sum bs) ~~ (xmult as bs) -- | \( \forall M,N \geq 1; a_1 \dots a_M, b_1 \dots b_N \in R: (\sum_{i=1}^M a_i) * (\sum_{j=1}^N b_j) = \sum_{i=1 j=1}^{i=M j=N} a_i * b_j \) -- --- If /R/ is also left-distributive then it supports (non-empty) cross-multiplication. --- -distributive_cross1_on :: Presemiring r => Apply f => Foldable1 f => Rel r b -> f r -> f r -> b -distributive_cross1_on (~~) as bs = (sum1 as * sum1 bs) ~~ (cross1 as bs) - ------------------------------------------------------------------------------------- --- Commutative presemirings and semirings - --- | \( \forall a, b \in R: a * b \sim b * a \) --- --- This is a an /optional/ property for semigroups, and a /optional/ property for semirings. --- It is a /required/ property for rings. --- -commutative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> r -> b -commutative_multiplication_on (~~) = Prop.commutative_on (~~) (*) - ------------------------------------------------------------------------------------- --- Properties of cancellative semigroups - --- | \( \forall a, b, c \in R: b + a \sim c + a \Rightarrow b = c \) --- --- If /R/ is right-cancellative wrt addition then for all /a/ --- the section /(a +)/ is injective. --- --- See < https://en.wikipedia.org/wiki/Cancellation_property > --- -cancellative_addition_on :: (Additive-Semigroup) r => Rel r Bool -> r -> r -> r -> Bool -cancellative_addition_on (~~) a = Prop.injective_on (~~) (+ a) - --- | \( \forall a, b, c \in R: b * a \sim c * a \Rightarrow b = c \) --- --- If /R/ is right-cancellative wrt multiplication then for all /a/ --- the section /(a *)/ is injective. --- -cancellative_multiplication_on :: (Multiplicative-Semigroup) r => Rel r Bool -> r -> r -> r -> Bool -cancellative_multiplication_on (~~) a = Prop.injective_on (~~) (* a) - --- | Idempotency property for additive semigroups. --- --- @ 'idempotent_addition' = 'absorbative_addition' 'one' @ --- --- See < https://en.wikipedia.org/wiki/Band_(mathematics) >. --- --- This is a required property for lattices. --- -idempotent_addition_on :: (Additive-Semigroup) r => Rel r b -> r -> b -idempotent_addition_on (~~) r = (r + r) ~~ r - --- | Idempotency property for semigroups. --- --- @ 'idempotent_multiplication' = 'absorbative_multiplication' 'zero' @ --- --- See < https://en.wikipedia.org/wiki/Band_(mathematics) >. --- --- This is a an /optional/ property for semigroups, and a /optional/ property for semirings. --- --- This is a /required/ property for lattices. +-- If /R/ is also left-distributive then it supports (non-empty) xmult-multiplication. -- -idempotent_multiplication_on :: (Multiplicative-Semigroup) r => Rel r b -> r -> b -idempotent_multiplication_on (~~) r = (r * r) ~~ r +distributive_xmult1_on :: Presemiring r => Apply f => Foldable1 f => Rel r b -> f r -> f r -> b +distributive_xmult1_on (~~) as bs = (sum1 as * sum1 bs) ~~ (xmult1 as bs) diff --git a/test/Test/Data/Dioid/Signed.hs b/test/Test/Data/Dioid/Signed.hs --- a/test/Test/Data/Dioid/Signed.hs +++ b/test/Test/Data/Dioid/Signed.hs @@ -146,7 +146,7 @@ x = f 0.37794903 y = f 0.3269925 -cosunit (residl x) y +cosaempty (residl x) y residl x = Conn (<>x) . (//x) $ y