    -- Copyright (c) David Amos, 2008-2015. All rights reserved.

    -- |Constructions of the finite geometries AG(n,Fq) and PG(n,Fq), their points, lines and flats,

    -- together with the incidence graphs between points and lines.

    module Math.Combinatorics.FiniteGeometry where

    import Prelude hiding ( (*>) )

    import Data.List as L
    import qualified Data.Set as S

    import Math.Common.ListSet (toListSet)
    import Math.Core.Utils

    import Math.Core.Field
    import Math.Algebra.LinearAlgebra -- hiding ( det )

    import Math.Combinatorics.Graph
    import Math.Combinatorics.GraphAuts -- for use in GHCi

    import Math.Algebra.Group.PermutationGroup hiding (elts) -- for use in GHCi

    import Math.Algebra.Group.SchreierSims as SS hiding (elts) -- for use in GHCi

    -- !! The following two functions previously required (FiniteField a) as context

    -- but this has been temporarily removed to enable them to work with Math.Core.Field

    -- |ptsAG n fq returns the points of the affine geometry AG(n,Fq), where fq are the elements of Fq

    ptsAG :: Int -> [a] -> [[a]]
    ptsAG :: Int -> [a] -> [[a]]
    ptsAG 0 fq :: [a]
    fq = [[]]
    ptsAG n :: Int
    n fq :: [a]
    fq = [a
    xa -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :[a]
    xs | a
    x <- [a]
    fq, [a]
    xs <- Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG (Int
    nInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) [a]
    fq]

    -- |ptsPG n fq returns the points of the projective geometry PG(n,Fq), where fq are the elements of Fq

    ptsPG :: Num a => Int -> [a] -> [[a]]
    ptsPG :: Int -> [a] -> [[a]]
    ptsPG 0 _ = [[1]]
    ptsPG n :: Int
    n fq :: [a]
    fq = ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map (0a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :) (Int -> [a] -> [[a]]
    forall a. Num a => Int -> [a] -> [[a]]
    ptsPG (Int
    nInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) [a]
    fq) [[a]] -> [[a]] -> [[a]]
    forall a. [a] -> [a] -> [a]
    ++ ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map (1a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :) (Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG Int
    n [a]
    fq)

    -- "projective normal form"

    pnf :: [a] -> [a]
    pnf (0:xs :: [a]
    xs) = 0 a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    : [a] -> [a]
    pnf [a]
    xs
    pnf (1:xs :: [a]
    xs) = 1 a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    : [a]
    xs
    pnf (x :: a
    x:xs :: [a]
    xs) = 1 a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    : (a -> a) -> [a] -> [a]
    forall a b. (a -> b) -> [a] -> [b]
    map (a -> a -> a
    forall a. Num a => a -> a -> a
    * a
    x') [a]
    xs where x' :: a
    x' = a -> a
    forall a. Fractional a => a -> a
    recip a
    x

    ispnf :: [a] -> Bool
    ispnf (0:xs :: [a]
    xs) = [a] -> Bool
    ispnf [a]
    xs
    ispnf (1:xs :: [a]
    xs) = Bool
    True
    ispnf _ = Bool
    False

    -- closure of points in AG(n,Fq)

    -- given p1, .., pk, we're looking for all a1 p1 + ... + ak pk, s.t. a1 + ... + ak = 1

    -- if m is the matrix with p1, .., pk as rows, and vs are the vectors [a1, .., ak]

    -- then this is the same as [v <*>> m | v <- vs] == [m' <<*> v | v <- vs]

    -- |Given a list of points in AG(n,Fq), return their closure, the smallest flat containing them

    closureAG :: (Num a, Ord a, FinSet a) => [[a]] -> [[a]]
    closureAG :: [[a]] -> [[a]]
    closureAG ps :: [[a]]
    ps =
        let vs :: [[a]]
    vs = [ (1 a -> a -> a
    forall a. Num a => a -> a -> a
    - [a] -> a
    forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
    sum [a]
    xs) a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    : [a]
    xs | [a]
    xs <- Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG (Int
    kInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) [a]
    fq ] -- k-vectors over fq whose sum is 1

        in [[a]] -> [[a]]
    forall b. Ord b => [b] -> [b]
    toListSet [[[a]]
    m' [[a]] -> [a] -> [a]
    forall a. Num a => [[a]] -> [a] -> [a]
    <<*> [a]
    v | [a]
    v <- [[a]]
    vs]
        where k :: Int
    k = [[a]] -> Int
    forall (t :: * -> *) a. Foldable t => t a -> Int
    length [[a]]
    ps -- the dimension of the flat (assuming ps are independent)

              m' :: [[a]]
    m' = [[a]] -> [[a]]
    forall a. [[a]] -> [[a]]
    L.transpose [[a]]
    ps
              fq :: [a]
    fq = [a]
    forall x. FinSet x => [x]
    elts
    -- toListSet call sorts the result, and also removes duplicates in case the points weren't independent

    {-
    closureAG ps =
        let vs = [ (1 - sum xs) : xs | xs <- ptsAG (k-1) fq ] -- k-vectors over fq whose sum is 1
        in toListSet [foldl1 (<+>) $ zipWith (*>) v ps | v <- vs]
        where k = length ps        -- the dimension of the flat
              fq = eltsFq undefined
    -}

    lineAG :: [[a]] -> [[a]]
    lineAG [p1 :: [a]
    p1,p2 :: [a]
    p2] = [[a]] -> [[a]]
    forall b. Ord b => [b] -> [b]
    L.sort [ [a]
    p1 [a] -> [a] -> [a]
    forall a. Num a => [a] -> [a] -> [a]
    <+> (a
    c a -> [a] -> [a]
    forall a. Num a => a -> [a] -> [a]
    *> [a]
    dp) | a
    c <- [a]
    fq ] where
        dp :: [a]
    dp = [a]
    p2 [a] -> [a] -> [a]
    forall a. Num a => [a] -> [a] -> [a]
    <-> [a]
    p1
        fq :: [a]
    fq = [a]
    forall x. FinSet x => [x]
    elts

    -- closure of points in PG(n,Fq)

    -- take all linear combinations of the points (ie the subspace generated by the points, considered as points in Fq ^(n+1) )

    -- then discard all which aren't in PNF (thus dropping back into PG(n,Fq))

    -- |Given a set of points in PG(n,Fq), return their closure, the smallest flat containing them

    closurePG :: (Num a, Ord a, FinSet a) => [[a]] -> [[a]]
    closurePG :: [[a]] -> [[a]]
    closurePG ps :: [[a]]
    ps = [[a]] -> [[a]]
    forall b. Ord b => [b] -> [b]
    toListSet ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ ([a] -> Bool) -> [[a]] -> [[a]]
    forall a. (a -> Bool) -> [a] -> [a]
    filter [a] -> Bool
    forall a. (Eq a, Num a) => [a] -> Bool
    ispnf ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map ([a] -> [[a]] -> [a]
    forall a. Num a => [a] -> [[a]] -> [a]
    <*>> [[a]]
    ps) ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG Int
    k [a]
    fq where
        k :: Int
    k = [[a]] -> Int
    forall (t :: * -> *) a. Foldable t => t a -> Int
    length [[a]]
    ps
        fq :: [a]
    fq = [a]
    forall x. FinSet x => [x]
    elts
    -- toListSet call sorts the result, and also removes duplicates in case the points weren't independent

    linePG :: [[a]] -> [[a]]
    linePG [p1 :: [a]
    p1,p2 :: [a]
    p2] = [[a]] -> [[a]]
    forall b. Ord b => [b] -> [b]
    toListSet ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ ([a] -> Bool) -> [[a]] -> [[a]]
    forall a. (a -> Bool) -> [a] -> [a]
    filter [a] -> Bool
    forall a. (Eq a, Num a) => [a] -> Bool
    ispnf [(a
    a a -> [a] -> [a]
    forall a. Num a => a -> [a] -> [a]
    *> [a]
    p1) [a] -> [a] -> [a]
    forall a. Num a => [a] -> [a] -> [a]
    <+> (a
    b a -> [a] -> [a]
    forall a. Num a => a -> [a] -> [a]
    *> [a]
    p2) | a
    a <- [a]
    fq, a
    b <- [a]
    fq]
        where fq :: [a]
    fq = [a]
    forall x. FinSet x => [x]
    elts

    -- van Lint & Wilson, p325, 332

    qtorial :: b -> a -> a
    qtorial n :: b
    n q :: a
    q | b
    n b -> b -> Bool
    forall a. Ord a => a -> a -> Bool
    >= 0 = [a] -> a
    forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
    product [(a
    qa -> b -> a
    forall a b. (Num a, Integral b) => a -> b -> a
    ^b
    i a -> a -> a
    forall a. Num a => a -> a -> a
    - 1) a -> a -> a
    forall a. Integral a => a -> a -> a
    `div` (a
    qa -> a -> a
    forall a. Num a => a -> a -> a
    -1) | b
    i <- [1..b
    n]]

    -- van Lint & Wilson, p326

    qnomial :: b -> b -> a -> a
    qnomial n :: b
    n k :: b
    k q :: a
    q = (b
    n b -> a -> a
    forall b a. (Integral a, Integral b) => b -> a -> a
    `qtorial` a
    q) a -> a -> a
    forall a. Integral a => a -> a -> a
    `div` ( (b
    k b -> a -> a
    forall b a. (Integral a, Integral b) => b -> a -> a
    `qtorial` a
    q) a -> a -> a
    forall a. Num a => a -> a -> a
    * ((b
    nb -> b -> b
    forall a. Num a => a -> a -> a
    -b
    k) b -> a -> a
    forall b a. (Integral a, Integral b) => b -> a -> a
    `qtorial` a
    q) )

    -- Cameron, p129

    numFlatsPG :: b -> a -> b -> a
    numFlatsPG n :: b
    n q :: a
    q k :: b
    k = b -> b -> a -> a
    forall a b. (Integral a, Integral b) => b -> b -> a -> a
    qnomial (b
    nb -> b -> b
    forall a. Num a => a -> a -> a
    +1) (b
    kb -> b -> b
    forall a. Num a => a -> a -> a
    +1) a
    q -- because it's the number of subspaces in AG n+1

    -- Cameron, p137

    numFlatsAG :: b -> a -> b -> a
    numFlatsAG n :: b
    n q :: a
    q k :: b
    k = a
    qa -> b -> a
    forall a b. (Num a, Integral b) => a -> b -> a
    ^(b
    nb -> b -> b
    forall a. Num a => a -> a -> a
    -b
    k) a -> a -> a
    forall a. Num a => a -> a -> a
    * b -> b -> a -> a
    forall a b. (Integral a, Integral b) => b -> b -> a -> a
    qnomial b
    n b
    k a
    q

    qtorials :: b -> [b]
    qtorials q :: b
    q = (b -> b -> b) -> b -> [b] -> [b]
    forall b a. (b -> a -> b) -> b -> [a] -> [b]
    scanl b -> b -> b
    forall a. Num a => a -> a -> a
    (*) 1 [(b
    qb -> Integer -> b
    forall a b. (Num a, Integral b) => a -> b -> a
    ^Integer
    i b -> b -> b
    forall a. Num a => a -> a -> a
    - 1) b -> b -> b
    forall a. Integral a => a -> a -> a
    `div` (b
    qb -> b -> b
    forall a. Num a => a -> a -> a
    -1) | Integer
    i <- [1..]]

    qnomials :: a -> [[a]]
    qnomials q :: a
    q = ([a] -> [a]) -> [a] -> [[a]]
    forall a. (a -> a) -> a -> [a]
    iterate [a] -> [a]
    succ [1] where
        succ :: [a] -> [a]
    succ xs :: [a]
    xs = (a -> a -> a -> a) -> [a] -> [a] -> [a] -> [a]
    forall a b c d. (a -> b -> c -> d) -> [a] -> [b] -> [c] -> [d]
    L.zipWith3 (\l :: a
    l r :: a
    r c :: a
    c -> a
    la -> a -> a
    forall a. Num a => a -> a -> a
    +a
    ca -> a -> a
    forall a. Num a => a -> a -> a
    *a
    r) (0a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :[a]
    xs) ([a]
    xs[a] -> [a] -> [a]
    forall a. [a] -> [a] -> [a]
    ++[0]) ((a -> a) -> a -> [a]
    forall a. (a -> a) -> a -> [a]
    iterate (a -> a -> a
    forall a. Num a => a -> a -> a
    *a
    q) 1)
        -- succ xs = zipWith (+) (0:xs) $ zipWith (*) (xs++[0]) $ iterate (*q) 1

    -- This amounts to saying

    -- [n+1,k]_q = [n,k-1]_q + q^k [n,k]_q

    -- Cameron, Combinatorics, p126

    -- FLATS VIA REDUCED ROW ECHELON FORMS

    -- Suggested by Cameron p125

    data ZeroOneStar = Zero | One | Star deriving (ZeroOneStar -> ZeroOneStar -> Bool
    (ZeroOneStar -> ZeroOneStar -> Bool)
    -> (ZeroOneStar -> ZeroOneStar -> Bool) -> Eq ZeroOneStar
    forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
    /= :: ZeroOneStar -> ZeroOneStar -> Bool
    $c/= :: ZeroOneStar -> ZeroOneStar -> Bool
    == :: ZeroOneStar -> ZeroOneStar -> Bool
    $c== :: ZeroOneStar -> ZeroOneStar -> Bool
    Eq)

    instance Show ZeroOneStar where
        show :: ZeroOneStar -> String
    show Zero = "0"
        show One  = "1"
        show Star = "*"

    -- reduced row echelon forms

    rrefs :: Int -> Int -> [[[ZeroOneStar]]]
    rrefs n :: Int
    n k :: Int
    k = ([Int] -> [[ZeroOneStar]]) -> [[Int]] -> [[[ZeroOneStar]]]
    forall a b. (a -> b) -> [a] -> [b]
    map (Int -> Int -> [Int] -> [[ZeroOneStar]]
    rref 1 1) (Int -> [Int] -> [[Int]]
    forall a. Int -> [a] -> [[a]]
    combinationsOf Int
    k [1..Int
    n]) where
        rref :: Int -> Int -> [Int] -> [[ZeroOneStar]]
    rref r :: Int
    r c :: Int
    c (x :: Int
    x:xs :: [Int]
    xs) =
            if Int
    c Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    x
            then (ZeroOneStar -> [ZeroOneStar] -> [ZeroOneStar])
    -> [ZeroOneStar] -> [[ZeroOneStar]] -> [[ZeroOneStar]]
    forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
    zipWith (:) (Int -> [ZeroOneStar]
    oneColumn Int
    r) (Int -> Int -> [Int] -> [[ZeroOneStar]]
    rref (Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1) (Int
    cInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1) [Int]
    xs)
            else (ZeroOneStar -> [ZeroOneStar] -> [ZeroOneStar])
    -> [ZeroOneStar] -> [[ZeroOneStar]] -> [[ZeroOneStar]]
    forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
    zipWith (:) (Int -> [ZeroOneStar]
    starColumn Int
    r) (Int -> Int -> [Int] -> [[ZeroOneStar]]
    rref Int
    r (Int
    cInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1) (Int
    xInt -> [Int] -> [Int]
    forall a. a -> [a] -> [a]
    :[Int]
    xs))
        rref _ c :: Int
    c [] = Int -> [ZeroOneStar] -> [[ZeroOneStar]]
    forall a. Int -> a -> [a]
    replicate Int
    k (Int -> ZeroOneStar -> [ZeroOneStar]
    forall a. Int -> a -> [a]
    replicate (Int
    nInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1Int -> Int -> Int
    forall a. Num a => a -> a -> a
    -Int
    c) ZeroOneStar
    Star)
        oneColumn :: Int -> [ZeroOneStar]
    oneColumn r :: Int
    r = Int -> ZeroOneStar -> [ZeroOneStar]
    forall a. Int -> a -> [a]
    replicate (Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) ZeroOneStar
    Zero [ZeroOneStar] -> [ZeroOneStar] -> [ZeroOneStar]
    forall a. [a] -> [a] -> [a]
    ++ ZeroOneStar
    One ZeroOneStar -> [ZeroOneStar] -> [ZeroOneStar]
    forall a. a -> [a] -> [a]
    : Int -> ZeroOneStar -> [ZeroOneStar]
    forall a. Int -> a -> [a]
    replicate (Int
    kInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -Int
    r) ZeroOneStar
    Zero
        starColumn :: Int -> [ZeroOneStar]
    starColumn r :: Int
    r = Int -> ZeroOneStar -> [ZeroOneStar]
    forall a. Int -> a -> [a]
    replicate (Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) ZeroOneStar
    Star [ZeroOneStar] -> [ZeroOneStar] -> [ZeroOneStar]
    forall a. [a] -> [a] -> [a]
    ++ Int -> ZeroOneStar -> [ZeroOneStar]
    forall a. Int -> a -> [a]
    replicate (Int
    kInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1Int -> Int -> Int
    forall a. Num a => a -> a -> a
    -Int
    r) ZeroOneStar
    Zero

    -- flatsPG :: (FiniteField a) => Int -> [a] -> Int -> [[[a]]]

    -- |@flatsPG n fq k@ returns the k-flats in PG(n,Fq), where fq are the elements of Fq.

    -- The returned flats are represented as matrices in reduced row echelon form,

    -- the rows of which are the points that generate the flat.

    -- The full set of points in the flat can be recovered by calling 'closurePG'

    flatsPG :: (Eq a, Num a) => Int -> [a] -> Int -> [[[a]]]
    flatsPG :: Int -> [a] -> Int -> [[[a]]]
    flatsPG n :: Int
    n fq :: [a]
    fq k :: Int
    k = ([[ZeroOneStar]] -> [[[a]]]) -> [[[ZeroOneStar]]] -> [[[a]]]
    forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
    concatMap [[ZeroOneStar]] -> [[[a]]]
    substStars ([[[ZeroOneStar]]] -> [[[a]]]) -> [[[ZeroOneStar]]] -> [[[a]]]
    forall a b. (a -> b) -> a -> b
    $ Int -> Int -> [[[ZeroOneStar]]]
    rrefs (Int
    nInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1) (Int
    kInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +1) where
        substStars :: [[ZeroOneStar]] -> [[[a]]]
    substStars (r :: [ZeroOneStar]
    r:rs :: [[ZeroOneStar]]
    rs) = [[a]
    r'[a] -> [[a]] -> [[a]]
    forall a. a -> [a] -> [a]
    :[[a]]
    rs' | [a]
    r' <- [ZeroOneStar] -> [[a]]
    substStars' [ZeroOneStar]
    r, [[a]]
    rs' <- [[ZeroOneStar]] -> [[[a]]]
    substStars [[ZeroOneStar]]
    rs]
        substStars [] = [[]]
        substStars' :: [ZeroOneStar] -> [[a]]
    substStars' (Star:xs :: [ZeroOneStar]
    xs) = [a
    x'a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :[a]
    xs' | a
    x' <- [a]
    fq, [a]
    xs' <- [ZeroOneStar] -> [[a]]
    substStars' [ZeroOneStar]
    xs]
        substStars' (Zero:xs :: [ZeroOneStar]
    xs) = ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map (0a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :) ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ [ZeroOneStar] -> [[a]]
    substStars' [ZeroOneStar]
    xs
        substStars' (One:xs :: [ZeroOneStar]
    xs) = ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map (1a -> [a] -> [a]
    forall a. a -> [a] -> [a]
    :) ([[a]] -> [[a]]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> a -> b
    $ [ZeroOneStar] -> [[a]]
    substStars' [ZeroOneStar]
    xs
        substStars' [] = [[]]

    -- Flats in AG(n,Fq) are just the flats in PG(n,Fq) which are not "at infinity"

    -- flatsAG :: (FiniteField a) => Int -> [a] -> Int -> [[[a]]]

    -- |flatsAG n fq k returns the k-flats in AG(n,Fq), where fq are the elements of Fq.

    flatsAG :: (Eq a, Num a) => Int -> [a] -> Int -> [[[a]]]
    flatsAG :: Int -> [a] -> Int -> [[[a]]]
    flatsAG n :: Int
    n fq :: [a]
    fq k :: Int
    k = [([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map [a] -> [a]
    forall a. [a] -> [a]
    tail ([a]
    r [a] -> [[a]] -> [[a]]
    forall a. a -> [a] -> [a]
    : ([a] -> [a]) -> [[a]] -> [[a]]
    forall a b. (a -> b) -> [a] -> [b]
    map ([a]
    r [a] -> [a] -> [a]
    forall a. Num a => [a] -> [a] -> [a]
    <+>) [[a]]
    rs) | r :: [a]
    r:rs :: [[a]]
    rs <- Int -> [a] -> Int -> [[[a]]]
    forall a. (Eq a, Num a) => Int -> [a] -> Int -> [[[a]]]
    flatsPG Int
    n [a]
    fq Int
    k, [a] -> a
    forall a. [a] -> a
    head [a]
    r a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == 1]
    -- The head r == 1 condition is saying that we want points which are in the "finite" part of PG(n,Fq), not points at infinity

    -- The reason we add r to each of the rs is to bring them into the "finite" part

    -- (If you don't do this, it can lead to incorrect results, for example some of the flats having the same closure)

    -- linesPG :: (FiniteField a) => Int -> [a] -> [[[a]]]

    -- |The lines (1-flats) in PG(n,fq)

    linesPG :: (Eq a, Num a) => Int -> [a] -> [[[a]]]
    linesPG :: Int -> [a] -> [[[a]]]
    linesPG n :: Int
    n fq :: [a]
    fq = Int -> [a] -> Int -> [[[a]]]
    forall a. (Eq a, Num a) => Int -> [a] -> Int -> [[[a]]]
    flatsPG Int
    n [a]
    fq 1

    -- linesAG :: (FiniteField a) => Int -> [a] -> [[[a]]]

    -- |The lines (1-flats) in AG(n,fq)

    linesAG :: (Eq a, Num a) => Int -> [a] -> [[[a]]]
    linesAG :: Int -> [a] -> [[[a]]]
    linesAG n :: Int
    n fq :: [a]
    fq = Int -> [a] -> Int -> [[[a]]]
    forall a. (Eq a, Num a) => Int -> [a] -> Int -> [[[a]]]
    flatsAG Int
    n [a]
    fq 1

    -- almost certainly not as efficient as linesAG, because requires lineAG/closureAG call

    -- among all pairs of distinct points, select those which are the first two in the line they generate

    linesAG1 :: Int -> [a] -> [[[a]]]
    linesAG1 n :: Int
    n fq :: [a]
    fq = [ [[a]
    x,[a]
    y] | [x :: [a]
    x,y :: [a]
    y] <- Int -> [[a]] -> [[[a]]]
    forall a. Int -> [a] -> [[a]]
    combinationsOf 2 (Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG Int
    n [a]
    fq),
                              [[a]
    x,[a]
    y] [[a]] -> [[a]] -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int -> [[a]] -> [[a]]
    forall a. Int -> [a] -> [a]
    take 2 ([[a]] -> [[a]]
    forall a. (Ord a, FinSet a, Num a) => [[a]] -> [[a]]
    lineAG [[a]
    x,[a]
    y]) ]
    -- the point of the condition at the end is to avoid listing the same line more than once

    -- almost certainly not as efficient as linesAG, because requires lineAG/closureAG call

    -- a line in AG(n,fq) is a translation (x) of a line through the origin (y)

    linesAG2 :: Int -> [a] -> [[[a]]]
    linesAG2 n :: Int
    n fq :: [a]
    fq = [ [[a]
    x,[a]
    z] | [a]
    x <- Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG Int
    n [a]
    fq, [a]
    y <- Int -> [a] -> [[a]]
    forall a. Num a => Int -> [a] -> [[a]]
    ptsPG (Int
    nInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -1) [a]
    fq,
                              [a]
    z <- [[a]
    x [a] -> [a] -> [a]
    forall a. Num a => [a] -> [a] -> [a]
    <+> [a]
    y], [[a]
    x,[a]
    z] [[a]] -> [[a]] -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int -> [[a]] -> [[a]]
    forall a. Int -> [a] -> [a]
    take 2 ([[a]] -> [[a]]
    forall a. (Num a, Ord a, FinSet a) => [[a]] -> [[a]]
    closureAG [[a]
    x,[a]
    z]) ]

    -- INCIDENCE GRAPH

    -- |Incidence graph of PG(n,fq), considered as an incidence structure between points and lines

    incidenceGraphPG :: (Num a, Ord a, FinSet a) => Int -> [a] -> Graph (Either [a] [[a]])
    incidenceGraphPG :: Int -> [a] -> Graph (Either [a] [[a]])
    incidenceGraphPG n :: Int
    n fq :: [a]
    fq = [Either [a] [[a]]]
    -> [[Either [a] [[a]]]] -> Graph (Either [a] [[a]])
    forall a. [a] -> [[a]] -> Graph a
    G [Either [a] [[a]]]
    vs [[Either [a] [[a]]]]
    es where
        points :: [[a]]
    points = Int -> [a] -> [[a]]
    forall a. Num a => Int -> [a] -> [[a]]
    ptsPG Int
    n [a]
    fq
        lines :: [[[a]]]
    lines = Int -> [a] -> [[[a]]]
    forall a. (Eq a, Num a) => Int -> [a] -> [[[a]]]
    linesPG Int
    n [a]
    fq
        vs :: [Either [a] [[a]]]
    vs = [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall b. Ord b => [b] -> [b]
    L.sort ([Either [a] [[a]]] -> [Either [a] [[a]]])
    -> [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> a -> b
    $ ([a] -> Either [a] [[a]]) -> [[a]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> [a] -> [b]
    map [a] -> Either [a] [[a]]
    forall a b. a -> Either a b
    Left [[a]]
    points [Either [a] [[a]]] -> [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall a. [a] -> [a] -> [a]
    ++ ([[a]] -> Either [a] [[a]]) -> [[[a]]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> [a] -> [b]
    map [[a]] -> Either [a] [[a]]
    forall a b. b -> Either a b
    Right [[[a]]]
    lines
        es :: [[Either [a] [[a]]]]
    es = [[Either [a] [[a]]]] -> [[Either [a] [[a]]]]
    forall b. Ord b => [b] -> [b]
    L.sort [ [[a] -> Either [a] [[a]]
    forall a b. a -> Either a b
    Left [a]
    x, [[a]] -> Either [a] [[a]]
    forall a b. b -> Either a b
    Right [[a]]
    b] | [[a]]
    b <- [[[a]]]
    lines, [a]
    x <- [[a]] -> [[a]]
    forall a. (Num a, Ord a, FinSet a) => [[a]] -> [[a]]
    closurePG [[a]]
    b]
    -- Could also consider incidence structure between points and planes, etc

    -- incidenceAuts (incidenceGraphPG n fq) == PGL(n+1,fq) * auts fq

    -- For example, incidenceAuts (incidenceGraphPG 2 f4) =

    -- PGL(3,f4) * auts f4

    -- where PGL(3,f4)/PSL(3,f4) == f4* (multiplicative group of f4),

    -- and auts f4 == { 1, x -> x^2 } (the Frobenius aut)

    -- The full group is called PGammaL(3,f4)

    -- |Incidence graph of AG(n,fq), considered as an incidence structure between points and lines

    incidenceGraphAG :: (Num a, Ord a, FinSet a) => Int -> [a] -> Graph (Either [a] [[a]])
    incidenceGraphAG :: Int -> [a] -> Graph (Either [a] [[a]])
    incidenceGraphAG n :: Int
    n fq :: [a]
    fq = [Either [a] [[a]]]
    -> [[Either [a] [[a]]]] -> Graph (Either [a] [[a]])
    forall a. [a] -> [[a]] -> Graph a
    G [Either [a] [[a]]]
    vs [[Either [a] [[a]]]]
    es where
        points :: [[a]]
    points = Int -> [a] -> [[a]]
    forall a. Int -> [a] -> [[a]]
    ptsAG Int
    n [a]
    fq
        lines :: [[[a]]]
    lines = Int -> [a] -> [[[a]]]
    forall a. (Eq a, Num a) => Int -> [a] -> [[[a]]]
    linesAG Int
    n [a]
    fq
        vs :: [Either [a] [[a]]]
    vs = [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall b. Ord b => [b] -> [b]
    L.sort ([Either [a] [[a]]] -> [Either [a] [[a]]])
    -> [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> a -> b
    $ ([a] -> Either [a] [[a]]) -> [[a]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> [a] -> [b]
    map [a] -> Either [a] [[a]]
    forall a b. a -> Either a b
    Left [[a]]
    points [Either [a] [[a]]] -> [Either [a] [[a]]] -> [Either [a] [[a]]]
    forall a. [a] -> [a] -> [a]
    ++ ([[a]] -> Either [a] [[a]]) -> [[[a]]] -> [Either [a] [[a]]]
    forall a b. (a -> b) -> [a] -> [b]
    map [[a]] -> Either [a] [[a]]
    forall a b. b -> Either a b
    Right [[[a]]]
    lines
        es :: [[Either [a] [[a]]]]
    es = [[Either [a] [[a]]]] -> [[Either [a] [[a]]]]
    forall b. Ord b => [b] -> [b]
    L.sort [ [[a] -> Either [a] [[a]]
    forall a b. a -> Either a b
    Left [a]
    x, [[a]] -> Either [a] [[a]]
    forall a b. b -> Either a b
    Right [[a]]
    b] | [[a]]
    b <- [[[a]]]
    lines, [a]
    x <- [[a]] -> [[a]]
    forall a. (Num a, Ord a, FinSet a) => [[a]] -> [[a]]
    closureAG [[a]]
    b]

    -- incidenceAuts (incidenceGraphAG n fq) == Aff(n,fq) * auts fq

    -- where Aff(n,fq), the affine group, is the semi-direct product GL(n,fq) * (fq^n)+

    -- where (fq^n)+ is the additive group of translations

    -- Each elt of Aff(n,fq) is of the form x -> ax + b, where a <- GL(n,fq), b <- (fq^n)+

    orderGL :: b -> a -> a
    orderGL n :: b
    n q :: a
    q = [a] -> a
    forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
    product [a
    qa -> b -> a
    forall a b. (Num a, Integral b) => a -> b -> a
    ^b
    n a -> a -> a
    forall a. Num a => a -> a -> a
    - a
    qa -> b -> a
    forall a b. (Num a, Integral b) => a -> b -> a
    ^b
    i | b
    i <- [0..b
    nb -> b -> b
    forall a. Num a => a -> a -> a
    -1] ]
    -- for the first row, we can choose any vector except zero, hence q^n-1

    -- for the second row, we can choose any vector except a multiple of the first, hence q^n-q

    -- etc

    orderAff :: b -> a -> a
    orderAff n :: b
    n q :: a
    q = a
    qa -> b -> a
    forall a b. (Num a, Integral b) => a -> b -> a
    ^b
    n a -> a -> a
    forall a. Num a => a -> a -> a
    * b -> a -> a
    forall a b. (Num a, Integral b) => b -> a -> a
    orderGL b
    n a
    q

    orderPGL :: b -> a -> a
    orderPGL n :: b
    n q :: a
    q = b -> a -> a
    forall a b. (Num a, Integral b) => b -> a -> a
    orderGL b
    n a
    q a -> a -> a
    forall a. Integral a => a -> a -> a
    `div` (a
    qa -> a -> a
    forall a. Num a => a -> a -> a
    -1)

    -- NOTE:

    -- AG(n,F2) is degenerate:

    -- Every pair of points is a line, so it is the complete graph on 2^n points

    -- And as such has aut group S(2^n)

    -- Heawood graph = incidence graph of Fano plane

    heawood :: Graph Integer
    heawood = Graph (Either [F2] [[F2]]) -> Graph Integer
    forall t a. (Ord t, Num t, Enum t, Ord a) => Graph a -> Graph t
    to1n (Graph (Either [F2] [[F2]]) -> Graph Integer)
    -> Graph (Either [F2] [[F2]]) -> Graph Integer
    forall a b. (a -> b) -> a -> b
    $ Int -> [F2] -> Graph (Either [F2] [[F2]])
    forall a.
    (Num a, Ord a, FinSet a) =>
    Int -> [a] -> Graph (Either [a] [[a]])
    incidenceGraphPG 2 [F2]
    f2