Cube stack problem How many distinct patterns are possible if you omit (a) 1 piece, (b) 2 pieces and (c) 3 pieces from a cube originally consisting of 27 smaller and equally sized cubes? 
 A: For those who are interested here is a version of the Maple program that has hard-coded vertices, edges and faces. It is not portable to higher dimension hypercubes but as the code is reduced considerably and readability profits, with the focus being on the essential.


with(numtheory);
with(group):
with(combinat):

pet_autom2cycles :=
proc(src, aut)
        local numa, numsubs;
        local marks, pos, cycs, cpos, clen;

        numsubs := [seq(src[k]=k, k=1..nops(src))];
        numa := subs(numsubs, aut);

        marks := [seq(true, pos=1..nops(aut))];

        cycs := []; pos := 1;

        while pos <= nops(aut) do
              if marks[pos] then
                 clen := 0; cpos := pos;

                 while marks[cpos] do
                       marks[cpos] := false;
                       cpos := numa[cpos];
                       clen := clen+1;
                 od;

                 cycs := [op(cycs), clen];
              fi;

              pos := pos+1;
        od;

        return mul(a[cycs[k]], k=1..nops(cycs));
end;


pet_varinto_cind :=
proc(poly, ind)
           local subs1, subs2, polyvars, indvars, v, pot, res;

           res := ind;

           polyvars := indets(poly);
           indvars := indets(ind);

           for v in indvars do
               pot := op(1, v);

               subs1 :=
               [seq(polyvars[k]=polyvars[k]^pot,
               k=1..nops(polyvars))];

               subs2 := [v=subs(subs1, poly)];

               res := subs(subs2, res);
           od;

           res;
end;

cube_vsyms :=
proc()
         option remember;
         local aut, adj, p, bitstr, bits, b, vsyms, vsym, flip,
         doflip, doperm, bits2ind, v1, v2, idx, idxterm,
         corners, edgemids, centers;

         bitstr := [];
         for b from 0 to 7 do
             bits := convert(8+b, base, 2);
             bitstr := [op(bitstr), [seq(bits[k], k=1..3)]];
         od;

         bits2ind := b -> 1+b[1]+2*b[2]+4*b[3];

         vsyms := [];

         for flip from 0 to 7 do
             bits := convert(8+flip, base, 2);

             doflip :=
             proc(v)
                 local res, bpos;
                 res := [];

                 for bpos to 3 do
                     if bits[bpos] = 1 then
                         res := [op(res), 1-v[bpos]];
                     else
                         res := [op(res), v[bpos]];
                     fi;
                 od;
             end;

             for p in permute(3) do
                 doperm :=
                 proc(v)
                     return [v[p[1]], v[p[2]], v[p[3]]];
                 end;

                 vsym := map(doflip, bitstr);
                 vsym := map(doperm, vsym);

                 vsyms := [op(vsyms), map(bits2ind, vsym)];
             od;
         od;

         corners := [1, 2, 3, 4, 5, 6, 7, 8];
         edgemids := [{1, 2}, {1, 3}, {1, 5}, {2, 4}, {2, 6},
                      {3, 4}, {3, 7}, {4, 8}, {5, 6}, {5, 7},
                      {6, 8}, {7, 8}];
         centers := [{1, 3, 5, 7}, {2, 4, 6, 8}, {1, 2, 5, 6},
                     {3, 4, 7, 8}, {1, 2, 3, 4}, {5, 6, 7, 8}];

         idx := 0;

         for aut in vsyms do
             idxterm := 1;

             p := [seq(k=aut[k], k=1..nops(aut))];

             idxterm :=
             idxterm *
             pet_autom2cycles(corners, subs(p, corners));

             idxterm :=
             idxterm *
             pet_autom2cycles(edgemids, subs(p, edgemids));

             idxterm :=
             idxterm *
             pet_autom2cycles(centers, subs(p, centers));

             idx := idx+idxterm*a[1];
         od;

         idx/nops(vsyms);
end;


v :=
proc(n)
         option remember;
         local gf, rep;

         rep := add(cat('X', k), k=1..n);
         gf := expand(pet_varinto_cind(rep, cube_vsyms()));

         subs([seq(cat('X',k)=1, k=1..n)], gf);
end;

Here is another interesting MSE cycle index computation I. This MSE cycle index computation II could be qualified as rather exotic.
