{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE TypeSynonymInstances #-}
{-# OPTIONS_GHC -Wall -fno-warn-orphans #-}
{-# OPTIONS_HADDOCK show-extensions #-}
module ToySolver.Data.Polynomial.Factorization.FiniteField
( factor
, sqfree
, berlekamp
, basisOfBerlekampSubalgebra
) where
import Control.Exception (assert)
import Data.FiniteField
import Data.List
import Data.Ord
import Data.Set (Set)
import qualified Data.Set as Set
import GHC.TypeLits
import ToySolver.Data.Polynomial.Base (Polynomial, UPolynomial, X (..))
import qualified ToySolver.Data.Polynomial.Base as P
import qualified ToySolver.Data.Polynomial.GroebnerBasis as GB
instance KnownNat p => P.Factor (UPolynomial (PrimeField p)) where
factor :: UPolynomial (PrimeField p)
-> [(UPolynomial (PrimeField p), Integer)]
factor = forall k.
(Ord k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
factor
instance KnownNat p => P.SQFree (UPolynomial (PrimeField p)) where
sqfree :: UPolynomial (PrimeField p)
-> [(UPolynomial (PrimeField p), Integer)]
sqfree = forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree
factor :: forall k. (Ord k, FiniteField k) => UPolynomial k -> [(UPolynomial k, Integer)]
factor :: forall k.
(Ord k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
factor UPolynomial k
f = do
(UPolynomial k
g,Integer
n) <- forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree UPolynomial k
f
if forall t. Degree t => t -> Integer
P.deg UPolynomial k
g forall a. Ord a => a -> a -> Bool
> Integer
0
then do
UPolynomial k
h <- forall k.
(Eq k, Ord k, FiniteField k) =>
UPolynomial k -> [UPolynomial k]
berlekamp UPolynomial k
g
forall (m :: * -> *) a. Monad m => a -> m a
return (UPolynomial k
h,Integer
n)
else do
forall (m :: * -> *) a. Monad m => a -> m a
return (UPolynomial k
g,Integer
n)
sqfree :: forall k. (Eq k, FiniteField k) => UPolynomial k -> [(UPolynomial k, Integer)]
sqfree :: forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree UPolynomial k
f
| k
c forall a. Eq a => a -> a -> Bool
== k
1 = forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' UPolynomial k
f
| Bool
otherwise = (forall k v. (Eq k, Num k) => k -> Polynomial k v
P.constant k
c, Integer
1) forall a. a -> [a] -> [a]
: forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' (forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff (forall a. Fractional a => a -> a -> a
/k
c) UPolynomial k
f)
where
c :: k
c = forall k v. Num k => MonomialOrder v -> Polynomial k v -> k
P.lc MonomialOrder X
P.nat UPolynomial k
f
sqfree' :: forall k. (Eq k, FiniteField k) => UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' :: forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' UPolynomial k
0 = []
sqfree' UPolynomial k
f
| UPolynomial k
g forall a. Eq a => a -> a -> Bool
== UPolynomial k
0 = [(UPolynomial k
h, Integer
nforall a. Num a => a -> a -> a
*Integer
p) | (UPolynomial k
h,Integer
n) <- forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' (forall k. (Eq k, FiniteField k) => UPolynomial k -> UPolynomial k
polyPthRoot UPolynomial k
f)]
| Bool
otherwise = forall {k}.
(FiniteField k, Eq k) =>
Integer
-> UPolynomial k
-> UPolynomial k
-> [(UPolynomial k, Integer)]
-> [(UPolynomial k, Integer)]
go Integer
1 UPolynomial k
c0 UPolynomial k
w0 []
where
p :: Integer
p = forall k. FiniteField k => k -> Integer
char (forall a. HasCallStack => a
undefined :: k)
g :: UPolynomial k
g = forall k v.
(Eq k, Num k, Ord v) =>
Polynomial k v -> v -> Polynomial k v
P.deriv UPolynomial k
f X
X
c0 :: UPolynomial k
c0 = forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.gcd UPolynomial k
f UPolynomial k
g
w0 :: UPolynomial k
w0 = forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.div UPolynomial k
f UPolynomial k
c0
go :: Integer
-> UPolynomial k
-> UPolynomial k
-> [(UPolynomial k, Integer)]
-> [(UPolynomial k, Integer)]
go !Integer
i UPolynomial k
c UPolynomial k
w ![(UPolynomial k, Integer)]
result
| UPolynomial k
w forall a. Eq a => a -> a -> Bool
== UPolynomial k
1 =
if UPolynomial k
c forall a. Eq a => a -> a -> Bool
== UPolynomial k
1
then [(UPolynomial k, Integer)]
result
else [(UPolynomial k, Integer)]
result forall a. [a] -> [a] -> [a]
++ [(UPolynomial k
h, Integer
nforall a. Num a => a -> a -> a
*Integer
p) | (UPolynomial k
h,Integer
n) <- forall k.
(Eq k, FiniteField k) =>
UPolynomial k -> [(UPolynomial k, Integer)]
sqfree' (forall k. (Eq k, FiniteField k) => UPolynomial k -> UPolynomial k
polyPthRoot UPolynomial k
c)]
| Bool
otherwise = Integer
-> UPolynomial k
-> UPolynomial k
-> [(UPolynomial k, Integer)]
-> [(UPolynomial k, Integer)]
go (Integer
iforall a. Num a => a -> a -> a
+Integer
1) UPolynomial k
c' UPolynomial k
w' [(UPolynomial k, Integer)]
result'
where
y :: UPolynomial k
y = forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.gcd UPolynomial k
w UPolynomial k
c
z :: UPolynomial k
z = UPolynomial k
w forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` UPolynomial k
y
c' :: UPolynomial k
c' = UPolynomial k
c forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` UPolynomial k
y
w' :: UPolynomial k
w' = UPolynomial k
y
result' :: [(UPolynomial k, Integer)]
result' = [(UPolynomial k
z,Integer
i) | UPolynomial k
z forall a. Eq a => a -> a -> Bool
/= UPolynomial k
1] forall a. [a] -> [a] -> [a]
++ [(UPolynomial k, Integer)]
result
polyPthRoot :: forall k. (Eq k, FiniteField k) => UPolynomial k -> UPolynomial k
polyPthRoot :: forall k. (Eq k, FiniteField k) => UPolynomial k -> UPolynomial k
polyPthRoot Polynomial k X
f = forall a. HasCallStack => Bool -> a -> a
assert (forall k v.
(Eq k, Num k, Ord v) =>
Polynomial k v -> v -> Polynomial k v
P.deriv Polynomial k X
f X
X forall a. Eq a => a -> a -> Bool
== Polynomial k X
0) forall a b. (a -> b) -> a -> b
$
forall k v. (Eq k, Num k, Ord v) => [Term k v] -> Polynomial k v
P.fromTerms [(forall k. FiniteField k => k -> k
pthRoot k
c, forall {t}. Degree t => t -> Monomial X
g Monomial X
mm) | (k
c,Monomial X
mm) <- forall k v. Polynomial k v -> [Term k v]
P.terms Polynomial k X
f]
where
p :: Integer
p = forall k. FiniteField k => k -> Integer
char (forall a. HasCallStack => a
undefined :: k)
g :: t -> Monomial X
g t
mm = forall a v. Var a v => v -> a
P.var X
X forall v. Ord v => Monomial v -> Integer -> Monomial v
`P.mpow` (forall t. Degree t => t -> Integer
P.deg t
mm forall a. Integral a => a -> a -> a
`div` Integer
p)
berlekamp :: forall k. (Eq k, Ord k, FiniteField k) => UPolynomial k -> [UPolynomial k]
berlekamp :: forall k.
(Eq k, Ord k, FiniteField k) =>
UPolynomial k -> [UPolynomial k]
berlekamp UPolynomial k
f = Set (UPolynomial k) -> [UPolynomial k] -> [UPolynomial k]
go (forall a. a -> Set a
Set.singleton UPolynomial k
f) [UPolynomial k]
basis
where
go :: Set (UPolynomial k) -> [UPolynomial k] -> [UPolynomial k]
go :: Set (UPolynomial k) -> [UPolynomial k] -> [UPolynomial k]
go Set (UPolynomial k)
_ [] = forall a. HasCallStack => [Char] -> a
error forall a b. (a -> b) -> a -> b
$ [Char]
"ToySolver.Data.Polynomial.Factorization.FiniteField.berlekamp: should not happen"
go Set (UPolynomial k)
fs (UPolynomial k
b:[UPolynomial k]
bs)
| forall a. Set a -> Int
Set.size Set (UPolynomial k)
fs forall a. Eq a => a -> a -> Bool
== Int
r = forall a. Set a -> [a]
Set.toList Set (UPolynomial k)
fs
| Bool
otherwise = Set (UPolynomial k) -> [UPolynomial k] -> [UPolynomial k]
go (forall (f :: * -> *) a. (Foldable f, Ord a) => f (Set a) -> Set a
Set.unions [UPolynomial k -> Set (UPolynomial k)
func UPolynomial k
fi | UPolynomial k
fi <- forall a. Set a -> [a]
Set.toList Set (UPolynomial k)
fs]) [UPolynomial k]
bs
where
func :: UPolynomial k -> Set (UPolynomial k)
func UPolynomial k
fi = forall a. Ord a => [a] -> Set a
Set.fromList forall a b. (a -> b) -> a -> b
$ [UPolynomial k]
hs2 forall a. [a] -> [a] -> [a]
++ [UPolynomial k]
hs1
where
hs1 :: [UPolynomial k]
hs1 = [UPolynomial k
h | k
k <- forall k. FiniteField k => [k]
allValues, let h :: UPolynomial k
h = forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.gcd UPolynomial k
fi (UPolynomial k
b forall a. Num a => a -> a -> a
- forall k v. (Eq k, Num k) => k -> Polynomial k v
P.constant k
k), forall t. Degree t => t -> Integer
P.deg UPolynomial k
h forall a. Ord a => a -> a -> Bool
> Integer
0]
hs2 :: [UPolynomial k]
hs2 = if forall t. Degree t => t -> Integer
P.deg UPolynomial k
g forall a. Ord a => a -> a -> Bool
> Integer
0 then [UPolynomial k
g] else []
where
g :: UPolynomial k
g = UPolynomial k
fi forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial k]
hs1
basis :: [UPolynomial k]
basis = forall k.
(Ord k, FiniteField k) =>
UPolynomial k -> [UPolynomial k]
basisOfBerlekampSubalgebra UPolynomial k
f
r :: Int
r = forall (t :: * -> *) a. Foldable t => t a -> Int
length [UPolynomial k]
basis
basisOfBerlekampSubalgebra :: forall k. (Ord k, FiniteField k) => UPolynomial k -> [UPolynomial k]
basisOfBerlekampSubalgebra :: forall k.
(Ord k, FiniteField k) =>
UPolynomial k -> [UPolynomial k]
basisOfBerlekampSubalgebra UPolynomial k
f =
forall b a. Ord b => (a -> b) -> [a] -> [a]
sortOn (forall a. a -> Down a
Down forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall t. Degree t => t -> Integer
P.deg) forall a b. (a -> b) -> a -> b
$
forall a b. (a -> b) -> [a] -> [b]
map (forall r v.
(Eq r, Fractional r) =>
MonomialOrder v -> Polynomial r v -> Polynomial r v
P.toMonic MonomialOrder X
P.nat) forall a b. (a -> b) -> a -> b
$
[UPolynomial k]
basis
where
q :: Integer
q = forall k. FiniteField k => k -> Integer
order (forall a. HasCallStack => a
undefined :: k)
d :: Integer
d = forall t. Degree t => t -> Integer
P.deg UPolynomial k
f
x :: UPolynomial k
x = forall a v. Var a v => v -> a
P.var X
X
qs :: [UPolynomial k]
qs :: [UPolynomial k]
qs = [(UPolynomial k
xforall a b. (Num a, Integral b) => a -> b -> a
^(Integer
qforall a. Num a => a -> a -> a
*Integer
i)) forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.mod` UPolynomial k
f | Integer
i <- [Integer
0 .. Integer
d forall a. Num a => a -> a -> a
- Integer
1]]
gb :: [Polynomial k Int]
gb :: [Polynomial k Int]
gb = forall k v.
(Eq k, Fractional k, Ord k, Ord v) =>
MonomialOrder v -> [Polynomial k v] -> [Polynomial k v]
GB.basis forall v. Ord v => MonomialOrder v
P.grlex [Polynomial k Int
p3 | (Polynomial k Int
p3,Monomial X
_) <- forall k v. Polynomial k v -> [Term k v]
P.terms UPolynomial (Polynomial k Int)
p2]
p1 :: Polynomial k Int
p1 :: Polynomial k Int
p1 = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [forall a v. Var a v => v -> a
P.var Int
i forall a. Num a => a -> a -> a
* (forall k v2 v1.
(Eq k, Num k, Ord v2) =>
Polynomial k v1 -> (v1 -> Polynomial k v2) -> Polynomial k v2
P.subst UPolynomial k
qi (\X
X -> forall a v. Var a v => v -> a
P.var (-Int
1)) forall a. Num a => a -> a -> a
- (forall a v. Var a v => v -> a
P.var (-Int
1) forall a b. (Num a, Integral b) => a -> b -> a
^ Int
i)) | (Int
i, UPolynomial k
qi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [Int
0..] [UPolynomial k]
qs]
p2 :: UPolynomial (Polynomial k Int)
p2 :: UPolynomial (Polynomial k Int)
p2 = forall k v.
(Ord k, Num k, Ord v) =>
Polynomial k v -> v -> UPolynomial (Polynomial k v)
P.toUPolynomialOf Polynomial k Int
p1 (-Int
1)
es :: [(Int, Polynomial k Int)]
es = [(Int
i, forall k v.
(Eq k, Fractional k, Ord v) =>
MonomialOrder v
-> Polynomial k v -> [Polynomial k v] -> Polynomial k v
P.reduce forall v. Ord v => MonomialOrder v
P.grlex (forall a v. Var a v => v -> a
P.var Int
i) [Polynomial k Int]
gb) | Int
i <- [Int
0 .. forall a b. (Integral a, Num b) => a -> b
fromIntegral Integer
d forall a. Num a => a -> a -> a
- Int
1]]
vs1 :: [Int]
vs1 = [Int
i | (Int
i, Polynomial k Int
gi_def) <- [(Int, Polynomial k Int)]
es, Polynomial k Int
gi_def forall a. Eq a => a -> a -> Bool
== forall a v. Var a v => v -> a
P.var Int
i]
vs2 :: [(Int, Polynomial k Int)]
vs2 = [(Int
i, Polynomial k Int
gi_def) | (Int
i, Polynomial k Int
gi_def) <- [(Int, Polynomial k Int)]
es, Polynomial k Int
gi_def forall a. Eq a => a -> a -> Bool
/= forall a v. Var a v => v -> a
P.var Int
i]
basis :: [UPolynomial k]
basis :: [UPolynomial k]
basis = [ UPolynomial k
xforall a b. (Num a, Integral b) => a -> b -> a
^Int
i forall a. Num a => a -> a -> a
+ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [forall k v. (Eq k, Num k) => k -> Polynomial k v
P.constant (forall k v. Num k => (v -> k) -> Polynomial k v -> k
P.eval (forall {a} {a}. (Eq a, Num a) => a -> a -> a
delta Int
i) Polynomial k Int
gj_def) forall a. Num a => a -> a -> a
* UPolynomial k
xforall a b. (Num a, Integral b) => a -> b -> a
^Int
j | (Int
j, Polynomial k Int
gj_def) <- [(Int, Polynomial k Int)]
vs2] | Int
i <- [Int]
vs1 ]
where
delta :: a -> a -> a
delta a
i a
k
| a
kforall a. Eq a => a -> a -> Bool
==a
i = a
1
| Bool
otherwise = a
0