{- |
module: Factor.Nfs
description: The general number field sieve
license: MIT

maintainer: Joe Leslie-Hurd <joe@gilith.com>
stability: provisional
portability: portable
-}
module Factor.Nfs
where

import Control.Monad (guard)
import qualified Data.IntSet as IntSet
import qualified Data.List as List
import Data.Maybe (isNothing,mapMaybe)

import qualified Factor.Bz as Bz
import qualified Factor.Gfpx as Gfpx
import Factor.Nfzw (Ideal,Nfzw(..))
import qualified Factor.Nfzw as Nfzw
import Factor.Prime (Prime)
import qualified Factor.Prime as Prime
import Factor.Util
import Factor.Zx (Zx)
import qualified Factor.Zx as Zx

-------------------------------------------------------------------------------
-- Polynomial selection
--
-- Given n, return m and f such that f is monic, irreducible, and f(m) == n
-------------------------------------------------------------------------------

data PolynomialDegree =
    FixedPolynomialDegree Int
  | OptimalPolynomialDegree
  deriving (PolynomialDegree -> PolynomialDegree -> Bool
(PolynomialDegree -> PolynomialDegree -> Bool)
-> (PolynomialDegree -> PolynomialDegree -> Bool)
-> Eq PolynomialDegree
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PolynomialDegree -> PolynomialDegree -> Bool
$c/= :: PolynomialDegree -> PolynomialDegree -> Bool
== :: PolynomialDegree -> PolynomialDegree -> Bool
$c== :: PolynomialDegree -> PolynomialDegree -> Bool
Eq,Eq PolynomialDegree
Eq PolynomialDegree
-> (PolynomialDegree -> PolynomialDegree -> Ordering)
-> (PolynomialDegree -> PolynomialDegree -> Bool)
-> (PolynomialDegree -> PolynomialDegree -> Bool)
-> (PolynomialDegree -> PolynomialDegree -> Bool)
-> (PolynomialDegree -> PolynomialDegree -> Bool)
-> (PolynomialDegree -> PolynomialDegree -> PolynomialDegree)
-> (PolynomialDegree -> PolynomialDegree -> PolynomialDegree)
-> Ord PolynomialDegree
PolynomialDegree -> PolynomialDegree -> Bool
PolynomialDegree -> PolynomialDegree -> Ordering
PolynomialDegree -> PolynomialDegree -> PolynomialDegree
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: PolynomialDegree -> PolynomialDegree -> PolynomialDegree
$cmin :: PolynomialDegree -> PolynomialDegree -> PolynomialDegree
max :: PolynomialDegree -> PolynomialDegree -> PolynomialDegree
$cmax :: PolynomialDegree -> PolynomialDegree -> PolynomialDegree
>= :: PolynomialDegree -> PolynomialDegree -> Bool
$c>= :: PolynomialDegree -> PolynomialDegree -> Bool
> :: PolynomialDegree -> PolynomialDegree -> Bool
$c> :: PolynomialDegree -> PolynomialDegree -> Bool
<= :: PolynomialDegree -> PolynomialDegree -> Bool
$c<= :: PolynomialDegree -> PolynomialDegree -> Bool
< :: PolynomialDegree -> PolynomialDegree -> Bool
$c< :: PolynomialDegree -> PolynomialDegree -> Bool
compare :: PolynomialDegree -> PolynomialDegree -> Ordering
$ccompare :: PolynomialDegree -> PolynomialDegree -> Ordering
$cp1Ord :: Eq PolynomialDegree
Ord,Int -> PolynomialDegree -> ShowS
[PolynomialDegree] -> ShowS
PolynomialDegree -> String
(Int -> PolynomialDegree -> ShowS)
-> (PolynomialDegree -> String)
-> ([PolynomialDegree] -> ShowS)
-> Show PolynomialDegree
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [PolynomialDegree] -> ShowS
$cshowList :: [PolynomialDegree] -> ShowS
show :: PolynomialDegree -> String
$cshow :: PolynomialDegree -> String
showsPrec :: Int -> PolynomialDegree -> ShowS
$cshowsPrec :: Int -> PolynomialDegree -> ShowS
Show)

data PolynomialBase =
    FixedPolynomialBase Integer
  | ClosestPolynomialBase
  | FloorPolynomialBase
  deriving (PolynomialBase -> PolynomialBase -> Bool
(PolynomialBase -> PolynomialBase -> Bool)
-> (PolynomialBase -> PolynomialBase -> Bool) -> Eq PolynomialBase
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PolynomialBase -> PolynomialBase -> Bool
$c/= :: PolynomialBase -> PolynomialBase -> Bool
== :: PolynomialBase -> PolynomialBase -> Bool
$c== :: PolynomialBase -> PolynomialBase -> Bool
Eq,Eq PolynomialBase
Eq PolynomialBase
-> (PolynomialBase -> PolynomialBase -> Ordering)
-> (PolynomialBase -> PolynomialBase -> Bool)
-> (PolynomialBase -> PolynomialBase -> Bool)
-> (PolynomialBase -> PolynomialBase -> Bool)
-> (PolynomialBase -> PolynomialBase -> Bool)
-> (PolynomialBase -> PolynomialBase -> PolynomialBase)
-> (PolynomialBase -> PolynomialBase -> PolynomialBase)
-> Ord PolynomialBase
PolynomialBase -> PolynomialBase -> Bool
PolynomialBase -> PolynomialBase -> Ordering
PolynomialBase -> PolynomialBase -> PolynomialBase
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: PolynomialBase -> PolynomialBase -> PolynomialBase
$cmin :: PolynomialBase -> PolynomialBase -> PolynomialBase
max :: PolynomialBase -> PolynomialBase -> PolynomialBase
$cmax :: PolynomialBase -> PolynomialBase -> PolynomialBase
>= :: PolynomialBase -> PolynomialBase -> Bool
$c>= :: PolynomialBase -> PolynomialBase -> Bool
> :: PolynomialBase -> PolynomialBase -> Bool
$c> :: PolynomialBase -> PolynomialBase -> Bool
<= :: PolynomialBase -> PolynomialBase -> Bool
$c<= :: PolynomialBase -> PolynomialBase -> Bool
< :: PolynomialBase -> PolynomialBase -> Bool
$c< :: PolynomialBase -> PolynomialBase -> Bool
compare :: PolynomialBase -> PolynomialBase -> Ordering
$ccompare :: PolynomialBase -> PolynomialBase -> Ordering
$cp1Ord :: Eq PolynomialBase
Ord,Int -> PolynomialBase -> ShowS
[PolynomialBase] -> ShowS
PolynomialBase -> String
(Int -> PolynomialBase -> ShowS)
-> (PolynomialBase -> String)
-> ([PolynomialBase] -> ShowS)
-> Show PolynomialBase
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [PolynomialBase] -> ShowS
$cshowList :: [PolynomialBase] -> ShowS
show :: PolynomialBase -> String
$cshow :: PolynomialBase -> String
showsPrec :: Int -> PolynomialBase -> ShowS
$cshowsPrec :: Int -> PolynomialBase -> ShowS
Show)

data PolynomialCoeff =
    FixedPolynomialCoeff [Integer]
  | SmallestPolynomialCoeff
  | PositivePolynomialCoeff
  deriving (PolynomialCoeff -> PolynomialCoeff -> Bool
(PolynomialCoeff -> PolynomialCoeff -> Bool)
-> (PolynomialCoeff -> PolynomialCoeff -> Bool)
-> Eq PolynomialCoeff
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c/= :: PolynomialCoeff -> PolynomialCoeff -> Bool
== :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c== :: PolynomialCoeff -> PolynomialCoeff -> Bool
Eq,Eq PolynomialCoeff
Eq PolynomialCoeff
-> (PolynomialCoeff -> PolynomialCoeff -> Ordering)
-> (PolynomialCoeff -> PolynomialCoeff -> Bool)
-> (PolynomialCoeff -> PolynomialCoeff -> Bool)
-> (PolynomialCoeff -> PolynomialCoeff -> Bool)
-> (PolynomialCoeff -> PolynomialCoeff -> Bool)
-> (PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff)
-> (PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff)
-> Ord PolynomialCoeff
PolynomialCoeff -> PolynomialCoeff -> Bool
PolynomialCoeff -> PolynomialCoeff -> Ordering
PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff
$cmin :: PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff
max :: PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff
$cmax :: PolynomialCoeff -> PolynomialCoeff -> PolynomialCoeff
>= :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c>= :: PolynomialCoeff -> PolynomialCoeff -> Bool
> :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c> :: PolynomialCoeff -> PolynomialCoeff -> Bool
<= :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c<= :: PolynomialCoeff -> PolynomialCoeff -> Bool
< :: PolynomialCoeff -> PolynomialCoeff -> Bool
$c< :: PolynomialCoeff -> PolynomialCoeff -> Bool
compare :: PolynomialCoeff -> PolynomialCoeff -> Ordering
$ccompare :: PolynomialCoeff -> PolynomialCoeff -> Ordering
$cp1Ord :: Eq PolynomialCoeff
Ord,Int -> PolynomialCoeff -> ShowS
[PolynomialCoeff] -> ShowS
PolynomialCoeff -> String
(Int -> PolynomialCoeff -> ShowS)
-> (PolynomialCoeff -> String)
-> ([PolynomialCoeff] -> ShowS)
-> Show PolynomialCoeff
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [PolynomialCoeff] -> ShowS
$cshowList :: [PolynomialCoeff] -> ShowS
show :: PolynomialCoeff -> String
$cshow :: PolynomialCoeff -> String
showsPrec :: Int -> PolynomialCoeff -> ShowS
$cshowsPrec :: Int -> PolynomialCoeff -> ShowS
Show)

data PolynomialConfig =
    PolynomialConfig
      {PolynomialConfig -> PolynomialDegree
polynomialDegree :: PolynomialDegree,
       PolynomialConfig -> PolynomialBase
polynomialBase :: PolynomialBase,
       PolynomialConfig -> PolynomialCoeff
polynomialCoeff :: PolynomialCoeff}
  deriving (PolynomialConfig -> PolynomialConfig -> Bool
(PolynomialConfig -> PolynomialConfig -> Bool)
-> (PolynomialConfig -> PolynomialConfig -> Bool)
-> Eq PolynomialConfig
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PolynomialConfig -> PolynomialConfig -> Bool
$c/= :: PolynomialConfig -> PolynomialConfig -> Bool
== :: PolynomialConfig -> PolynomialConfig -> Bool
$c== :: PolynomialConfig -> PolynomialConfig -> Bool
Eq,Eq PolynomialConfig
Eq PolynomialConfig
-> (PolynomialConfig -> PolynomialConfig -> Ordering)
-> (PolynomialConfig -> PolynomialConfig -> Bool)
-> (PolynomialConfig -> PolynomialConfig -> Bool)
-> (PolynomialConfig -> PolynomialConfig -> Bool)
-> (PolynomialConfig -> PolynomialConfig -> Bool)
-> (PolynomialConfig -> PolynomialConfig -> PolynomialConfig)
-> (PolynomialConfig -> PolynomialConfig -> PolynomialConfig)
-> Ord PolynomialConfig
PolynomialConfig -> PolynomialConfig -> Bool
PolynomialConfig -> PolynomialConfig -> Ordering
PolynomialConfig -> PolynomialConfig -> PolynomialConfig
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: PolynomialConfig -> PolynomialConfig -> PolynomialConfig
$cmin :: PolynomialConfig -> PolynomialConfig -> PolynomialConfig
max :: PolynomialConfig -> PolynomialConfig -> PolynomialConfig
$cmax :: PolynomialConfig -> PolynomialConfig -> PolynomialConfig
>= :: PolynomialConfig -> PolynomialConfig -> Bool
$c>= :: PolynomialConfig -> PolynomialConfig -> Bool
> :: PolynomialConfig -> PolynomialConfig -> Bool
$c> :: PolynomialConfig -> PolynomialConfig -> Bool
<= :: PolynomialConfig -> PolynomialConfig -> Bool
$c<= :: PolynomialConfig -> PolynomialConfig -> Bool
< :: PolynomialConfig -> PolynomialConfig -> Bool
$c< :: PolynomialConfig -> PolynomialConfig -> Bool
compare :: PolynomialConfig -> PolynomialConfig -> Ordering
$ccompare :: PolynomialConfig -> PolynomialConfig -> Ordering
$cp1Ord :: Eq PolynomialConfig
Ord,Int -> PolynomialConfig -> ShowS
[PolynomialConfig] -> ShowS
PolynomialConfig -> String
(Int -> PolynomialConfig -> ShowS)
-> (PolynomialConfig -> String)
-> ([PolynomialConfig] -> ShowS)
-> Show PolynomialConfig
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [PolynomialConfig] -> ShowS
$cshowList :: [PolynomialConfig] -> ShowS
show :: PolynomialConfig -> String
$cshow :: PolynomialConfig -> String
showsPrec :: Int -> PolynomialConfig -> ShowS
$cshowsPrec :: Int -> PolynomialConfig -> ShowS
Show)

defaultPolynomialConfig :: PolynomialConfig
defaultPolynomialConfig :: PolynomialConfig
defaultPolynomialConfig =
    PolynomialConfig :: PolynomialDegree
-> PolynomialBase -> PolynomialCoeff -> PolynomialConfig
PolynomialConfig
      {polynomialDegree :: PolynomialDegree
polynomialDegree = PolynomialDegree
OptimalPolynomialDegree,
       polynomialBase :: PolynomialBase
polynomialBase = PolynomialBase
ClosestPolynomialBase,
       polynomialCoeff :: PolynomialCoeff
polynomialCoeff = PolynomialCoeff
SmallestPolynomialCoeff}

fixedPolynomialConfig :: Zx -> Integer -> PolynomialConfig
fixedPolynomialConfig :: Zx -> Integer -> PolynomialConfig
fixedPolynomialConfig Zx
f Integer
m =
    PolynomialConfig :: PolynomialDegree
-> PolynomialBase -> PolynomialCoeff -> PolynomialConfig
PolynomialConfig
      {polynomialDegree :: PolynomialDegree
polynomialDegree = Int -> PolynomialDegree
FixedPolynomialDegree (Zx -> Int
Zx.degree Zx
f),
       polynomialBase :: PolynomialBase
polynomialBase = Integer -> PolynomialBase
FixedPolynomialBase Integer
m,
       polynomialCoeff :: PolynomialCoeff
polynomialCoeff = [Integer] -> PolynomialCoeff
FixedPolynomialCoeff (Zx -> [Integer]
Zx.toCoeff Zx
f)}

selectPolynomialDegree :: PolynomialDegree -> Integer -> Int
selectPolynomialDegree :: PolynomialDegree -> Integer -> Int
selectPolynomialDegree (FixedPolynomialDegree Int
d) Integer
_ = Int
d
selectPolynomialDegree PolynomialDegree
OptimalPolynomialDegree Integer
n =
    Double -> Int
forall a b. (RealFrac a, Integral b) => a -> b
round (Double -> Int) -> Double -> Int
forall a b. (a -> b) -> a -> b
$ (((Double
3.0 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
l) Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double -> Double
forall a. Floating a => a -> a
log Double
l) Double -> Double -> Double
forall a. Floating a => a -> a -> a
** (Double
1.0Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
3.0))
  where
    l :: Double
l = Integer -> Double
logInteger Integer
n

selectPolynomialBase :: PolynomialBase -> Integer -> Int -> Integer
selectPolynomialBase :: PolynomialBase -> Integer -> Int -> Integer
selectPolynomialBase (FixedPolynomialBase Integer
m) Integer
_ Int
_ = Integer
m
selectPolynomialBase PolynomialBase
ClosestPolynomialBase Integer
n Int
d = Integer -> Integer -> Integer
nthRootClosest (Int -> Integer
forall a. Integral a => a -> Integer
toInteger Int
d) Integer
n
selectPolynomialBase PolynomialBase
FloorPolynomialBase Integer
n Int
d = Integer -> Integer -> Integer
nthRoot (Int -> Integer
forall a. Integral a => a -> Integer
toInteger Int
d) Integer
n

selectPolynomialCoeff :: PolynomialCoeff -> Integer -> Int -> Integer -> [Integer]
selectPolynomialCoeff :: PolynomialCoeff -> Integer -> Int -> Integer -> [Integer]
selectPolynomialCoeff (FixedPolynomialCoeff [Integer]
c) Integer
_ Int
_ Integer
_ = [Integer]
c
selectPolynomialCoeff PolynomialCoeff
cfg Integer
n Int
d Integer
m = Integer
x Integer -> [Integer] -> [Integer]
forall a. a -> [a] -> [a]
: [Integer]
c
  where
    (Integer
x,[Integer]
c) = (Integer -> (Integer, [Integer]) -> (Integer, [Integer]))
-> (Integer, [Integer]) -> [Integer] -> (Integer, [Integer])
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr Integer -> (Integer, [Integer]) -> (Integer, [Integer])
reduce (Integer
n Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- [Integer] -> Integer
forall a. [a] -> a
last [Integer]
ms, [Integer
1]) ([Integer] -> [Integer]
forall a. [a] -> [a]
init [Integer]
ms)

    reduce :: Integer -> (Integer, [Integer]) -> (Integer, [Integer])
reduce Integer
mi (Integer
t,[Integer]
l) = (Integer
r, Integer
q Integer -> [Integer] -> [Integer]
forall a. a -> [a] -> [a]
: [Integer]
l)
      where (Integer
q,Integer
r) = Integer -> Integer -> (Integer, Integer)
idiv Integer
t Integer
mi

    ms :: [Integer]
ms = Int -> [Integer] -> [Integer]
forall a. Int -> [a] -> [a]
take Int
d ([Integer] -> [Integer]) -> [Integer] -> [Integer]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> Integer -> [Integer]
forall a. (a -> a) -> a -> [a]
iterate (Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
(*) Integer
m) Integer
m

    idiv :: Integer -> Integer -> (Integer, Integer)
idiv = if PolynomialCoeff
cfg PolynomialCoeff -> PolynomialCoeff -> Bool
forall a. Eq a => a -> a -> Bool
== PolynomialCoeff
PositivePolynomialCoeff then Integer -> Integer -> (Integer, Integer)
division else Integer -> Integer -> (Integer, Integer)
divisionClosest

selectPolynomial :: PolynomialConfig -> Integer -> (Zx,Integer)
selectPolynomial :: PolynomialConfig -> Integer -> (Zx, Integer)
selectPolynomial PolynomialConfig
cfg Integer
n = ([Integer] -> Zx
Zx.fromCoeff [Integer]
c, Integer
m)
  where
    d :: Int
d = PolynomialDegree -> Integer -> Int
selectPolynomialDegree (PolynomialConfig -> PolynomialDegree
polynomialDegree PolynomialConfig
cfg) Integer
n
    m :: Integer
m = PolynomialBase -> Integer -> Int -> Integer
selectPolynomialBase (PolynomialConfig -> PolynomialBase
polynomialBase PolynomialConfig
cfg) Integer
n Int
d
    c :: [Integer]
c = PolynomialCoeff -> Integer -> Int -> Integer -> [Integer]
selectPolynomialCoeff (PolynomialConfig -> PolynomialCoeff
polynomialCoeff PolynomialConfig
cfg) Integer
n Int
d Integer
m

-------------------------------------------------------------------------------
-- Factor bases
-------------------------------------------------------------------------------

type FactorBase = [Prime]

data FactorBaseConfig =
    FixedFactorBase Integer
  | OptimalFactorBase Double
  deriving (FactorBaseConfig -> FactorBaseConfig -> Bool
(FactorBaseConfig -> FactorBaseConfig -> Bool)
-> (FactorBaseConfig -> FactorBaseConfig -> Bool)
-> Eq FactorBaseConfig
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c/= :: FactorBaseConfig -> FactorBaseConfig -> Bool
== :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c== :: FactorBaseConfig -> FactorBaseConfig -> Bool
Eq,Eq FactorBaseConfig
Eq FactorBaseConfig
-> (FactorBaseConfig -> FactorBaseConfig -> Ordering)
-> (FactorBaseConfig -> FactorBaseConfig -> Bool)
-> (FactorBaseConfig -> FactorBaseConfig -> Bool)
-> (FactorBaseConfig -> FactorBaseConfig -> Bool)
-> (FactorBaseConfig -> FactorBaseConfig -> Bool)
-> (FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig)
-> (FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig)
-> Ord FactorBaseConfig
FactorBaseConfig -> FactorBaseConfig -> Bool
FactorBaseConfig -> FactorBaseConfig -> Ordering
FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig
$cmin :: FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig
max :: FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig
$cmax :: FactorBaseConfig -> FactorBaseConfig -> FactorBaseConfig
>= :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c>= :: FactorBaseConfig -> FactorBaseConfig -> Bool
> :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c> :: FactorBaseConfig -> FactorBaseConfig -> Bool
<= :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c<= :: FactorBaseConfig -> FactorBaseConfig -> Bool
< :: FactorBaseConfig -> FactorBaseConfig -> Bool
$c< :: FactorBaseConfig -> FactorBaseConfig -> Bool
compare :: FactorBaseConfig -> FactorBaseConfig -> Ordering
$ccompare :: FactorBaseConfig -> FactorBaseConfig -> Ordering
$cp1Ord :: Eq FactorBaseConfig
Ord,Int -> FactorBaseConfig -> ShowS
[FactorBaseConfig] -> ShowS
FactorBaseConfig -> String
(Int -> FactorBaseConfig -> ShowS)
-> (FactorBaseConfig -> String)
-> ([FactorBaseConfig] -> ShowS)
-> Show FactorBaseConfig
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [FactorBaseConfig] -> ShowS
$cshowList :: [FactorBaseConfig] -> ShowS
show :: FactorBaseConfig -> String
$cshow :: FactorBaseConfig -> String
showsPrec :: Int -> FactorBaseConfig -> ShowS
$cshowsPrec :: Int -> FactorBaseConfig -> ShowS
Show)

defaultFactorBaseConfig :: FactorBaseConfig
defaultFactorBaseConfig :: FactorBaseConfig
defaultFactorBaseConfig = Double -> FactorBaseConfig
OptimalFactorBase Double
3.0

maxFactorBase :: FactorBaseConfig -> Integer -> Integer
maxFactorBase :: FactorBaseConfig -> Integer -> Integer
maxFactorBase (FixedFactorBase Integer
b) Integer
_ = Integer
b
maxFactorBase (OptimalFactorBase Double
c) Integer
n =
    Double -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Double -> Integer) -> Double -> Integer
forall a b. (a -> b) -> a -> b
$ Double
c Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double -> Double
forall a. Floating a => a -> a
exp (Double
t2Double -> Double -> Double
forall a. Floating a => a -> a -> a
**Double
t2 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
lDouble -> Double -> Double
forall a. Floating a => a -> a -> a
**Double
t1 Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
llDouble -> Double -> Double
forall a. Floating a => a -> a -> a
**Double
t2)
  where
    t1 :: Double
t1 = Double
1.0 Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
3.0
    t2 :: Double
t2 = Double
2.0 Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
3.0
    l :: Double
l = Integer -> Double
logInteger Integer
n
    ll :: Double
ll = Double -> Double
forall a. Floating a => a -> a
log Double
l

-------------------------------------------------------------------------------
-- Finding smooth elements of Z[w]
-------------------------------------------------------------------------------

destSmoothInteger :: FactorBase -> Integer -> Maybe Integer
destSmoothInteger :: [Integer] -> Integer -> Maybe Integer
destSmoothInteger [Integer]
_ Integer
n | Integer -> Integer
forall a. Num a => a -> a
abs Integer
n Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
1 = Maybe Integer
forall a. Maybe a
Nothing
destSmoothInteger [] Integer
n = if Integer
n Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
/= Integer
0 Bool -> Bool -> Bool
&& Integer -> Bool
isSquare (Integer -> Integer
forall a. Num a => a -> a
abs Integer
n) then Maybe Integer
forall a. Maybe a
Nothing else Integer -> Maybe Integer
forall a. a -> Maybe a
Just Integer
n
destSmoothInteger (Integer
p : [Integer]
ps) Integer
n = [Integer] -> Integer -> Maybe Integer
destSmoothInteger [Integer]
ps ((Integer, Integer) -> Integer
forall a b. (a, b) -> b
snd (Integer -> Integer -> (Integer, Integer)
divPower Integer
p Integer
n))

isSmoothInteger :: FactorBase -> Integer -> Bool
isSmoothInteger :: [Integer] -> Integer -> Bool
isSmoothInteger [Integer]
fb Integer
n = Maybe Integer -> Bool
forall a. Maybe a -> Bool
isNothing ([Integer] -> Integer -> Maybe Integer
destSmoothInteger [Integer]
fb Integer
n)

notSmoothInteger :: FactorBase -> Integer -> Bool
notSmoothInteger :: [Integer] -> Integer -> Bool
notSmoothInteger [Integer]
fb = Bool -> Bool
not (Bool -> Bool) -> (Integer -> Bool) -> Integer -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Integer] -> Integer -> Bool
isSmoothInteger [Integer]
fb

rationalNorm :: Integer -> Nfzw -> Integer
rationalNorm :: Integer -> Nfzw -> Integer
rationalNorm = Integer -> Nfzw -> Integer
Nfzw.toInteger

algebraicNorm :: Zx -> Nfzw -> Integer
algebraicNorm :: Zx -> Nfzw -> Integer
algebraicNorm = Zx -> Nfzw -> Integer
Nfzw.norm

isSmoothNfzw :: Zx -> Integer -> FactorBase -> FactorBase -> Nfzw -> Bool
isSmoothNfzw :: Zx -> Integer -> [Integer] -> [Integer] -> Nfzw -> Bool
isSmoothNfzw Zx
f Integer
m [Integer]
rfb [Integer]
afb Nfzw
x =
    [Integer] -> Integer -> Bool
isSmoothInteger [Integer]
rfb (Integer -> Nfzw -> Integer
rationalNorm Integer
m Nfzw
x) Bool -> Bool -> Bool
&&
    [Integer] -> Integer -> Bool
isSmoothInteger [Integer]
afb (Zx -> Nfzw -> Integer
algebraicNorm Zx
f Nfzw
x)

smoothNfzw :: Zx -> Integer -> FactorBase -> FactorBase -> [Nfzw]
smoothNfzw :: Zx -> Integer -> [Integer] -> [Integer] -> [Nfzw]
smoothNfzw Zx
f Integer
m [Integer]
rfb [Integer]
afb = (Nfzw -> Bool) -> [Nfzw]
Nfzw.filterValid (Zx -> Integer -> [Integer] -> [Integer] -> Nfzw -> Bool
isSmoothNfzw Zx
f Integer
m [Integer]
rfb [Integer]
afb)

-------------------------------------------------------------------------------
-- Quadratic characters
-------------------------------------------------------------------------------

data QuadraticCharacterConfig =
    FixedQuadraticCharacterConfig Int
  | LinearQuadraticCharacterConfig Int Double
  deriving (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
(QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> Eq QuadraticCharacterConfig
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c/= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
== :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c== :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
Eq,Eq QuadraticCharacterConfig
Eq QuadraticCharacterConfig
-> (QuadraticCharacterConfig
    -> QuadraticCharacterConfig -> Ordering)
-> (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> (QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool)
-> (QuadraticCharacterConfig
    -> QuadraticCharacterConfig -> QuadraticCharacterConfig)
-> (QuadraticCharacterConfig
    -> QuadraticCharacterConfig -> QuadraticCharacterConfig)
-> Ord QuadraticCharacterConfig
QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
QuadraticCharacterConfig -> QuadraticCharacterConfig -> Ordering
QuadraticCharacterConfig
-> QuadraticCharacterConfig -> QuadraticCharacterConfig
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: QuadraticCharacterConfig
-> QuadraticCharacterConfig -> QuadraticCharacterConfig
$cmin :: QuadraticCharacterConfig
-> QuadraticCharacterConfig -> QuadraticCharacterConfig
max :: QuadraticCharacterConfig
-> QuadraticCharacterConfig -> QuadraticCharacterConfig
$cmax :: QuadraticCharacterConfig
-> QuadraticCharacterConfig -> QuadraticCharacterConfig
>= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c>= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
> :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c> :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
<= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c<= :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
< :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
$c< :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Bool
compare :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Ordering
$ccompare :: QuadraticCharacterConfig -> QuadraticCharacterConfig -> Ordering
$cp1Ord :: Eq QuadraticCharacterConfig
Ord,Int -> QuadraticCharacterConfig -> ShowS
[QuadraticCharacterConfig] -> ShowS
QuadraticCharacterConfig -> String
(Int -> QuadraticCharacterConfig -> ShowS)
-> (QuadraticCharacterConfig -> String)
-> ([QuadraticCharacterConfig] -> ShowS)
-> Show QuadraticCharacterConfig
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [QuadraticCharacterConfig] -> ShowS
$cshowList :: [QuadraticCharacterConfig] -> ShowS
show :: QuadraticCharacterConfig -> String
$cshow :: QuadraticCharacterConfig -> String
showsPrec :: Int -> QuadraticCharacterConfig -> ShowS
$cshowsPrec :: Int -> QuadraticCharacterConfig -> ShowS
Show)

defaultQuadraticCharacterConfig :: QuadraticCharacterConfig
defaultQuadraticCharacterConfig :: QuadraticCharacterConfig
defaultQuadraticCharacterConfig = Int -> Double -> QuadraticCharacterConfig
LinearQuadraticCharacterConfig Int
20 Double
0.5

quadraticCharacters :: QuadraticCharacterConfig -> Integer -> Int
quadraticCharacters :: QuadraticCharacterConfig -> Integer -> Int
quadraticCharacters (FixedQuadraticCharacterConfig Int
q) Integer
_ = Int
q
quadraticCharacters (LinearQuadraticCharacterConfig Int
m Double
r) Integer
n =
    Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
m (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Double -> Int
forall a b. (RealFrac a, Integral b) => a -> b
floor (Double
r Double -> Double -> Double
forall a. Num a => a -> a -> a
* Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Integer -> Int
widthInteger Integer
n))

isQuadraticCharacter :: Zx -> [Nfzw] -> Ideal -> Bool
isQuadraticCharacter :: Zx -> [Nfzw] -> (Integer, Integer) -> Bool
isQuadraticCharacter Zx
f' [Nfzw]
xs (Integer, Integer)
i =
    Integer -> Gfpx -> Integer -> Integer
Gfpx.evaluate Integer
q Gfpx
fq' Integer
s Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
/= Integer
0 Bool -> Bool -> Bool
&&
    (Nfzw -> Bool) -> [Nfzw] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (\Nfzw
x -> Bool -> Bool
not (Nfzw -> (Integer, Integer) -> Bool
Nfzw.inIdeal Nfzw
x (Integer, Integer)
i)) [Nfzw]
xs
  where
    fq' :: Gfpx
fq' = Integer -> Zx -> Gfpx
Gfpx.fromZx Integer
q Zx
f'
    (Integer
s,Integer
q) = (Integer, Integer)
i

nextQuadraticCharacter :: Zx -> [Nfzw] -> [Ideal] -> (Ideal,[Ideal])
nextQuadraticCharacter :: Zx
-> [Nfzw]
-> [(Integer, Integer)]
-> ((Integer, Integer), [(Integer, Integer)])
nextQuadraticCharacter Zx
f' [Nfzw]
xs ((Integer, Integer)
i : [(Integer, Integer)]
il) | Zx -> [Nfzw] -> (Integer, Integer) -> Bool
isQuadraticCharacter Zx
f' [Nfzw]
xs (Integer, Integer)
i = ((Integer, Integer)
i,[(Integer, Integer)]
il)
nextQuadraticCharacter Zx
f' [Nfzw]
xs [(Integer, Integer)]
il = Zx
-> [Nfzw]
-> [(Integer, Integer)]
-> ((Integer, Integer), [(Integer, Integer)])
nextQuadraticCharacter Zx
f' [Nfzw]
xs ([(Integer, Integer)] -> [(Integer, Integer)]
forall a. [a] -> [a]
tail [(Integer, Integer)]
il)

takeQuadraticCharacters :: Zx -> [Nfzw] -> Int -> [Ideal] -> ([Ideal],[Ideal])
takeQuadraticCharacters :: Zx
-> [Nfzw]
-> Int
-> [(Integer, Integer)]
-> ([(Integer, Integer)], [(Integer, Integer)])
takeQuadraticCharacters Zx
f' [Nfzw]
xs = ([(Integer, Integer)]
 -> ((Integer, Integer), [(Integer, Integer)]))
-> Int
-> [(Integer, Integer)]
-> ([(Integer, Integer)], [(Integer, Integer)])
forall b a. (b -> (a, b)) -> Int -> b -> ([a], b)
unfoldlN (([(Integer, Integer)]
  -> ((Integer, Integer), [(Integer, Integer)]))
 -> Int
 -> [(Integer, Integer)]
 -> ([(Integer, Integer)], [(Integer, Integer)]))
-> ([(Integer, Integer)]
    -> ((Integer, Integer), [(Integer, Integer)]))
-> Int
-> [(Integer, Integer)]
-> ([(Integer, Integer)], [(Integer, Integer)])
forall a b. (a -> b) -> a -> b
$ Zx
-> [Nfzw]
-> [(Integer, Integer)]
-> ((Integer, Integer), [(Integer, Integer)])
nextQuadraticCharacter Zx
f' [Nfzw]
xs

isQuadraticResidue :: Ideal -> Nfzw -> Bool
isQuadraticResidue :: (Integer, Integer) -> Nfzw -> Bool
isQuadraticResidue (Integer
s,Integer
q) Nfzw
x =
    case Integer -> Integer -> Residue
jacobiSymbol (Integer -> Nfzw -> Integer
Nfzw.toInteger Integer
s Nfzw
x) Integer
q of
      Residue
Residue -> Bool
True
      Residue
ZeroResidue -> String -> Bool
forall a. HasCallStack => String -> a
error String
"zero residue"
      Residue
NonResidue -> Bool
False

notQuadraticResidue :: Ideal -> Nfzw -> Bool
notQuadraticResidue :: (Integer, Integer) -> Nfzw -> Bool
notQuadraticResidue (Integer, Integer)
i = Bool -> Bool
not (Bool -> Bool) -> (Nfzw -> Bool) -> Nfzw -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer, Integer) -> Nfzw -> Bool
isQuadraticResidue (Integer, Integer)
i

productIsQuadraticResidue :: Ideal -> [Nfzw] -> Bool
productIsQuadraticResidue :: (Integer, Integer) -> [Nfzw] -> Bool
productIsQuadraticResidue (Integer, Integer)
i = (Nfzw -> Bool -> Bool) -> Bool -> [Nfzw] -> Bool
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr (\Nfzw
x Bool
b -> (Integer, Integer) -> Nfzw -> Bool
isQuadraticResidue (Integer, Integer)
i Nfzw
x Bool -> Bool -> Bool
forall a. Eq a => a -> a -> Bool
== Bool
b) Bool
True

-------------------------------------------------------------------------------
-- Creating the matrix
-------------------------------------------------------------------------------

type Row = [Bool]

formRow :: Zx -> Integer -> [Prime] -> [Ideal] -> [Ideal] -> Nfzw -> Row
formRow :: Zx
-> Integer
-> [Integer]
-> [(Integer, Integer)]
-> [(Integer, Integer)]
-> Nfzw
-> Row
formRow Zx
f Integer
m [Integer]
rps [(Integer, Integer)]
aps [(Integer, Integer)]
qcs Nfzw
x =
    [Integer] -> Integer -> Row
rationalRow [Integer]
rps (Integer -> Nfzw -> Integer
Nfzw.toInteger Integer
m Nfzw
x) Row -> Row -> Row
forall a. [a] -> [a] -> [a]
++
    Zx -> [(Integer, Integer)] -> [(Integer, Integer)] -> Nfzw -> Row
algebraicRow Zx
f [(Integer, Integer)]
aps [(Integer, Integer)]
qcs Nfzw
x

rationalRow :: [Prime] -> Integer -> Row
rationalRow :: [Integer] -> Integer -> Row
rationalRow [Integer]
rps Integer
n =
    (Integer
n Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
0) Bool -> Row -> Row
forall a. a -> [a] -> [a]
: ([(Integer, Integer)], Row) -> Row
forall a b. (a, b) -> b
snd (([(Integer, Integer)] -> Integer -> ([(Integer, Integer)], Bool))
-> [(Integer, Integer)] -> [Integer] -> ([(Integer, Integer)], Row)
forall (t :: * -> *) a b c.
Traversable t =>
(a -> b -> (a, c)) -> a -> t b -> (a, t c)
List.mapAccumL [(Integer, Integer)] -> Integer -> ([(Integer, Integer)], Bool)
forall a. Eq a => [(a, Integer)] -> a -> ([(a, Integer)], Bool)
oddPower [(Integer, Integer)]
pks [Integer]
rps)
  where
    ([(Integer, Integer)]
pks,Integer
_) = [Integer] -> Integer -> ([(Integer, Integer)], Integer)
Prime.factor [Integer]
rps (Integer -> Integer
forall a. Num a => a -> a
abs Integer
n)

algebraicRow :: Zx -> [Ideal] -> [Ideal] -> Nfzw -> Row
algebraicRow :: Zx -> [(Integer, Integer)] -> [(Integer, Integer)] -> Nfzw -> Row
algebraicRow Zx
f [(Integer, Integer)]
aps [(Integer, Integer)]
qcs Nfzw
x =
    ([((Integer, Integer), Integer)], Row) -> Row
forall a b. (a, b) -> b
snd (([((Integer, Integer), Integer)]
 -> (Integer, Integer) -> ([((Integer, Integer), Integer)], Bool))
-> [((Integer, Integer), Integer)]
-> [(Integer, Integer)]
-> ([((Integer, Integer), Integer)], Row)
forall (t :: * -> *) a b c.
Traversable t =>
(a -> b -> (a, c)) -> a -> t b -> (a, t c)
List.mapAccumL [((Integer, Integer), Integer)]
-> (Integer, Integer) -> ([((Integer, Integer), Integer)], Bool)
forall a. Eq a => [(a, Integer)] -> a -> ([(a, Integer)], Bool)
oddPower [((Integer, Integer), Integer)]
iks [(Integer, Integer)]
aps) Row -> Row -> Row
forall a. [a] -> [a] -> [a]
++
    ((Integer, Integer) -> Bool) -> [(Integer, Integer)] -> Row
forall a b. (a -> b) -> [a] -> [b]
map (((Integer, Integer) -> Nfzw -> Bool)
-> Nfzw -> (Integer, Integer) -> Bool
forall a b c. (a -> b -> c) -> b -> a -> c
flip (Integer, Integer) -> Nfzw -> Bool
notQuadraticResidue Nfzw
x) [(Integer, Integer)]
qcs
  where
    ([((Integer, Integer), Integer)]
iks,Integer
_) = Zx
-> [(Integer, Integer)]
-> Nfzw
-> ([((Integer, Integer), Integer)], Integer)
Nfzw.factor Zx
f [(Integer, Integer)]
aps Nfzw
x

oddPower :: Eq a => [(a,Integer)] -> a -> ([(a,Integer)],Bool)
oddPower :: [(a, Integer)] -> a -> ([(a, Integer)], Bool)
oddPower ((a
p,Integer
k) : [(a, Integer)]
pks) a
q | a
p a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
q = ([(a, Integer)]
pks, Integer -> Bool
forall a. Integral a => a -> Bool
odd Integer
k)
oddPower [(a, Integer)]
pks a
_ = ([(a, Integer)]
pks,Bool
False)

-------------------------------------------------------------------------------
-- Gaussian elimination
-------------------------------------------------------------------------------

gaussianElimination :: [Row] -> [[Int]]
gaussianElimination :: [Row] -> [[Int]]
gaussianElimination [Row]
rows =
    (IntSet -> [Int]) -> [IntSet] -> [[Int]]
forall a b. (a -> b) -> [a] -> [b]
map IntSet -> [Int]
IntSet.toList ([IntSet] -> [[Int]]) -> [IntSet] -> [[Int]]
forall a b. (a -> b) -> a -> b
$
    (IntSet -> Int) -> [IntSet] -> [IntSet]
forall b a. Ord b => (a -> b) -> [a] -> [a]
List.sortOn IntSet -> Int
IntSet.size ([IntSet] -> [IntSet]) -> [IntSet] -> [IntSet]
forall a b. (a -> b) -> a -> b
$
    ([IntSet] -> IntSet -> [IntSet])
-> [IntSet] -> [IntSet] -> [IntSet]
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
List.foldl' [IntSet] -> IntSet -> [IntSet]
processColumn [IntSet]
ortho [IntSet]
columns
  where
    ortho :: [IntSet]
ortho = (Int -> IntSet) -> [Int] -> [IntSet]
forall a b. (a -> b) -> [a] -> [b]
map Int -> IntSet
IntSet.singleton [Int
0..([Row] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Row]
rows Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)]

    columns :: [IntSet]
columns = (Row -> IntSet) -> [Row] -> [IntSet]
forall a b. (a -> b) -> [a] -> [b]
map Row -> IntSet
mkVec ([Row] -> [IntSet]) -> [Row] -> [IntSet]
forall a b. (a -> b) -> a -> b
$ [Row] -> [Row]
forall a. [[a]] -> [[a]]
List.transpose [Row]
rows

    processColumn :: [IntSet] -> IntSet -> [IntSet]
processColumn [IntSet]
basis IntSet
column =
        case (IntSet -> Bool) -> [IntSet] -> ([IntSet], [IntSet])
forall a. (a -> Bool) -> [a] -> ([a], [a])
List.partition (IntSet -> IntSet -> Bool
evalVec IntSet
column) [IntSet]
basis of
          ([],[IntSet]
_) -> [IntSet]
basis
          (IntSet
u : [IntSet]
us, [IntSet]
vs) -> (IntSet -> IntSet) -> [IntSet] -> [IntSet]
forall a b. (a -> b) -> [a] -> [b]
map (IntSet -> IntSet -> IntSet
addVec IntSet
u) [IntSet]
us [IntSet] -> [IntSet] -> [IntSet]
forall a. [a] -> [a] -> [a]
++ [IntSet]
vs

    mkVec :: Row -> IntSet
mkVec = [Int] -> IntSet
IntSet.fromDistinctAscList ([Int] -> IntSet) -> (Row -> [Int]) -> Row -> IntSet
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((Int, Bool) -> Int) -> [(Int, Bool)] -> [Int]
forall a b. (a -> b) -> [a] -> [b]
map (Int, Bool) -> Int
forall a b. (a, b) -> a
fst ([(Int, Bool)] -> [Int]) -> (Row -> [(Int, Bool)]) -> Row -> [Int]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((Int, Bool) -> Bool) -> [(Int, Bool)] -> [(Int, Bool)]
forall a. (a -> Bool) -> [a] -> [a]
filter (Int, Bool) -> Bool
forall a b. (a, b) -> b
snd ([(Int, Bool)] -> [(Int, Bool)])
-> (Row -> [(Int, Bool)]) -> Row -> [(Int, Bool)]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Int] -> Row -> [(Int, Bool)]
forall a b. [a] -> [b] -> [(a, b)]
zip [Int
0..]

    evalVec :: IntSet -> IntSet -> Bool
evalVec IntSet
u IntSet
v  = Int -> Bool
forall a. Integral a => a -> Bool
odd (IntSet -> Int
IntSet.size (IntSet -> IntSet -> IntSet
IntSet.intersection IntSet
u IntSet
v))

    addVec :: IntSet -> IntSet -> IntSet
addVec IntSet
u IntSet
v = IntSet -> IntSet -> IntSet
IntSet.difference (IntSet -> IntSet -> IntSet
IntSet.union IntSet
u IntSet
v) (IntSet -> IntSet -> IntSet
IntSet.intersection IntSet
u IntSet
v)

-------------------------------------------------------------------------------
-- Square roots
--
-- https://hal.inria.fr/hal-00756838/en
-------------------------------------------------------------------------------

rationalSquareRoot :: Integer -> Integer -> Zx -> [Prime] -> [Nfzw] -> Integer
rationalSquareRoot :: Integer -> Integer -> Zx -> [Integer] -> [Nfzw] -> Integer
rationalSquareRoot Integer
n Integer
m Zx
f' [Integer]
rps [Nfzw]
xs =
    Integer -> [Integer] -> Integer
Prime.product Integer
n (Integer
fm' Integer -> [Integer] -> [Integer]
forall a. a -> [a] -> [a]
: ((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> Integer
forall a. Integral a => (Integer, a) -> Integer
sqrtPk [(Integer, Integer)]
pks [Integer] -> [Integer] -> [Integer]
forall a. [a] -> [a] -> [a]
++ (Integer -> Integer) -> [Integer] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map Integer -> Integer
sqrtS [Integer]
sl)
  where
    fm' :: Integer
fm' = Integer -> Integer -> Integer
Prime.fromInteger Integer
n (Zx -> Integer -> Integer
Zx.evaluate Zx
f' Integer
m)

    ([(Integer, Integer)]
pks,[Integer]
sl) = [Integer] -> [Integer] -> ([(Integer, Integer)], [Integer])
Prime.factorProduct [Integer]
rps ((Nfzw -> Integer) -> [Nfzw] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> Nfzw -> Integer
Nfzw.toInteger Integer
m) [Nfzw]
xs)

    sqrtPk :: (Integer, a) -> Integer
sqrtPk (Integer
p,a
k) =
        if a -> Bool
forall a. Integral a => a -> Bool
odd a
k then String -> Integer
forall a. HasCallStack => String -> a
error String
"prime power is not a square"
        else Integer -> Integer -> Integer -> Integer
Prime.exp Integer
n (Integer -> Integer -> Integer
Prime.fromInteger Integer
n Integer
p) (a -> Integer
forall a. Integral a => a -> Integer
toInteger (a
k a -> a -> a
forall a. Integral a => a -> a -> a
`div` a
2))

    sqrtS :: Integer -> Integer
sqrtS Integer
s = case Integer -> Maybe Integer
destSquare (Integer -> Integer
forall a. Num a => a -> a
abs Integer
s) of
                Maybe Integer
Nothing -> String -> Integer
forall a. HasCallStack => String -> a
error String
"smooth remainder is not a square"
                Just Integer
r -> Integer
r

algebraicSquareRoot ::
    Integer -> Zx -> Integer -> Zx -> [Ideal] -> [Nfzw] ->
    (Integer -> Bool) -> Verbose (Maybe Integer)
algebraicSquareRoot :: Integer
-> Zx
-> Integer
-> Zx
-> [(Integer, Integer)]
-> [Nfzw]
-> (Integer -> Bool)
-> Verbose (Maybe Integer)
algebraicSquareRoot Integer
n Zx
f Integer
m Zx
f' [(Integer, Integer)]
qil0 [Nfzw]
sq Integer -> Bool
sameSquare =
    case (Integer -> Maybe (Integer, Maybe [(Integer, Integer)]))
-> [Integer] -> [(Integer, Maybe [(Integer, Integer)])]
forall a b. (a -> Maybe b) -> [a] -> [b]
mapMaybe Integer -> Maybe (Integer, Maybe [(Integer, Integer)])
splits ([Integer] -> [Integer]
forall a. [a] -> [a]
tail [Integer]
Prime.primes) of
      [] ->
          String -> Verbose (Maybe Integer)
forall a. HasCallStack => String -> a
error String
"ran out of primes"
      (Integer
p,Maybe [(Integer, Integer)]
Nothing) : [(Integer, Maybe [(Integer, Integer)])]
_ -> do
          String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Working modulo " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
p String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    String
" shows that no algebraic square root exists"
          Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure Maybe Integer
forall a. Maybe a
Nothing
      (Integer
p, Just [(Integer, Integer)]
rcs) : [(Integer, Maybe [(Integer, Integer)])]
_ -> do
          String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Reducing modulo prime " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
p String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    String
"\n  totally splits f(x) as " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    ((Integer, Integer) -> String) -> [(Integer, Integer)] -> String
forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\(Integer
r,Integer
_) -> String
"(x - " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
r String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
")") [(Integer, Integer)]
rcs String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    String
"\n  and algebraic square root is " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    [Integer] -> String
forall a. Show a => a -> String
show (((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> Integer -> Integer
Prime.toSmallestInteger Integer
p (Integer -> Integer)
-> ((Integer, Integer) -> Integer) -> (Integer, Integer) -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer, Integer) -> Integer
forall a b. (a, b) -> b
snd) [(Integer, Integer)]
rcs)
          Integer
-> Int
-> Integer
-> [(Integer, Integer)]
-> [[(Integer, Integer)]]
-> Verbose (Maybe Integer)
forall a.
(Num a, Show a) =>
Integer
-> a
-> Integer
-> [(Integer, Integer)]
-> [[(Integer, Integer)]]
-> Verbose (Maybe Integer)
checkSqrt Integer
p (Int
1 :: Int) Integer
p [(Integer, Integer)]
qil0 ([[(Integer, Integer)]] -> Verbose (Maybe Integer))
-> [[(Integer, Integer)]] -> Verbose (Maybe Integer)
forall a b. (a -> b) -> a -> b
$ [[(Integer, Integer)]] -> [[(Integer, Integer)]]
forall a. [[a]] -> [[a]]
List.transpose ([[(Integer, Integer)]] -> [[(Integer, Integer)]])
-> [[(Integer, Integer)]] -> [[(Integer, Integer)]]
forall a b. (a -> b) -> a -> b
$ ((Integer, Integer) -> [(Integer, Integer)])
-> [(Integer, Integer)] -> [[(Integer, Integer)]]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> (Integer, Integer) -> [(Integer, Integer)]
liftSqrt Integer
p) [(Integer, Integer)]
rcs
  where
    evalSquare :: Integer -> Integer -> Integer
evalSquare Integer
pk Integer
x = Integer -> [Integer] -> Integer
Prime.product Integer
pk (Integer
fx2' Integer -> [Integer] -> [Integer]
forall a. a -> [a] -> [a]
: [Integer]
sqx)
      where
        fx2' :: Integer
fx2' = Integer -> Integer -> Integer
Prime.square Integer
pk (Integer -> Gfpx -> Integer -> Integer
Gfpx.evaluate Integer
pk (Integer -> Zx -> Gfpx
Gfpx.fromZx Integer
pk Zx
f') Integer
x)
        sqx :: [Integer]
sqx = (Nfzw -> Integer) -> [Nfzw] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> Integer -> Integer
Prime.fromInteger Integer
pk (Integer -> Integer) -> (Nfzw -> Integer) -> Nfzw -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Integer -> Nfzw -> Integer
Nfzw.toInteger Integer
x) [Nfzw]
sq

    splits :: Integer -> Maybe (Integer, Maybe [(Integer, Integer)])
splits Integer
p = do
        [Integer]
rs <- Zx -> Integer -> Maybe [Integer]
Gfpx.totallySplits Zx
f Integer
p
        [Integer]
ss <- [Integer] -> Maybe [Integer]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Integer] -> Maybe [Integer]) -> [Integer] -> Maybe [Integer]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> [Integer] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> Integer -> Integer
evalSquare Integer
p) [Integer]
rs
        if (Integer -> Bool) -> [Integer] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
any (Integer -> Integer -> Bool
Prime.nonResidue Integer
p) [Integer]
ss then (Integer, Maybe [(Integer, Integer)])
-> Maybe (Integer, Maybe [(Integer, Integer)])
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer
p,Maybe [(Integer, Integer)]
forall a. Maybe a
Nothing) else do
          [Integer]
cs <- [Integer] -> Maybe [Integer]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Integer] -> Maybe [Integer]) -> [Integer] -> Maybe [Integer]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> [Integer] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer -> Integer -> Integer
Prime.sqrt Integer
p) [Integer]
ss
          Bool -> Maybe ()
forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Bool -> Maybe ()) -> Bool -> Maybe ()
forall a b. (a -> b) -> a -> b
$ (Integer -> Bool) -> [Integer] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
(/=) Integer
0) [Integer]
cs
          (Integer, Maybe [(Integer, Integer)])
-> Maybe (Integer, Maybe [(Integer, Integer)])
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer
p, [(Integer, Integer)] -> Maybe [(Integer, Integer)]
forall a. a -> Maybe a
Just ([Integer] -> [Integer] -> [(Integer, Integer)]
forall a b. [a] -> [b] -> [(a, b)]
zip [Integer]
rs [Integer]
cs))

    liftSqrt :: Integer -> (Integer, Integer) -> [(Integer, Integer)]
liftSqrt Integer
p (Integer
r,Integer
c) =
        [Integer] -> [Integer] -> [(Integer, Integer)]
forall a b. [a] -> [b] -> [(a, b)]
zip (((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> Integer
forall a b. (a, b) -> b
snd [(Integer, Integer)]
rks) ((Integer, [Integer]) -> [Integer]
forall a b. (a, b) -> b
snd ((Integer, [Integer]) -> [Integer])
-> (Integer, [Integer]) -> [Integer]
forall a b. (a -> b) -> a -> b
$ (Integer -> (Integer, Integer) -> (Integer, Integer))
-> Integer -> [(Integer, Integer)] -> (Integer, [Integer])
forall (t :: * -> *) a b c.
Traversable t =>
(a -> b -> (a, c)) -> a -> t b -> (a, t c)
List.mapAccumL Integer -> (Integer, Integer) -> (Integer, Integer)
lift Integer
c ([(Integer, Integer)] -> [(Integer, Integer)]
forall a. [a] -> [a]
tail [(Integer, Integer)]
rks))
      where
        lift :: Integer -> (Integer, Integer) -> (Integer, Integer)
lift Integer
ck (Integer
pk,Integer
rk) = ((Integer
ck Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- ((Integer
ckInteger -> Integer -> Integer
forall a. Num a => a -> a -> a
*Integer
ck Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- Integer -> Integer -> Integer
evalSquare Integer
pk Integer
rk) Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
a)) Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
`mod` Integer
pk, Integer
ck)
        rks :: [(Integer, Integer)]
rks = Zx -> Integer -> Integer -> [(Integer, Integer)]
Gfpx.liftRoot Zx
f Integer
p Integer
r
        a :: Integer
a = Integer -> Integer -> Integer
Prime.invert Integer
p (Integer -> Integer -> Integer -> Integer
Prime.multiply Integer
p Integer
2 Integer
c)

    checkSqrt :: Integer
-> a
-> Integer
-> [(Integer, Integer)]
-> [[(Integer, Integer)]]
-> Verbose (Maybe Integer)
checkSqrt Integer
_ a
_ Integer
_ [(Integer, Integer)]
_ [] = String -> Verbose (Maybe Integer)
forall a. HasCallStack => String -> a
error String
"ran out of lifted square roots"
    checkSqrt Integer
p a
k Integer
pk [(Integer, Integer)]
qil ([(Integer, Integer)]
rcks : [[(Integer, Integer)]]
rcksl) =
        case (Integer -> Bool) -> [Integer] -> [Integer]
forall a. (a -> Bool) -> [a] -> [a]
filter Integer -> Bool
sameSquare ([Integer] -> [Integer]) -> [Integer] -> [Integer]
forall a b. (a -> b) -> a -> b
$ ((Integer, [(Integer, Integer)]) -> Integer)
-> [(Integer, [(Integer, Integer)])] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, [(Integer, Integer)]) -> Integer
crtSqrt ([(Integer, [(Integer, Integer)])] -> [Integer])
-> [(Integer, [(Integer, Integer)])] -> [Integer]
forall a b. (a -> b) -> a -> b
$ Integer
-> [(Integer, Integer)] -> [(Integer, [(Integer, Integer)])]
forall b a. Num b => b -> [(a, b)] -> [(b, [(a, b)])]
allSigns Integer
pk [(Integer, Integer)]
rcks of
          [] -> if Bool
ok then Integer
-> a
-> Integer
-> [(Integer, Integer)]
-> [[(Integer, Integer)]]
-> Verbose (Maybe Integer)
checkSqrt Integer
p (a
k a -> a -> a
forall a. Num a => a -> a -> a
+ a
1) (Integer
p Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
pk) [(Integer, Integer)]
qil' [[(Integer, Integer)]]
rcksl else do
              String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Quadratic character " String -> ShowS
forall a. [a] -> [a] -> [a]
++ (Integer, Integer) -> String
forall a. Show a => a -> String
show (Integer, Integer)
qc String -> ShowS
forall a. [a] -> [a] -> [a]
++
                        String
" shows that no algebraic square root exists"
              Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure Maybe Integer
forall a. Maybe a
Nothing
            where
              ((Integer, Integer)
qc,[(Integer, Integer)]
qil') = Zx
-> [Nfzw]
-> [(Integer, Integer)]
-> ((Integer, Integer), [(Integer, Integer)])
nextQuadraticCharacter Zx
f' [Nfzw]
sq [(Integer, Integer)]
qil
              ok :: Bool
ok = (Integer, Integer) -> [Nfzw] -> Bool
productIsQuadraticResidue (Integer, Integer)
qc [Nfzw]
sq
          Integer
asq : [Integer]
_ -> do
              String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Lifted algebraic square root modulo " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                        Integer -> String
forall a. Show a => a -> String
show Integer
p String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
"^" String -> ShowS
forall a. [a] -> [a] -> [a]
++ a -> String
forall a. Show a => a -> String
show a
k String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" has same square modulo n"
              Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Maybe Integer -> Verbose (Maybe Integer))
-> Maybe Integer -> Verbose (Maybe Integer)
forall a b. (a -> b) -> a -> b
$ Integer -> Maybe Integer
forall a. a -> Maybe a
Just Integer
asq

    allSigns :: b -> [(a, b)] -> [(b, [(a, b)])]
allSigns b
_ [] = String -> [(b, [(a, b)])]
forall a. HasCallStack => String -> a
error String
"no roots"
    allSigns b
pk ((a, b)
rck : [(a, b)]
rcks) = ([(a, b)] -> (b, [(a, b)])) -> [[(a, b)]] -> [(b, [(a, b)])]
forall a b. (a -> b) -> [a] -> [b]
map (\[(a, b)]
z -> (b
pk, (a, b)
rck (a, b) -> [(a, b)] -> [(a, b)]
forall a. a -> [a] -> [a]
: [(a, b)]
z)) ([[(a, b)]] -> [(b, [(a, b)])]) -> [[(a, b)]] -> [(b, [(a, b)])]
forall a b. (a -> b) -> a -> b
$ ((a, b) -> [[(a, b)]] -> [[(a, b)]])
-> [[(a, b)]] -> [(a, b)] -> [[(a, b)]]
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr (a, b) -> [[(a, b)]] -> [[(a, b)]]
forall a. (a, b) -> [[(a, b)]] -> [[(a, b)]]
pm [[]] [(a, b)]
rcks
      where pm :: (a, b) -> [[(a, b)]] -> [[(a, b)]]
pm (a
rk,b
ck) [[(a, b)]]
z = ([(a, b)] -> [(a, b)]) -> [[(a, b)]] -> [[(a, b)]]
forall a b. (a -> b) -> [a] -> [b]
map ((:) (a
rk,b
ck)) [[(a, b)]]
z [[(a, b)]] -> [[(a, b)]] -> [[(a, b)]]
forall a. [a] -> [a] -> [a]
++ ([(a, b)] -> [(a, b)]) -> [[(a, b)]] -> [[(a, b)]]
forall a b. (a -> b) -> [a] -> [b]
map ((:) (a
rk, b
pk b -> b -> b
forall a. Num a => a -> a -> a
- b
ck)) [[(a, b)]]
z

    crtSqrt :: (Integer, [(Integer, Integer)]) -> Integer
crtSqrt (Integer
pk,[(Integer, Integer)]
rcks) =
        Integer -> Integer -> Integer
Prime.fromInteger Integer
n (Integer -> Integer) -> Integer -> Integer
forall a b. (a -> b) -> a -> b
$
        Integer -> Integer -> Integer
Prime.toSmallestInteger Integer
pk (Integer -> Integer) -> Integer -> Integer
forall a b. (a -> b) -> a -> b
$
        Integer -> Integer -> Integer -> Integer
Prime.multiply Integer
pk (Integer -> Integer -> Integer
Prime.fromInteger Integer
pk Integer
n) (Integer -> Integer) -> Integer -> Integer
forall a b. (a -> b) -> a -> b
$
        Integer -> [Integer] -> Integer
Prime.sum Integer
pk (((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> Integer
evalM [(Integer, Integer)]
rcks)
      where
        evalM :: (Integer, Integer) -> Integer
evalM (Integer
rk,Integer
ck) =
            Integer -> Integer -> Integer -> Integer
Prime.multiply Integer
pk Integer
ck (Integer -> Integer) -> Integer -> Integer
forall a b. (a -> b) -> a -> b
$
            Integer -> Integer -> Integer
Prime.invert Integer
pk (Integer -> Integer -> Integer -> Integer
Prime.multiply Integer
pk Integer
m_rk Integer
frk')
          where
            m_rk :: Integer
m_rk = Integer -> Integer -> Integer
Prime.fromInteger Integer
pk (Integer
m Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- Integer
rk)
            frk' :: Integer
frk' = Integer -> Gfpx -> Integer -> Integer
Gfpx.evaluate Integer
pk (Integer -> Zx -> Gfpx
Gfpx.fromZx Integer
pk Zx
f') Integer
rk

-------------------------------------------------------------------------------
-- The number field sieve
--
-- Input: An odd integer greater than 5
-- Output: Either a nontrivial factor of the input or failure
-------------------------------------------------------------------------------

data Config =
    Config
      {Config -> PolynomialConfig
polynomialConfig :: PolynomialConfig,
       Config -> FactorBaseConfig
rationalFactorBaseConfig :: FactorBaseConfig,
       Config -> FactorBaseConfig
algebraicFactorBaseConfig :: FactorBaseConfig,
       Config -> QuadraticCharacterConfig
quadraticCharacterConfig :: QuadraticCharacterConfig,
       Config -> Int
extraRankConfig :: Int,
       Config -> Bool
verboseConfig :: Bool}
  deriving (Config -> Config -> Bool
(Config -> Config -> Bool)
-> (Config -> Config -> Bool) -> Eq Config
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Config -> Config -> Bool
$c/= :: Config -> Config -> Bool
== :: Config -> Config -> Bool
$c== :: Config -> Config -> Bool
Eq,Eq Config
Eq Config
-> (Config -> Config -> Ordering)
-> (Config -> Config -> Bool)
-> (Config -> Config -> Bool)
-> (Config -> Config -> Bool)
-> (Config -> Config -> Bool)
-> (Config -> Config -> Config)
-> (Config -> Config -> Config)
-> Ord Config
Config -> Config -> Bool
Config -> Config -> Ordering
Config -> Config -> Config
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: Config -> Config -> Config
$cmin :: Config -> Config -> Config
max :: Config -> Config -> Config
$cmax :: Config -> Config -> Config
>= :: Config -> Config -> Bool
$c>= :: Config -> Config -> Bool
> :: Config -> Config -> Bool
$c> :: Config -> Config -> Bool
<= :: Config -> Config -> Bool
$c<= :: Config -> Config -> Bool
< :: Config -> Config -> Bool
$c< :: Config -> Config -> Bool
compare :: Config -> Config -> Ordering
$ccompare :: Config -> Config -> Ordering
$cp1Ord :: Eq Config
Ord,Int -> Config -> ShowS
[Config] -> ShowS
Config -> String
(Int -> Config -> ShowS)
-> (Config -> String) -> ([Config] -> ShowS) -> Show Config
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Config] -> ShowS
$cshowList :: [Config] -> ShowS
show :: Config -> String
$cshow :: Config -> String
showsPrec :: Int -> Config -> ShowS
$cshowsPrec :: Int -> Config -> ShowS
Show)

defaultConfig :: Config
defaultConfig :: Config
defaultConfig =
    Config :: PolynomialConfig
-> FactorBaseConfig
-> FactorBaseConfig
-> QuadraticCharacterConfig
-> Int
-> Bool
-> Config
Config
      {polynomialConfig :: PolynomialConfig
polynomialConfig = PolynomialConfig
defaultPolynomialConfig,
       rationalFactorBaseConfig :: FactorBaseConfig
rationalFactorBaseConfig = Double -> FactorBaseConfig
OptimalFactorBase Double
3.0,
       algebraicFactorBaseConfig :: FactorBaseConfig
algebraicFactorBaseConfig = Double -> FactorBaseConfig
OptimalFactorBase Double
10.0,
       quadraticCharacterConfig :: QuadraticCharacterConfig
quadraticCharacterConfig = QuadraticCharacterConfig
defaultQuadraticCharacterConfig,
       extraRankConfig :: Int
extraRankConfig = Int
5,
       verboseConfig :: Bool
verboseConfig = Bool
False}

setVerboseConfig :: Bool -> Config -> Config
setVerboseConfig :: Bool -> Config -> Config
setVerboseConfig Bool
v Config
cfg = Config
cfg {verboseConfig :: Bool
verboseConfig = Bool
v}

setQuadraticCharacterConfig :: Maybe Int -> Config -> Config
setQuadraticCharacterConfig :: Maybe Int -> Config -> Config
setQuadraticCharacterConfig Maybe Int
Nothing Config
cfg = Config
cfg
setQuadraticCharacterConfig (Just Int
q) Config
cfg =
    Config
cfg {quadraticCharacterConfig :: QuadraticCharacterConfig
quadraticCharacterConfig = Int -> QuadraticCharacterConfig
FixedQuadraticCharacterConfig Int
q}

verboseList :: Config -> String -> [String] -> String
verboseList :: Config -> String -> [String] -> String
verboseList Config
cfg String
s = if Config -> Bool
verboseConfig Config
cfg then [String] -> String
unabbrevList else String -> [String] -> String
abbrevList String
s

factorSquareRoots ::
    Config -> Integer -> Zx -> Integer -> FactorBase -> Zx -> [Ideal] ->
    [[Nfzw]] -> Verbose (Maybe Integer)
factorSquareRoots :: Config
-> Integer
-> Zx
-> Integer
-> [Integer]
-> Zx
-> [(Integer, Integer)]
-> [[Nfzw]]
-> Verbose (Maybe Integer)
factorSquareRoots Config
_ Integer
_ Zx
_ Integer
_ [Integer]
_ Zx
_ [(Integer, Integer)]
_ [] = do
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"No more square products, NFS factorization failed"
    Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure Maybe Integer
forall a. Maybe a
Nothing
factorSquareRoots Config
cfg Integer
n Zx
f Integer
m [Integer]
rfb Zx
f' [(Integer, Integer)]
qil ([Nfzw]
sq : [[Nfzw]]
sqs) = do
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Considering square product " String -> ShowS
forall a. [a] -> [a] -> [a]
++
              (if [Nfzw] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Nfzw]
sq Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 then String
"consisting of a single element of Z[w]"
               else String
"of " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show ([Nfzw] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Nfzw]
sq) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" elements of Z[w]") String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
":" String -> ShowS
forall a. [a] -> [a] -> [a]
++ Config -> String -> [String] -> String
verboseList Config
cfg String
"elements" ((Nfzw -> String) -> [Nfzw] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map Nfzw -> String
forall a. Show a => a -> String
show [Nfzw]
sq)
    Integer
rsq <- Integer -> Verbose Integer
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer -> Verbose Integer) -> Integer -> Verbose Integer
forall a b. (a -> b) -> a -> b
$ Integer -> Integer -> Zx -> [Integer] -> [Nfzw] -> Integer
rationalSquareRoot Integer
n Integer
m Zx
f' [Integer]
rfb [Nfzw]
sq
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Rational square root modulo n is " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
rsq
    Integer -> Bool
sameSquare <- (Integer -> Bool) -> Verbose (Integer -> Bool)
forall (f :: * -> *) a. Applicative f => a -> f a
pure ((Integer -> Bool) -> Verbose (Integer -> Bool))
-> (Integer -> Bool) -> Verbose (Integer -> Bool)
forall a b. (a -> b) -> a -> b
$ Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
(==) (Integer -> Integer -> Integer
Prime.square Integer
n Integer
rsq) (Integer -> Bool) -> (Integer -> Integer) -> Integer -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Integer -> Integer -> Integer
Prime.square Integer
n
    Maybe Integer
asqm <- Integer
-> Zx
-> Integer
-> Zx
-> [(Integer, Integer)]
-> [Nfzw]
-> (Integer -> Bool)
-> Verbose (Maybe Integer)
algebraicSquareRoot Integer
n Zx
f Integer
m Zx
f' [(Integer, Integer)]
qil [Nfzw]
sq Integer -> Bool
sameSquare
    case Maybe Integer
asqm of
      Maybe Integer
Nothing -> Config
-> Integer
-> Zx
-> Integer
-> [Integer]
-> Zx
-> [(Integer, Integer)]
-> [[Nfzw]]
-> Verbose (Maybe Integer)
factorSquareRoots Config
cfg Integer
n Zx
f Integer
m [Integer]
rfb Zx
f' [(Integer, Integer)]
qil [[Nfzw]]
sqs
      Just Integer
asq -> do
        String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Algebraic square root modulo n is " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
asq
        Integer
g <- Integer -> Verbose Integer
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer -> Verbose Integer) -> Integer -> Verbose Integer
forall a b. (a -> b) -> a -> b
$ Integer -> Integer -> Integer
forall a. Integral a => a -> a -> a
gcd Integer
n (Integer
rsq Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
asq)
        Bool
s <- Bool -> Verbose Bool
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Bool -> Verbose Bool) -> Bool -> Verbose Bool
forall a b. (a -> b) -> a -> b
$ Integer
1 Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
g Bool -> Bool -> Bool
&& Integer
g Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
n
        String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Greatest common divisor of n and " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                  String
"sum of square roots is " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                  (if Integer
g Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
n then String
"n" else Integer -> String
forall a. Show a => a -> String
show Integer
g) String -> ShowS
forall a. [a] -> [a] -> [a]
++
                  (if Bool
s then String
"" else String
" (bad luck)")
        if Bool
s then Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer -> Maybe Integer
forall a. a -> Maybe a
Just Integer
g) else Config
-> Integer
-> Zx
-> Integer
-> [Integer]
-> Zx
-> [(Integer, Integer)]
-> [[Nfzw]]
-> Verbose (Maybe Integer)
factorSquareRoots Config
cfg Integer
n Zx
f Integer
m [Integer]
rfb Zx
f' [(Integer, Integer)]
qil [[Nfzw]]
sqs

factorWithPolynomial ::
    Config -> Integer -> Zx -> Integer -> Verbose (Maybe Integer)
factorWithPolynomial :: Config -> Integer -> Zx -> Integer -> Verbose (Maybe Integer)
factorWithPolynomial Config
cfg Integer
n Zx
f Integer
m = do
    Integer
rfm <- Integer -> Verbose Integer
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer -> Verbose Integer) -> Integer -> Verbose Integer
forall a b. (a -> b) -> a -> b
$ FactorBaseConfig -> Integer -> Integer
maxFactorBase (Config -> FactorBaseConfig
rationalFactorBaseConfig Config
cfg) Integer
n
    [Integer]
rfb <- [Integer] -> Verbose [Integer]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Integer] -> Verbose [Integer]) -> [Integer] -> Verbose [Integer]
forall a b. (a -> b) -> a -> b
$ (Integer -> Bool) -> [Integer] -> [Integer]
forall a. (a -> Bool) -> [a] -> [a]
takeWhile (Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
(>=) Integer
rfm) [Integer]
Prime.primes
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Rational factor base contains " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show ([Integer] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Integer]
rfb) String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
" prime integers:" String -> ShowS
forall a. [a] -> [a] -> [a]
++ Config -> String -> [String] -> String
verboseList Config
cfg String
"primes" ((Integer -> String) -> [Integer] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map Integer -> String
forall a. Show a => a -> String
show [Integer]
rfb)
    Integer
afm <- Integer -> Verbose Integer
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Integer -> Verbose Integer) -> Integer -> Verbose Integer
forall a b. (a -> b) -> a -> b
$ FactorBaseConfig -> Integer -> Integer
maxFactorBase (Config -> FactorBaseConfig
algebraicFactorBaseConfig Config
cfg) Integer
n
    ([(Integer, Integer)]
ail,[(Integer, Integer)]
qil) <- ([(Integer, Integer)], [(Integer, Integer)])
-> Verbose ([(Integer, Integer)], [(Integer, Integer)])
forall (f :: * -> *) a. Applicative f => a -> f a
pure (([(Integer, Integer)], [(Integer, Integer)])
 -> Verbose ([(Integer, Integer)], [(Integer, Integer)]))
-> ([(Integer, Integer)], [(Integer, Integer)])
-> Verbose ([(Integer, Integer)], [(Integer, Integer)])
forall a b. (a -> b) -> a -> b
$ ((Integer, Integer) -> Bool)
-> [(Integer, Integer)]
-> ([(Integer, Integer)], [(Integer, Integer)])
forall a. (a -> Bool) -> [a] -> ([a], [a])
span (Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
(>=) Integer
afm (Integer -> Bool)
-> ((Integer, Integer) -> Integer) -> (Integer, Integer) -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Integer, Integer) -> Integer
forall a b. (a, b) -> b
snd) (Zx -> [(Integer, Integer)]
Nfzw.ideals Zx
f)
    [Integer]
afb <- [Integer] -> Verbose [Integer]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Integer] -> Verbose [Integer]) -> [Integer] -> Verbose [Integer]
forall a b. (a -> b) -> a -> b
$ ([Integer] -> Integer) -> [[Integer]] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map [Integer] -> Integer
forall a. [a] -> a
head ([[Integer]] -> [Integer]) -> [[Integer]] -> [Integer]
forall a b. (a -> b) -> a -> b
$ [Integer] -> [[Integer]]
forall a. Eq a => [a] -> [[a]]
List.group ([Integer] -> [[Integer]]) -> [Integer] -> [[Integer]]
forall a b. (a -> b) -> a -> b
$ ((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> Integer
forall a b. (a, b) -> b
snd [(Integer, Integer)]
ail
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Algebraic factor base contains " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show ([(Integer, Integer)] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [(Integer, Integer)]
ail) String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
" first degree prime ideals:\n" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
"  (r,p) such that f(r) = 0 (mod p)" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Config -> String -> [String] -> String
verboseList Config
cfg String
"prime ideals" (((Integer, Integer) -> String) -> [(Integer, Integer)] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> String
forall a. Show a => a -> String
show [(Integer, Integer)]
ail)
    Int
qcs <- Int -> Verbose Int
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Int -> Verbose Int) -> Int -> Verbose Int
forall a b. (a -> b) -> a -> b
$ QuadraticCharacterConfig -> Integer -> Int
quadraticCharacters (Config -> QuadraticCharacterConfig
quadraticCharacterConfig Config
cfg) Integer
n
    Int
extra <- Int -> Verbose Int
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Int -> Verbose Int) -> Int -> Verbose Int
forall a b. (a -> b) -> a -> b
$ Config -> Int
extraRankConfig Config
cfg
    Int
cols <- Int -> Verbose Int
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Int -> Verbose Int) -> Int -> Verbose Int
forall a b. (a -> b) -> a -> b
$ Int
1 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ [Integer] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Integer]
rfb Int -> Int -> Int
forall a. Num a => a -> a -> a
+ [(Integer, Integer)] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [(Integer, Integer)]
ail Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
qcs
    [Nfzw]
xs <- [Nfzw] -> Verbose [Nfzw]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Nfzw] -> Verbose [Nfzw]) -> [Nfzw] -> Verbose [Nfzw]
forall a b. (a -> b) -> a -> b
$ Int -> [Nfzw] -> [Nfzw]
forall a. Int -> [a] -> [a]
take (Int
cols Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
extra) (Zx -> Integer -> [Integer] -> [Integer] -> [Nfzw]
smoothNfzw Zx
f Integer
m [Integer]
rfb [Integer]
afb)
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Searching for 1+" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Int -> String
forall a. Show a => a -> String
show ([Integer] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [Integer]
rfb) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
"+" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Int -> String
forall a. Show a => a -> String
show ([(Integer, Integer)] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [(Integer, Integer)]
ail) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
"+" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Int -> String
forall a. Show a => a -> String
show Int
qcs String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
"+" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Int -> String
forall a. Show a => a -> String
show Int
extra String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" = " String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Int -> String
forall a. Show a => a -> String
show (Int
cols Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
extra) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" smooth elements of Z[w]:\n" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
"  a + bw |-> (a + bm, (-b)^degree(f) * f(a/(-b)))" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Config -> String -> [String] -> String
verboseList Config
cfg String
"smooth elements"
                ((Nfzw -> String) -> [Nfzw] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map (\Nfzw
x -> Nfzw -> String
forall a. Show a => a -> String
show Nfzw
x String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" |-> (" String -> ShowS
forall a. [a] -> [a] -> [a]
++
                            Integer -> String
forall a. Show a => a -> String
show (Integer -> Nfzw -> Integer
rationalNorm Integer
m Nfzw
x) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
", " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                            Integer -> String
forall a. Show a => a -> String
show (Zx -> Nfzw -> Integer
algebraicNorm Zx
f Nfzw
x) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
")") [Nfzw]
xs)
    Zx
f' <- Zx -> Verbose Zx
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Zx -> Verbose Zx) -> Zx -> Verbose Zx
forall a b. (a -> b) -> a -> b
$ Zx -> Zx
Zx.derivative Zx
f
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Derivative of f is f'(x) = " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Zx -> String
forall a. Show a => a -> String
show Zx
f'
    ([(Integer, Integer)]
qcl,[(Integer, Integer)]
qil') <- ([(Integer, Integer)], [(Integer, Integer)])
-> Verbose ([(Integer, Integer)], [(Integer, Integer)])
forall (f :: * -> *) a. Applicative f => a -> f a
pure (([(Integer, Integer)], [(Integer, Integer)])
 -> Verbose ([(Integer, Integer)], [(Integer, Integer)]))
-> ([(Integer, Integer)], [(Integer, Integer)])
-> Verbose ([(Integer, Integer)], [(Integer, Integer)])
forall a b. (a -> b) -> a -> b
$ Zx
-> [Nfzw]
-> Int
-> [(Integer, Integer)]
-> ([(Integer, Integer)], [(Integer, Integer)])
takeQuadraticCharacters Zx
f' [Nfzw]
xs Int
qcs [(Integer, Integer)]
qil
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Generated " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show Int
qcs String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" quadratic characters " String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
"nonzero for f' and all smooth elements:" String -> ShowS
forall a. [a] -> [a] -> [a]
++
              Config -> String -> [String] -> String
verboseList Config
cfg String
"quadratic characters" (((Integer, Integer) -> String) -> [(Integer, Integer)] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map (Integer, Integer) -> String
forall a. Show a => a -> String
show [(Integer, Integer)]
qcl)
    [Row]
rows <- [Row] -> Verbose [Row]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([Row] -> Verbose [Row]) -> [Row] -> Verbose [Row]
forall a b. (a -> b) -> a -> b
$ (Nfzw -> Row) -> [Nfzw] -> [Row]
forall a b. (a -> b) -> [a] -> [b]
map (Zx
-> Integer
-> [Integer]
-> [(Integer, Integer)]
-> [(Integer, Integer)]
-> Nfzw
-> Row
formRow Zx
f Integer
m [Integer]
rfb [(Integer, Integer)]
ail [(Integer, Integer)]
qcl) [Nfzw]
xs
    [[Nfzw]]
sql <- [[Nfzw]] -> Verbose [[Nfzw]]
forall (f :: * -> *) a. Applicative f => a -> f a
pure ([[Nfzw]] -> Verbose [[Nfzw]]) -> [[Nfzw]] -> Verbose [[Nfzw]]
forall a b. (a -> b) -> a -> b
$ ([Int] -> [Nfzw]) -> [[Int]] -> [[Nfzw]]
forall a b. (a -> b) -> [a] -> [b]
map ((Int -> Nfzw) -> [Int] -> [Nfzw]
forall a b. (a -> b) -> [a] -> [b]
map (\Int
i -> [Nfzw]
xs [Nfzw] -> Int -> Nfzw
forall a. [a] -> Int -> a
!! Int
i)) ([Row] -> [[Int]]
gaussianElimination [Row]
rows)
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Gaussian elimination resulted in " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show ([[Nfzw]] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [[Nfzw]]
sql) String -> ShowS
forall a. [a] -> [a] -> [a]
++
              String
" square products"
    Config
-> Integer
-> Zx
-> Integer
-> [Integer]
-> Zx
-> [(Integer, Integer)]
-> [[Nfzw]]
-> Verbose (Maybe Integer)
factorSquareRoots Config
cfg Integer
n Zx
f Integer
m [Integer]
rfb Zx
f' [(Integer, Integer)]
qil' [[Nfzw]]
sql

factor :: Config -> Integer -> Verbose (Maybe Integer)
factor :: Config -> Integer -> Verbose (Maybe Integer)
factor Config
cfg Integer
n = do
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"NFS configuration = " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Config -> String
forall a. Show a => a -> String
show Config
cfg
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Attempting to factor n = " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
n
    (Zx
f,Integer
m) <- (Zx, Integer) -> Verbose (Zx, Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure ((Zx, Integer) -> Verbose (Zx, Integer))
-> (Zx, Integer) -> Verbose (Zx, Integer)
forall a b. (a -> b) -> a -> b
$ PolynomialConfig -> Integer -> (Zx, Integer)
selectPolynomial (Config -> PolynomialConfig
polynomialConfig Config
cfg) Integer
n
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Working in Z[w] where w is a complex root of f and f(m) = n"
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"  where f(x) = " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Zx -> String
forall a. Show a => a -> String
show Zx
f
    String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"  and m = " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Integer -> String
forall a. Show a => a -> String
show Integer
m
    case Zx -> (Integer, [(Zx, Integer)])
Bz.factor Zx
f of
      (Integer
1,[(Zx
_,Integer
1)]) -> do
          String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Verified that f(x) is irreducible in Z[x]"
          Config -> Integer -> Zx -> Integer -> Verbose (Maybe Integer)
factorWithPolynomial Config
cfg Integer
n Zx
f Integer
m
      (Integer
_,[(Zx, Integer)]
gks) -> do
          String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Factored f(x) in Z[x] as the product of:" String -> ShowS
forall a. [a] -> [a] -> [a]
++ [(Zx, Integer)] -> String
showProd [(Zx, Integer)]
gks
          String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Evaluating at m factors n in Z as the product of:" String -> ShowS
forall a. [a] -> [a] -> [a]
++
                    [(Zx, Integer)] -> String
showProd (((Zx, Integer) -> (Zx, Integer))
-> [(Zx, Integer)] -> [(Zx, Integer)]
forall a b. (a -> b) -> [a] -> [b]
map (Zx, Integer) -> (Zx, Integer)
forall b. (Zx, b) -> (Zx, b)
evalFacAtM [(Zx, Integer)]
gks)
          if Integer
x Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
n then Maybe Integer -> Verbose (Maybe Integer)
forall (f :: * -> *) a. Applicative f => a -> f a
pure (Maybe Integer -> Verbose (Maybe Integer))
-> Maybe Integer -> Verbose (Maybe Integer)
forall a b. (a -> b) -> a -> b
$ Integer -> Maybe Integer
forall a. a -> Maybe a
Just Integer
x else do
              String -> Verbose ()
comment (String -> Verbose ()) -> String -> Verbose ()
forall a b. (a -> b) -> a -> b
$ String
"Adjust f(x) to " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Zx -> String
forall a. Show a => a -> String
show Zx
h
              Config -> Integer -> Zx -> Integer -> Verbose (Maybe Integer)
factorWithPolynomial Config
cfg Integer
n Zx
h Integer
m
        where
          (Integer
x,Zx
h) = [(Integer, Zx)] -> (Integer, Zx)
forall (t :: * -> *) a. (Foldable t, Ord a) => t a -> a
maximum ([(Integer, Zx)] -> (Integer, Zx))
-> [(Integer, Zx)] -> (Integer, Zx)
forall a b. (a -> b) -> a -> b
$ ((Zx, Integer) -> (Integer, Zx))
-> [(Zx, Integer)] -> [(Integer, Zx)]
forall a b. (a -> b) -> [a] -> [b]
map (Zx, Integer) -> (Integer, Zx)
forall b. (Zx, b) -> (Integer, Zx)
absEvalAtM [(Zx, Integer)]
gks
          evalAtM :: Zx -> Integer
evalAtM Zx
g = Zx -> Integer -> Integer
Zx.evaluate Zx
g Integer
m
          absEvalAtM :: (Zx, b) -> (Integer, Zx)
absEvalAtM (Zx
g,b
_) = (Integer -> Integer
forall a. Num a => a -> a
abs (Zx -> Integer
evalAtM Zx
g), Zx
g)
          evalFacAtM :: (Zx, b) -> (Zx, b)
evalFacAtM (Zx
g,b
k) = (Integer -> Zx
Zx.constant (Zx -> Integer
evalAtM Zx
g), b
k)
          showProd :: [(Zx, Integer)] -> String
showProd = ((Zx, Integer) -> String) -> [(Zx, Integer)] -> String
forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\(Zx, Integer)
gk -> String
"\n  " String -> ShowS
forall a. [a] -> [a] -> [a]
++ (Zx, Integer) -> String
forall a a. (Eq a, Num a, Show a, Show a) => (a, a) -> String
showFac (Zx, Integer)
gk)
          showFac :: (a, a) -> String
showFac (a
g,a
1) = a -> String
forall a. Show a => a -> String
show a
g
          showFac (a
g,a
k) = String
"(" String -> ShowS
forall a. [a] -> [a] -> [a]
++ a -> String
forall a. Show a => a -> String
show a
g String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
") ^ " String -> ShowS
forall a. [a] -> [a] -> [a]
++ a -> String
forall a. Show a => a -> String
show a
k

{-
ghci Nfs.hs
:break algebraicFactorBase
factor defaultConfig (3119 * 4261)
factor defaultConfig (48029 * 57641)
factor defaultConfig (3119 * 48029 * 57641)
factor defaultConfig (10712543 * 13777067)
factor defaultConfig (2616707501 * 3620408573)
-}