Polynomial factorization in Haskell

Question

With hammar's help I have made a template Haskell bit which compiles

$(zModP 5)

to

newtype Z5 = Z5 Int
instance Additive.C Z5 where
  (Z5 x) + (Z5 y) = Z5 $ (x + y) `mod` 5
...

I'm now facing a problem that I don't think I can solve this way.

A remarkable fact about polynomials is that they are irreducible in the rationals if they are irreducible modulo some prime p. I already have a method which brute-force attempts to factor polynomials over a given (finite) field.

I want to try running this function for multiple fields. Here's kind of what I want:

isIrreducible :: (FiniteField.C a) => Poly.T a -> Bool
isIrreducible p = ...

intPolyIrreducible :: Poly.T Int -> Bool
intPolyIrreducible p = isIrreducible (p :: Poly.T Z2) ||
                       isIrreducible (p :: Poly.T Z3) ||
                       isIrreducible (p :: Poly.T Z5) ||
                       ...

Basically I want to try running my factoring algorithm for a large number of definitions of "division".

I think this is possible to do with TH, but it seems like it would take forever. I'm wondering if it would be easier to just pass my arithmetical operations in as a parameter to isIrreducible?

Alternatively it seems like this might be something the Newtype module could help with, but I can't think of how it would work without using TH in a way which would be just as hard...

Anyone have any thoughts on how best to accomplish this?

Answer 1

You can do computations in finite fields using type-level numerics, for example with the type-level package:

{-# LANGUAGE ScopedTypeVariables #-}
module Mod where
import Data.TypeLevel.Num (Nat,toNum, reifyIntegral)

data Z p = Z Integer

instance Eq (Z p) where Z x == Z y = x == y
instance Ord (Z p) where -- dummy instance
instance Show (Z p) where show (Z n) = show n

instance Nat p => Num (Z p) where
    Z x + Z y = Z $ (x + y) `mod` toNum (undefined :: p)
    Z x - Z y = Z $ (x - y) `mod` toNum (undefined :: p)
    Z x * Z y = Z $ (x * y) `mod` toNum (undefined :: p)
    fromInteger n = Z (n `mod` toNum (undefined :: p))
    -- etc

-- Compute x^2-6 (mod p)
poly :: Nat p => Z p -> Z p
poly x = x*x-6

-- Computes whether x^2-6==0 (mod p), for x=3
checkPoly :: Integer -> Bool
checkPoly n = reifyIntegral n test
  where
    test :: forall p . Nat p => p -> Bool
    test _ = poly (3::Z p) == 0

test1 = map checkPoly [2,3,5]
-- Result: [False,True,False]

This approach has the advantage of not requiring a new template haskell instance for each numeric type. The disadvantage is that it's probably slower than the template haskell solution, since each operation passes the size of the finite field around via a class dictionary.