computational-algebra 0.1.3.6 → 0.2.0.0
raw patch · 10 files changed
+133/−671 lines, 10 filesdep +equational-reasoningdep +singletonsdep +sized-vectordep ~monad-loops
Dependencies added: equational-reasoning, singletons, sized-vector, type-natural
Dependency ranges changed: monad-loops
Files
- Algebra/Algorithms/Groebner.hs +31/−27
- Algebra/Algorithms/Groebner/Monomorphic.hs +14/−8
- Algebra/Internal.hs +4/−508
- Algebra/Ring/Noetherian.hs +9/−6
- Algebra/Ring/Polynomial.hs +51/−45
- Algebra/Ring/Polynomial/Monomorphic.hs +9/−11
- Algebra/Ring/Polynomial/Parser.hs +1/−1
- Monomorphic.hs +0/−56
- computational-algebra.cabal +13/−9
- examples/polymorphic.hs +1/−0
Algebra/Algorithms/Groebner.hs view
@@ -1,8 +1,7 @@-{-# LANGUAGE ConstraintKinds, DataKinds, FlexibleContexts, FlexibleInstances #-}-{-# LANGUAGE GADTs, MultiParamTypeClasses, NoImplicitPrelude #-}-{-# LANGUAGE ParallelListComp, RankNTypes, ScopedTypeVariables #-}-{-# LANGUAGE StandaloneDeriving, TemplateHaskell, TypeFamilies #-}-{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE ConstraintKinds, DataKinds, FlexibleContexts, FlexibleInstances #-}+{-# LANGUAGE GADTs, MultiParamTypeClasses, NoImplicitPrelude #-}+{-# LANGUAGE ParallelListComp, RankNTypes, ScopedTypeVariables #-}+{-# LANGUAGE StandaloneDeriving, TemplateHaskell, TypeFamilies, TypeOperators #-} {-# OPTIONS_GHC -fno-warn-type-defaults -fno-warn-orphans #-} module Algebra.Algorithms.Groebner ( -- * Polynomial division@@ -36,8 +35,14 @@ import Data.Maybe import Data.Proxy import Data.STRef+import Data.Type.Monomorphic+import Data.Type.Natural hiding (max, one, zero)+import Data.Vector.Sized hiding (all, drop, foldr, head, map,+ take, zipWith)+import qualified Data.Vector.Sized as V import Numeric.Algebra import Prelude hiding (Num (..), recip)+import Proof.Equational -- | Calculate a polynomial quotient and remainder w.r.t. second argument. divModPolynomial :: (IsMonomialOrder order, IsPolynomial r n, Field r)@@ -95,7 +100,7 @@ where calc acc = [ q | f <- acc, g <- acc, f /= g , let f0 = leadingMonomial f, let g0 = leadingMonomial g- , lcmMonomial f0 g0 /= zipWithV (+) f0 g0+ , lcmMonomial f0 g0 /= V.zipWithSame (+) f0 g0 , let q = sPolynomial f g `modPolynomial` acc, q /= zero ] @@ -112,12 +117,12 @@ (xs', ys') -> (x :- xs', y :- ys') padVec def (x :- xs) Nil = case padVec def xs Nil of- (xs', ys') -> case maxZR (sLengthV xs) of- Eql -> (x :- xs', def :- ys')+ (xs', ys') -> case maxZR (sLength xs) of+ Refl -> (x :- xs', def :- ys') padVec def Nil (y :- ys) = case padVec def Nil ys of- (xs', ys') -> case maxZL (sLengthV ys) of- Eql -> (def :- xs', y :- ys')+ (xs', ys') -> case maxZL (sLength ys) of+ Refl -> (def :- xs', y :- ys') combinations :: [a] -> [(a, a)] combinations xs = concat $ zipWith (map . (,)) xs $ drop 1 $ tails xs@@ -170,7 +175,7 @@ && (all (\k -> H.all ((/=k) . H.payload) rest) [(f, h), (g, h), (h, f), (h, g)]) && leadingMonomial h `divs` l) gs0- when (l /= zipWithV (+) f0 g0 && not redundant) $ do+ when (l /= V.zipWithSame (+) f0 g0 && not redundant) $ do len0 <- readSTRef len let qs = (H.toList gs0) s = sPolynomial f g `modPolynomial` map H.payload qs@@ -312,9 +317,9 @@ -> Ideal (OrderedPolynomial k ord (m :-: n)) thEliminationIdealWith ord n ideal = case singInstance n of- SingInstance -> toIdeal $ [ transformMonomial (dropV n) f+ SingInstance -> toIdeal $ [ transformMonomial (V.drop n) f | f <- calcGroebnerBasisWith ord ideal- , all (all (== 0) . take (toInt n) . toList . snd) $ getTerms f+ , all (all (== 0) . take (sNatToInt n) . toList . snd) $ getTerms f ] -- | Calculate n-th elimination ideal using the specified n-th elimination type order.@@ -328,9 +333,9 @@ -> Ideal (OrderedPolynomial k ord (m :-: n)) unsafeThEliminationIdealWith ord n ideal = case singInstance n of- SingInstance -> toIdeal $ [ transformMonomial (dropV n) f+ SingInstance -> toIdeal $ [ transformMonomial (V.drop n) f | f <- calcGroebnerBasisWith ord ideal- , all (all (== 0) . take (toInt n) . toList . snd) $ getTerms f+ , all (all (== 0) . take (sNatToInt n) . toList . snd) $ getTerms f ] -- | An intersection ideal of given ideals (using 'WeightedEliminationOrder').@@ -340,15 +345,15 @@ ) => Vector (Ideal (OrderedPolynomial r ord n)) k -> Ideal (OrderedPolynomial r ord n)-intersection Nil = Ideal $ singletonV one+intersection Nil = Ideal $ singleton one intersection idsv@(_ :- _) =- let sk = sLengthV idsv+ let sk = sLength idsv sn = sing :: SNat n ts = genVars (sk %+ sn) tis = zipWith (\ideal t -> mapIdeal ((t *) . shiftR sk) ideal) (toList idsv) ts j = foldr appendIdeal (principalIdeal (one - foldr (+) zero ts)) tis in case plusMinusEqR sn sk of- Eql -> case propToBoolLeq (plusLeqL sk sn) of+ Refl -> case propToBoolLeq (plusLeqL sk sn) of LeqTrueInstance -> thEliminationIdeal sk j -- | Ideal quotient by a principal ideals.@@ -357,8 +362,8 @@ -> OrderedPolynomial k ord n -> Ideal (OrderedPolynomial k ord n) quotByPrincipalIdeal i g =- case intersection (i :- (Ideal $ singletonV g) :- Nil) of- Ideal gs -> Ideal $ mapV (snd . head . (`divPolynomial` [g])) gs+ case intersection (i :- (Ideal $ singleton g) :- Nil) of+ Ideal gs -> Ideal $ V.map (snd . head . (`divPolynomial` [g])) gs -- | Ideal quotient by the given ideal. quotIdeal :: forall k ord n. (IsPolynomial k n, Field k, IsMonomialOrder ord)@@ -366,10 +371,10 @@ -> Ideal (OrderedPolynomial k ord n) -> Ideal (OrderedPolynomial k ord n) quotIdeal i (Ideal g) =- case singInstance (sLengthV g) of+ case singInstance (sLength g) of SingInstance ->- case singInstance (sLengthV g %+ (sing :: SNat n)) of- SingInstance -> intersection $ mapV (i `quotByPrincipalIdeal`) g+ case singInstance (sLength g %+ (sing :: SNat n)) of+ SingInstance -> intersection $ V.map (i `quotByPrincipalIdeal`) g -- | Saturation by a principal ideal. saturationByPrincipalIdeal :: (Field k, IsPolynomial k n, IsMonomialOrder ord)@@ -385,8 +390,7 @@ -> Ideal (OrderedPolynomial k ord n) -> Ideal (OrderedPolynomial k ord n) saturationIdeal i (Ideal g) =- case singInstance (sLengthV g) of+ case singInstance (sLength g) of SingInstance ->- case singInstance (sLengthV g %+ (sing :: SNat n)) of- SingInstance -> intersection $ mapV (i `saturationByPrincipalIdeal`) g-+ case singInstance (sLength g %+ (sing :: SNat n)) of+ SingInstance -> intersection $ V.map (i `saturationByPrincipalIdeal`) g
Algebra/Algorithms/Groebner/Monomorphic.hs view
@@ -1,6 +1,7 @@-{-# LANGUAGE ConstraintKinds, FlexibleContexts, FlexibleInstances, GADTs #-}-{-# LANGUAGE PolyKinds, RecordWildCards, ScopedTypeVariables, TypeFamilies #-}-{-# LANGUAGE TypeOperators, UndecidableInstances #-}+{-# LANGUAGE ConstraintKinds, FlexibleContexts, FlexibleInstances, GADTs #-}+{-# LANGUAGE IncoherentInstances, OverlappingInstances, PolyKinds #-}+{-# LANGUAGE RecordWildCards, ScopedTypeVariables, TypeFamilies #-}+{-# LANGUAGE TypeOperators, UndecidableInstances #-} -- | Monomorphic interface for Groenber basis. module Algebra.Algorithms.Groebner.Monomorphic ( Groebnerable@@ -38,6 +39,13 @@ import Control.Arrow import Data.List import qualified Data.Map as M+import Data.Singletons hiding (demote, promote)+import Data.Type.Monomorphic+import Data.Type.Natural hiding (demote, one,+ promote, zero)+import Data.Vector.Sized (Vector (..), sLength,+ toList)+import qualified Data.Vector.Sized as V import Numeric.Algebra import Prelude hiding (Num (..)) @@ -57,10 +65,10 @@ SingInstance -> case promote ps :: Monomorphic (Vector [Polynomial r]) of Monomorphic vec ->- let slen = sLengthV vec+ let slen = sLength vec in case singInstance slen of SingInstance ->- let ids = mapV (toIdeal . map (flip orderedBy Lex . Poly.polynomial . M.mapKeys (Poly.OrderedMonomial . Poly.fromList sdim . encodeMonomList vars) . unPolynomial)) vec+ let ids = V.map (toIdeal . map (flip orderedBy Lex . Poly.polynomial . M.mapKeys (Poly.OrderedMonomial . Poly.fromList sdim . encodeMonomList vars) . unPolynomial)) vec in case singInstance (slen %+ sdim) of SingInstance -> demoteComposed $ Gr.intersection ids @@ -220,7 +228,7 @@ let fs' = map ((flip Poly.orderedBy Poly.Lex) . Poly.scastPolynomial newDim) fs in case propToBoolLeq $ maxLeqL sk sdim of LeqTrueInstance ->- case singInstance (newDim %- sk) of+ case singInstance (newDim %:- sk) of SingInstance -> map (renameVars rest) $ demoteComposed $ Gr.unsafeThEliminationIdealWith ord sk (toIdeal fs') where@@ -230,8 +238,6 @@ eliminate :: forall r. Groebnerable r => [Variable] -> [Polynomial r] -> [Polynomial r] eliminate vs j = eliminateWith Lex vs j- where- vars = nub $ sort $ concatMap buildVarsList j -- | Computes nth elimination ideal. thEliminationIdeal :: Groebnerable r => Int -> [Polynomial r] -> [Polynomial r]
Algebra/Internal.hs view
@@ -1,514 +1,10 @@-{-# LANGUAGE RankNTypes #-}-{-# LANGUAGE DataKinds, FlexibleContexts, FlexibleInstances, GADTs #-}-{-# LANGUAGE MultiParamTypeClasses, PolyKinds, StandaloneDeriving #-}-{-# LANGUAGE TypeFamilies, TypeOperators #-}+{-# LANGUAGE DataKinds, FlexibleContexts, FlexibleInstances, GADTs #-}+{-# LANGUAGE MultiParamTypeClasses, PolyKinds, RankNTypes, StandaloneDeriving #-}+{-# LANGUAGE TypeFamilies, TypeOperators #-} {-# OPTIONS_GHC -fwarn-incomplete-patterns #-}-module Algebra.Internal ( toProxy, Nat(..), SNat(..), Vector(..), Sing(..)- , SingInstance(..), singInstance, toInt- , Min, Max, sMin, sMax, sZ, sS, (:+:), (%+), (:-:), (%-)- , sZero, sOne, sTwo, sThree, Zero, One, Two, Three- , SZero, SOne, STwo, SThree- , lengthV, sLengthV, takeV, dropV, splitAtV, appendV- , foldrV, foldlV, singletonV, zipWithV, toList, allV- , mapV, headV, tailV, splitAtLess, takeVAtMost- , Leq(..), (:<<=), (:<=), LeqInstance(..)- , LeqTrueInstance(..), boolToPropLeq, boolToClassLeq- , propToClassLeq, propToBoolLeq- , leqRefl, leqSucc, Eql(..), eqlRefl, eqlSymm- , eqlTrans, plusZR, plusZL, eqPreservesS, plusAssociative- , sAndPlusOne, plusCommutative, minusCongEq, minusNilpotent- , eqSuccMinus, plusMinusEqL, plusMinusEqR, plusLeqL, plusLeqR- , zAbsorbsMinR, zAbsorbsMinL, minLeqL, minLeqR, plusSR- , leqRhs, leqLhs, leqTrans, minComm, leqAnitsymmetric- , maxZL, maxComm, maxZR, maxLeqL, maxLeqR, plusMonotone- , module Monomorphic+module Algebra.Internal ( toProxy, module Data.Proxy ) where import Data.Proxy-import Monomorphic toProxy :: a -> Proxy a toProxy _ = Proxy--data Nat = Z | S Nat-data Vector (a :: *) (n :: Nat) where- Nil :: Vector a Z- (:-) :: a -> Vector a n -> Vector a (S n)--infixr 5 :--deriving instance Show a => Show (Vector a n)--type family Min (n :: Nat) (m :: Nat) :: Nat-type instance Min Z Z = Z-type instance Min Z (S n) = Z-type instance Min (S m) Z = Z-type instance Min (S n) (S m) = S (Min n m)--type family Max (n :: Nat) (m :: Nat) :: Nat-type instance Max Z Z = Z-type instance Max Z (S n) = S n-type instance Max (S n) Z = S n-type instance Max (S n) (S m) = S (Max n m)---- | The smart constructor for @SZ@.-sZ :: SNat Z-sZ = case singInstance SZ of- SingInstance -> SZ---- | The smart constructor for @SS n@.-sS :: SNat n -> SNat (S n)-sS n = case singInstance n of- SingInstance -> SS n--type Zero = Z-type One = S Z-type Two = S (S Z)-type Three = S (S (S Z))--type SZero = SNat Zero-type SOne = SNat One-type STwo = SNat Two-type SThree = SNat Three--sZero :: SZero-sZero = SZ--sOne :: SOne-sOne = SS sZero--sTwo :: STwo-sTwo = SS sOne--sThree :: SThree-sThree = SS sTwo--sMin :: SNat n -> SNat m -> SNat (Min n m)-sMin SZ SZ = SZ-sMin (SS _) SZ = SZ-sMin SZ (SS _) = SZ-sMin (SS n) (SS m) = SS (sMin n m)--sMax :: SNat n -> SNat m -> SNat (Max n m)-sMax SZ SZ = SZ-sMax (SS n) SZ = SS n-sMax SZ (SS n) = SS n-sMax (SS n) (SS m) = SS (sMax n m)--class Sing (n :: Nat) where- sing :: SNat n--instance Sing Z where- sing = SZ--instance Sing n => Sing (S n) where- sing = SS sing--{-# SPECIALISE INLINE sing :: SNat Z #-}-{-# SPECIALISE INLINE sing :: SNat (S Z) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S Z)) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S (S Z))) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S (S (S Z)))) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S (S (S (S Z))))) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S (S (S (S (S Z)))))) #-}-{-# SPECIALISE INLINE sing :: SNat (S (S (S (S (S (S (S Z))))))) #-}--class (n :: Nat) :<= (m :: Nat)-instance Zero :<= n-instance (n :<= m) => S n :<= S m--data SNat (n :: Nat) where- SZ :: SNat Z- SS :: SNat n -> SNat (S n)--deriving instance Show (SNat n)--instance Monomorphicable SNat where- type MonomorphicRep SNat = Int- demote (Monomorphic sn) = toInt sn- promote n- | n < 0 = error "negative integer!"- | n == 0 = Monomorphic SZ- | otherwise = withPolymorhic (n - 1) $ \sn -> Monomorphic $ SS sn--instance Monomorphicable (Vector a) where- type MonomorphicRep (Vector a) = [a]- demote (Monomorphic vec) = toList vec- promote [] = Monomorphic Nil- promote (x:xs) =- case promote xs of- Monomorphic vec -> Monomorphic $ x :- vec--lengthV :: Vector a n -> Int-lengthV Nil = 0-lengthV (_ :- xs) = 1 + lengthV xs--type family (n :: Nat) :+: (m :: Nat) :: Nat-type instance Z :+: n = n-type instance S m :+: n = S (m :+: n)--(%+) :: SNat n -> SNat m -> SNat (n :+: m)-SZ %+ n = n-SS n %+ m = SS (n %+ m)--type family (n :: Nat) :-: (m :: Nat) :: Nat-type instance n :-: Z = n-type instance Z :-: m = Z-type instance S n :-: S m = n :-: m--(%-) :: (m :<<= n) ~ True => SNat n -> SNat m -> SNat (n :-: m)-n %- SZ = n-SS n %- SS m = n %- m-_ %- _ = error "impossible!"---- | Comparison function-type family (n :: Nat) :<<= (m :: Nat) :: Bool-type instance Z :<<= n = True-type instance S n :<<= Z = False-type instance S n :<<= S m = n :<<= m---- | Comparison witness via GADTs.-data Leq (n :: Nat) (m :: Nat) where- ZeroLeq :: SNat m -> Leq Zero m- SuccLeqSucc :: Leq n m -> Leq (S n) (S m)--data LeqInstance n m where- LeqInstance :: (n :<= m) => LeqInstance n m--boolToPropLeq :: (n :<<= m) ~ True => SNat n -> SNat m -> Leq n m-boolToPropLeq SZ m = ZeroLeq m-boolToPropLeq (SS n) (SS m) = SuccLeqSucc $ boolToPropLeq n m-boolToPropLeq _ _ = error "impossible happend!"--boolToClassLeq :: (n :<<= m) ~ True => SNat n -> SNat m -> LeqInstance n m-boolToClassLeq SZ _ = LeqInstance-boolToClassLeq (SS n) (SS m) =- case boolToClassLeq n m of- LeqInstance -> LeqInstance-boolToClassLeq _ _ = error "impossible!"--propToClassLeq :: Leq n m -> LeqInstance n m-propToClassLeq (ZeroLeq _) = LeqInstance-propToClassLeq (SuccLeqSucc leq) =- case propToClassLeq leq of- LeqInstance -> LeqInstance--appendV :: Vector a n -> Vector a m -> Vector a (n :+: m)-appendV (x :- xs) ys = x :- appendV xs ys-appendV Nil ys = ys--foldrV :: (a -> b -> b) -> b -> Vector a n -> b-foldrV _ b Nil = b-foldrV f a (x :- xs) = f x (foldrV f a xs)--foldlV :: (a -> b -> a) -> a -> Vector b n -> a-foldlV _ a Nil = a-foldlV f a (b :- bs) = foldlV f (f a b) bs--singletonV :: a -> Vector a (S Z)-singletonV = (:- Nil)--zipWithV :: (a -> b -> c) -> Vector a n -> Vector b n -> Vector c n-zipWithV _ Nil Nil = Nil-zipWithV f (x :- xs) (y :- ys) = f x y :- zipWithV f xs ys-zipWithV _ _ _ = error "cannot happen"--toList :: Vector a n -> [a]-toList = foldrV (:) []--instance (Eq a) => Eq (Vector a n) where- Nil == Nil = True- (x :- xs) == (y :- ys) = x == y && xs == ys- _ == _ = error "impossible!"--allV :: (a -> Bool) -> Vector a n -> Bool-allV p = foldrV ((&&) . p) False--dropV :: (n :<<= m) ~ True => SNat n -> Vector a m -> Vector a (m :-: n)-dropV n = snd . splitAtV n--takeV :: (n :<<= m) ~ True => SNat n -> Vector a m -> Vector a n-takeV n = fst . splitAtV n--takeVAtMost :: SNat n -> Vector a m -> Vector a (Min n m)-takeVAtMost = (fst .) . splitAtLess--splitAtLess :: SNat n -> Vector a m -> (Vector a (Min n m), Vector a (m :-: n))-splitAtLess SZ Nil = (Nil, Nil)-splitAtLess SZ (x :- xs) = (Nil, x :- xs)-splitAtLess (SS _) Nil = (Nil, Nil)-splitAtLess (SS n) (x :- xs) =- case splitAtLess n xs of- (ys, zs) -> (x :- ys, zs)--toInt :: SNat n -> Int-toInt SZ = 0-toInt (SS n) = 1 + toInt n--splitAtV :: (n :<<= m) ~ True- => SNat n -> Vector a m -> (Vector a n, Vector a (m :-: n))-splitAtV SZ xs = (Nil, xs)-splitAtV (SS n) (x :- xs) =- case splitAtV n xs of- (xs', ys') -> (x :- xs', ys')-splitAtV _ _ = error "could not happen!"--sLengthV :: Vector a n -> SNat n-sLengthV Nil = SZ-sLengthV (_ :- xs) = sOne %+ sLengthV xs-{-# RULES-"sLengthV/zero" forall (v :: Vector a Z). sLengthV v = sZ-"sLengthV/one" forall (v :: Vector a (S Z)). sLengthV v = sS sZ-"sLengthV/two" forall (v :: Vector a (S (S Z))). sLengthV v = sS (sS sZ)-"sLengthV/three" forall (v :: Vector a (S (S (S Z)))). sLengthV v = sS (sS (sS sZ))-"sLengthV/four" forall (v :: Vector a (S (S (S (S Z))))). sLengthV v = sS (sS (sS (sS sZ)))-"sLengthV/five" forall (v :: Vector a (S (S (S (S (S Z)))))). sLengthV v = sS (sS (sS (sS (sS sZ))))-"sLengthV/sing" forall (v :: Sing n => Vector a n). sLengthV (v :: Vector a n) = sing :: SNat n- #-}--mapV :: (a -> b) -> Vector a n -> Vector b n-mapV _ Nil = Nil-mapV f (x :- xs) = f x :- mapV f xs--headV :: Vector a (S n) -> a-headV (x :- _) = x--tailV :: Vector a (S n) -> Vector a n-tailV (_ :- xs) = xs--data SingInstance a where- SingInstance :: Sing a => SingInstance a--singInstance :: SNat n -> SingInstance n-singInstance SZ = SingInstance-singInstance (SS n) =- case singInstance n of- SingInstance -> SingInstance--data Reason x y where- Because :: SNat y -> Eql x y -> Reason x y--because :: SNat y -> Eql x y -> Reason x y-because = Because--infixl 4 ===, =~=-infix 5 `Because`-infix 5 `because`---(===) :: Eql x y -> Reason y z -> Eql x z-eq === (_ `Because` eq') = eqlTrans eq eq'--(=~=) :: Eql x y -> SNat y -> Eql x y-eq =~= _ = eq--start :: SNat a -> Eql a a-start = eqlRefl--definition, byDefinition :: Sing a => Eql a a-byDefinition = eqlRefl sing-definition = eqlRefl sing--admitted :: Reason x y-admitted = undefined-{-# WARNING admitted "There are some goals left yet unproven." #-}--infix 4 :=:-type a :=: b = Eql a b--cong' :: (SNat m -> SNat (f m)) -> a :=: b -> f a :=: f b-cong' _ Eql = Eql---leqRefl :: SNat n -> Leq n n-leqRefl SZ = ZeroLeq sZ-leqRefl (SS n) = SuccLeqSucc $ leqRefl n--leqSucc :: SNat n -> Leq n (S n)-leqSucc SZ = ZeroLeq sOne-leqSucc (SS n) = SuccLeqSucc $ leqSucc n--data Eql a b where- Eql :: Eql a a--eqlRefl :: SNat a -> Eql a a-eqlRefl _ = Eql--eqlSymm :: Eql a b -> Eql b a-eqlSymm Eql = Eql--eqlTrans :: Eql a b -> Eql b c -> Eql a c-eqlTrans Eql Eql = Eql--plusZR :: SNat n -> Eql (n :+: Z) n-plusZR SZ = Eql-plusZR (SS n) =- start (sS n %+ sZ)- =~= sS (n %+ sZ)- === sS n `because` cong' sS (plusZR n)--plusZL :: SNat n -> Eql (Z :+: n) n-plusZL _ = Eql--eqPreservesS :: Eql n m -> Eql (S n) (S m)-eqPreservesS Eql = Eql--plusAssociative :: SNat n -> SNat m -> SNat l- -> Eql (n :+: (m :+: l)) ((n :+: m) :+: l)-plusAssociative SZ _ _ = Eql-plusAssociative (SS n) m l =- start (sS n %+ (m %+ l))- =~= sS (n %+ (m %+ l))- === sS ((n %+ m) %+ l) `because` cong' sS (plusAssociative n m l)- =~= sS (n %+ m) %+ l- =~= (sS n %+ m) %+ l--sAndPlusOne :: SNat n -> Eql (S n) (n :+: One)-sAndPlusOne SZ = Eql-sAndPlusOne (SS n) =- start (sS (sS n))- === sS (n %+ sOne) `because` cong' sS (sAndPlusOne n)- =~= sS n %+ sOne--plusCongL :: SNat n -> m :=: m' -> n :+: m :=: n :+: m'-plusCongL _ Eql = Eql--plusCongR :: SNat n -> m :=: m' -> m :+: n :=: m' :+: n-plusCongR _ Eql = Eql--plusCommutative :: SNat n -> SNat m -> Eql (n :+: m) (m :+: n)-plusCommutative SZ SZ = Eql-plusCommutative SZ (SS m) =- start (sZ %+ sS m)- =~= sS m- === sS (m %+ sZ) `because` cong' sS (plusCommutative SZ m)- =~= sS m %+ sZ-plusCommutative (SS n) m =- start (sS n %+ m)- =~= sS (n %+ m)- === sS (m %+ n) `because` cong' sS (plusCommutative n m)- === (m %+ n) %+ sOne `because` sAndPlusOne (m %+ n)- === m %+ (n %+ sOne) `because` eqlSymm (plusAssociative m n sOne)- === m %+ sS n `because` plusCongL m (eqlSymm $ sAndPlusOne n)--minusCongEq :: Eql n m -> SNat l -> Eql (n :-: l) (m :-: l)-minusCongEq Eql _ = Eql--minusNilpotent :: SNat n -> Eql (n :-: n) Zero-minusNilpotent SZ = Eql-minusNilpotent (SS n) =- case minusNilpotent n of- Eql -> Eql--eqSuccMinus :: ((m :<<= n) ~ True)- => SNat n -> SNat m -> Eql (S n :-: m) (S (n :-: m))-eqSuccMinus _ SZ = Eql-eqSuccMinus (SS n) (SS m) = case eqSuccMinus n m of Eql -> Eql-eqSuccMinus _ _ = error "impossible!"--plusMinusEqL :: SNat n -> SNat m -> Eql ((n :+: m) :-: m) n-plusMinusEqL SZ m = minusNilpotent m-plusMinusEqL (SS n) m =- case propToBoolLeq (plusLeqR n m) of- LeqTrueInstance -> eqlTrans (eqSuccMinus (n %+ m) m) (eqPreservesS $ plusMinusEqL n m)--plusMinusEqR :: SNat n -> SNat m -> Eql ((m :+: n) :-: m) n-plusMinusEqR n m = eqlTrans (minusCongEq (plusCommutative m n) m) (plusMinusEqL n m)--data LeqTrueInstance a b where- LeqTrueInstance :: (a :<<= b) ~ True => LeqTrueInstance a b--propToBoolLeq :: Leq n m -> LeqTrueInstance n m-propToBoolLeq (ZeroLeq _) = LeqTrueInstance-propToBoolLeq (SuccLeqSucc leq) =- case propToBoolLeq leq of- LeqTrueInstance -> LeqTrueInstance---plusLeqL :: SNat n -> SNat m -> Leq n (n :+: m)-plusLeqL SZ m = case plusZR m of Eql -> ZeroLeq m-plusLeqL (SS n) m = SuccLeqSucc $ plusLeqL n m--plusLeqR :: SNat n -> SNat m -> Leq m (n :+: m)-plusLeqR n m =- case plusCommutative n m of- Eql -> plusLeqL m n--zAbsorbsMinR :: SNat n -> Eql (Min n Z) Z-zAbsorbsMinR SZ = Eql-zAbsorbsMinR (SS n) =- case zAbsorbsMinR n of- Eql -> Eql--zAbsorbsMinL :: SNat n -> Eql (Min Z n) Z-zAbsorbsMinL SZ = Eql-zAbsorbsMinL (SS n) = case zAbsorbsMinL n of Eql -> Eql--minLeqL :: SNat n -> SNat m -> Leq (Min n m) n-minLeqL SZ m = case zAbsorbsMinL m of Eql -> ZeroLeq sZ-minLeqL n SZ = case zAbsorbsMinR n of Eql -> ZeroLeq n-minLeqL (SS n) (SS m) = SuccLeqSucc (minLeqL n m)--minLeqR :: SNat n -> SNat m -> Leq (Min n m) m-minLeqR n m = case minComm n m of Eql -> minLeqL m n--leqRhs :: Leq n m -> SNat m-leqRhs (ZeroLeq m) = m-leqRhs (SuccLeqSucc leq) = SS $ leqRhs leq--leqLhs :: Leq n m -> SNat n-leqLhs (ZeroLeq _) = SZ-leqLhs (SuccLeqSucc leq) = SS $ leqLhs leq--leqTrans :: Leq n m -> Leq m l -> Leq n l-leqTrans (ZeroLeq _) leq = ZeroLeq $ leqRhs leq-leqTrans (SuccLeqSucc nLeqm) (SuccLeqSucc mLeql) = SuccLeqSucc $ leqTrans nLeqm mLeql-leqTrans _ _ = error "impossible!"--minComm :: SNat n -> SNat m -> Eql (Min n m) (Min m n)-minComm SZ SZ = Eql-minComm SZ (SS _) = Eql-minComm (SS _) SZ = Eql-minComm (SS n) (SS m) = case minComm n m of Eql -> Eql--leqAnitsymmetric :: Leq n m -> Leq m n -> Eql n m-leqAnitsymmetric (ZeroLeq _) (ZeroLeq _) = Eql-leqAnitsymmetric (SuccLeqSucc leq1) (SuccLeqSucc leq2) = eqPreservesS $ leqAnitsymmetric leq1 leq2-leqAnitsymmetric _ _ = error "impossible"--maxZL :: SNat n -> Eql (Max Z n) n-maxZL SZ = Eql-maxZL (SS _) = Eql--maxComm :: SNat n -> SNat m -> Eql (Max n m) (Max m n)-maxComm SZ SZ = Eql-maxComm SZ (SS _) = Eql-maxComm (SS _) SZ = Eql-maxComm (SS n) (SS m) = case maxComm n m of Eql -> Eql--maxZR :: SNat n -> Eql (Max n Z) n-maxZR n = eqlTrans (maxComm n sZ) (maxZL n)--maxLeqL :: SNat n -> SNat m -> Leq n (Max n m)-maxLeqL SZ m = ZeroLeq (sMax sZ m)-maxLeqL n SZ = case maxZR n of- Eql -> leqRefl n-maxLeqL (SS n) (SS m) = SuccLeqSucc $ maxLeqL n m--maxLeqR :: SNat n -> SNat m -> Leq m (Max n m)-maxLeqR n m = case maxComm n m of- Eql -> maxLeqL m n--plusSR :: SNat n -> SNat m -> Eql (S (n :+: m)) (n :+: S m)-plusSR n m =- start (sS (n %+ m))- === (n %+ m) %+ sOne `because` sAndPlusOne (n %+ m)- === n %+ (m %+ sOne) `because` eqlSymm (plusAssociative n m sOne)- === n %+ sS m `because` plusCongL n (eqlSymm $ sAndPlusOne m)--plusMonotone :: Leq n m -> Leq l k -> Leq (n :+: l) (m :+: k)-plusMonotone (ZeroLeq m) (ZeroLeq k) = ZeroLeq (m %+ k)-plusMonotone (ZeroLeq m) (SuccLeqSucc leq) =- case plusSR m (leqRhs leq) of- Eql -> SuccLeqSucc $ plusMonotone (ZeroLeq m) leq-plusMonotone (SuccLeqSucc leq) leq' = SuccLeqSucc $ plusMonotone leq leq'---- (m + S n) - m = S (m + n) - m
Algebra/Ring/Noetherian.hs view
@@ -15,6 +15,9 @@ import Prelude hiding (negate, subtract, (*), (+), (-)) import qualified Prelude as P+import qualified Data.Vector.Sized as V+import Data.Vector.Sized (Vector(..))+import Data.Type.Natural class (Commutative r, Ring r) => NoetherianRing r where @@ -79,7 +82,7 @@ (*) = (P.*) pow1p r n = r ^^ (n P.+ 1) -data Ideal r = forall n. Ideal (Vector r n)+data Ideal r = forall n. Ideal (V.Vector r n) instance Eq r => Eq (Ideal r) where (==) = (==) `on` generators@@ -97,16 +100,16 @@ toIdeal = foldr addToIdeal (Ideal Nil) appendIdeal :: Ideal r -> Ideal r -> Ideal r-appendIdeal (Ideal is) (Ideal js) = Ideal (is `appendV` js)+appendIdeal (Ideal is) (Ideal js) = Ideal (is `V.append` js) generators :: Ideal r -> [r]-generators (Ideal is) = toList is+generators (Ideal is) = V.toList is filterIdeal :: NoetherianRing r => (r -> Bool) -> Ideal r -> Ideal r-filterIdeal p (Ideal i) = foldrV (\h -> if p h then addToIdeal h else id) (toIdeal []) i+filterIdeal p (Ideal i) = V.foldr (\h -> if p h then addToIdeal h else id) (toIdeal []) i principalIdeal :: r -> Ideal r-principalIdeal = Ideal . singletonV+principalIdeal = Ideal . V.singleton mapIdeal :: (r -> r') -> Ideal r -> Ideal r'-mapIdeal fun (Ideal xs) = Ideal $ mapV fun xs+mapIdeal fun (Ideal xs) = Ideal $ V.map fun xs
Algebra/Ring/Polynomial.hs view
@@ -1,8 +1,9 @@+{-# LANGUAGE ViewPatterns #-} {-# LANGUAGE ConstraintKinds, DataKinds, FlexibleContexts, FlexibleInstances #-} {-# LANGUAGE GADTs, GeneralizedNewtypeDeriving, LiberalTypeSynonyms #-} {-# LANGUAGE MultiParamTypeClasses, OverlappingInstances, PolyKinds #-} {-# LANGUAGE RankNTypes, ScopedTypeVariables, StandaloneDeriving #-}-{-# LANGUAGE TypeFamilies, TypeOperators, UndecidableInstances, ViewPatterns #-}+{-# LANGUAGE TypeFamilies, TypeOperators, UndecidableInstances #-} {-# OPTIONS_GHC -fno-warn-orphans -fno-warn-type-defaults #-} module Algebra.Ring.Polynomial ( Polynomial, Monomial, MonomialOrder, EliminationType, EliminationOrder@@ -23,6 +24,7 @@ import Algebra.Ring.Noetherian import Control.Arrow import Control.Lens+import Data.Function import Data.List (intercalate) import Data.Map (Map) import qualified Data.Map as M@@ -31,10 +33,14 @@ import Data.Ord import Data.Proxy import Data.Ratio-import Numeric.Algebra hiding (Order(..), sum)-import Prelude hiding (lex, (*), (+), (-), (^), (^^), recip, negate)+import Data.Type.Monomorphic+import Data.Type.Natural hiding (max, one, promote, zero)+import Data.Vector.Sized (Vector (..))+import qualified Data.Vector.Sized as V+import Numeric.Algebra hiding (Order (..), sum)+import Prelude hiding (lex, negate, recip, (*), (+),+ (-), (^), (^^)) import qualified Prelude as P-import Data.Function -- | N-ary Monomial. IntMap contains degrees for each x_i. type Monomial (n :: Nat) = Vector Int n@@ -71,7 +77,7 @@ -- | apply monomial ordering to monomorphic monomials. cmpMonomial' :: IsMonomialOrder ord => Proxy ord -> [Int] -> [Int] -> Ordering cmpMonomial' pxy xs ys =- withPolymorhic (max (length xs) (length ys)) $ \n ->+ withPolymorhic (P.max (length xs) (length ys)) $ \n -> cmpMonomial pxy (fromList n xs) (fromList n ys) -- | Monomial order (of degree n). This should satisfy following laws:@@ -81,7 +87,7 @@ type MonomialOrder = forall n. Monomial n -> Monomial n -> Ordering totalDegree :: Monomial n -> Int-totalDegree = foldlV (+) 0+totalDegree = V.foldl (+) 0 {-# INLINE totalDegree #-} -- | Lexicographical order. This *is* a monomial order.@@ -156,11 +162,11 @@ data ProductOrder (n :: Nat) (a :: *) (b :: *) where ProductOrder :: SNat n -> ord -> ord' -> ProductOrder n ord ord' -productOrder :: forall ord ord' n m. (IsOrder ord, IsOrder ord', Sing n)+productOrder :: forall ord ord' n m. (IsOrder ord, IsOrder ord', SingRep n) => Proxy (ProductOrder n ord ord') -> Monomial m -> Monomial m -> Ordering productOrder _ m m' = case sing :: SNat n of- n -> case (splitAtLess n m, splitAtLess n m') of+ n -> case (V.splitAtMost n m, V.splitAtMost n m') of ((xs, xs'), (ys, ys')) -> cmpMonomial (Proxy :: Proxy ord) xs ys <> cmpMonomial (Proxy :: Proxy ord') xs' ys' productOrder' :: forall n ord ord' m.(IsOrder ord, IsOrder ord')@@ -185,9 +191,9 @@ instance ToWeightVector '[] where calcOrderWeight Proxy' _ = 0 -instance (Sing n, ToWeightVector ns) => ToWeightVector (n ': ns) where+instance (SingRep n, ToWeightVector ns) => ToWeightVector (n ': ns) where calcOrderWeight Proxy' Nil = 0- calcOrderWeight Proxy' (x :- xs) = x * toInt (sing :: SNat n) + calcOrderWeight (Proxy' :: Proxy' ns) xs+ calcOrderWeight Proxy' (x :- xs) = x * sNatToInt (sing :: SNat n) + calcOrderWeight (Proxy' :: Proxy' ns) xs weightOrder :: forall ns ord m. (ToWeightVector ns, IsOrder ord) => Proxy (WeightOrder ns ord) -> Monomial m -> Monomial m -> Ordering@@ -197,7 +203,7 @@ instance (ToWeightVector ws, IsOrder ord) => IsOrder (WeightOrder ws ord) where cmpMonomial p = weightOrder p -instance (IsOrder ord, IsOrder ord', Sing n) => IsOrder (ProductOrder n ord ord') where+instance (IsOrder ord, IsOrder ord', SingRep n) => IsOrder (ProductOrder n ord ord') where cmpMonomial p = productOrder p -- They're all total orderings.@@ -221,14 +227,14 @@ instance IsMonomialOrder Grlex instance IsMonomialOrder Grevlex instance IsMonomialOrder Lex-instance (Sing n, IsMonomialOrder o, IsMonomialOrder o') => IsMonomialOrder (ProductOrder n o o')+instance (SingRep n, IsMonomialOrder o, IsMonomialOrder o') => IsMonomialOrder (ProductOrder n o o') instance (ToWeightVector ws, IsMonomialOrder ord) => IsMonomialOrder (WeightOrder ws ord) -- | Monomial order which can be use to calculate n-th elimination ideal. -- This should judge it as bigger that contains variables to eliminate.-class (IsMonomialOrder ord, Sing n) => EliminationType n ord-instance Sing n => EliminationType n Lex-instance (Sing n, IsMonomialOrder ord, IsMonomialOrder ord') => EliminationType n (ProductOrder n ord ord')+class (IsMonomialOrder ord, SingRep n) => EliminationType n ord+instance SingRep n => EliminationType n Lex+instance (SingRep n, IsMonomialOrder ord, IsMonomialOrder ord') => EliminationType n (ProductOrder n ord ord') instance (IsMonomialOrder ord) => EliminationType Z (WeightOrder '[] ord) instance (IsMonomialOrder ord, ToWeightVector ns, EliminationType n (WeightOrder ns ord)) => EliminationType (S n) (WeightOrder (One ': ns) ord)@@ -243,14 +249,14 @@ weightedEliminationOrder :: SNat n -> WeightedEliminationOrder n Grevlex weightedEliminationOrder n = WEOrder n (Proxy :: Proxy Grevlex) -type family EWeight n :: [Nat]+type family EWeight (n :: Nat) :: [Nat] type instance EWeight Z = '[] type instance EWeight (S n) = One ': EWeight n data WeightedEliminationOrder (n :: Nat) (ord :: *) where WEOrder :: SNat n -> Proxy ord -> WeightedEliminationOrder n ord -instance (Sing n, IsMonomialOrder ord) => IsOrder (WeightedEliminationOrder n ord) where+instance (SingRep n, IsMonomialOrder ord) => IsOrder (WeightedEliminationOrder n ord) where cmpMonomial Proxy m m' = comparing (calc (sing :: SNat n)) m m' <> cmpMonomial (Proxy :: Proxy ord) m m' where calc :: SNat l -> Vector Int m -> Int@@ -258,9 +264,9 @@ calc SZ _ = 0 calc (SS l) (x :- xs)= x + calc l xs -instance (Sing n, IsMonomialOrder ord) => IsMonomialOrder (WeightedEliminationOrder n ord)+instance (SingRep n, IsMonomialOrder ord) => IsMonomialOrder (WeightedEliminationOrder n ord) -instance (Sing n, IsMonomialOrder ord) => EliminationType n (WeightedEliminationOrder n ord) where+instance (SingRep n, IsMonomialOrder ord) => EliminationType n (WeightedEliminationOrder n ord) where -- | Special ordering for ordered-monomials. instance (Eq (Monomial n), IsOrder name) => Ord (OrderedMonomial name n) where@@ -275,7 +281,7 @@ type Polynomial r = OrderedPolynomial r Grevlex -- | Type-level constraint to check whether it forms polynomial ring or not.-type IsPolynomial r n = (NoetherianRing r, Sing n, Eq r)+type IsPolynomial r n = (NoetherianRing r, SingRep n, Eq r) -- | coefficient for a degree. coeff :: (IsOrder order, IsPolynomial r n) => Monomial n -> OrderedPolynomial r order n -> r@@ -285,18 +291,18 @@ (OrderedPolynomial r order n) (OrderedPolynomial q order' m) where wrapped = iso Polynomial terms -castMonomial :: (IsOrder o, IsOrder o', Sing m, n :<= m) => OrderedMonomial o n -> OrderedMonomial o' m-castMonomial = unwrapped %~ fromList sing . toList+castMonomial :: (IsOrder o, IsOrder o', SingRep m, n :<= m) => OrderedMonomial o n -> OrderedMonomial o' m+castMonomial = unwrapped %~ fromList sing . V.toList scastMonomial :: (n :<= m) => SNat m -> OrderedMonomial o n -> OrderedMonomial o m-scastMonomial snat = unwrapped %~ fromList snat . toList+scastMonomial snat = unwrapped %~ fromList snat . V.toList -castPolynomial :: (IsPolynomial r n, IsPolynomial r m, Sing m, IsOrder o, IsOrder o', n :<= m)+castPolynomial :: (IsPolynomial r n, IsPolynomial r m, SingRep m, IsOrder o, IsOrder o', n :<= m) => OrderedPolynomial r o n -> OrderedPolynomial r o' m castPolynomial = unwrapped %~ M.mapKeys castMonomial -scastPolynomial :: (IsOrder o, IsOrder o', IsPolynomial r n, IsPolynomial r m, n :<= m, Sing m)+scastPolynomial :: (IsOrder o, IsOrder o', IsPolynomial r n, IsPolynomial r m, n :<= m, SingRep m) => SNat m -> OrderedPolynomial r o n -> OrderedPolynomial r o' m scastPolynomial _ = castPolynomial @@ -332,7 +338,7 @@ one = injectCoeff one instance (IsOrder order, IsPolynomial r n) => Multiplicative (OrderedPolynomial r order n) where Polynomial (M.toList -> d1) * Polynomial (M.toList -> d2) =- let dic = [ (OrderedMonomial $ zipWithV (+) a b, r * r') | (getMonomial -> a, r) <- d1, (getMonomial -> b, r') <- d2 ]+ let dic = [ (OrderedMonomial $ V.zipWithSame (+) a b, r * r') | (getMonomial -> a, r) <- d1, (getMonomial -> b, r') <- d2 ] in normalize $ Polynomial $ M.fromListWith (+) dic instance (IsOrder order, IsPolynomial r n) => Semiring (OrderedPolynomial r order n) where instance (IsOrder order, IsPolynomial r n) => Commutative (OrderedPolynomial r order n) where@@ -341,10 +347,10 @@ instance (Eq r, IsPolynomial r n, IsOrder order, Show r) => Show (OrderedPolynomial r order n) where show = showPolynomialWithVars [(n, "X_"++ show n) | n <- [1..]] -instance (Sing n, IsOrder order) => Show (OrderedPolynomial Rational order n) where+instance (SingRep n, IsOrder order) => Show (OrderedPolynomial Rational order n) where show = showPolynomialWith [(n, "X_"++ show n) | n <- [1..]] showRational -showPolynomialWithVars :: (Eq a, Show a, Sing n, NoetherianRing a, IsOrder ordering)+showPolynomialWithVars :: (Eq a, Show a, SingRep n, NoetherianRing a, IsOrder ordering) => [(Int, String)] -> OrderedPolynomial a ordering n -> String showPolynomialWithVars dic p0@(Polynomial d) | p0 == zero = "0"@@ -354,11 +360,11 @@ | c == zero = Nothing | otherwise = let cstr = if (c == zero - one)- then if any (/= zero) (toList deg) then "-" else "-1"+ then if any (/= zero) (V.toList deg) then "-" else "-1" else if (c /= one || isConstantMonomial deg) then show c ++ " " else ""- in Just $ cstr ++ unwords (mapMaybe showDeg (zip [1..] $ toList deg))+ in Just $ cstr ++ unwords (mapMaybe showDeg (zip [1..] $ V.toList deg)) showDeg (n, p) | p == 0 = Nothing | p == 1 = Just $ showVar n | otherwise = Just $ showVar n ++ "^" ++ show p@@ -375,7 +381,7 @@ formatRat q | denominator q == 1 = show $ numerator q | otherwise = show (numerator q) ++ "/" ++ show (denominator q) ++ " " -showPolynomialWith :: (Eq a, Show a, Sing n, NoetherianRing a, IsOrder ordering)+showPolynomialWith :: (Eq a, Show a, SingRep n, NoetherianRing a, IsOrder ordering) => [(Int, String)] -> (a -> Coefficient) -> OrderedPolynomial a ordering n -> String showPolynomialWith vDic showCoeff p0@(Polynomial d) | p0 == zero = "0"@@ -396,14 +402,14 @@ | c == one = Positive "" | c == negate one = Negative "" | otherwise = cKind- in Just $ (cff, unwords (mapMaybe showDeg (zip [1..] $ toList deg)))+ in Just $ (cff, unwords (mapMaybe showDeg (zip [1..] $ V.toList deg))) showDeg (n, p) | p == 0 = Nothing | p == 1 = Just $ showVar n | otherwise = Just $ showVar n ++ "^" ++ show p showVar n = fromMaybe ("X_" ++ show n) $ lookup n vDic isConstantMonomial :: (Eq a, Num a) => Vector a n -> Bool-isConstantMonomial v = all (== 0) $ toList v+isConstantMonomial v = all (== 0) $ V.toList v -- | We provide Num instance to use trivial injection R into R[X]. -- Do not use signum or abs.@@ -415,16 +421,16 @@ abs = id negate = ((P.negate 1 :: Integer) .*) -varX :: (NoetherianRing r, Sing n, One :<= n) => OrderedPolynomial r order n+varX :: (NoetherianRing r, SingRep n, One :<= n) => OrderedPolynomial r order n varX = Polynomial $ M.singleton (OrderedMonomial $ fromList sing [1]) one -var :: (NoetherianRing r, Sing m, S n :<= m) => SNat (S n) -> OrderedPolynomial r order m+var :: (NoetherianRing r, SingRep m, S n :<= m) => SNat (S n) -> OrderedPolynomial r order m var vIndex = Polynomial $ M.singleton (OrderedMonomial $ fromList sing (buildIndex vIndex)) one toPolynomial :: (IsOrder order, IsPolynomial r n) => (r, Monomial n) -> OrderedPolynomial r order n toPolynomial (c, deg) = Polynomial $ M.singleton (OrderedMonomial deg) c -polynomial :: (Sing n, Eq r, NoetherianRing r, IsOrder order) => Map (OrderedMonomial order n) r -> OrderedPolynomial r order n+polynomial :: (SingRep n, Eq r, NoetherianRing r, IsOrder order) => Map (OrderedMonomial order n) r -> OrderedPolynomial r order n polynomial dic = normalize $ Polynomial dic buildIndex :: SNat (S n) -> [Int]@@ -449,15 +455,15 @@ leadingCoeff = fst . leadingTerm divs :: Monomial n -> Monomial n -> Bool-xs `divs` ys = and $ toList $ zipWithV (<=) xs ys+xs `divs` ys = and $ V.toList $ V.zipWith (<=) xs ys tryDiv :: Field r => (r, Monomial n) -> (r, Monomial n) -> (r, Monomial n) tryDiv (a, f) (b, g)- | g `divs` f = (a * recip b, zipWithV (-) f g)+ | g `divs` f = (a * recip b, V.zipWithSame (-) f g) | otherwise = error "cannot divide." lcmMonomial :: Monomial n -> Monomial n -> Monomial n-lcmMonomial = zipWithV max+lcmMonomial = V.zipWithSame max sPolynomial :: (IsPolynomial k n, Field k, IsOrder order) => OrderedPolynomial k order n@@ -466,11 +472,11 @@ let h = (one, lcmMonomial (leadingMonomial f) (leadingMonomial g)) in toPolynomial (h `tryDiv` leadingTerm f) * f - toPolynomial (h `tryDiv` leadingTerm g) * g -changeOrder :: (Eq (Monomial n), IsOrder o, IsOrder o', Sing n)+changeOrder :: (Eq (Monomial n), IsOrder o, IsOrder o', SingRep n) => o' -> OrderedPolynomial k o n -> OrderedPolynomial k o' n changeOrder _ = unwrapped %~ M.mapKeys (OrderedMonomial . getMonomial) -changeOrderProxy :: (Eq (Monomial n), IsOrder o, IsOrder o', Sing n)+changeOrderProxy :: (Eq (Monomial n), IsOrder o, IsOrder o', SingRep n) => Proxy o' -> OrderedPolynomial k o n -> OrderedPolynomial k o' n changeOrderProxy _ = unwrapped %~ M.mapKeys (OrderedMonomial . getMonomial) @@ -488,17 +494,17 @@ => SNat k -> OrderedPolynomial r ord n -> OrderedPolynomial r ord (k :+: n) shiftR k = case singInstance k of- SingInstance -> transformMonomial (appendV (fromList k []))+ SingInstance -> transformMonomial (V.append (fromList k [])) genVars :: forall k o n. (IsPolynomial k (S n), IsOrder o) => SNat (S n) -> [OrderedPolynomial k o (S n)] genVars sn =- let n = toInt sn+ let n = sNatToInt sn seed = cycle $ 1 : replicate (n - 1) 0 in map (\m -> Polynomial $ M.singleton (OrderedMonomial $ fromList sn $ take n (drop (n-m) seed)) one) [0..n-1] sDegree :: OrderedPolynomial k ord n -> SNat n-sDegree (Polynomial dic) = sLengthV $ getMonomial $ fst $ M.findMin dic+sDegree (Polynomial dic) = V.sLength $ getMonomial $ fst $ M.findMin dic {-# RULES "sDegree/zero" forall (v :: OrderedPolynomial k ord Z). sDegree v = SZ "sDegree/one" forall (v :: OrderedPolynomial k ord (S Z)). sDegree v = SS SZ@@ -506,5 +512,5 @@ "sDegree/three" forall (v :: OrderedPolynomial k ord (S (S (S Z)))). sDegree v = SS (SS (sS SZ)) "sDegree/four" forall (v :: OrderedPolynomial k ord (S (S (S (S Z))))). sDegree v = SS (SS (SS (SS SZ))) "sDegree/five" forall (v :: OrderedPolynomial k ord (S (S (S (S (S Z)))))). sDegree v = SS (SS (SS (SS (SS SZ))))-"sDegree/sing" forall (v :: Sing n => OrderedPolynomial k ord n). sDegree (v :: OrderedPolynomial k ord n) = sing :: SNat n+"sDegree/sing" forall (v :: SingRep n => OrderedPolynomial k ord n). sDegree (v :: OrderedPolynomial k ord n) = sing :: SNat n #-}
Algebra/Ring/Polynomial/Monomorphic.hs view
@@ -10,8 +10,12 @@ import Data.List import qualified Data.Map as M import Data.Maybe+import Data.Type.Natural hiding (one, zero, promote)+import Data.Type.Monomorphic import qualified Numeric.Algebra as NA import Data.Ratio+import Data.Singletons hiding (promote)+import qualified Data.Vector.Sized as V data Variable = Variable { varName :: Char , varIndex :: Maybe Int@@ -42,12 +46,6 @@ normalizeMonom :: Monomial -> Monomial normalizeMonom = M.filter (/= 0) -instance (Eq r, NoetherianRing r) => NA.LeftModule r (Polynomial r) where- c .* Polynomial d = normalize $ Polynomial $ fmap (c NA.*) d--instance (Eq r, NoetherianRing r) => NA.RightModule r (Polynomial r) where- (*.) = flip (NA..*)- instance (Eq r, NoetherianRing r) => NoetherianRing (Polynomial r) instance (Eq r, NoetherianRing r) => NA.Commutative (Polynomial r) instance (Eq r, NoetherianRing r) => NA.Multiplicative (Polynomial r) where@@ -82,7 +80,7 @@ encodeMonomList :: [Variable] -> Monomial -> [Int] encodeMonomList vars mono = map (maybe 0 fromInteger . flip M.lookup mono) vars -encodeMonomial :: [Variable] -> Monomial -> Monomorphic (Vector Int)+encodeMonomial :: [Variable] -> Monomial -> Monomorphic (V.Vector Int) encodeMonomial vars mono = promote $ encodeMonomList vars mono encodePolynomial :: (Monomorphicable (Poly.Polynomial r))@@ -95,7 +93,7 @@ , dimension = promote $ length $ buildVarsList p } -data PolynomialSetting r = PolySetting { dimension :: Monomorphic SNat+data PolynomialSetting r = PolySetting { dimension :: Monomorphic (Sing :: Nat -> *) , polyn :: Polynomial r } @@ -144,7 +142,7 @@ => Monomorphicable (Ideal :.: Poly.OrderedPolynomial r ord) where type MonomorphicRep (Ideal :.: Poly.OrderedPolynomial r ord) = [Polynomial r] promote = uniformlyPromote- demote (Monomorphic (Comp (Ideal v))) = map (polyn . demote . Monomorphic) $ toList v+ demote (Monomorphic (Comp (Ideal v))) = map (polyn . demote . Monomorphic) $ V.toList v promoteList :: (Eq r, NoetherianRing r, Poly.IsMonomialOrder ord) => [Polynomial r] -> Monomorphic ([] :.: Poly.OrderedPolynomial r ord)@@ -188,7 +186,7 @@ case singInstance (Poly.sDegree f') of SingInstance -> Poly.showPolynomialWithVars dic f' where- dic = zip [1..] $ map show $ buildVarsList f+ dic = zip [1 :: Int ..] $ map show $ buildVarsList f showRatPolynomial :: (Integral a, Show a) => Polynomial (Ratio a) -> String showRatPolynomial f =@@ -197,7 +195,7 @@ case singInstance (Poly.sDegree f') of SingInstance -> Poly.showPolynomialWith dic Poly.showRational f' where- dic = zip [1..] $ map show $ buildVarsList f+ dic = zip [1 :: Int ..] $ map show $ buildVarsList f injectVar :: NA.Unital r => Variable -> Polynomial r injectVar var = Polynomial $ M.singleton (M.singleton var 1) NA.one
Algebra/Ring/Polynomial/Parser.hs view
@@ -29,7 +29,7 @@ / term term :: Polynomial Rational- = number space* monoms { $1 NA..* $3 }+ = number space* monoms { injectCoeff $1 * $3 } / number { injectCoeff $1 } / monoms
Monomorphic.hs view
@@ -1,61 +1,5 @@ {-# LANGUAGE DataKinds, ExistentialQuantification, FlexibleContexts, GADTs #-} {-# LANGUAGE PolyKinds, RankNTypes, TypeFamilies, TypeOperators #-} {-# LANGUAGE UndecidableInstances, CPP #-}-#if __GLASGOW_HASKELL__ >= 761 module Monomorphic (module Data.Type.Monomorphic) where import Data.Type.Monomorphic-#else-module Monomorphic ( Monomorphic (..), Monomorphicable(..)- , demote', demoteComposed, monomorphicCompose- , withPolymorhic, liftPoly, viaPoly, (:.:)(..)- ) where-import Control.Arrow--newtype (:.:) f g a = Comp (f (g a))---- | A wrapper type for polymophic types.-data Monomorphic k = forall a. Monomorphic (k a)---- | A types which have the monomorphic representation.-class Monomorphicable k where- -- | Monomorphic representation- type MonomorphicRep k :: *- -- | Promote the monomorphic value to the polymophic one.- promote :: MonomorphicRep k -> Monomorphic k- -- | Demote the polymorphic value to the monomorphic representation.- demote :: Monomorphic k -> MonomorphicRep k---- | Convinience function to demote polymorphic types into monomorphic one directly.-demote' :: Monomorphicable k => k a -> MonomorphicRep k-demote' = demote . Monomorphic---- | Demote polymorphic nested types directly into monomorphic representation.-demoteComposed :: Monomorphicable (f :.: g) => f (g a) -> MonomorphicRep (f :.: g)-demoteComposed = demote . Monomorphic . Comp--monomorphicCompose :: f (g a) -> Monomorphic (f :.: g)-monomorphicCompose = Monomorphic . Comp---- | Apply dependently-typed function to the monomorphic representation.-withPolymorhic :: Monomorphicable k- => MonomorphicRep k -> (forall a. k a -> b) -> b-withPolymorhic k trans =- case promote k of- Monomorphic a -> trans a---- | Flipped version of 'withPolymorhic'.-liftPoly :: Monomorphicable k- => (forall a. k a -> b) -> MonomorphicRep k -> b-liftPoly = flip withPolymorhic---- | Demote the function between polymorphic types into the one between monomorphic one.-viaPoly :: (Monomorphicable k, Monomorphicable k')- => (forall x y. k x -> k' y) -> MonomorphicRep k -> MonomorphicRep k'-viaPoly f a = demote $ Monomorphic $ liftPoly f a--instance (Show (MonomorphicRep k), Monomorphicable k) => Show (Monomorphic k) where- showsPrec d x = showString "Polymorphic " . showsPrec (d + 1) (demote x)--instance (Read (MonomorphicRep k), Monomorphicable k) => Read (Monomorphic k) where- readsPrec i = map (first promote) . readsPrec i-#endif
computational-algebra.cabal view
@@ -2,7 +2,7 @@ -- further documentation, see http://haskell.org/cabal/users-guide/ name: computational-algebra-version: 0.1.3.6+version: 0.2.0.0 synopsis: Well-kinded computational algebra library, currently supporting Groebner basis. description: Dependently-typed computational algebra libray for Groebner basis. homepage: https://github.com/konn/computational-algebra@@ -28,13 +28,17 @@ , Algebra.Ring.Polynomial.Parser , Algebra.Internal other-modules: Monomorphic- build-depends: base >= 2.0 && < 5- , algebra == 3.*- , tagged >= 0.4 && < 1- , lens == 3.*- , containers >= 0.4 && < 0.6- , peggy == 0.3.*- , monad-loops == 0.3.*- , heaps == 0.2.*+ build-depends: base >= 2.0 && < 5+ , algebra == 3.*+ , tagged >= 0.4 && < 1+ , lens == 3.*+ , containers >= 0.4 && < 0.6+ , peggy == 0.3.*+ , monad-loops >= 0.3 && <0.5+ , heaps == 0.2.*+ , type-natural >= 0.0.2.1+ , sized-vector == 0.0.*+ , singletons >= 0.8+ , equational-reasoning == 0.0.* if impl(ghc >= 7.6.1) build-depends: monomorphic == 0.0.*
examples/polymorphic.hs view
@@ -6,6 +6,7 @@ import Algebra.Ring.Noetherian import Algebra.Ring.Polynomial import Data.Ratio+import Data.Type.Natural hiding (one, zero) import Numeric.Algebra import Prelude hiding (Fractional (..), Integral (..), Num (..), (^), (^^))