diff --git a/src/Data/AdditiveGroup.hs b/src/Data/AdditiveGroup.hs --- a/src/Data/AdditiveGroup.hs +++ b/src/Data/AdditiveGroup.hs @@ -1,50 +1,74 @@ -{-# LANGUAGE TypeOperators #-} +{-# LANGUAGE TypeOperators, CPP #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE ScopedTypeVariables #-} ---------------------------------------------------------------------- -- | -- Module : Data.AdditiveGroup -- Copyright : (c) Conal Elliott and Andy J Gill 2008 -- License : BSD3 --- +-- -- Maintainer : conal@conal.net, andygill@ku.edu -- Stability : experimental --- +-- -- Groups: zero, addition, and negation (additive inverse) ---------------------------------------------------------------------- module Data.AdditiveGroup - ( - AdditiveGroup(..), (^-^), sumV + ( + AdditiveGroup(..), sumV , Sum(..), inSum, inSum2 ) where import Prelude hiding (foldr) import Control.Applicative +#if !(MIN_VERSION_base(4,8,0)) import Data.Monoid (Monoid(..)) -import Data.Foldable (Foldable,foldr) +import Data.Foldable (Foldable) +#endif +import Data.Foldable (foldr) import Data.Complex hiding (magnitude) import Data.Ratio +#if !(MIN_VERSION_base(4,11,0)) +import Data.Semigroup (Semigroup(..)) +#endif +import Foreign.C.Types (CSChar, CInt, CShort, CLong, CLLong, CIntMax, CFloat, CDouble) import Data.MemoTrie +import Data.VectorSpace.Generic +import qualified GHC.Generics as Gnrx +import GHC.Generics (Generic, (:*:)(..)) + infixl 6 ^+^, ^-^ -- | Additive group @v@. class AdditiveGroup v where -- | The zero element: identity for '(^+^)' zeroV :: v + default zeroV :: (Generic v, AdditiveGroup (VRep v)) => v + zeroV = Gnrx.to (zeroV :: VRep v) + {-# INLINE zeroV #-} -- | Add vectors (^+^) :: v -> v -> v + default (^+^) :: (Generic v, AdditiveGroup (VRep v)) => v -> v -> v + v ^+^ v' = Gnrx.to (Gnrx.from v ^+^ Gnrx.from v' :: VRep v) + {-# INLINE (^+^) #-} -- | Additive inverse negateV :: v -> v - --- | Group subtraction -(^-^) :: AdditiveGroup v => v -> v -> v -v ^-^ v' = v ^+^ negateV v' + default negateV :: (Generic v, AdditiveGroup (VRep v)) => v -> v + negateV v = Gnrx.to (negateV $ Gnrx.from v :: VRep v) + {-# INLINE negateV #-} + -- | Group subtraction + (^-^) :: v -> v -> v + v ^-^ v' = v ^+^ negateV v' -- | Sum over several vectors sumV :: (Foldable f, AdditiveGroup v) => f v -> v sumV = foldr (^+^) zeroV +{-# INLINE sumV #-} instance AdditiveGroup () where zeroV = () @@ -52,16 +76,28 @@ negateV = id -- For 'Num' types: --- +-- -- instance AdditiveGroup n where {zeroV=0; (^+^) = (+); negateV = negate} -instance AdditiveGroup Int where {zeroV=0; (^+^) = (+); negateV = negate} -instance AdditiveGroup Integer where {zeroV=0; (^+^) = (+); negateV = negate} -instance AdditiveGroup Float where {zeroV=0; (^+^) = (+); negateV = negate} -instance AdditiveGroup Double where {zeroV=0; (^+^) = (+); negateV = negate} -instance Integral a => AdditiveGroup (Ratio a) where - {zeroV=0; (^+^) = (+); negateV = negate} +#define ScalarTypeCon(con,t) \ + instance con => AdditiveGroup (t) where {zeroV=0; (^+^) = (+); negateV = negate} +#define ScalarType(t) ScalarTypeCon((),t) + +ScalarType(Int) +ScalarType(Integer) +ScalarType(Float) +ScalarType(Double) +ScalarType(CSChar) +ScalarType(CInt) +ScalarType(CShort) +ScalarType(CLong) +ScalarType(CLLong) +ScalarType(CIntMax) +ScalarType(CFloat) +ScalarType(CDouble) +ScalarTypeCon(Integral a,Ratio a) + instance (RealFloat v, AdditiveGroup v) => AdditiveGroup (Complex v) where zeroV = zeroV :+ zeroV (^+^) = (+) @@ -149,6 +185,7 @@ instance Functor Sum where fmap f (Sum a) = Sum (f a) + {-# INLINE fmap #-} -- instance Applicative Sum where -- pure a = Sum a @@ -156,32 +193,40 @@ instance Applicative Sum where pure = Sum + {-# INLINE pure #-} (<*>) = inSum2 ($) + {-# INLINE (<*>) #-} +instance AdditiveGroup a => Semigroup (Sum a) where + (<>) = liftA2 (^+^) + {-# INLINE (<>) #-} + instance AdditiveGroup a => Monoid (Sum a) where mempty = Sum zeroV - mappend = liftA2 (^+^) - +#if !(MIN_VERSION_base(4,11,0)) + mappend = (<>) +#endif -- | Application a unary function inside a 'Sum' inSum :: (a -> b) -> (Sum a -> Sum b) inSum = getSum ~> Sum +{-# INLINE inSum #-} -- | Application a binary function inside a 'Sum' inSum2 :: (a -> b -> c) -> (Sum a -> Sum b -> Sum c) inSum2 = getSum ~> inSum - +{-# INLINE inSum2 #-} instance AdditiveGroup a => AdditiveGroup (Sum a) where - zeroV = mempty - (^+^) = mappend + zeroV = Sum zeroV + (^+^) = (<>) negateV = inSum negateV - ---- to go elsewhere (~>) :: (a' -> a) -> (b -> b') -> ((a -> b) -> (a' -> b')) (i ~> o) f = o . f . i +{-# INLINE (~>) #-} -- result :: (b -> b') -> ((a -> b) -> (a -> b')) -- result = (.) @@ -190,3 +235,33 @@ -- argument = flip (.) -- g ~> f = result g . argument f + + + +instance AdditiveGroup a => AdditiveGroup (Gnrx.Rec0 a s) where + zeroV = Gnrx.K1 zeroV + {-# INLINE zeroV #-} + negateV (Gnrx.K1 v) = Gnrx.K1 $ negateV v + {-# INLINE negateV #-} + Gnrx.K1 v ^+^ Gnrx.K1 w = Gnrx.K1 $ v ^+^ w + {-# INLINE (^+^) #-} + Gnrx.K1 v ^-^ Gnrx.K1 w = Gnrx.K1 $ v ^-^ w + {-# INLINE (^-^) #-} +instance AdditiveGroup (f p) => AdditiveGroup (Gnrx.M1 i c f p) where + zeroV = Gnrx.M1 zeroV + {-# INLINE zeroV #-} + negateV (Gnrx.M1 v) = Gnrx.M1 $ negateV v + {-# INLINE negateV #-} + Gnrx.M1 v ^+^ Gnrx.M1 w = Gnrx.M1 $ v ^+^ w + {-# INLINE (^+^) #-} + Gnrx.M1 v ^-^ Gnrx.M1 w = Gnrx.M1 $ v ^-^ w + {-# INLINE (^-^) #-} +instance (AdditiveGroup (f p), AdditiveGroup (g p)) => AdditiveGroup ((f :*: g) p) where + zeroV = zeroV :*: zeroV + {-# INLINE zeroV #-} + negateV (x:*:y) = negateV x :*: negateV y + {-# INLINE negateV #-} + (x:*:y) ^+^ (ξ:*:υ) = (x^+^ξ) :*: (y^+^υ) + {-# INLINE (^+^) #-} + (x:*:y) ^-^ (ξ:*:υ) = (x^-^ξ) :*: (y^-^υ) + {-# INLINE (^-^) #-} diff --git a/src/Data/AffineSpace.hs b/src/Data/AffineSpace.hs --- a/src/Data/AffineSpace.hs +++ b/src/Data/AffineSpace.hs @@ -1,4 +1,10 @@ -{-# LANGUAGE MultiParamTypeClasses, FlexibleContexts, TypeFamilies #-} +{-# LANGUAGE MultiParamTypeClasses, FlexibleContexts, TypeFamilies, CPP #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE ScopedTypeVariables #-} +{-# LANGUAGE TypeOperators #-} +{-# LANGUAGE UndecidableInstances #-} +{-# LANGUAGE DeriveGeneric #-} ---------------------------------------------------------------------- -- | -- Module : Data.AffineSpace @@ -13,26 +19,51 @@ module Data.AffineSpace ( - AffineSpace(..), (.-^), distanceSq, distance, alerp + AffineSpace(..), (.-^), distanceSq, distance, alerp, affineCombo ) where - +#if !MIN_VERSION_base(4,10,0) import Control.Applicative (liftA2) +#endif import Data.Ratio +import Foreign.C.Types (CSChar, CInt, CShort, CLong, CLLong, CIntMax, CFloat, CDouble) +import Control.Arrow(first) import Data.VectorSpace +import Data.Basis -infix 4 .+^, .-^, .-. +import Data.VectorSpace.Generic +import qualified GHC.Generics as Gnrx +import GHC.Generics (Generic, (:*:)(..)) +-- Through 0.8.4, I used the following fixities. +-- +-- infix 4 .+^, .-^, .-. +-- +-- Changed in 0.8.5 to match precedence of + and -, and to associate usefully. +-- Thanks to Ben Gamari for suggesting left-associativity. + +infixl 6 .+^, .-^ +infix 6 .-. + + -- TODO: Convert AffineSpace from fundep to associated type, and eliminate -- FunctionalDependencies above. class AdditiveGroup (Diff p) => AffineSpace p where -- | Associated vector space type Diff p + type Diff p = GenericDiff p -- | Subtract points (.-.) :: p -> p -> Diff p + default (.-.) :: ( Generic p, Diff p ~ GenericDiff p, AffineSpace (VRep p) ) + => p -> p -> Diff p + p .-. q = GenericDiff + $ (Gnrx.from p .-. (Gnrx.from q :: VRep p)) -- | Point plus vector (.+^) :: p -> Diff p -> p + default (.+^) :: ( Generic p, Diff p ~ GenericDiff p, AffineSpace (VRep p) ) + => p -> Diff p -> p + p .+^ GenericDiff q = Gnrx.to (Gnrx.from p .+^ q :: VRep p) -- | Point minus vector (.-^) :: AffineSpace p => p -> Diff p -> p @@ -56,21 +87,46 @@ p -> p -> Scalar (Diff p) -> p alerp p p' s = p .+^ (s *^ (p' .-. p)) -instance AffineSpace Double where - type Diff Double = Double - (.-.) = (-) - (.+^) = (+) +-- | Compute an affine combination (weighted average) of points. +-- The first element is used as origin and is weighted +-- such that all coefficients sum to 1. For example, +-- +-- > affineCombo a [(0.3,b), (0.2,c)] +-- +-- is equal to +-- +-- > a .+^ (0.3 *^ (b .-. a) ^+^ 0.2 *^ (c .-. a)) +-- +-- and if @a@, @b@, and @c@ were in a vector space would also be equal to +-- +-- > 0.5 *^ a ^+^ 0.3 *^ b ^+^ 0.2 *^ c +-- +-- See also 'linearCombo' (on vector spaces). +affineCombo :: (AffineSpace p, v ~ Diff p, VectorSpace v) => p -> [(p,Scalar v)] -> p +affineCombo z l = z .+^ linearCombo (map (first (.-. z)) l) -instance AffineSpace Float where - type Diff Float = Float - (.-.) = (-) - (.+^) = (+) +#define ScalarTypeCon(con,t) \ + instance con => AffineSpace (t) where \ + { type Diff (t) = t \ + ; (.-.) = (-) \ + ; (.+^) = (+) } -instance Integral a => AffineSpace (Ratio a) where - type Diff (Ratio a) = Ratio a - (.-.) = (-) - (.+^) = (+) +#define ScalarType(t) ScalarTypeCon((),t) +ScalarType(Int) +ScalarType(Integer) +ScalarType(Double) +ScalarType(Float) +ScalarType(CSChar) +ScalarType(CInt) +ScalarType(CShort) +ScalarType(CLong) +ScalarType(CLLong) +ScalarType(CIntMax) +ScalarType(CDouble) +ScalarType(CFloat) +ScalarTypeCon(Integral a,Ratio a) + instance (AffineSpace p, AffineSpace q) => AffineSpace (p,q) where type Diff (p,q) = (Diff p, Diff q) (p,q) .-. (p',q') = (p .-. p', q .-. q') @@ -86,3 +142,59 @@ type Diff (a -> p) = a -> Diff p (.-.) = liftA2 (.-.) (.+^) = liftA2 (.+^) + + + +newtype GenericDiff p = GenericDiff (Diff (VRep p)) + deriving (Generic) + +instance AdditiveGroup (Diff (VRep p)) => AdditiveGroup (GenericDiff p) +instance VectorSpace (Diff (VRep p)) => VectorSpace (GenericDiff p) +instance (AdditiveGroup (Scalar (Diff (VRep p))), InnerSpace (Diff (VRep p))) => InnerSpace (GenericDiff p) +instance HasBasis (Diff (VRep p)) => HasBasis (GenericDiff p) + +data AffineDiffProductSpace f g p = AffineDiffProductSpace + !(Diff (f p)) !(Diff (g p)) deriving (Generic) +instance (AffineSpace (f p), AffineSpace (g p)) + => AdditiveGroup (AffineDiffProductSpace f g p) +instance ( AffineSpace (f p), AffineSpace (g p) + , VectorSpace (Diff (f p)), VectorSpace (Diff (g p)) + , Scalar (Diff (f p)) ~ Scalar (Diff (g p)) ) + => VectorSpace (AffineDiffProductSpace f g p) +instance ( AdditiveGroup (Scalar (Diff (g p))) + , AffineSpace (f p), AffineSpace (g p) + , InnerSpace (Diff (f p)), InnerSpace (Diff (g p)) + , Scalar (Diff (f p)) ~ Scalar (Diff (g p)) + , Num (Scalar (Diff (f p))) ) + => InnerSpace (AffineDiffProductSpace f g p) +instance (AffineSpace (f p), AffineSpace (g p)) + => AffineSpace (AffineDiffProductSpace f g p) where + type Diff (AffineDiffProductSpace f g p) = AffineDiffProductSpace f g p + (.+^) = (^+^) + (.-.) = (^-^) +instance ( AffineSpace (f p), AffineSpace (g p) + , HasBasis (Diff (f p)), HasBasis (Diff (g p)) + , Scalar (Diff (f p)) ~ Scalar (Diff (g p)) ) + => HasBasis (AffineDiffProductSpace f g p) where + type Basis (AffineDiffProductSpace f g p) = Either (Basis (Diff (f p))) + (Basis (Diff (g p))) + basisValue (Left bf) = AffineDiffProductSpace (basisValue bf) zeroV + basisValue (Right bg) = AffineDiffProductSpace zeroV (basisValue bg) + decompose (AffineDiffProductSpace vf vg) + = map (first Left) (decompose vf) ++ map (first Right) (decompose vg) + decompose' (AffineDiffProductSpace vf _) (Left bf) = decompose' vf bf + decompose' (AffineDiffProductSpace _ vg) (Right bg) = decompose' vg bg + + +instance AffineSpace a => AffineSpace (Gnrx.Rec0 a s) where + type Diff (Gnrx.Rec0 a s) = Diff a + Gnrx.K1 v .+^ w = Gnrx.K1 $ v .+^ w + Gnrx.K1 v .-. Gnrx.K1 w = v .-. w +instance AffineSpace (f p) => AffineSpace (Gnrx.M1 i c f p) where + type Diff (Gnrx.M1 i c f p) = Diff (f p) + Gnrx.M1 v .+^ w = Gnrx.M1 $ v .+^ w + Gnrx.M1 v .-. Gnrx.M1 w = v .-. w +instance (AffineSpace (f p), AffineSpace (g p)) => AffineSpace ((f :*: g) p) where + type Diff ((f:*:g) p) = AffineDiffProductSpace f g p + (x:*:y) .+^ AffineDiffProductSpace ξ υ = (x.+^ξ) :*: (y.+^υ) + (x:*:y) .-. (ξ:*:υ) = AffineDiffProductSpace (x.-.ξ) (y.-.υ) diff --git a/src/Data/Basis.hs b/src/Data/Basis.hs --- a/src/Data/Basis.hs +++ b/src/Data/Basis.hs @@ -1,10 +1,8 @@ --- WARNING: this module depends on type families working fairly well, and --- requires ghc version at least 6.9. I didn't find a way to specify that --- dependency in the .cabal. --- {-# LANGUAGE TypeOperators, TypeFamilies, UndecidableInstances - , FlexibleInstances, MultiParamTypeClasses - #-} + , FlexibleInstances, MultiParamTypeClasses, CPP #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE ScopedTypeVariables #-} {-# OPTIONS_GHC -Wall -fno-warn-orphans #-} ---------------------------------------------------------------------- -- | @@ -24,37 +22,50 @@ -- import Control.Applicative ((<$>)) import Control.Arrow (first) import Data.Ratio +import Foreign.C.Types (CFloat, CDouble) +import Data.Kind -- import Data.Either import Data.VectorSpace +import Data.VectorSpace.Generic +import qualified GHC.Generics as Gnrx +import GHC.Generics (Generic, (:*:)(..)) + -- using associated data type instead of associated type synonym to work -- around ghc bug class VectorSpace v => HasBasis v where -- | Representation of the canonical basis for @v@ - type Basis v :: * + type Basis v :: Type + type Basis v = Basis (VRep v) -- | Interpret basis rep as a vector basisValue :: Basis v -> v + default basisValue :: (Generic v, HasBasis (VRep v), Basis (VRep v) ~ Basis v) + => Basis v -> v + basisValue b = Gnrx.to (basisValue b :: VRep v) -- | Extract coordinates decompose :: v -> [(Basis v, Scalar v)] + default decompose :: ( Generic v, HasBasis (VRep v) + , Scalar (VRep v) ~ Scalar v, Basis (VRep v) ~ Basis v ) + => v -> [(Basis v, Scalar v)] + decompose v = decompose (Gnrx.from v :: VRep v) -- | Experimental version. More elegant definitions, and friendly to -- infinite-dimensional vector spaces. decompose' :: v -> (Basis v -> Scalar v) + default decompose' :: ( Generic v, HasBasis (VRep v) + , Scalar (VRep v) ~ Scalar v, Basis (VRep v) ~ Basis v ) + => v -> Basis v -> Scalar v + decompose' v = decompose' (Gnrx.from v :: VRep v) -- Defining property: recompose . decompose == id --- | Linear combination -linearCombo :: VectorSpace v => [(v,Scalar v)] -> v -linearCombo ps = sumV [s *^ v | (v,s) <- ps] - -- Turn a basis decomposition back into a vector. recompose :: HasBasis v => [(Basis v, Scalar v)] -> v recompose = linearCombo . fmap (first basisValue) -- recompose ps = linearCombo (first basisValue <$> ps) - -- I don't know how to define -- -- recompose' :: HasBasis v => (Basis v -> Scalar v) -> v @@ -62,23 +73,20 @@ -- However, I don't seem to use recompose anywhere. -- I don't even use basisValue or decompose. -instance HasBasis Float where - type Basis Float = () - basisValue () = 1 - decompose s = [((),s)] - decompose' s = const s +#define ScalarTypeCon(con,t) \ + instance con => HasBasis (t) where \ + { type Basis (t) = () \ + ; basisValue () = 1 \ + ; decompose s = [((),s)] \ + ; decompose' s = const s } -instance HasBasis Double where - type Basis Double = () - basisValue () = 1 - decompose s = [((),s)] - decompose' s = const s +#define ScalarType(t) ScalarTypeCon((),t) -instance Integral a => HasBasis (Ratio a) where - type Basis (Ratio a) = () - basisValue () = 1 - decompose s = [((),s)] - decompose' s = const s +ScalarType(Float) +ScalarType(CFloat) +ScalarType(Double) +ScalarType(CDouble) +ScalarTypeCon(Integral a, Ratio a) instance ( HasBasis u, s ~ Scalar u , HasBasis v, s ~ Scalar v ) @@ -143,3 +151,21 @@ t4 = basisValue (Right (Left ())) :: (Float,Double,Float) -} + +instance HasBasis a => HasBasis (Gnrx.Rec0 a s) where + type Basis (Gnrx.Rec0 a s) = Basis a + basisValue = Gnrx.K1 . basisValue + decompose = decompose . Gnrx.unK1 + decompose' = decompose' . Gnrx.unK1 +instance HasBasis (f p) => HasBasis (Gnrx.M1 i c f p) where + type Basis (Gnrx.M1 i c f p) = Basis (f p) + basisValue = Gnrx.M1 . basisValue + decompose = decompose . Gnrx.unM1 + decompose' = decompose' . Gnrx.unM1 +instance (HasBasis (f p), HasBasis (g p), Scalar (f p) ~ Scalar (g p)) + => HasBasis ((f :*: g) p) where + type Basis ((f:*:g) p) = Either (Basis (f p)) (Basis (g p)) + basisValue (Left bf) = basisValue bf :*: zeroV + basisValue (Right bg) = zeroV :*: basisValue bg + decompose (u:*:v) = decomp2 Left u ++ decomp2 Right v + decompose' (u:*:v) = decompose' u `either` decompose' v diff --git a/src/Data/Cross.hs b/src/Data/Cross.hs --- a/src/Data/Cross.hs +++ b/src/Data/Cross.hs @@ -1,6 +1,6 @@ {-# LANGUAGE FlexibleInstances, FlexibleContexts, TypeOperators - , TypeFamilies, TypeSynonymInstances - #-} + , TypeFamilies, TypeSynonymInstances + , UndecidableInstances #-} {-# OPTIONS_GHC -Wall #-} ---------------------------------------------------------------------- -- | @@ -49,8 +49,7 @@ instance AdditiveGroup u => HasCross2 (u,u) where cross2 (x,y) = (negateV y,x) -- or @(y,-x)@? -instance ( HasBasis a, HasTrie (Basis a) - , VectorSpace v, HasCross2 v) => HasCross2 (a:>v) where +instance (HasTrie (Basis a), HasCross2 v) => HasCross2 (a:>v) where -- 2d cross-product is linear cross2 = fmapD cross2 @@ -74,8 +73,7 @@ -- l `atB` b = maybe zeroV (`untrie` b) l -instance ( Num s, VectorSpace s - , HasBasis s, HasTrie (Basis s), Basis s ~ ()) +instance (VectorSpace s, HasBasis s, HasTrie (Basis s), Basis s ~ ()) => HasNormal (Two (One s :> s)) where normalVec = unpairD . normalVec . pairD @@ -102,7 +100,7 @@ where d = derivAtBasis v -instance ( Num s, VectorSpace s, HasBasis s, HasTrie (Basis s) - , HasNormal (Two s :> Three s)) +instance ( VectorSpace s, HasBasis s, HasTrie (Basis s) + , HasNormal (Two s :> Three s) ) => HasNormal (Three (Two s :> s)) where normalVec = untripleD . normalVec . tripleD diff --git a/src/Data/LinearMap.hs b/src/Data/LinearMap.hs --- a/src/Data/LinearMap.hs +++ b/src/Data/LinearMap.hs @@ -1,9 +1,11 @@ -{-# LANGUAGE TypeOperators, FlexibleContexts, TypeFamilies, GeneralizedNewtypeDeriving, StandaloneDeriving #-} +{-# LANGUAGE TupleSections #-} +{-# LANGUAGE CPP, TypeOperators, FlexibleContexts, TypeFamilies + , GeneralizedNewtypeDeriving, StandaloneDeriving, UndecidableInstances #-} {-# OPTIONS_GHC -Wall -fno-warn-orphans #-} ---------------------------------------------------------------------- -- | -- Module : Data.LinearMap --- Copyright : (c) Conal Elliott 2008 +-- Copyright : (c) Conal Elliott 2008-2016 -- License : BSD3 -- -- Maintainer : conal@conal.net @@ -17,19 +19,22 @@ , inLMap, inLMap2, inLMap3 , liftMS, liftMS2, liftMS3 , liftL, liftL2, liftL3 - , firstL + , exlL, exrL, forkL, firstL, secondL + , inlL, inrL, joinL -- , leftL, rightL ) where -import Control.Applicative (Applicative,liftA2,liftA3) -import Control.Arrow (first) +#if !(MIN_VERSION_base(4,8,0)) +import Control.Applicative (Applicative, liftA2) +#endif +import Control.Applicative (liftA3) +import Control.Arrow (first,second) -import Data.MemoTrie ((:->:)(..)) +import Data.MemoTrie (HasTrie(..),(:->:)) import Data.AdditiveGroup (Sum(..), AdditiveGroup(..)) import Data.VectorSpace (VectorSpace(..)) import Data.Basis (HasBasis(..), linearCombo) - -- Linear maps are almost but not quite a Control.Category. The type -- class constraints interfere. They're almost an Arrow also, but for the -- constraints and the generality of arr. @@ -42,6 +47,7 @@ type LMap' u v = MSum (Basis u :->: v) +infixr 1 :-* -- | Linear map, represented as an optional memo-trie from basis to -- values, where 'Nothing' means the zero map (an optimization). newtype u :-* v = LMap { unLMap :: LMap' u v } @@ -53,13 +59,39 @@ type Scalar (u :-* v) = Scalar v (*^) s = (inLMap.liftMS.fmap) (s *^) -firstL :: ( HasBasis u, HasBasis u', HasBasis v - , HasTrie (Basis u), HasTrie (Basis v) - , Scalar u ~ Scalar v, Scalar u ~ Scalar u' - ) => - (u :-* u') -> ((u,v) :-* (u',v)) -firstL = linear.first.lapply +-- In GHC 7.10: +-- Constraint is no smaller than the instance head +-- in the constraint: HasTrie (Basis u) +-- (Use UndecidableInstances to permit this) +exlL :: ( HasBasis a, HasTrie (Basis a), HasBasis b, HasTrie (Basis b) + , Scalar a ~ Scalar b ) + => (a,b) :-* a +exlL = linear fst + +exrL :: ( HasBasis a, HasTrie (Basis a), HasBasis b, HasTrie (Basis b) + , Scalar a ~ Scalar b ) + => (a,b) :-* b +exrL = linear snd + +forkL :: (HasTrie (Basis a), HasBasis c, HasBasis d) + => (a :-* c) -> (a :-* d) -> (a :-* (c,d)) +forkL = (inLMap2.liftL2) (,) + +firstL :: ( HasBasis u, HasBasis u', HasBasis v + , HasTrie (Basis u), HasTrie (Basis v) + , Scalar u ~ Scalar v, Scalar u ~ Scalar u' + ) => + (u :-* u') -> ((u,v) :-* (u',v)) +firstL = linear.first.lapply + +secondL :: ( HasBasis u, HasBasis v, HasBasis v' + , HasTrie (Basis u), HasTrie (Basis v) + , Scalar u ~ Scalar v, Scalar u ~ Scalar v' + ) => + (v :-* v') -> ((u,v) :-* (u,v')) +secondL = linear.second.lapply + -- TODO: more efficient firstL -- liftMS :: (AdditiveGroup a) => (a -> b) -> (MSum a -> MSum b) @@ -70,6 +102,21 @@ -- (inLMap.liftMS.fmap) (s *^) :: (u :-* v) -> (u :-* v) +inlL :: (HasBasis a, HasTrie (Basis a), HasBasis b) + => a :-* (a,b) +inlL = linear (,zeroV) + +inrL :: (HasBasis a, HasBasis b, HasTrie (Basis b)) + => b :-* (a,b) +inrL = linear (zeroV,) + +joinL :: ( HasBasis a, HasTrie (Basis a) + , HasBasis b, HasTrie (Basis b) + , Scalar a ~ Scalar b, Scalar a ~ Scalar c + , VectorSpace c ) + => (a :-* c) -> (b :-* c) -> ((a,b) :-* c) +f `joinL` g = linear (\ (a,b) -> lapply f a ^+^ lapply g b) + -- Before version 0.7, u :-* v was a type synonym, resulting in a subtle -- ambiguity: u:-*v == u':-*v' does not imply that u==u', since Basis -- might map different types to the same basis (e.g., Float & Double). @@ -129,7 +176,7 @@ infixr 9 *.* -- | Compose linear maps -(*.*) :: ( HasBasis u, HasTrie (Basis u) +(*.*) :: ( HasTrie (Basis u) , HasBasis v, HasTrie (Basis v) , VectorSpace w , Scalar v ~ Scalar w ) => @@ -174,9 +221,7 @@ -- to values and then decomposed, followed by recombination of the -- results. -liftMS :: (AdditiveGroup a) => - (a -> b) - -> (MSum a -> MSum b) +liftMS :: (a -> b) -> (MSum a -> MSum b) -- liftMS _ Nothing = Nothing -- liftMS h ma = Just (Sum (h (z ma))) @@ -201,8 +246,7 @@ -- | Apply a linear function to each element of a linear map. -- @liftL f l == linear f *.* l@, but works more efficiently. -liftL :: (Functor f, AdditiveGroup (f a)) => - (a -> b) -> MSum (f a) -> MSum (f b) +liftL :: Functor f => (a -> b) -> MSum (f a) -> MSum (f b) liftL = liftMS . fmap -- | Apply a linear binary function (not to be confused with a bilinear @@ -253,7 +297,6 @@ -} - ----- diff --git a/src/Data/Maclaurin.hs b/src/Data/Maclaurin.hs --- a/src/Data/Maclaurin.hs +++ b/src/Data/Maclaurin.hs @@ -1,8 +1,7 @@ {-# LANGUAGE TypeOperators, MultiParamTypeClasses, UndecidableInstances , TypeSynonymInstances, FlexibleInstances , FlexibleContexts, TypeFamilies - , ScopedTypeVariables - #-} + , ScopedTypeVariables, CPP #-} -- The ScopedTypeVariables is there just as a bug work-around. Without it -- I get a bogus error about context mismatch for mutually recursive @@ -25,7 +24,7 @@ -- Stability : experimental -- -- Infinite derivative towers via linear maps, using the Maclaurin --- representation. See blog posts . +-- representation. See blog posts . ---------------------------------------------------------------------- module Data.Maclaurin @@ -53,6 +52,10 @@ import Data.Boolean +#if MIN_VERSION_base(4,8,0) +import Prelude hiding ((<*)) +#endif + infixr 9 `D` -- | Tower of derivatives. data a :> b = D { powVal :: b, derivative :: a :-* (a :> b) } @@ -71,8 +74,7 @@ infixl 4 <$>> -- | Map a /linear/ function over a derivative tower. -fmapD, (<$>>) :: (HasBasis a, HasTrie (Basis a), AdditiveGroup b) => - (b -> c) -> (a :> b) -> (a :> c) +fmapD, (<$>>) :: HasTrie (Basis a) => (b -> c) -> (a :> b) -> (a :> c) fmapD f = lf where lf (D b0 b') = D (f b0) ((inLMap.liftL) lf b') @@ -108,9 +110,7 @@ -- | Differentiable identity function. Sometimes called "the -- derivation variable" or similar, but it's not really a variable. -idD :: ( VectorSpace u, s ~ Scalar u - , VectorSpace (u :> u), VectorSpace s - , HasBasis u, HasTrie (Basis u)) => +idD :: (VectorSpace u , HasBasis u, HasTrie (Basis u)) => u :~> u idD = linearD id @@ -165,9 +165,7 @@ -- | Derivative tower for applying a binary function that distributes over -- addition, such as multiplication. A bit weaker assumption than -- bilinearity. Is bilinearity necessary for correctness here? -distrib :: forall a b c u. - ( HasBasis a, HasTrie (Basis a) - , AdditiveGroup b, AdditiveGroup c, AdditiveGroup u) => +distrib :: forall a b c u. (HasBasis a, HasTrie (Basis a) , AdditiveGroup u) => (b -> c -> u) -> (a :> b) -> (a :> c) -> (a :> u) distrib op = (#) @@ -187,18 +185,19 @@ instance Show b => Show (a :> b) where show (D b0 _) = "D " ++ show b0 ++ " ..." -instance Eq b => Eq (a :> b) where (==) = noOv "(==)" +instance Eq (a :> b) where (==) = noOv "(==)" -instance (AdditiveGroup v, HasBasis u, HasTrie (Basis u), IfB b v) => - IfB b (u :> v) where +type instance BooleanOf (a :> b) = BooleanOf b + +instance (AdditiveGroup v, HasBasis u, HasTrie (Basis u), IfB v) => + IfB (u :> v) where ifB = liftD2 . ifB -instance (AdditiveGroup v, HasBasis u, HasTrie (Basis u), OrdB b v) => - OrdB b (u :> v) where +instance OrdB v => OrdB (u :> v) where (<*) = (<*) `on` powVal instance ( AdditiveGroup b, HasBasis a, HasTrie (Basis a) - , OrdB bool b, IfB bool b, Ord b) => Ord (a :> b) where + , OrdB b, IfB b, Ord b) => Ord (a :> b) where compare = compare `on` powVal min = minB max = maxB @@ -213,8 +212,7 @@ -- Less efficient: adds zero -- (^+^) = liftD2 (^+^) -instance ( HasBasis a, HasTrie (Basis a) - , VectorSpace u, AdditiveGroup (Scalar u) ) +instance (HasBasis a, HasTrie (Basis a), VectorSpace u) => VectorSpace (a :> u) where type Scalar (a :> u) = (a :> Scalar u) (*^) = distrib (*^) @@ -236,8 +234,7 @@ infix 0 >-< -- | Specialized chain rule. See also '(\@.)' -(>-<) :: ( HasBasis a, HasTrie (Basis a), VectorSpace u - , AdditiveGroup (Scalar u)) => +(>-<) :: (HasBasis a, HasTrie (Basis a), VectorSpace u) => (u -> u) -> ((a :> u) -> (a :> Scalar u)) -> (a :> u) -> (a :> u) f >-< f' = \ u@(D u0 u') -> D (f u0) ((inLMap.liftMS) (f' u *^) u') @@ -293,31 +290,21 @@ ---- Misc -pairD :: ( HasBasis a, HasTrie (Basis a) - , VectorSpace b, VectorSpace c - , Scalar b ~ Scalar c - ) => (a:>b,a:>c) -> a:>(b,c) +pairD :: (HasBasis a, HasTrie (Basis a), VectorSpace b, VectorSpace c) + => (a:>b,a:>c) -> a:>(b,c) pairD (u,v) = liftD2 (,) u v -unpairD :: ( HasBasis a, HasTrie (Basis a) - , VectorSpace a, VectorSpace b, VectorSpace c - , Scalar b ~ Scalar c - ) => (a :> (b,c)) -> (a:>b, a:>c) +unpairD :: HasTrie (Basis a) => (a :> (b,c)) -> (a:>b, a:>c) unpairD d = (fst <$>> d, snd <$>> d) tripleD :: ( HasBasis a, HasTrie (Basis a) , VectorSpace b, VectorSpace c, VectorSpace d - , Scalar b ~ Scalar c, Scalar c ~ Scalar d ) => (a:>b,a:>c,a:>d) -> a:>(b,c,d) tripleD (u,v,w) = liftD3 (,,) u v w -untripleD :: ( HasBasis a, HasTrie (Basis a) - , VectorSpace a, VectorSpace b, VectorSpace c, VectorSpace d - , Scalar b ~ Scalar c, Scalar c ~ Scalar d - ) => - (a :> (b,c,d)) -> (a:>b, a:>c, a:>d) +untripleD :: HasTrie (Basis a) => (a :> (b,c,d)) -> (a:>b, a:>c, a:>d) untripleD d = ((\ (a,_,_) -> a) <$>> d, (\ (_,b,_) -> b) <$>> d, (\ (_,_,c) -> c) <$>> d) diff --git a/src/Data/VectorSpace.hs b/src/Data/VectorSpace.hs --- a/src/Data/VectorSpace.hs +++ b/src/Data/VectorSpace.hs @@ -1,18 +1,21 @@ {-# LANGUAGE MultiParamTypeClasses, TypeOperators - , TypeFamilies, UndecidableInstances - #-} + , TypeFamilies, UndecidableInstances, CPP + , FlexibleContexts #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE DefaultSignatures #-} +{-# LANGUAGE ScopedTypeVariables #-} {-# OPTIONS_GHC -Wall #-} ---------------------------------------------------------------------- -- | -- Module : Data.VectorSpace -- Copyright : (c) Conal Elliott and Andy J Gill 2008 -- License : BSD3 --- +-- -- Maintainer : conal@conal.net, andygill@ku.edu -- Stability : experimental --- +-- -- Vector spaces --- +-- -- This version uses associated types instead of fundeps and -- requires ghc-6.10 or later ---------------------------------------------------------------------- @@ -20,83 +23,119 @@ -- NB: I'm attempting to replace fundeps with associated types. See -- NewVectorSpace.hs. Ran into trouble with type equality constraints not -- getting propagated. Manuel Ch is looking into it. --- +-- -- Blocking bug: http://hackage.haskell.org/trac/ghc/ticket/2448 module Data.VectorSpace ( module Data.AdditiveGroup , VectorSpace(..), (^/), (^*) , InnerSpace(..) - , lerp, magnitudeSq, magnitude, normalized, project + , lerp, linearCombo, magnitudeSq, magnitude, normalized, project ) where - +#if !(MIN_VERSION_base(4,8,0)) import Control.Applicative (liftA2) +#endif import Data.Complex hiding (magnitude) +import Foreign.C.Types (CSChar, CInt, CShort, CLong, CLLong, CIntMax, CFloat, CDouble) import Data.Ratio +import Data.Kind import Data.AdditiveGroup import Data.MemoTrie +import Data.VectorSpace.Generic +import qualified GHC.Generics as Gnrx +import GHC.Generics (Generic, (:*:)(..)) + infixr 7 *^ -- | Vector space @v@. class AdditiveGroup v => VectorSpace v where - type Scalar v :: * + type Scalar v :: Type + type Scalar v = Scalar (VRep v) -- | Scale a vector (*^) :: Scalar v -> v -> v + default (*^) :: (Generic v, VectorSpace (VRep v), Scalar (VRep v) ~ Scalar v) + => Scalar v -> v -> v + μ *^ v = Gnrx.to (μ *^ Gnrx.from v :: VRep v) + {-# INLINE (*^) #-} infixr 7 <.> -- | Adds inner (dot) products. -class VectorSpace v => InnerSpace v where +class (VectorSpace v, AdditiveGroup (Scalar v)) => InnerSpace v where -- | Inner/dot product (<.>) :: v -> v -> Scalar v + default (<.>) :: (Generic v, InnerSpace (VRep v), Scalar (VRep v) ~ Scalar v) + => v -> v -> Scalar v + v<.>w = (Gnrx.from v :: VRep v) <.> Gnrx.from w + {-# INLINE (<.>) #-} infixr 7 ^/ infixl 7 ^* -- | Vector divided by scalar (^/) :: (VectorSpace v, s ~ Scalar v, Fractional s) => v -> s -> v -v ^/ s = (1/s) *^ v +v ^/ s = recip s *^ v +{-# INLINE (^/) #-} -- | Vector multiplied by scalar (^*) :: (VectorSpace v, s ~ Scalar v) => v -> s -> v (^*) = flip (*^) +{-# INLINE (^*) #-} -- | Linear interpolation between @a@ (when @t==0@) and @b@ (when @t==1@). -- lerp :: (VectorSpace v, s ~ Scalar v, Num s) => v -> v -> s -> v lerp :: VectorSpace v => v -> v -> Scalar v -> v lerp a b t = a ^+^ t *^ (b ^-^ a) +{-# INLINE lerp #-} +-- | Linear combination of vectors +linearCombo :: VectorSpace v => [(v,Scalar v)] -> v +linearCombo ps = sumV [v ^* s | (v,s) <- ps] +{-# INLINE linearCombo #-} + -- | Square of the length of a vector. Sometimes useful for efficiency. -- See also 'magnitude'. magnitudeSq :: (InnerSpace v, s ~ Scalar v) => v -> s magnitudeSq v = v <.> v +{-# INLINE magnitudeSq #-} -- | Length of a vector. See also 'magnitudeSq'. magnitude :: (InnerSpace v, s ~ Scalar v, Floating s) => v -> s magnitude = sqrt . magnitudeSq +{-# INLINE magnitude #-} --- | Vector in same direction as given one but with length of one. If --- given the zero vector, then return it. +-- | Vector in same direction as given one but with length of one. +-- Divides by zero for the zero vector. normalized :: (InnerSpace v, s ~ Scalar v, Floating s) => v -> v normalized v = v ^/ magnitude v +{-# INLINE normalized #-} -- | @project u v@ computes the projection of @v@ onto @u@. -project :: (InnerSpace v, s ~ Scalar v, Floating s) => v -> v -> v -project u v = (v <.> u') *^ u' - where u' = normalized u +project :: (InnerSpace v, s ~ Scalar v, Fractional s) => v -> v -> v +project u v = ((v <.> u) / magnitudeSq u) *^ u +{-# INLINE project #-} -instance VectorSpace Double where - type Scalar Double = Double - (*^) = (*) -instance InnerSpace Double where (<.>) = (*) +#define ScalarType(t) \ + instance VectorSpace (t) where \ + { type Scalar (t) = (t) \ + ; (*^) = (*) } ; \ + instance InnerSpace (t) where (<.>) = (*) -instance VectorSpace Float where - type Scalar Float = Float - (*^) = (*) -instance InnerSpace Float where (<.>) = (*) +ScalarType(Int) +ScalarType(Integer) +ScalarType(Double) +ScalarType(Float) +ScalarType(CSChar) +ScalarType(CInt) +ScalarType(CShort) +ScalarType(CLong) +ScalarType(CLLong) +ScalarType(CIntMax) +ScalarType(CDouble) +ScalarType(CFloat) instance Integral a => VectorSpace (Ratio a) where type Scalar (Ratio a) = Ratio a @@ -107,7 +146,7 @@ type Scalar (Complex v) = Scalar v s*^(u :+ v) = s*^u :+ s*^v -instance (RealFloat v, InnerSpace v, s ~ Scalar v, AdditiveGroup s) +instance (RealFloat v, InnerSpace v) => InnerSpace (Complex v) where (u :+ v) <.> (u' :+ v') = (u <.> u') ^+^ (v <.> v') @@ -127,8 +166,7 @@ s *^ (u,v) = (s*^u,s*^v) instance ( InnerSpace u, s ~ Scalar u - , InnerSpace v, s ~ Scalar v - , AdditiveGroup (Scalar v) ) + , InnerSpace v, s ~ Scalar v ) => InnerSpace (u,v) where (u,v) <.> (u',v') = (u <.> u') ^+^ (v <.> v') @@ -141,8 +179,7 @@ instance ( InnerSpace u, s ~ Scalar u , InnerSpace v, s ~ Scalar v - , InnerSpace w, s ~ Scalar w - , AdditiveGroup s ) + , InnerSpace w, s ~ Scalar w ) => InnerSpace (u,v,w) where (u,v,w) <.> (u',v',w') = u<.>u' ^+^ v<.>v' ^+^ w<.>w' @@ -157,8 +194,7 @@ instance ( InnerSpace u, s ~ Scalar u , InnerSpace v, s ~ Scalar v , InnerSpace w, s ~ Scalar w - , InnerSpace x, s ~ Scalar x - , AdditiveGroup s ) + , InnerSpace x, s ~ Scalar x ) => InnerSpace (u,v,w,x) where (u,v,w,x) <.> (u',v',w',x') = u<.>u' ^+^ v<.>v' ^+^ w<.>w' ^+^ x<.>x' @@ -196,7 +232,7 @@ -instance (InnerSpace a, AdditiveGroup (Scalar a)) => InnerSpace (Maybe a) where +instance InnerSpace a => InnerSpace (Maybe a) where -- dotting with zero (vector) yields zero (scalar) Nothing <.> _ = zeroV _ <.> Nothing = zeroV @@ -205,3 +241,30 @@ -- mu <.> mv = fromMaybe zeroV (liftA2 (<.>) mu mv) -- (<.>) = (fmap.fmap) (fromMaybe zeroV) (liftA2 (<.>)) + + +instance VectorSpace a => VectorSpace (Gnrx.Rec0 a s) where + type Scalar (Gnrx.Rec0 a s) = Scalar a + μ *^ Gnrx.K1 v = Gnrx.K1 $ μ*^v + {-# INLINE (*^) #-} +instance VectorSpace (f p) => VectorSpace (Gnrx.M1 i c f p) where + type Scalar (Gnrx.M1 i c f p) = Scalar (f p) + μ *^ Gnrx.M1 v = Gnrx.M1 $ μ*^v + {-# INLINE (*^) #-} +instance (VectorSpace (f p), VectorSpace (g p), Scalar (f p) ~ Scalar (g p)) + => VectorSpace ((f :*: g) p) where + type Scalar ((f:*:g) p) = Scalar (f p) + μ *^ (x:*:y) = μ*^x :*: μ*^y + {-# INLINE (*^) #-} + +instance InnerSpace a => InnerSpace (Gnrx.Rec0 a s) where + Gnrx.K1 v <.> Gnrx.K1 w = v<.>w + {-# INLINE (<.>) #-} +instance InnerSpace (f p) => InnerSpace (Gnrx.M1 i c f p) where + Gnrx.M1 v <.> Gnrx.M1 w = v<.>w + {-# INLINE (<.>) #-} +instance ( InnerSpace (f p), InnerSpace (g p) + , Scalar (f p) ~ Scalar (g p), Num (Scalar (f p)) ) + => InnerSpace ((f :*: g) p) where + (x:*:y) <.> (ξ:*:υ) = x<.>ξ + y<.>υ + {-# INLINE (<.>) #-} diff --git a/src/Data/VectorSpace/Generic.hs b/src/Data/VectorSpace/Generic.hs new file mode 100644 --- /dev/null +++ b/src/Data/VectorSpace/Generic.hs @@ -0,0 +1,20 @@ +-- | +-- Module : Data.VectorSpace.Generic +-- Copyright : (c) Conal Elliott and Justus Sagemüller 2017 +-- License : BSD3 +-- +-- Maintainer : conal@conal.net, (@) jsagemue $ uni-koeln.de +-- Stability : experimental +-- +-- Underpinnings of the type that represents vector / affine / etc. spaces +-- with GHC generics + +module Data.VectorSpace.Generic where + + +import qualified GHC.Generics as Gnrx + +import Data.Void + + +type VRep v = Gnrx.Rep v Void diff --git a/vector-space.cabal b/vector-space.cabal --- a/vector-space.cabal +++ b/vector-space.cabal @@ -1,6 +1,6 @@ Name: vector-space -Version: 0.8.0 -Cabal-Version: >= 1.2 +Version: 0.19 +Cabal-Version: >= 1.10 Synopsis: Vector & affine spaces, linear maps, and derivatives Category: math Description: @@ -11,26 +11,30 @@ (scalars, vectors, matrices, ...). . /Warning/: this package depends on type families working fairly well, - and requires ghc version at least 6.9. + requiring GHC version at least 6.9. . Project wiki page: . - © 2008-2011 by Conal Elliott; BSD3 license. + © 2008-2012 by Conal Elliott; BSD3 license. Author: Conal Elliott Maintainer: conal@conal.net -Homepage: http://haskell.org/haskellwiki/vector-space -Package-Url: http://code.haskell.org/vector-space -Copyright: (c) 2008-2011 by Conal Elliott +Copyright: (c) 2008-2012 by Conal Elliott License: BSD3 License-File: COPYING Stability: experimental build-type: Simple +source-repository head + type: git + location: git://github.com/conal/vector-space.git + Library + default-language: Haskell2010 hs-Source-Dirs: src Extensions: - Build-Depends: base<5, MemoTrie >= 0.4.2, Boolean >= 0.0.1, - NumInstances >= 1.0 + Build-Depends: base<5, MemoTrie >= 0.5 + , Boolean >= 0.1.0 + , NumInstances >= 1.0 Exposed-Modules: Data.AdditiveGroup Data.VectorSpace @@ -41,13 +45,11 @@ Data.Derivative Data.Cross Data.AffineSpace - + Other-Modules: + Data.VectorSpace.Generic -- This library relies on type families working as well as in 6.10. - if impl(ghc < 6.10) { - buildable: False - } - ghc-options: -Wall -O2 --- ghc-prof-options: -prof -auto-all - --- For ghc-options: -ddump-simpl-stats -ddump-rules -ddump-simpl -ddump-simpl-phases + if impl(ghc < 6.10) { buildable: False } + if !impl(ghc >= 7.6) { Build-Depends: ghc-prim >= 0.2 } + if !impl(ghc >= 7.9) { Build-Depends: void >= 0.4 } + if !impl(ghc >= 8.0) { Build-Depends: semigroups >= 0.16 }