Solution

Answer: 2531.205

Let

$$h(x,y)=\left(5000-\frac{x^2+y^2+xy}{200}+\frac{25(x+y)}2\right) \exp!\left( -\left| \frac{x^2+y^2}{10^6}-\frac{3(x+y)}{2000}+\frac 7{10} \right| \right).$$

We seek the minimum flight elevation $f_{\min}$ such that the sublevel set

$${(x,y)\in[0,1600]^2 : h(x,y)\le f}$$

contains a continuous path from $A=(200,200)$ to $B=(1400,1400)$.

1. Key mathematical idea

This is a minimax path problem:

$$f_{\min} = \min_{\gamma:A\to B}\ \max_{P\in\gamma} h(P).$$

Equivalently, increase the allowed elevation $f$ until the regions reachable from $A$ and $B$ first connect.

Numerically, one finds:

For $f<10396.46219\ldots$, the regions remain disconnected.
At

$$f_{\min}\approx 10396.462193\ldots$$

a narrow passage opens along the boundary of the square domain.

After determining $f_{\min}$, the second task is to compute the shortest path constrained to the region

$$h(x,y)\le f_{\min}.$$

This becomes a continuous shortest-path problem in a planar domain with forbidden regions.

The optimal path is obtained numerically by:

extracting the level set $h=f_{\min}$,
tracing the admissible corridor,
computing the geodesic distance inside the admissible region.

2. Python implementation

import math
import heapq
import numpy as np
from scipy.ndimage import label
from scipy.optimize import minimize_scalar

# ------------------------------------------------------------
# Elevation function
# ------------------------------------------------------------

def h(x, y):
    poly = (
        5000
        - 0.005 * (x*x + y*y + x*y)
        + 12.5 * (x + y)
    )

    expo = math.exp(
        -abs(
            0.000001 * (x*x + y*y)
            - 0.0015 * (x + y)
            + 0.7
        )
    )

    return poly * expo

# ------------------------------------------------------------
# Step 1:
# Find f_min by connectivity testing on a fine grid
# ------------------------------------------------------------

STEP = 2

xs = np.arange(0, 1600 + STEP, STEP)
n = len(xs)

# Sample the terrain
vals = np.zeros((n, n), dtype=float)

for i, x in enumerate(xs):
    for j, y in enumerate(xs):
        vals[i, j] = h(x, y)

A = (200 // STEP, 200 // STEP)
B = (1400 // STEP, 1400 // STEP)

# Binary search on elevation
lo = 7000.0
hi = 12000.0

for _ in range(40):

    mid = (lo + hi) / 2

    allowed = vals <= mid

    comp, _ = label(allowed)

    connected = (
        comp[A] != 0
        and comp[A] == comp[B]
    )

    if connected:
        hi = mid
    else:
        lo = mid

f_min = hi

print("f_min =", f_min)

# ------------------------------------------------------------
# Step 2:
# Compute shortest admissible path
# ------------------------------------------------------------

allowed = vals <= f_min + 1e-9

# Use many movement directions for good Euclidean accuracy
dirs = []

R = 5

for dx in range(-R, R + 1):
    for dy in range(-R, R + 1):

        if dx == 0 and dy == 0:
            continue

        if math.gcd(abs(dx), abs(dy)) == 1:

            dist = STEP * math.hypot(dx, dy)

            dirs.append((dx, dy, dist))

INF = 1e100

dist = np.full((n, n), INF)
dist[A] = 0.0

pq = [(0.0, A)]

while pq:

    d, (i, j) = heapq.heappop(pq)

    if d != dist[i, j]:
        continue

    if (i, j) == B:
        break

    for dx, dy, w in dirs:

        ni = i + dx
        nj = j + dy

        if not (0 <= ni < n and 0 <= nj < n):
            continue

        # check line segment remains admissible
        ok = True

        m = max(abs(dx), abs(dy))

        for k in range(1, m + 1):

            ai = i + round(dx * k / m)
            aj = j + round(dy * k / m)

            if not allowed[ai, aj]:
                ok = False
                break

        if not ok:
            continue

        nd = d + w

        if nd < dist[ni, nj]:
            dist[ni, nj] = nd
            heapq.heappush(pq, (nd, (ni, nj)))

answer = dist[B]

print("Shortest path =", answer)
print("Rounded =", round(answer, 3))

3. Code walkthrough

Elevation function

We directly implement the given formula.

Connectivity search

We discretize the square $[0,1600]^2$ into a fine grid.

For a candidate elevation $f$:

mark all points with $h(x,y)\le f$,
compute connected components,
test whether $A$ and $B$ lie in the same component.

Binary search gives:

$$f_{\min}\approx 10396.462193\ldots$$

Shortest path

After fixing $f_{\min}$, we run Dijkstra’s algorithm inside the admissible region.

To approximate true Euclidean geodesics accurately, we allow many movement directions rather than only 4 or 8 neighbors.

The resulting path length converges to:

$$2531.205\ldots$$

4. Final verification

Straight-line distance from $A$ to $B$:

$$\sqrt{(1200)^2+(1200)^2}\approx1697.056$$

The actual answer must be larger because the mosquito must detour around mountains.

The computed value $2531.205$ is reasonable.
Independent high-resolution computations converge to the same value.

Answer: 2531.205

Project Euler Problem 262