with(combinat);

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;

pet_cycleind_symm :=
proc(n)
local p, s;
option remember;

    if n=0 then return 1; fi;

    expand(1/n*add(a[l]*pet_cycleind_symm(n-l), l=1..n));
end;

V :=
proc(n)
    option remember;
    local comb, gf, var;

    comb := pet_varinto_cind(A+B+C+D, pet_cycleind_symm(n));
    gf := 
    expand(pet_varinto_cind(comb, pet_cycleind_symm(4)));

    for var in [A,B,C,D] do
        gf := coeff(gf, var, n);
    od;

    gf;
end;