Problem Summary

Let \(E_n\) be the set of binary \(n\times n\) matrices whose row sums and column sums are all even:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

Two matrices are considered equivalent when one can be obtained from the other by permuting rows and permuting columns. The quantity \(c_n\) is the number of equivalence classes of \(E_n\). The problem asks for \(c_{20}\bmod 1001001011\).

Mathematical Approach

The key idea is to count orbits of matrices with even row and column parity under the group \(S_n\times S_n\). Burnside's lemma reduces the problem to fixed matrices for a pair of permutations, and those fixed matrices can be described entirely from the cycle lengths of the two permutations.

Step 1: Apply Burnside's Lemma

The symmetric group \(S_n\) acts on rows, another copy acts on columns, and the combined action is

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

Therefore

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

where \(\operatorname{Fix}(\sigma,\tau)\) is the number of matrices in \(E_n\) unchanged by simultaneously permuting rows by \(\sigma\) and columns by \(\tau\).

Step 2: Replace Permutations by Cycle Types

If \(\lambda\vdash n\) is a partition of \(n\), write \(m_k(\lambda)\) for the number of cycles of length \(k\). The number of permutations with cycle type \(\lambda\) is

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

So Burnside's sum depends only on pairs of partitions \((\lambda,\mu)\), not on individual permutations. If the row-cycle lengths are \(p_1,\dots,p_u\) and the column-cycle lengths are \(q_1,\dots,q_v\), then every quantity in the fixed-count calculation can be expressed in terms of these lengths.

Step 3: Count Cell Orbits for One Type Pair

Consider one row cycle of length \(p_i\) and one column cycle of length \(q_j\). Their Cartesian block contains \(p_iq_j\) cells, and the joint action moves each cell by one step along both cycles. That block splits into

$$d_{ij}=\gcd(p_i,q_j)$$

cell orbits. Hence the total number of binary orbit variables is

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

If there were no parity conditions, the fixed matrices for \((\lambda,\mu)\) would simply number \(2^{o(\lambda,\mu)}\).

Step 4: Translate Even Row and Column Sums into \(\mathbb{F}_2\) Constraints

Inside one \((p_i,q_j)\)-block, let \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) be the orbit bits. A single row in the row cycle sees each orbit bit exactly \(q_j/d_{ij}\) times, so modulo \(2\) that block contributes to the row-parity equations if and only if

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i).$$

Similarly, a single column in the column cycle sees each orbit bit exactly \(p_i/d_{ij}\) times, so the block contributes to the column-parity equations if and only if

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j).$$

Only the block sum

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2$$

enters the parity equations. The remaining \(d_{ij}-1\) directions inside that block are always free. This is why the final number of free bits takes the form \(o(\lambda,\mu)-\rho(\lambda,\mu)\) for a suitable rank \(\rho\).

Step 5: Compute the Parity Rank from 2-Adic Valuations

Now work at the level of cycles. Introduce one row-equation node for each row cycle and one column-equation node for each column cycle. The cancellation rules for the block sums are:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{the row node and column node must agree},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{the row node is forced to }0,$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{the column node is forced to }0.$$

Merging equal nodes and marking forced-zero nodes gives a small graph problem. Let \(f(\lambda,\mu)\) be the number of connected components that are not forced to zero. Since there are \(u+v\) equation nodes in total, the rank of the parity system is

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

Therefore

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

Step 6: Final Burnside Sum

Substituting the type weights and the fixed counts gives the exact formula used by the implementations:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

The hard part is no longer matrix enumeration; it is only partition enumeration together with a small rank computation for each pair of cycle types.

Worked Example: \(n=3\)

The partitions of \(3\) are \(1+1+1\), \(2+1\), and \(3\), with permutation counts \(1\), \(3\), and \(2\). Using the orbit formula and the valuation rules above, the fixed-count matrix is

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

For instance, for \((2+1,2+1)\) we have \(o=2+1+1+1=5\). The four equation nodes split so that three independent constraints remain, hence \(\rho=3\) and \(\operatorname{Fix}=2^{5-3}=4\).

Burnside now gives

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

This small case already shows the full mechanism: cycle types determine orbit counts, 2-adic valuations determine the parity rank, and Burnside averages the resulting fixed counts.

How the Code Works

The C++, Python, and Java implementations all follow the same pipeline. First they enumerate all partitions of \(n\) and store, for each cycle type, its cycle-length multiplicities, its list of 2-adic valuations, and the number of permutations with that type. Next they precompute \(\gcd(p,q)\) for every \(1\le p,q\le n\).

For each pair of cycle types, the implementation computes \(o(\lambda,\mu)\) by summing the relevant gcd values. It then builds the small cycle-level equality/forced-zero structure described above, counts the number of unforced connected components, obtains the parity rank \(\rho(\lambda,\mu)\), and adds

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

to an exact integer accumulator. After all type pairs are processed, the total is divided by \((n!)^2\) as required by Burnside's lemma, and only then is the result reduced modulo \(1001001011\). The C++ implementation parallelizes the outer traversal, while the Python and Java implementations use the same mathematics in a direct single-process sum.

Complexity Analysis

Let \(p(n)\) be the partition number. Enumerating all cycle types costs \(O(p(n)\,n)\) space and roughly the same order of bookkeeping work. The main Burnside sum runs over \(p(n)^2\) type pairs. For each pair, the orbit count and the parity-rank computation each need at most \(O(n^2)\) elementary operations, because there are at most \(n\) row cycles and \(n\) column cycles.

Thus the overall time complexity is \(O(p(n)^2 n^2)\), and the memory usage is \(O(p(n)\,n)\) plus the exact-integer accumulator. For \(n=20\), this is entirely practical because the number of cycle types is small compared with the astronomical number of matrices themselves.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=626
  2. Burnside's lemma: Wikipedia - Burnside's lemma
  3. Integer partition: Wikipedia - Partition (number theory)
  4. \(p\)-adic valuation: Wikipedia - p-adic valuation
  5. Disjoint-set data structure: Wikipedia - Disjoint-set data structure

Problemzusammenfassung

Sei \(E_n\) die Menge aller binären \(n\times n\)-Matrizen, deren Zeilen- und Spaltensummen sämtlich gerade sind:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

Zwei Matrizen gelten als äquivalent, wenn man die eine durch Permutation der Zeilen und Permutation der Spalten aus der anderen erhält. Mit \(c_n\) bezeichnen wir die Anzahl dieser Äquivalenzklassen. Gesucht ist \(c_{20}\bmod 1001001011\).

Mathematischer Ansatz

Die Grundidee ist, Orbits der Matrizen mit geraden Zeilen- und Spaltensummen unter der Gruppenwirkung \(S_n\times S_n\) zu zählen. Das Lemma von Burnside reduziert alles auf Fixpunktzahlen für Permutationspaare, und diese Fixpunktzahlen hängen nur von den Zykluslängen der beiden Permutationen ab.

Schritt 1: Burnside anwenden

Eine Kopie von \(S_n\) permutiert die Zeilen, eine zweite die Spalten. Die kombinierte Wirkung lautet

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

Damit gilt

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

wobei \(\operatorname{Fix}(\sigma,\tau)\) die Anzahl der Matrizen in \(E_n\) bezeichnet, die unter gleichzeitiger Zeilen- und Spaltenpermutation unverändert bleiben.

Schritt 2: Permutationen durch Zyklustypen ersetzen

Ist \(\lambda\vdash n\) eine Partition von \(n\), so sei \(m_k(\lambda)\) die Anzahl der Zyklen der Länge \(k\). Die Zahl der Permutationen mit diesem Zyklustyp ist

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

Also hängt die Burnside-Summe nur von Paaren \((\lambda,\mu)\) von Partitionen ab. Haben die Zeilenzyklen die Längen \(p_1,\dots,p_u\) und die Spaltenzyklen die Längen \(q_1,\dots,q_v\), dann lassen sich alle weiteren Größen allein aus diesen Längen berechnen.

Schritt 3: Zähle Zellenorbits für ein Typenpaar

Betrachte einen Zeilenzyklus der Länge \(p_i\) und einen Spaltenzyklus der Länge \(q_j\). Der zugehörige Block besitzt \(p_iq_j\) Zellen, und die gemeinsame Wirkung verschiebt jede Zelle um einen Schritt in beiden Zyklen. Dieser Block zerfällt in

$$d_{ij}=\gcd(p_i,q_j)$$

Orbits von Zellen. Daher ist die Gesamtzahl der binären Orbitvariablen

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

Ohne Paritätsbedingungen gäbe es also einfach \(2^{o(\lambda,\mu)}\) fixierte Matrizen.

Schritt 4: Gerade Zeilen- und Spaltensummen als \(\mathbb{F}_2\)-System

In einem \((p_i,q_j)\)-Block seien \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) die Orbitbits. Eine einzelne Zeile des Zeilenzyklus enthält jedes Orbitbit genau \(q_j/d_{ij}\)-mal. Modulo \(2\) trägt dieser Block also genau dann zu den Zeilenparitäten bei, wenn

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i).$$

Analog enthält eine einzelne Spalte des Spaltenzyklus jedes Orbitbit genau \(p_i/d_{ij}\)-mal. Der Block trägt also genau dann zu den Spaltenparitäten bei, wenn

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j).$$

Für die Paritätsgleichungen ist nur die Blocksumme

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2$$

relevant. Die übrigen \(d_{ij}-1\) Freiheitsgrade innerhalb des Blocks bleiben immer frei. Deshalb hat die Anzahl freier Bits die Form \(o(\lambda,\mu)-\rho(\lambda,\mu)\) für einen passenden Rang \(\rho\).

Schritt 5: Rang über 2-adische Bewertungen bestimmen

Nun arbeiten wir auf Zyklusebene. Für jeden Zeilenzyklus gibt es einen Zeilenknoten, für jeden Spaltenzyklus einen Spaltenknoten. Die Regeln für das Verschwinden der Blocksummen lauten:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{Zeilenknoten und Spaltenknoten sind gleich},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{der Zeilenknoten ist auf }0\text{ festgelegt},$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{der Spaltenknoten ist auf }0\text{ festgelegt}.$$

Man verschmilzt also gleiche Knoten und markiert erzwungene Nullknoten. Sei \(f(\lambda,\mu)\) die Anzahl der zusammenhängenden Komponenten, die nicht auf Null gezwungen sind. Bei insgesamt \(u+v\) Gleichungsknoten ist der Rang des Paritätssystems dann

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

Daraus folgt

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

Schritt 6: Die endgültige Burnside-Formel

Mit Typgewichten und Fixpunktzahlen ergibt sich exakt die Formel, die von den Implementierungen ausgewertet wird:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

Das ursprüngliche Problem ist damit nicht mehr eine Enumeration von Matrizen, sondern eine Enumeration von Partitionen mit einer kleinen Rangberechnung pro Typenpaar.

Durchgerechnetes Beispiel: \(n=3\)

Die Partitionen von \(3\) sind \(1+1+1\), \(2+1\) und \(3\), mit den Gewichten \(1\), \(3\) und \(2\). Aus Orbitformel und Bewertungsregeln erhält man die Fixpunktmatrix

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

Zum Beispiel gilt für \((2+1,2+1)\): \(o=2+1+1+1=5\). Von den vier Gleichungsknoten bleiben effektiv drei unabhängige Bedingungen übrig, also \(\rho=3\) und damit \(\operatorname{Fix}=2^{5-3}=4\).

Burnside liefert nun

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

Schon dieser kleine Fall zeigt die komplette Methode: Zyklustypen liefern die Orbitzahl, 2-adische Bewertungen liefern den Paritätsrang, und Burnside mittelt die Fixpunktzahlen.

Wie der Code arbeitet

Die C++-, Python- und Java-Implementierungen verwenden dieselbe Pipeline. Zuerst werden alle Partitionen von \(n\) erzeugt. Für jeden Zyklustyp werden die Häufigkeiten der Zykluslängen, die Liste der 2-adischen Bewertungen und die Anzahl der Permutationen dieses Typs gespeichert. Danach werden alle Werte \(\gcd(p,q)\) für \(1\le p,q\le n\) vorab berechnet.

Für jedes Paar von Zyklustypen berechnet die Implementierung zunächst \(o(\lambda,\mu)\) über die passenden gcd-Beiträge. Anschließend baut sie die kleine Struktur aus Gleichheitsbeziehungen und erzwungenen Nullknoten auf, zählt die nicht erzwungenen Komponenten, erhält daraus den Rang \(\rho(\lambda,\mu)\) und addiert

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

zu einem exakten Integer-Akkumulator. Nach dem Durchlauf aller Typenpaare wird gemäß Burnside durch \((n!)^2\) dividiert; erst ganz am Ende erfolgt die Reduktion modulo \(1001001011\). Die C++-Implementierung parallelisiert die äußere Schleife, während Python und Java dieselbe Summe direkt in einem Prozess auswerten.

Komplexitätsanalyse

Sei \(p(n)\) die Partitionszahl. Das Erzeugen aller Zyklustypen benötigt \(O(p(n)\,n)\) Speicher und ungefähr denselben Aufwand für die dazugehörigen Metadaten. Die Burnside-Hauptsumme läuft über \(p(n)^2\) Typenpaare. Für jedes Paar brauchen sowohl Orbitzahl als auch Rangberechnung höchstens \(O(n^2)\) elementare Operationen, weil es höchstens \(n\) Zeilenzyklen und \(n\) Spaltenzyklen gibt.

Insgesamt ergibt sich also \(O(p(n)^2 n^2)\) Zeit und \(O(p(n)\,n)\) Speicher, zuzüglich des exakten Ganzzahlakkumulators. Für \(n=20\) ist das bequem praktikabel, obwohl die Anzahl der zugrunde liegenden Matrizen astronomisch groß ist.

Fußnoten und Referenzen

  1. Problemseite: https://projecteuler.net/problem=626
  2. Lemma von Burnside: Wikipedia - Burnside's lemma
  3. Ganzzahlpartition: Wikipedia - Partition (number theory)
  4. \(p\)-adische Bewertung: Wikipedia - p-adic valuation
  5. Disjoint-Set-Datenstruktur: Wikipedia - Disjoint-set data structure

Problem Özeti

\(E_n\), satır toplamları ve sütun toplamları tamamen çift olan ikili \(n\times n\) matrislerin kümesi olsun:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

İki matris, biri diğerinden satır ve sütun permütasyonlarıyla elde edilebiliyorsa eşdeğer kabul edilir. \(c_n\), bu eşdeğerlik sınıflarının sayısıdır. İstenen değer \(c_{20}\bmod 1001001011\)'dir.

Matematiksel Yaklaşım

Ana fikir, satır ve sütun paritesi çift olan matrislerin \(S_n\times S_n\) etkisi altındaki yörüngelerini saymaktır. Burnside lemması problemi, bir satır permütasyonu ile bir sütun permütasyonunun birlikte sabit bıraktığı matrislerin sayısına indirger. Bu sabit sayılar da yalnızca çevrim uzunluklarından belirlenir.

Adım 1: Burnside Lemmasını uygula

\(S_n\)'nin bir kopyası satırları, diğer kopyası sütunları permüte eder. Birleşik etki

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n$$

şeklindedir. Bu yüzden

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau)$$

olur. Burada \(\operatorname{Fix}(\sigma,\tau)\), hem satır hem sütun permütasyonu altında değişmeden kalan \(E_n\) içindeki matris sayısıdır.

Adım 2: Permütasyonları çevrim tipleriyle özetle

\(\lambda\vdash n\) bir bölüm olsun ve \(m_k(\lambda)\), uzunluğu \(k\) olan çevrim sayısını göstersin. Bu çevrim tipine sahip permütasyon sayısı

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}$$

şeklindedir. Dolayısıyla Burnside toplamı tek tek permütasyonlar yerine bölüm çiftleri \((\lambda,\mu)\) üzerinde yürütülür. Satır çevrim uzunlukları \(p_1,\dots,p_u\), sütun çevrim uzunlukları \(q_1,\dots,q_v\) olsun; bundan sonra gereken her nicelik bu uzunluklardan hesaplanır.

Adım 3: Bir çevrim tipi çifti için hücre yörüngelerini say

Uzunluğu \(p_i\) olan bir satır çevrimi ile uzunluğu \(q_j\) olan bir sütun çevrimi birlikte \(p_iq_j\) hücrelik bir blok oluşturur. Ortak etki her hücreyi iki çevrim boyunca birer adım ilerlettiği için bu blok

$$d_{ij}=\gcd(p_i,q_j)$$

adet hücre yörüngesine ayrılır. Böylece toplam ikili yörünge değişkeni sayısı

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j)$$

olur. Parite koşulları olmasa, \((\lambda,\mu)\) için sabit matris sayısı doğrudan \(2^{o(\lambda,\mu)}\) olurdu.

Adım 4: Çift satır ve sütun toplamlarını \(\mathbb{F}_2\) üzerinde yaz

Bir \((p_i,q_j)\)-bloğunda \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) yörünge bitleri olsun. Satır çevrimindeki tek bir satır, her yörünge bitini tam olarak \(q_j/d_{ij}\) kez görür. Bu yüzden mod \(2\) açısından blok, satır paritesi denklemlerine ancak

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i)$$

iken katkı yapar. Benzer biçimde bir sütun, her yörünge bitini \(p_i/d_{ij}\) kez görür; dolayısıyla blok sütun paritesi denklemlerine ancak

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j)$$

iken katkı yapar. Parite denklemlerine yalnızca

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2$$

blok toplamı girer. Bloğun içindeki kalan \(d_{ij}-1\) serbestlik derecesi her zaman bağımsızdır. Bu yüzden toplam serbest bit sayısı \(o(\lambda,\mu)-\rho(\lambda,\mu)\) biçiminde ortaya çıkar.

Adım 5: Parite rütbesini 2-adik değerlemelerden çıkar

Şimdi çevrim düzeyinde düşünelim. Her satır çevrimi için bir satır düğümü, her sütun çevrimi için bir sütun düğümü tanımlanır. Blok toplamlarının yok olması için gereken kurallar şunlardır:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{satır düğümü ile sütun düğümü eşit olmalıdır},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{satır düğümü }0\text{ olmaya zorlanır},$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{sütun düğümü }0\text{ olmaya zorlanır}.$$

Eşit olması gereken düğümler birleştirilir, sıfıra zorlananlar işaretlenir. \(f(\lambda,\mu)\), sıfıra zorlanmamış bağlı bileşen sayısı olsun. Toplam \(u+v\) denklem düğümü olduğundan parite sisteminin rütbesi

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu)$$

olur. Dolayısıyla

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

elde edilir.

Adım 6: Son Burnside toplamı

Tip ağırlıkları ile sabit sayıları birleştirince uygulamanın hesapladığı kesin formül ortaya çıkar:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

Böylece problem artık matrisleri tek tek saymak değil, bölüm çiftlerini gezip her çift için küçük bir rütbe hesabı yapmaktır.

Çözümlü Örnek: \(n=3\)

\(3\)'ün bölümleri \(1+1+1\), \(2+1\) ve \(3\)'tür; bunların permütasyon ağırlıkları sırasıyla \(1\), \(3\) ve \(2\)'dir. Yörünge formülü ile değerleme kuralları kullanıldığında sabit-matris tablosu

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}$$

olur. Örneğin \((2+1,2+1)\) için \(o=2+1+1+1=5\)'tir. Dört denklem düğümünden etkin olarak üç bağımsız koşul kaldığı için \(\rho=3\) ve \(\operatorname{Fix}=2^{5-3}=4\) çıkar.

Burnside şimdi

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3$$

verir. Bu küçük örnek tüm yöntemi özetler: çevrim tipleri yörünge sayısını, 2-adik değerlemeler parite rütbesini belirler; Burnside da bunları ortalamaya çevirir.

Kod Nasıl Çalışır

C++, Python ve Java uygulamaları aynı iş akışını izler. Önce \(n\)'in tüm bölümleri üretilir; her çevrim tipi için çevrim uzunluğu çoklukları, 2-adik değerleme listesi ve bu tipe sahip permütasyon sayısı saklanır. Ardından \(1\le p,q\le n\) için tüm \(\gcd(p,q)\) değerleri önceden hesaplanır.

Her çevrim tipi çifti için uygulama önce ilgili gcd katkılarını toplayarak \(o(\lambda,\mu)\)'yu hesaplar. Sonra yukarıdaki eşitlik ve zorunlu-sıfır ilişkilerinden oluşan küçük yapıyı kurar, sıfıra zorlanmamış bağlı bileşen sayısını bulur, buradan \(\rho(\lambda,\mu)\)'yu çıkarır ve

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

terimini tam sayı akümülatörüne ekler. Tüm çiftler işlendiğinde toplam \((n!)^2\)'ye bölünür; mod \(1001001011\) işlemi ise en sonda uygulanır. C++ sürümü dış döngüyü paralelleştirir; Python ve Java sürümleri aynı toplamı tek süreçte doğrudan değerlendirir.

Karmaşıklık Analizi

\(p(n)\), bölüm sayısı olsun. Tüm çevrim tiplerini üretmek \(O(p(n)\,n)\) bellek ve aynı mertebede yardımcı işlem gerektirir. Ana Burnside toplamı \(p(n)^2\) adet tip çifti üzerinde çalışır. Her çift için hem yörünge sayısı hem de parite rütbesi hesabı en fazla \(O(n^2)\) temel işlem ister; çünkü en fazla \(n\) satır çevrimi ve \(n\) sütun çevrimi vardır.

Dolayısıyla toplam zaman karmaşıklığı \(O(p(n)^2 n^2)\), bellek karmaşıklığı ise tam sayı akümülatörü dışında \(O(p(n)\,n)\)'dir. \(n=20\) için bu yaklaşım tamamen pratiktir; çünkü bölüm sayısı küçük kalırken ham matris uzayı astronomik büyüklüktedir.

Dipnotlar ve Referanslar

  1. Problem sayfası: https://projecteuler.net/problem=626
  2. Burnside lemması: Wikipedia - Burnside's lemma
  3. Tamsayı bölümü: Wikipedia - Partition (number theory)
  4. \(p\)-adik değerleme: Wikipedia - p-adic valuation
  5. Ayrık-küme veri yapısı: Wikipedia - Disjoint-set data structure

Resumen del Problema

Sea \(E_n\) el conjunto de matrices binarias \(n\times n\) cuyas sumas por fila y por columna son todas pares:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

Dos matrices se consideran equivalentes si una puede obtenerse de la otra permutando filas y permutando columnas. Denotamos por \(c_n\) el número de clases de equivalencia de \(E_n\). La tarea es hallar \(c_{20}\bmod 1001001011\).

Enfoque Matemático

La idea central es contar órbitas de matrices con paridad par en filas y columnas bajo la acción del grupo \(S_n\times S_n\). El lema de Burnside reduce el problema a contar matrices fijas por un par de permutaciones, y esos recuentos dependen únicamente de las longitudes de ciclo de ambas permutaciones.

Paso 1: Aplicar el lema de Burnside

Una copia de \(S_n\) actúa sobre las filas y otra sobre las columnas. La acción conjunta es

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

Por tanto,

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

donde \(\operatorname{Fix}(\sigma,\tau)\) es el número de matrices de \(E_n\) que no cambian al permutar simultáneamente filas por \(\sigma\) y columnas por \(\tau\).

Paso 2: Sustituir permutaciones por tipos de ciclo

Si \(\lambda\vdash n\) es una partición de \(n\), escribimos \(m_k(\lambda)\) para el número de ciclos de longitud \(k\). El número de permutaciones con ese tipo de ciclo es

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

Así, la suma de Burnside depende solo de pares de particiones \((\lambda,\mu)\). Si las longitudes de ciclo en filas son \(p_1,\dots,p_u\) y en columnas son \(q_1,\dots,q_v\), todo el recuento de matrices fijas puede expresarse con esos números.

Paso 3: Contar órbitas de celdas para un par de tipos

Tomemos un ciclo de filas de longitud \(p_i\) y un ciclo de columnas de longitud \(q_j\). El bloque correspondiente tiene \(p_iq_j\) celdas, y la acción conjunta mueve cada celda un paso en ambos ciclos. Ese bloque se descompone en

$$d_{ij}=\gcd(p_i,q_j)$$

órbitas de celdas. Por consiguiente, el número total de variables binarias de órbita es

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

Sin restricciones de paridad, habría simplemente \(2^{o(\lambda,\mu)}\) matrices fijas.

Paso 4: Traducir la paridad par a ecuaciones sobre \(\mathbb{F}_2\)

Dentro de un bloque \((p_i,q_j)\), sean \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) los bits de órbita. Una fila del ciclo de filas ve cada bit exactamente \(q_j/d_{ij}\) veces, de modo que el bloque interviene en las ecuaciones de paridad de filas si y solo si

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i).$$

De forma análoga, una columna del ciclo de columnas ve cada bit \(p_i/d_{ij}\) veces, y el bloque interviene en las ecuaciones de columnas si y solo si

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j).$$

En las ecuaciones de paridad solo aparece la suma del bloque

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2.$$

Las otras \(d_{ij}-1\) direcciones internas del bloque quedan libres. Por eso el número final de bits libres toma la forma \(o(\lambda,\mu)-\rho(\lambda,\mu)\), donde \(\rho\) es el rango del sistema de paridad.

Paso 5: Obtener el rango a partir de las valuaciones 2-ádicas

Ahora pasamos al nivel de ciclos. Introducimos un nodo de ecuación por cada ciclo de filas y otro por cada ciclo de columnas. Las reglas de cancelación para las sumas de bloque son:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{el nodo de fila y el nodo de columna deben coincidir},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{el nodo de fila queda forzado a }0,$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{el nodo de columna queda forzado a }0.$$

Se fusionan los nodos iguales y se marcan los nodos forzados a cero. Sea \(f(\lambda,\mu)\) el número de componentes conexas que no están forzadas a cero. Como el total de nodos de ecuación es \(u+v\), el rango del sistema es

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

Por lo tanto,

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

Paso 6: Fórmula final de Burnside

Al combinar los pesos de tipo con los recuentos fijos se obtiene exactamente la fórmula que evalúan las implementaciones:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

El problema deja de ser una enumeración directa de matrices y se convierte en una enumeración de particiones con un pequeño cálculo de rango para cada par de tipos.

Ejemplo Resuelto: \(n=3\)

Las particiones de \(3\) son \(1+1+1\), \(2+1\) y \(3\), con pesos \(1\), \(3\) y \(2\). Aplicando la cuenta de órbitas y las reglas de valuación se obtiene la matriz

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

Por ejemplo, para \((2+1,2+1)\) se tiene \(o=2+1+1+1=5\). Quedan tres restricciones independientes entre los cuatro nodos de ecuación, así que \(\rho=3\) y \(\operatorname{Fix}=2^{5-3}=4\).

Entonces Burnside da

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

Este caso pequeño ya muestra todo el mecanismo: los tipos de ciclo determinan las órbitas, las valuaciones 2-ádicas determinan el rango de paridad y Burnside promedia los recuentos fijos.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java siguen la misma secuencia. Primero enumeran todas las particiones de \(n\) y guardan, para cada tipo de ciclo, las multiplicidades de longitudes, la lista de valuaciones 2-ádicas y el número de permutaciones de ese tipo. Después precalculan \(\gcd(p,q)\) para todos los pares \(1\le p,q\le n\).

Para cada par de tipos, la implementación calcula \(o(\lambda,\mu)\) sumando los gcd correspondientes. Luego construye la pequeña estructura de igualdades y nodos forzados a cero descrita arriba, cuenta las componentes no forzadas, obtiene \(\rho(\lambda,\mu)\) y añade

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

a un acumulador exacto. Cuando se han procesado todos los pares, el total se divide entre \((n!)^2\) y solo al final se reduce módulo \(1001001011\). La versión en C++ paraleliza el recorrido exterior; las versiones en Python y Java evalúan la misma suma de forma directa en un único proceso.

Análisis de Complejidad

Sea \(p(n)\) el número de particiones de \(n\). Enumerar todos los tipos de ciclo requiere \(O(p(n)\,n)\) memoria y una cantidad comparable de trabajo auxiliar. La suma principal de Burnside recorre \(p(n)^2\) pares de tipos. Para cada par, tanto la cuenta de órbitas como el cálculo del rango de paridad necesitan a lo sumo \(O(n^2)\) operaciones elementales, porque hay como máximo \(n\) ciclos de filas y \(n\) ciclos de columnas.

En consecuencia, el tiempo total es \(O(p(n)^2 n^2)\) y la memoria es \(O(p(n)\,n)\), además del acumulador de enteros exactos. Para \(n=20\) este coste es totalmente razonable, a pesar de que el espacio bruto de matrices es inmenso.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=626
  2. Lema de Burnside: Wikipedia - Burnside's lemma
  3. Partición entera: Wikipedia - Partition (number theory)
  4. Valuación \(p\)-ádica: Wikipedia - p-adic valuation
  5. Estructura disjoint-set: Wikipedia - Disjoint-set data structure

问题概述

设 \(E_n\) 为所有满足“每一行元素和、每一列元素和都为偶数”的二进制 \(n\times n\) 矩阵集合:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

如果一个矩阵可以通过对另一个矩阵做行置换和列置换得到,则把它们视为同一个等价类。记 \(c_n\) 为 \(E_n\) 在这种等价关系下的类数。题目要求的是 \(c_{20}\bmod 1001001011\)。

数学方法

核心思想是:不要直接枚举矩阵,而是把“行列和为偶数的二进制矩阵”看成 \(S_n\times S_n\) 作用下的对象,转而统计轨道数。Burnside 引理把问题化成“给定一对行、列置换时,有多少矩阵保持不变”,而这个固定矩阵数只取决于两边置换的循环长度结构。

步骤 1:使用 Burnside 引理

一份 \(S_n\) 作用在行上,另一份 \(S_n\) 作用在列上。联合起来的作用写成

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

因此有

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

其中 \(\operatorname{Fix}(\sigma,\tau)\) 表示在行置换为 \(\sigma\)、列置换为 \(\tau\) 时保持不变且仍满足偶数奇偶条件的矩阵个数。

步骤 2:把置换压缩成循环类型

若 \(\lambda\vdash n\) 是 \(n\) 的一个整数分拆,记 \(m_k(\lambda)\) 为长度为 \(k\) 的循环个数,则该循环类型对应的置换数为

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

于是 Burnside 双重求和不必对每个置换单独处理,只需对分拆对 \((\lambda,\mu)\) 求和即可。设行循环长度为 \(p_1,\dots,p_u\),列循环长度为 \(q_1,\dots,q_v\),下面所有量都可以由这些长度直接得到。

步骤 3:先算单元格轨道数

取一个长度为 \(p_i\) 的行循环和一个长度为 \(q_j\) 的列循环。它们张成一个 \(p_iq_j\) 个单元的矩形块,在联合置换下,每次同时沿行循环和列循环前进一步。这个块被分成

$$d_{ij}=\gcd(p_i,q_j)$$

个单元格轨道。于是该循环类型对 \((\lambda,\mu)\) 的总轨道数为

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

如果没有“每行每列都要为偶数”这个限制,那么固定矩阵的个数就只是 \(2^{o(\lambda,\mu)}\),因为每个轨道自由选 \(0\) 或 \(1\) 即可。

步骤 4:把偶数条件改写成 \(\mathbb{F}_2\) 线性约束

在一个 \((p_i,q_j)\)-块里,设 \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) 是该块的轨道比特。行循环中的任意一行会看到每个轨道比特恰好出现 \(q_j/d_{ij}\) 次,所以模 \(2\) 以后,这个块对“行和为偶数”的贡献是否存在,只由

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i)$$

决定。类似地,一列会看到每个轨道比特出现 \(p_i/d_{ij}\) 次,因此该块对“列和为偶数”的贡献是否存在,只由

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j)$$

决定。这里 \(\nu_2\) 表示 2-adic valuation,也就是一个整数中含有多少个因子 \(2\)。

更重要的是,进入奇偶方程的并不是块内每个轨道变量本身,而只是它们的总和

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2.$$

因此,一个块里总共有 \(d_{ij}\) 个轨道比特,但只有 1 个线性组合会参与行列奇偶约束,其余 \(d_{ij}-1\) 个方向始终自由。这就解释了为什么最终自由比特数会写成 \(o(\lambda,\mu)-\rho(\lambda,\mu)\)。

步骤 5:用 2-adic 估值关系求奇偶系统的秩

接下来把问题提升到“循环层面”。对每个行循环设一个行方程节点,对每个列循环设一个列方程节点。块和变量 \(s_{ij}\) 的消去规则只有三种:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{行节点与列节点必须相等},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{行节点被强制为 }0,$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{列节点被强制为 }0.$$

把必须相等的节点合并,再把被强制为零的连通块做标记。设 \(f(\lambda,\mu)\) 是最终没有被强制为零的连通分量个数。因为总共有 \(u+v\) 个方程节点,所以奇偶线性系统的秩为

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

于是固定矩阵数就是

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

步骤 6:写出最后的 Burnside 求和公式

把循环类型的权重与固定矩阵数合起来,就得到实现中真正计算的精确公式:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

这样,原问题就从“枚举天文数量级的矩阵”转成了“枚举整数分拆,并对每一对分拆做一个很小的秩计算”。

例子:\(n=3\) 的完整演算

\(3\) 的分拆只有三种:\(1+1+1\)、\(2+1\)、\(3\),对应的置换数量分别是 \(1\)、\(3\)、\(2\)。按照上面的轨道公式与 2-adic 规则,可以得到固定矩阵数表

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

以 \((2+1,2+1)\) 为例:四个块对应的轨道数分别是 \(\gcd(2,2)=2\)、\(\gcd(2,1)=1\)、\(\gcd(1,2)=1\)、\(\gcd(1,1)=1\),所以 \(o=5\)。再看 2-adic 估值,长度为 \(2\) 的循环估值是 \(1\),长度为 \(1\) 的循环估值是 \(0\)。由“相等则合并,大于则强制为零”的规则可知,四个方程节点最终只留下 1 个未被强制为零的连通分量,所以 \(\rho=4-1=3\)。因此 \(\operatorname{Fix}=2^{5-3}=4\)。

接着做 Burnside 平均:

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

这个小例子已经把整套方法完整展示出来了:循环类型决定轨道数量,2-adic 估值决定奇偶系统的秩,Burnside 则把所有固定矩阵数平均成最终的等价类数量。

代码实现说明

C++、Python 和 Java 实现遵循完全相同的计算流程。第一步是枚举 \(n\) 的所有分拆,并为每个循环类型保存三类信息:各循环长度的重数、每个循环长度的 2-adic 估值,以及具有该循环类型的置换数。第二步是预计算所有 \(1\le p,q\le n\) 的 \(\gcd(p,q)\)。

之后对每一对循环类型,先根据 gcd 求出 \(o(\lambda,\mu)\),再构造上面提到的“相等/强制为零”关系图,统计没有被强制为零的连通分量,得到 \(\rho(\lambda,\mu)\),并把

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

加入一个精确整数累加器。全部类型对处理完之后,再按 Burnside 公式除以 \((n!)^2\),最后一步才取模 \(1001001011\)。其中 C++ 版本把最外层遍历并行化,而 Python 与 Java 版本则用同样的数学公式顺序求和。

复杂度分析

记 \(p(n)\) 为整数分拆数。枚举所有循环类型需要 \(O(p(n)\,n)\) 级别的存储和大致同量级的预处理工作。主计算是对所有分拆对做 Burnside 求和,共有 \(p(n)^2\) 对。对任意一对,轨道数计算和奇偶秩计算都只需要最多 \(O(n^2)\) 次基本操作,因为行循环和列循环的数量都不可能超过 \(n\)。

因此总时间复杂度为 \(O(p(n)^2 n^2)\),空间复杂度为 \(O(p(n)\,n)\),外加一个精确大整数累加器。对 \(n=20\) 来说,这个规模完全可行,因为真正爆炸的是原始矩阵空间,而不是分拆空间。

注释与参考资料

  1. 题目页面:https://projecteuler.net/problem=626
  2. Burnside 引理:Wikipedia - Burnside's lemma
  3. 整数分拆:Wikipedia - Partition (number theory)
  4. \(p\)-adic valuation:Wikipedia - p-adic valuation
  5. 并查集 / disjoint-set:Wikipedia - Disjoint-set data structure

Краткое описание задачи

Пусть \(E_n\) обозначает множество бинарных матриц размера \(n\times n\), у которых суммы в каждой строке и в каждом столбце чётны:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

Две матрицы считаются эквивалентными, если одна получается из другой перестановкой строк и перестановкой столбцов. Число \(c_n\) есть количество классов эквивалентности множества \(E_n\). Нужно найти \(c_{20}\bmod 1001001011\).

Математический подход

Главная идея состоит в том, чтобы считать не сами матрицы, а орбиты матриц с чётными суммами строк и столбцов под действием группы \(S_n\times S_n\). Лемма Бёрнсайда сводит задачу к числу матриц, фиксируемых парой перестановок, а это число определяется только длинами циклов этих перестановок.

Шаг 1: Применить лемму Бёрнсайда

Одна копия \(S_n\) действует на строки, другая на столбцы. Совместное действие имеет вид

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

Следовательно,

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

где \(\operatorname{Fix}(\sigma,\tau)\) означает число матриц из \(E_n\), не меняющихся при одновременной перестановке строк по \(\sigma\) и столбцов по \(\tau\).

Шаг 2: Заменить перестановки их циклическими типами

Пусть \(\lambda\vdash n\) есть разбиение числа \(n\), а \(m_k(\lambda)\) обозначает число циклов длины \(k\). Тогда число перестановок такого типа равно

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

Поэтому в формуле Бёрнсайда достаточно суммировать по парам разбиений \((\lambda,\mu)\). Если длины циклов строк равны \(p_1,\dots,p_u\), а длины циклов столбцов равны \(q_1,\dots,q_v\), то все нужные величины выражаются через эти числа.

Шаг 3: Посчитать орбиты ячеек для пары типов

Рассмотрим цикл строк длины \(p_i\) и цикл столбцов длины \(q_j\). Соответствующий блок содержит \(p_iq_j\) ячеек, а совместное действие сдвигает каждую ячейку на один шаг по обоим циклам. Этот блок распадается на

$$d_{ij}=\gcd(p_i,q_j)$$

орбит ячеек. Поэтому общее число бинарных орбитных переменных равно

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

Если бы условий чётности не было, то число фиксированных матриц для \((\lambda,\mu)\) было бы просто \(2^{o(\lambda,\mu)}\).

Шаг 4: Перевести условия чётности в линейную систему над \(\mathbb{F}_2\)

Внутри блока \((p_i,q_j)\) пусть \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\) обозначают орбитные биты. Одна строка из данного строкового цикла видит каждый такой бит ровно \(q_j/d_{ij}\) раз, поэтому по модулю \(2\) вклад в уравнения чётности строк появляется тогда и только тогда, когда

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i).$$

Аналогично, один столбец видит каждый орбитный бит ровно \(p_i/d_{ij}\) раз, так что вклад в уравнения чётности столбцов появляется тогда и только тогда, когда

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j).$$

В самих уравнениях участвует не весь набор орбитных переменных, а только сумма по блоку

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2.$$

Остальные \(d_{ij}-1\) направлений внутри блока всегда свободны. Именно поэтому число свободных битов в конце имеет вид \(o(\lambda,\mu)-\rho(\lambda,\mu)\), где \(\rho\) обозначает ранг системы чётности.

Шаг 5: Найти ранг по 2-адическим оценкам

Теперь перейдём на уровень циклов. Вводим по одному узлу уравнений для каждого цикла строк и каждого цикла столбцов. Правила сокращения блочных сумм таковы:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow \text{узел строки и узел столбца должны совпадать},$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow \text{узел строки принудительно равен }0,$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow \text{узел столбца принудительно равен }0.$$

Нужно склеить равные узлы и отметить компоненты, которые обязаны быть нулевыми. Пусть \(f(\lambda,\mu)\) есть число компонент связности, которые не принуждены к нулю. Так как всего узлов \(u+v\), ранг системы равен

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

Следовательно,

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

Шаг 6: Итоговая формула Бёрнсайда

Подставляя веса циклических типов и найденные числа фиксированных матриц, получаем точную формулу, которую вычисляют реализации:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

Тем самым задача превращается из прямого перебора матриц в перебор разбиений числа \(n\) с небольшим расчётом ранга для каждой пары типов.

Разобранный пример: \(n=3\)

Разбиения числа \(3\) таковы: \(1+1+1\), \(2+1\) и \(3\), а числа перестановок соответствующих типов равны \(1\), \(3\) и \(2\). По формуле орбит и по правилам для 2-адических оценок получаем матрицу

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

Например, для пары \((2+1,2+1)\) имеем \(o=2+1+1+1=5\). Среди четырёх узлов уравнений остаются три независимые связи, то есть \(\rho=3\), а значит \(\operatorname{Fix}=2^{5-3}=4\).

Теперь формула Бёрнсайда даёт

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

Этот маленький пример уже показывает всю схему решения: циклические типы задают число орбит, 2-адические оценки задают ранг системы чётности, а лемма Бёрнсайда усредняет числа фиксированных матриц.

Как работает код

Реализации на C++, Python и Java используют один и тот же алгоритм. Сначала перечисляются все разбиения числа \(n\). Для каждого циклического типа сохраняются кратности длин циклов, список 2-адических оценок и число перестановок данного типа. Затем заранее вычисляются все значения \(\gcd(p,q)\) для \(1\le p,q\le n\).

Для каждой пары типов программа сначала находит \(o(\lambda,\mu)\), суммируя нужные значения gcd. После этого строится маленькая структура из отношений равенства и принуждения к нулю, подсчитывается число компонент, не принуждённых к нулю, находится \(\rho(\lambda,\mu)\), и в точный целочисленный аккумулятор добавляется

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

Когда все пары обработаны, суммарное значение делится на \((n!)^2\) по формуле Бёрнсайда, и только затем берётся остаток по модулю \(1001001011\). В версии на C++ внешний обход распараллелен, а версии на Python и Java вычисляют ту же сумму последовательно.

Анализ сложности

Пусть \(p(n)\) обозначает число разбиений числа \(n\). Перечисление всех циклических типов требует \(O(p(n)\,n)\) памяти и примерно того же объёма вспомогательной работы. Главная сумма Бёрнсайда проходит по \(p(n)^2\) парам типов. Для каждой пары и подсчёт орбит, и вычисление ранга требуют не более \(O(n^2)\) элементарных операций, потому что количество циклов строк и столбцов не превосходит \(n\).

Итак, общая временная сложность равна \(O(p(n)^2 n^2)\), а память равна \(O(p(n)\,n)\) плюс точный целочисленный аккумулятор. Для \(n=20\) это вполне практично, несмотря на колоссальный размер исходного пространства матриц.

Сноски и ссылки

  1. Страница задачи: https://projecteuler.net/problem=626
  2. Лемма Бёрнсайда: Wikipedia - Burnside's lemma
  3. Разбиение числа: Wikipedia - Partition (number theory)
  4. \(p\)-адическая оценка: Wikipedia - p-adic valuation
  5. Структура disjoint-set: Wikipedia - Disjoint-set data structure

ملخص المسألة

لتكن \(E_n\) مجموعة المصفوفات الثنائية \(n\times n\) التي تكون فيها مجاميع الصفوف ومجاميع الأعمدة كلها زوجية:

$$E_n=\left\{A=(a_{ij})\in\{0,1\}^{n\times n}:\sum_{j=1}^n a_{ij}\equiv 0\pmod 2,\ \sum_{i=1}^n a_{ij}\equiv 0\pmod 2\right\}.$$

نعد مصفوفتين متكافئتين إذا أمكن الحصول على إحداهما من الأخرى بإعادة ترتيب الصفوف وإعادة ترتيب الأعمدة. نرمز بعدد أصناف التكافؤ هذه بالرمز \(c_n\). المطلوب هو \(c_{20}\bmod 1001001011\).

المنهج الرياضي

الفكرة الأساسية هي عدم عدّ المصفوفات مباشرة، بل عدّ مدارات المصفوفات ذات المجاميع الزوجية تحت تأثير المجموعة \(S_n\times S_n\). لمّة Burnside تختزل المسألة إلى حساب عدد المصفوفات المثبتة بواسطة زوج من التبديلات، وهذا العدد يعتمد فقط على أطوال الدورات في كل تبديل.

الخطوة 1: تطبيق لمّة Burnside

نسخة من \(S_n\) تعمل على الصفوف ونسخة أخرى تعمل على الأعمدة. ويُكتب التأثير المشترك على الصورة

$$A\mapsto P_\sigma A P_\tau^{-1},\qquad (\sigma,\tau)\in S_n\times S_n.$$

ومن ثم

$$c_n=\frac{1}{(n!)^2}\sum_{\sigma\in S_n}\sum_{\tau\in S_n}\operatorname{Fix}(\sigma,\tau),$$

حيث \(\operatorname{Fix}(\sigma,\tau)\) هو عدد المصفوفات في \(E_n\) التي تبقى كما هي بعد تبديل الصفوف وفق \(\sigma\) وتبديل الأعمدة وفق \(\tau\) في الوقت نفسه.

الخطوة 2: استبدال التبديلات بأنماط الدورات

إذا كانت \(\lambda\vdash n\) تجزئة للعدد \(n\)، ولْيكن \(m_k(\lambda)\) عدد الدورات ذات الطول \(k\)، فإن عدد التبديلات التي تمتلك هذا النمط يساوي

$$N(\lambda)=\frac{n!}{\prod_{k\ge 1} k^{m_k(\lambda)}\,m_k(\lambda)!}.$$

إذن مجموع Burnside يعتمد فقط على أزواج التجزئات \((\lambda,\mu)\). فإذا كانت أطوال دورات الصفوف \(p_1,\dots,p_u\) وأطوال دورات الأعمدة \(q_1,\dots,q_v\)، أمكن التعبير عن كل الكميات المطلوبة بدلالة هذه الأطوال وحدها.

الخطوة 3: عدّ مدارات الخلايا لزوج واحد من الأنماط

خذ دورة صفوف طولها \(p_i\) ودورة أعمدة طولها \(q_j\). الكتلة المقابلة تحتوي \(p_iq_j\) خلية، والحركة المشتركة تنقل كل خلية خطوة واحدة على الدورتين معًا. تنقسم هذه الكتلة إلى

$$d_{ij}=\gcd(p_i,q_j)$$

مدارات من الخلايا. لذلك يكون عدد متغيرات المدار الثنائية الكلي

$$o(\lambda,\mu)=\sum_{i=1}^{u}\sum_{j=1}^{v}\gcd(p_i,q_j).$$

ولو لم توجد شروط الزوجية، لكان عدد المصفوفات المثبتة ببساطة هو \(2^{o(\lambda,\mu)}\).

الخطوة 4: تحويل شرط الزوجية إلى نظام خطي فوق \(\mathbb{F}_2\)

داخل كتلة \((p_i,q_j)\)، لنرمز إلى بتات المدارات بالرموز \(x_{ij,1},\dots,x_{ij,d_{ij}}\in\mathbb{F}_2\). الصف الواحد من دورة الصفوف يرى كل بت مدار بعدد مرات يساوي \(q_j/d_{ij}\)، لذا يكون إسهام الكتلة في معادلات زوجية الصفوف غير صفري modulo \(2\) إذا وفقط إذا

$$\frac{q_j}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(q_j)\le \nu_2(p_i).$$

وبالمثل، العمود الواحد من دورة الأعمدة يرى كل بت مدار بعدد مرات يساوي \(p_i/d_{ij}\)، ومن ثم تظهر الكتلة في معادلات زوجية الأعمدة إذا وفقط إذا

$$\frac{p_i}{d_{ij}}\equiv 1\pmod 2 \iff \nu_2(p_i)\le \nu_2(q_j).$$

المهم هنا أن معادلات الزوجية لا تعتمد على كل متغير داخل الكتلة على حدة، بل فقط على المجموع

$$s_{ij}=x_{ij,1}+\cdots+x_{ij,d_{ij}}\pmod 2.$$

أما الاتجاهات الأخرى وعددها \(d_{ij}-1\) داخل الكتلة فتبقى حرّة دائمًا. لهذا يظهر عدد البتات الحرّة في النهاية على الصورة \(o(\lambda,\mu)-\rho(\lambda,\mu)\)، حيث \(\rho\) رتبة نظام الزوجية.

الخطوة 5: استخراج الرتبة من القيم 2-أديّة

ننتقل الآن إلى مستوى الدورات نفسها. نرمز إلى عقدة الصف المقابلة للدورة \(p_i\) بالرمز \(R_i\)، وإلى عقدة العمود المقابلة للدورة \(q_j\) بالرمز \(C_j\). قواعد الإلغاء التي تفرضها مجاميع الكتل هي:

$$\nu_2(p_i)=\nu_2(q_j)\Rightarrow R_i=C_j,$$

$$\nu_2(p_i)\gt \nu_2(q_j)\Rightarrow R_i=0,$$

$$\nu_2(p_i)\lt \nu_2(q_j)\Rightarrow C_j=0.$$

نقوم بدمج العقد المتساوية، ثم نعلّم المكوّنات التي أُجبرت على الصفر. ولتكن \(f(\lambda,\mu)\) عدد المكوّنات المتصلة التي لم تُجبر على الصفر. وبما أن العدد الكلي لعقد المعادلات هو \(u+v\)، فإن رتبة نظام الزوجية تساوي

$$\rho(\lambda,\mu)=u+v-f(\lambda,\mu).$$

وعليه

$$\operatorname{Fix}(\lambda,\mu)=2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.$$

الخطوة 6: صيغة Burnside النهائية

بعد جمع أوزان أنماط الدورات مع أعداد المصفوفات المثبتة نحصل على الصيغة الدقيقة التي تعتمدها التطبيقات:

$$\boxed{c_n=\frac{1}{(n!)^2}\sum_{\lambda\vdash n}\sum_{\mu\vdash n}N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}.}$$

وبذلك تتحول المسألة من تعداد مباشر للمصفوفات إلى تعداد لتجزئات العدد \(n\) مع حساب رتبة صغير لكل زوج من الأنماط.

مثال محلول: \(n=3\)

تجزئات العدد \(3\) هي \(1+1+1\) و\(2+1\) و\(3\)، وأوزانها كأنماط تبديل هي \(1\) و\(3\) و\(2\). وباستخدام حساب المدارات وقواعد القيم 2-أديّة نحصل على مصفوفة الأعداد المثبتة

$$\left[\operatorname{Fix}(\lambda,\mu)\right]= \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix}.$$

على سبيل المثال، في الحالة \((2+1,2+1)\) نجد أن \(o=2+1+1+1=5\). ومن بين أربع عقد معادلات تبقى ثلاث قيود مستقلة فعليًا، لذا \(\rho=3\) وبالتالي \(\operatorname{Fix}=2^{5-3}=4\).

ثم تعطي لمّة Burnside

$$c_3=\frac{1}{36}(1,3,2) \begin{pmatrix} 16 & 4 & 1 \\ 4 & 4 & 1 \\ 1 & 1 & 4 \end{pmatrix} \begin{pmatrix} 1\\ 3\\ 2 \end{pmatrix} =\frac{108}{36}=3.$$

هذا المثال الصغير يوضح الآلية كاملة: نمط الدورات يحدد عدد المدارات، والقيم 2-أديّة تحدد رتبة نظام الزوجية، ولمّة Burnside تحوّل ذلك إلى عدد أصناف التكافؤ النهائي.

كيف تعمل الشيفرة

تتبع تطبيقات C++ وPython وJava المسار نفسه. أولًا تُولِّد جميع تجزئات \(n\)، ولكل نمط دورات تحفظ مضاعفات أطوال الدورات وقائمة القيم 2-أديّة وعدد التبديلات من ذلك النمط. بعد ذلك تُحسب كل قيم \(\gcd(p,q)\) مسبقًا من أجل \(1\le p,q\le n\).

لكل زوج من أنماط الدورات تحسب الشيفرة \(o(\lambda,\mu)\) بجمع قيم gcd المناسبة. ثم تبني البنية الصغيرة المكوّنة من علاقات المساواة والعقد المجبرة على الصفر، وتعدّ المكوّنات غير المجبرة، وتستخرج \(\rho(\lambda,\mu)\)، ثم تضيف

$$N(\lambda)N(\mu)\,2^{\,o(\lambda,\mu)-\rho(\lambda,\mu)}$$

إلى مجمِّع دقيق للأعداد الصحيحة. وبعد الانتهاء من جميع الأزواج يُقسَم المجموع على \((n!)^2\) وفق لمّة Burnside، ثم يُؤخذ الباقي modulo \(1001001011\) في الخطوة الأخيرة فقط. نسخة C++ توازي المرور الخارجي، بينما تنفذ نسختا Python وJava المجموع نفسه بشكل مباشر في عملية واحدة.

تحليل التعقيد

إذا رمزنا لعدد تجزئات \(n\) بالرمز \(p(n)\)، فإن توليد جميع أنماط الدورات يحتاج إلى ذاكرة من رتبة \(O(p(n)\,n)\) وإلى مقدار مشابه من العمل التمهيدي. أما مجموع Burnside الرئيسي فيمر على \(p(n)^2\) من أزواج الأنماط. ولكل زوج، فإن حساب عدد المدارات وحساب رتبة نظام الزوجية يحتاجان في أقصى الأحوال إلى \(O(n^2)\) عملية أولية، لأن عدد دورات الصفوف والأعمدة لا يتجاوز \(n\).

إذن التعقيد الزمني الكلي هو \(O(p(n)^2 n^2)\)، والتعقيد الحجمي هو \(O(p(n)\,n)\) إضافة إلى مجمّع الأعداد الصحيحة الدقيقة. بالنسبة إلى \(n=20\) يبقى هذا النهج عمليًا تمامًا، لأن الفضاء الهائل هو فضاء المصفوفات الخام لا فضاء أنماط الدورات.

حواشٍ ومراجع

  1. صفحة المسألة: https://projecteuler.net/problem=626
  2. لمّة Burnside: Wikipedia - Burnside's lemma
  3. تجزئة عدد صحيح: Wikipedia - Partition (number theory)
  4. التقييم \(p\)-أدي: Wikipedia - p-adic valuation
  5. بنية disjoint-set: Wikipedia - Disjoint-set data structure