Introductory remark. The following discussion uses two different
definition of spanning subgraphs of $K_{n,m}$, one being subgraphs
with the same vertex set and the second subgraphs with the same vertex
set where there are no isolated vertices. The first of these is
equivalent to coloring the graph with two colors. Call these two model
$Q$ and model $P$ respectively.
The goal here is to enumerate spanning subgraphs of $K_{n,m}$ where we
will treat the simple case where even if $n=m$ there are no flips
above a central vertical axis i.e. no reflections. We can do much
better than NAUTY as we are only counting these graphs as opposed to
enumerating them. We use the Polya Enumeration Theorem (PET) to obtain
the count of all non-isomorphic subgraphs of $K_{n,m}$ (model $Q$) and
the principle of inclusion-exclusion (PIE) to extract the count of the
spanning subgraphs (model $P$).
We will consult NAUTY just the same to get sample data to match
against in our mathematical analysis. The following Perl script was
used.
#! /usr/bin/perl -w
#
sub decode_graph {
my ($str) = @_;
sub R {
my (@args) = map {
sprintf "%06b", $_;
} @_;
join '', @args;
}
my (@ents) = map {
ord($_) - 63
} split //, $str;
my $n = shift @ents;
my @adj_data = split //, R(@ents);
my $adj = []; my $pos = 0;
for(my $ind2 = 1; $ind2 < $n; $ind2++){
for(my $ind1 = 0; $ind1 < $ind2; $ind1++){
$adj->[$ind1]->[$ind2] = $adj_data[$pos];
$adj->[$ind2]->[$ind1] = $adj_data[$pos];
$pos++;
}
}
return $adj;
}
MAIN: {
my $mx = shift || 3;
die "out of range for GENBG: $mx"
if 2*$mx < 2 || 2*$mx > 32;
for(my $comp_a = 1; $comp_a <= $mx; $comp_a++){
for(my $comp_b = 1; $comp_b <= $mx; $comp_b++){
my $vcount = $comp_a + $comp_b;
my $cmd = sprintf "./genbg %d %d",
$comp_a, $comp_b;
my $count = 0;
open GENBG, "$cmd 2>/dev/null|";
while(my $bp = <GENBG>){
chomp $bp; my $adj = decode_graph $bp;
for($v = 0; $v < $vcount; $v++){
my $deg = 0;
for(my $w = 0; $w < $vcount; $w++){
my $ent = $adj->[$v]->[$w];
$deg++ if defined($ent) && $ent == 1;
}
last if $deg == 0;
}
$count++ if $v == $vcount;
}
close GENBG;
printf " " if $comp_b > 1;
printf "%06d", $count;
}
printf "\n";
};
}
This gave the following table:
$ ./scripts/bipartite.pl 6
000001 000001 000001 000001 000001 000001
000001 000003 000005 000008 000011 000015
000001 000005 000017 000042 000091 000180
000001 000008 000042 000179 000633 002001
000001 000011 000091 000633 003835 020755
000001 000015 000180 002001 020755 200082
Now for the mathematics. We use the Polya Enumeration Theorem as
conjectured by the OP. To do this we need the cycle index of the
action on the edges of the group that permutes the vertices in
partition $A$ of size $n$ according to the symmetric group $S_n$ and
the vertices in partition $B$ of size $m$ according to $S_m.$
These cycle indices are easy to compute and we do not need to iterate
over all $n!\times m!$ pairs of permutations (acting on $A$ and $B$)
but instead it is sufficient to iterate over pairs of terms from the
cycle indices $Z(S_n)$ and $Z(S_m)$ of the symmetric groups $S_n$ and
$S_m$ according to their multiplicities to obtain the cycle index
$Z(Q_{n,m})$ of the combined action on $A$ and $B$. The number of
terms here is the much better count of the number of partitions of $n$
and $m$ (upper bound).
The classic approach to the calculation of these cycle indices is
based on the simple observation that for two cycles, one of length
$l_1$ from a permutation $\alpha$ of $A$ and another of length $l_2$
from a column permutation $\beta$ of $B$ their contribution to the
disjoint cycle decomposition product for $(\alpha,\beta)$ in the cycle
index $Z(Q_{n,m})$ is by inspection
$$a_{\mathrm{lcm}(l_1, l_2)}^{l_1 l_2 / \mathrm{lcm}(l_1, l_2)}
= a_{\mathrm{lcm}(l_1, l_2)}^{\gcd(l_1, l_2)}.$$
Once we have the cycle indices we evaluate
$$Z(Q_{n,m})(1+z)$$
which is the standard substitution to produce the OGF. If we are only
looking to obtain the count, we may use
$$Z(Q_{n,m})_{\Large a_1=a_2=a_3=\cdots=2}.$$
Here is an example:
$$Z(Q_{3,4}) =
{\frac {{a_{{1}}}^{12}}{144}}+1/24\,{a_{{2}}}^{3}{a_{{1}}}^{6}+1/18\,{a_{{3
}}}^{3}{a_{{1}}}^{3}+1/12\,{a_{{2}}}^{6}+1/6\,{a_{{4}}}^{3}
\\+1/48\,{a_{{2}}}
^{4}{a_{{1}}}^{4}+1/8\,{a_{{1}}}^{2}{a_{{2}}}^{5}+1/6\,a_{{1}}a_{{2}}a_{{3}
}a_{{6}}+1/8\,{a_{{3}}}^{4}
\\+1/12\,{a_{{3}}}^{2}a_{{6}}+1/24\,{a_{{6}}}^{2}+
1/12\,a_{{12}}.$$
The substituted cycle index becomes
$$Z(Q_{3,4})(1+z) =
{z}^{12}+{z}^{11}+3\,{z}^{10}+6\,{z}^{9}+11\,{z}^{8}
\\ +13\,{z}^{7}+17\,{z}^{6
}+13\,{z}^{5}+11\,{z}^{4}+6\,{z}^{3}+3\,{z}^{2}+z+1.$$
At this point we have just about everything we need, the only problem
as is evident from the substituted cycle index (indexed by edge count)
is that we are counting all subgraphs including those that obviously
cannot span $K_{3,4}.$ Observe however that $Z(Q_{3,4})$ includes the
count from all graphs $K_{a,b}$ where $1\le a \le 3$ and $1\le b \le
4$ and this observation holds for the $Z(Q_{a,b})$ as well. The way to
identify a spanning subgraph of $K_{3,4}$ is that every vertex in the
vertex set has degree at least one, which means these are just the
graphs that cannot possibly be counted by $Z(Q_{a,b})$ with $(a,b)\ne
(3,4)$ because of the missing vertices. Therefore we apply PIE to the
poset where the nodes corresponding to the $(a,b)$ is the set of
graphs counted by the corresponding substituted cycle index. We have
by inspection that this poset is isomorphic to the divisor poset of
$2^{n-1}\times 3^{m-1}$ so that we may use the ordinary Möbius
function from number theory as our Möbius function for the PIE
computation. (We are not including the empty graph in the sets at the
nodes from the poset.)
This is implemented in the following Maple program.
with(numtheory);
pet_cycleind_symm :=
proc(n)
option remember;
if n=0 then return 1; fi;
expand(1/n*add(a[l]*pet_cycleind_symm(n-l), l=1..n));
end;
pet_varinto_cind :=
proc(poly, ind)
local subs1, subsl, polyvars, indvars, v, pot;
polyvars := indets(poly);
indvars := indets(ind);
subsl := [];
for v in indvars do
pot := op(1, v);
subs1 :=
[seq(polyvars[k]=polyvars[k]^pot,
k=1..nops(polyvars))];
subsl := [op(subsl), v=subs(subs1, poly)];
od;
subs(subsl, ind);
end;
pet_cycleind_knm :=
proc(n, m)
option remember;
local cind, sind1, sind2, t1, t2, q,
v1, v2, len, len1, len2;
cind := 0;
if n=1 then
sind1 := [a[1]];
else
sind1 := pet_cycleind_symm(n);
fi;
if m=1 then
sind2 := [a[1]];
else
sind2 := pet_cycleind_symm(m);
fi;
for t1 in sind1 do
for t2 in sind2 do
q := 1;
for v1 in indets(t1) do
len1 := op(1, v1);
for v2 in indets(t2) do
len2 := op(1, v2);
len := lcm(len1, len2);
q := q *
a[len]^((len1*len2/len) *
degree(t1, v1)*degree(t2, v2));
od;
od;
cind := cind +
lcoeff(t1)*lcoeff(t2)*q;
od;
od;
cind;
end;
v_pre_pie :=
proc(n, m)
option remember;
local cind;
cind := pet_cycleind_knm(n, m);
subs([seq(a[v]=2, v=1..n*m)], cind);
end;
v :=
proc(n, m)
local q, a, b, res;
q := 2^(n-1)*3^(m-1);
res := 0;
for a to n do
for b to m do
res := res +
mobius(q/2^(a-1)/3^(b-1))*
(v_pre_pie(a, b)-1);
od;
od;
res;
end;
print_table :=
proc(mx)
local n, m;
for n to mx do
for m to mx do
if m>1 then printf(" ") fi;
printf("%06d", v(n, m));
od;
printf("\n");
od;
end;
The above Maple code produces the following table:
> print_table(6);
000001 000001 000001 000001 000001 000001
000001 000003 000005 000008 000011 000015
000001 000005 000017 000042 000091 000180
000001 000008 000042 000179 000633 002001
000001 000011 000091 000633 003835 020755
000001 000015 000180 002001 020755 200082
It can calculate values that are completely out of reach for NAUTY
like the sequence of non-isomorphic spanning subgraphs of $K_{n,n}$
which is
$$1, 3, 17, 179, 3835, 200082, 29610804, 13702979132, \\
20677458750966, 103609939177198046, 1745061194503344181714, \\
99860890306900024150675406,\ldots$$
which points us to OEIS A054976
where we find confirmation of the above calculation and a slightly
different interpretation of the problem statement.
The function v in the Maple program implements model $P$ and the
function v_pre_pie implements model $Q.$
For bicolored versions of $K_{n,n}$ model $Q$ gives the sequence
$$2, 7, 36, 317, 5624, 251610, 33642660, 14685630688, \\
21467043671008, 105735224248507784, 1764356230257807614296, \\
100455994644460412263071692,\ldots$$
which points us to OEIS A002724,
where the calculation is confirmed.
This MSE Meta Link
has many more PET computations by various users.
Thanks go to the authors of the NAUTY package.
