{-# LANGUAGE NoMonomorphismRestriction, PartialTypeSignatures #-}
{-# LANGUAGE PatternSynonyms, ScopedTypeVariables #-}
module Algebra.Algorithms.Groebner.F4
(f4', f4, f4WithStrategy', f4WithStrategy, normalStrategy) where
import Algebra.Matrix.DataMatrix (DMatrix)
import Algebra.Matrix.Generic
import Algebra.Prelude.Core hiding (Min)
import Control.Lens hiding ((.=))
import Control.Monad.Loops
import Control.Monad.ST
import Control.Monad.ST.Combinators
import qualified Data.Foldable as F
import qualified Data.Heap as H
import Data.Monoid (First (..))
import qualified Data.Set as S
import qualified Data.Vector as V
import qualified Data.Vector.Generic as GV
import qualified Data.Vector.Mutable as MV
type Strategy f w = f -> f -> w
normalStrategy :: IsOrderedPolynomial poly => poly -> poly -> Int
normalStrategy :: poly -> poly -> Int
normalStrategy = poly -> poly -> Int
forall poly. IsOrderedPolynomial poly => poly -> poly -> Int
degPair
degPair :: (IsOrderedPolynomial poly) => poly -> poly -> Int
degPair :: poly -> poly -> Int
degPair poly
f poly
g = OrderedMonomial (MOrder poly) (Arity poly) -> Int
forall k (ord :: k) (n :: Nat). OrderedMonomial ord n -> Int
totalDegree (OrderedMonomial (MOrder poly) (Arity poly) -> Int)
-> OrderedMonomial (MOrder poly) (Arity poly) -> Int
forall a b. (a -> b) -> a -> b
$ poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
lcmLeadingMonomial poly
f poly
g
lcmLeadingMonomial :: (IsOrderedPolynomial poly)
=> poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
lcmLeadingMonomial :: poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
lcmLeadingMonomial poly
f poly
g = OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
forall k (n :: Nat) (ord :: k).
KnownNat n =>
OrderedMonomial ord n
-> OrderedMonomial ord n -> OrderedMonomial ord n
lcmMonomial (poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial poly
f) (poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial poly
g)
viewMins :: Ord w => H.Heap w -> Maybe (H.Heap w, H.Heap w)
viewMins :: Heap w -> Maybe (Heap w, Heap w)
viewMins Heap w
h = do
(w
a, Heap w
_) <- Heap w -> Maybe (w, Heap w)
forall a. Heap a -> Maybe (a, Heap a)
H.viewMin Heap w
h
(Heap w, Heap w) -> Maybe (Heap w, Heap w)
forall (m :: * -> *) a. Monad m => a -> m a
return ((Heap w, Heap w) -> Maybe (Heap w, Heap w))
-> (Heap w, Heap w) -> Maybe (Heap w, Heap w)
forall a b. (a -> b) -> a -> b
$ (w -> Bool) -> Heap w -> (Heap w, Heap w)
forall a. (a -> Bool) -> Heap a -> (Heap a, Heap a)
H.span (w -> w -> Bool
forall a. Eq a => a -> a -> Bool
== w
a) Heap w
h
f4 :: (Field (Coefficient poly), Num (Coefficient poly),
IsOrderedPolynomial poly, Normed (Coefficient poly))
=> Ideal poly -> [poly]
f4 :: Ideal poly -> [poly]
f4 = Proxy DMatrix -> Ideal poly -> [poly]
forall poly (mat :: * -> *) (proxy :: (* -> *) -> *).
(Normed (Coefficient poly), IsOrderedPolynomial poly,
Field (Coefficient poly), Matrix mat (Coefficient poly)) =>
proxy mat -> Ideal poly -> [poly]
f4' (Proxy DMatrix
forall k (t :: k). Proxy t
Proxy :: Proxy DMatrix)
f4' ::(Normed (Coefficient poly),
IsOrderedPolynomial poly,
Field (Coefficient poly),
Matrix mat (Coefficient poly))
=> proxy mat -> Ideal poly -> [poly]
f4' :: proxy mat -> Ideal poly -> [poly]
f4' proxy mat
pxy = proxy mat -> Strategy poly Int -> Ideal poly -> [poly]
forall poly w (mat :: * -> *) (proxy :: (* -> *) -> *).
(Normed (Coefficient poly), Ord w, IsOrderedPolynomial poly,
Field (Coefficient poly), Matrix mat (Coefficient poly)) =>
proxy mat -> Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy' proxy mat
pxy Strategy poly Int
forall poly. IsOrderedPolynomial poly => poly -> poly -> Int
normalStrategy
f4WithStrategy' :: (Normed (Coefficient poly),
Ord w, IsOrderedPolynomial poly,
Field (Coefficient poly),
Matrix mat (Coefficient poly))
=> proxy mat -> Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy' :: proxy mat -> Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy' proxy mat
mrep Strategy poly w
select Ideal poly
ideal = (forall s. ST s [poly]) -> [poly]
forall a. (forall s. ST s a) -> a
runST ((forall s. ST s [poly]) -> [poly])
-> (forall s. ST s [poly]) -> [poly]
forall a b. (a -> b) -> a -> b
$ do
let is0 :: Vector poly
is0 = [poly] -> Vector poly
forall a. [a] -> Vector a
V.fromList ([poly] -> Vector poly) -> [poly] -> Vector poly
forall a b. (a -> b) -> a -> b
$ Ideal poly -> [poly]
forall r. Ideal r -> [r]
generators Ideal poly
ideal
buildHeap :: Vector poly -> [Int] -> Heap (Entry w (Int, Int))
buildHeap Vector poly
fs [Int]
is =
[Entry w (Int, Int)] -> Heap (Entry w (Int, Int))
forall a. Ord a => [a] -> Heap a
H.fromList [ w -> (Int, Int) -> Entry w (Int, Int)
forall p a. p -> a -> Entry p a
H.Entry (Strategy poly w
select (Vector poly
fs Vector poly -> Int -> poly
forall a. Vector a -> Int -> a
V.! Int
i) (Vector poly
fs Vector poly -> Int -> poly
forall a. Vector a -> Int -> a
V.! Int
j)) (Int
i, Int
j)
| Int
j <- [Int]
is, Int
i <- [Int
0..Int
jInt -> Int -> Int
forall r. Group r => r -> r -> r
-Int
1]
]
STRef s (MVector s poly)
gs <- MVector s poly -> ST s (STRef s (MVector s poly))
forall a s. a -> ST s (STRef s a)
newSTRef (MVector s poly -> ST s (STRef s (MVector s poly)))
-> ST s (MVector s poly) -> ST s (STRef s (MVector s poly))
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Vector poly -> ST s (MVector (PrimState (ST s)) poly)
forall (m :: * -> *) a.
PrimMonad m =>
Vector a -> m (MVector (PrimState m) a)
V.unsafeThaw Vector poly
is0
STRef s (Heap (Entry w (Int, Int)))
bs <- Heap (Entry w (Int, Int))
-> ST s (STRef s (Heap (Entry w (Int, Int))))
forall a s. a -> ST s (STRef s a)
newSTRef (Heap (Entry w (Int, Int))
-> ST s (STRef s (Heap (Entry w (Int, Int)))))
-> Heap (Entry w (Int, Int))
-> ST s (STRef s (Heap (Entry w (Int, Int))))
forall a b. (a -> b) -> a -> b
$ Vector poly -> [Int] -> Heap (Entry w (Int, Int))
buildHeap Vector poly
is0 [Int
0.. Vector poly -> Int
forall a. Vector a -> Int
V.length Vector poly
is0 Int -> Int -> Int
forall r. Group r => r -> r -> r
- Int
1]
let cancel :: Int -> Int -> ST s [poly]
cancel Int
i Int
j = do
(poly
fi, poly
fj) <- (,) (poly -> poly -> (poly, poly))
-> ST s poly -> ST s (poly -> (poly, poly))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> STRef s (MVector s poly)
gs STRef s (MVector s poly) -> Int -> ST s poly
forall s a. ArrayRef s a -> Int -> ST s a
%! Int
i ST s (poly -> (poly, poly)) -> ST s poly -> ST s (poly, poly)
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> STRef s (MVector s poly)
gs STRef s (MVector s poly) -> Int -> ST s poly
forall s a. ArrayRef s a -> Int -> ST s a
%! Int
j
let l :: OrderedMonomial (MOrder poly) (Arity poly)
l = poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> poly -> OrderedMonomial (MOrder poly) (Arity poly)
lcmLeadingMonomial poly
fi poly
fj
[poly] -> ST s [poly]
forall (m :: * -> *) a. Monad m => a -> m a
return ([poly] -> ST s [poly]) -> [poly] -> ST s [poly]
forall a b. (a -> b) -> a -> b
$ (poly -> poly) -> [poly] -> [poly]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (\poly
g -> poly -> poly
forall poly.
(Field (Coefficient poly), IsOrderedPolynomial poly) =>
poly -> poly
monoize (poly -> poly) -> poly -> poly
forall a b. (a -> b) -> a -> b
$ OrderedMonomial (MOrder poly) (Arity poly)
l OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
forall r. Division r => r -> r -> r
/ poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial poly
g OrderedMonomial (MOrder poly) (Arity poly) -> poly -> poly
forall poly.
IsOrderedPolynomial poly =>
OrderedMonomial (MOrder poly) (Arity poly) -> poly -> poly
>* poly
g) [poly
fi, poly
fj]
ST s (Maybe (Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int))))
-> ((Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int)))
-> ST s ())
-> ST s ()
forall (m :: * -> *) a b.
Monad m =>
m (Maybe a) -> (a -> m b) -> m ()
whileJust_ (Heap (Entry w (Int, Int))
-> Maybe (Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int)))
forall w. Ord w => Heap w -> Maybe (Heap w, Heap w)
viewMins (Heap (Entry w (Int, Int))
-> Maybe (Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int))))
-> ST s (Heap (Entry w (Int, Int)))
-> ST
s (Maybe (Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int))))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> STRef s (Heap (Entry w (Int, Int)))
-> ST s (Heap (Entry w (Int, Int)))
forall s a. STRef s a -> ST s a
readSTRef STRef s (Heap (Entry w (Int, Int)))
bs) (((Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int)))
-> ST s ())
-> ST s ())
-> ((Heap (Entry w (Int, Int)), Heap (Entry w (Int, Int)))
-> ST s ())
-> ST s ()
forall a b. (a -> b) -> a -> b
$ \(Heap (Entry w (Int, Int))
cur, Heap (Entry w (Int, Int))
bs') -> do
STRef s (Heap (Entry w (Int, Int)))
bs STRef s (Heap (Entry w (Int, Int)))
-> Heap (Entry w (Int, Int)) -> ST s ()
forall s a. STRef s a -> a -> ST s ()
.= Heap (Entry w (Int, Int))
bs'
[poly]
g <- STRef s (MVector s poly) -> ST s [poly]
forall s a. ArrayRef s a -> ST s [a]
arrayToList STRef s (MVector s poly)
gs
[poly]
ls <- [[poly]] -> [poly]
forall w. Monoid w => [w] -> w
concat ([[poly]] -> [poly]) -> ST s [[poly]] -> ST s [poly]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> (Entry w (Int, Int) -> ST s [poly])
-> [Entry w (Int, Int)] -> ST s [[poly]]
forall (t :: * -> *) (m :: * -> *) a b.
(Traversable t, Monad m) =>
(a -> m b) -> t a -> m (t b)
mapM ((Int -> Int -> ST s [poly]) -> (Int, Int) -> ST s [poly]
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry Int -> Int -> ST s [poly]
cancel ((Int, Int) -> ST s [poly])
-> (Entry w (Int, Int) -> (Int, Int))
-> Entry w (Int, Int)
-> ST s [poly]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Entry w (Int, Int) -> (Int, Int)
forall p a. Entry p a -> a
H.payload) (Heap (Entry w (Int, Int)) -> [Entry w (Int, Int)]
forall (t :: * -> *) a. Foldable t => t a -> [a]
F.toList Heap (Entry w (Int, Int))
cur)
let (Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs, mat (Coefficient poly)
mat) = proxy mat
-> [poly]
-> [poly]
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
forall poly (mat :: * -> *) (proxy :: (* -> *) -> *).
(IsOrderedPolynomial poly, Field (Coefficient poly),
Matrix mat (Coefficient poly)) =>
proxy mat
-> [poly]
-> [poly]
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
computeMatrix proxy mat
mrep [poly]
ls [poly]
g
ms :: [OrderedMonomial (MOrder poly) (Arity poly)]
ms = (Row (Mutable mat) (Coefficient poly)
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly)))
-> [Row (Mutable mat) (Coefficient poly)]
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a b. (a -> Maybe b) -> [a] -> [b]
mapMaybe (\Row (Mutable mat) (Coefficient poly)
xs -> First (OrderedMonomial (MOrder poly) (Arity poly))
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly))
forall a. First a -> Maybe a
getFirst (First (OrderedMonomial (MOrder poly) (Arity poly))
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly)))
-> First (OrderedMonomial (MOrder poly) (Arity poly))
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly))
forall a b. (a -> b) -> a -> b
$ (Int
-> OrderedMonomial (MOrder poly) (Arity poly)
-> First (OrderedMonomial (MOrder poly) (Arity poly)))
-> Vector (OrderedMonomial (MOrder poly) (Arity poly))
-> First (OrderedMonomial (MOrder poly) (Arity poly))
forall i (f :: * -> *) m a.
(FoldableWithIndex i f, Monoid m) =>
(i -> a -> m) -> f a -> m
ifoldMap (\Int
i OrderedMonomial (MOrder poly) (Arity poly)
m -> Maybe (OrderedMonomial (MOrder poly) (Arity poly))
-> First (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Maybe a -> First a
First (Maybe (OrderedMonomial (MOrder poly) (Arity poly))
-> First (OrderedMonomial (MOrder poly) (Arity poly)))
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly))
-> First (OrderedMonomial (MOrder poly) (Arity poly))
forall a b. (a -> b) -> a -> b
$
if Coefficient poly -> Bool
forall r. DecidableZero r => r -> Bool
isZero (Row (Mutable mat) (Coefficient poly)
xs Row (Mutable mat) (Coefficient poly) -> Int -> Coefficient poly
forall (v :: * -> *) a. Vector v a => v a -> Int -> a
GV.! Int
i)
then Maybe (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Maybe a
Nothing
else OrderedMonomial (MOrder poly) (Arity poly)
-> Maybe (OrderedMonomial (MOrder poly) (Arity poly))
forall a. a -> Maybe a
Just OrderedMonomial (MOrder poly) (Arity poly)
m)
Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs) ([Row (Mutable mat) (Coefficient poly)]
-> [OrderedMonomial (MOrder poly) (Arity poly)])
-> [Row (Mutable mat) (Coefficient poly)]
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a b. (a -> b) -> a -> b
$ mat (Coefficient poly) -> [Row mat (Coefficient poly)]
forall (mat :: * -> *) a. Matrix mat a => mat a -> [Row mat a]
toRows mat (Coefficient poly)
mat
(mat (Coefficient poly)
red, mat (Coefficient poly)
_, Coefficient poly
_) = mat (Coefficient poly)
-> (mat (Coefficient poly), mat (Coefficient poly),
Coefficient poly)
forall (mat :: * -> *) a.
(Matrix mat a, Eq a, Normed a, Field a) =>
mat a -> (mat a, mat a, a)
gaussReduction mat (Coefficient poly)
mat
ps :: [poly]
ps = (poly -> Bool) -> [poly] -> [poly]
forall a. (a -> Bool) -> [a] -> [a]
filter (\poly
f -> (OrderedMonomial (MOrder poly) (Arity poly) -> Bool)
-> [OrderedMonomial (MOrder poly) (Arity poly)] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (Bool -> Bool
not (Bool -> Bool)
-> (OrderedMonomial (MOrder poly) (Arity poly) -> Bool)
-> OrderedMonomial (MOrder poly) (Arity poly)
-> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly) -> Bool
forall k (n :: Nat) (ord :: k).
KnownNat n =>
OrderedMonomial ord n -> OrderedMonomial ord n -> Bool
`divs` poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial poly
f)) [OrderedMonomial (MOrder poly) (Arity poly)]
ms) ([poly] -> [poly]) -> [poly] -> [poly]
forall a b. (a -> b) -> a -> b
$
Vector (OrderedMonomial (MOrder poly) (Arity poly))
-> mat (Coefficient poly) -> [poly]
forall (mat :: * -> *) poly.
(Matrix mat (Coefficient poly), IsOrderedPolynomial poly) =>
Vector (OMonom poly) -> mat (Coefficient poly) -> [poly]
decodeMatrix Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs mat (Coefficient poly)
red
Bool -> ST s () -> ST s ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
unless ([poly] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [poly]
ps) (ST s () -> ST s ()) -> ST s () -> ST s ()
forall a b. (a -> b) -> a -> b
$ do
MVector s poly
g0 <- STRef s (MVector s poly) -> ST s (MVector s poly)
forall s a. STRef s a -> ST s a
readSTRef STRef s (MVector s poly)
gs
let len0 :: Int
len0 = MVector s poly -> Int
forall s a. MVector s a -> Int
MV.length MVector s poly
g0
ps' :: Vector poly
ps' = [poly] -> Vector poly
forall a. [a] -> Vector a
V.fromList [poly]
ps
size :: Int
size = Vector poly -> Int
forall a. Vector a -> Int
V.length Vector poly
ps'
MVector s poly
mv' <- MVector (PrimState (ST s)) poly
-> Int -> ST s (MVector (PrimState (ST s)) poly)
forall (m :: * -> *) a.
PrimMonad m =>
MVector (PrimState m) a -> Int -> m (MVector (PrimState m) a)
MV.unsafeGrow MVector s poly
MVector (PrimState (ST s)) poly
g0 Int
size
STRef s (MVector s poly)
gs STRef s (MVector s poly) -> MVector s poly -> ST s ()
forall s a. STRef s a -> a -> ST s ()
.= MVector s poly
mv'
MVector (PrimState (ST s)) poly
-> MVector (PrimState (ST s)) poly -> ST s ()
forall (m :: * -> *) a.
PrimMonad m =>
MVector (PrimState m) a -> MVector (PrimState m) a -> m ()
MV.copy (Int -> Int -> MVector s poly -> MVector s poly
forall s a. Int -> Int -> MVector s a -> MVector s a
MV.slice Int
len0 Int
size MVector s poly
mv') (MVector s poly -> ST s ()) -> ST s (MVector s poly) -> ST s ()
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Vector poly -> ST s (MVector (PrimState (ST s)) poly)
forall (m :: * -> *) a.
PrimMonad m =>
Vector a -> m (MVector (PrimState m) a)
V.unsafeThaw Vector poly
ps'
Vector poly
fs <- MVector (PrimState (ST s)) poly -> ST s (Vector poly)
forall (m :: * -> *) a.
PrimMonad m =>
MVector (PrimState m) a -> m (Vector a)
V.unsafeFreeze MVector s poly
MVector (PrimState (ST s)) poly
mv'
STRef s (Heap (Entry w (Int, Int)))
bs STRef s (Heap (Entry w (Int, Int)))
-> (Heap (Entry w (Int, Int)) -> Heap (Entry w (Int, Int)))
-> ST s ()
forall s a. STRef s a -> (a -> a) -> ST s ()
.%= Heap (Entry w (Int, Int))
-> Heap (Entry w (Int, Int)) -> Heap (Entry w (Int, Int))
forall a. Heap a -> Heap a -> Heap a
H.union (Vector poly -> [Int] -> Heap (Entry w (Int, Int))
buildHeap Vector poly
fs [Int
len0..Int
len0 Int -> Int -> Int
forall r. Additive r => r -> r -> r
+ Int
size Int -> Int -> Int
forall r. Group r => r -> r -> r
- Int
1])
STRef s (MVector s poly) -> ST s [poly]
forall s a. ArrayRef s a -> ST s [a]
arrayToList STRef s (MVector s poly)
gs
decodeMatrix :: (Matrix mat (Coefficient poly), IsOrderedPolynomial poly)
=> Vector (OMonom poly)
-> mat (Coefficient poly)
-> [poly]
decodeMatrix :: Vector (OMonom poly) -> mat (Coefficient poly) -> [poly]
decodeMatrix Vector (OMonom poly)
labs mat (Coefficient poly)
mat =
(poly -> Bool) -> [poly] -> [poly]
forall a. (a -> Bool) -> [a] -> [a]
filter (Bool -> Bool
not (Bool -> Bool) -> (poly -> Bool) -> poly -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. poly -> Bool
forall r. DecidableZero r => r -> Bool
isZero) ([poly] -> [poly]) -> [poly] -> [poly]
forall a b. (a -> b) -> a -> b
$
(Row (Mutable mat) (Coefficient poly) -> poly)
-> [Row (Mutable mat) (Coefficient poly)] -> [poly]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map ((Int -> Coefficient poly -> poly -> poly)
-> poly -> Row (Mutable mat) (Coefficient poly) -> poly
forall (v :: * -> *) a b.
Vector v a =>
(Int -> a -> b -> b) -> b -> v a -> b
GV.ifoldr (\Int
j Coefficient poly
d poly
c -> (Coefficient poly, OMonom poly) -> poly
forall poly.
IsOrderedPolynomial poly =>
(Coefficient poly, OrderedMonomial (MOrder poly) (Arity poly))
-> poly
toPolynomial (Coefficient poly
d, Vector (OMonom poly)
labs Vector (OMonom poly) -> Int -> OMonom poly
forall a. Vector a -> Int -> a
V.! Int
j) poly -> poly -> poly
forall r. Additive r => r -> r -> r
+ poly
c) poly
forall m. Monoidal m => m
zero) ([Row (Mutable mat) (Coefficient poly)] -> [poly])
-> [Row (Mutable mat) (Coefficient poly)] -> [poly]
forall a b. (a -> b) -> a -> b
$
mat (Coefficient poly) -> [Row mat (Coefficient poly)]
forall (mat :: * -> *) a. Matrix mat a => mat a -> [Row mat a]
toRows mat (Coefficient poly)
mat
f4WithStrategy :: (Field (Coefficient poly),
IsOrderedPolynomial poly,
Normed (Coefficient poly),
Num (Coefficient poly),
Ord w)
=> Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy :: Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy = Proxy DMatrix -> Strategy poly w -> Ideal poly -> [poly]
forall poly w (mat :: * -> *) (proxy :: (* -> *) -> *).
(Normed (Coefficient poly), Ord w, IsOrderedPolynomial poly,
Field (Coefficient poly), Matrix mat (Coefficient poly)) =>
proxy mat -> Strategy poly w -> Ideal poly -> [poly]
f4WithStrategy' (Proxy DMatrix
forall k (t :: k). Proxy t
Proxy :: Proxy DMatrix)
computeMatrix :: (IsOrderedPolynomial poly,
Field (Coefficient poly),
Matrix mat (Coefficient poly)
)
=> proxy mat -> [poly] -> [poly]
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)), mat (Coefficient poly))
computeMatrix :: proxy mat
-> [poly]
-> [poly]
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
computeMatrix proxy mat
_ [poly]
ls [poly]
gs = (forall s.
ST
s
(Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly)))
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
forall a. (forall s. ST s a) -> a
runST ((forall s.
ST
s
(Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly)))
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly)))
-> (forall s.
ST
s
(Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly)))
-> (Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
forall a b. (a -> b) -> a -> b
$ do
STRef s [poly]
hs <- [poly] -> ST s (STRef s [poly])
forall a s. a -> ST s (STRef s a)
newSTRef [poly]
ls
STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
done <- Set (OrderedMonomial (MOrder poly) (Arity poly))
-> ST
s (STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly))))
forall a s. a -> ST s (STRef s a)
newSTRef (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> ST
s (STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))))
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
-> ST
s (STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly))))
forall a b. (a -> b) -> a -> b
$ [OrderedMonomial (MOrder poly) (Arity poly)]
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Ord a => [a] -> Set a
S.fromList ([OrderedMonomial (MOrder poly) (Arity poly)]
-> Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall a b. (a -> b) -> a -> b
$ (poly -> OrderedMonomial (MOrder poly) (Arity poly))
-> [poly] -> [OrderedMonomial (MOrder poly) (Arity poly)]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial [poly]
ls
STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
monHs <- Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST
s
(STRef
s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))))
forall a s. a -> ST s (STRef s a)
newSTRef (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST
s
(STRef
s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST
s
(STRef
s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))))
forall a b. (a -> b) -> a -> b
$ [Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a. Ord a => [a] -> Heap a
H.fromList ([Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a b. (a -> b) -> a -> b
$ (OrderedMonomial (MOrder poly) (Arity poly)
-> Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map OrderedMonomial (MOrder poly) (Arity poly)
-> Down (OrderedMonomial (MOrder poly) (Arity poly))
forall a. a -> Down a
Down ([OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))])
-> [OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
forall a b. (a -> b) -> a -> b
$ Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a. Set a -> [a]
S.toList (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)])
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a b. (a -> b) -> a -> b
$
(poly -> Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> [poly] -> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall (t :: * -> *) m a.
(Foldable t, Monoid m) =>
(a -> m) -> t a -> m
foldMap (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Set a -> Set a
S.deleteMax (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> (poly -> Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> poly
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall b c a. (b -> c) -> (a -> b) -> a -> c
. poly -> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall poly.
IsOrderedPolynomial poly =>
poly -> Set (OrderedMonomial (MOrder poly) (Arity poly))
orderedMonomials) [poly]
ls
let ins :: poly -> ST s ()
ins poly
f = do
STRef s [poly] -> ([poly] -> [poly]) -> ST s ()
forall s a. STRef s a -> (a -> a) -> ST s ()
modifySTRef' STRef s [poly]
hs (poly
fpoly -> [poly] -> [poly]
forall a. a -> [a] -> [a]
:)
Set (OrderedMonomial (MOrder poly) (Arity poly))
del <- STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
forall s a. STRef s a -> ST s a
readSTRef STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
done
let ms :: Set (OrderedMonomial (MOrder poly) (Arity poly))
ms = poly -> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall poly.
IsOrderedPolynomial poly =>
poly -> Set (OrderedMonomial (MOrder poly) (Arity poly))
orderedMonomials poly
f Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Ord a => Set a -> Set a -> Set a
S.\\ Set (OrderedMonomial (MOrder poly) (Arity poly))
del
STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ()
forall s a. STRef s a -> (a -> a) -> ST s ()
modifySTRef' STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
monHs ((Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ())
-> (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ()
forall a b. (a -> b) -> a -> b
$ Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a. Heap a -> Heap a -> Heap a
H.union (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a b. (a -> b) -> a -> b
$ [Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a. Ord a => [a] -> Heap a
H.fromList ([Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
forall a b. (a -> b) -> a -> b
$ (OrderedMonomial (MOrder poly) (Arity poly)
-> Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map OrderedMonomial (MOrder poly) (Arity poly)
-> Down (OrderedMonomial (MOrder poly) (Arity poly))
forall a. a -> Down a
Down ([OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))])
-> [OrderedMonomial (MOrder poly) (Arity poly)]
-> [Down (OrderedMonomial (MOrder poly) (Arity poly))]
forall a b. (a -> b) -> a -> b
$ Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a. Set a -> [a]
S.toList Set (OrderedMonomial (MOrder poly) (Arity poly))
ms
ST
s
(Maybe
(Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))))
-> ((Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ())
-> ST s ()
forall (m :: * -> *) a b.
Monad m =>
m (Maybe a) -> (a -> m b) -> m ()
whileJust_ (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Maybe
(Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
forall a. Heap a -> Maybe (a, Heap a)
H.viewMin (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> Maybe
(Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))))
-> ST s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST
s
(Maybe
(Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
forall s a. STRef s a -> ST s a
readSTRef STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
monHs) (((Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ())
-> ST s ())
-> ((Down (OrderedMonomial (MOrder poly) (Arity poly)),
Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> ST s ())
-> ST s ()
forall a b. (a -> b) -> a -> b
$ \(Down OrderedMonomial (MOrder poly) (Arity poly)
m, Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
mHS') -> do
STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
monHs STRef s (Heap (Down (OrderedMonomial (MOrder poly) (Arity poly))))
-> Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s ()
forall s a. STRef s a -> a -> ST s ()
.= Heap (Down (OrderedMonomial (MOrder poly) (Arity poly)))
mHS'
STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
done STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s ()
forall s a. STRef s a -> (a -> a) -> ST s ()
.%= OrderedMonomial (MOrder poly) (Arity poly)
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
forall a. Ord a => a -> Set a -> Set a
S.insert OrderedMonomial (MOrder poly) (Arity poly)
m
Maybe poly -> (poly -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ ((poly -> Bool) -> [poly] -> Maybe poly
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Maybe a
find ((OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly) -> Bool
forall k (n :: Nat) (ord :: k).
KnownNat n =>
OrderedMonomial ord n -> OrderedMonomial ord n -> Bool
`divs` OrderedMonomial (MOrder poly) (Arity poly)
m)(OrderedMonomial (MOrder poly) (Arity poly) -> Bool)
-> (poly -> OrderedMonomial (MOrder poly) (Arity poly))
-> poly
-> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial) [poly]
gs) ((poly -> ST s ()) -> ST s ()) -> (poly -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \poly
g ->
poly -> ST s ()
ins (poly -> ST s ()) -> poly -> ST s ()
forall a b. (a -> b) -> a -> b
$ OrderedMonomial (MOrder poly) (Arity poly)
m OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
-> OrderedMonomial (MOrder poly) (Arity poly)
forall r. Division r => r -> r -> r
/ poly -> OrderedMonomial (MOrder poly) (Arity poly)
forall poly.
IsOrderedPolynomial poly =>
poly -> OrderedMonomial (MOrder poly) (Arity poly)
leadingMonomial poly
g OrderedMonomial (MOrder poly) (Arity poly) -> poly -> poly
forall poly.
IsOrderedPolynomial poly =>
OrderedMonomial (MOrder poly) (Arity poly) -> poly -> poly
>* poly
g
Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs <- [OrderedMonomial (MOrder poly) (Arity poly)]
-> Vector (OrderedMonomial (MOrder poly) (Arity poly))
forall a. [a] -> Vector a
V.fromList ([OrderedMonomial (MOrder poly) (Arity poly)]
-> Vector (OrderedMonomial (MOrder poly) (Arity poly)))
-> (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)])
-> Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Vector (OrderedMonomial (MOrder poly) (Arity poly))
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Set (OrderedMonomial (MOrder poly) (Arity poly))
-> [OrderedMonomial (MOrder poly) (Arity poly)]
forall a. Set a -> [a]
S.toDescList (Set (OrderedMonomial (MOrder poly) (Arity poly))
-> Vector (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s (Vector (OrderedMonomial (MOrder poly) (Arity poly)))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
-> ST s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
forall s a. STRef s a -> ST s a
readSTRef STRef s (Set (OrderedMonomial (MOrder poly) (Arity poly)))
done
[poly]
fs <- STRef s [poly] -> ST s [poly]
forall s a. STRef s a -> ST s a
readSTRef STRef s [poly]
hs
let mat :: mat (Coefficient poly)
mat = [Row mat (Coefficient poly)] -> mat (Coefficient poly)
forall (mat :: * -> *) a. Matrix mat a => [Row mat a] -> mat a
fromRows ([Row mat (Coefficient poly)] -> mat (Coefficient poly))
-> [Row mat (Coefficient poly)] -> mat (Coefficient poly)
forall a b. (a -> b) -> a -> b
$
(poly -> Row (Mutable mat) (Coefficient poly))
-> [poly] -> [Row (Mutable mat) (Coefficient poly)]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (\poly
f -> Int
-> (Int -> Coefficient poly)
-> Row (Mutable mat) (Coefficient poly)
forall (v :: * -> *) a. Vector v a => Int -> (Int -> a) -> v a
GV.generate (Vector (OrderedMonomial (MOrder poly) (Arity poly)) -> Int
forall a. Vector a -> Int
V.length Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs) ((Int -> Coefficient poly) -> Row (Mutable mat) (Coefficient poly))
-> (Int -> Coefficient poly)
-> Row (Mutable mat) (Coefficient poly)
forall a b. (a -> b) -> a -> b
$ \Int
i -> OrderedMonomial (MOrder poly) (Arity poly)
-> poly -> Coefficient poly
forall poly.
IsOrderedPolynomial poly =>
OrderedMonomial (MOrder poly) (Arity poly)
-> poly -> Coefficient poly
coeff (Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs Vector (OrderedMonomial (MOrder poly) (Arity poly))
-> Int -> OrderedMonomial (MOrder poly) (Arity poly)
forall (v :: * -> *) a. Vector v a => v a -> Int -> a
GV.! Int
i) poly
f) [poly]
fs
(Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
-> ST
s
(Vector (OrderedMonomial (MOrder poly) (Arity poly)),
mat (Coefficient poly))
forall (m :: * -> *) a. Monad m => a -> m a
return (Vector (OrderedMonomial (MOrder poly) (Arity poly))
labs, mat (Coefficient poly)
mat)