functorial composition

[mantepse]

Provide an implementation of functorial composition.

This is the first implementation that can do the molecular expansion. The algorithm is entirely naive, except for using gap to the best of my knowledge. Thus, even rather small examples cannot be computed. I doubt that a significantly better algorithm exists, however.

Dependencies: #38544

[github-actions]

Documentation preview for this PR (built with commit a3db032; changes) is ready! 🎉
This preview will update shortly after each push to this PR.

[DaveWitteMorris]

It's great if GAP actions will solve the problem, but I don't know anything about them. Here is a sketch of a different method I think you could use to have sagemath and GAP do the calculation if the stabilizers are small, meaning that they each have only a few dozen subgroups.

A group G is acting transitively on a set X, with stabilizer G_x, and we have a subgroup H of G. For each subgroup H' of H, here is a way to calculate the number of points x in X, such that H' is the stabilizer of x in the subgroup H. (Dividing this number by |H/H'| then tells you the coefficient of the action H/H' in the restriction to H.)

For simplicity, let's first take the case where H' = H (so we are looking for the number of fixed points of H). Let n_0 be the number of subgroups of G_x that are conjugate to H'. (For each subgroup K of G_x, GAP can check whether H' is conjugate to K. If it is, then |N_G(H')| is the number of elements of G that conjugate H' to K, where N_G(H') is the normalizer of H' in G.) Then
the number of fixed points is n_0 |N_G(H')| / |G_x|.

The next step is to look at the case where H' is a maximal subgroup of H. Let n_1 be the number of subgroups of G_x that are conjugate to H'. The number of elements of G that conjugate H' into a subgroup of G_x is n_1 |N_G(H')|. However, some of these elements conjugate all of H into G_x, and therefore correspond to fixed points of the action of H, rather than orbits isomorphic to H/H'. Therefore,
the number of copies of H/H' is (n_1 - n_0) |N_G(H')| / (|G_x| |H/H'|)

Continuing in this way (essentially using inclusion-exclusion) to subtract off elements of G that give smaller orbits, sage can find the coefficient of H/H' for every subgroup H' of H. (Only take one representative of each conjugacy class of subgroups of H, though.)

[dimpase]

this seems like a huge overkill - what's unclear about GAP's actions? An action on a domain is specified by a function a(g,x) which takes a group element g and a domain element x, and computes g(x).

There are also few functions which actually don't adhere to this generic standard, but compute actions directly, without using a(g,x), e.g. the one for actions on cosets of a subgroup

[DaveWitteMorris]

I know almost nothing about GAP actions, so I don't know whether they can do the necessary calculations (such as finding the number of points that have a given stabilizer) for an action on a set with too many elements to list.

[dimpase]

finding the pointwise stabiliser of a subset of the domain in a permutation action is a standard thing in algorithms for permutation groups (involving "bases"- minimal wrt inclusion subsets B={B_1,..., B_k} fixed by every element, and "strong stabiliser chains" - subgroups G_j fixing B_1,...B_j, for all 0<j<k)

[DaveWitteMorris]

I don't think that's what we're looking for here. Basically, we're given a subgroup of the acting group (not a subset of the domain), and we want to know how many points it fixes. (More precisely, we want to know how many points have it as the exact stabilizer.) My understanding is that there are too many points (say 10^30) to actually find the fixed points, so we just want to know how many there are.

[DaveWitteMorris]

There are also few functions which actually don't adhere to this generic standard, but compute actions directly, without using a(g,x), e.g. the one for actions on cosets of a subgroup

An action on G/H is what we're looking at. For K_1 < K < G, we want to know how many points of G/H have K_1 as the stabilizer for the action of K on G/H.

[dimpase]

OK, then it's probably already in GAP - tables of marks: https://docs.gap-system.org/doc/ref/chap70.html

[DaveWitteMorris]

I think G might be too big to calculate the subgroup lattice, say S_20 or S_30, so we can't get its table of marks. But H and K are small enough (maybe S_4 or S_5) that we can compute their subgroup lattices. It would be great if GAP could compute the conjugacy information just for subgroups of H and K. (That's what I was trying to do by hand in my first answer.)

[dimpase]

On Thu, Jun 19, 2025 at 04:25:08PM -0700, DaveWitteMorris wrote: DaveWitteMorris left a comment (sagemath/sage#40186) I think G might be too big to calculate the subgroup lattice, say S_20 or S_30, so we can't get its table of marks. But H and K are small enough (maybe S_4 or S_5) that we can compute their subgroup lattices. It would be great if GAP could compute the conjugacy information just for subgroups of H and K. (That's what I was trying to do by hand in my first answer.)

GAP can compute the double cosets (H,K), i.e. {HgK | g in G}. Will it do?

[mantepse]

I don't think that's what we're looking for here. Basically, we're given a subgroup of the acting group (not a subset of the domain), and we want to know how many points it fixes. (More precisely, we want to know how many points have it as the exact stabilizer.) My understanding is that there are too many points (say 10^30) to actually find the fixed points, so we just want to know how many there are.

I am very interested, but I haven't digested your idea completely yet. Do you think you could walk me through? For example, following my notation from the sage-devel post:

• if F is the action whose stabilizer is the cyclic group generated by (1,2,3,4,5,6) and G is the action with trivial stabilizer group as a subgroup of S_3, we should get

  • 3 times the stabilizer group (1,2,3)
  • 8 times the stabilizer group generated by (1,2)
  • 15 times the trivial stabilizer.

• if F is the action with stabilizer generated by (1,2,3,4,5,6,7,8) and G is the action with one stabilizer subgroup being trivial as subgroup of S_3 and another generated by (1,2,3) we should get

  • 48 times the stabilizer group generated by (1,2)
  • 816 times the trivial stabilizer

• if F is the action with stabilizer generated by (1,2,...,24) and G is the action with trivial stabilizer as a subgroup of S_4, I cannot compute the stabilizers anymore. The size of the set S_4 is acting on is 25852016738884976640000 = 23!.

[DaveWitteMorris]

@dimpase: I think most of your suggestions would improve both the efficiency and the simplicity of the existing code, but the improvements would be somewhat incremental, because they will rely on making a list of the orbits, and then counting the number of orbits. (That is what the code currently does.) The OP wants a way to calculate the number of orbits (actually, the number of orbits with each possible stabilizer) when there are too many of them to make a list. E.g., he gives an example where the number of orbits is on the order of 23!. My comment above makes the simple observation that this counting is not difficult to do if there are only a few possibilities for the stabilizers. We do some small calculations, and then multiply by |N_G(H')| (which may be huge) to get the number of orbits of a given type.

@mantepse: I will try to respond soon with a more concrete explanation of how to use my suggestion to do the calculations for the examples you mentioned. It should not be difficult to automate in sagemath (using GAP), but may require some thought to get the inclusion-exclusion right.

[LudovicSchwob]

I have trouble when computing the functorial composition with cycles. When I do

L. = LazyCombinatorialSpecies(QQ)
E = L.Sets()
C = L.Cycles()
It takes forever to compute
E.functorial_composition(C)[6].

Also with

C.functorial_composition(C)[5]

I get
GAPError: Error, reached the pre-set memory limit
(change it with the -o command line option)

Small question: wouldn't it be possible to do lazy computing for the restriction of species ?
Here all coefficients are computed and it can take some time.

E.restrict(1,25)
X + E_2 + E_3 + E_4 + E_5 + E_6 + E_7 + E_8 + E_9 + E_10 + E_11 + E_12 + E_13 + E_14 + E_15 + E_16 + E_17 + E_18 + E_19 + E_20 + E_21 + E_22 + E_23 + E_24 + E_25

[mantepse]

Thank you for this comment!

Most of the time was spent in the renaming, so I made this optional - I will push this to #41569 because it has nothing to do with functorial composition.

However, the main problem remains.

I have trouble when computing the functorial composition with cycles. When I do

L. = LazyCombinatorialSpecies(QQ)
E = L.Sets()
C = L.Cycles()
It takes forever to compute
E.functorial_composition(C)[6].

Yes, this takes forever, and it would be really nice to have an algorithm that allows the computation if the output is small.

In the case at hand, what happens is as follows. Let us do it for n=8:

We call coefficient(8):

        def coefficient(n):
            S_n = SymmetricGroup(n)
            G_n = G[n]
            g_n = factorial(n) * G.generating_series()[n]

g_n is the number of G[n] structures, so 5040 = 7! in our case.

            result = R.zero()
            for f, c in left[g_n]:
...

Even though left is just the species of sets, this will take forever, because g_n is so big. The computation involved is to compute the decomposition of E_n into a direct product via PermutationGroup_generic.disjoint_direct_product_decomposition, i.e., a product of atomic species. In the case at hand, we know that E_n is atomic itself, so we could avoid it, but not in general. Almost all the time spent in this decomposition is used in Gap to obtain a stabilzer chain, so it is unlikely that we can improve it.

I have to check whether I can tell species easily that their n-th part is atomic.

[mantepse]

With #41569 we can do:

sage: L.Sets().functorial_composition(L.Cycles())
1 + X + E_2 + E_3 + E_4 + E_5 + E_6 + O^7

quickly.

[mantepse]

#41569 is actually quite a lot better than this one now, so maybe it is better to review that one instead.

[mantepse]

There is a hidden bug here. For example

sage: L.<X> = LazyCombinatorialSpecies(QQ)
sage: E = L.Sets()
sage: H = (E^2).functorial_composition(E^2); H.generating_series()
2 + 4*X + 8*X^2 + 128/3*X^3 + 8192/3*X^4 + 536870912/15*X^5 + 1152921504606846976/45*X^6 + O(X^7)
sage: H[2]
8*E_2 + 4*X^2
sage: H[1]
...
GAPError: Error, no method found! Error, no 1st choice method found for `FactorCosetAction' on 2 arguments

or, slightly simpler (the result should be 2X):

sage: (X^2).functorial_composition(2*X)
...
GAPError: Error, no method found! Error, no 1st choice method found for `FactorCosetAction' on 2 arguments

The problem is that l_G may be empty. The more important problem is that I forgot what my code does.

[mantepse]

@LudovicSchwob, I think it is correct now. See also #41569!

[mantepse]

Concerning performance, the following situation would be nice to improve. Consider

sage:  F.functorial_composition(G)

where, without loss of generality, F is molecular of degree n (i.e., corresponds to a transitive action of $\mathfrak S_n$) and produces only very few structures (i.e., the stabilizer subgroup has only very few cosets).

The extreme case would be the species of sets on n labels, in which case there is only a single structure, which would be particularly easy to handle. Currently, in this case, we are computing the action of a possibly huge symmetric group:

    f_act = libgap.FactorCosetAction(SymmetricGroup(g_count), f.permutation_group()[0])
    f_images = [libgap.Image(f_act, image) for image in images]

For example, if $E$ is the species of sets and suppose we want to compute E.functorial_composition(E^3). Of course this is silly for humans - we know that the result must be E itself. But in the two lines above, for the $n$-th term we have g_count = 3^n.

I wonder whether we can only do the special case of the species of sets (i.e., the full symmetric group as stabilizer subgroup), or whether we can somehow deal with slightly smaller stabilizer subgroups, too.

[mantepse]

Superseded by #41655.