-- |
-- Module:      Math.NumberTheory.DirichletCharacters
-- Copyright:   (c) 2018 Bhavik Mehta
-- Licence:     MIT
-- Maintainer:  Bhavik Mehta <bhavikmehta8@gmail.com>
--
-- Implementation and enumeration of Dirichlet characters.
--

{-# LANGUAGE CPP                        #-}
{-# LANGUAGE DataKinds                  #-}
{-# LANGUAGE GADTs                      #-}
{-# LANGUAGE GeneralizedNewtypeDeriving #-}
{-# LANGUAGE KindSignatures             #-}
{-# LANGUAGE LambdaCase                 #-}
{-# LANGUAGE PatternSynonyms            #-}
{-# LANGUAGE ScopedTypeVariables        #-}
{-# LANGUAGE ViewPatterns               #-}

module Math.NumberTheory.DirichletCharacters
  (
  -- * An absorbing semigroup
  OrZero, pattern Zero, pattern NonZero
  , orZeroToNum
  -- * Dirichlet characters
  , DirichletCharacter
  -- ** Construction
  , indexToChar
  , indicesToChars
  , characterNumber
  , allChars
  , fromTable
  -- ** Evaluation
  , eval
  , evalGeneral
  , evalAll
  -- ** Special Dirichlet characters
  , principalChar
  , isPrincipal
  , orderChar
  -- ** Real Dirichlet characters
  , RealCharacter
  , isRealCharacter
  , getRealChar
  , toRealFunction
  , jacobiCharacter
  -- ** Primitive characters
  , PrimitiveCharacter
  , isPrimitive
  , getPrimitiveChar
  , induced
  , makePrimitive
  , WithNat(..)
  -- * Roots of unity
  , RootOfUnity(..)
  , toRootOfUnity
  , toComplex
  -- * Debugging
  , validChar
  ) where

#if !MIN_VERSION_base(4,12,0)
import Control.Applicative                                 (liftA2)
#endif
import Data.Bits                                           (Bits(..))
import Data.Constraint
import Data.Foldable
import Data.Functor.Identity                               (Identity(..))
import Data.Kind
import Data.List                                           (sort, unfoldr)
import Data.Maybe                                          (mapMaybe, fromJust, fromMaybe)
import Data.Mod
#if MIN_VERSION_base(4,12,0)
import Data.Monoid                                         (Ap(..))
#endif
import Data.Proxy                                          (Proxy(..))
import Data.Ratio                                          ((%), numerator, denominator)
import Data.Semigroup                                      (Semigroup(..),Product(..))
import Data.Traversable
import qualified Data.Vector as V
import qualified Data.Vector.Mutable as MV
import Data.Vector                                         (Vector, (!))
import GHC.TypeNats                                        (KnownNat, Nat, SomeNat(..), natVal, someNatVal)
import Numeric.Natural                                     (Natural)

import Math.NumberTheory.ArithmeticFunctions               (totient)
import Math.NumberTheory.Moduli.Chinese
import Math.NumberTheory.Moduli.Internal                   (discreteLogarithmPP)
import Math.NumberTheory.Moduli.Multiplicative
import Math.NumberTheory.Moduli.Singleton
import Math.NumberTheory.Primes
import Math.NumberTheory.RootsOfUnity
import Math.NumberTheory.Utils
import Math.NumberTheory.Utils.FromIntegral

-- | A Dirichlet character mod \(n\) is a group homomorphism from \((\mathbb{Z}/n\mathbb{Z})^*\)
-- to \(\mathbb{C}^*\), represented abstractly by `DirichletCharacter`. In particular, they take
-- values at roots of unity and can be evaluated using `eval`.
-- A Dirichlet character can be extended to a completely multiplicative function on \(\mathbb{Z}\)
-- by assigning the value 0 for \(a\) sharing a common factor with \(n\), using `evalGeneral`.
--
-- There are finitely many possible Dirichlet characters for a given modulus, in particular there
-- are \(\phi(n)\) characters modulo \(n\), where \(\phi\) refers to Euler's `totient` function.
-- This gives rise to `Enum` and `Bounded` instances.
newtype DirichletCharacter (n :: Nat) = Generated [DirichletFactor]

-- | The group (Z/nZ)^* decomposes to a product (Z/2^k0 Z)^* x (Z/p1^k1 Z)^* x ... x (Z/pi^ki Z)^*
-- where n = 2^k0 p1^k1 ... pi^ki, and the pj are odd primes, k0 possibly 0. Thus, a group
-- homomorphism from (Z/nZ)^* is characterised by group homomorphisms from each of these factor
-- groups. Furthermore, for odd p, we have (Z/p^k Z)^* isomorphic to Z / p^(k-1)*(p-1) Z, an
-- additive group, where an isomorphism is specified by a choice of primitive root.
-- Similarly, for k >= 2, (Z/2^k Z)^* is isomorphic to Z/2Z * (Z / 2^(k-2) Z) (and for k < 2
-- it is trivial).  (See @lambda@ for this isomorphism).
-- Thus, to specify a Dirichlet character, it suffices to specify the value of generators
-- of each of these cyclic groups, when primitive roots are given. This data is given by a
-- DirichletFactor.
-- We have the invariant that the factors must be given in strictly increasing order, and the
-- generator is as given by @generator@, and are each non-trivial. These conditions are verified
-- using `validChar`.
data DirichletFactor = OddPrime { DirichletFactor -> Prime Natural
_getPrime :: Prime Natural
                                , DirichletFactor -> Word
_getPower :: Word
                                , DirichletFactor -> Natural
_getGenerator :: Natural
                                , DirichletFactor -> RootOfUnity
_getValue :: RootOfUnity
                                }
                     | TwoPower { DirichletFactor -> Int
_getPower2 :: Int -- this ought to be Word, but many applications
                                                    -- needed to use wordToInt, so Int is cleaner
                                                    -- Required to be >= 2
                                , DirichletFactor -> RootOfUnity
_getFirstValue :: RootOfUnity
                                , DirichletFactor -> RootOfUnity
_getSecondValue :: RootOfUnity
                                }
                     | Two

instance Eq (DirichletCharacter n) where
  Generated [DirichletFactor]
a == :: DirichletCharacter n -> DirichletCharacter n -> Bool
== Generated [DirichletFactor]
b = [DirichletFactor]
a forall a. Eq a => a -> a -> Bool
== [DirichletFactor]
b

instance Eq DirichletFactor where
  TwoPower Int
_ RootOfUnity
x1 RootOfUnity
x2 == :: DirichletFactor -> DirichletFactor -> Bool
== TwoPower Int
_ RootOfUnity
y1 RootOfUnity
y2 = RootOfUnity
x1 forall a. Eq a => a -> a -> Bool
== RootOfUnity
y1 Bool -> Bool -> Bool
&& RootOfUnity
x2 forall a. Eq a => a -> a -> Bool
== RootOfUnity
y2
  OddPrime Prime Natural
_ Word
_ Natural
_ RootOfUnity
x == OddPrime Prime Natural
_ Word
_ Natural
_ RootOfUnity
y = RootOfUnity
x forall a. Eq a => a -> a -> Bool
== RootOfUnity
y
  DirichletFactor
Two              == DirichletFactor
Two              = Bool
True
  DirichletFactor
_ == DirichletFactor
_ = Bool
False

-- | For primes, define the canonical primitive root as the smallest such.
generator :: Prime Natural -> Word -> Natural
generator :: Prime Natural -> Word -> Natural
generator Prime Natural
p Word
k = case forall a.
(Eq a, Num a) =>
[(Prime a, Word)] -> Maybe (Some (CyclicGroup a))
cyclicGroupFromFactors [(Prime Natural
p, Word
k)] of
  Maybe (Some (CyclicGroup Natural))
Nothing -> forall a. HasCallStack => [Char] -> a
error [Char]
"illegal"
  Just (Some CyclicGroup Natural m
cg)
    | Sub Dict (KnownNat m)
(() :: Constraint) => Dict (KnownNat m)
Dict <- forall a (m :: Natural).
Integral a =>
CyclicGroup a m -> (() :: Constraint) :- KnownNat m
proofFromCyclicGroup CyclicGroup Natural m
cg ->
      forall (m :: Natural). Mod m -> Natural
unMod forall a b. (a -> b) -> a -> b
$ forall (m :: Natural). MultMod m -> Mod m
multElement forall a b. (a -> b) -> a -> b
$ forall (m :: Natural). PrimitiveRoot m -> MultMod m
unPrimitiveRoot forall a b. (a -> b) -> a -> b
$ forall a. [a] -> a
head forall a b. (a -> b) -> a -> b
$
        forall a b. (a -> Maybe b) -> [a] -> [b]
mapMaybe (forall a (m :: Natural).
(Integral a, UniqueFactorisation a) =>
CyclicGroup a m -> Mod m -> Maybe (PrimitiveRoot m)
isPrimitiveRoot CyclicGroup Natural m
cg) [Mod m
2..forall a. Bounded a => a
maxBound]

-- | Implement the function \(\lambda\) from page 5 of
-- https://www2.eecs.berkeley.edu/Pubs/TechRpts/1984/CSD-84-186.pdf
lambda :: Integer -> Int -> Integer
lambda :: Integer -> Int -> Integer
lambda Integer
x Int
e = ((Integer
xPower forall a. Num a => a -> a -> a
- Integer
1) forall a. Bits a => a -> Int -> a
`shiftR` (Int
eforall a. Num a => a -> a -> a
+Int
1)) forall a. Bits a => a -> a -> a
.&. (Integer
modulus forall a. Num a => a -> a -> a
- Integer
1)
  where
    modulus :: Integer
modulus  = Integer
1 forall a. Bits a => a -> Int -> a
`shiftL` (Int
e forall a. Num a => a -> a -> a
- Int
2)
    largeMod :: Natural
largeMod = Natural
1 forall a. Bits a => a -> Int -> a
`shiftL` (Int
2 forall a. Num a => a -> a -> a
* Int
e forall a. Num a => a -> a -> a
- Int
1)
    xPower :: Integer
xPower = case Natural -> SomeNat
someNatVal Natural
largeMod of
      SomeNat (Proxy n
_ :: Proxy largeMod) ->
        forall a. Integral a => a -> Integer
toInteger (forall (m :: Natural). Mod m -> Natural
unMod (forall a. Num a => Integer -> a
fromInteger Integer
x forall a b. (Num a, Integral b) => a -> b -> a
^ (Integer
2 forall a. Num a => a -> a -> a
* Integer
modulus) :: Mod largeMod))


-- | For elements of the multiplicative group \((\mathbb{Z}/n\mathbb{Z})^*\), a Dirichlet
-- character evaluates to a root of unity.
eval :: DirichletCharacter n -> MultMod n -> RootOfUnity
eval :: forall (n :: Natural).
DirichletCharacter n -> MultMod n -> RootOfUnity
eval (Generated [DirichletFactor]
ds) MultMod n
m = forall (t :: * -> *) m a.
(Foldable t, Monoid m) =>
(a -> m) -> t a -> m
foldMap (Integer -> DirichletFactor -> RootOfUnity
evalFactor Integer
m') [DirichletFactor]
ds
  where
    m' :: Integer
m' = forall a. Integral a => a -> Integer
toInteger forall a b. (a -> b) -> a -> b
$ forall (m :: Natural). Mod m -> Natural
unMod forall a b. (a -> b) -> a -> b
$ forall (m :: Natural). MultMod m -> Mod m
multElement MultMod n
m

-- | Evaluate each factor of the Dirichlet character.
evalFactor :: Integer -> DirichletFactor -> RootOfUnity
evalFactor :: Integer -> DirichletFactor -> RootOfUnity
evalFactor Integer
m =
  \case
    OddPrime (forall a. Integral a => a -> Integer
toInteger forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Prime a -> a
unPrime -> Integer
p) Word
k (forall a. Integral a => a -> Integer
toInteger -> Integer
a) RootOfUnity
b ->
      Integer -> Word -> Integer -> Integer -> Natural
discreteLogarithmPP Integer
p Word
k Integer
a (Integer
m forall a. Integral a => a -> a -> a
`rem` Integer
pforall a b. (Num a, Integral b) => a -> b -> a
^Word
k) forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b
    TwoPower Int
k RootOfUnity
s RootOfUnity
b -> (if forall a. Bits a => a -> Int -> Bool
testBit Integer
m Int
1 then RootOfUnity
s else forall a. Monoid a => a
mempty)
                   forall a. Semigroup a => a -> a -> a
<> Integer -> Int -> Integer
lambda (forall p. (Bits p, Num p) => Int -> p -> p
thingy Int
k Integer
m) Int
k forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b
    DirichletFactor
Two -> forall a. Monoid a => a
mempty

thingy :: (Bits p, Num p) => Int -> p -> p
thingy :: forall p. (Bits p, Num p) => Int -> p -> p
thingy Int
k p
m = if forall a. Bits a => a -> Int -> Bool
testBit p
m Int
1
                then forall a. Bits a => Int -> a
bit Int
k forall a. Num a => a -> a -> a
- p
m'
                else p
m'
  where m' :: p
m' = p
m forall a. Bits a => a -> a -> a
.&. (forall a. Bits a => Int -> a
bit Int
k forall a. Num a => a -> a -> a
- p
1)

-- | A character can evaluate to a root of unity or zero: represented by @Nothing@.
evalGeneral :: KnownNat n => DirichletCharacter n -> Mod n -> OrZero RootOfUnity
evalGeneral :: forall (n :: Natural).
KnownNat n =>
DirichletCharacter n -> Mod n -> OrZero RootOfUnity
evalGeneral DirichletCharacter n
chi Mod n
t = case forall (m :: Natural). KnownNat m => Mod m -> Maybe (MultMod m)
isMultElement Mod n
t of
                      Maybe (MultMod n)
Nothing -> forall a. OrZero a
Zero
                      Just MultMod n
x -> forall a. a -> OrZero a
NonZero forall a b. (a -> b) -> a -> b
$ forall (n :: Natural).
DirichletCharacter n -> MultMod n -> RootOfUnity
eval DirichletCharacter n
chi MultMod n
x

-- | Give the principal character for this modulus: a principal character mod \(n\) is 1 for
-- \(a\) coprime to \(n\), and 0 otherwise.
principalChar :: KnownNat n => DirichletCharacter n
principalChar :: forall (n :: Natural). KnownNat n => DirichletCharacter n
principalChar = forall a. Bounded a => a
minBound

mulChars :: DirichletCharacter n -> DirichletCharacter n -> DirichletCharacter n
mulChars :: forall (n :: Natural).
DirichletCharacter n
-> DirichletCharacter n -> DirichletCharacter n
mulChars (Generated [DirichletFactor]
x) (Generated [DirichletFactor]
y) = forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated (forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith DirichletFactor -> DirichletFactor -> DirichletFactor
combine [DirichletFactor]
x [DirichletFactor]
y)
  where combine :: DirichletFactor -> DirichletFactor -> DirichletFactor
        combine :: DirichletFactor -> DirichletFactor -> DirichletFactor
combine DirichletFactor
Two DirichletFactor
Two = DirichletFactor
Two
        combine (OddPrime Prime Natural
p Word
k Natural
g RootOfUnity
n) (OddPrime Prime Natural
_ Word
_ Natural
_ RootOfUnity
m) =
          Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g (RootOfUnity
n forall a. Semigroup a => a -> a -> a
<> RootOfUnity
m)
        combine (TwoPower Int
k RootOfUnity
a RootOfUnity
n) (TwoPower Int
_ RootOfUnity
b RootOfUnity
m) =
          Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k (RootOfUnity
a forall a. Semigroup a => a -> a -> a
<> RootOfUnity
b) (RootOfUnity
n forall a. Semigroup a => a -> a -> a
<> RootOfUnity
m)
        combine DirichletFactor
_ DirichletFactor
_ = forall a. HasCallStack => [Char] -> a
error [Char]
"internal error: malformed DirichletCharacter"

-- | This Semigroup is in fact a group, so @stimes@ can be called with a negative first argument.
instance Semigroup (DirichletCharacter n) where
  <> :: DirichletCharacter n
-> DirichletCharacter n -> DirichletCharacter n
(<>) = forall (n :: Natural).
DirichletCharacter n
-> DirichletCharacter n -> DirichletCharacter n
mulChars
  stimes :: forall b.
Integral b =>
b -> DirichletCharacter n -> DirichletCharacter n
stimes = forall a (n :: Natural).
Integral a =>
a -> DirichletCharacter n -> DirichletCharacter n
stimesChar

instance KnownNat n => Monoid (DirichletCharacter n) where
  mempty :: DirichletCharacter n
mempty = forall (n :: Natural). KnownNat n => DirichletCharacter n
principalChar
  mappend :: DirichletCharacter n
-> DirichletCharacter n -> DirichletCharacter n
mappend = forall a. Semigroup a => a -> a -> a
(<>)

stimesChar :: Integral a => a -> DirichletCharacter n -> DirichletCharacter n
stimesChar :: forall a (n :: Natural).
Integral a =>
a -> DirichletCharacter n -> DirichletCharacter n
stimesChar a
s (Generated [DirichletFactor]
xs) = forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated (forall a b. (a -> b) -> [a] -> [b]
map DirichletFactor -> DirichletFactor
mult [DirichletFactor]
xs)
  where mult :: DirichletFactor -> DirichletFactor
        mult :: DirichletFactor -> DirichletFactor
mult (OddPrime Prime Natural
p Word
k Natural
g RootOfUnity
n) = Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g (a
s forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
n)
        mult (TwoPower Int
k RootOfUnity
a RootOfUnity
b) = Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k (a
s forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
a) (a
s forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b)
        mult DirichletFactor
Two = DirichletFactor
Two

-- | We define `succ` and `pred` with more efficient implementations than
-- @`toEnum` . (+1) . `fromEnum`@.
instance KnownNat n => Enum (DirichletCharacter n) where
  toEnum :: Int -> DirichletCharacter n
toEnum = forall (n :: Natural).
KnownNat n =>
Natural -> DirichletCharacter n
indexToChar forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Natural
intToNatural
  fromEnum :: DirichletCharacter n -> Int
fromEnum = Integer -> Int
integerToInt forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (n :: Natural). DirichletCharacter n -> Integer
characterNumber
  succ :: DirichletCharacter n -> DirichletCharacter n
succ DirichletCharacter n
x = forall a (n :: Natural).
Integral a =>
DirichletCharacter n -> a -> DirichletCharacter n
makeChar DirichletCharacter n
x (forall (n :: Natural). DirichletCharacter n -> Integer
characterNumber DirichletCharacter n
x forall a. Num a => a -> a -> a
+ Integer
1)
  pred :: DirichletCharacter n -> DirichletCharacter n
pred DirichletCharacter n
x = forall a (n :: Natural).
Integral a =>
DirichletCharacter n -> a -> DirichletCharacter n
makeChar DirichletCharacter n
x (forall (n :: Natural). DirichletCharacter n -> Integer
characterNumber DirichletCharacter n
x forall a. Num a => a -> a -> a
- Integer
1)

  enumFromTo :: DirichletCharacter n
-> DirichletCharacter n -> [DirichletCharacter n]
enumFromTo DirichletCharacter n
x DirichletCharacter n
y       = forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x [forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
x..forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
y]
  enumFrom :: DirichletCharacter n -> [DirichletCharacter n]
enumFrom DirichletCharacter n
x           = forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x [forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
x..]
  enumFromThenTo :: DirichletCharacter n
-> DirichletCharacter n
-> DirichletCharacter n
-> [DirichletCharacter n]
enumFromThenTo DirichletCharacter n
x DirichletCharacter n
y DirichletCharacter n
z = forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x [forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
x, forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
y..forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
z]
  enumFromThen :: DirichletCharacter n
-> DirichletCharacter n -> [DirichletCharacter n]
enumFromThen DirichletCharacter n
x DirichletCharacter n
y     = forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x [forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
x, forall a. Enum a => a -> Int
fromEnum DirichletCharacter n
y..]

instance KnownNat n => Bounded (DirichletCharacter n) where
  minBound :: DirichletCharacter n
minBound = forall (n :: Natural).
KnownNat n =>
Natural -> DirichletCharacter n
indexToChar Natural
0
  maxBound :: DirichletCharacter n
maxBound = forall (n :: Natural).
KnownNat n =>
Natural -> DirichletCharacter n
indexToChar (forall n. UniqueFactorisation n => n -> n
totient Natural
n forall a. Num a => a -> a -> a
- Natural
1)
    where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)

-- | We have a (non-canonical) enumeration of dirichlet characters.
characterNumber :: DirichletCharacter n -> Integer
characterNumber :: forall (n :: Natural). DirichletCharacter n -> Integer
characterNumber (Generated [DirichletFactor]
y) = forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' Integer -> DirichletFactor -> Integer
go Integer
0 [DirichletFactor]
y
  where go :: Integer -> DirichletFactor -> Integer
go Integer
x (OddPrime Prime Natural
p Word
k Natural
_ RootOfUnity
a) = Integer
x forall a. Num a => a -> a -> a
* Integer
m forall a. Num a => a -> a -> a
+ forall a. Ratio a -> a
numerator (RootOfUnity -> Ratio Integer
fromRootOfUnity RootOfUnity
a forall a. Num a => a -> a -> a
* (Integer
m forall a. Integral a => a -> a -> Ratio a
% Integer
1))
          where p' :: Integer
p' = Natural -> Integer
naturalToInteger (forall a. Prime a -> a
unPrime Prime Natural
p)
                m :: Integer
m = Integer
p'forall a b. (Num a, Integral b) => a -> b -> a
^(Word
kforall a. Num a => a -> a -> a
-Word
1)forall a. Num a => a -> a -> a
*(Integer
p'forall a. Num a => a -> a -> a
-Integer
1)
        go Integer
x (TwoPower Int
k RootOfUnity
a RootOfUnity
b)   = Integer
x' forall a. Num a => a -> a -> a
* Integer
2 forall a. Num a => a -> a -> a
+ forall a. Ratio a -> a
numerator (RootOfUnity -> Ratio Integer
fromRootOfUnity RootOfUnity
a forall a. Num a => a -> a -> a
* Ratio Integer
2)
          where m :: Integer
m = forall a. Bits a => Int -> a
bit (Int
kforall a. Num a => a -> a -> a
-Int
2) :: Integer
                x' :: Integer
x' = Integer
x forall a. Bits a => a -> Int -> a
`shiftL` (Int
kforall a. Num a => a -> a -> a
-Int
2) forall a. Num a => a -> a -> a
+ forall a. Ratio a -> a
numerator (RootOfUnity -> Ratio Integer
fromRootOfUnity RootOfUnity
b forall a. Num a => a -> a -> a
* (Integer
m forall a. Integral a => a -> a -> Ratio a
% Integer
1))
        go Integer
x DirichletFactor
Two = Integer
x

-- | Give the dirichlet character from its number.
-- Inverse of `characterNumber`.
indexToChar :: forall n. KnownNat n => Natural -> DirichletCharacter n
indexToChar :: forall (n :: Natural).
KnownNat n =>
Natural -> DirichletCharacter n
indexToChar = forall a. Identity a -> a
runIdentity forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (n :: Natural) (f :: * -> *).
(KnownNat n, Functor f) =>
f Natural -> f (DirichletCharacter n)
indicesToChars forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. a -> Identity a
Identity

-- | Give a collection of dirichlet characters from their numbers. This may be more efficient than
-- `indexToChar` for multiple characters, as it prevents some internal recalculations.
indicesToChars :: forall n f. (KnownNat n, Functor f) => f Natural -> f (DirichletCharacter n)
indicesToChars :: forall (n :: Natural) (f :: * -> *).
(KnownNat n, Functor f) =>
f Natural -> f (DirichletCharacter n)
indicesToChars = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Template] -> Natural -> [DirichletFactor]
unroll [Template]
t forall b c a. (b -> c) -> (a -> b) -> a -> c
. (forall a. Integral a => a -> a -> a
`mod` Natural
m))
  where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        (Product Natural
m, [Template]
t) = Natural -> (Product Natural, [Template])
mkTemplate Natural
n

-- | List all characters for the modulus. This is preferred to using @[minBound..maxBound]@.
allChars :: forall n. KnownNat n => [DirichletCharacter n]
allChars :: forall (n :: Natural). KnownNat n => [DirichletCharacter n]
allChars = forall (n :: Natural) (f :: * -> *).
(KnownNat n, Functor f) =>
f Natural -> f (DirichletCharacter n)
indicesToChars [Natural
0..Natural
mforall a. Num a => a -> a -> a
-Natural
1]
  where m :: Natural
m = forall n. UniqueFactorisation n => n -> n
totient forall a b. (a -> b) -> a -> b
$ forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)

-- | The same as `indexToChar`, but if we're given a character we can create others more efficiently.
makeChar :: Integral a => DirichletCharacter n -> a -> DirichletCharacter n
makeChar :: forall a (n :: Natural).
Integral a =>
DirichletCharacter n -> a -> DirichletCharacter n
makeChar DirichletCharacter n
x = forall a. Identity a -> a
runIdentity forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. a -> Identity a
Identity

-- | Use one character to make many more: better than indicesToChars since it avoids recalculating
-- some primitive roots
bulkMakeChars :: (Integral a, Functor f) => DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars :: forall a (f :: * -> *) (n :: Natural).
(Integral a, Functor f) =>
DirichletCharacter n -> f a -> f (DirichletCharacter n)
bulkMakeChars DirichletCharacter n
x = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Template] -> Natural -> [DirichletFactor]
unroll [Template]
t forall b c a. (b -> c) -> (a -> b) -> a -> c
. (forall a. Integral a => a -> a -> a
`mod` Natural
m) forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a b. (Integral a, Num b) => a -> b
fromIntegral')
  where (Product Natural
m, [Template]
t) = forall (n :: Natural).
DirichletCharacter n -> (Product Natural, [Template])
templateFromCharacter DirichletCharacter n
x

-- We assign each natural a unique Template, which can be decorated (eg in `unroll`) to
-- form a DirichletCharacter. A Template effectively holds the information carried around
-- in a DirichletFactor which depends only on the modulus of the character.
data Template = OddTemplate { Template -> Prime Natural
_getPrime'     :: Prime Natural
                            , Template -> Word
_getPower'     :: Word
                            , Template -> Natural
_getGenerator' :: !Natural
                            , Template -> Natural
_getModulus'   :: !Natural
                            }
              | TwoPTemplate { Template -> Int
_getPower2'    :: Int
                             , _getModulus'   :: !Natural
                             } -- the modulus is derivable from the other values, but calculation
                               -- may be expensive, so we pre-calculate it
                               -- morally getModulus should be a prefactored but seems to be
                               -- pointless here
              | TwoTemplate

templateFromCharacter :: DirichletCharacter n -> (Product Natural, [Template])
templateFromCharacter :: forall (n :: Natural).
DirichletCharacter n -> (Product Natural, [Template])
templateFromCharacter (Generated [DirichletFactor]
t) = forall (t :: * -> *) (f :: * -> *) a b.
(Traversable t, Applicative f) =>
(a -> f b) -> t a -> f (t b)
traverse DirichletFactor -> (Product Natural, Template)
go [DirichletFactor]
t
  where go :: DirichletFactor -> (Product Natural, Template)
go (OddPrime Prime Natural
p Word
k Natural
g RootOfUnity
_) = (forall a. a -> Product a
Product Natural
m, Prime Natural -> Word -> Natural -> Natural -> Template
OddTemplate Prime Natural
p Word
k Natural
g Natural
m)
          where p' :: Natural
p' = forall a. Prime a -> a
unPrime Prime Natural
p
                m :: Natural
m = Natural
p'forall a b. (Num a, Integral b) => a -> b -> a
^(Word
kforall a. Num a => a -> a -> a
-Word
1)forall a. Num a => a -> a -> a
*(Natural
p'forall a. Num a => a -> a -> a
-Natural
1)
        go (TwoPower Int
k RootOfUnity
_ RootOfUnity
_) = (forall a. a -> Product a
Product (Natural
2forall a. Num a => a -> a -> a
*Natural
m), Int -> Natural -> Template
TwoPTemplate Int
k Natural
m)
          where m :: Natural
m = forall a. Bits a => Int -> a
bit (Int
kforall a. Num a => a -> a -> a
-Int
2)
        go DirichletFactor
Two = (forall a. a -> Product a
Product Natural
1, Template
TwoTemplate)

mkTemplate :: Natural -> (Product Natural, [Template])
mkTemplate :: Natural -> (Product Natural, [Template])
mkTemplate = [(Prime Natural, Word)] -> (Product Natural, [Template])
go forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Ord a => [a] -> [a]
sort forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. UniqueFactorisation a => a -> [(Prime a, Word)]
factorise
  where go :: [(Prime Natural, Word)] -> (Product Natural, [Template])
        go :: [(Prime Natural, Word)] -> (Product Natural, [Template])
go ((forall a. Prime a -> a
unPrime -> Natural
2, Word
1): [(Prime Natural, Word)]
xs) = (forall a. a -> Product a
Product Natural
1, [Template
TwoTemplate]) forall a. Semigroup a => a -> a -> a
<> forall (t :: * -> *) (f :: * -> *) a b.
(Traversable t, Applicative f) =>
(a -> f b) -> t a -> f (t b)
traverse (Prime Natural, Word) -> (Product Natural, Template)
odds [(Prime Natural, Word)]
xs
        go ((forall a. Prime a -> a
unPrime -> Natural
2, Word -> Int
wordToInt -> Int
k): [(Prime Natural, Word)]
xs) = (forall a. a -> Product a
Product (Natural
2forall a. Num a => a -> a -> a
*Natural
m), [Int -> Natural -> Template
TwoPTemplate Int
k Natural
m]) forall a. Semigroup a => a -> a -> a
<> forall (t :: * -> *) (f :: * -> *) a b.
(Traversable t, Applicative f) =>
(a -> f b) -> t a -> f (t b)
traverse (Prime Natural, Word) -> (Product Natural, Template)
odds [(Prime Natural, Word)]
xs
          where m :: Natural
m = forall a. Bits a => Int -> a
bit (Int
kforall a. Num a => a -> a -> a
-Int
2)
        go [(Prime Natural, Word)]
xs = forall (t :: * -> *) (f :: * -> *) a b.
(Traversable t, Applicative f) =>
(a -> f b) -> t a -> f (t b)
traverse (Prime Natural, Word) -> (Product Natural, Template)
odds [(Prime Natural, Word)]
xs
        odds :: (Prime Natural, Word) -> (Product Natural, Template)
        odds :: (Prime Natural, Word) -> (Product Natural, Template)
odds (Prime Natural
p, Word
k) = (forall a. a -> Product a
Product Natural
m, Prime Natural -> Word -> Natural -> Natural -> Template
OddTemplate Prime Natural
p Word
k (Prime Natural -> Word -> Natural
generator Prime Natural
p Word
k) Natural
m)
          where p' :: Natural
p' = forall a. Prime a -> a
unPrime Prime Natural
p
                m :: Natural
m = Natural
p'forall a b. (Num a, Integral b) => a -> b -> a
^(Word
kforall a. Num a => a -> a -> a
-Word
1)forall a. Num a => a -> a -> a
*(Natural
p'forall a. Num a => a -> a -> a
-Natural
1)

-- the validity of the producted dirichletfactor list here requires the template to be valid
unroll :: [Template] -> Natural -> [DirichletFactor]
unroll :: [Template] -> Natural -> [DirichletFactor]
unroll [Template]
t Natural
m = forall a b. (a, b) -> b
snd (forall (t :: * -> *) s a b.
Traversable t =>
(s -> a -> (s, b)) -> s -> t a -> (s, t b)
mapAccumL Natural -> Template -> (Natural, DirichletFactor)
func Natural
m [Template]
t)
  where func :: Natural -> Template -> (Natural, DirichletFactor)
        func :: Natural -> Template -> (Natural, DirichletFactor)
func Natural
a (OddTemplate Prime Natural
p Word
k Natural
g Natural
n) = (Natural
a1, Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g (Ratio Integer -> RootOfUnity
toRootOfUnity forall a b. (a -> b) -> a -> b
$ forall a. Integral a => a -> Integer
toInteger Natural
a2 forall a. Integral a => a -> a -> Ratio a
% forall a. Integral a => a -> Integer
toInteger Natural
n))
          where (Natural
a1,Natural
a2) = forall a. Integral a => a -> a -> (a, a)
quotRem Natural
a Natural
n
        func Natural
a (TwoPTemplate Int
k Natural
n) = (Natural
b1, Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k (Ratio Integer -> RootOfUnity
toRootOfUnity forall a b. (a -> b) -> a -> b
$ forall a. Integral a => a -> Integer
toInteger Natural
a2 forall a. Integral a => a -> a -> Ratio a
% Integer
2) (Ratio Integer -> RootOfUnity
toRootOfUnity forall a b. (a -> b) -> a -> b
$ forall a. Integral a => a -> Integer
toInteger Natural
b2 forall a. Integral a => a -> a -> Ratio a
% forall a. Integral a => a -> Integer
toInteger Natural
n))
          where (Natural
a1,Natural
a2) = forall a. Integral a => a -> a -> (a, a)
quotRem Natural
a Natural
2
                (Natural
b1,Natural
b2) = forall a. Integral a => a -> a -> (a, a)
quotRem Natural
a1 Natural
n
        func Natural
a Template
TwoTemplate = (Natural
a, DirichletFactor
Two)

-- | Test if a given Dirichlet character is prinicpal for its modulus: a principal character mod
-- \(n\) is 1 for \(a\) coprime to \(n\), and 0 otherwise.
isPrincipal :: DirichletCharacter n -> Bool
isPrincipal :: forall (n :: Natural). DirichletCharacter n -> Bool
isPrincipal DirichletCharacter n
chi = forall (n :: Natural). DirichletCharacter n -> Integer
characterNumber DirichletCharacter n
chi forall a. Eq a => a -> a -> Bool
== Integer
0

-- | Induce a Dirichlet character to a higher modulus. If \(d \mid n\), then \(a \bmod{n}\) can be
-- reduced to \(a \bmod{d}\). Thus, the multiplicative function on \(\mathbb{Z}/d\mathbb{Z}\)
-- induces a multiplicative function on \(\mathbb{Z}/n\mathbb{Z}\).
--
-- >>> :set -XTypeApplications -XDataKinds
-- >>> chi = indexToChar 5 :: DirichletCharacter 45
-- >>> chi2 = induced @135 chi :: Maybe (DirichletCharacter 135)
induced :: forall n d. (KnownNat d, KnownNat n) => DirichletCharacter d -> Maybe (DirichletCharacter n)
induced :: forall (n :: Natural) (d :: Natural).
(KnownNat d, KnownNat n) =>
DirichletCharacter d -> Maybe (DirichletCharacter n)
induced (Generated [DirichletFactor]
start) = if Natural
n forall a. Integral a => a -> a -> a
`rem` Natural
d forall a. Eq a => a -> a -> Bool
== Natural
0
                            then forall a. a -> Maybe a
Just (forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated ([Template] -> [DirichletFactor] -> [DirichletFactor]
combine (forall a b. (a, b) -> b
snd forall a b. (a -> b) -> a -> b
$ Natural -> (Product Natural, [Template])
mkTemplate Natural
n) [DirichletFactor]
start))
                            else forall a. Maybe a
Nothing
  where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        d :: Natural
d = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy d)
        combine :: [Template] -> [DirichletFactor] -> [DirichletFactor]
        combine :: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [] [DirichletFactor]
_ = []
        combine [Template]
ts [] = forall a b. (a -> b) -> [a] -> [b]
map Template -> DirichletFactor
newFactor [Template]
ts
        combine (Template
t:[Template]
xs) (DirichletFactor
y:[DirichletFactor]
ys) = case (Template
t,DirichletFactor
y) of
                                  (Template
TwoTemplate, DirichletFactor
Two) -> DirichletFactor
Twoforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs [DirichletFactor]
ys
                                  (Template
TwoTemplate, DirichletFactor
_) -> DirichletFactor
Twoforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs (DirichletFactor
yforall a. a -> [a] -> [a]
:[DirichletFactor]
ys)
                                  (TwoPTemplate Int
k Natural
_, DirichletFactor
Two) -> Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k forall a. Monoid a => a
mempty forall a. Monoid a => a
memptyforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs [DirichletFactor]
ys
                                  (TwoPTemplate Int
k Natural
_, TwoPower Int
_ RootOfUnity
a RootOfUnity
b) -> Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k RootOfUnity
a RootOfUnity
bforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs [DirichletFactor]
ys
                                  (TwoPTemplate Int
k Natural
_, DirichletFactor
_) -> Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k forall a. Monoid a => a
mempty forall a. Monoid a => a
memptyforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs (DirichletFactor
yforall a. a -> [a] -> [a]
:[DirichletFactor]
ys)
                                  (OddTemplate Prime Natural
p Word
k Natural
_ Natural
_, OddPrime Prime Natural
q Word
_ Natural
g RootOfUnity
a) | Prime Natural
p forall a. Eq a => a -> a -> Bool
== Prime Natural
q -> Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g RootOfUnity
aforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs [DirichletFactor]
ys
                                  (OddTemplate Prime Natural
p Word
k Natural
g Natural
_, OddPrime Prime Natural
q Word
_ Natural
_ RootOfUnity
_) | Prime Natural
p forall a. Ord a => a -> a -> Bool
< Prime Natural
q -> Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g forall a. Monoid a => a
memptyforall a. a -> [a] -> [a]
: [Template] -> [DirichletFactor] -> [DirichletFactor]
combine [Template]
xs (DirichletFactor
yforall a. a -> [a] -> [a]
:[DirichletFactor]
ys)
                                  (Template, DirichletFactor)
_ -> forall a. HasCallStack => [Char] -> a
error [Char]
"internal error in induced: please report this as a bug"
        newFactor :: Template -> DirichletFactor
        newFactor :: Template -> DirichletFactor
newFactor Template
TwoTemplate = DirichletFactor
Two
        newFactor (TwoPTemplate Int
k Natural
_) = Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k forall a. Monoid a => a
mempty forall a. Monoid a => a
mempty
        newFactor (OddTemplate Prime Natural
p Word
k Natural
g Natural
_) = Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g forall a. Monoid a => a
mempty

-- | The <https://en.wikipedia.org/wiki/Jacobi_symbol Jacobi symbol> gives a real Dirichlet
-- character for odd moduli.
jacobiCharacter :: forall n. KnownNat n => Maybe (RealCharacter n)
jacobiCharacter :: forall (n :: Natural). KnownNat n => Maybe (RealCharacter n)
jacobiCharacter = if forall a. Integral a => a -> Bool
odd Natural
n
                     then forall a. a -> Maybe a
Just forall a b. (a -> b) -> a -> b
$ forall (n :: Natural). DirichletCharacter n -> RealCharacter n
RealChar forall a b. (a -> b) -> a -> b
$ forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map Template -> DirichletFactor
go forall a b. (a -> b) -> a -> b
$ forall a b. (a, b) -> b
snd forall a b. (a -> b) -> a -> b
$ Natural -> (Product Natural, [Template])
mkTemplate Natural
n
                     else forall a. Maybe a
Nothing
  where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        go :: Template -> DirichletFactor
        go :: Template -> DirichletFactor
go (OddTemplate Prime Natural
p Word
k Natural
g Natural
_) = Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g forall a b. (a -> b) -> a -> b
$ Ratio Integer -> RootOfUnity
toRootOfUnity (forall a. Integral a => a -> Integer
toInteger Word
k forall a. Integral a => a -> a -> Ratio a
% Integer
2)
          -- jacobi symbol of a primitive root mod p over p is always -1
        go Template
_ = forall a. HasCallStack => [Char] -> a
error [Char]
"internal error in jacobiCharacter: please report this as a bug"
          -- every factor of n should be odd

-- | A Dirichlet character is real if it is real-valued.
newtype RealCharacter n = RealChar { -- | Extract the character itself from a `RealCharacter`.
                                     forall (n :: Natural). RealCharacter n -> DirichletCharacter n
getRealChar :: DirichletCharacter n
                                   }
                                   deriving RealCharacter n -> RealCharacter n -> Bool
forall (n :: Natural). RealCharacter n -> RealCharacter n -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: RealCharacter n -> RealCharacter n -> Bool
$c/= :: forall (n :: Natural). RealCharacter n -> RealCharacter n -> Bool
== :: RealCharacter n -> RealCharacter n -> Bool
$c== :: forall (n :: Natural). RealCharacter n -> RealCharacter n -> Bool
Eq

-- | Test if a given `DirichletCharacter` is real, and if so give a `RealCharacter`.
isRealCharacter :: DirichletCharacter n -> Maybe (RealCharacter n)
isRealCharacter :: forall (n :: Natural).
DirichletCharacter n -> Maybe (RealCharacter n)
isRealCharacter t :: DirichletCharacter n
t@(Generated [DirichletFactor]
xs) = if forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all DirichletFactor -> Bool
real [DirichletFactor]
xs then forall a. a -> Maybe a
Just (forall (n :: Natural). DirichletCharacter n -> RealCharacter n
RealChar DirichletCharacter n
t) else forall a. Maybe a
Nothing
  where real :: DirichletFactor -> Bool
        real :: DirichletFactor -> Bool
real (OddPrime Prime Natural
_ Word
_ Natural
_ RootOfUnity
a) = RootOfUnity
a forall a. Semigroup a => a -> a -> a
<> RootOfUnity
a forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty
        real (TwoPower Int
_ RootOfUnity
_ RootOfUnity
b) = RootOfUnity
b forall a. Semigroup a => a -> a -> a
<> RootOfUnity
b forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty
        real DirichletFactor
Two = Bool
True

-- TODO: it should be possible to calculate this without eval/evalGeneral
-- and thus avoid using discrete log calculations: consider the order of m
-- inside each of the factor groups?
-- | Evaluate a real Dirichlet character, which can only take values \(-1,0,1\).
toRealFunction :: KnownNat n => RealCharacter n -> Mod n -> Int
toRealFunction :: forall (n :: Natural).
KnownNat n =>
RealCharacter n -> Mod n -> Int
toRealFunction (RealChar DirichletCharacter n
chi) Mod n
m = case forall (n :: Natural).
KnownNat n =>
DirichletCharacter n -> Mod n -> OrZero RootOfUnity
evalGeneral DirichletCharacter n
chi Mod n
m of
                                    OrZero RootOfUnity
Zero -> Int
0
                                    NonZero RootOfUnity
t | RootOfUnity
t forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty -> Int
1
                                    NonZero RootOfUnity
t | RootOfUnity
t forall a. Eq a => a -> a -> Bool
== Ratio Integer -> RootOfUnity
RootOfUnity (Integer
1 forall a. Integral a => a -> a -> Ratio a
% Integer
2) -> -Int
1
                                    OrZero RootOfUnity
_ -> forall a. HasCallStack => [Char] -> a
error [Char]
"internal error in toRealFunction: please report this as a bug"
                                      -- A real character should not be able to evaluate to
                                      -- anything other than {-1,0,1}, so should not reach this branch

-- | Test if the internal DirichletCharacter structure is valid.
validChar :: forall n. KnownNat n => DirichletCharacter n -> Bool
validChar :: forall (n :: Natural). KnownNat n => DirichletCharacter n -> Bool
validChar (Generated [DirichletFactor]
xs) = Bool
correctDecomposition Bool -> Bool -> Bool
&& forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all DirichletFactor -> Bool
correctPrimitiveRoot [DirichletFactor]
xs Bool -> Bool -> Bool
&& forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all DirichletFactor -> Bool
validValued [DirichletFactor]
xs
  where correctDecomposition :: Bool
correctDecomposition = forall a. Ord a => [a] -> [a]
sort (forall a. UniqueFactorisation a => a -> [(Prime a, Word)]
factorise Natural
n) forall a. Eq a => a -> a -> Bool
== forall a b. (a -> b) -> [a] -> [b]
map DirichletFactor -> (Prime Natural, Word)
getPP [DirichletFactor]
xs
        getPP :: DirichletFactor -> (Prime Natural, Word)
getPP (TwoPower Int
k RootOfUnity
_ RootOfUnity
_) = (Prime Natural
two, Int -> Word
intToWord Int
k)
        getPP (OddPrime Prime Natural
p Word
k Natural
_ RootOfUnity
_) = (Prime Natural
p, Word
k)
        getPP DirichletFactor
Two = (Prime Natural
two,Word
1)
        correctPrimitiveRoot :: DirichletFactor -> Bool
correctPrimitiveRoot (OddPrime Prime Natural
p Word
k Natural
g RootOfUnity
_) = Natural
g forall a. Eq a => a -> a -> Bool
== Prime Natural -> Word -> Natural
generator Prime Natural
p Word
k
        correctPrimitiveRoot DirichletFactor
_ = Bool
True
        validValued :: DirichletFactor -> Bool
validValued (TwoPower Int
k RootOfUnity
a RootOfUnity
b) = RootOfUnity
a forall a. Semigroup a => a -> a -> a
<> RootOfUnity
a forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty Bool -> Bool -> Bool
&& (forall a. Bits a => Int -> a
bit (Int
kforall a. Num a => a -> a -> a
-Int
2) :: Integer) forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty
        validValued (OddPrime (forall a. Prime a -> a
unPrime -> Natural
p) Word
k Natural
_ RootOfUnity
a) = (Natural
pforall a b. (Num a, Integral b) => a -> b -> a
^(Word
kforall a. Num a => a -> a -> a
-Word
1)forall a. Num a => a -> a -> a
*(Natural
pforall a. Num a => a -> a -> a
-Natural
1)) forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
a forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty
        validValued DirichletFactor
Two = Bool
True
        n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        two :: Prime Natural
two = forall a.
(Bits a, Integral a, UniqueFactorisation a) =>
a -> Prime a
nextPrime Natural
2

-- | Get the order of the Dirichlet Character.
orderChar :: DirichletCharacter n -> Integer
orderChar :: forall (n :: Natural). DirichletCharacter n -> Integer
orderChar (Generated [DirichletFactor]
xs) = forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' forall a. Integral a => a -> a -> a
lcm Integer
1 forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map DirichletFactor -> Integer
orderFactor [DirichletFactor]
xs
  where orderFactor :: DirichletFactor -> Integer
orderFactor (TwoPower Int
_ (RootOfUnity Ratio Integer
a) (RootOfUnity Ratio Integer
b)) = forall a. Ratio a -> a
denominator Ratio Integer
a forall a. Integral a => a -> a -> a
`lcm` forall a. Ratio a -> a
denominator Ratio Integer
b
        orderFactor (OddPrime Prime Natural
_ Word
_ Natural
_ (RootOfUnity Ratio Integer
a)) = forall a. Ratio a -> a
denominator Ratio Integer
a
        orderFactor DirichletFactor
Two = Integer
1

-- | Test if a Dirichlet character is <https://en.wikipedia.org/wiki/Dirichlet_character#Primitive_characters_and_conductor primitive>.
isPrimitive :: DirichletCharacter n -> Maybe (PrimitiveCharacter n)
isPrimitive :: forall (n :: Natural).
DirichletCharacter n -> Maybe (PrimitiveCharacter n)
isPrimitive t :: DirichletCharacter n
t@(Generated [DirichletFactor]
xs) = if forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all DirichletFactor -> Bool
primitive [DirichletFactor]
xs then forall a. a -> Maybe a
Just (forall (n :: Natural). DirichletCharacter n -> PrimitiveCharacter n
PrimitiveCharacter DirichletCharacter n
t) else forall a. Maybe a
Nothing
  where primitive :: DirichletFactor -> Bool
        primitive :: DirichletFactor -> Bool
primitive DirichletFactor
Two = Bool
False
        -- for odd p, we're testing if phi(p^(k-1)) `stimes` a is 1, since this means the
        -- character can come from some the smaller modulus p^(k-1)
        primitive (OddPrime Prime Natural
_ Word
1 Natural
_ RootOfUnity
a) = RootOfUnity
a forall a. Eq a => a -> a -> Bool
/= forall a. Monoid a => a
mempty
        primitive (OddPrime (forall a. Prime a -> a
unPrime -> Natural
p) Word
k Natural
_ RootOfUnity
a) = (Natural
pforall a b. (Num a, Integral b) => a -> b -> a
^(Word
kforall a. Num a => a -> a -> a
-Word
2)forall a. Num a => a -> a -> a
*(Natural
pforall a. Num a => a -> a -> a
-Natural
1)) forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
a forall a. Eq a => a -> a -> Bool
/= forall a. Monoid a => a
mempty
        primitive (TwoPower Int
2 RootOfUnity
a RootOfUnity
_) = RootOfUnity
a forall a. Eq a => a -> a -> Bool
/= forall a. Monoid a => a
mempty
        primitive (TwoPower Int
k RootOfUnity
_ RootOfUnity
b) = (forall a. Bits a => Int -> a
bit (Int
kforall a. Num a => a -> a -> a
-Int
3) :: Integer) forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b forall a. Eq a => a -> a -> Bool
/= forall a. Monoid a => a
mempty

-- | A Dirichlet character is primitive if cannot be 'induced' from any character with
-- strictly smaller modulus.
newtype PrimitiveCharacter n = PrimitiveCharacter { -- | Extract the character itself from a `PrimitiveCharacter`.
                                                    forall (n :: Natural). PrimitiveCharacter n -> DirichletCharacter n
getPrimitiveChar :: DirichletCharacter n
                                                    }
                                                    deriving PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
forall (n :: Natural).
PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
$c/= :: forall (n :: Natural).
PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
== :: PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
$c== :: forall (n :: Natural).
PrimitiveCharacter n -> PrimitiveCharacter n -> Bool
Eq

-- | Wrapper to hide an unknown type-level natural.
data WithNat (a :: Nat -> Type) where
  WithNat :: KnownNat m => a m -> WithNat a

-- | This function also provides access to the new modulus on type level, with a KnownNat instance
makePrimitive :: DirichletCharacter n -> WithNat PrimitiveCharacter
makePrimitive :: forall (n :: Natural).
DirichletCharacter n -> WithNat PrimitiveCharacter
makePrimitive (Generated [DirichletFactor]
xs) =
  case Natural -> SomeNat
someNatVal (forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Natural]
mods) of
    SomeNat (Proxy n
Proxy :: Proxy m) -> forall (m :: Natural) (a :: Natural -> *).
KnownNat m =>
a m -> WithNat a
WithNat (forall (n :: Natural). DirichletCharacter n -> PrimitiveCharacter n
PrimitiveCharacter (forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated [DirichletFactor]
ys) :: PrimitiveCharacter m)
  where ([Natural]
mods,[DirichletFactor]
ys) = forall a b. [(a, b)] -> ([a], [b])
unzip (forall a b. (a -> Maybe b) -> [a] -> [b]
mapMaybe DirichletFactor -> Maybe (Natural, DirichletFactor)
prim [DirichletFactor]
xs)
        prim :: DirichletFactor -> Maybe (Natural, DirichletFactor)
        prim :: DirichletFactor -> Maybe (Natural, DirichletFactor)
prim DirichletFactor
Two = forall a. Maybe a
Nothing
        prim (OddPrime Prime Natural
p' Word
k Natural
g RootOfUnity
a) = case forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Maybe a
find (Word, Natural) -> Bool
works [(Word, Natural)]
options of
                                     Maybe (Word, Natural)
Nothing -> forall a. HasCallStack => [Char] -> a
error [Char]
"invalid character"
                                     Just (Word
0,Natural
_) -> forall a. Maybe a
Nothing
                                     Just (Word
i,Natural
_) -> forall a. a -> Maybe a
Just (Natural
pforall a b. (Num a, Integral b) => a -> b -> a
^Word
i, Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p' Word
i Natural
g RootOfUnity
a)
          where options :: [(Word, Natural)]
options = (Word
0,Natural
1)forall a. a -> [a] -> [a]
: [(Word
i,Natural
pforall a b. (Num a, Integral b) => a -> b -> a
^(Word
iforall a. Num a => a -> a -> a
-Word
1)forall a. Num a => a -> a -> a
*(Natural
pforall a. Num a => a -> a -> a
-Natural
1)) | Word
i <- [Word
1..Word
k]]
                works :: (Word, Natural) -> Bool
works (Word
_,Natural
phi) = Natural
phi forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
a forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty
                p :: Natural
p = forall a. Prime a -> a
unPrime Prime Natural
p'
        prim (TwoPower Int
k RootOfUnity
a RootOfUnity
b) = case forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Maybe a
find (Int, Natural) -> Bool
worksb [(Int, Natural)]
options of
                                  Maybe (Int, Natural)
Nothing -> forall a. HasCallStack => [Char] -> a
error [Char]
"invalid character"
                                  Just (Int
2,Natural
_) | RootOfUnity
a forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty -> forall a. Maybe a
Nothing
                                  Just (Int
i,Natural
_) -> forall a. a -> Maybe a
Just (forall a. Bits a => Int -> a
bit Int
i :: Natural, Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
i RootOfUnity
a RootOfUnity
b)
          where options :: [(Int, Natural)]
options = [(Int
i, forall a. Bits a => Int -> a
bit (Int
iforall a. Num a => a -> a -> a
-Int
2) :: Natural) | Int
i <- [Int
2..Int
k]]
                worksb :: (Int, Natural) -> Bool
worksb (Int
_,Natural
phi) = Natural
phi forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b forall a. Eq a => a -> a -> Bool
== forall a. Monoid a => a
mempty

#if !MIN_VERSION_base(4,12,0)
newtype Ap f a = Ap { getAp :: f a }
  deriving (Eq, Functor, Applicative, Monad)

instance (Applicative f, Semigroup a) => Semigroup (Ap f a) where
  (<>) = liftA2 (<>)

instance (Applicative f, Semigroup a, Monoid a) => Monoid (Ap f a) where
  mempty = pure mempty
  mappend = (<>)
#endif

-- | Similar to Maybe, but with different Semigroup and Monoid instances.
type OrZero a = Ap Maybe a

-- | 'Ap' 'Nothing'
pattern Zero :: OrZero a
pattern $bZero :: forall a. OrZero a
$mZero :: forall {r} {a}. OrZero a -> ((# #) -> r) -> ((# #) -> r) -> r
Zero = Ap Nothing

-- | 'Ap' ('Just' x)
pattern NonZero :: a -> OrZero a
pattern $bNonZero :: forall a. a -> OrZero a
$mNonZero :: forall {r} {a}. OrZero a -> (a -> r) -> ((# #) -> r) -> r
NonZero x = Ap (Just x)

{-# COMPLETE Zero, NonZero #-}

-- | Interpret an `OrZero` as a number, taking the `Zero` case to be 0.
orZeroToNum :: Num a => (b -> a) -> OrZero b -> a
orZeroToNum :: forall a b. Num a => (b -> a) -> OrZero b -> a
orZeroToNum b -> a
_ OrZero b
Zero = a
0
orZeroToNum b -> a
f (NonZero b
x) = b -> a
f b
x

-- | In general, evaluating a DirichletCharacter at a point involves solving the discrete logarithm
-- problem, which can be hard: the implementations here are around O(sqrt n).
-- However, evaluating a dirichlet character at every point amounts to solving the discrete
-- logarithm problem at every point also, which can be done together in O(n) time, better than
-- using a complex algorithm at each point separately. Thus, if a large number of evaluations
-- of a dirichlet character are required, `evalAll` will be better than `evalGeneral`, since
-- computations can be shared.
evalAll :: forall n. KnownNat n => DirichletCharacter n -> Vector (OrZero RootOfUnity)
evalAll :: forall (n :: Natural).
KnownNat n =>
DirichletCharacter n -> Vector (OrZero RootOfUnity)
evalAll (Generated [DirichletFactor]
xs) = forall a. Int -> (Int -> a) -> Vector a
V.generate (Natural -> Int
naturalToInt Natural
n) Int -> OrZero RootOfUnity
func
  where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        vectors :: [(Int, Vector (OrZero RootOfUnity))]
vectors = forall a b. (a -> b) -> [a] -> [b]
map DirichletFactor -> (Int, Vector (OrZero RootOfUnity))
mkVector [DirichletFactor]
xs
        func :: Int -> OrZero RootOfUnity
        func :: Int -> OrZero RootOfUnity
func Int
m = forall (t :: * -> *) m a.
(Foldable t, Monoid m) =>
(a -> m) -> t a -> m
foldMap (Int, Vector (OrZero RootOfUnity)) -> OrZero RootOfUnity
go [(Int, Vector (OrZero RootOfUnity))]
vectors
          where go :: (Int, Vector (OrZero RootOfUnity)) -> OrZero RootOfUnity
                go :: (Int, Vector (OrZero RootOfUnity)) -> OrZero RootOfUnity
go (Int
modulus,Vector (OrZero RootOfUnity)
v) = Vector (OrZero RootOfUnity)
v forall a. Vector a -> Int -> a
! (Int
m forall a. Integral a => a -> a -> a
`mod` Int
modulus)
        mkVector :: DirichletFactor -> (Int, Vector (OrZero RootOfUnity))
        mkVector :: DirichletFactor -> (Int, Vector (OrZero RootOfUnity))
mkVector DirichletFactor
Two = (Int
2, forall a. [a] -> Vector a
V.fromList [forall a. OrZero a
Zero, forall a. Monoid a => a
mempty])
        mkVector (OddPrime Prime Natural
p Word
k (Natural -> Int
naturalToInt -> Int
g) RootOfUnity
a) = (Int
modulus, Vector (OrZero RootOfUnity)
w)
          where
            p' :: Natural
p' = forall a. Prime a -> a
unPrime Prime Natural
p
            modulus :: Int
modulus = Natural -> Int
naturalToInt (Natural
p'forall a b. (Num a, Integral b) => a -> b -> a
^Word
k) :: Int
            w :: Vector (OrZero RootOfUnity)
w = forall a. (forall s. ST s (MVector s a)) -> Vector a
V.create forall a b. (a -> b) -> a -> b
$ do
              MVector s (OrZero RootOfUnity)
v <- forall (m :: * -> *) a.
PrimMonad m =>
Int -> a -> m (MVector (PrimState m) a)
MV.replicate Int
modulus forall a. OrZero a
Zero
              -- TODO: we're in the ST monad here anyway, could be better to use STRefs to manage
              -- this loop, the current implementation probably doesn't fuse well
              let powers :: [(Int, RootOfUnity)]
powers = forall a. (a -> Maybe a) -> a -> [a]
iterateMaybe (Int, RootOfUnity) -> Maybe (Int, RootOfUnity)
go (Int
1,forall a. Monoid a => a
mempty)
                  go :: (Int, RootOfUnity) -> Maybe (Int, RootOfUnity)
go (Int
m,RootOfUnity
x) = if Int
m' forall a. Ord a => a -> a -> Bool
> Int
1
                                then forall a. a -> Maybe a
Just (Int
m', RootOfUnity
xforall a. Semigroup a => a -> a -> a
<>RootOfUnity
a)
                                else forall a. Maybe a
Nothing
                    where m' :: Int
m' = Int
mforall a. Num a => a -> a -> a
*Int
g forall a. Integral a => a -> a -> a
`mod` Int
modulus
              forall (t :: * -> *) (f :: * -> *) a b.
(Foldable t, Applicative f) =>
t a -> (a -> f b) -> f ()
for_ [(Int, RootOfUnity)]
powers forall a b. (a -> b) -> a -> b
$ \(Int
m,RootOfUnity
x) -> forall (m :: * -> *) a.
PrimMonad m =>
MVector (PrimState m) a -> Int -> a -> m ()
MV.unsafeWrite MVector s (OrZero RootOfUnity)
v Int
m (forall a. a -> OrZero a
NonZero RootOfUnity
x)
              -- don't bother with bounds check since m was reduced mod p^k
              forall (m :: * -> *) a. Monad m => a -> m a
return MVector s (OrZero RootOfUnity)
v
        -- for powers of two we use lambda directly instead, since the generators of the cyclic
        -- groups aren't obvious; it's possible to get them though:
        -- 5^(lambda(5)^{-1} mod 2^(p-2)) mod 2^p
        mkVector (TwoPower Int
k RootOfUnity
a RootOfUnity
b) = (Int
modulus, Vector (OrZero RootOfUnity)
w)
          where
            modulus :: Int
modulus = forall a. Bits a => Int -> a
bit Int
k
            w :: Vector (OrZero RootOfUnity)
w = forall a. Int -> (Int -> a) -> Vector a
V.generate Int
modulus Int -> OrZero RootOfUnity
f
            f :: Int -> OrZero RootOfUnity
f Int
m
              | forall a. Integral a => a -> Bool
even Int
m = forall a. OrZero a
Zero
              | Bool
otherwise = forall a. a -> OrZero a
NonZero ((if forall a. Bits a => a -> Int -> Bool
testBit Int
m Int
1 then RootOfUnity
a else forall a. Monoid a => a
mempty) forall a. Semigroup a => a -> a -> a
<> Integer -> Int -> Integer
lambda (forall a. Integral a => a -> Integer
toInteger Int
m'') Int
k forall a b. (Semigroup a, Integral b) => b -> a -> a
`stimes` RootOfUnity
b)
              where m'' :: Int
m'' = forall p. (Bits p, Num p) => Int -> p -> p
thingy Int
k Int
m

-- somewhere between unfoldr and iterate
iterateMaybe :: (a -> Maybe a) -> a -> [a]
iterateMaybe :: forall a. (a -> Maybe a) -> a -> [a]
iterateMaybe a -> Maybe a
f a
x = forall b a. (b -> Maybe (a, b)) -> b -> [a]
unfoldr (forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (\a
t -> (a
t, a -> Maybe a
f a
t))) (forall a. a -> Maybe a
Just a
x)

-- | Attempt to construct a character from its table of values.
-- An inverse to `evalAll`, defined only on its image.
fromTable :: forall n. KnownNat n => Vector (OrZero RootOfUnity) -> Maybe (DirichletCharacter n)
fromTable :: forall (n :: Natural).
KnownNat n =>
Vector (OrZero RootOfUnity) -> Maybe (DirichletCharacter n)
fromTable Vector (OrZero RootOfUnity)
v = if forall (t :: * -> *) a. Foldable t => t a -> Int
length Vector (OrZero RootOfUnity)
v forall a. Eq a => a -> a -> Bool
== Natural -> Int
naturalToInt Natural
n
                 then forall (t :: * -> *) (f :: * -> *) a b.
(Traversable t, Applicative f) =>
(a -> f b) -> t a -> f (t b)
traverse Template -> Maybe DirichletFactor
makeFactor [Template]
tmpl forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= DirichletCharacter n -> Maybe (DirichletCharacter n)
check forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (n :: Natural). [DirichletFactor] -> DirichletCharacter n
Generated
                 else forall a. Maybe a
Nothing
  where n :: Natural
n = forall (n :: Natural) (proxy :: Natural -> *).
KnownNat n =>
proxy n -> Natural
natVal (forall {k} (t :: k). Proxy t
Proxy :: Proxy n)
        n' :: Integer
n' = Natural -> Integer
naturalToInteger Natural
n :: Integer
        tmpl :: [Template]
tmpl = forall a b. (a, b) -> b
snd (Natural -> (Product Natural, [Template])
mkTemplate Natural
n)
        check :: DirichletCharacter n -> Maybe (DirichletCharacter n)
        check :: DirichletCharacter n -> Maybe (DirichletCharacter n)
check DirichletCharacter n
chi = if forall (n :: Natural).
KnownNat n =>
DirichletCharacter n -> Vector (OrZero RootOfUnity)
evalAll DirichletCharacter n
chi forall a. Eq a => a -> a -> Bool
== Vector (OrZero RootOfUnity)
v then forall a. a -> Maybe a
Just DirichletCharacter n
chi else forall a. Maybe a
Nothing
        makeFactor :: Template -> Maybe DirichletFactor
        makeFactor :: Template -> Maybe DirichletFactor
makeFactor Template
TwoTemplate = forall a. a -> Maybe a
Just DirichletFactor
Two
        makeFactor (TwoPTemplate Int
k Natural
_) = Int -> RootOfUnity -> RootOfUnity -> DirichletFactor
TwoPower Int
k forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> (Integer, Integer) -> Maybe RootOfUnity
getValue (-Integer
1,forall a. Bits a => Int -> a
bit Int
k) forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> (Integer, Integer) -> Maybe RootOfUnity
getValue (Int -> Integer
exp4 Int
k, forall a. Bits a => Int -> a
bit Int
k)
        makeFactor (OddTemplate Prime Natural
p Word
k Natural
g Natural
_) = Prime Natural -> Word -> Natural -> RootOfUnity -> DirichletFactor
OddPrime Prime Natural
p Word
k Natural
g forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> (Integer, Integer) -> Maybe RootOfUnity
getValue (forall a. Integral a => a -> Integer
toInteger Natural
g, forall a. Integral a => a -> Integer
toInteger (forall a. Prime a -> a
unPrime Prime Natural
p)forall a b. (Num a, Integral b) => a -> b -> a
^Word
k)
        getValue :: (Integer, Integer) -> Maybe RootOfUnity
        getValue :: (Integer, Integer) -> Maybe RootOfUnity
getValue (Integer
g, Integer
m) = forall {k} (f :: k -> *) (a :: k). Ap f a -> f a
getAp (Vector (OrZero RootOfUnity)
v forall a. Vector a -> Int -> a
! forall a. Num a => Integer -> a
fromInteger (forall a b. (a, b) -> a
fst (forall a. HasCallStack => Maybe a -> a
fromJust (forall a.
(Eq a, Ring a, Euclidean a) =>
(a, a) -> (a, a) -> Maybe (a, a)
chinese (Integer
g, Integer
m) (Integer
1, Integer
n' forall a. Integral a => a -> a -> a
`quot` Integer
m))) forall a. Integral a => a -> a -> a
`mod` Integer
n'))

exp4terms :: [Rational]
exp4terms :: [Ratio Integer]
exp4terms = [Integer
4forall a b. (Num a, Integral b) => a -> b -> a
^Integer
k forall a. Integral a => a -> a -> Ratio a
% forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
product [Integer
1..Integer
k] | Integer
k <- [Integer
0..]]

-- For reasons that aren't clear to me, `exp4` gives the inverse of 1 under lambda, so it gives the generator
-- This is the same as https://oeis.org/A320814
-- In particular, lambda (exp4 n) n == 1 (for n >= 3)
-- I've verified this for 3 <= n <= 2000, so the reasoning in fromTable should be accurate for moduli below 2^2000
exp4 :: Int -> Integer
exp4 :: Int -> Integer
exp4 Int
n
  = (forall a. Integral a => a -> a -> a
`mod` forall a. Bits a => Int -> a
bit Int
n)
  forall a b. (a -> b) -> a -> b
$ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum
  forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map (\Ratio Integer
q -> (forall a. Ratio a -> a
numerator Ratio Integer
q forall a. Num a => a -> a -> a
* forall a. a -> Maybe a -> a
fromMaybe (forall a. HasCallStack => [Char] -> a
error [Char]
"error in exp4") (Integer -> Integer -> Maybe Integer
recipMod (forall a. Ratio a -> a
denominator Ratio Integer
q) (forall a. Bits a => Int -> a
bit Int
n))) forall a. Integral a => a -> a -> a
`mod` forall a. Bits a => Int -> a
bit Int
n)
  forall a b. (a -> b) -> a -> b
$ forall a. Int -> [a] -> [a]
take Int
n [Ratio Integer]
exp4terms