Michi’s blog » Computational Group Theory in Haskell (1 in a series)

## Computational Group Theory in Haskell (1 in a series)

• October 18th, 2006
• 3:37 pm

This term, I’m listening to a lecture course on Computational Group Theory. As a good exercise, I plan to implement everything we talk about as Haskell functions as well.

The first lecture was mainly an introduction to the area, ending with a very naïve algorithm to generate a permutation group from a set of generators. Next week will bring less naïve algorithms with not quite as horrible complexity.

Before the algorithm can be brought, we’d want some undergrowth: we’d want to be able to work with permutations at all. So, we’ll start with the basic group theory and permutation implementations. A lot of this is stolen or rewritten from this permutation group code.

Our code will make use of two libraries, so if you collate code snippets while reading this, you’ll want to use

import Array
import List

If you don’t want to bother with that, the code is available here.

One of the things I love the most with Haskell is the way that you just write down what you think it means, and end up with good code. So, what is a group? It is this type class:

class Group a where
identity :: Int -> a
isIdentity :: a -> Bool
(*>) :: a -> a -> a     – right multiplication
inverse :: a -> a

with some conditions on the functions given. This is geared towards implementing permutation groups in the first line — so the idea is that you give your group as a subgroup of Sn, and you can retrieve the identity in Sn by using the identity function.

We tend to want to be able to conjugate elements and multiply a whole list of elements, so let’s add a few functions:

conjugate g h = inverse h *> g *> h
prod gs = foldl1 (*>) gs

And now we can start looking at permutations. We will store a permutation as the list of images, in order. So for instance, the permutation given by (1 2)(1 3) will be stored as the list [2,3,1]. (note that in the good tradition of computational group theory, we act from the right – thus motivating this result and not [2,1,3])

The interesting bit about permutations, though, is that they form a group; so we instantiate our type class as well.

newtype Permutation = PL [Int] deriving (Eq, Ord, Read, Show)
x .^ (PL g) = g !! (x-1)
instance Group Permutation where
identity n = PL [1..n]
isIdentity (PL xs) = all (\(i,j) -> i==j) (zip [1..] xs)
(PL g) *> h = PL (map (.^ h) g)
inverse (PL g) = PL (map snd (sort (zip g [1..])))

We are used to being able to write permutations in cycle notation, so we’d want to be able to convert to and from that, too.

fromCycles cs = PL (map snd (sort pointActions))
where
pointActions = concat (map fromCycle cs)                                                                                                                    fromCycle is = zip is (rotateL is)                                                                                                                              fromCycles’ n cs = PL (elems (array (1,n) [(i,i) | i <- [1..n]] // (concat (map fromCycle cs))))
findOrbit (PL permutation) elt =
let buildOrbit os ps e =
if (elem e os)
then os
else let ei = elemIndex e ps in
case ei of
Nothing -> os
Just e’ -> buildOrbit ([e] ++ os) ps (e’+1)
in buildOrbit [] permutation elt                                                                                                                              toCycles (PL permutation) =
let
getOrbits out [] per = out
getOrbits out els per =
let
orbit = findOrbit per (head els)
els’ = els \\ orbit
in getOrbits (out ++ [orbit]) els’ per
in getOrbits [] permutation (PL permutation)
rotateL (x:xs) = xs ++ [x]
saturateCycles :: [[Int]] -> Int -> [[Int]]                                     saturateCycles cs n = cs ++ ( group ( [1..n] \\ foldl1 (++) cs ) )
desaturateCycles cs = filter (\x -> length x > 1) cs

We can use saturateCycles and desaturateCycles in order to convert to and from the normal way of writing permutations. The toCycles and fromCycles functions both expect fixpoints to be listed as well. We need, however, to tell the saturateCycles which symmetric group we want to live in.

Now, the most naïve algorithm to list all elements in a group starts with the one element we know will live in the group: the identity. And for each iteration, all generators are multiplied onto the first known candidate in a list of group element candidates, and all results that are not yet known are added to this list.

generateGroup generators =
i = g0 *> inverse g0
addElts group [] gens = group
let cand = take 1 cands
group’ = group ++ cand
newcands = map (*> (head cand)) gens
cands’ = (tail cands) ++ (newcands \\ (group’ ++ cands))
in

We keep on adding tested candidates to the first list, and anything not yet seen that we get from our multiplication with generators gets added to the middle list. Once the middle list is empty, we’re done.

### 4 People had this to say...

Are you sure you want `identity :: n -> a` rather than `identity :: a`

• Michi
• October 19th, 2006
• 14:56

Actually, I do, since I view every group as a subgroup of some symmetric group. The identity is a function taking as argument which symmetric group I operate in, and returning the appropriate identity element.

Sure, translating group axioms straight off would make `identity :: a` be much saner; but on the other hand, I write the group type class with a specific use in mind.

Hang on, I’m sure there was meant to be more of that comment, to the effect that what you really want is for each symmetric group to be a type in its own right. But then you’d need types parametrized by integers, which you can do in C++ or Ada, but not (IIRC) in Haskell. I remember running into this problem a while ago while implementing the algorithms from my Computational Methods in Finite Fields course.

• Michi
• October 19th, 2006
• 15:38

Parametrizing the type “Subgroup of a symmetric group” would make a lot of things much easier; handling the cycle representation would look different, as would the multiplication code. But I’m not very much inclined to start writing things like
`newtype Permutation3 = ...`, and Haskell does not, at the current stage, seem to allow parametrized types.

There may possibly be some extension floating around in an alpha version capable of it though…