<p>You are given two integers <code>m</code> and <code>n</code>. Consider an <code>m x n</code> grid where each cell is initially white. You can paint each cell <strong>red</strong>, <strong>green</strong>, or <strong>blue</strong>. All cells <strong>must</strong> be painted.</p>

<p>Return<em> the number of ways to color the grid with <strong>no two adjacent cells having the same color</strong></em>. Since the answer can be very large, return it <strong>modulo</strong> <code>10<sup>9</sup> + 7</code>.</p>

<p>&nbsp;</p>
<p><strong class="example">Example 1:</strong></p>
<img alt="" src="https://assets.leetcode.com/uploads/2021/06/22/colorthegrid.png" style="width: 200px; height: 50px;" />
<pre>
<strong>Input:</strong> m = 1, n = 1
<strong>Output:</strong> 3
<strong>Explanation:</strong> The three possible colorings are shown in the image above.
</pre>

<p><strong class="example">Example 2:</strong></p>
<img alt="" src="https://assets.leetcode.com/uploads/2021/06/22/copy-of-colorthegrid.png" style="width: 321px; height: 121px;" />
<pre>
<strong>Input:</strong> m = 1, n = 2
<strong>Output:</strong> 6
<strong>Explanation:</strong> The six possible colorings are shown in the image above.
</pre>

<p><strong class="example">Example 3:</strong></p>

<pre>
<strong>Input:</strong> m = 5, n = 5
<strong>Output:</strong> 580986
</pre>

<p>&nbsp;</p>
<p><strong>Constraints:</strong></p>

<ul>
	<li><code>1 &lt;= m &lt;= 5</code></li>
	<li><code>1 &lt;= n &lt;= 1000</code></li>
</ul>

import numpy as np

class Solution:
    def colorTheGrid(self, m: int, n: int) -> int:
        # m x n grid
        #
        # each cell is initially white
        #   > paint each cell r, g, b
        #   > all must be painted
        # Return the number of ways to color the grid with no two adj cells having the same color

        # neat. combinatorics, row/col optimization probably, etc.

        # wlog have to color the first row, which is special b/c there are no constraints on it
        # other than within the same row

        # 3 ways to pick the first cell
        # then 2 diff colors exist for the next
        # and 2 diff after that for the next, etc.
        # for a total of 3*2**(n-1) ways

        # c0  c1  c2  c3  c4  ...
        # c1  ?                   # if we start with c1, we have two choices b/c there's one unique color
        #                         # if we start with the third color, then ? has only one choice
        # this gets complicated, probably there's an optimization but hard to say

        ## so something a bit more brute force: notice that m <= 5 so rows are very short, cheap
        #
        # idea: we can form all of the 3*2**4 = 48 column configurations
        # and the columns they can attach to
        # then we can use a transition matrix approach

        BIG = 10**9+7

        if m == 1:
            return 3 * pow(2, n-1, BIG) % BIG # 3 choices for first, each later has 2 different colors

        # all configs: DFS + bt
        configs = []
        curr = []
        def dfs(r: int) -> None:
            if r == 0:
                for c in 'rgb':
                    curr.append(c)
                    dfs(1)
                    curr.pop()
            elif r == m:
                configs.append(''.join(curr))
            else:
                for c in 'rgb':
                    if c == curr[-1]: continue

                    curr.append(c)
                    dfs(r+1)
                    curr.pop()

        dfs(0)

        # build transition matrix
        C = len(configs)
        T = [[0]*C for _ in range(C)]
        for c1 in range(C):
            for c2 in range(c1):
                if all(a != b for a, b, in zip(configs[c1], configs[c2])):
                    # no conflicts: we can append a column of config c1 to c2 and vice-versa
                    T[c1][c2] = 1
                    T[c2][c1] = 1

        # TODO: use Lapack, etc. to exploit the symmetric nature of the transition matrix
        
        # we advance to 2, 3, etc. ways, n-1 total. If we represent ways_n, the ways to end with
        # each column pattern after n steps as a vector of length C, then
        #    > ways_n = T*ways_{n-1}
        #             = T*T*ways{n-2}
        #             = ...
        #             = T**(n-1) * ways_1
        #
        # We know ways_1 is the ones vector; there is one way to have a single column in each of
        # the single columns
        #
        # So we can use scaling-and-squaring mod BIG to get T**(n-1) in log time
        
        def mult_mod(A, B):
            """Returns A.dot(B) mod BIG"""
            # we need this because A.dot(B) in numpy doesn't have a mod argument
            # as a result we have have up to 48 entries of size BIG-1 multiplied together
            # and 48*(BIG-1)**2 is above the max uint64 value!!

            out = [[0]*C for _ in range(C)]
            for i in range(C):
                for j in range(C):
                    for k in range(C):
                        out[i][j] = (out[i][j] + A[i][k]*B[k][j]) % BIG

            return out

        # repeated scaling and squaring
        n -= 1
        Tn = [[0]*C for _ in range(C)]
        for i in range(C): Tn[i][i] = 1
        P = [t[:] for t in T]
        while n:
            if n & 1:
                Tn = mult_mod(Tn, P)

            # move on to next bit of n
            #     lsb will be 2*lsb, and thus P**(2*lsb) == (P**lsb) ** 2
            P = mult_mod(P, P)
            n >>= 1

        return sum(sum(row) for row in Tn) % BIG
1240