diff --git a/HaskellForMaths.cabal b/HaskellForMaths.cabal
--- a/HaskellForMaths.cabal
+++ b/HaskellForMaths.cabal
@@ -1,7 +1,8 @@
    Name:                HaskellForMaths
-   Version:             0.1.7
+   Version:             0.1.8
    Category:            Math
    Description:         Math library - combinatorics, group theory, commutative algebra, non-commutative algebra
+   Synopsis:            Combinatorics, group theory, commutative algebra, non-commutative algebra
    License:             BSD3
    License-file:        license.txt
    Author:              David Amos
@@ -9,7 +10,17 @@
    Homepage:            http://www.polyomino.f2s.com/haskellformathsv2/HaskellForMathsv2.html
    Build-Type:          Simple
    Cabal-Version:       >=1.2
-   
+
+   Extra-source-files:
+        Math/Test/TCommutativeAlgebra.hs,
+        Math/Test/TDesign.hs,
+        Math/Test/TField.hs,
+        Math/Test/TFiniteGeometry.hs,
+        Math/Test/TGraph.hs,
+        Math/Test/TNonCommutativeAlgebra.hs,
+        Math/Test/TPermutationGroup.hs,
+        Math/Test/TestAll.hs
+
    Library
      Build-Depends:     base >=3 && < 4, containers
      Exposed-modules:
@@ -23,7 +34,9 @@
         Math.Combinatorics.Design, Math.Combinatorics.FiniteGeometry, Math.Combinatorics.Hypergraph,
         Math.Common.IntegerAsType, Math.Common.ListSet,
         Math.Projects.RootSystem,
+        Math.Projects.Rubik,
         Math.Projects.ChevalleyGroup.Classical, Math.Projects.ChevalleyGroup.Exceptional,
         Math.Projects.KnotTheory.Braid,
         Math.Projects.KnotTheory.LaurentMPoly, Math.Projects.KnotTheory.TemperleyLieb, Math.Projects.KnotTheory.IwahoriHecke
+
      ghc-options:       -w
diff --git a/Math/Algebra/Field/Extension.hs b/Math/Algebra/Field/Extension.hs
--- a/Math/Algebra/Field/Extension.hs
+++ b/Math/Algebra/Field/Extension.hs
@@ -4,6 +4,7 @@
 
 module Math.Algebra.Field.Extension where
 
+import Data.Ratio
 import Data.List as L (elemIndex)
 
 import Math.Common.IntegerAsType
@@ -108,7 +109,7 @@
     show (Ext (UP as)) = showUP "a" as
 
 instance (Num k, Fractional k, PolynomialAsType k poly) => Num (ExtensionField k poly) where
-    Ext x + Ext y = Ext $ (x+y) `modUP` pvalue (undefined :: (k,poly))
+    Ext x + Ext y = Ext $ (x+y) -- `modUP` pvalue (undefined :: (k,poly))
     Ext x * Ext y = Ext $ (x*y) `modUP` pvalue (undefined :: (k,poly))
     negate (Ext x) = Ext $ negate x
     fromInteger x = Ext $ fromInteger x
@@ -116,9 +117,15 @@
 instance (Num k, Fractional k, PolynomialAsType k poly) => Fractional (ExtensionField k poly) where
     recip 0 = error "ExtensionField.recip 0"
     recip (Ext f) = let g = pvalue (undefined :: (k,poly))
-                        (u,v,1) = extendedEuclidUP f g -- so u*f + v*g == 1. (We know the gcd is 1 as g is irreducible)
-                    in Ext $ u `modUP` g
+                        (u,v,d@(UP [c])) = extendedEuclidUP f g
+                        -- so u*f + v*g == d. We know the d is a unit, ie field element, since g is irreducible
+                    in Ext $ (c /> u) `modUP` g
+    fromRational q = fromInteger a / fromInteger b where a = numerator q; b = denominator q
 
+-- divide through
+c /> f@(UP as) | c == 1 = f
+               | c /= 0 = UP (map (c' *) as) where c' = recip c
+
 instance (FiniteField k, PolynomialAsType k poly) => FiniteField (ExtensionField k poly) where
     eltsFq _ = map Ext (polys (d-1) fp) where
         fp = eltsFq (undefined :: k)
@@ -133,7 +140,8 @@
 
 polys d fp = map toUPoly $ polys' d where
     polys' 0 = [[]]
-    polys' d = [x:xs | x <- fp, xs <- polys' (d-1)]
+    polys' d = [x:xs | x <- fp, xs <- polys' (d-1)] -- return in ascending order
+    -- polys' d = [x:xs | xs <- polys' (d-1), x <- fp] -- return with elts of fp first
 
 -- Conway polynomials from Holt, Handbook of Computational Group Theory, p60
 
diff --git a/Math/Algebra/LinearAlgebra.hs b/Math/Algebra/LinearAlgebra.hs
--- a/Math/Algebra/LinearAlgebra.hs
+++ b/Math/Algebra/LinearAlgebra.hs
@@ -2,6 +2,14 @@
 
 {-# OPTIONS_GHC -fglasgow-exts #-}
 
+-- |A module providing elementary operations involving scalars, vectors, and matrices
+-- over a ring or field. Vectors are represented as [a], matrices as [[a]].
+-- (No distinction is made between row and column vectors.)
+-- It is the caller's responsibility to ensure that the lists have the correct number of elements.
+--
+-- The mnemonic for many of the arithmetic operations is that the number of angle brackets
+-- on each side indicates the dimension of the argument on that side. For example,
+-- v <*>> m is multiplication of a vector on the left by a matrix on the right.
 module Math.Algebra.LinearAlgebra where
 
 import qualified Data.List as L
@@ -15,34 +23,54 @@
 
 -- The mnemonic for these operations is that the number of angle brackets on each side indicates the dimension of the argument on that side
 
+
+-- vector operations
+
+-- |u <+> v returns the sum u+v of vectors
+(<+>) :: (Num a) => [a] -> [a] -> [a]
 u <+> v = zipWith (+) u v
 
+-- |u <-> v returns the difference u-v of vectors
+(<->) :: (Num a) => [a] -> [a] -> [a]
 u <-> v = zipWith (-) u v
 
--- scalar multiplication
+-- |k *> v returns the product k*v of the scalar k and the vector v
+(*>) :: (Num a) => a -> [a] -> [a]
 k *> v = map (k*) v
 
-k *>> m = (map . map) (k*) m
-
--- dot product of vectors (also called inner or scalar product)
+-- |u <.> v returns the dot product of vectors (also called inner or scalar product)
+(<.>) :: (Num a) => [a] -> [a] -> a
 u <.> v = sum (zipWith (*) u v)
 
--- tensor product of vectors (also called outer or matrix product)
+-- |u <*> v returns the tensor product of vectors (also called outer or matrix product)
+(<*>) :: (Num a) => [a] -> [a] -> [[a]]
 u <*> v = [ [a*b | b <- v] | a <- u]
 
 
 -- matrix operations
 
+-- |a <<+>> b returns the sum a+b of matrices
+(<<+>>) :: (Num a) => [[a]] -> [[a]] -> [[a]]
 a <<+>> b = (zipWith . zipWith) (+) a b
 
+-- |a <<->> b returns the difference a-b of matrices
+(<<->>) :: (Num a) => [[a]] -> [[a]] -> [[a]]
 a <<->> b = (zipWith . zipWith) (-) a b
 
+-- |a <<*>> b returns the product a*b of matrices
+(<<*>>) :: (Num a) => [[a]] -> [[a]] -> [[a]]
 a <<*>> b = [ [u <.> v | v <- L.transpose b] | u <- a]
  
--- action on the left
+-- |k *> m returns the product k*m of the scalar k and the matrix m
+(*>>) :: (Num a) => a -> [[a]] -> [[a]]
+k *>> m = (map . map) (k*) m
+
+-- |m <<*> v is multiplication of a vector by a matrix on the left
+(<<*>) :: (Num a) => [[a]] -> [a] -> [a]
 m <<*> v = map (<.> v) m
 
--- action on the right
+-- |v <*>> m is multiplication of a vector by a matrix on the right
+(<*>>) :: (Num a) => [a] -> [[a]] -> [a]
 v <*>> m = map (v <.>) (L.transpose m)
 
 
@@ -58,8 +86,16 @@
     idMxs = map snd $ iterate next (0,[])
     next (j,m) = (j+1, (1 : replicate j 0) : map (0:) m)
 
+-- |iMx n is the n*n identity matrix
+iMx :: (Num t) => Int -> [[t]]
+iMx n = idMx n
+
+-- |jMx n is the n*n matrix of all 1s
+jMx :: (Num t) => Int -> [[t]]
 jMx n = replicate n (replicate n 1)
 
+-- |zMx n is the n*n matrix of all 0s
+zMx :: (Num t) => Int -> [[t]]
 zMx n = replicate n (replicate n 0)
 
 {-
@@ -100,6 +136,9 @@
 		Just a' -> M a'
 -}
 
+
+-- |The inverse of a matrix (over a field), if it exists
+inverse :: (Fractional a) => [[a]] -> Maybe [[a]]
 inverse m =
     let d = length m -- the dimension
         i = idMx d
@@ -140,6 +179,7 @@
          r:_ -> rowEchelonForm (((x:xs) <+> r) : rs)
 rowEchelonForm zs@([]:_) = zs
 
+reducedRowEchelonForm :: (Fractional a) => [[a]] -> [[a]]
 reducedRowEchelonForm m = reverse $ reduce $ reverse $ rowEchelonForm m where
     reduce (r:rs) = let r':rs' = reduceStep (r:rs) in r' : reduce rs' -- is this scanl or similar?
     reduce [] = []
@@ -168,6 +208,8 @@
 
 -- t (M m) = M (L.transpose m)
 
+-- |The determinant of a matrix (over a field)
+det :: (Fractional a) => [[a]] -> a
 det [[x]] = x
 det ((x:xs):rs) =
     if x /= 0
diff --git a/Math/Combinatorics/FiniteGeometry.hs b/Math/Combinatorics/FiniteGeometry.hs
--- a/Math/Combinatorics/FiniteGeometry.hs
+++ b/Math/Combinatorics/FiniteGeometry.hs
@@ -1,15 +1,19 @@
 -- Copyright (c) David Amos, 2008-2009. 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 Data.List as L
 import qualified Data.Set as S
 
 import Math.Algebra.Field.Base
-import Math.Algebra.Field.Extension hiding ( (<+>) )
+import Math.Algebra.Field.Extension hiding ( (<+>) ) -- , (*>) )
 import Math.Algebra.LinearAlgebra -- hiding ( det )
 
 import Math.Combinatorics.Graph
+import Math.Combinatorics.GraphAuts -- for use in GHCi
+import Math.Algebra.Group.PermutationGroup -- for use in GHCi
 
 -- |ptsAG n fq returns the points of the affine geometry AG(n,Fq), where fq are the elements of Fq
 ptsAG :: (FiniteField a) => Int -> [a] -> [[a]]
@@ -40,6 +44,9 @@
           k = length ps        -- the dimension of the flat
           fq = eltsFq undefined
 
+lineAG [p1,p2] = L.sort [ p1 <+> (c *> dp) | c <- fq ] where
+    dp = p2 <-> p1
+    fq = eltsFq undefined
 
 -- 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) )
@@ -125,9 +132,15 @@
 linesAG n fq = flatsAG n fq 1
 
 
--- less efficient but perhaps more intuitive
+-- 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 n fq = [ [x,y] | [x,y] <- combinationsOf 2 (ptsAG n fq),
+                          [x,y] == take 2 (lineAG [x,y]) ]
+
+
+-- 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)
-linesAG1 n fq = [ [x,z] | x <- ptsAG n fq, y <- ptsPG (n-1) fq,
+linesAG2 n fq = [ [x,z] | x <- ptsAG n fq, y <- ptsPG (n-1) fq,
                           z <- [x <+> y], [x,z] == take 2 (closureAG [x,z]) ]
 -- the point of the condition at the end is to avoid listing the same line more than once
 
@@ -176,4 +189,4 @@
 -- 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 as aut group S(2^n)
+-- And as such has aut group S(2^n)
diff --git a/Math/Combinatorics/Graph.hs b/Math/Combinatorics/Graph.hs
--- a/Math/Combinatorics/Graph.hs
+++ b/Math/Combinatorics/Graph.hs
@@ -1,5 +1,8 @@
 -- Copyright (c) David Amos, 2008. All rights reserved.
 
+-- |A module defining a polymorphic data type for (simple, undirected) graphs,
+-- together with constructions of some common families of graphs,
+-- new from old constructions, and calculation of simple properties of graphs.
 module Math.Combinatorics.Graph where
 
 import qualified Data.List as L
@@ -238,7 +241,8 @@
     [] -> -1 -- infinite
     p:ps -> length p - 1
 
--- diameter of a graph is maximum distance between two distinct vertices
+-- |The diameter of a graph is maximum distance between two distinct vertices
+diameter :: (Ord t) => Graph t -> Int
 diameter g@(G vs es)
     | isConnected g = maximum $ map maxDistance vs
     | otherwise = -1
@@ -252,8 +256,9 @@
     bfs ((z:zs) : nodes) = (z:zs) : bfs (nodes ++ [ w:z:zs | w <- nbrs g z, w `notElem` zs])
     bfs [] = []
 
--- girth of a graph is the size of the smallest cycle it contains
--- Note: If graph contains no cycles, we return -1, representing infinity
+-- |The girth of a graph is the size of the smallest cycle that it contains.
+-- Note: If the graph contains no cycles, we return -1, representing infinity.
+girth :: (Eq t) => Graph t -> Int
 girth g@(G vs es) = minimum' $ map minCycle vs where
     minimum' xs = let (zs,nzs) = L.partition (==0) xs in if null nzs then -1 else minimum nzs
     minCycle v = case findCycles g v of
@@ -289,6 +294,9 @@
 -- j v k i is isomorphic to j v (v-k) (v-2k+i), so may as well have v >= 2k
 
 -- kneser v k | v >= 2*k = j v k 0
+-- |kneser n k returns the kneser graph KG n,k -
+-- whose vertices are the k-element subsets of [1..n], with edges joining disjoint subsets
+kneser :: (Integral t) => t -> t -> Graph [t]
 kneser n k | 2*k <= n = graph (vs,es) where
     vs = combinationsOf k [1..n]
     es = [ [v1,v2] | [v1,v2] <- combinationsOf 2 vs, disjoint v1 v2]
diff --git a/Math/Common/ListSet.hs b/Math/Common/ListSet.hs
--- a/Math/Common/ListSet.hs
+++ b/Math/Common/ListSet.hs
@@ -1,5 +1,5 @@
--- Copyright (c) David Amos, 2008. All rights reserved.
 
+
 module Math.Common.ListSet where
 
 import Data.List (group,sort)
@@ -16,7 +16,8 @@
     LT -> x : union xs (y:ys)
     EQ -> x : union xs ys
     GT -> y : union (x:xs) ys
-union xs ys = xs ++ ys -- one of them is null
+union [] ys = ys
+union xs [] = xs
 
 intersect (x:xs) (y:ys) =
     case compare x y of
@@ -38,7 +39,8 @@
     LT -> x : symDiff xs (y:ys)
     EQ -> symDiff xs ys
     GT -> y : symDiff (x:xs) ys
-symDiff xs ys = xs ++ ys -- one of them is null
+symDiff [] ys = ys
+symDiff xs [] = xs
 
 disjoint xs ys = null (intersect xs ys)
 
diff --git a/Math/Projects/Rubik.hs b/Math/Projects/Rubik.hs
new file mode 100644
--- /dev/null
+++ b/Math/Projects/Rubik.hs
@@ -0,0 +1,26 @@
+-- Copyright (c) David Amos, 2009. All rights reserved.
+
+module Math.Projects.Rubik where
+
+import Math.Algebra.Group.PermutationGroup
+import Math.Algebra.Group.SchreierSims
+
+
+--           11 12 13
+--           14  U 16
+--           17 18 19
+-- 21 22 23   1  2  3  41 42 43  51 52 53
+-- 24  L 26   4  F  6  44  R 46  54  B 56
+-- 27 28 29   7  8  9  47 48 49  57 58 59
+--           31 32 33
+--           34  D 36
+--           37 38 39
+
+f = p [[ 1, 3, 9, 7],[ 2, 6, 8, 4],[17,41,33,29],[18,44,32,26],[19,47,31,23]]
+b = p [[51,53,59,57],[52,56,58,54],[11,27,39,43],[12,24,38,46],[13,21,37,49]]
+l = p [[21,23,29,27],[22,26,28,24],[ 1,31,59,11],[ 4,34,56,14],[ 7,37,53,17]]
+r = p [[41,43,49,47],[42,46,48,44],[ 3,13,57,33],[ 6,16,54,36],[ 9,19,51,39]]
+u = p [[11,13,19,17],[12,16,18,14],[ 1,21,51,41],[ 2,22,52,42],[ 3,23,53,43]]
+d = p [[31,33,39,37],[32,36,38,34],[ 7,47,57,27],[ 8,48,58,28],[ 9,49,59,29]]
+
+-- In Singmaster notation these would be capital letters.
diff --git a/Math/Test/TCommutativeAlgebra.hs b/Math/Test/TCommutativeAlgebra.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TCommutativeAlgebra.hs
@@ -0,0 +1,95 @@
+-- Copyright (c) David Amos, 2008. All rights reserved.
+
+{-# LANGUAGE FlexibleInstances #-}
+
+module Math.Test.TCommutativeAlgebra where
+
+import Math.Algebra.Field.Base
+import Math.Algebra.Commutative.Monomial
+import Math.Algebra.Commutative.MPoly
+import Math.Algebra.Commutative.GBasis
+
+import Test.QuickCheck
+
+-- > quickCheck prop_CommRingMPoly
+-- > verboseCheck prop_ComRingMPoly -- to see what input data is being used
+
+-- Commutative Ring (with 1)
+prop_CommRing (a,b,c) =
+    a+(b+c) == (a+b)+c   &&  -- addition is associative
+    a+b == b+a           &&  -- addition is commutative
+    a+0 == a             &&  -- additive identity
+    a+(-a) == 0          &&  -- additive inverse
+    a*(b*c) == (a*b)*c   &&  -- multiplication is associative
+    a*b == b*a           &&  -- multiplication is commutative
+    a*1 == a             &&  -- multiplicative identity
+    a*(b+c) == a*b + a*c     -- distributivity
+
+monomial is = product $ zipWith (^) (map x_ [1..]) (map (max 0) is)
+
+-- mpoly :: [(Integer,[Int])] -> MPoly Grevlex Q
+mpoly ais = sum [fromInteger a * monomial is | (a,is) <- ais]
+
+{-
+-- can take a long time to run, probably because of the test for associativity of multiplication
+prop_CommRingMPoly (ais,bjs,cks) = prop_CommRing (f,g,h) where
+    f = mpoly ais
+    g = mpoly bjs
+    h = mpoly cks
+    types = (ais,bjs,cks) :: ( [(Integer,[Int])], [(Integer,[Int])], [(Integer,[Int])] )
+-}
+
+instance Arbitrary (MPoly Grevlex Q) where
+    -- arbitrary = do ais <- arbitrary :: Gen [(Integer,[Int])]
+    arbitrary = do ais <- sized $ \n -> resize (n `div` 2) arbitrary :: Gen [(Integer,[Int])]
+                   return (mpoly ais)
+    coarbitrary = undefined -- !! only required if we want to test functions over the type
+
+prop_CommRingMPoly (f,g,h) = prop_CommRing (f,g,h) where
+    types = (f,g,h) :: (MPoly Grevlex Q, MPoly Grevlex Q, MPoly Grevlex Q)
+
+
+-- Sources for tests:
+-- [IVA] - Cox, Little, O'Shea: Ideals, Varieties and Algorithms
+-- [UAG] - Cox, Little, O'Shea: Using Algebraic Geometry
+
+
+test = and [
+    gb (map toGlex [x*z-y^2,x^3-z^2]) == map toGlex [y^6-z^5,x*y^4-z^4,x^2*y^2-z^3,x^3-z^2,x*z-y^2], -- IVA p93
+    gb (map toLex [x^2+y^2+z^2-1,x^2+z^2-y,x-z]) == map toLex [x-z,y-2*z^2,z^4+1/2*z^2-1/4], -- IVA p94
+    gb (map toLex [x^2+y^2+z^2-1,x*y*z-1]) == map toLex [x+y^3*z+y*z^3-y*z,y^4*z^2+y^2*z^4-y^2*z^2+1], -- IVA p116
+    gb [x*y+z-x*z,x^2-z,2*x^3-x^2*y*z-1] == [z^4-3*z^3-4*y*z+2*z^2-y+2*z-2,y*z^2+2*y*z-2*z^2+1,y^2-2*y*z+z^2-z,x+y-z] -- Grevlex, UAG p50-1
+    ]
+
+
+
+{-
+http://www.cs.amherst.edu/~dac/iva.html
+states that IVA, 2nd ed, 5th printing (the one I have) has a production error causing many +s and -s to appear incorrectly
+
+This explains the following misprints I've found:
+p117:
+gb (map toLex [x*y-4,y^2-(x^3-1)])
+-> [x-1/16y^4-1/16y^2,y^5+y^3-64]
+IVA p117 claims it should be -y^3 in the second poly
+But my answer is clearly correct, by looking at the reduction sequence for x*y-4
+x*y-4 -> 1/16(y^5+y^3)-4 -> 0
+  x-1/16(y^4+y^2)  y^5+y^3-64
+By contrast, reducing over their set clearly stops at 1/8y^3
+
+gb (map toLex [x-t-u,y-t^2-2*t*u,z-t^3-3*t^2*u])
+The answer I get has some sign differences compared to IVA p127
+-}
+
+{-
+The code has no trouble chomping through some of the examples that took a long time in the Sugar paper, eg
+gb [x+y+z+t+u, x*y+y*z+z*t+t*u+u*x, x*y*z+y*z*t+z*t*u+t*u*x+u*x*y, x*y*z*t+y*z*t*u+z*t*u*x+t*u*x*y+u*x*y*z, x*y*z*t*u-1]
+gb $ map toLex [x+y+z+t+u, x*y+y*z+z*t+t*u+u*x, x*y*z+y*z*t+z*t*u+t*u*x+u*x*y, x*y*z*t+y*z*t*u+z*t*u*x+t*u*x*y+u*x*y*z, x*y*z*t*u-1]
+gb [w^31-w^6-w-x, w^8-y, w^10-z]
+gb $ map toLex [w^31-w^6-w-x, w^8-y, w^10-z]
+
+However, for some reason, the code gets indigestion on the following
+gb $ map toLex [y*(1+x^2)^4 - 2*(5+19*x^2-45*x^4+x^6-4*x^8), z*(1+x^2)^4-2*(x+51*x^3+3*x^5+17*x^7)]
+
+(For comparison, the v1 implementation of gbasis can manage, even though its performance on the sugar examples is only comparable)
+-}
diff --git a/Math/Test/TDesign.hs b/Math/Test/TDesign.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TDesign.hs
@@ -0,0 +1,42 @@
+-- Copyright (c) David Amos, 2008. All rights reserved.
+
+{-# LANGUAGE FlexibleInstances #-}
+
+module Math.Test.TDesign where
+
+import qualified Data.List as L
+
+-- import PermutationGroup as P
+import Math.Algebra.Field.Base
+import Math.Algebra.Field.Extension
+-- import Graph as G
+-- import StronglyRegularGraph as SRG
+import Math.Combinatorics.Design as D
+import Math.Algebra.Group.SchreierSims as SS
+-- import LinearAlgebra hiding ( (^-) )
+
+
+factorial n = product [1..n]
+
+choose n m | m <= n = product [m+1..n] `div` product [1..n-m]
+
+test = and [designParamsTest, designAutTest]
+
+
+designParamsTest = and
+    [designParams (ag2 f2) == Just (2,(4,2,1))
+    ,designParams (ag2 f3) == Just (2,(9,3,1))
+    ,designParams (ag2 f4) == Just (2,(16,4,1))
+    ,designParams (pg2 f2) == Just (2,(7,3,1))
+    ,designParams (pg2 f3) == Just (2,(13,4,1))
+    ,designParams (pg2 f4) == Just (2,(21,5,1))
+    ]
+
+designAutTest = all (uncurry (==)) designAutTests
+
+designAutTests =
+    [(SS.order $ designAuts $ pg2 f2, 168) -- this is L3(2), see Atlas
+    ,(SS.order $ designAuts $ pg2 f3, 5616) -- this is L3(3)
+    ,(SS.order $ designAuts $ pg2 f4, 120960) -- this is S3.L3(4)
+--    ,(SS.order $ designAuts $ pg2 f5, 372000) -- this is L3(5)
+    ]
diff --git a/Math/Test/TField.hs b/Math/Test/TField.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TField.hs
@@ -0,0 +1,9 @@
+module Math.Test.TField where
+
+import Math.Algebra.Field.Base
+import Math.Algebra.Field.Extension
+
+test = and
+    [ (1/5 :: QSqrt3) * 5 == 1 -- regression test for defect
+    , (1/4 :: F25) * 4 == 1
+    ]
diff --git a/Math/Test/TFiniteGeometry.hs b/Math/Test/TFiniteGeometry.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TFiniteGeometry.hs
@@ -0,0 +1,29 @@
+module Math.Test.TFiniteGeometry where
+
+import Math.Combinatorics.FiniteGeometry
+import Math.Algebra.Field.Base
+import Math.Algebra.Field.Extension
+import Math.Combinatorics.GraphAuts
+import Math.Algebra.Group.PermutationGroup
+
+test = and
+    [numFlatsAG 2 2 0 == length (flatsAG 2 f2 0)
+    ,numFlatsAG 2 2 1 == length (flatsAG 2 f2 1)
+    ,numFlatsAG 2 2 2 == length (flatsAG 2 f2 2)
+    ,numFlatsAG 3 2 1 == length (flatsAG 3 f2 1)
+    ,numFlatsAG 3 3 1 == length (flatsAG 3 f3 1)
+    ,numFlatsAG 3 4 1 == length (flatsAG 3 f4 1)
+    ,numFlatsAG 3 4 2 == length (flatsAG 3 f4 2)
+    ,numFlatsAG 3 4 3 == length (flatsAG 3 f4 3)
+    ,numFlatsPG 2 2 0 == length (flatsPG 2 f2 0)
+    ,numFlatsPG 2 2 1 == length (flatsPG 2 f2 1)
+    ,numFlatsPG 2 2 2 == length (flatsPG 2 f2 2)
+    ,numFlatsPG 3 2 1 == length (flatsPG 3 f2 1)
+    ,numFlatsPG 3 3 1 == length (flatsPG 3 f3 1)
+    ,numFlatsPG 3 4 1 == length (flatsPG 3 f4 1)
+    ,numFlatsPG 3 4 2 == length (flatsPG 3 f4 2)
+    ,numFlatsPG 3 4 3 == length (flatsPG 3 f4 3)
+    ,(orderSGS $ incidenceAuts $ incidenceGraphAG 2 f2) == orderAff 2 2 * degree f2 
+    ,(orderSGS $ incidenceAuts $ incidenceGraphAG 2 f3) == orderAff 2 3 * degree f3 
+    ,(orderSGS $ incidenceAuts $ incidenceGraphAG 2 f4) == orderAff 2 4 * degree f4 
+    ]
diff --git a/Math/Test/TGraph.hs b/Math/Test/TGraph.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TGraph.hs
@@ -0,0 +1,92 @@
+-- Copyright (c) David Amos, 2008. All rights reserved.
+
+{-# LANGUAGE FlexibleInstances #-}
+
+module Math.Test.TGraph where
+
+import qualified Data.List as L
+
+import Math.Combinatorics.Graph as G
+import Math.Combinatorics.StronglyRegularGraph as SRG
+import Math.Combinatorics.Hypergraph as H
+import Math.Combinatorics.GraphAuts
+import Math.Algebra.Group.PermutationGroup as P -- not used
+import Math.Algebra.Group.SchreierSims as SS
+
+import Math.Algebra.Group.StringRewriting
+
+
+-- Sources
+-- [AGT] - Godsil and Royle, Algebraic Graph Theory
+
+factorial n = product [1..n]
+
+choose n m | m <= n = product [m+1..n] `div` product [1..n-m]
+
+test = and [graphPropsTest, graphTransitivityTest, srgParamTest, graphAutTest]
+
+
+graphPropsTest = all (uncurry (==)) graphPropsTestsBool
+              && all (uncurry (==)) graphPropsTestsInt
+
+graphPropsTestsBool =
+    [(isConnected nullGraph, True)] ++
+    [(isConnected (c n), True) | n <- [3..8] ] ++
+    [(isConnected $ complement $ k n, False) | n <- [3..6] ]
+
+graphPropsTestsInt =
+    [(diameter (c n), n `div` 2) | n <- [3..8] ] ++
+    [(girth (c n), n) | n <- [3..8] ] ++
+    [(girth (kb m n), 4) | m <- [2..4], n <- [2..4] ] ++
+    [(girth petersen, 5), (girth heawoodGraph, 6), (girth coxeterGraph, 7), (girth tutteCoxeterGraph, 8)]
+
+graphTransitivityTest = and graphTransitivityTests
+
+graphTransitivityTests =
+    [(not . isVertexTransitive) (kb m n) | n <- [1..3], m <- [1..3], m /= n] ++
+    [isEdgeTransitive (kb m n) | n <- [1..3], m <- [1..3]] ++
+    map isArcTransitive [k 4, kb 3 3, q 3, dodecahedron, G.to1n heawoodGraph, G.to1n coxeterGraph, G.to1n tutteCoxeterGraph] ++
+    map is2ArcTransitive [c 7, q 3, G.to1n coxeterGraph] ++
+    map is3ArcTransitive [c 7, G.to1n petersen] ++
+    map (not . is3ArcTransitive) [q 3] ++
+    [isArcTransitive (j v k i) | v <- [3..5], k <- [1..v `div` 2], i <- [0..k] ] ++ -- [AGT] p60
+    [is2ArcTransitive (j (2*k+1) k 0) | k <- [1..2] ] ++
+    [isDistanceTransitive (j v k (k-1)) | v <- [3..5], k <- [1..v `div` 2] ] ++ -- [AGT] p75
+    [isDistanceTransitive (j (2*k+1) k 0) | k <- [1..2] ] ++
+    [p doyleGraph | p <- [isVertexTransitive, isEdgeTransitive, not . isArcTransitive, not . isDistanceTransitive] ]
+
+-- Most of the graphs we construct are highly symmetric, and turn out to be arc- and distance-transitive
+-- On the other hand, those which aren't arc- or distance-transitive are often trivially not so,
+-- by virtue of not even being vertex- or edge-transitive
+-- It is actually rather hard to find graphs which are vertex- and edge-transitive but not arc-transitive, but here is one
+-- Doyle, "A 27-vertex graph that is vertex-transitive and edge-transitive but not 1-transitive"
+doyleGraph = G gs es where
+    relations = knuthBendix [("aaaaaaaaa",""), ("ccccccccc",""), ("aaaaaa","ccc"), ("cccccc","aaa"), ("ccccccccac","aaaa"), ("aaaaaaaaca","cccc")]
+    gs = L.sort $ nfs ("ac",relations) -- all elements of the group, reduced to normal form
+    hs = ["a","c","aaaaaaaa","cccccccc"] -- a, c, a^-1, c^-1
+    es = [ [v,v'] | v <- gs, v' <- L.sort [rewrite relations (v ++ h) | h <- hs], v < v']
+    -- so the edges join g to ga, gc, ga^-1, and gc^-1
+
+
+srgParamTest = all (uncurry (==)) srgParamTests
+
+-- van Lint & Wilson 262
+srgParamTests =
+    [(srgParams $ SRG.t m, Just (m `choose` 2, 2*(m-2), m-2, 4) ) | m <- [4..7] ]
+    ++ [(srgParams $ l2 m, Just (m^2, 2*(m-1), m-2, 2) ) | m <- [2..6] ]
+--     ++ [(srgParams $ paleyGraph fq, Just (q, (q-1) `div` 2, (q-5) `div` 4, (q-1) `div` 4) ) | (q,fq) <- [(5,f5),(9,f9),(13,f13),(17,f17)] ]
+    ++ [(srgParams $ G.petersen, Just (10,3,0,1) ) ]
+    ++ [(srgParams $ clebsch, Just (16,5,0,2) ) ]
+    ++ [(srgParams $ hoffmanSingleton, Just (50,7,0,1) ) ]
+    ++ [(srgParams $ higmanSimsGraph, Just (100,22,0,6) ) ]
+    ++ [(srgParams $ sp (2*r), Just (2^(2*r)-1,2^(2*r-1),2^(2*r-2),2^(2*r-2))) | r <- [2..3] ]
+
+graphAutTest = all (uncurry (==)) graphAutTests
+
+graphAutTests =
+    [(SS.order $ graphAuts $ c n, 2*n) | n <- [3..6] ] -- Aut(C n) = _D2 n
+    ++ [(SS.order $ graphAuts $ k n, factorial n) | n <- [3..6] ] -- Aut(K n) = S n
+    ++ [(SS.order $ graphAuts $ kb m n, factorial m * factorial n) | m <- [1..4], n <- [m+1..5] ] -- Aut(K m n) = S m * S n (m /= n)
+    ++ [(SS.order $ graphAuts $ kb n n, 2 * (factorial n)^2 ) | n <- [1..5] ] -- Aut(K n n) = S n * S n * C2 (m == n)
+    ++ [(SS.order $ graphAuts $ l2 n, 2 * (factorial n)^2 ) | n <- [2..5] ] -- Aut(L2 n) = S m * S m * C2
+
diff --git a/Math/Test/TNonCommutativeAlgebra.hs b/Math/Test/TNonCommutativeAlgebra.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TNonCommutativeAlgebra.hs
@@ -0,0 +1,39 @@
+-- Copyright (c) David Amos, 2008. All rights reserved.
+
+{-# LANGUAGE FlexibleInstances #-}
+
+module Math.Test.TNonCommutativeAlgebra where
+
+import Math.Algebra.Field.Base
+import Math.Algebra.NonCommutative.NCPoly
+import Math.Algebra.NonCommutative.TensorAlgebra
+
+import Test.QuickCheck
+
+-- > quickCheck prop_NonCommRingNPoly
+
+-- Non-Commutative Ring (with 1)
+prop_NonCommRing (a,b,c) =
+    a+(b+c) == (a+b)+c   &&  -- addition is associative
+    a+b == b+a           &&  -- addition is commutative
+    a+0 == a             &&  -- additive identity
+    a+(-a) == 0          &&  -- additive inverse
+    a*(b*c) == (a*b)*c   &&  -- multiplication is associative
+    a*1 == a && 1*a == a &&  -- multiplicative identity
+    a*(b+c) == a*b + a*c &&  -- left distributivity
+    (a+b)*c == a*c + b*c     -- left distributivity
+
+monomial is = product $ map (e_ . abs) is
+
+-- npoly :: [(Integer,[Int])] -> NPoly Q Basis
+npoly ais = sum [fromInteger a * monomial is | (a,is) <- ais]
+
+instance Arbitrary (NPoly Q Basis) where
+    -- arbitrary = do ais <- arbitrary :: Gen [(Integer,[Int])]
+    arbitrary = do ais <- sized $ \n -> resize (n `div` 2) arbitrary :: Gen [(Integer,[Int])]
+                   return (npoly ais)
+    coarbitrary = undefined -- !! only required if we want to test functions over the type
+
+prop_NonCommRingNPoly (f,g,h) = prop_NonCommRing (f,g,h) where
+    types = (f,g,h) :: (NPoly Q Basis, NPoly Q Basis, NPoly Q Basis)
+
diff --git a/Math/Test/TPermutationGroup.hs b/Math/Test/TPermutationGroup.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TPermutationGroup.hs
@@ -0,0 +1,119 @@
+-- Copyright (c) David Amos, 2008. All rights reserved.
+
+{-# LANGUAGE FlexibleInstances #-}
+
+module Math.Test.TPermutationGroup where
+
+import qualified Data.List as L
+
+import Math.Algebra.Group.PermutationGroup as P
+import Math.Algebra.Group.SchreierSims as SS
+import Math.Combinatorics.Graph
+import Math.Combinatorics.GraphAuts
+
+import Test.QuickCheck hiding (choose)
+
+factorials = scanl (*) 1 [1..] :: [Integer]
+
+-- factorial representation
+-- express n as a sum [ai * i! | i <- [1..]]
+facRep n = facRep' [] facs n where
+    facs = reverse $ takeWhile (<= n) $ tail factorials -- [i!, ..., 2!, 1!] where n >= i!
+    facRep' as (f:fs) n = let (q,r) = n `quotRem` f
+                          in facRep' (q : as) fs r
+    facRep' as [] 0     = as
+
+
+-- Unrank a permutation in Sn, using lexicographic order
+-- eg for S3, we have
+-- r    facRep r     unrankSn 3 r
+-- 0    0.1!+0.2!    [1,2,3]
+-- 1    1.1!+0.2!    [1,3,2]
+-- 2    0.1!+1.2!    [2,1,3]
+-- 3    1.1!+1.2!    [2,3,1]
+-- 4    0.1!+2.2!    [3,1,2]
+-- 5    1.1!+2.2!    [3,2,1]
+-- So the image of 1 is determined by the most significant digit of the facRep, and then recurse on the remaining
+unrankSn n r | r < factorial (toInteger n) =
+    let ds = reverse $ take (n-1) $ facRep r ++ repeat 0
+    in unrank' ds [1..n]
+    where unrank' (d:ds) xs = let x = xs !! fromIntegral d in x : unrank' ds (L.delete x xs)
+          unrank' [] [x] = [x]
+
+-- unrank permutations of S(N) (where N is the natural numbers from 1)
+-- doesn't use lexicographic order any more, but still a 1-1 mapping from N to permutations
+unrankSN r | r >= 0 = let ds = reverse (facRep r) in reverse (unrank' ds [1..length ds + 1]) where
+    unrank' (d:ds) xs = let x = if d==0 then last xs else xs !! (fromIntegral d-1)
+                        in x : unrank' ds (L.delete x xs)
+    unrank' [] [x] = [x]
+
+-- perm r = fromPairs $ zip [1..] $ unrankSN r
+
+instance Arbitrary (Permutation Int) where
+    arbitrary = do r <- arbitrary -- :: Gen Integer
+                   return (fromList $ unrankSN $ abs r)
+                   -- return (perm (abs r))
+                   -- return (perm (r^2)) -- to get some larger perms
+    coarbitrary = undefined
+
+
+prop_Group (g,h,k) =
+    g*(h*k)==(g*h)*k           && -- associativity
+    1*g == g && g*1 == g       && -- identity
+    g*g^-1 == 1 && g^-1*g == 1    -- inverse
+
+prop_GroupPerm (g,h,k) = prop_Group (g,h,k)
+    where types = (g,h,k) :: (Permutation Int, Permutation Int, Permutation Int)
+
+
+-- Could do more, like taking arbitrary lists of perms as generators of a group,
+-- and checking that the centre has the required property, etc
+
+
+factorial n = product [1..n]
+
+choose n m | m <= n = product [m+1..n] `div` product [1..n-m]
+
+test = and [sgsTest, ssTest, ccTest]
+
+
+sgsTest = all (uncurry (==)) sgsTests
+
+sgsTests =
+    [(sgsOrder $ _A n, SS.order $ _A n) | n <- [4..7] ] ++
+    [(sgsOrder $ _S n, SS.order $ _S n) | n <- [4..7] ] ++
+    [(sgsOrder $ _D2 n, SS.order $ _D2 n) | n <- [4..10] ] ++
+    [let _G = toSn (_S 3 `dp` _S 3) in (sgsOrder _G, SS.order _G) ] ++
+    [let _G = toSn (_C 3 `wr` _S 3) in (sgsOrder _G, SS.order _G) ] ++
+    [let _G = toSn (_S 3 `wr` _C 3) in (sgsOrder _G, SS.order _G) ]
+    where sgsOrder = orderTGS . tgsFromSgs . sgs
+
+
+ssTest = all (uncurry (==)) ssTests
+
+ssTests =
+    [(L.sort $ P.elts $ _C n, L.sort $ SS.elts $ _C n) | n <- [2..6] ]
+    ++ [(L.sort $ P.elts $ _D2 n, L.sort $ SS.elts $ _D2 n) | n <- [3..6] ]
+    ++ [(L.sort $ P.elts $ _S n, L.sort $ SS.elts $ _S n) | n <- [3..5] ]
+    ++ [(L.sort $ P.elts $ _A n, L.sort $ SS.elts $ _A n) | n <- [3..5] ]
+    ++ [let _G = toSn (_S 3 `dp` _S 3) in (L.sort $ P.elts _G, L.sort $ SS.elts _G) ]
+    ++ [let _G = toSn (_C 3 `wr` _S 3) in (L.sort $ P.elts _G, L.sort $ SS.elts _G) ]
+    ++ [let _G = toSn (_S 3 `wr` _C 3) in (L.sort $ P.elts _G, L.sort $ SS.elts _G) ]
+
+ccTest = and ccTests
+
+ccTests =
+    [conjClassReps (graphAuts2 (c 5)) == [(p [],1),(p [[1,2],[3,5]],5),(p [[1,2,3,4,5]],2),(p [[1,3,5,2,4]],2)]
+    ,conjClassReps (graphAuts2 (q 3)) ==
+        [(p [],1)
+        ,(p [[0,1],[2,3],[4,5],[6,7]],3)
+        ,(p [[0,1],[2,5],[3,4],[6,7]],6)
+        ,(p [[0,1,3,2],[4,5,7,6]],6)
+        ,(p [[0,1,3,7,6,4],[2,5]],8)
+        ,(p [[0,3],[1,2],[4,7],[5,6]],3)
+        ,(p [[0,3,6,5],[1,2,7,4]],6)
+        ,(p [[0,3,6],[1,7,4]],8)
+        ,(p [[0,3],[4,7]],6)
+        ,(p [[0,7],[1,6],[2,5],[3,4]],1)
+        ]
+    ]
diff --git a/Math/Test/TestAll.hs b/Math/Test/TestAll.hs
new file mode 100644
--- /dev/null
+++ b/Math/Test/TestAll.hs
@@ -0,0 +1,26 @@
+module Math.Test.TestAll where
+
+import Math.Test.TGraph
+import Math.Test.TDesign
+import Math.Test.TPermutationGroup
+import Math.Test.TFiniteGeometry
+import Math.Test.TCommutativeAlgebra
+import Math.Test.TNonCommutativeAlgebra
+import Math.Test.TField
+
+import Test.QuickCheck
+
+testall = and
+    [Math.Test.TGraph.test
+    ,Math.Test.TDesign.test
+    ,Math.Test.TPermutationGroup.test
+    ,Math.Test.TFiniteGeometry.test
+    ,Math.Test.TCommutativeAlgebra.test
+    ,Math.Test.TField.test
+    ]
+
+quickCheckAll =
+    do
+    quickCheck prop_CommRingMPoly
+    quickCheck prop_NonCommRingNPoly
+    quickCheck prop_GroupPerm