{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# OPTIONS_GHC -Wall #-}
{-# OPTIONS_HADDOCK show-extensions #-}
-----------------------------------------------------------------------------
-- |
-- Module      :  ToySolver.Data.Polynomial.Factorization.Hensel.Internal
-- Copyright   :  (c) Masahiro Sakai 2013-2014
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  non-portable
--
-- References:
--
-- * <http://www.math.kobe-u.ac.jp/Asir/ca.pdf>
--
-- * <http://www14.in.tum.de/konferenzen/Jass07/courses/1/Bulwahn/Buhlwahn_Paper.pdf>
--
-----------------------------------------------------------------------------
module ToySolver.Data.Polynomial.Factorization.Hensel.Internal
  ( hensel
  , cabook_proposition_5_10
  , cabook_proposition_5_11
  ) where

import Control.Exception (assert)
import Data.FiniteField
import Data.Proxy
import GHC.TypeLits

import ToySolver.Data.Polynomial.Base (UPolynomial)
import qualified ToySolver.Data.Polynomial.Base as P

-- import Text.PrettyPrint.HughesPJClass

hensel :: forall p. KnownNat p => UPolynomial Integer -> [UPolynomial (PrimeField p)] -> Integer -> [UPolynomial Integer]
hensel :: forall (p :: Nat).
KnownNat p =>
UPolynomial Integer
-> [UPolynomial (PrimeField p)] -> Integer -> [UPolynomial Integer]
hensel UPolynomial Integer
f [UPolynomial (PrimeField p)]
fs1 Integer
k
  | Integer
k forall a. Ord a => a -> a -> Bool
<= Integer
0    = forall a. (?callStack::CallStack) => [Char] -> a
error [Char]
"ToySolver.Data.Polynomial.Factorization.Hensel.hensel: k <= 0"
  | Bool
otherwise = forall a. (?callStack::CallStack) => Bool -> a -> a
assert Bool
precondition forall a b. (a -> b) -> a -> b
$ Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
go Integer
1 (forall a b. (a -> b) -> [a] -> [b]
map (forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall (p :: Nat). PrimeField p -> Integer
Data.FiniteField.toInteger) [UPolynomial (PrimeField p)]
fs1)
  where
    precondition :: Bool
precondition =
      forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall a. Num a => Integer -> a
fromInteger UPolynomial Integer
f forall a. Eq a => a -> a -> Bool
== forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial (PrimeField p)]
fs1 Bool -> Bool -> Bool
&&
      forall t. Degree t => t -> Integer
P.deg UPolynomial Integer
f forall a. Eq a => a -> a -> Bool
== forall t. Degree t => t -> Integer
P.deg (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial (PrimeField p)]
fs1)

    p :: Integer
    p :: Integer
p = forall (n :: Nat) (proxy :: Nat -> *).
KnownNat n =>
proxy n -> Integer
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy p)

    go :: Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
    go :: Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
go !Integer
i [UPolynomial Integer]
fs
      | Integer
iforall a. Eq a => a -> a -> Bool
==Integer
k      = forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer -> [UPolynomial Integer] -> Bool
check Integer
i [UPolynomial Integer]
fs) forall a b. (a -> b) -> a -> b
$ [UPolynomial Integer]
fs
      | Bool
otherwise = forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer -> [UPolynomial Integer] -> Bool
check Integer
i [UPolynomial Integer]
fs) forall a b. (a -> b) -> a -> b
$ Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
go (Integer
iforall a. Num a => a -> a -> a
+Integer
1) (Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
lift Integer
i [UPolynomial Integer]
fs)

    check :: Integer -> [UPolynomial Integer] -> Bool
    check :: Integer -> [UPolynomial Integer] -> Bool
check Integer
k' [UPolynomial Integer]
fs =
        forall (t :: * -> *). Foldable t => t Bool -> Bool
and
        [ forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff (forall a. Integral a => a -> a -> a
`mod` Integer
pk) UPolynomial Integer
f forall a. Eq a => a -> a -> Bool
== forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff (forall a. Integral a => a -> a -> a
`mod` Integer
pk) (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial Integer]
fs)
        , [UPolynomial (PrimeField p)]
fs1 forall a. Eq a => a -> a -> Bool
== forall a b. (a -> b) -> [a] -> [b]
map (forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall a. Num a => Integer -> a
fromInteger) [UPolynomial Integer]
fs
        , forall (t :: * -> *). Foldable t => t Bool -> Bool
and [forall t. Degree t => t -> Integer
P.deg UPolynomial (PrimeField p)
fi1 forall a. Eq a => a -> a -> Bool
== forall t. Degree t => t -> Integer
P.deg UPolynomial Integer
fik | (UPolynomial (PrimeField p)
fi1, UPolynomial Integer
fik) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial (PrimeField p)]
fs1 [UPolynomial Integer]
fs]
        ]
      where
        pk :: Integer
pk = Integer
p forall a b. (Num a, Integral b) => a -> b -> a
^ Integer
k'

    lift :: Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
    lift :: Integer -> [UPolynomial Integer] -> [UPolynomial Integer]
lift Integer
k' [UPolynomial Integer]
fs = [UPolynomial Integer]
fs'
      where
        pk :: Integer
pk  = Integer
pforall a b. (Num a, Integral b) => a -> b -> a
^Integer
k'
        pk1 :: Integer
pk1 = Integer
pforall a b. (Num a, Integral b) => a -> b -> a
^(Integer
k'forall a. Num a => a -> a -> a
+Integer
1)

        -- f ≡ product fs + p^k h  (mod p^(k+1))
        h :: UPolynomial Integer
        h :: UPolynomial Integer
h = forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff (\Integer
c -> (Integer
c forall a. Integral a => a -> a -> a
`mod` Integer
pk1) forall a. Integral a => a -> a -> a
`div` Integer
pk) (UPolynomial Integer
f forall a. Num a => a -> a -> a
- forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial Integer]
fs)

        hs :: [UPolynomial (PrimeField p)]
        hs :: [UPolynomial (PrimeField p)]
hs = forall k.
(Fractional k, Ord k) =>
[UPolynomial k] -> UPolynomial k -> [UPolynomial k]
cabook_proposition_5_11 (forall a b. (a -> b) -> [a] -> [b]
map (forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall a. Num a => Integer -> a
fromInteger) [UPolynomial Integer]
fs) (forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall a. Num a => Integer -> a
fromInteger UPolynomial Integer
h)

        fs' :: [UPolynomial Integer]
        fs' :: [UPolynomial Integer]
fs' = [ forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff (forall a. Integral a => a -> a -> a
`mod` Integer
pk1) (UPolynomial Integer
fi forall a. Num a => a -> a -> a
+ forall k v. (Eq k, Num k) => k -> Polynomial k v
P.constant Integer
pk forall a. Num a => a -> a -> a
* forall k1 k v.
(Eq k1, Num k1) =>
(k -> k1) -> Polynomial k v -> Polynomial k1 v
P.mapCoeff forall (p :: Nat). PrimeField p -> Integer
Data.FiniteField.toInteger UPolynomial (PrimeField p)
hi)
              | (UPolynomial Integer
fi, UPolynomial (PrimeField p)
hi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial Integer]
fs [UPolynomial (PrimeField p)]
hs ]

-- http://www.math.kobe-u.ac.jp/Asir/ca.pdf
cabook_proposition_5_10
  :: forall k. (Num k, Fractional k, Eq k)
  => [UPolynomial k] -> [UPolynomial k]
cabook_proposition_5_10 :: forall k.
(Num k, Fractional k, Eq k) =>
[UPolynomial k] -> [UPolynomial k]
cabook_proposition_5_10 [UPolynomial k]
fs = [UPolynomial k] -> [UPolynomial k]
normalize ([UPolynomial k] -> [UPolynomial k]
go [UPolynomial k]
fs)
  where
    check :: [UPolynomial k] -> [UPolynomial k] -> Bool
    check :: [UPolynomial k] -> [UPolynomial k] -> Bool
check [UPolynomial k]
es [UPolynomial k]
fs' = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [UPolynomial k
ei forall a. Num a => a -> a -> a
* (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial k]
fs' forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` UPolynomial k
fi) | (UPolynomial k
ei,UPolynomial k
fi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial k]
es [UPolynomial k]
fs'] forall a. Eq a => a -> a -> Bool
== UPolynomial k
1

    go :: [UPolynomial k] -> [UPolynomial k]
    go :: [UPolynomial k] -> [UPolynomial k]
go [] = forall a. (?callStack::CallStack) => [Char] -> a
error [Char]
"cabook_proposition_5_10: empty list"
    go [UPolynomial k
fi] = forall a. (?callStack::CallStack) => Bool -> a -> a
assert ([UPolynomial k] -> [UPolynomial k] -> Bool
check [UPolynomial k
1] [UPolynomial k
fi]) [UPolynomial k
1]
    go (UPolynomial k
fi : [UPolynomial k]
fs') =
      case forall k.
(Eq k, Fractional k) =>
UPolynomial k
-> UPolynomial k -> (UPolynomial k, UPolynomial k, UPolynomial k)
P.exgcd (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial k]
fs') UPolynomial k
fi of
        (UPolynomial k
g,UPolynomial k
ei,UPolynomial k
v) ->
           forall a. (?callStack::CallStack) => Bool -> a -> a
assert (UPolynomial k
g forall a. Eq a => a -> a -> Bool
== UPolynomial k
1) forall a b. (a -> b) -> a -> b
$
             let es' :: [UPolynomial k]
es' = [UPolynomial k] -> [UPolynomial k]
go [UPolynomial k]
fs'
                 es :: [UPolynomial k]
es  = UPolynomial k
ei forall a. a -> [a] -> [a]
: forall a b. (a -> b) -> [a] -> [b]
map (UPolynomial k
vforall a. Num a => a -> a -> a
*) [UPolynomial k]
es'
             in forall a. (?callStack::CallStack) => Bool -> a -> a
assert ([UPolynomial k] -> [UPolynomial k] -> Bool
check [UPolynomial k]
es (UPolynomial k
fi forall a. a -> [a] -> [a]
: [UPolynomial k]
fs')) [UPolynomial k]
es

    normalize :: [UPolynomial k] -> [UPolynomial k]
    normalize :: [UPolynomial k] -> [UPolynomial k]
normalize [UPolynomial k]
es = forall a. (?callStack::CallStack) => Bool -> a -> a
assert ([UPolynomial k] -> [UPolynomial k] -> Bool
check [UPolynomial k]
es2 [UPolynomial k]
fs) [UPolynomial k]
es2
      where
        es2 :: [UPolynomial k]
es2 = forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.mod [UPolynomial k]
es [UPolynomial k]
fs

-- http://www.math.kobe-u.ac.jp/Asir/ca.pdf
cabook_proposition_5_11
  :: forall k. (Fractional k, Ord k)
  => [UPolynomial k] -> UPolynomial k -> [UPolynomial k]
cabook_proposition_5_11 :: forall k.
(Fractional k, Ord k) =>
[UPolynomial k] -> UPolynomial k -> [UPolynomial k]
cabook_proposition_5_11 [UPolynomial k]
fs UPolynomial k
g =
  forall a. (?callStack::CallStack) => Bool -> a -> a
assert (forall t. Degree t => t -> Integer
P.deg UPolynomial k
g forall a. Ord a => a -> a -> Bool
<= forall t. Degree t => t -> Integer
P.deg (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial k]
fs)) forall a b. (a -> b) -> a -> b
$
  forall a. (?callStack::CallStack) => Bool -> a -> a
assert (forall t. Degree t => t -> Integer
P.deg UPolynomial k
c forall a. Ord a => a -> a -> Bool
<= Integer
0) forall a b. (a -> b) -> a -> b
$
  forall a. (?callStack::CallStack) => Bool -> a -> a
assert ([UPolynomial k] -> Bool
check [UPolynomial k]
es2) forall a b. (a -> b) -> a -> b
$
    [UPolynomial k]
es2
  where
    es :: [UPolynomial k]
es  = forall a b. (a -> b) -> [a] -> [b]
map (UPolynomial k
gforall a. Num a => a -> a -> a
*) forall a b. (a -> b) -> a -> b
$ forall k.
(Num k, Fractional k, Eq k) =>
[UPolynomial k] -> [UPolynomial k]
cabook_proposition_5_10 [UPolynomial k]
fs
    c :: UPolynomial k
c   = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [UPolynomial k
ei forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` UPolynomial k
fi | (UPolynomial k
ei,UPolynomial k
fi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial k]
es [UPolynomial k]
fs]
    es2 :: [UPolynomial k]
es2 = case forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
P.mod [UPolynomial k]
es [UPolynomial k]
fs of
            UPolynomial k
e2' : [UPolynomial k]
es2' -> UPolynomial k
e2' forall a. Num a => a -> a -> a
+ UPolynomial k
c forall a. Num a => a -> a -> a
* forall a. [a] -> a
head [UPolynomial k]
fs forall a. a -> [a] -> [a]
: [UPolynomial k]
es2'

    check :: [UPolynomial k] -> Bool
    check :: [UPolynomial k] -> Bool
check [UPolynomial k]
es' =
      forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [UPolynomial k
ei forall a. Num a => a -> a -> a
* (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [UPolynomial k]
fs forall k.
(Eq k, Fractional k) =>
UPolynomial k -> UPolynomial k -> UPolynomial k
`P.div` UPolynomial k
fi) | (UPolynomial k
ei,UPolynomial k
fi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial k]
es' [UPolynomial k]
fs] forall a. Eq a => a -> a -> Bool
== UPolynomial k
g Bool -> Bool -> Bool
&&
      forall (t :: * -> *). Foldable t => t Bool -> Bool
and [forall t. Degree t => t -> Integer
P.deg UPolynomial k
ei forall a. Ord a => a -> a -> Bool
<= forall t. Degree t => t -> Integer
P.deg UPolynomial k
fi | (UPolynomial k
ei,UPolynomial k
fi) <- forall a b. [a] -> [b] -> [(a, b)]
zip [UPolynomial k]
es' [UPolynomial k]
fs]