Introduction to GAP

This chapter provides an introduction to working with permutation groups. This will building towards Stabilizer Chains, the most important data structure when computing with permutation groups.

Throughout this book, we will assume we are working with finite permutation groups in GAP. Before reading this post, please go and install GAP so you can play along!


Let’s begin with a simple permutation group:

G := Group((1,3,2)(5,11,6,9,7,10)(8,12)(13,15,14),
Group([ (1,3,2)(5,11,6,9,7,10)(8,12)(13,15,14), (1,9,5)(2,11,6,3,10,7)(4,12,8)  (14,15) ])

There are many obvious questions we might want to answer about G, here are some simple examples:

What is the orbit of 1 in G?

Orbit(G, 1);
[ 1, 3, 9, 2, 10, 7, 5, 11, 6 ]

What is the orbit of 4 in G?

Orbit(G, 4);
[ 4, 12, 8 ]

Does G contain a particular permutation?

(5,9)(6,10)(7,11)(8,12) in G;

What is the size of G?


What is the subgroup of G which fixes 4?

Stabilizer(G, 4);
Group([ (2,3)(6,7)(10,11)(14,15), (5,9)(6,10)(7,11)(8,12), (1,2,3)(5,6,7)  (9,10,11)(13,14,15) ])

How can we calculate these values efficently? All of these could be done fairly easily if we checked every element of G in turn, to see where it maps 1, or if it fizes 4. However, there are many groups we are interested in studying which cannot be completely enumerated. Therefore we need to find ways of solving group problems without complete enumeration of the group.

We will start with some easier calculations, and build our way towards the most difficult to solve problems.

Orbit Calculation

We begin by calculating the orbit of a pont. This is the easiest thing to calculate, and a vital building block towards a structure known as a stabilizer chain, which help us solve many important problems.

The standard algorithm for calculating the orbit of an integer uses the idea that all elements of a (finite) group can be made by combining the the generators. Therefore if we want to know all the possible values of \(\{a^g | g \in G\}\), we can instead calculate all the places we can get to using a generating set \(s\) where \(\langle s \rangle = G\).

Except, if we calculate \(\{a^s | s \in S\}\), we can miss some points in the orbit, because we might need to apply elements multiple times. I think this is easier to understand in code, so lets write some GAP! This function is a little inefficient, but will be OK as a first attempt.

CalculateOrbitSlow := function(G, point)
    local knownOrbit, p, g, gens;
    gens := GeneratorsOfGroup(G);

1    knownOrbit := [point];
    for p in knownOrbit do
        for g in gens do
            if not(p^g in knownOrbit) then
                Add(knownOrbit, p^g);
    return knownOrbit;
Start knownOrbit with the only point we know is in the orbit – point
function( G, point ) ... end

Let’s write a couple of simple tests, to convince ourselves our function works.

CalculateOrbitSlow(G, 1);
# [ 1, 3, 9, 2, 10, 7, 5, 11, 6 ]
CalculateOrbitSlow(G, 4);
# [ 4, 12, 8 ]

This function has two limitations: Firstly, it is a inefficient because the line p^g in knownOrbit is slow – it has to search through the whole orbit we have seen so far. Secondly, in practice we often want to know the permutation which will get us to each point in the orbit, but we throw that information away!

So, let’s do an improved version. Here we will start with a base value, and then track whenever we find a new value in the orbit, which permutation got us there. By following the permutations we will be able to get back to the base.

Enough chat, I find this data structure is easier to understand once you have it:

CalculateSchreier := function(G, point)
    local knownOrbit, vec, img, p, g, gens;
    gens := GeneratorsOfGroup(G);

    vec := [];
    knownOrbit := [point];
    vec[point] := ();

    for p in knownOrbit do
        for g in gens do
            img := p/g; # Note this line!
            if not(IsBound(vec[img])) then
                vec[img] := g;
                Add(knownOrbit, img);

    # Return everything we found
    return rec(generators := gens,
               orbit := knownOrbit,
               transveral := vec);

Notice that line img := p/g. You might not have seen /. This returns the value p maps to g. It’s the same as applying g to the inverse of p, but is more efficient.

Why do we do this? Well, we want to be able to get back to the base. We could instead store the inverse of p, but that would require inverting the permutations.

So, what is this vector useful for? It’s main use is to let us find a permutation which maps any integer back to our ``base’’ point. I would play with this function for a while, to be sure you really understand why it works!

RepresentativePerm := function(schvec, val)
    local ret, gen;
    if not(IsBound(schvec.transveral[val])) then
        return fail;

    ret := ();

    while val <> schvec.orbit[1] do
        gen := schvec.transveral[val];
        val := val^gen;
        ret := ret*gen;

    return ret;

Let’s see how this function works in practice. It gives a permutation which maps it’s first argument back to the base point, or returns fail if no such permutation exists.

sch := CalculateSchreier(G, 1);;
RepresentativePerm(sch, 2);
# (1,3,2)(5,11,6,9,7,10)(8,12)(13,15,14)
RepresentativePerm(sch, 6);
# (1,10,5,2,9,6)(3,11,7)(4,12,8)(13,14)
RepresentativePerm(sch, 4);
# fail

Given these two together, we can take a permutation \(p\) and a group \(G\) and find another permutation \(q\), where \(q\) is in \(G\) if and only if \(p\) is, and \(q\) fixes the base point of our schreier vector. This function looks like this:

MapToBase := function(schvec, perm)
    local map;
    map := RepresentativePerm(schvec, schvec.orbit[1]^perm);
    if map = fail then return fail; fi;
    return perm * map;

This function takes a permutation and multiplies it by a permutation in our group, to fix the base point (or returns fail if no such permutation exists).

sch := CalculateSchreier(G, 1);;
MapToBase(sch, (1,2,3,4));
# (3,4)(5,11,6,9,7,10)(8,12)(13,15,14)
MapToBase(sch, (1,4));
# fail

So, we have turned the problem of finding if a permutation is in our group to the problem of finding if a different permutation is in the stabilizer of a single point in that group.

So how do we solve this problem? Well, build Stabilizer(G,1), build a schreier vector for that, and recurse until no group is left!

StabilizerChain := function(G)
    local root, pnt, Gstab;
    pnt := SmallestMovedPoint(G);
    root := CalculateSchreier(G, pnt);
    Gstab := Stabilizer(G, pnt);
    if not IsTrivial(Gstab) then
        root.stabilizer := StabilizerChain(Gstab);
    return root;

Now we have a stabilizer chain, we are (finally) in a position to implement checking if a permutation is in a group. This function is quite long, but follows a basic recursive design:

  • Look where the permutation (p) maps the base point. Is that outside the orbit of the base point? Then it is not in our group (G).
  • Calculate a permutation (q) in (G) such that (p.q) maps the base point to itself, then search in the stabilizer.

The final terminating case is to detect when we reach the bottom of the stabilizer chain, in which case we must have the identity permutation left.

PermInGroup := function(chain, perm)
    local basemap;
    basemap := MapToBase(chain, perm);
    if basemap = fail then return false; fi;
    if not IsBound(chain.stabilizer) then
        return basemap = ();
    return PermInGroup(chain.stabilizer, basemap);

Let’s test this new function:

chain := StabilizerChain(G);;
PermInGroup(chain, (1,5)(2,6)(3,7)(4,8));
# true
PermInGroup(chain, (1,5)(2,6)(3,7)(4,9));
# false

So, what else can we do with a stabilizer chain? Lots of things, almost anything which explores a group. Here is one easy thing we can do – find the size of the group using the orbit - stabilizer theorem.

GroupSize := function(stabchain)
    if not IsBound(stabchain.stabilizer) then
        return Length(stabchain.orbit);
    return Length(stabchain.orbit) *

# 36

We can do much more interesting calculations, like begin to build a group intersection algorithm.