Htykuut

I'm interested in the number of binary matrices of a given size that are distinct with regard to row and column permutations.

If $sim$ is the equivalence relation on $ntimes m$ binary matrices such that $A sim B$ iff one can obtain B from applying a permutation matrix to A, I'm interested in the number of $sim$-equivalence classes over all $ntimes m$ binary matrices.

I know there are $2^{nm}$ binary matrices of size $ntimes m$, and $n!m!$ possible permutations, but somehow I fail to get an intuition on what this implies for the equivalence classes.

This is solved here using Pólya enumeration theory. For the square case ($n=m$), see this sequence.

Comment: I found these by searching for $1,2,7$ in the OEIS.

Here is a computational contribution that treats the case of a square
matrix. As pointed out this problem can be solved using the Polya
Enumeration Theorem. In fact if we are interested only in counting
these matrices, then the Burnside lemma will suffice. We just need to
compute the cycle index of the group acting on the slots of the
matrix.

These cycle indices are easy to compute and we do not need to iterate
over all $(n!)^2$ pairs of permutations (acting on rows and columns)
but instead it is sufficient to iterate over pairs of terms from the
cycle index $Z(S_n)$ of the symmetric group $S_n$ according to their
multiplicities to obtain the cycle index $Z(Q_n)$ of the combined
action on rows and columns. The number of terms here is the much
better count of the number of partitions of $n$ squared (upper bound).

Now for a pair of cycles, one of length $l_1$ from a row permutation
$alpha$ and another of length $l_2$ from a column permutation $beta$
their contribution to the disjoint cycle decomposition product for
$(alpha,beta)$ in the cycle index $Z(Q_n)$ 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)}.$$

The algorithm now becomes very simple -- iterate over pairs of terms
as described above, collect the contribution from each pair of cycles
and add it to the cycle index being computed.

This gives the following cycle indices (only the first four are
shown):

$$Z(Q_2) = 1/4,{a_{{1}}}^{4}+3/4,{a_{{2}}}^{2},$$

$$Z(Q_3) =
1/36,{a_{{1}}}^{9}+1/6,{a_{{1}}}^{3}{a_{{2}}}^{3}+1/4,a_{{
1}}{a_{{2}}}^{4}+2/9,{a_{{3}}}^{3}+1/3,a_{{3}}a_{{6}},$$

$$Z(Q_4) =
{frac {{a_{{1}}}^{16}}{576}}+1/48,{a_{{1}}}^{8}{a_{{2}}}^{4
}+1/16,{a_{{1}}}^{4}{a_{{2}}}^{6}+1/36,{a_{{1}}}^{4}{a_{{3}
}}^{4}+{frac {17,{a_{{2}}}^{8}}{192}}\+1/6,{a_{{1}}}^{2}a_{
{2}}{a_{{3}}}^{2}a_{{6}}+1/9,a_{{1}}{a_{{3}}}^{5}+1/12,{a_{
{2}}}^{2}{a_{{6}}}^{2}+{frac {13,{a_{{4}}}^{4}}{48}}+1/6,a
_{{4}}a_{{12}},$$

and

$$Z(Q_5) =
{frac {{a_{{1}}}^{25}}{14400}}+{frac {{a_{{1}}}^{15}{a_{{2}
}}^{5}}{720}}+{frac {{a_{{1}}}^{9}{a_{{2}}}^{8}}{144}}+{
frac {{a_{{1}}}^{10}{a_{{3}}}^{5}}{360}}+{frac {{a_{{1}}}^{
5}{a_{{2}}}^{10}}{480}}\+1/48,{a_{{1}}}^{3}{a_{{2}}}^{11}+{
frac {a_{{1}}{a_{{2}}}^{12}}{64}}+1/36,{a_{{1}}}^{6}{a_{{2}
}}^{2}{a_{{3}}}^{3}a_{{6}}+1/36,{a_{{1}}}^{4}{a_{{3}}}^{7}+{
frac {{a_{{1}}}^{5}{a_{{4}}}^{5}}{240}}\+{frac {{a_{{2}}}^{5
}{a_{{3}}}^{5}}{360}}+1/24,{a_{{1}}}^{3}a_{{2}}{a_{{4}}}^{5}
+1/24,{a_{{1}}}^{2}{a_{{2}}}^{4}a_{{3}}{a_{{6}}}^{2}+1/36,{
a_{{2}}}^{5}{a_{{3}}}^{3}a_{{6}}\+1/16,a_{{1}}{a_{{2}}}^{2}{a
_{{4}}}^{5}+1/24,{a_{{2}}}^{5}a_{{3}}{a_{{6}}}^{2}+1/18,{a_
{{2}}}^{2}{a_{{3}}}^{5}a_{{6}}+1/16,a_{{1}}{a_{{4}}}^{6}\+1/
36,{a_{{2}}}^{2}{a_{{3}}}^{3}{a_{{6}}}^{2}+1/12,{a_{{1}}}^{
2}a_{{3}}{a_{{4}}}^{2}a_{{12}}+1/12,a_{{2}}a_{{3}}{a_{{4}}}^
{2}a_{{12}}+{frac {13,{a_{{5}}}^{5}}{300}}\+1/30,{a_{{5}}}^
{3}a_{{10}}+1/15,{a_{{5}}}^{2}a_{{15}}+1/20,a_{{5}}{a_{{10}
}}^{2}+1/10,a_{{5}}a_{{20}}+1/15,a_{{10}}a_{{15}}.$$

Evaluating these cycle indices with the variables set to two we
quickly obtain the sequence

$$2, 7, 36, 317, 5624, 251610, 33642660, 14685630688,\
21467043671008, 105735224248507784,1764356230257807614296,\
100455994644460412263071692,19674097197480928600253198363072,\
13363679231028322645152300040033513414,\
31735555932041230032311939400670284689732948,ldots$$
which is indeed OEIS A002724.

Note that the cycle indices make it possible to enumerate
configurations with more than two possible entries or with entries
having different weights. For example, with a 3x3 square and three
colors $A,B$ and $C$ we get the generating function

$$Z(Q_3)(A+B+C) = 1/36, left( A+B+C right) ^{9}+1/6,
left( A+B+C right) ^{3} left( {A}^{2}+{B}^{2}+{C}^{2}
right) ^{3}\+2/9, left( {A}^{3}+{B}^{3}+{C}^{3} right)^{3}
+1/4, left( A+B+C right) left( {A}^{2}+{B}^{2}+{C}^{2
} right) ^{4}\+1/3, left( {A}^{3}+{B}^{3}+{C}^{3}
right) left( {A}^{6}+{B}^{6}+{C}^{6} right)$$

which expands to

$${A}^{9}+{A}^{8}B+{A}^{8}C+3,{A}^{7}{B}^{2}+3,{A}^{7}B
C+3,{A}^{7}{C}^{2}+6,{A}^{6}{B}^{3}+10,{A}^{6}{B}^{2
}C\+10,{A}^{6}B{C}^{2}+6,{A}^{6}{C}^{3}+7,{A}^{5}{B}^
{4}+17,{A}^{5}{B}^{3}C+28,{A}^{5}{B}^{2}{C}^{2}\+17,{
A}^{5}B{C}^{3}+7,{A}^{5}{C}^{4}+7,{A}^{4}{B}^{5}+22,
{A}^{4}{B}^{4}C+43,{A}^{4}{B}^{3}{C}^{2}+43,{A}^{4}{B
}^{2}{C}^{3}\+22,{A}^{4}B{C}^{4}+7,{A}^{4}{C}^{5}+6,{
A}^{3}{B}^{6}+17,{A}^{3}{B}^{5}C+43,{A}^{3}{B}^{4}{C}
^{2}+54,{A}^{3}{B}^{3}{C}^{3}\+43,{A}^{3}{B}^{2}{C}^{4
}+17,{A}^{3}B{C}^{5}+6,{A}^{3}{C}^{6}+3,{A}^{2}{B}^{
7}+10,{A}^{2}{B}^{6}C+28,{A}^{2}{B}^{5}{C}^{2}\+43,{A
}^{2}{B}^{4}{C}^{3}+43,{A}^{2}{B}^{3}{C}^{4}+28,{A}^{
2}{B}^{2}{C}^{5}+10,{A}^{2}B{C}^{6}+3,{A}^{2}{C}^{7}\+
A{B}^{8}+3,A{B}^{7}C+10,A{B}^{6}{C}^{2}+17,A{B}^{5}{
C}^{3}+22,A{B}^{4}{C}^{4}+17,A{B}^{3}{C}^{5}\+10,A{B}
^{2}{C}^{6}+3,AB{C}^{7}+A{C}^{8}+{B}^{9}+{B}^{8}C+3,{
B}^{7}{C}^{2}+6,{B}^{6}{C}^{3}+7,{B}^{5}{C}^{4}\+7,{B
}^{4}{C}^{5}+6,{B}^{3}{C}^{6}+3,{B}^{2}{C}^{7}+B{C}^{
8}+{C}^{9}.$$

This is the Maple code for this computation. Here we have two slightly
different ways of evaluating the count, the first by substituting into
the cycle index and the second by skipping the cycle index altogether
and evaluating all variables at two during processing. The latter
should be used when we are interested ony in the count as opposed to
classifying configurations according to the number of each color /
value that are present.

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, 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_sqmat :=

proc(n)

option remember;

local sind, cind, term_a, term_b, v_a, v_b,

    len_a, len_b, inst_a, inst_b, p;



    cind := 0;



    if n=1 then

        sind := [a[1]];

    else

        sind := pet_cycleind_symm(n);

    fi;



    for term_a in sind do

        for term_b in sind do

            p := 1;

            for v_a in indets(term_a) do

                len_a := op(1, v_a);

                inst_a := degree(term_a, v_a);



                for v_b in indets(term_b) do

                    len_b := op(1, v_b);

                    inst_b := degree(term_b, v_b);



                    p := p*a[lcm(len_a, len_b)]

                    ^(gcd(len_a, len_b)*inst_a*inst_b);

                od;

            od;



            cind := cind +

            lcoeff(term_a)*lcoeff(term_b)*p;

        od;

    od;



    cind;

end;



v :=

proc(n)

    option remember;

    local cind, vars, sbl;



    cind := pet_cycleind_sqmat(n);



    vars := indets(cind);

    sbl := [seq(v=2, v in vars)];



    subs(sbl, cind);

end;



w := 

proc(n)

option remember;

local sind, count, term_a, term_b, v_a, v_b,

    len_a, len_b, inst_a, inst_b, p;



    count := 0;



    if n=1 then

        sind := [a[1]];

    else

        sind := pet_cycleind_symm(n);

    fi;



    for term_a in sind do

        for term_b in sind do

            p := 1;

            for v_a in indets(term_a) do

                len_a := op(1, v_a);

                inst_a := degree(term_a, v_a);



                for v_b in indets(term_b) do

                    len_b := op(1, v_b);

                    inst_b := degree(term_b, v_b);



                    p := p*

                    2^(gcd(len_a, len_b)*inst_a*inst_b);

                od;

            od;



            count := count +

            lcoeff(term_a)*lcoeff(term_b)*p;

        od;

    od;



    count;

end;

This MSE Meta Link
has many more PET computations by various users.