Problem 380: Amazing Mazes!

Detailed mathematical approach

Problem Summary

An \(m \times n\) maze is a rectangular grid of \(m \times n\) square cells with walls between neighbouring cells placed so that there is exactly one path from the top-left cell to every other cell. Let \(C(m,n)\) count the distinct \(m \times n\) mazes; two that differ only by a rotation or a reflection still count as two. ^[1] The statement gives four values to test an implementation against: \(C(1,1)=1\), \(C(2,2)=4\), \(C(3,4)=2415\) and \(C(9,12)=2.5720\times10^{46}\) to five significant digits. We must find \(C(100,500)\) in the same form, with a lowercase e between mantissa and exponent.

The whole solution in one sentence: a maze is exactly a spanning tree of the grid graph whose vertices are the cells, Kirchhoff's matrix–tree theorem turns the number of spanning trees into a product of \(mn-1\) explicit cosine expressions, and since that product is far too large for a machine number we add the base-ten logarithms of its factors instead.

Mathematical Approach

1. A maze is a spanning tree

Turn the picture into a graph. The vertices are the \(mn\) cells, and we draw an edge between two neighbouring cells whenever the wall between them is absent, so one can walk across. There are \(m(n-1)\) wall positions between horizontal neighbours and \((m-1)n\) between vertical ones, so a maze is a choice of which of the

$$E=m(n-1)+(m-1)n=2mn-m-n$$

wall positions to open. Call that graph \(T\).

Every cell must be reachable from the top-left cell \(r\) by exactly one path. "At least one path" means \(T\) is connected. "At most one path" forbids cycles: let \(v\) be a vertex of a cycle at minimum distance from \(r\). A shortest path \(P\) from \(r\) to \(v\) cannot meet the cycle again, since such a point would be nearer to \(r\) than \(v\). With \(a,b\) the neighbours of \(v\) on the cycle, \(P\) plus the edge \(va\) and \(P\) plus the arc \(v,b,\dots,a\) are two different paths to \(a\) — a contradiction. Conversely, in a connected acyclic graph any two vertices are joined by exactly one path, as two would enclose a cycle.

So the mazes are precisely the spanning trees of the grid graph: connected cycle-free subgraphs using all \(mn\) cells, hence exactly \(mn-1\) open passages out of \(E\).

2. The smallest cases by hand

A single cell has no internal wall positions, so \(C(1,1)=1\). A row of \(n\) cells has \(n-1\) wall positions, all of which must be open, so \(C(1,n)=1\). For \(m=n=2\) the four cells and \(E=4\) wall positions form one square cycle; a tree keeps \(3\) passages open, so we merely choose the wall that stays, and \(C(2,2)=4\), as the statement says.

For \(2\times3\) we have \(E=7\) and a tree uses \(5\) edges, giving \(\binom75=21\) candidates. A cycle inside five edges must be one of the two unit squares, since the outer cycle needs six edges. Each square uses four edges and the fifth may be any of the three others; no five edges hold both squares, as together they already use all seven. So exactly \(2\cdot3=6\) candidates have a cycle and \(C(2,3)=21-6=15\). Enumeration gives \(C(3,3)=192\) and, over the \(\binom{17}{11}=12376\) candidates of the \(3\times4\) grid, \(C(3,4)=2415\) — the third checkpoint.

3. Why enumeration stops here

For \(9\times12\) one would have to find the trees among \(\binom{195}{107}\) subsets, and \(100\times500\) has \(E=99400\) wall positions and \(49999\) open passages. Since \(C(9,12)\) alone is a \(47\)-digit number, no method that touches one maze at a time can work. We need a closed formula.

4. Kirchhoff's matrix–tree theorem

The Laplacian of a graph on \(N\) vertices is \(L=D-A\), with \(D\) diagonal carrying the vertex degrees and \(A\) the adjacency matrix. ^[3] Every row and column of \(L\) sums to zero, so \(L\) is singular and the all-ones vector lies in its kernel. Kirchhoff's matrix–tree theorem says the number \(\tau(G)\) of spanning trees of a connected \(G\) equals any cofactor of \(L\). ^[2] Deleting row \(a\) and column \(b\) leaves a smaller matrix whose determinant must still be multiplied by \((-1)^{a+b}\); deleting the row and column of the same vertex makes that sign \(+1\), so there the count is the determinant itself.

Check it on the \(2\times2\) maze. All four cells have degree \(2\); deleting the row and column of the bottom-right cell leaves

$$\begin{pmatrix}2&-1&-1\\-1&2&0\\-1&0&2\end{pmatrix},$$

whose determinant is \(2\cdot4-(-1)\cdot(-2)+(-1)\cdot2=4\), in agreement with \(C(2,2)=4\).

5. From the cofactor to the eigenvalues

A determinant of size \(mn-1\) is no use for \(m=100\), \(n=500\), but eigenvalues help. An eigenvector of \(L\) is a nonzero \(v\) that \(L\) merely rescales, \(Lv=\lambda v\); \(\lambda\) is its eigenvalue. For the two-cell path \(L_2\) sends \((1,1)\) to \((0,0)\) and \((1,-1)\) to \((2,-2)\), so its eigenvalues are \(0\) and \(2\). In general \(L\) is symmetric, so all \(N\) eigenvalues are real, and \(x^{\mathsf T}Lx\) is the sum of \((x_u-x_v)^2\) over the edges, so none is negative; on a connected graph that sum vanishes only for constant vectors, so exactly one is \(0\). Order them \(0=\lambda_0\lt\lambda_1\le\dots\le\lambda_{N-1}\) and consider

$$\begin{aligned}p(t)&=\det(tI+L)\\&=t\prod_{k=1}^{N-1}(t+\lambda_k).\end{aligned}$$

Expanding the determinant, the coefficient of \(t\) collects the terms taking \(t\) from one diagonal entry and \(L\) elsewhere, so it is the sum of the \(N\) principal minors of order \(N-1\) — the \(N\) principal cofactors, each equal to \(\tau(G)\). On the right that coefficient is \(\lambda_1\cdots\lambda_{N-1}\). Equating them, ^[2]

$$\tau(G)=\frac{1}{N}\prod_{k=1}^{N-1}\lambda_k.$$

For the \(2\times2\) maze \(L\) is that of a four-cycle, with eigenvalues \(0,2,2,4\), and indeed \(2\cdot2\cdot4/4=4\). It remains to find all \(mn\) grid eigenvalues in closed form.

6. The grid is a product of two paths

Index the cells by \((r,c)\) with \(0\le r\lt m\) and \(0\le c\lt n\). Two cells are neighbours exactly when one coordinate agrees and the other differs by \(1\); that is the Cartesian product \(P_m\,\square\,P_n\) of the paths on \(m\) and on \(n\) vertices, and the Cartesian product of two paths is precisely a grid graph. ^[4]

The degree of \((r,c)\) is the degree of \(r\) in \(P_m\) plus that of \(c\) in \(P_n\), and each neighbour changes exactly one coordinate, so \(L\) is the Kronecker sum ^[6]

$$L=L_m\otimes I_n+I_m\otimes L_n,$$

with \(L_m,L_n\) the path Laplacians and \(I_k\) the identity of size \(k\). For vectors \(u,w\) let \(u\otimes w\) be the vector giving cell \((r,c)\) the value \(u_rw_c\). If \(L_mu=\alpha u\) and \(L_nw=\beta w\), the first term acts only on the row index and returns \(\alpha\,u\otimes w\), the second only on the column index and returns \(\beta\,u\otimes w\), so

$$L(u\otimes w)=(\alpha+\beta)\,(u\otimes w).$$

Taking eigenvector bases \(u_0,\dots,u_{m-1}\) and \(w_0,\dots,w_{n-1}\), the \(mn\) products \(u_i\otimes w_j\) form a basis of the whole space, so the grid eigenvalues are exactly the sums \(\alpha_i+\beta_j\). We confirmed this for the \(2\times3\), \(3\times4\) and \(4\times5\) grids: the sorted spectra match the sorted sums to within \(4\cdot10^{-15}\).

7. The spectrum of a path

For \(m=1\) the path has no edges, \(L_1=(0)\) and the only eigenvalue is \(0\). For \(m\ge2\), \(L_m\) is tridiagonal with \(-1\) next to the diagonal and endpoint degrees \(1\), interior degrees \(2\) on it. Try the shifted cosine

$$v(x)=\cos\bigl(\theta(x+\tfrac12)\bigr),\quad 0\le x\lt m.$$

An interior row computes \(2v(x)-v(x-1)-v(x+1)\), and \(\cos(A-\theta)+\cos(A+\theta)=2\cos\theta\cos A\) with \(A=\theta(x+\tfrac12)\) gives \(v(x-1)+v(x+1)=2\cos\theta\,v(x)\), so that row returns \((2-2\cos\theta)v(x)\). The first row returns only \(v(0)-v(1)\), and the half-step shift makes that work too: the same identity with \(A=\theta/2\) and \(\cos(-\theta/2)=\cos(\theta/2)\) gives \(v(0)+v(1)=2\cos\theta\,v(0)\), as if a mirror cell with \(v(-1)=v(0)\) sat outside the left end. On the right the mirror needs \(v(m)=v(m-1)\), which holds exactly when \(\theta m\) is a multiple of \(\pi\), as \(\cos\) is even about such points. So \(\theta=\pi k/m\) and

$$\alpha_k=2-2\cos\frac{\pi k}{m}=4\sin^2\frac{\pi k}{2m}$$

for \(k=0,1,\dots,m-1\). These \(m\) numbers are distinct, since \(\cos\) strictly decreases on \([0,\pi]\) and the arguments lie in \([0,\pi)\); an \(m \times m\) matrix with \(m\) distinct eigenvalues has no others, so this is the full spectrum, with \(\alpha_0=0\) for the constant vector. The same eigenvalues and shifted cosines are the standard ones for the discrete second derivative with Neumann boundary conditions. ^[5]

8. The closed formula

With \(\beta_j=2-2\cos(\pi j/n)\), Sections 5 to 7 give the grid eigenvalues

$$\lambda_{i,j}=4-2\cos\frac{\pi i}{m}-2\cos\frac{\pi j}{n},$$

for \(0\le i\lt m\) and \(0\le j\lt n\). Since \(\alpha_i=0\) only at \(i=0\) and \(\beta_j=0\) only at \(j=0\), the pair \((0,0)\) gives the single zero eigenvalue — which also proves the grid connected — and the rest are positive. With \(N=mn\),

$$C(m,n)=\frac{1}{mn}\prod_{(i,j)\ne(0,0)}\lambda_{i,j}.$$

For \(m=n=1\) the product is empty, hence \(1\), and the formula still returns \(C(1,1)=1\). It is symmetric in \(m\) and \(n\), as it must be, since transposing a maze gives a maze. For square grids it is the product recorded in OEIS A007341, whose terms \(1,4,192,100352,557568000\) and \(32565539635200\) we reproduced exactly. ^[7]

9. A worked example: \(C(3,4)=2415\)

With \(m=3\) the path eigenvalues are \(\alpha=(0,1,3)\), because \(2-2\cos(\pi/3)=1\) and \(2-2\cos(2\pi/3)=3\). With \(n=4\) they are \(\beta=(0,\,2-\sqrt2,\,2,\,2+\sqrt2)\). Multiplying the four entries of each row of the table \(\alpha_i+\beta_j\), and dropping the one zero factor:

• \(i=0\): \((2-\sqrt2)\cdot2\cdot(2+\sqrt2)\)\({}=2(4-2)=4\);
• \(i=1\): \(1\cdot(3-\sqrt2)\cdot3\cdot(3+\sqrt2)\)\({}=3(9-2)=21\);
• \(i=2\): \(3\cdot(5-\sqrt2)\cdot5\cdot(5+\sqrt2)\)\({}=15(25-2)=345\).

The product is \(4\cdot21\cdot345=28980\), and dividing by \(mn=12\) gives \(2415\), the checkpoint of the statement, from nothing but conjugate pairs that multiply out to integers. That cancellation is why an irrational-looking product lands on a whole number.

10. Working in logarithms

The product cannot be formed directly. \(C(100,500)\) has more than twenty thousand decimal digits — choosing a maze independently in each of the \(500\) disjoint \(10\times10\) blocks and joining them by fixed passages already gives at least \(C(10,10)^{500}\) — so it overflows the floating-point types these programs use. Nor is it a product of integers: for \(m=100\), \(n=500\) all but \(105\) of the \(49999\) factors are irrational, so an exact evaluation would have to return to the cofactor of Section 4. Every factor is positive, so we add logarithms:

$$\begin{aligned}\mathcal L&=\log_{10}C(m,n)\\&=\sum_{(i,j)\ne(0,0)}\log_{10}\lambda_{i,j}\\&\quad-\log_{10}(mn).\end{aligned}$$

The requested format then follows: put \(e=\lfloor\mathcal L\rfloor\) and \(\mu=10^{\mathcal L-e}\), so \(1\le\mu\lt10\) and \(C=\mu\cdot10^{e}\); rounding \(\mu\) to four decimals gives five significant digits. One case needs care: if the rounding carries \(\mu\) to \(10.0000\), the mantissa must be divided by ten and the exponent raised, which is why the code tests for it. For \(9\times12\) the sum of \(107\) logarithms yields mantissa \(2.5720\) and exponent \(46\).

11. How accurate must the sum be?

If \(\widehat{\mathcal L}\) approximates \(\mathcal L\) with the same floor, the mantissas obey \(|\widehat\mu-\mu|=\mu\,|10^{\delta}-1|\) with \(\delta=\widehat{\mathcal L}-\mathcal L\), hence at most \(10(10^{|\delta|}-1)\). Taking

$$|\delta|\le2\cdot10^{-6}$$

holds the mantissa error below \(4.61\cdot10^{-5}\), just under half a unit in the fifth digit, whatever the mantissa turns out to be. That is still not the same as printing the correctly rounded five digits: two mantissas on opposite sides of a rounding midpoint can be far closer than that, as \(1.23454999\) and \(1.23455001\) are. A budget must therefore be paired with a margin to the nearest midpoint, and the margin measured below is itself smaller than \(4.61\cdot10^{-5}\), so this loose budget alone would not settle the last digit. Two things can spoil \(\mathcal L\).

First, the factors. When \(i\) and \(j\) are small both cosines are close to \(1\), so the form of Section 8 subtracts two quantities close to \(2\) from \(4\) and digits cancel. The smallest factor used for \(m=100\), \(n=500\) is \(\lambda_{0,1}\approx3.95\cdot10^{-5}\); in double precision it carries a relative error of about \(7\cdot10^{-13}\) and contributes some \(3\cdot10^{-13}\) to \(\mathcal L\). The form \(4\sin^2(\pi i/2m)+4\sin^2(\pi j/2n)\) of Section 7 avoids the cancellation entirely, and summing with it moves \(\mathcal L\) by under \(10^{-10}\): the cancellation is harmless here.

Second, the accumulation. Let \(A\) be the sum of the absolute values of the \(T\) terms. Sequential addition has an absolute error of order \(TuA\), with \(u\approx1.1\cdot10^{-16}\) the unit roundoff, while Kahan's compensated summation cuts the leading part to order \(uA\), leaving a term of order \(Tu^2A\). ^[9] A relative error also brings in the ratio of \(A\) to the exact sum, so cancellation matters; here \(T=49999\), exactly \(4361\) terms are negative, \(A\) is of order \(10^{4}\) and that ratio is only \(1.15\). Against a \(50\)-digit recomputation the plain double sum is off by \(1.2\cdot10^{-10}\), the compensated one by \(3.8\cdot10^{-12}\).

12. A sanity check on the size

A coarse independent check on the scale comes from the square-grid asymptotic: with \(\log\) the natural logarithm,

$$\lim_{k\to\infty}\frac{\log C(k,k)}{k^2}=\frac{4G}{\pi},$$

where \(G=0.9159655942\dots\) is Catalan's constant, so the limit is \(1.166243616\dots\), the spanning-tree constant of the square lattice. ^{[7] [8]} Our formula gives per-cell logarithms \(0.9391\), \(1.0530\), \(1.1101\) and \(1.1383\) for \(k=8,16,32,64\): a slow approach, matching the boundary correction in that asymptotic, since in our samples \(\bigl(1.1662436\,mn-\log C\bigr)\)\({}/(m+n)\) stayed between \(0.892\) and \(0.906\) from \(16\times16\) to \(20\times100\). This fixes the order of magnitude, but certifies no digits.

How the Code Works

The C++ program keeps \(m=100\), \(n=500\) and a flag in an Options record. parse_arguments, with the helper parse_int_after_prefix, accepts --m=, --n= and --skip-checkpoints, rejects anything else and requires \(m,n\ge1\) after parsing, so one binary reproduces any checkpoint. The helper accumulates digits into an int with no overflow test, so a dimension must fit that type; an oversized argument is not safely rejected.

compute_log10_spanning_trees is Sections 8 and 10 in code. It fills cos_m and cos_n with the \(m+n\) needed cosines, loops over all pairs \((i,j)\), skips \((0,0)\), forms eigenvalue as \(4-2c_i-2c_j\) and adds its logarithm to sum_log10. The lines with y, t and compensation are Kahan's summation from Section 11: compensation carries the part of each term the running sum could not hold and feeds it back. Subtracting \(\log_{10}(mn)\) supplies the \(1/(mn)\) of the theorem.

format_scientific_from_log10 is the last step of Section 10: exponent is the floor of \(\mathcal L\), mantissa is \(10^{\mathcal L-e}\), scale is \(10^{4}\), rounded is the rounded scaled mantissa, and the test rounded >= 10 * scale handles the carry. An ostringstream with std::fixed and precision \(4\) prints mantissa, e and exponent.

Unless --skip-checkpoints is given, main calls run_checkpoints first: \(C(1,1)\) must have logarithm \(0\); exponentiating the logarithms of \(C(2,2)\) and \(C(3,4)\) and applying llround must give \(4\) and \(2415\); \(C(9,12)\) must format to the string of the statement; and the logarithms of \(C(3,4)\) and \(C(4,3)\) must agree under is_nearly_equal, a relative comparison. A failed check exits with status \(2\) and prints no result. Passing them certifies nothing about the requested grid or its rounding.

The cosine tables, logarithms, compensated sum and mantissa are long double, while dimensions and indices are int, the exponent long long, and the rounded mantissa is converted to double before streaming. How much precision long double carries is a property of the compiler and platform, not of the source: on the build we measured, std::numeric_limits reports the same \(53\)-bit significand as double, so there the accuracy edge over the other two comes from the Kahan compensation alone.

Python runs the same loops with both cosine tables cached, accumulates with a plain +=, and formats with round and an f-string. Java keeps ci outside the inner loop, computing it \(m\) times, but calls Math.cos(pi * j / n) for each of the \(mn-1\) counted pairs: \(mn+m-1\) cosine calls and no table. Its rounded mantissa is a long, and the carry rounded /= 10 is integer division, exact because the carry can only fire at \(100000\). Both fix \(m,n=100,500\) with no option to change them and run no checkpoints, so neither reproduces a smaller grid unedited. One caveat: Java's String.format uses the default formatting locale, so where the decimal separator is a comma the mantissa gets a comma, while the required format needs a point.

Complexity and Verification

The loop costs one logarithm per counted pair, so the time is \(\Theta(mn)\): for \(100\times500\) that is \(49999\) logarithms plus one for \(mn\), with \(m+n\) cosines in C++ and Python or \(mn+m-1\) in Java. Memory is \(O(m+n)\) for the cosine tables, or \(O(1)\) in Java. Measured on one machine, the loop takes well under a millisecond in C++ and about four milliseconds in Python. The cofactor route of Section 4 would face a matrix of order \(N-1=49999\): dense fraction-free elimination costs \(O(N^3)\) integer operations and \(O(N^2)\) stored entries, and since the entries grow that is not even a bit-complexity bound. Sparse methods do better, but the explicit spectrum avoids the matrix entirely.

Beyond those checkpoints, the following was verified for this article with code outside the published solutions. A literal brute force over all \(2^{E}\) wall configurations, keeping those in which every cell is reached from the top-left cell by exactly one path, gives \(1,1,4,15\) and \(192\) for \(1\times1\), \(1\times4\), \(2\times2\), \(2\times3\) and \(3\times3\): Section 1 is thus tested against the problem's wording, not assumed. Exact integer cofactors by fraction-free elimination give \(1,4,15,56,192,2415\) and \(100352\) for \(1\times1\), \(2\times2\), \(2\times3\), \(2\times4\), \(3\times3\), \(3\times4\) and \(4\times4\), and the cosine formula at \(60\) digits rounds to the right integer for all \(24\) pairs with \(m\le4\), \(n\le6\).

For larger grids the exact cofactors are \(C(5,7)=10021992194369\), the \(47\)-digit \(C(9,12)=25720403635745\dots\) — that is \(2.5720\times10^{46}\), as the statement says — and \(C(10,10)\), \(C(12,12)\) with \(43\) and \(64\) digits; the formula matched every digit compared. ^[7] With --m= and --n= the C++ program prints \(2.4150\mathrm{e}3\), \(1.0022\mathrm{e}13\) and \(5.6943\mathrm{e}42\) for \(4\times3\), \(5\times7\) and \(10\times10\); Python and Java carry no such options, so their loops were checked in separate small-parameter harnesses, which reproduce those strings under a decimal-point locale.

Finally, recomputing the \(49999\)-term sum in \(50\)-digit arithmetic confirms the error figures of Section 11. In units of the fifth digit, the plain double sum's error of \(1.2\cdot10^{-10}\) in \(\mathcal L\) moves the mantissa by at most \(2.8\cdot10^{-5}\) of a unit, while the true mantissa sits \(0.13\) of a unit from the nearest rounding midpoint — nearly \(5\cdot10^{3}\) times as far. That, and not the loose budget, is what certifies the five printed digits.

Footnotes and References

The references below supply the background results used above. The maze–tree correspondence, the spectral product and the path eigenvalues are derived in this article; Section 12 uses a cited asymptotic as a sanity check.

1. Project Euler 380 — Amazing Mazes!. The official statement defines an m × n maze as a grid of cells with walls between cells such that there is exactly one path from the top-left square to any other square, counts rotations and reflections of a maze as distinct, records C(1,1) = 1, C(2,2) = 4, C(3,4) = 2415 and C(9,12) = 2.5720e46, and asks for C(100,500) to five significant digits.
2. Wikipedia — Kirchhoff's theorem. Gives the matrix–tree theorem as "the number of spanning trees equals any cofactor of the Laplacian", with the cofactor sign, and the spectral form (1/n)·λ₁λ₂⋯λₙ₋₁ for Laplacian eigenvalues 0 = λ₀ ≤ … ≤ λₙ₋₁. Sections 4 and 5 use both.
3. Wikipedia — Laplacian matrix. Defines the Laplacian as degree matrix minus adjacency matrix, notes the zero row and column sums, and states that the multiplicity of the eigenvalue 0 is the number of connected components — the fact that turns the single zero factor of Section 8 into a proof of connectedness.
4. Wikipedia — Cartesian product of graphs. Defines the product on vertex pairs, adjacent when one coordinate is equal and the other adjacent, states that the Cartesian product of two path graphs is a grid graph, and gives its adjacency matrix as the Kronecker sum A₁ ⊗ I + I ⊗ A₂. Section 6 rests on these.
5. Wikipedia — Eigenvalues and eigenvectors of the second derivative. For the discrete second derivative with pure Neumann boundary conditions on n points it lists eigenvalues −(4/h²)·sin²(π(j−1)/(2n)) and eigenvectors cos(π(j−1)(i−0.5)/n): with h = 1 and the opposite sign, the path Laplacian spectrum and shifted cosines of Section 7.
6. Wikipedia — Kronecker product. Defines the Kronecker product and the Kronecker sum A ⊗ Iₘ + Iₙ ⊗ B, distinguishing it from the direct sum. It gives the spectrum only of the plain product, so the eigenvalue claim of Section 6 is proved in the article.
7. OEIS A007341 — Number of spanning trees in n × n grid. Lists 1, 4, 192, 100352, 557568000, 32565539635200, … and a product formula with prefactor 2^(n²−1)/n² over the factors 2 − cos(π·n₁/n) − cos(π·n₂/n) for (n₁, n₂) ≠ (0,0) — the formula of Section 8 after pulling a 2 out of each of the n² − 1 terms. Its asymptotic entry contains exp(4Gn²/π), the limit used in Section 12.
8. OEIS A218387 — Spanning tree constant of the square lattice. Gives the expansion 1.16624361612327512… and identifies the constant as Catalan's constant times 4/π, the per-cell growth rate of Section 12.
9. Wikipedia — Kahan summation algorithm. Describes the compensated summation used by the C++ solution and bounds its error by the sum of the absolute values of the inputs; sequential summation has a worst-case error growing proportionally to n, and compensated summation a relative error effectively independent of n for a fixed condition number. Section 11 uses these bounds.