{-# LANGUAGE BangPatterns, DataKinds, ExtendedDefaultRules #-}
{-# LANGUAGE FlexibleContexts, GADTs, MultiParamTypeClasses #-}
{-# LANGUAGE NoImplicitPrelude, OverloadedLabels, OverloadedStrings #-}
{-# LANGUAGE PatternSynonyms, PolyKinds, RecordWildCards #-}
{-# LANGUAGE ScopedTypeVariables, TupleSections, TypeApplications #-}
{-# LANGUAGE ViewPatterns #-}
{-# OPTIONS_GHC -fno-warn-name-shadowing #-}
{-# OPTIONS_GHC -fno-warn-type-defaults #-}
module Algebra.Ring.Polynomial.Factorise
(
factorise, factorQBigPrime, factorHensel,
isReducible,
distinctDegFactor,
equalDegreeSplitM, equalDegreeFactorM,
henselStep, multiHensel, clearDenom,
squareFreePart, pthRoot,
yun, squareFreeDecompFiniteField
) where
import Algebra.Arithmetic hiding (modPow)
import Algebra.Field.Prime
import Algebra.Prelude.Core
import Algebra.Ring.Euclidean.Quotient
import Algebra.Ring.Polynomial.Univariate
import Control.Applicative ((<|>))
import Control.Lens (both, ifoldMap, (%~), (&))
import Control.Monad (guard, replicateM)
import Control.Monad.Loops (untilJust)
import Control.Monad.Random (MonadRandom, Random,
getRandom, getRandomR)
import qualified Data.DList as DL
import qualified Data.FMList as FML
import Data.IntMap (IntMap)
import qualified Data.IntMap.Strict as IM
import Data.Maybe (fromJust)
import Data.Monoid ((<>))
import Data.Proxy (Proxy (..))
import qualified Data.Set as S
import qualified Data.Traversable as F
import qualified Data.Vector as V
import Math.NumberTheory.Primes (Prime, precPrime)
import Math.NumberTheory.Primes (Prime (unPrime))
import qualified Math.NumberTheory.Primes as PRIMES
import Numeric.Decidable.Zero (isZero)
import Numeric.Domain.GCD (gcd, lcm)
import qualified Numeric.Field.Fraction as F
import qualified Prelude as P
distinctDegFactor :: forall k. (Eq k, FiniteField k)
=> Unipol k
-> [(Natural, Unipol k)]
distinctDegFactor :: Unipol k -> [(Natural, Unipol k)]
distinctDegFactor Unipol k
f0 = [Natural] -> [Unipol k] -> [(Natural, Unipol k)]
forall a b. [a] -> [b] -> [(a, b)]
zip [Natural
1..] ([Unipol k] -> [(Natural, Unipol k)])
-> [Unipol k] -> [(Natural, Unipol k)]
forall a b. (a -> b) -> a -> b
$ ([Unipol k] -> [Unipol k])
-> Unipol k -> Unipol k -> [Unipol k] -> [Unipol k]
go [Unipol k] -> [Unipol k]
forall a. a -> a
id (Ordinal (Arity (Unipol k)) -> Unipol k
forall poly. IsPolynomial poly => Ordinal (Arity poly) -> poly
var Ordinal (Arity (Unipol k))
forall (n :: Nat). (0 < n) => Ordinal n
OZ :: Unipol k) Unipol k
f0 []
where
go :: ([Unipol k] -> [Unipol k])
-> Unipol k -> Unipol k -> [Unipol k] -> [Unipol k]
go [Unipol k] -> [Unipol k]
gs Unipol k
h Unipol k
f =
let h' :: Unipol k
h' = Unipol k -> Natural -> Unipol k -> Unipol k
forall poly. Euclidean poly => poly -> Natural -> poly -> poly
modPow Unipol k
h (Proxy k -> Natural
forall k (proxy :: * -> *). FiniteField k => proxy k -> Natural
order (Proxy k
forall k (t :: k). Proxy t
Proxy :: Proxy k)) Unipol k
f
g' :: Unipol k
g' = Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd (Unipol k
h' Unipol k -> Unipol k -> Unipol k
forall r. Group r => r -> r -> r
- Ordinal (Arity (Unipol k)) -> Unipol k
forall poly. IsPolynomial poly => Ordinal (Arity poly) -> poly
var Ordinal (Arity (Unipol k))
0) Unipol k
f
f' :: Unipol k
f' = Unipol k
f Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
g'
gs' :: [Unipol k] -> [Unipol k]
gs' = [Unipol k] -> [Unipol k]
gs ([Unipol k] -> [Unipol k])
-> ([Unipol k] -> [Unipol k]) -> [Unipol k] -> [Unipol k]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Unipol k
g' Unipol k -> [Unipol k] -> [Unipol k]
forall a. a -> [a] -> [a]
:)
in if Unipol k
f' Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
== Unipol k
forall r. Unital r => r
one
then [Unipol k] -> [Unipol k]
gs'
else ([Unipol k] -> [Unipol k])
-> Unipol k -> Unipol k -> [Unipol k] -> [Unipol k]
go [Unipol k] -> [Unipol k]
gs' Unipol k
h' Unipol k
f'
modPow :: (Euclidean poly)
=> poly -> Natural -> poly -> poly
modPow :: poly -> Natural -> poly -> poly
modPow poly
a Natural
p poly
f = poly -> (forall r. Reifies r poly => Quotient r poly) -> poly
forall a.
Euclidean a =>
a -> (forall r. Reifies r a => Quotient r a) -> a
withQuotient poly
f ((forall r. Reifies r poly => Quotient r poly) -> poly)
-> (forall r. Reifies r poly => Quotient r poly) -> poly
forall a b. (a -> b) -> a -> b
$
Quotient r poly -> Natural -> Quotient r poly
forall r. Multiplicative r => r -> Natural -> r
repeatedSquare (poly -> Quotient r poly
forall k (r :: k) a.
(Reifies r a, Euclidean a) =>
a -> Quotient r a
quotient poly
a) Natural
p
traceCharTwo :: (Unital m, Monoidal m) => Natural -> m -> m
traceCharTwo :: Natural -> m -> m
traceCharTwo Natural
m m
a = (m -> m -> m) -> m -> [m] -> m
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' m -> m -> m
forall r. Additive r => r -> r -> r
(+) m
forall m. Monoidal m => m
zero ([m] -> m) -> [m] -> m
forall a b. (a -> b) -> a -> b
$ Int -> [m] -> [m]
forall a. Int -> [a] -> [a]
take (Natural -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Natural
m) ([m] -> [m]) -> [m] -> [m]
forall a b. (a -> b) -> a -> b
$ (m -> m) -> m -> [m]
forall a. (a -> a) -> a -> [a]
iterate (\(!m
x) ->m
xm -> m -> m
forall r. Multiplicative r => r -> r -> r
*m
x) m
a
equalDegreeSplitM
:: forall k m. (MonadRandom m, Random k, CoeffRing k, FiniteField k)
=> Unipol k
-> Natural
-> m (Maybe (Unipol k))
equalDegreeSplitM :: Unipol k -> Natural -> m (Maybe (Unipol k))
equalDegreeSplitM Unipol k
f Natural
d
| Int
n Int -> Int -> Int
forall a. Integral a => a -> a -> a
`mod` Natural -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Natural
d Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
0 = Maybe (Unipol k) -> m (Maybe (Unipol k))
forall (m :: * -> *) a. Monad m => a -> m a
return Maybe (Unipol k)
forall a. Maybe a
Nothing
| Bool
otherwise = do
let q :: Natural
q = Natural -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Natural -> Natural) -> Natural -> Natural
forall a b. (a -> b) -> a -> b
$ Proxy k -> Natural
forall k (proxy :: * -> *). FiniteField k => proxy k -> Natural
order (Proxy k
forall k (t :: k). Proxy t
Proxy :: Proxy k)
Int
e <- (Int, Int) -> m Int
forall (m :: * -> *) a. (MonadRandom m, Random a) => (a, a) -> m a
getRandomR (Int
1, Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
P.- Int
1)
[k]
cs <- Int -> m k -> m [k]
forall (m :: * -> *) a. Applicative m => Int -> m a -> m [a]
replicateM (Int -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
e) m k
forall (m :: * -> *) a. (MonadRandom m, Random a) => m a
getRandom
let a :: Unipol k
a = Ordinal (Arity (Unipol k)) -> Unipol k
forall poly. IsPolynomial poly => Ordinal (Arity poly) -> poly
var Ordinal (Arity (Unipol k))
0 Unipol k -> Natural -> Unipol k
forall r. Unital r => r -> Natural -> r
^ Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
e Unipol k -> Unipol k -> Unipol k
forall r. Additive r => r -> r -> r
+
[Unipol k] -> Unipol k
forall (f :: * -> *) m. (Foldable f, Monoidal m) => f m -> m
sum ((Unipol k -> Unipol k -> Unipol k)
-> [Unipol k] -> [Unipol k] -> [Unipol k]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Unipol k -> Unipol k -> Unipol k
forall r. Multiplicative r => r -> r -> r
(*) ((k -> Unipol k) -> [k] -> [Unipol k]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map k -> Unipol k
forall poly. IsPolynomial poly => Coefficient poly -> poly
injectCoeff [k]
cs) [Ordinal (Arity (Unipol k)) -> Unipol k
forall poly. IsPolynomial poly => Ordinal (Arity poly) -> poly
var Ordinal (Arity (Unipol k))
0 Unipol k -> Natural -> Unipol k
forall r. Unital r => r -> Natural -> r
^ Natural
l | Natural
l <-[Natural
0..]])
g1 :: Unipol k
g1 = Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd Unipol k
a Unipol k
f
Maybe (Unipol k) -> m (Maybe (Unipol k))
forall (m :: * -> *) a. Monad m => a -> m a
return (Maybe (Unipol k) -> m (Maybe (Unipol k)))
-> Maybe (Unipol k) -> m (Maybe (Unipol k))
forall a b. (a -> b) -> a -> b
$ (Bool -> Maybe ()
forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Unipol k
g1 Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= Unipol k
forall r. Unital r => r
one) Maybe () -> Maybe (Unipol k) -> Maybe (Unipol k)
forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> Unipol k -> Maybe (Unipol k)
forall (m :: * -> *) a. Monad m => a -> m a
return Unipol k
g1)
Maybe (Unipol k) -> Maybe (Unipol k) -> Maybe (Unipol k)
forall (f :: * -> *) a. Alternative f => f a -> f a -> f a
<|> do let b :: Unipol k
b | Unipol k -> Natural
forall r. Characteristic r => Unipol r -> Natural
charUnipol Unipol k
f Natural -> Natural -> Bool
forall a. Eq a => a -> a -> Bool
== Natural
2 =
Unipol k
-> (forall r. Reifies r (Unipol k) => Quotient r (Unipol k))
-> Unipol k
forall a.
Euclidean a =>
a -> (forall r. Reifies r a => Quotient r a) -> a
withQuotient Unipol k
f ((forall r. Reifies r (Unipol k) => Quotient r (Unipol k))
-> Unipol k)
-> (forall r. Reifies r (Unipol k) => Quotient r (Unipol k))
-> Unipol k
forall a b. (a -> b) -> a -> b
$
Natural -> Quotient r (Unipol k) -> Quotient r (Unipol k)
forall m. (Unital m, Monoidal m) => Natural -> m -> m
traceCharTwo (Unipol k -> Natural
forall r. FiniteField r => Unipol r -> Natural
powerUnipol Unipol k
fNatural -> Natural -> Natural
forall r. Multiplicative r => r -> r -> r
*Natural
d) (Unipol k -> Quotient r (Unipol k)
forall k (r :: k) a.
(Reifies r a, Euclidean a) =>
a -> Quotient r a
quotient Unipol k
a)
| Bool
otherwise = Unipol k -> Natural -> Unipol k -> Unipol k
forall poly. Euclidean poly => poly -> Natural -> poly -> poly
modPow Unipol k
a (Natural -> Natural
forall a. Enum a => a -> a
pred (Natural
qNatural -> Natural -> Natural
forall r. Unital r => r -> Natural -> r
^Natural
d) Natural -> Natural -> Natural
forall a. Integral a => a -> a -> a
`div`Natural
2) Unipol k
f
g2 :: Unipol k
g2 = Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd (Unipol k
b Unipol k -> Unipol k -> Unipol k
forall r. Group r => r -> r -> r
- Unipol k
forall r. Unital r => r
one) Unipol k
f
Bool -> Maybe ()
forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Unipol k
g2 Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= Unipol k
forall r. Unital r => r
one Bool -> Bool -> Bool
&& Unipol k
g2 Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= Unipol k
f)
Unipol k -> Maybe (Unipol k)
forall (m :: * -> *) a. Monad m => a -> m a
return Unipol k
g2
where n :: Int
n = Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol k
f
equalDegreeFactorM :: (Eq k, FiniteField k, MonadRandom m, Random k)
=> Unipol k -> Natural -> m [Unipol k]
equalDegreeFactorM :: Unipol k -> Natural -> m [Unipol k]
equalDegreeFactorM Unipol k
f Natural
d = Unipol k -> m ([Unipol k] -> [Unipol k])
go Unipol k
f m ([Unipol k] -> [Unipol k])
-> (([Unipol k] -> [Unipol k]) -> m [Unipol k]) -> m [Unipol k]
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= \[Unipol k] -> [Unipol k]
a -> [Unipol k] -> m [Unipol k]
forall (m :: * -> *) a. Monad m => a -> m a
return ([Unipol k] -> [Unipol k]
a [])
where
go :: Unipol k -> m ([Unipol k] -> [Unipol k])
go Unipol k
h
| Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol k
h Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 = ([Unipol k] -> [Unipol k]) -> m ([Unipol k] -> [Unipol k])
forall (m :: * -> *) a. Monad m => a -> m a
return [Unipol k] -> [Unipol k]
forall a. a -> a
id
| Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol k
h) Natural -> Natural -> Bool
forall a. Eq a => a -> a -> Bool
== Natural
d = ([Unipol k] -> [Unipol k]) -> m ([Unipol k] -> [Unipol k])
forall (m :: * -> *) a. Monad m => a -> m a
return (Unipol k
hUnipol k -> [Unipol k] -> [Unipol k]
forall a. a -> [a] -> [a]
:)
| Bool
otherwise = do
Unipol k
g <- m (Maybe (Unipol k)) -> m (Unipol k)
forall (m :: * -> *) a. Monad m => m (Maybe a) -> m a
untilJust (Unipol k -> Natural -> m (Maybe (Unipol k))
forall k (m :: * -> *).
(MonadRandom m, Random k, CoeffRing k, FiniteField k) =>
Unipol k -> Natural -> m (Maybe (Unipol k))
equalDegreeSplitM Unipol k
h Natural
d)
[Unipol k] -> [Unipol k]
l <- Unipol k -> m ([Unipol k] -> [Unipol k])
go Unipol k
g
[Unipol k] -> [Unipol k]
r <- Unipol k -> m ([Unipol k] -> [Unipol k])
go (Unipol k
h Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
g)
([Unipol k] -> [Unipol k]) -> m ([Unipol k] -> [Unipol k])
forall (m :: * -> *) a. Monad m => a -> m a
return (([Unipol k] -> [Unipol k]) -> m ([Unipol k] -> [Unipol k]))
-> ([Unipol k] -> [Unipol k]) -> m ([Unipol k] -> [Unipol k])
forall a b. (a -> b) -> a -> b
$ [Unipol k] -> [Unipol k]
l ([Unipol k] -> [Unipol k])
-> ([Unipol k] -> [Unipol k]) -> [Unipol k] -> [Unipol k]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Unipol k] -> [Unipol k]
r
factorSquareFree :: (Eq k, Random k, FiniteField k, MonadRandom m)
=> Unipol k -> m [Unipol k]
factorSquareFree :: Unipol k -> m [Unipol k]
factorSquareFree Unipol k
f =
[[Unipol k]] -> [Unipol k]
forall w. Monoid w => [w] -> w
concat ([[Unipol k]] -> [Unipol k]) -> m [[Unipol k]] -> m [Unipol k]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> ((Natural, Unipol k) -> m [Unipol k])
-> [(Natural, Unipol k)] -> m [[Unipol k]]
forall (t :: * -> *) (m :: * -> *) a b.
(Traversable t, Monad m) =>
(a -> m b) -> t a -> m (t b)
mapM ((Natural -> Unipol k -> m [Unipol k])
-> (Natural, Unipol k) -> m [Unipol k]
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry ((Natural -> Unipol k -> m [Unipol k])
-> (Natural, Unipol k) -> m [Unipol k])
-> (Natural -> Unipol k -> m [Unipol k])
-> (Natural, Unipol k)
-> m [Unipol k]
forall a b. (a -> b) -> a -> b
$ (Unipol k -> Natural -> m [Unipol k])
-> Natural -> Unipol k -> m [Unipol k]
forall a b c. (a -> b -> c) -> b -> a -> c
flip Unipol k -> Natural -> m [Unipol k]
forall k (m :: * -> *).
(Eq k, FiniteField k, MonadRandom m, Random k) =>
Unipol k -> Natural -> m [Unipol k]
equalDegreeFactorM) (((Natural, Unipol k) -> Bool)
-> [(Natural, Unipol k)] -> [(Natural, Unipol k)]
forall a. (a -> Bool) -> [a] -> [a]
filter ((Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= Unipol k
forall r. Unital r => r
one) (Unipol k -> Bool)
-> ((Natural, Unipol k) -> Unipol k) -> (Natural, Unipol k) -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Natural, Unipol k) -> Unipol k
forall a b. (a, b) -> b
snd) ([(Natural, Unipol k)] -> [(Natural, Unipol k)])
-> [(Natural, Unipol k)] -> [(Natural, Unipol k)]
forall a b. (a -> b) -> a -> b
$ Unipol k -> [(Natural, Unipol k)]
forall k.
(Eq k, FiniteField k) =>
Unipol k -> [(Natural, Unipol k)]
distinctDegFactor Unipol k
f)
squareFreePart :: (Eq k, FiniteField k)
=> Unipol k -> Unipol k
squareFreePart :: Unipol k -> Unipol k
squareFreePart Unipol k
f =
let !n :: Natural
n = Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Natural) -> Int -> Natural
forall a b. (a -> b) -> a -> b
$ Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol k
f
u :: Unipol k
u = Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd Unipol k
f (Ordinal (Arity (Unipol k)) -> Unipol k -> Unipol k
forall poly.
IsOrderedPolynomial poly =>
Ordinal (Arity poly) -> poly -> poly
diff Ordinal (Arity (Unipol k))
0 Unipol k
f)
v :: Unipol k
v = Unipol k
f Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
u
f' :: Unipol k
f' = Unipol k
u Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd Unipol k
u (Unipol k
v Unipol k -> Natural -> Unipol k
forall r. Unital r => r -> Natural -> r
^ Natural
n)
in if Unipol k
f' Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
== Unipol k
forall r. Unital r => r
one
then Unipol k
v
else Unipol k
v Unipol k -> Unipol k -> Unipol k
forall r. Multiplicative r => r -> r -> r
* Unipol k -> Unipol k
forall k. (Eq k, FiniteField k) => Unipol k -> Unipol k
squareFreePart (Unipol k -> Unipol k
forall r. (CoeffRing r, FiniteField r) => Unipol r -> Unipol r
pthRoot Unipol k
f')
yun :: (CoeffRing r, Field r)
=> Unipol r -> IntMap (Unipol r)
yun :: Unipol r -> IntMap (Unipol r)
yun Unipol r
f = let f' :: Unipol r
f' = Ordinal (Arity (Unipol r)) -> Unipol r -> Unipol r
forall poly.
IsOrderedPolynomial poly =>
Ordinal (Arity poly) -> poly -> poly
diff Ordinal (Arity (Unipol r))
forall (n :: Nat). (0 < n) => Ordinal n
OZ Unipol r
f
u :: Unipol r
u = Unipol r -> Unipol r -> Unipol r
forall d. GCDDomain d => d -> d -> d
gcd Unipol r
f Unipol r
f'
in Int
-> IntMap (Unipol r) -> Unipol r -> Unipol r -> IntMap (Unipol r)
forall t.
(IsOrderedPolynomial t, Euclidean t, (1 <=? Arity t) ~ 'True) =>
Int -> IntMap t -> t -> t -> IntMap t
go Int
1 IntMap (Unipol r)
forall a. IntMap a
IM.empty (Unipol r
f Unipol r -> Unipol r -> Unipol r
forall d. Euclidean d => d -> d -> d
`quot` Unipol r
u) (Unipol r
f' Unipol r -> Unipol r -> Unipol r
forall d. Euclidean d => d -> d -> d
`quot` Unipol r
u)
where
go :: Int -> IntMap t -> t -> t -> IntMap t
go !Int
i IntMap t
dic t
v t
w =
let t :: t
t = t
w t -> t -> t
forall r. Group r => r -> r -> r
- Ordinal (Arity t) -> t -> t
forall poly.
IsOrderedPolynomial poly =>
Ordinal (Arity poly) -> poly -> poly
diff Ordinal (Arity t)
forall (n :: Nat). (0 < n) => Ordinal n
OZ t
v
h :: t
h = t -> t -> t
forall d. GCDDomain d => d -> d -> d
gcd t
v t
t
v' :: t
v' = t
v t -> t -> t
forall d. Euclidean d => d -> d -> d
`quot` t
h
w' :: t
w' = t
t t -> t -> t
forall d. Euclidean d => d -> d -> d
`quot` t
h
dic' :: IntMap t
dic' = Int -> t -> IntMap t -> IntMap t
forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
i t
h IntMap t
dic
in if t
v' t -> t -> Bool
forall a. Eq a => a -> a -> Bool
== t
forall r. Unital r => r
one
then IntMap t
dic'
else Int -> IntMap t -> t -> t -> IntMap t
go (Int
iInt -> Int -> Int
forall r. Additive r => r -> r -> r
+Int
1) IntMap t
dic' t
v' t
w'
charUnipol :: forall r. Characteristic r => Unipol r -> Natural
charUnipol :: Unipol r -> Natural
charUnipol Unipol r
_ = Proxy r -> Natural
forall r (proxy :: * -> *). Characteristic r => proxy r -> Natural
char (Proxy r
forall k (t :: k). Proxy t
Proxy :: Proxy r)
powerUnipol :: forall r. FiniteField r => Unipol r -> Natural
powerUnipol :: Unipol r -> Natural
powerUnipol Unipol r
_ = Proxy r -> Natural
forall k (proxy :: * -> *). FiniteField k => proxy k -> Natural
power (Proxy r
forall k (t :: k). Proxy t
Proxy :: Proxy r)
pthRoot :: forall r. (CoeffRing r, FiniteField r) => Unipol r -> Unipol r
pthRoot :: Unipol r -> Unipol r
pthRoot Unipol r
f =
let !p :: Natural
p = Unipol r -> Natural
forall r. Characteristic r => Unipol r -> Natural
charUnipol Unipol r
f
!q :: Natural
q = Proxy r -> Natural
forall k (proxy :: * -> *). FiniteField k => proxy k -> Natural
order (Proxy r -> Natural) -> Proxy r -> Natural
forall a b. (a -> b) -> a -> b
$ Proxy r
forall k (t :: k). Proxy t
Proxy @r
!s :: Natural
s = Natural
q Natural -> Natural -> Natural
forall a. Integral a => a -> a -> a
`div` Natural
p
in [Unipol r] -> Unipol r
forall (f :: * -> *) m. (Foldable f, Monoidal m) => f m -> m
sum [ (r
ar -> Natural -> r
forall r. Unital r => r -> Natural -> r
^Natural
s) r -> Unipol r -> Unipol r
forall r m. Module (Scalar r) m => r -> m -> m
.*. IsLabel "x" (Unipol r)
Unipol r
#xUnipol r -> Natural -> Unipol r
forall r. Unital r => r -> Natural -> r
^Natural
i
| (Int
ip, r
a) <- IntMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList (IntMap r -> [(Int, r)]) -> IntMap r -> [(Int, r)]
forall a b. (a -> b) -> a -> b
$ Unipol r -> IntMap r
forall k. Unipol k -> IntMap k
termsUnipol Unipol r
f
, (Natural
i, Natural
0) <- [Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
ip Natural -> Natural -> (Natural, Natural)
forall a. Integral a => a -> a -> (a, a)
`P.divMod` Natural
p]
]
squareFreeDecompFiniteField
:: (Eq k, FiniteField k)
=> Unipol k -> IntMap (Unipol k)
squareFreeDecompFiniteField :: Unipol k -> IntMap (Unipol k)
squareFreeDecompFiniteField Unipol k
f =
let dcmp :: IntMap (Unipol k)
dcmp = Unipol k -> IntMap (Unipol k)
forall r. (CoeffRing r, Field r) => Unipol r -> IntMap (Unipol r)
yun Unipol k
f
h' :: Unipol k
h' = Mult (Unipol k) -> Unipol k
forall a. Mult a -> a
runMult (Mult (Unipol k) -> Unipol k) -> Mult (Unipol k) -> Unipol k
forall a b. (a -> b) -> a -> b
$ (Int -> Unipol k -> Mult (Unipol k))
-> IntMap (Unipol k) -> Mult (Unipol k)
forall i (f :: * -> *) m a.
(FoldableWithIndex i f, Monoid m) =>
(i -> a -> m) -> f a -> m
ifoldMap (\Int
i Unipol k
g -> Unipol k -> Mult (Unipol k)
forall a. a -> Mult a
Mult (Unipol k -> Mult (Unipol k)) -> Unipol k -> Mult (Unipol k)
forall a b. (a -> b) -> a -> b
$ Unipol k
g Unipol k -> Natural -> Unipol k
forall r. Unital r => r -> Natural -> r
^ Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
i) IntMap (Unipol k)
dcmp
fhInv :: Unipol k
fhInv = Unipol k
f Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
h'
p :: Int
p = Natural -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Natural -> Int) -> Natural -> Int
forall a b. (a -> b) -> a -> b
$ Unipol k -> Natural
forall r. Characteristic r => Unipol r -> Natural
charUnipol Unipol k
f
powers :: IntMap (Unipol k)
powers =
(Int -> Int) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. (Int -> Int) -> IntMap a -> IntMap a
IM.mapKeysMonotonic (Int
pInt -> Int -> Int
forall r. Multiplicative r => r -> r -> r
*) (IntMap (Unipol k) -> IntMap (Unipol k))
-> IntMap (Unipol k) -> IntMap (Unipol k)
forall a b. (a -> b) -> a -> b
$
Unipol k -> IntMap (Unipol k)
forall k. (Eq k, FiniteField k) => Unipol k -> IntMap (Unipol k)
squareFreeDecompFiniteField (Unipol k -> IntMap (Unipol k)) -> Unipol k -> IntMap (Unipol k)
forall a b. (a -> b) -> a -> b
$ Unipol k -> Unipol k
forall r. (CoeffRing r, FiniteField r) => Unipol r -> Unipol r
pthRoot Unipol k
fhInv
((IntMap (Unipol k)
pns, IntMap (Unipol k)
pnjs), IntMap (Unipol k)
js) =
((IntMap (Unipol k), IntMap (Unipol k))
-> Int
-> Unipol k
-> ((IntMap (Unipol k), IntMap (Unipol k)), Unipol k))
-> (IntMap (Unipol k), IntMap (Unipol k))
-> IntMap (Unipol k)
-> ((IntMap (Unipol k), IntMap (Unipol k)), IntMap (Unipol k))
forall a b c.
(a -> Int -> b -> (a, c)) -> a -> IntMap b -> (a, IntMap c)
IM.mapAccumWithKey
(\(IntMap (Unipol k)
ps', IntMap (Unipol k)
pnjs') Int
j Unipol k
g0 ->
let (Unipol k
g', IntMap (Unipol k, Unipol k)
pnjsps'') =
(Unipol k -> Unipol k -> (Unipol k, (Unipol k, Unipol k)))
-> Unipol k
-> IntMap (Unipol k)
-> (Unipol k, IntMap (Unipol k, Unipol k))
forall a b c. (a -> b -> (a, c)) -> a -> IntMap b -> (a, IntMap c)
IM.mapAccum
(\Unipol k
g Unipol k
h ->
let u :: Unipol k
u = Unipol k
g Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
`gcd` Unipol k
h
h' :: Unipol k
h' = Unipol k
h Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
u
in (Unipol k
g Unipol k -> Unipol k -> Unipol k
forall d. Euclidean d => d -> d -> d
`quot` Unipol k
u, (Unipol k
u, Unipol k
h'))
)
Unipol k
g0 IntMap (Unipol k)
ps'
pnjs'' :: IntMap (Unipol k)
pnjs'' = (Unipol k -> Bool) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. (a -> Bool) -> IntMap a -> IntMap a
IM.filter ((Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>Int
0) (Int -> Bool) -> (Unipol k -> Int) -> Unipol k -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree')
(IntMap (Unipol k) -> IntMap (Unipol k))
-> IntMap (Unipol k) -> IntMap (Unipol k)
forall a b. (a -> b) -> a -> b
$ (Int -> Int) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. (Int -> Int) -> IntMap a -> IntMap a
IM.mapKeysMonotonic (Int -> Int -> Int
forall r. Additive r => r -> r -> r
+ Int
j)
(IntMap (Unipol k) -> IntMap (Unipol k))
-> IntMap (Unipol k) -> IntMap (Unipol k)
forall a b. (a -> b) -> a -> b
$ (Unipol k, Unipol k) -> Unipol k
forall a b. (a, b) -> a
fst ((Unipol k, Unipol k) -> Unipol k)
-> IntMap (Unipol k, Unipol k) -> IntMap (Unipol k)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> IntMap (Unipol k, Unipol k)
pnjsps''
ps'' :: IntMap (Unipol k)
ps'' = (Unipol k -> Bool) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. (a -> Bool) -> IntMap a -> IntMap a
IM.filter ((Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>Int
0) (Int -> Bool) -> (Unipol k -> Int) -> Unipol k -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree')
(IntMap (Unipol k) -> IntMap (Unipol k))
-> IntMap (Unipol k) -> IntMap (Unipol k)
forall a b. (a -> b) -> a -> b
$ (Unipol k, Unipol k) -> Unipol k
forall a b. (a, b) -> b
snd ((Unipol k, Unipol k) -> Unipol k)
-> IntMap (Unipol k, Unipol k) -> IntMap (Unipol k)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> IntMap (Unipol k, Unipol k)
pnjsps''
in (( IntMap (Unipol k)
ps'', IntMap (Unipol k)
pnjs' IntMap (Unipol k) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. Semigroup a => a -> a -> a
<> IntMap (Unipol k)
pnjs''), Unipol k
g')
)
(IntMap (Unipol k)
powers, IntMap (Unipol k)
forall a. IntMap a
IM.empty)
(IntMap (Unipol k)
-> ((IntMap (Unipol k), IntMap (Unipol k)), IntMap (Unipol k)))
-> IntMap (Unipol k)
-> ((IntMap (Unipol k), IntMap (Unipol k)), IntMap (Unipol k))
forall a b. (a -> b) -> a -> b
$ (IntMap (Unipol k), IntMap (Unipol k)) -> IntMap (Unipol k)
forall a b. (a, b) -> a
fst (Int -> IntMap (Unipol k) -> (IntMap (Unipol k), IntMap (Unipol k))
forall a. Int -> IntMap a -> (IntMap a, IntMap a)
IM.split Int
p IntMap (Unipol k)
dcmp)
in if (Int, Unipol k) -> Int
forall a b. (a, b) -> a
fst (IntMap (Unipol k) -> (Int, Unipol k)
forall a. IntMap a -> (Int, a)
IM.findMax IntMap (Unipol k)
dcmp) Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
p Bool -> Bool -> Bool
&& Unipol k
fhInv Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
== Unipol k
forall r. Unital r => r
one
then IntMap (Unipol k)
dcmp
else [IntMap (Unipol k)] -> IntMap (Unipol k)
forall (f :: * -> *) a. Foldable f => f (IntMap a) -> IntMap a
IM.unions [(Unipol k -> Bool) -> IntMap (Unipol k) -> IntMap (Unipol k)
forall a. (a -> Bool) -> IntMap a -> IntMap a
IM.filter ((Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>Int
0) (Int -> Bool) -> (Unipol k -> Int) -> Unipol k -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree') IntMap (Unipol k)
js, IntMap (Unipol k)
pnjs, IntMap (Unipol k)
pns]
newtype GCD a = GCD { GCD a -> a
runGCD :: a }
instance Euclidean a => P.Semigroup (GCD a) where
<> :: GCD a -> GCD a -> GCD a
(<>) = (a -> GCD a
forall a. a -> GCD a
GCD (a -> GCD a) -> (GCD a -> a) -> GCD a -> GCD a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.) ((GCD a -> a) -> GCD a -> GCD a)
-> (GCD a -> GCD a -> a) -> GCD a -> GCD a -> GCD a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (a -> a -> a
forall d. GCDDomain d => d -> d -> d
gcd (a -> a -> a) -> (GCD a -> a) -> GCD a -> GCD a -> a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
`on` GCD a -> a
forall a. GCD a -> a
runGCD)
instance Euclidean a => Monoid (GCD a) where
mappend :: GCD a -> GCD a -> GCD a
mappend = GCD a -> GCD a -> GCD a
forall a. Semigroup a => a -> a -> a
(<>)
mempty :: GCD a
mempty = a -> GCD a
forall a. a -> GCD a
GCD a
forall m. Monoidal m => m
zero
newtype LCM a = LCM { LCM a -> a
runLCM :: a }
instance Euclidean a => P.Semigroup (LCM a) where
<> :: LCM a -> LCM a -> LCM a
(<>) = (a -> LCM a
forall a. a -> LCM a
LCM (a -> LCM a) -> (LCM a -> a) -> LCM a -> LCM a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.) ((LCM a -> a) -> LCM a -> LCM a)
-> (LCM a -> LCM a -> a) -> LCM a -> LCM a -> LCM a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (a -> a -> a
forall d. GCDDomain d => d -> d -> d
lcm (a -> a -> a) -> (LCM a -> a) -> LCM a -> LCM a -> a
forall b c a. (b -> b -> c) -> (a -> b) -> a -> a -> c
`on` LCM a -> a
forall a. LCM a -> a
runLCM)
instance Euclidean a => Monoid (LCM a) where
mappend :: LCM a -> LCM a -> LCM a
mappend = LCM a -> LCM a -> LCM a
forall a. Semigroup a => a -> a -> a
(<>)
mempty :: LCM a
mempty = a -> LCM a
forall a. a -> LCM a
LCM a
forall r. Unital r => r
one
factorise :: (MonadRandom m, CoeffRing k, FiniteField k, Random k)
=> Unipol k -> m [(Unipol k, Natural)]
factorise :: Unipol k -> m [(Unipol k, Natural)]
factorise Unipol k
f0
| Unipol k -> Bool
forall r. DecidableZero r => r -> Bool
isZero Unipol k
f0 = [(Unipol k, Natural)] -> m [(Unipol k, Natural)]
forall (f :: * -> *) a. Applicative f => a -> f a
pure [(Unipol k
forall m. Monoidal m => m
zero, Natural
1)]
| Bool
otherwise = do
let f :: Unipol k
f = Unipol k -> Unipol k
forall poly.
(Field (Coefficient poly), IsOrderedPolynomial poly) =>
poly -> poly
monoize Unipol k
f0
c :: Coefficient (Unipol k)
c = Unipol k -> Coefficient (Unipol k)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol k
f0
[(Unipol k, Natural)]
monoFacts <- [[(Unipol k, Natural)]] -> [(Unipol k, Natural)]
forall w. Monoid w => [w] -> w
concat
([[(Unipol k, Natural)]] -> [(Unipol k, Natural)])
-> m [[(Unipol k, Natural)]] -> m [(Unipol k, Natural)]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> ((Int, Unipol k) -> m [(Unipol k, Natural)])
-> [(Int, Unipol k)] -> m [[(Unipol k, Natural)]]
forall (t :: * -> *) (m :: * -> *) a b.
(Traversable t, Monad m) =>
(a -> m b) -> t a -> m (t b)
mapM (\(Int
r, Unipol k
h) -> (Unipol k -> (Unipol k, Natural))
-> [Unipol k] -> [(Unipol k, Natural)]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (,Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
r) ([Unipol k] -> [(Unipol k, Natural)])
-> m [Unipol k] -> m [(Unipol k, Natural)]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Unipol k -> m [Unipol k]
forall k (m :: * -> *).
(Eq k, Random k, FiniteField k, MonadRandom m) =>
Unipol k -> m [Unipol k]
factorSquareFree Unipol k
h) (IntMap (Unipol k) -> [(Int, Unipol k)]
forall a. IntMap a -> [(Int, a)]
IM.toList (IntMap (Unipol k) -> [(Int, Unipol k)])
-> IntMap (Unipol k) -> [(Int, Unipol k)]
forall a b. (a -> b) -> a -> b
$ Unipol k -> IntMap (Unipol k)
forall k. (Eq k, FiniteField k) => Unipol k -> IntMap (Unipol k)
squareFreeDecompFiniteField Unipol k
f)
[(Unipol k, Natural)] -> m [(Unipol k, Natural)]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([(Unipol k, Natural)] -> m [(Unipol k, Natural)])
-> [(Unipol k, Natural)] -> m [(Unipol k, Natural)]
forall a b. (a -> b) -> a -> b
$ if k
Coefficient (Unipol k)
c k -> k -> Bool
forall a. Eq a => a -> a -> Bool
== k
forall r. Unital r => r
one
then [(Unipol k, Natural)]
monoFacts
else (Coefficient (Unipol k) -> Unipol k
forall poly. IsPolynomial poly => Coefficient poly -> poly
injectCoeff Coefficient (Unipol k)
c, Natural
1) (Unipol k, Natural)
-> [(Unipol k, Natural)] -> [(Unipol k, Natural)]
forall a. a -> [a] -> [a]
: [(Unipol k, Natural)]
monoFacts
clearDenom :: (CoeffRing a, Euclidean a)
=> Unipol (Fraction a) -> (a, Unipol a)
clearDenom :: Unipol (Fraction a) -> (a, Unipol a)
clearDenom Unipol (Fraction a)
f =
let g :: a
g = (Fraction a -> a -> a) -> a -> Map (Monomial 1) (Fraction a) -> a
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr (a -> a -> a
forall d. GCDDomain d => d -> d -> d
lcm (a -> a -> a) -> (Fraction a -> a) -> Fraction a -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Fraction a -> a
forall t. Fraction t -> t
denominator) a
forall r. Unital r => r
one (Map (Monomial 1) (Fraction a) -> a)
-> Map (Monomial 1) (Fraction a) -> a
forall a b. (a -> b) -> a -> b
$ Unipol (Fraction a)
-> Map
(Monomial (Arity (Unipol (Fraction a))))
(Coefficient (Unipol (Fraction a)))
forall poly.
IsPolynomial poly =>
poly -> Map (Monomial (Arity poly)) (Coefficient poly)
terms' Unipol (Fraction a)
f
in (a
g, (Fraction a -> a) -> Unipol (Fraction a) -> Unipol a
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Fraction a -> a
forall t. Fraction t -> t
numerator (Fraction a -> a) -> (Fraction a -> Fraction a) -> Fraction a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((a
g a -> a -> Fraction a
forall d. GCDDomain d => d -> d -> Fraction d
F.% a
forall r. Unital r => r
one)Fraction a -> Fraction a -> Fraction a
forall r. Multiplicative r => r -> r -> r
*)) Unipol (Fraction a)
f)
factorQBigPrime :: (MonadRandom m)
=> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
factorQBigPrime :: Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
factorQBigPrime = (Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
forall (m :: * -> *).
MonadRandom m =>
(Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
wrapSQFFactor Unipol Integer -> m [Unipol Integer]
forall (m :: * -> *).
MonadRandom m =>
Unipol Integer -> m [Unipol Integer]
factorSqFreeQBP
factorHensel :: (MonadRandom m)
=> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
factorHensel :: Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
factorHensel = (Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
forall (m :: * -> *).
MonadRandom m =>
(Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
wrapSQFFactor Unipol Integer -> m [Unipol Integer]
forall (m :: * -> *).
MonadRandom m =>
Unipol Integer -> m [Unipol Integer]
factorHenselSqFree
wrapSQFFactor :: (MonadRandom m)
=> (Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
wrapSQFFactor :: (Unipol Integer -> m [Unipol Integer])
-> Unipol Integer -> m (Integer, IntMap (Set (Unipol Integer)))
wrapSQFFactor Unipol Integer -> m [Unipol Integer]
fac Unipol Integer
f0
| Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 = (Integer, IntMap (Set (Unipol Integer)))
-> m (Integer, IntMap (Set (Unipol Integer)))
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f0, IntMap (Set (Unipol Integer))
forall a. IntMap a
IM.empty)
| Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 = (Integer, IntMap (Set (Unipol Integer)))
-> m (Integer, IntMap (Set (Unipol Integer)))
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Unipol Integer -> Coefficient (Unipol Integer)
forall poly.
(IsPolynomial poly, Euclidean (Coefficient poly)) =>
poly -> Coefficient poly
content Unipol Integer
f0, Int -> Set (Unipol Integer) -> IntMap (Set (Unipol Integer))
forall a. Int -> a -> IntMap a
IM.singleton Int
1 (Set (Unipol Integer) -> IntMap (Set (Unipol Integer)))
-> Set (Unipol Integer) -> IntMap (Set (Unipol Integer))
forall a b. (a -> b) -> a -> b
$ Unipol Integer -> Set (Unipol Integer)
forall a. a -> Set a
S.singleton (Unipol Integer -> Set (Unipol Integer))
-> Unipol Integer -> Set (Unipol Integer)
forall a b. (a -> b) -> a -> b
$ Unipol Integer -> Unipol Integer
forall poly.
(Euclidean (Coefficient poly), IsPolynomial poly) =>
poly -> poly
pp Unipol Integer
f0)
| Bool
otherwise = do
let (Unipol Integer
g, Integer
c) = (Unipol Integer -> Unipol Integer
forall poly.
(Euclidean (Coefficient poly), IsPolynomial poly) =>
poly -> poly
pp Unipol Integer
f0, Unipol Integer -> Coefficient (Unipol Integer)
forall poly.
(IsPolynomial poly, Euclidean (Coefficient poly)) =>
poly -> Coefficient poly
content Unipol Integer
f0)
IntMap (Integer, [Unipol Integer])
ts0 <- (Unipol (Fraction Integer) -> m (Integer, [Unipol Integer]))
-> IntMap (Unipol (Fraction Integer))
-> m (IntMap (Integer, [Unipol Integer]))
forall (t :: * -> *) (m :: * -> *) a b.
(Traversable t, Monad m) =>
(a -> m b) -> t a -> m (t b)
F.mapM ((Unipol Integer -> m [Unipol Integer])
-> (Integer, Unipol Integer) -> m (Integer, [Unipol Integer])
forall (f :: * -> *) t a t1.
Functor f =>
(t -> f a) -> (t1, t) -> f (t1, a)
secondM Unipol Integer -> m [Unipol Integer]
fac ((Integer, Unipol Integer) -> m (Integer, [Unipol Integer]))
-> (Unipol (Fraction Integer) -> (Integer, Unipol Integer))
-> Unipol (Fraction Integer)
-> m (Integer, [Unipol Integer])
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Unipol (Fraction Integer) -> (Integer, Unipol Integer)
forall a.
(CoeffRing a, Euclidean a) =>
Unipol (Fraction a) -> (a, Unipol a)
clearDenom) (Unipol (Fraction Integer) -> IntMap (Unipol (Fraction Integer))
forall r. (CoeffRing r, Field r) => Unipol r -> IntMap (Unipol r)
yun (Unipol (Fraction Integer) -> IntMap (Unipol (Fraction Integer)))
-> Unipol (Fraction Integer) -> IntMap (Unipol (Fraction Integer))
forall a b. (a -> b) -> a -> b
$ Unipol (Fraction Integer) -> Unipol (Fraction Integer)
forall poly.
(Field (Coefficient poly), IsOrderedPolynomial poly) =>
poly -> poly
monoize (Unipol (Fraction Integer) -> Unipol (Fraction Integer))
-> Unipol (Fraction Integer) -> Unipol (Fraction Integer)
forall a b. (a -> b) -> a -> b
$ (Integer -> Fraction Integer)
-> Unipol Integer -> Unipol (Fraction Integer)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Integer -> Integer -> Fraction Integer
forall d. GCDDomain d => d -> d -> Fraction d
F.% Integer
1) Unipol Integer
g)
let anss :: [(Int, (Integer, [Unipol Integer]))]
anss = IntMap (Integer, [Unipol Integer])
-> [(Int, (Integer, [Unipol Integer]))]
forall a. IntMap a -> [(Int, a)]
IM.toList IntMap (Integer, [Unipol Integer])
ts0
k :: Integer
k = Integer
c Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
g Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`div` [Integer] -> Integer
forall (f :: * -> *) r. (Foldable f, Unital r) => f r -> r
product (((Int, (Integer, [Unipol Integer])) -> Integer)
-> [(Int, (Integer, [Unipol Integer]))] -> [Integer]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (\(Int
p,(Integer
i,[Unipol Integer]
_)) -> Integer
i Integer -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
P.^ Int
p) [(Int, (Integer, [Unipol Integer]))]
anss)
(Integer, IntMap (Set (Unipol Integer)))
-> m (Integer, IntMap (Set (Unipol Integer)))
forall (m :: * -> *) a. Monad m => a -> m a
return (Integer
k, [(Int, Set (Unipol Integer))] -> IntMap (Set (Unipol Integer))
forall a. [(Int, a)] -> IntMap a
IM.fromList ([(Int, Set (Unipol Integer))] -> IntMap (Set (Unipol Integer)))
-> [(Int, Set (Unipol Integer))] -> IntMap (Set (Unipol Integer))
forall a b. (a -> b) -> a -> b
$ ((Int, (Integer, [Unipol Integer])) -> (Int, Set (Unipol Integer)))
-> [(Int, (Integer, [Unipol Integer]))]
-> [(Int, Set (Unipol Integer))]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (((Integer, [Unipol Integer]) -> Set (Unipol Integer))
-> (Int, (Integer, [Unipol Integer]))
-> (Int, Set (Unipol Integer))
forall (a :: * -> * -> *) b c d.
Arrow a =>
a b c -> a (d, b) (d, c)
second (((Integer, [Unipol Integer]) -> Set (Unipol Integer))
-> (Int, (Integer, [Unipol Integer]))
-> (Int, Set (Unipol Integer)))
-> ((Integer, [Unipol Integer]) -> Set (Unipol Integer))
-> (Int, (Integer, [Unipol Integer]))
-> (Int, Set (Unipol Integer))
forall a b. (a -> b) -> a -> b
$ [Unipol Integer] -> Set (Unipol Integer)
forall a. Ord a => [a] -> Set a
S.fromList ([Unipol Integer] -> Set (Unipol Integer))
-> ((Integer, [Unipol Integer]) -> [Unipol Integer])
-> (Integer, [Unipol Integer])
-> Set (Unipol Integer)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer, [Unipol Integer]) -> [Unipol Integer]
forall a b. (a, b) -> b
snd) [(Int, (Integer, [Unipol Integer]))]
anss)
where
n :: Int
n = Unipol Integer -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol Integer
f0
secondM :: Functor f => (t -> f a) -> (t1, t) -> f (t1, a)
secondM :: (t -> f a) -> (t1, t) -> f (t1, a)
secondM t -> f a
f (t1
a, t
b)= (t1
a,) (a -> (t1, a)) -> f a -> f (t1, a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> t -> f a
f t
b
factorSqFreeQBP :: (MonadRandom m)
=> Unipol Integer -> m [Unipol Integer]
{-# INLINE factorSqFreeQBP #-}
factorSqFreeQBP :: Unipol Integer -> m [Unipol Integer]
factorSqFreeQBP =
(Prime Integer -> Unipol Integer -> BasicConsts)
-> (BasicConsts -> Unipol Integer -> Maybe (Integer, Integer))
-> (BasicConsts -> Prime Integer -> Unipol Integer -> Integer)
-> (BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
forall (m :: * -> *) a.
MonadRandom m =>
(Prime Integer -> Unipol Integer -> a)
-> (a -> Unipol Integer -> Maybe (Integer, Integer))
-> (a -> Prime Integer -> Unipol Integer -> Integer)
-> (a
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
factorPrimitiveIntegerPolWith
((Unipol Integer -> BasicConsts)
-> Prime Integer -> Unipol Integer -> BasicConsts
forall a b. a -> b -> a
const Unipol Integer -> BasicConsts
calcBasicConsts)
(\BasicConsts{Double
Integer
lc :: BasicConsts -> Integer
kB :: BasicConsts -> Double
kA :: BasicConsts -> Integer
lc :: Integer
kB :: Double
kA :: Integer
..} -> Maybe (Integer, Integer)
-> Unipol Integer -> Maybe (Integer, Integer)
forall a b. a -> b -> a
const (Maybe (Integer, Integer)
-> Unipol Integer -> Maybe (Integer, Integer))
-> Maybe (Integer, Integer)
-> Unipol Integer
-> Maybe (Integer, Integer)
forall a b. (a -> b) -> a -> b
$
(Integer, Integer) -> Maybe (Integer, Integer)
forall a. a -> Maybe a
Just (Double -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
floor (Double -> Integer) -> Double -> Integer
forall a b. (a -> b) -> a -> b
$ Double
2Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
*Double
kB, Double -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
ceiling (Double -> Integer) -> Double -> Integer
forall a b. (a -> b) -> a -> b
$ Double
4Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
*Double
kB)
)
((Prime Integer -> Unipol Integer -> Integer)
-> BasicConsts -> Prime Integer -> Unipol Integer -> Integer
forall a b. a -> b -> a
const ((Prime Integer -> Unipol Integer -> Integer)
-> BasicConsts -> Prime Integer -> Unipol Integer -> Integer)
-> (Prime Integer -> Unipol Integer -> Integer)
-> BasicConsts
-> Prime Integer
-> Unipol Integer
-> Integer
forall a b. (a -> b) -> a -> b
$ Integer -> Unipol Integer -> Integer
forall a b. a -> b -> a
const (Integer -> Unipol Integer -> Integer)
-> (Prime Integer -> Integer)
-> Prime Integer
-> Unipol Integer
-> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Prime Integer -> Integer
forall a. Prime a -> a
unPrime)
((BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer -> m [Unipol Integer])
-> (BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ (Prime Integer
-> Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
forall a b. a -> b -> a
const ((Prime Integer
-> Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> (Prime Integer
-> Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> BasicConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ (Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
forall a b. a -> b -> a
const ((Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> (Unipol Integer -> [Unipol Integer] -> [Unipol Integer])
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ ([Unipol Integer] -> [Unipol Integer])
-> Unipol Integer -> [Unipol Integer] -> [Unipol Integer]
forall a b. a -> b -> a
const [Unipol Integer] -> [Unipol Integer]
forall a. a -> a
id
data BasicConsts =
BasicConsts
{ BasicConsts -> Integer
kA :: !Integer
, BasicConsts -> Double
kB :: !Double
, BasicConsts -> Integer
lc :: !Integer
}
data HenselConsts = HenselConsts
{ HenselConsts -> BasicConsts
basicConsts :: !BasicConsts
, HenselConsts -> Natural
primePower :: Natural
}
factorHenselSqFree
:: MonadRandom m
=> Unipol Integer -> m [Unipol Integer]
{-# INLINE factorHenselSqFree #-}
factorHenselSqFree :: Unipol Integer -> m [Unipol Integer]
factorHenselSqFree =
(Prime Integer -> Unipol Integer -> HenselConsts)
-> (HenselConsts -> Unipol Integer -> Maybe (Integer, Integer))
-> (HenselConsts -> Prime Integer -> Unipol Integer -> Integer)
-> (HenselConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
forall (m :: * -> *) a.
MonadRandom m =>
(Prime Integer -> Unipol Integer -> a)
-> (a -> Unipol Integer -> Maybe (Integer, Integer))
-> (a -> Prime Integer -> Unipol Integer -> Integer)
-> (a
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
factorPrimitiveIntegerPolWith
Prime Integer -> Unipol Integer -> HenselConsts
calcHenselConsts
((Unipol Integer -> Maybe (Integer, Integer))
-> HenselConsts -> Unipol Integer -> Maybe (Integer, Integer)
forall a b. a -> b -> a
const ((Unipol Integer -> Maybe (Integer, Integer))
-> HenselConsts -> Unipol Integer -> Maybe (Integer, Integer))
-> (Unipol Integer -> Maybe (Integer, Integer))
-> HenselConsts
-> Unipol Integer
-> Maybe (Integer, Integer)
forall a b. (a -> b) -> a -> b
$ Maybe (Integer, Integer)
-> Unipol Integer -> Maybe (Integer, Integer)
forall a b. a -> b -> a
const Maybe (Integer, Integer)
forall a. Maybe a
Nothing)
(\HenselConsts{Natural
BasicConsts
primePower :: Natural
basicConsts :: BasicConsts
primePower :: HenselConsts -> Natural
basicConsts :: HenselConsts -> BasicConsts
..} (Prime Integer -> Integer
forall a. Prime a -> a
unPrime -> Integer
p) Unipol Integer
_ ->
Integer
p Integer -> Natural -> Integer
forall r. Unital r => r -> Natural -> r
^ Natural
primePower
)
((HenselConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer -> m [Unipol Integer])
-> (HenselConsts
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ \HenselConsts{basicConsts :: HenselConsts -> BasicConsts
basicConsts = BasicConsts{Double
Integer
lc :: Integer
kB :: Double
kA :: Integer
lc :: BasicConsts -> Integer
kB :: BasicConsts -> Double
kA :: BasicConsts -> Integer
..}, Natural
primePower :: Natural
primePower :: HenselConsts -> Natural
..} (Prime Integer -> Integer
forall a. Prime a -> a
unPrime -> Integer
p) Unipol Integer
f [Unipol Integer]
gs ->
Natural
-> Int -> Unipol Integer -> [Unipol Integer] -> [Unipol Integer]
multiHensel (Integer -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral Integer
p) (Natural -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Natural
primePower) Unipol Integer
f [Unipol Integer]
gs
{-# INLINE calcHenselConsts #-}
calcHenselConsts :: Prime Integer -> Unipol Integer -> HenselConsts
calcHenselConsts :: Prime Integer -> Unipol Integer -> HenselConsts
calcHenselConsts Prime Integer
p Unipol Integer
f =
let primePower :: Natural
primePower = Double -> Natural
forall a b. (RealFrac a, Integral b) => a -> b
ceiling (Double -> Natural) -> Double -> Natural
forall a b. (a -> b) -> a -> b
$ Double -> Double -> Double
forall a. Floating a => a -> a -> a
logBase (Integer -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Integer -> Double) -> Integer -> Double
forall a b. (a -> b) -> a -> b
$ Prime Integer -> Integer
forall a. Prime a -> a
unPrime Prime Integer
p) (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Double
2 Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
* Double
kB Double -> Double -> Double
forall r. Additive r => r -> r -> r
+ Double
1
basicConsts :: BasicConsts
basicConsts@BasicConsts{Double
Integer
lc :: Integer
kA :: Integer
kB :: Double
lc :: BasicConsts -> Integer
kB :: BasicConsts -> Double
kA :: BasicConsts -> Integer
..} = Unipol Integer -> BasicConsts
calcBasicConsts Unipol Integer
f
in HenselConsts :: BasicConsts -> Natural -> HenselConsts
HenselConsts{Natural
BasicConsts
basicConsts :: BasicConsts
primePower :: Natural
primePower :: Natural
basicConsts :: BasicConsts
..}
calcBasicConsts :: Unipol Integer -> BasicConsts
calcBasicConsts :: Unipol Integer -> BasicConsts
calcBasicConsts Unipol Integer
f =
let lc :: Coefficient (Unipol Integer)
lc = Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f
kA :: Norm (Coefficient (Unipol Integer))
kA = Unipol Integer -> Norm (Coefficient (Unipol Integer))
forall poly.
(IsPolynomial poly, Normed (Coefficient poly)) =>
poly -> Norm (Coefficient poly)
maxNorm Unipol Integer
f
n :: Int
n = Unipol Integer -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol Integer
f
kB :: Double
kB = Double -> Double
forall a. Floating a => a -> a
sqrt (Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n Double -> Double -> Double
forall r. Additive r => r -> r -> r
+ Double
1) Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
* Double
2 Double -> Int -> Double
forall a b. (Num a, Integral b) => a -> b -> a
P.^ Int
n Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
* Integer -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Integer
Norm (Coefficient (Unipol Integer))
kA Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Integer
Coefficient (Unipol Integer)
lc)
in BasicConsts :: Integer -> Double -> Integer -> BasicConsts
BasicConsts{Double
Integer
Coefficient (Unipol Integer)
Norm (Coefficient (Unipol Integer))
kB :: Double
kA :: Norm (Coefficient (Unipol Integer))
lc :: Coefficient (Unipol Integer)
lc :: Integer
kB :: Double
kA :: Integer
..}
factorPrimitiveIntegerPolWith
:: MonadRandom m
=> (Prime Integer -> Unipol Integer -> a)
-> (a -> Unipol Integer -> Maybe (Integer, Integer))
-> (a -> Prime Integer -> Unipol Integer -> Integer)
-> (a -> Prime Integer -> Unipol Integer
-> [Unipol Integer] -> [Unipol Integer])
-> Unipol Integer -> m [Unipol Integer]
factorPrimitiveIntegerPolWith :: (Prime Integer -> Unipol Integer -> a)
-> (a -> Unipol Integer -> Maybe (Integer, Integer))
-> (a -> Prime Integer -> Unipol Integer -> Integer)
-> (a
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer])
-> Unipol Integer
-> m [Unipol Integer]
factorPrimitiveIntegerPolWith Prime Integer -> Unipol Integer -> a
pre a -> Unipol Integer -> Maybe (Integer, Integer)
mRange a -> Prime Integer -> Unipol Integer -> Integer
modulus a
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
postProcess Unipol Integer
f =
let lc :: Coefficient (Unipol Integer)
lc = Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f
info :: a
info = Prime Integer -> Unipol Integer -> a
pre Prime Integer
p Unipol Integer
f
range :: [Prime Integer]
range = [Prime Integer]
-> ((Integer, Integer) -> [Prime Integer])
-> Maybe (Integer, Integer)
-> [Prime Integer]
forall b a. b -> (a -> b) -> Maybe a -> b
maybe [Prime Integer]
forall a. Integral a => [Prime a]
PRIMES.primes (\(Integer
l,Integer
r) -> [Integer -> Prime Integer
forall a.
(Bits a, Integral a, UniqueFactorisation a) =>
a -> Prime a
precPrime Integer
l .. Integer -> Prime Integer
forall a.
(Bits a, Integral a, UniqueFactorisation a) =>
a -> Prime a
precPrime Integer
r])
(Maybe (Integer, Integer) -> [Prime Integer])
-> Maybe (Integer, Integer) -> [Prime Integer]
forall a b. (a -> b) -> a -> b
$ a -> Unipol Integer -> Maybe (Integer, Integer)
mRange a
info Unipol Integer
f
Just Prime Integer
p = (Prime Integer -> Bool) -> [Prime Integer] -> Maybe (Prime Integer)
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Maybe a
find Prime Integer -> Bool
isGoodPrime [Prime Integer]
range
in Integer
-> (forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> m [Unipol Integer])
-> m [Unipol Integer]
forall a.
Integer -> (forall (p :: Nat). KnownNat p => Proxy (F p) -> a) -> a
reifyPrimeField (Prime Integer -> Integer
forall a. Prime a -> a
unPrime Prime Integer
p) ((forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> m [Unipol Integer])
-> m [Unipol Integer])
-> (forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> m [Unipol Integer])
-> m [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ \Proxy (F p)
fp -> do
let lc' :: F p
lc' = Proxy (F p) -> Integer -> F p
forall k (proxy :: * -> *) (p :: k).
Reifies p Integer =>
proxy (F p) -> Integer -> F p
modNat' Proxy (F p)
fp Integer
Coefficient (Unipol Integer)
lc
f0 :: Unipol (F p)
f0 = (Integer -> F p) -> Unipol Integer -> Unipol (F p)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol ((F p -> F p -> F p
forall r. Division r => r -> r -> r
/F p
lc') (F p -> F p) -> (Integer -> F p) -> Integer -> F p
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Proxy (F p) -> Integer -> F p
forall k (proxy :: * -> *) (p :: k).
Reifies p Integer =>
proxy (F p) -> Integer -> F p
modNat' Proxy (F p)
fp) Unipol Integer
f
[Unipol Integer]
fps <- ((Unipol (F p), Natural) -> Unipol Integer)
-> [(Unipol (F p), Natural)] -> [Unipol Integer]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (Integer -> Unipol Integer -> Unipol Integer
normalizeMod (Prime Integer -> Integer
forall a. Prime a -> a
unPrime Prime Integer
p)
(Unipol Integer -> Unipol Integer)
-> ((Unipol (F p), Natural) -> Unipol Integer)
-> (Unipol (F p), Natural)
-> Unipol Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (F p -> Integer) -> Unipol (F p) -> Unipol Integer
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol F p -> Integer
forall k (p :: k). F p -> Integer
naturalRepr (Unipol (F p) -> Unipol Integer)
-> ((Unipol (F p), Natural) -> Unipol (F p))
-> (Unipol (F p), Natural)
-> Unipol Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Unipol (F p), Natural) -> Unipol (F p)
forall a b. (a, b) -> a
fst)
([(Unipol (F p), Natural)] -> [Unipol Integer])
-> m [(Unipol (F p), Natural)] -> m [Unipol Integer]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Unipol (F p) -> m [(Unipol (F p), Natural)]
forall (m :: * -> *) k.
(MonadRandom m, CoeffRing k, FiniteField k, Random k) =>
Unipol k -> m [(Unipol k, Natural)]
factorise Unipol (F p)
f0
let gs :: Vector (Unipol Integer)
gs = [Unipol Integer] -> Vector (Unipol Integer)
forall a. [a] -> Vector a
V.fromList ([Unipol Integer] -> Vector (Unipol Integer))
-> [Unipol Integer] -> Vector (Unipol Integer)
forall a b. (a -> b) -> a -> b
$ a
-> Prime Integer
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
postProcess a
info Prime Integer
p Unipol Integer
f [Unipol Integer]
fps
count :: Int
count = Vector (Unipol Integer) -> Int
forall a. Vector a -> Int
V.length Vector (Unipol Integer)
gs
alives :: [Int]
alives = [Int
0 .. Int
count Int -> Int -> Int
forall r. Group r => r -> r -> r
- Int
1]
[Unipol Integer] -> m [Unipol Integer]
forall (m :: * -> *) a. Monad m => a -> m a
return ([Unipol Integer] -> m [Unipol Integer])
-> [Unipol Integer] -> m [Unipol Integer]
forall a b. (a -> b) -> a -> b
$ Integer
-> [Int]
-> Int
-> Int
-> Unipol Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
loop (a -> Prime Integer -> Unipol Integer -> Integer
modulus a
info Prime Integer
p Unipol Integer
f) [Int]
alives Int
count Int
1 Unipol Integer
f Vector (Unipol Integer)
gs []
where
lc :: Coefficient (Unipol Integer)
lc = Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f
isGoodPrime :: Prime Integer -> Bool
isGoodPrime Prime Integer
p = Integer
-> (forall (p :: Nat). KnownNat p => Proxy (F p) -> Bool) -> Bool
forall a.
Integer -> (forall (p :: Nat). KnownNat p => Proxy (F p) -> a) -> a
reifyPrimeField (Prime Integer -> Integer
forall a. Prime a -> a
unPrime Prime Integer
p) ((forall (p :: Nat). KnownNat p => Proxy (F p) -> Bool) -> Bool)
-> (forall (p :: Nat). KnownNat p => Proxy (F p) -> Bool) -> Bool
forall a b. (a -> b) -> a -> b
$ \Proxy (F p)
fp ->
Integer
Coefficient (Unipol Integer)
lc Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`mod` Prime Integer -> Integer
forall a. Prime a -> a
unPrime Prime Integer
p Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
/= Integer
0 Bool -> Bool -> Bool
&& Unipol (F p) -> Bool
forall poly.
(IsOrderedPolynomial poly, GCDDomain poly) =>
poly -> Bool
isSquareFree ((Integer -> F p) -> Unipol Integer -> Unipol (F p)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Proxy (F p) -> Integer -> F p
forall k (proxy :: * -> *) (p :: k).
Reifies p Integer =>
proxy (F p) -> Integer -> F p
modNat' Proxy (F p)
fp) Unipol Integer
f)
kA :: Norm (Coefficient (Unipol Integer))
kA = Unipol Integer -> Norm (Coefficient (Unipol Integer))
forall poly.
(IsPolynomial poly, Normed (Coefficient poly)) =>
poly -> Norm (Coefficient poly)
maxNorm Unipol Integer
f
n :: Int
n = Unipol Integer -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol Integer
f
kB :: Double
kB :: Double
kB = Double -> Double
forall a. Floating a => a -> a
sqrt (Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n Double -> Double -> Double
forall r. Additive r => r -> r -> r
+ Double
1) Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
* Double
2 Double -> Int -> Double
forall a b. (Num a, Integral b) => a -> b -> a
P.^ Int
n Double -> Double -> Double
forall r. Multiplicative r => r -> r -> r
* Integer -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Integer
Norm (Coefficient (Unipol Integer))
kA Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Integer
Coefficient (Unipol Integer)
lc)
loop :: Integer
-> [Int]
-> Int
-> Int
-> Unipol Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
loop !Integer
pk [Int]
alives !Int
count !Int
l !Unipol Integer
h Vector (Unipol Integer)
gs [Unipol Integer]
acc
| Int
2 Int -> Int -> Int
forall r. Multiplicative r => r -> r -> r
* Int
l Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
> Int
count =
if Unipol Integer
h Unipol Integer -> Unipol Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Unipol Integer
forall r. Unital r => r
one then [Unipol Integer]
acc else Unipol Integer
h Unipol Integer -> [Unipol Integer] -> [Unipol Integer]
forall a. a -> [a] -> [a]
: [Unipol Integer]
acc
| Bool
otherwise = Integer
-> [Int]
-> Int
-> Int
-> [([Int], [Int])]
-> Unipol Integer
-> Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
go Integer
pk [Int]
alives Int
count Int
l
(Int -> [Int] -> [([Int], [Int])]
combinations Int
l [Int]
alives) Unipol Integer
h (Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
h) Vector (Unipol Integer)
gs [Unipol Integer]
acc
go :: Integer
-> [Int]
-> Int
-> Int
-> [([Int], [Int])]
-> Unipol Integer
-> Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
go Integer
pk [Int]
alives !Int
count Int
l [] Unipol Integer
f' Integer
_ Vector (Unipol Integer)
pols [Unipol Integer]
acc = Integer
-> [Int]
-> Int
-> Int
-> Unipol Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
loop Integer
pk [Int]
alives Int
count (Int
l Int -> Int -> Int
forall r. Additive r => r -> r -> r
+ Int
1) Unipol Integer
f' Vector (Unipol Integer)
pols [Unipol Integer]
acc
go Integer
pk [Int]
alives !Int
count Int
l (([Int]
gs, [Int]
hs) : [([Int], [Int])]
cands) Unipol Integer
f' Integer
b Vector (Unipol Integer)
pols [Unipol Integer]
acc
| Unipol Integer -> Norm (Coefficient (Unipol Integer))
forall poly.
(IsPolynomial poly, Normed (Coefficient poly),
Monoidal (Norm (Coefficient poly))) =>
poly -> Norm (Coefficient poly)
oneNorm Unipol Integer
g Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Unipol Integer -> Norm (Coefficient (Unipol Integer))
forall poly.
(IsPolynomial poly, Normed (Coefficient poly),
Monoidal (Norm (Coefficient poly))) =>
poly -> Norm (Coefficient poly)
oneNorm Unipol Integer
h Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Double -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
floor Double
kB =
Integer
-> [Int]
-> Int
-> Int
-> Unipol Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
loop Integer
pk [Int]
hs (Int
count Int -> Int -> Int
forall r. Group r => r -> r -> r
- Int
l) Int
l (Unipol Integer -> Unipol Integer
forall poly.
(Euclidean (Coefficient poly), IsPolynomial poly) =>
poly -> poly
pp Unipol Integer
h) Vector (Unipol Integer)
pols (Unipol Integer -> Unipol Integer
forall poly.
(Euclidean (Coefficient poly), IsPolynomial poly) =>
poly -> poly
pp Unipol Integer
g Unipol Integer -> [Unipol Integer] -> [Unipol Integer]
forall a. a -> [a] -> [a]
: [Unipol Integer]
acc)
| Bool
otherwise = Integer
-> [Int]
-> Int
-> Int
-> [([Int], [Int])]
-> Unipol Integer
-> Integer
-> Vector (Unipol Integer)
-> [Unipol Integer]
-> [Unipol Integer]
go Integer
pk [Int]
alives Int
count Int
l [([Int], [Int])]
cands Unipol Integer
f' Integer
b Vector (Unipol Integer)
pols [Unipol Integer]
acc
where
g :: Unipol Integer
g = Integer -> Unipol Integer -> Unipol Integer
normalizeMod Integer
pk (Unipol Integer -> Unipol Integer)
-> Unipol Integer -> Unipol Integer
forall a b. (a -> b) -> a -> b
$ Integer
b Integer -> Unipol Integer -> Unipol Integer
forall r m. LeftModule r m => r -> m -> m
.* Mult (Unipol Integer) -> Unipol Integer
forall a. Mult a -> a
runMult ((Int -> Mult (Unipol Integer)) -> [Int] -> Mult (Unipol Integer)
forall (t :: * -> *) m a.
(Foldable t, Monoid m) =>
(a -> m) -> t a -> m
foldMap (Unipol Integer -> Mult (Unipol Integer)
forall a. a -> Mult a
Mult (Unipol Integer -> Mult (Unipol Integer))
-> (Int -> Unipol Integer) -> Int -> Mult (Unipol Integer)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector (Unipol Integer) -> Int -> Unipol Integer
forall a. Vector a -> Int -> a
V.unsafeIndex Vector (Unipol Integer)
pols) [Int]
gs)
h :: Unipol Integer
h = Integer -> Unipol Integer -> Unipol Integer
normalizeMod Integer
pk (Unipol Integer -> Unipol Integer)
-> Unipol Integer -> Unipol Integer
forall a b. (a -> b) -> a -> b
$ Integer
b Integer -> Unipol Integer -> Unipol Integer
forall r m. LeftModule r m => r -> m -> m
.* Mult (Unipol Integer) -> Unipol Integer
forall a. Mult a -> a
runMult ((Int -> Mult (Unipol Integer)) -> [Int] -> Mult (Unipol Integer)
forall (t :: * -> *) m a.
(Foldable t, Monoid m) =>
(a -> m) -> t a -> m
foldMap (Unipol Integer -> Mult (Unipol Integer)
forall a. a -> Mult a
Mult (Unipol Integer -> Mult (Unipol Integer))
-> (Int -> Unipol Integer) -> Int -> Mult (Unipol Integer)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector (Unipol Integer) -> Int -> Unipol Integer
forall a. Vector a -> Int -> a
V.unsafeIndex Vector (Unipol Integer)
pols) [Int]
hs)
henselStep
:: (Eq r, Euclidean r)
=> r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> (Unipol r, Unipol r, Unipol r, Unipol r)
henselStep :: r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> (Unipol r, Unipol r, Unipol r, Unipol r)
henselStep r
m Unipol r
f Unipol r
g Unipol r
h Unipol r
s Unipol r
t =
let modCoeff :: Unipol r -> Unipol r
modCoeff = (r -> r) -> Unipol r -> Unipol r
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (r -> r -> r
forall d. Euclidean d => d -> d -> d
`rem` r
mr -> Natural -> r
forall r. Unital r => r -> Natural -> r
^Natural
2)
divModSq :: Unipol r -> Unipol r -> (Unipol r, Unipol r)
divModSq Unipol r
u Unipol r
v =
(r -> Fraction r) -> Unipol r -> Unipol (Fraction r)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (r -> r -> Fraction r
forall d. GCDDomain d => d -> d -> Fraction d
F.% r
forall r. Unital r => r
one) Unipol r
u Unipol (Fraction r)
-> Unipol (Fraction r)
-> (Unipol (Fraction r), Unipol (Fraction r))
forall d. Euclidean d => d -> d -> (d, d)
`divide` (r -> Fraction r) -> Unipol r -> Unipol (Fraction r)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (r -> r -> Fraction r
forall d. GCDDomain d => d -> d -> Fraction d
F.% r
forall r. Unital r => r
one) Unipol r
v
(Unipol (Fraction r), Unipol (Fraction r))
-> ((Unipol (Fraction r), Unipol (Fraction r))
-> (Unipol r, Unipol r))
-> (Unipol r, Unipol r)
forall a b. a -> (a -> b) -> b
& (Unipol (Fraction r) -> Identity (Unipol r))
-> (Unipol (Fraction r), Unipol (Fraction r))
-> Identity (Unipol r, Unipol r)
forall (r :: * -> * -> *) a b.
Bitraversable r =>
Traversal (r a a) (r b b) a b
both ((Unipol (Fraction r) -> Identity (Unipol r))
-> (Unipol (Fraction r), Unipol (Fraction r))
-> Identity (Unipol r, Unipol r))
-> (Unipol (Fraction r) -> Unipol r)
-> (Unipol (Fraction r), Unipol (Fraction r))
-> (Unipol r, Unipol r)
forall s t a b. ASetter s t a b -> (a -> b) -> s -> t
%~ (Fraction r -> r) -> Unipol (Fraction r) -> Unipol r
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Maybe r -> r
forall a. HasCallStack => Maybe a -> a
fromJust (Maybe r -> r) -> (Fraction r -> Maybe r) -> Fraction r -> r
forall b c a. (b -> c) -> (a -> b) -> a -> c
. r -> Fraction r -> Maybe r
forall s. (Euclidean s, Eq s) => s -> Fraction s -> Maybe s
modFraction (r
mr -> Natural -> r
forall r. Unital r => r -> Natural -> r
^Natural
2))
e :: Unipol r
e = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
f Unipol r -> Unipol r -> Unipol r
forall r. Group r => r -> r -> r
- Unipol r
g Unipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
* Unipol r
h
(Unipol r
q, Unipol r
r) = Unipol r -> Unipol r -> (Unipol r, Unipol r)
divModSq (Unipol r
sUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
e) Unipol r
h
g' :: Unipol r
g' = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
g Unipol r -> Unipol r -> Unipol r
forall r. Additive r => r -> r -> r
+ Unipol r
t Unipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
* Unipol r
e Unipol r -> Unipol r -> Unipol r
forall r. Additive r => r -> r -> r
+ Unipol r
q Unipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
* Unipol r
g
h' :: Unipol r
h' = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
h Unipol r -> Unipol r -> Unipol r
forall r. Additive r => r -> r -> r
+ Unipol r
r
b :: Unipol r
b = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
sUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
g' Unipol r -> Unipol r -> Unipol r
forall r. Additive r => r -> r -> r
+ Unipol r
tUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
h' Unipol r -> Unipol r -> Unipol r
forall r. Group r => r -> r -> r
- Unipol r
forall r. Unital r => r
one
(Unipol r
c, Unipol r
d) = Unipol r -> Unipol r -> (Unipol r, Unipol r)
divModSq (Unipol r
sUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
b) Unipol r
h'
s' :: Unipol r
s' = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
s Unipol r -> Unipol r -> Unipol r
forall r. Group r => r -> r -> r
- Unipol r
d
t' :: Unipol r
t' = Unipol r -> Unipol r
modCoeff (Unipol r -> Unipol r) -> Unipol r -> Unipol r
forall a b. (a -> b) -> a -> b
$ Unipol r
t Unipol r -> Unipol r -> Unipol r
forall r. Group r => r -> r -> r
- Unipol r
tUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
b Unipol r -> Unipol r -> Unipol r
forall r. Group r => r -> r -> r
- Unipol r
cUnipol r -> Unipol r -> Unipol r
forall r. Multiplicative r => r -> r -> r
*Unipol r
g'
in (Unipol r
g', Unipol r
h', Unipol r
s', Unipol r
t')
repeatHensel
:: Integer -> Int
-> Unipol Integer
-> Unipol Integer -> Unipol Integer
-> Unipol Integer -> Unipol Integer
-> (Unipol Integer, Unipol Integer, Unipol Integer, Unipol Integer)
repeatHensel :: Integer
-> Int
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> (Unipol Integer, Unipol Integer, Unipol Integer, Unipol Integer)
repeatHensel !Integer
m Int
0 Unipol Integer
_ Unipol Integer
g Unipol Integer
h Unipol Integer
s Unipol Integer
t = (Integer -> Unipol Integer -> Unipol Integer
normalizeMod Integer
m Unipol Integer
g, Integer -> Unipol Integer -> Unipol Integer
normalizeMod Integer
m Unipol Integer
h, Unipol Integer
s, Unipol Integer
t)
repeatHensel !Integer
m Int
n Unipol Integer
f Unipol Integer
g Unipol Integer
h Unipol Integer
s Unipol Integer
t =
let (Unipol Integer
g', Unipol Integer
h', Unipol Integer
s', Unipol Integer
t') = Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> (Unipol Integer, Unipol Integer, Unipol Integer, Unipol Integer)
forall r.
(Eq r, Euclidean r) =>
r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> Unipol r
-> (Unipol r, Unipol r, Unipol r, Unipol r)
henselStep Integer
m Unipol Integer
f Unipol Integer
g Unipol Integer
h Unipol Integer
s Unipol Integer
t
in Integer
-> Int
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> (Unipol Integer, Unipol Integer, Unipol Integer, Unipol Integer)
repeatHensel (Integer
m Integer -> Natural -> Integer
forall r. Unital r => r -> Natural -> r
^ Natural
2) (Int
n Int -> Int -> Int
forall r. Group r => r -> r -> r
- Int
1) Unipol Integer
f Unipol Integer
g' Unipol Integer
h' Unipol Integer
s' Unipol Integer
t'
multiHensel :: Natural
-> Int
-> Unipol Integer
-> [Unipol Integer]
-> [Unipol Integer]
{-# INLINE multiHensel #-}
multiHensel :: Natural
-> Int -> Unipol Integer -> [Unipol Integer] -> [Unipol Integer]
multiHensel (Natural -> Integer
forall a. Integral a => a -> Integer
toInteger -> Integer
p) Int
n = \Unipol Integer
f -> DList (Unipol Integer) -> [Unipol Integer]
forall a. DList a -> [a]
DL.toList (DList (Unipol Integer) -> [Unipol Integer])
-> ([Unipol Integer] -> DList (Unipol Integer))
-> [Unipol Integer]
-> [Unipol Integer]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Unipol Integer -> Vector (Unipol Integer) -> DList (Unipol Integer)
go Unipol Integer
f (Vector (Unipol Integer) -> DList (Unipol Integer))
-> ([Unipol Integer] -> Vector (Unipol Integer))
-> [Unipol Integer]
-> DList (Unipol Integer)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Unipol Integer] -> Vector (Unipol Integer)
forall a. [a] -> Vector a
V.fromList
where
go :: Unipol Integer -> Vector (Unipol Integer) -> DList (Unipol Integer)
go Unipol Integer
f Vector (Unipol Integer)
gs
| Int
len Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 = DList (Unipol Integer)
forall a. DList a
DL.empty
| Int
len Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 = Unipol Integer -> DList (Unipol Integer)
forall a. a -> DList a
DL.singleton (Unipol Integer -> DList (Unipol Integer))
-> Unipol Integer -> DList (Unipol Integer)
forall a b. (a -> b) -> a -> b
$
let q :: Integer
q = Integer
p Integer -> Int -> Integer
forall a b. (Num a, Integral b) => a -> b -> a
P.^ Int
n
(Integer
u,Integer
_,Integer
bInv) = Integer -> Integer -> (Integer, Integer, Integer)
forall d. PID d => d -> d -> (d, d, d)
egcd Integer
q (Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f)
in (Integer -> Integer) -> Unipol Integer -> Unipol Integer
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`mod` Integer
q)
(Unipol Integer -> Unipol Integer)
-> Unipol Integer -> Unipol Integer
forall a b. (a -> b) -> a -> b
$ (Integer
bInv Integer -> Integer -> Integer
forall d. Euclidean d => d -> d -> d
`quot` Integer
u) Integer -> Unipol Integer -> Unipol Integer
forall r m. Module (Scalar r) m => r -> m -> m
.*. Unipol Integer
f
| Bool
otherwise = Integer
-> (forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> DList (Unipol Integer))
-> DList (Unipol Integer)
forall a.
Integer -> (forall (p :: Nat). KnownNat p => Proxy (F p) -> a) -> a
reifyPrimeField Integer
p ((forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> DList (Unipol Integer))
-> DList (Unipol Integer))
-> (forall (p :: Nat).
KnownNat p =>
Proxy (F p) -> DList (Unipol Integer))
-> DList (Unipol Integer)
forall a b. (a -> b) -> a -> b
$ \Proxy (F p)
fp ->
let k :: Int
k = Int
len Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2
d :: Int
d = Int -> Int
ceilingLogBase2 Int
n
(Vector (Unipol Integer)
ls, Vector (Unipol Integer)
rs) = Int
-> Vector (Unipol Integer)
-> (Vector (Unipol Integer), Vector (Unipol Integer))
forall a. Int -> Vector a -> (Vector a, Vector a)
V.splitAt Int
k Vector (Unipol Integer)
gs
(Unipol Integer
g0, Unipol Integer
h0) = (Unipol Integer -> Coefficient (Unipol Integer)
forall poly. IsOrderedPolynomial poly => poly -> Coefficient poly
leadingCoeff Unipol Integer
f Integer -> Unipol Integer -> Unipol Integer
forall r m. Module (Scalar r) m => r -> m -> m
.*. Vector (Unipol Integer) -> Unipol Integer
forall (f :: * -> *) r. (Foldable f, Unital r) => f r -> r
product Vector (Unipol Integer)
ls, Vector (Unipol Integer) -> Unipol Integer
forall (f :: * -> *) r. (Foldable f, Unital r) => f r -> r
product Vector (Unipol Integer)
rs)
(Unipol Integer, Unipol Integer)
-> ((Unipol Integer, Unipol Integer)
-> (Unipol Integer, Unipol Integer))
-> (Unipol Integer, Unipol Integer)
forall a b. a -> (a -> b) -> b
& (Unipol Integer -> Identity (Unipol Integer))
-> (Unipol Integer, Unipol Integer)
-> Identity (Unipol Integer, Unipol Integer)
forall (r :: * -> * -> *) a b.
Bitraversable r =>
Traversal (r a a) (r b b) a b
both ((Unipol Integer -> Identity (Unipol Integer))
-> (Unipol Integer, Unipol Integer)
-> Identity (Unipol Integer, Unipol Integer))
-> (Unipol Integer -> Unipol Integer)
-> (Unipol Integer, Unipol Integer)
-> (Unipol Integer, Unipol Integer)
forall s t a b. ASetter s t a b -> (a -> b) -> s -> t
%~ (Integer -> Integer) -> Unipol Integer -> Unipol Integer
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`mod` Integer
p)
(Unipol (F p)
_, Unipol (F p)
s0, Unipol (F p)
t0) = Unipol (F p)
-> Unipol (F p) -> (Unipol (F p), Unipol (F p), Unipol (F p))
forall d. PID d => d -> d -> (d, d, d)
egcd
((Integer -> F p) -> Unipol Integer -> Unipol (F p)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Proxy (F p) -> Integer -> F p
forall k (proxy :: * -> *) (p :: k).
Reifies p Integer =>
proxy (F p) -> Integer -> F p
modNat' Proxy (F p)
fp) Unipol Integer
g0)
((Integer -> F p) -> Unipol Integer -> Unipol (F p)
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Proxy (F p) -> Integer -> F p
forall k (proxy :: * -> *) (p :: k).
Reifies p Integer =>
proxy (F p) -> Integer -> F p
modNat' Proxy (F p)
fp) Unipol Integer
h0)
(Unipol Integer
s, Unipol Integer
t) = (Unipol (F p)
s0, Unipol (F p)
t0) (Unipol (F p), Unipol (F p))
-> ((Unipol (F p), Unipol (F p))
-> (Unipol Integer, Unipol Integer))
-> (Unipol Integer, Unipol Integer)
forall a b. a -> (a -> b) -> b
& (Unipol (F p) -> Identity (Unipol Integer))
-> (Unipol (F p), Unipol (F p))
-> Identity (Unipol Integer, Unipol Integer)
forall (r :: * -> * -> *) a b.
Bitraversable r =>
Traversal (r a a) (r b b) a b
both ((Unipol (F p) -> Identity (Unipol Integer))
-> (Unipol (F p), Unipol (F p))
-> Identity (Unipol Integer, Unipol Integer))
-> (Unipol (F p) -> Unipol Integer)
-> (Unipol (F p), Unipol (F p))
-> (Unipol Integer, Unipol Integer)
forall s t a b. ASetter s t a b -> (a -> b) -> s -> t
%~ (F p -> Integer) -> Unipol (F p) -> Unipol Integer
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol F p -> Integer
forall k (p :: k). F p -> Integer
naturalRepr
(Unipol Integer
g, Unipol Integer
h, Unipol Integer
_, Unipol Integer
_) = Integer
-> Int
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> Unipol Integer
-> (Unipol Integer, Unipol Integer, Unipol Integer, Unipol Integer)
repeatHensel Integer
p Int
d Unipol Integer
f Unipol Integer
g0 Unipol Integer
h0 Unipol Integer
s Unipol Integer
t
in Unipol Integer -> Vector (Unipol Integer) -> DList (Unipol Integer)
go Unipol Integer
g Vector (Unipol Integer)
ls DList (Unipol Integer)
-> DList (Unipol Integer) -> DList (Unipol Integer)
forall a. Semigroup a => a -> a -> a
<> Unipol Integer -> Vector (Unipol Integer) -> DList (Unipol Integer)
go Unipol Integer
h Vector (Unipol Integer)
rs
where len :: Int
len = Vector (Unipol Integer) -> Int
forall a. Vector a -> Int
V.length Vector (Unipol Integer)
gs
recipMod :: (Euclidean a, Eq a) => a -> a -> Maybe a
recipMod :: a -> a -> Maybe a
recipMod a
m a
u =
let (a
a, a
r, a
_) : [(a, a, a)]
_ = a -> a -> [(a, a, a)]
forall d. Euclidean d => d -> d -> [(d, d, d)]
euclid a
u a
m
in if a
a a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
forall r. Unital r => r
one
then a -> Maybe a
forall a. a -> Maybe a
Just a
r else Maybe a
forall a. Maybe a
Nothing
modFraction :: (Euclidean s, Eq s) => s -> Fraction s -> Maybe s
modFraction :: s -> Fraction s -> Maybe s
modFraction s
m Fraction s
d = ((Fraction s -> s
forall t. Fraction t -> t
numerator Fraction s
d s -> s -> s
forall d. Euclidean d => d -> d -> d
`rem` s
m) s -> s -> s
forall r. Multiplicative r => r -> r -> r
*) (s -> s) -> Maybe s -> Maybe s
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> s -> s -> Maybe s
forall a. (Euclidean a, Eq a) => a -> a -> Maybe a
recipMod s
m (Fraction s -> s
forall t. Fraction t -> t
denominator Fraction s
d)
combinations :: Int -> [Int] -> [([Int], [Int])]
combinations :: Int -> [Int] -> [([Int], [Int])]
combinations Int
n = FMList ([Int], [Int]) -> [([Int], [Int])]
forall (t :: * -> *) a. Foldable t => t a -> [a]
FML.toList (FMList ([Int], [Int]) -> [([Int], [Int])])
-> ([Int] -> FMList ([Int], [Int])) -> [Int] -> [([Int], [Int])]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> [Int] -> FMList ([Int], [Int])
forall r a. (Eq r, Num r, Group r) => r -> [a] -> FMList ([a], [a])
go Int
n
where
go :: r -> [a] -> FMList ([a], [a])
go r
0 [a]
xs = ([a], [a]) -> FMList ([a], [a])
forall a. a -> FMList a
FML.singleton ([], [a]
xs)
go r
_ [] = FMList ([a], [a])
forall (f :: * -> *) a. Alternative f => f a
FML.empty
go r
k (a
x:[a]
xs) =
let present :: FMList ([a], [a])
present = r -> [a] -> FMList ([a], [a])
go (r
k r -> r -> r
forall r. Group r => r -> r -> r
- r
1) [a]
xs
absent :: FMList ([a], [a])
absent = r -> [a] -> FMList ([a], [a])
go r
k [a]
xs
in (([a], [a]) -> ([a], [a]))
-> FMList ([a], [a]) -> FMList ([a], [a])
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (([a] -> [a]) -> ([a], [a]) -> ([a], [a])
forall (a :: * -> * -> *) b c d.
Arrow a =>
a b c -> a (b, d) (c, d)
first (([a] -> [a]) -> ([a], [a]) -> ([a], [a]))
-> ([a] -> [a]) -> ([a], [a]) -> ([a], [a])
forall a b. (a -> b) -> a -> b
$ (:) a
x) FMList ([a], [a])
present
FMList ([a], [a]) -> FMList ([a], [a]) -> FMList ([a], [a])
forall a. Semigroup a => a -> a -> a
<> (([a], [a]) -> ([a], [a]))
-> FMList ([a], [a]) -> FMList ([a], [a])
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (([a] -> [a]) -> ([a], [a]) -> ([a], [a])
forall (a :: * -> * -> *) b c d.
Arrow a =>
a b c -> a (d, b) (d, c)
second (([a] -> [a]) -> ([a], [a]) -> ([a], [a]))
-> ([a] -> [a]) -> ([a], [a]) -> ([a], [a])
forall a b. (a -> b) -> a -> b
$ (:) a
x) FMList ([a], [a])
absent
normalizeMod :: Integer -> Unipol Integer -> Unipol Integer
normalizeMod :: Integer -> Unipol Integer -> Unipol Integer
normalizeMod Integer
p = (Integer -> Integer) -> Unipol Integer -> Unipol Integer
forall b a. DecidableZero b => (a -> b) -> Unipol a -> Unipol b
mapCoeffUnipol (Integer -> Integer -> Integer
forall r. Group r => r -> r -> r
subtract Integer
half (Integer -> Integer) -> (Integer -> Integer) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`mod` Integer
p) (Integer -> Integer) -> (Integer -> Integer) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer -> Integer -> Integer
forall r. Additive r => r -> r -> r
+ Integer
half))
where half :: Integer
half = Integer
p Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`div` Integer
2
isSquareFree :: forall poly. (IsOrderedPolynomial poly, GCDDomain poly) => poly -> Bool
isSquareFree :: poly -> Bool
isSquareFree poly
f = (poly
f poly -> poly -> poly
forall d. GCDDomain d => d -> d -> d
`gcd` Ordinal (Arity poly) -> poly -> poly
forall poly.
IsOrderedPolynomial poly =>
Ordinal (Arity poly) -> poly -> poly
diff Ordinal (Arity poly)
0 poly
f) poly -> poly -> Bool
forall a. Eq a => a -> a -> Bool
== poly
forall r. Unital r => r
one
isReducible
:: forall k. (FiniteField k, Eq k)
=> Unipol k -> Bool
isReducible :: Unipol k -> Bool
isReducible Unipol k
f
| Natural
n Natural -> Natural -> Bool
forall a. Eq a => a -> a -> Bool
== Natural
0 = Bool
True
| Natural
n Natural -> Natural -> Bool
forall a. Eq a => a -> a -> Bool
== Natural
1 = Bool
False
| Bool
otherwise =
let divisors :: [Natural]
divisors = ((Prime Natural, Word) -> Natural)
-> [(Prime Natural, Word)] -> [Natural]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (Prime Natural -> Natural
forall a. Prime a -> a
PRIMES.unPrime (Prime Natural -> Natural)
-> ((Prime Natural, Word) -> Prime Natural)
-> (Prime Natural, Word)
-> Natural
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Prime Natural, Word) -> Prime Natural
forall a b. (a, b) -> a
fst) ([(Prime Natural, Word)] -> [Natural])
-> [(Prime Natural, Word)] -> [Natural]
forall a b. (a -> b) -> a -> b
$ Natural -> [(Prime Natural, Word)]
forall a. UniqueFactorisation a => a -> [(Prime a, Word)]
PRIMES.factorise Natural
n
q :: Natural
q = Proxy k -> Natural
forall k (proxy :: * -> *). FiniteField k => proxy k -> Natural
order (Proxy k -> Natural) -> Proxy k -> Natural
forall a b. (a -> b) -> a -> b
$ Proxy k
forall k (t :: k). Proxy t
Proxy @k
xqn :: Unipol k
xqn = Unipol k -> Natural -> Unipol k -> Unipol k
forall poly. Euclidean poly => poly -> Natural -> poly -> poly
modPow IsLabel "x" (Unipol k)
Unipol k
#x (Natural
qNatural -> Natural -> Natural
forall r. Unital r => r -> Natural -> r
^Natural
n) Unipol k
f
test :: Natural -> Bool
test Natural
t =
let b :: Unipol k
b = Unipol k -> Natural -> Unipol k -> Unipol k
forall poly. Euclidean poly => poly -> Natural -> poly -> poly
modPow IsLabel "x" (Unipol k)
Unipol k
#x (Natural
q Natural -> Natural -> Natural
forall r. Unital r => r -> Natural -> r
^ (Natural
n Natural -> Natural -> Natural
forall a. Integral a => a -> a -> a
`div` Natural
t)) Unipol k
f
in Unipol k -> Unipol k -> Unipol k
forall d. GCDDomain d => d -> d -> d
gcd (Unipol k
b Unipol k -> Unipol k -> Unipol k
forall r. Group r => r -> r -> r
- IsLabel "x" (Unipol k)
Unipol k
#x) Unipol k
f Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= Unipol k
forall r. Unital r => r
one
in Unipol k
xqn Unipol k -> Unipol k -> Bool
forall a. Eq a => a -> a -> Bool
/= IsLabel "x" (Unipol k)
Unipol k
#x Bool -> Bool -> Bool
|| (Natural -> Bool) -> [Natural] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
any Natural -> Bool
test [Natural]
divisors
where
n :: Natural
n = (Integral Int, Num Natural) => Int -> Natural
forall a b. (Integral a, Num b) => a -> b
fromIntegral @_ @Natural (Int -> Natural) -> Int -> Natural
forall a b. (a -> b) -> a -> b
$ Unipol k -> Int
forall poly. IsPolynomial poly => poly -> Int
totalDegree' Unipol k
f