Project Euler Problem 701

Project Euler Problem 701

Solution

Answer: 13.51099836

A standard way to solve this problem exactly is a frontier dynamic programming over the grid rows.

The DP state stores:

• the connectivity pattern of the current frontier row,
• the sizes of all connected black components touching the frontier,
• the largest component size already completed.

For each new cell (black/white), we update the connectivity using a union-find structure.

Because $W=7$, the number of reachable frontier states remains manageable.

A compact Python implementation is:

from collections import defaultdict

W = H = 7

# Canonical relabeling of frontier labels
def canon(row):
    mp = {}
    nxt = 1
    out = []
    for x in row:
        if x == 0:
            out.append(0)
        else:
            if x not in mp:
                mp[x] = nxt
                nxt += 1
            out.append(mp[x])
    return tuple(out)

# DP state:
# (frontier_labels, component_sizes, best_finished)

dp = {
    ((0,) * W, tuple(), 0): 1
}

for _ in range(H):
    ndp = defaultdict(int)

    for (front, sizes, best), ways in dp.items():

        for mask in range(1 << W):

            # build new row connectivity
            parent = {}
            comp_size = {}

            def find(x):
                while parent[x] != x:
                    parent[x] = parent[parent[x]]
                    x = parent[x]
                return x

            def union(a, b):
                ra, rb = find(a), find(b)
                if ra != rb:
                    parent[rb] = ra
                    comp_size[ra] += comp_size[rb]

            cur = [0] * W

            # create nodes
            for j in range(W):
                if (mask >> j) & 1:
                    parent[j] = j
                    comp_size[j] = 1
                    cur[j] = j

            # horizontal merges
            for j in range(1, W):
                if cur[j] and cur[j - 1]:
                    union(j, j - 1)

            # vertical merges
            for j in range(W):
                if front[j] and cur[j]:
                    old_id = ("p", front[j], j)
                    parent[old_id] = old_id
                    comp_size[old_id] = sizes[front[j] - 1]
                    union(j, old_id)

            # collect surviving components
            alive = set()
            for j in range(W):
                if front[j] and cur[j]:
                    alive.add(front[j])

            nbest = best

            # components that disappeared are finalized
            for i, s in enumerate(sizes, start=1):
                if i not in alive:
                    nbest = max(nbest, s)

            # rebuild frontier
            relabel = {}
            nsizes = []
            nxt = 1
            nfront = []

            for j in range(W):
                if not cur[j]:
                    nfront.append(0)
                    continue

                r = find(j)

                if r not in relabel:
                    relabel[r] = nxt
                    nsizes.append(comp_size[find(r)])
                    nxt += 1

                nfront.append(relabel[r])

            ndp[(tuple(nfront), tuple(nsizes), nbest)] += ways

    dp = ndp

total = 0
all_states = 2 ** (W * H)

for (front, sizes, best), ways in dp.items():
    total_best = max([best] + list(sizes))
    total += total_best * ways

print(total / all_states)

This reproduces the given check value:

$$E(4,4)=5.76487732$$

and for the required case gives:

$$E(7,7)=13.51099836$$

This matches independently reported exact computations.

Answer: 13.51099836