Domino Recurrence Generator

Question

Challenge

We once had a challenge to count domino tilings of m by n grid, and we all know that, for any fixed number of rows, the number of domino tilings by columns forms a linear recurrence. Then why not have a challenge to compute the linear recurrence?!

Let's define \$D_m(n)\$ as the number of domino tilings on a grid of \$m\$ rows and \$n\$ columns. Then the task is: given a single integer \$m \ge 1\$ as input, output the linear recurrence relation for \$D_m(n)\$.

If the relation has order \$k\$ (that is, \$D_m(n+k)\$ depends on \$k\$ previous terms), you need to output the coefficients \$a_i\$ of the recurrence relation

$$ D_m(n+k)=a_{k-1}D_m(n+k-1) + a_{k-2}D_m(n+k-2) + \cdots + a_0 D_m(n) $$

in the order of \$a_0\$ to \$a_{k-1}\$ or the reverse. There are infinitely many correct such relations; you don't need to minimize the order of the relation. But, to ensure that the result is at least minimally useful, the order \$k\$ cannot exceed \$2^m\$ for any input value of \$m\$.

(Side note: An actual sequence is defined only if the initial \$k\$ terms are given along with the recurrence equation. That part is omitted for simplicity of output, and to give incentive to approaches not using the brute-forced terms.)

Note that, for odd \$m\$, every odd-column term will be zero, so you will get a recurrence different from the OEIS entries which strip away zeroes (e.g. 3 rows, 5 rows, 7 rows).

Standard code-golf rules apply. The shortest code in bytes wins.

Examples

Here are the representations from the OEIS, adjusted for odd \$m\$. Initial terms start at \$D_m(0)\$, and the coefficients are presented from \$a_{k-1}\$ to \$a_0\$. Again, your program only needs to output the coefficients. To empirically check the correctness of your output of length \$k\$, plug in the \$k\$ initial terms from the respective OEIS entry, and see if the next \$k\$ terms agree.

m = 1
Initial terms [1, 0]  # D(0) = 1, D(1) = 0
Coefficients  [0, 1]  # D(n+2) = D(n)

m = 2
Initial terms [1, 1]
Coefficients  [1, 1]

m = 3
Initial terms [1, 0, 3, 0]
Coefficients  [0, 4, 0, -1]  # D(n+4) = 4D(n+2) - D(n)

m = 4
Initial terms [1, 1, 5, 11]
Coefficients  [1, 5, 1, -1]

m = 5
Initial terms [1, 0, 8, 0, 95, 0, 1183, 0]
Coefficients  [0, 15, 0, -32, 0, 15, 0, -1]

m = 6
Initial terms [1, 1, 13, 41, 281, 1183, 6728, 31529]
Coefficients  [1, 20, 10, -38, -10, 20, -1, -1]

Possible approaches

There is at least one way to find the recurrence without brute forcing the tilings, outlined below:

1. Compute the transition matrix \$A\$ of \$2^m\$ states, so that the target sequence is in the form of \$D_m(n) = u^T A^n v\$ for some column vectors \$u,v\$.
2. Find the characteristic polynomial or minimal polynomial of \$A\$ as

$$x^k - a_{k-1}x^{k-1} - a_{k-2}x^{k-2} - \cdots - a_0 $$

3. Then the corresponding recurrence relation is

$$s_{n+k} = a_{k-1}s_{n+k-1} + a_{k-2}s_{n+k-2} + \cdots + a_0s_n$$

An example algorithm of computing the minimal polynomial of a matrix can be found on this pdf.

(Of course, you can just brute force the domino tilings for small \$n\$ and plug into a recurrence finder.)

Answer 1

APL (Dyalog Unicode), 72 bytes

{⍵=1:⍳2⋄(2*⌈⍵÷2)(⌷⌹⍉⍤↑)⍉L↑⍉↑,¨+.×\(L←2*⍵)⍴⊂∘.((0∊+)⍱1∊2|×≢¨⍤⊆⊣)⍨⍸1⍴⍨⍵⍴2}

Try it online! (uses a polyfill for ⍤ since TIO is not updated to 18.0 yet)

Requires ⎕pp←2 (implicit rounding of output) and ⎕IO←0 (0-indexing).

We compute a transition matrix, then use the approach listed in S. Białas and M.Białas to find the minimum polynomial and hence the recurrence relation.

Defining the Transition Matrix

Each possible binary fill of a column of \$m\$ cells is one state, so there are \$2^m\$ states.

For m=3, one example state (1 0 0) is

█
▒
▒

The first cell is filled (it is the right side of a horizontal domino sticking out from the previous column), but the second and third cells are empty. This could occur e.g. as the second column in the following tiling (n=4, 3×4 grid):

━━┃┃
┃┃┃┃
┃┃━━

When considering state transitions, we have to be careful to avoid double counting.

My approach is to require full horizontal dominos to be placed whenever possible, then vertical dominos can optionally be placed in the next state's column.

For example, if the current state is 1 0 0:

█
▒
▒

then we force horizontal dominos on the bottom two rows

█▒
━━
━━

so the next state would have to be 0 1 1:

▒
█
█

This rule guarantees the current column to be filled in completely. In addition, it avoids double-counting transitions because it never places vertical dominos in the current column.

Vertical dominos go in the next column. There is no space for vertical dominos in the past example. As an example where vertical dominos can be placed, take the current state to be 1 1 1:

█▒
█▒
█▒

One possibility would be to place no vertical dominos at all, so 1 1 1 → 0 0 0 is a valid state transition.

In addition, a vertical domino can be placed in either of two positions:

█┃        █▒
█┃   or   █┃
█▒        █┃

so 1 1 1 → 1 1 0 and 1 1 1 → 0 1 1 are valid state transitions.

Obtaining the Recurrence from the Transition Matrix.

The paper describes the approach well, but I modified it a little while golfing.

As given, the problem is to find coefficients \$a_i\$ for a given recurrence order \$k\$ such that, for all \$n\$¸

$$a_0 D_m(n) + a_1 D_m(n+1) + \cdots + a_{k-1}D_m(n+k-1) = D_m(n+k)$$

With regards to powers of the transition matrix \$A\$, this can be rewritten as finding coefficients $c_i$ such that

$$c_1 A^1 + c_2 A^2 + \cdots + c_k A^k = A^{k+1}$$

(the paper starts with \$A^0=I_L\$, but that is expensive in terms of bytes)

Let \$L=2*m\$ be the width (and height) of the transition matrix \$A\$. Denoting the entries of \$A^i\$ as \$a_{11}^{(i)}, a_{12}^{(i)}, \ldots a_{LL}^{(i)}\$, we can rewrite the recurrence as \$L^2\$ equations

$$ \begin{align*} c_1 a_{11}^{(1)} + c_2 a_{11}^{(2)} + &\cdots + c_k a_{11}^{(k)} = a_{11}^{(k+1)} \\ c_1 a_{12}^{(1)} + c_2 a_{12}^{(2)} + &\cdots + c_k a_{12}^{(k)} = a_{12}^{(k+1)} \\ &\;\;\,\vdots \\ c_1 a_{LL}^{(1)} + c_2 a_{LL}^{(2)} + &\cdots + c_k a_{LL}^{(k)} = a_{LL}^{(k+1)} \end{align*} $$

Treating \$a_{hi}^{(j)}\$ as constants (since we know \$A\$), this is a system of \$L^2\$ equations in \$k\$ variables \$c_i\$.

Let \$B\$ be the augmented matrix for this massive system of equations for \$k=L\$. Solving the full \$B\$ would give a recurrence of order \$L=2^m\$, but we need a smaller recurrence.

To find a smaller recurrence, we simple use a smaller \$k\$. The bulk of the paper is in proving that the minimum possible \$k\$ is the rank of \$B\$. However, for this specific problem, the minimum \$k\$ is \$k_0=2^{\lceil \frac{m}{2} \rceil}\$ (Source --- floor since the row \$k\$ has \$m=k-1\$). Thus we can take the \$k_0 \times (k_0+1)\$ submatrix at the top left of \$B\$ and solve it to find the \$k_0\$ coefficients of a useful recurrence.

{⍵=1:⍳2⋄(2*⌈⍵÷2)(⌷⌹⍉⍤↑)⍉L↑⍉↑,¨+.×\(L←2*⍵)⍴⊂∘.((0∊+)⍱1∊2|×≢¨⍤⊆⊣)⍨⍸1⍴⍨⍵⍴2}
{...} ⍝ Dfn with right argument m
⍵=1:⍳2⋄ ⍝ Special case m=1: return [0 1]
⍝ Compute the transition matrix A:
 ⍸1⍴⍨⍵⍴2  ⍝ All 2^m states: cartesian m-th power of [0 1]
            ⍝ (m=1 yields a vector of scalars rather than vectors, which is why we need a special case for m=1)
 ∘.{...}⍨ ⍝ Outer product with itself (all current→next state pairs) using function:
  ⍱         ⍝ Neither of the following are true:
  (0∊+)       ⍝ 0→0 in mapping (invalid since we require a horizontal domino when the current state has a 0)
  1∊2|×≢¨⍤⊆⊣  ⍝ Some run of 1→1 has odd length (requires a half vertical domino, impossible)
⍝ Compute the minimal polynomial of A
 +.×\(L←2*⍵)⍴⊂ ⍝ Produce matrix powers of A: A, A*2, ... A*L, where L=2*m
 ↑,¨           ⍝ B: Ravel each (Vec A*k) and join into single (L×L) × L matrix
 ⍉L↑⍉          ⍝ B': Trim to first L rows (for numerical stability in later gauss-jordan elimination)
 (2*⌈⍵÷2) ⍝ Rank r
 ⌷⌹⍉⍤↑ ⍝ Compute recurrence coefficients α←first r entries of b˜÷B̃
  ⍉⍤↑   ⍝ B̃: columns 0 to r-1, inclusive, of B' (taller than B̃ in paper)
  ⌷     ⍝ b˜: r-th column of B' (taller than b˜ of paper)
  ⌹     ⍝ matrix divide b˜÷B̃ to get coefficients

Answer 2

Python 2,  327 ... 249  246 bytes

Saved 37 bytes thanks to fireflame241!

This is using a port of my JS answer to Number of domino tilings to feed SymPy's find_linear_recurrence() method.

import re,sympy,sympy.abc as s
L=[1];N=2**input()-1;a=[0]*N+L;R=range(N+1)
for _ in[0]+R:a=[sum(a[k]*(~k&~i&N<bool(re.match("0b(0*11)*0*$",bin(k&i))))for k in R)for i in R];L+=a[-1:]
print sympy.sequence(L,(s.n,1,N+3)).find_linear_recurrence(N+3)

Try it online! or Run a test suite

How?

State transitions

Given \$n-1\$ rows that are entirely filled and given an \$n\$^th row which is partially filled with state \$S_m(n)\$, we want to find out what are the compatible states \$S_m(n+1)\$ for the next row.

In the example below, we have \$m=5\$ and \$S_5(n)=7\$ (in blue). There are three valid ways of setting the next row while completing the \$n\$^th row. The compatible states \$S_5(n+1)\$ for the next row are \$24\$, \$27\$ and \$30\$.

As a rule of thumb, empty cells in the \$n\$^th row have to be filled with vertical dominoes (in yellow) and we may then insert horizontal dominoes (in green) in the remaining free spaces of the new row.

In the Python code, we use the variables k and i for \$S_m(n)\$ and \$S_m(n+1)\$ respectively.

For the vertical dominoes, we make sure that the bits that are cleared in k are not cleared in i by testing if the following expression evaluates to \$0\$:

~k & ~i & N

where N is a constant bit mask set to \$2^m-1\$.

For the horizontal dominoes, we make sure that the islands of bits that are set in both k and i all include an even number of bits. We do that with a regular expression:

re.match("0b(0*11)*0*$", bin(k & i))

Both tests are combined into:

~k & ~i & N < bool(re.match("0b(0*11)*0*$", bin(k & i)))

Number of valid tilings

The variable a holds a list of \$2^m\$ entries describing how many times each state appeared in the previous iteration. We update a by using the above tests: the new value for a[i] is the sum of all previous values a[k] for all pairs of compatible states (k,i):

a = [sum(a[k] * (...) for k in R) for i in R]

The total number of valid tilings is the number of times we reach the 'full' state (\$2^m-1\$) for the last row, which is a[-1].

Final solution

We use this method to compute the first \$2^m+2\$ terms of the sequence in the list L and inject it into find_linear_recurrence() to get the final solution.

Note: According to OEIS (and as already pointed out by fireflame241), computing \$2^{\lceil m/2\rceil}\$ terms would be enough and would make the code faster, but also a bit longer.

Answer 3

Python 3.8, 228 bytes

Similarly to Arnauld's answer, this uses my answer to Number of domino tilings to feed SymPy's find_linear_recurrence function.

from math import*
import sympy,sympy.abc as s
def f(m):N=2**m+2;return sympy.sequence([round(abs(prod(2*cos((i//k+1)*pi/-~m)+2j*cos((i%k+1)*pi/-~k)for i in range(m*k)))**.5)for k in range(N)],(s.n,1,N)).find_linear_recurrence(N)

Try it online!. TIO doesn't have sympy in its Python 3.8 installation, so the link includes a polyfill of math.prod, which is new to 3.8.

Since this multiplies (floating-point) complex numbers together, it loses precision for \$m\geq 5\$, leading to a completely incorrect recurrence.