1gap> # some random group
2gap> G := DirectProduct(SymmetricGroup(3), SymmetricGroup(3));;
3gap> irreps := IrreducibleRepresentations(G);;
4gap> # rho is irreps2 oplus irreps2 oplus irreps2 (rho has degree 3)
5gap> rho := DirectSumOfRepresentations([irreps[2], irreps[2], irreps[2]]);;
6gap> # so the canonical summand is just the whole space
7gap> summand := Cyclotomics^3;;
8gap> # and each axis is an irreducible G-invariant space
9gap> DecomposeCanonicalSummand@RepnDecomp(rho, irreps[2], summand); # should get 3 bits, 1 for each irrep
10[ rec( basis := [ [ 1, 0, 0 ] ] ),
11  rec( basis := [ [ 0, 1, 0 ] ] ),
12  rec( basis := [ [ 0, 0, 1 ] ] ) ]
13gap> # should also work if we use small-degree kronecker trick (same up to scaling)
14gap> DecomposeCanonicalSummand@RepnDecomp(rho, irreps[2], summand : use_kronecker);
15[ rec( basis := [ [ 36, 0, 0 ] ] ),
16  rec( basis := [ [ 0, 36, 0 ] ] ),
17  rec( basis := [ [ 0, 0, 36 ] ] ) ]