Project Euler

Convex Holes

Given 500 points in the plane generated by the linear congruential generator: s_0 = 290797 s_(n+1) = s_n^2 mod 50515093 t_n = (s_n mod 2000) - 1000 Points are P_k = (t_2k, t_(2k+1)) for k = 0, 1,.....

Source sync Apr 19, 2026

Problem #0252

Level Level 15

Solved By 1,022

Languages C++, Python

Answer 104924.0

Length 516 words

geometrydynamic_programminglinear_algebra

Given a set of points on a plane, we define a convex hole to be a convex polygon having as vertices any of the given points and not containing any of the given points in its interior (in addition to the vertices, other given points may lie on the perimeter of the polygon).

As an example, the image below shows a set of twenty points and a few such convex holes. The convex hole shown as a red heptagon has an area equal to $1049694.5$ square units, which is the highest possible area for a convex hole on the given set of points.

For our example, we used the first $20$ points $(T_{2k - 1}, T_{2k})$, for $k = 1,2,\dots ,20$, produced with the pseudo-random number generator: \begin {align*} S_0 &= 290797\\ S_{n + 1} &= S_n^2 \bmod 50515093\\ T_n &= (S_n \bmod 2000) - 1000 \end {align*}

i.e. $(527, 144), (-488, 732), (-454, -947), \dots $

What is the maximum area for a convex hole on the set containing the first $500$ points in the pseudo-random sequence?

Specify your answer including one digit after the decimal point.

Problem 252: Convex Holes

Formal Mathematical Development

Definition 1. Let $S \subset \mathbb{R}^2$ be a finite point set. A convex hole in $S$ is a convex polygon $P$ such that $\operatorname{Vert}(P) \subseteq S$ and $\operatorname{int}(P) \cap S = \emptyset$ , where $\operatorname{int}(P)$ denotes the open interior of $P$ .

Definition 2 (Signed Area). For three points $A, B, C \in \mathbb{R}^2$ , the signed area function is:

\sigma(A, B, C) = (B_x - A_x)(C_y - A_y) - (C_x - A_x)(B_y - A_y).

The value $\sigma(A, B, C) > 0$ if and only if the ordered triple $(A, B, C)$ is oriented counterclockwise.

Lemma 1 (Point-in-Triangle Test). A point $Q$ lies strictly inside triangle $\triangle P_1 P_2 P_3$ if and only if

\operatorname{sgn}\bigl(\sigma(P_1, P_2, Q)\bigr) = \operatorname{sgn}\bigl(\sigma(P_2, P_3, Q)\bigr) = \operatorname{sgn}\bigl(\sigma(P_3, P_1, Q)\bigr) \ne 0.

Proof. The barycentric coordinates of $Q$ with respect to $\triangle P_1 P_2 P_3$ are $\lambda_i = \sigma_i / \sigma$ , where $\sigma = \sigma(P_1, P_2, P_3)$ and $\sigma_1 = \sigma(P_2, P_3, Q)$ , $\sigma_2 = \sigma(P_3, P_1, Q)$ , $\sigma_3 = \sigma(P_1, P_2, Q)$ . The point $Q$ lies in the open interior if and only if all $\lambda_i > 0$ , which holds precisely when all three signed areas share the sign of $\sigma$ . $\blacksquare$

Theorem 1 (Fan Triangulation Emptiness). Let $P$ be a convex polygon with vertices $v_0, v_1, \ldots, v_{k-1}$ in counterclockwise order, and let $S$ be a finite point set with $\operatorname{Vert}(P) \subseteq S$ . Then $\operatorname{int}(P) \cap S = \emptyset$ if and only if for every $i \in \{1, \ldots, k-2\}$ , the triangle $\triangle v_0 v_i v_{i+1}$ contains no point of $S$ in its interior.

Proof. The fan triangulation $\{\triangle v_0 v_i v_{i+1}\}_{i=1}^{k-2}$ partitions the interior of $P$ (since $P$ is convex, the fan from any vertex is a valid triangulation). A point $Q \in \operatorname{int}(P)$ must lie in the interior of exactly one triangle of this fan. Hence $\operatorname{int}(P) \cap S = \emptyset$ if and only if $\operatorname{int}(\triangle v_0 v_i v_{i+1}) \cap S = \emptyset$ for all $i$ . $\blacksquare$

Theorem 2 (Anchor-Based Convex Hole DP). Fix an anchor point $p_0 \in S$ . Sort the remaining points $\{p_1, \ldots, p_{m}\}$ by polar angle from $p_0$ . For indices $j < k$ in this angular order, define $\mathrm{dp}[j][k]$ as the maximum double-area of a convex polygon with $p_0$ as anchor, first vertex $p_f$ (some $f \le j$ ), and last two vertices $p_j, p_k$ in counterclockwise order. Then:

\mathrm{dp}[k][\ell] = \max_{j: \, (j,k,\ell) \text{ valid}} \bigl(\mathrm{dp}[j][k] + \sigma(p_0, p_k, p_\ell)\bigr),

where “valid” requires: (i) the turn at $p_k$ is left (convexity), (ii) the fan triangle $\triangle p_0 p_k p_\ell$ is empty, (iii) the angular span from first vertex to $p_\ell$ is less than $\pi$ (closure feasibility).

Proof. The polygon is built incrementally from the anchor $p_0$ , appending vertices in counterclockwise angular order. Each extension from edge $(p_j, p_k)$ to vertex $p_\ell$ adds the triangle $\triangle p_0 p_k p_\ell$ to the polygon area. By Theorem 1, checking emptiness of each such fan triangle suffices for the emptiness of the entire polygon. The left-turn condition at $p_k$ (i.e., $\sigma(p_j, p_k, p_\ell) > 0$ ) maintains convexity. Tracking the first vertex $p_f$ and requiring $\sigma(p_\ell, p_0, p_f) > 0$ ensures the polygon can be closed convexly. $\blacksquare$

Corollary. The maximum-area convex hole in $S$ is obtained by maximizing over all anchors $p_0$ and all final edges $(k, \ell)$ for which closure is valid:

\max_{p_0 \in S} \max_{(k,\ell)} \frac{\mathrm{dp}[k][\ell]}{2}.

Editorial

A convex hole is a convex polygon whose vertices are from the point set and whose interior contains no other points from the set. Algorithm: For each anchor point p0, sort remaining points by polar angle, then DP on (last two vertices) tracking first vertex for closure checks. Empty-triangle tests use precomputed bitsets. We iterate over each anchor p0 in S.

Pseudocode

    n = |S|
    Precompute bitset B[i][j] for each directed edge (i, j):
        B[i][j] = {p in S : sigma(i, j, p) > 0} // points strictly left of i->j

        Return B[a][b] AND B[b][c] AND B[c][a] == 0 // bitwise intersection

    best = 0
    for each anchor p0 in S:
        Sort remaining points by polar angle from p0
        Initialize dp[j][k] from empty base triangles (p0, j, k)
        Extend DP: for each (j, k) with dp[j][k] > 0, try all l > k
            Check left-turn, fan-triangle emptiness, angular span
            Update dp[k][l] and check closure for new maximum
    Return best / 2

Complexity Analysis

Time: Precomputing all bitsets $B[i][j]$ takes $O(n^3)$ . The DP for each anchor runs in $O(n^3)$ worst case (over all triples $(j, k, \ell)$ ). Total: $O(n^4)$ . With $n = 500$ , this is approximately $1.56 \times 10^{10}$ , which requires significant optimization (bitset operations, pruning) but is feasible.

Space: $O(n^3 / 64)$ for bitset storage of $B[i][j]$ (using 64-bit words); $O(n^2)$ for the DP table per anchor.

Answer

\boxed{104924.0}

C++ project_euler/problem_252/solution.cpp

/*
 * Project Euler Problem 252: Convex Holes
 *
 * Find the largest area convex hole among 500 LCG-generated points.
 * A convex hole is a convex polygon whose vertices are from the point set
 * and whose interior contains no other points.
 *
 * Algorithm: For each anchor p0, sort remaining points by polar angle,
 * then DP on (last two vertices) tracking the first vertex for closure.
 * Empty-triangle tests use precomputed bitsets (64-bit words).
 *
 * Answer: 104924.0
 */

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef pair<ll, ll> pll;

ll cross(pll O, pll A, pll B) {
    return (A.first - O.first) * (B.second - O.second)
         - (A.second - O.second) * (B.first - O.first);
}

int main() {
    // Generate points
    vector<pll> pts;
    ll s = 290797;
    vector<ll> t;
    for (int i = 0; i < 1000; i++) {
        s = (s * s) % 50515093LL;
        t.push_back(s % 2000 - 1000);
    }
    for (int k = 0; k < 500; k++)
        pts.push_back({t[2 * k], t[2 * k + 1]});

    int n = pts.size();
    const int WORDS = (n + 63) / 64;

    // Precompute bitsets: B[i][j] = points strictly left of directed edge i->j
    vector<vector<vector<uint64_t>>> B(n, vector<vector<uint64_t>>(n, vector<uint64_t>(WORDS, 0)));
    for (int i = 0; i < n; i++)
        for (int j = 0; j < n; j++) {
            if (i == j) continue;
            for (int p = 0; p < n; p++) {
                if (p == i || p == j) continue;
                if (cross(pts[i], pts[j], pts[p]) > 0)
                    B[i][j][p / 64] |= (1ULL << (p % 64));
            }
        }

    auto triangle_empty = [&](int a, int b, int c) -> bool {
        for (int w = 0; w < WORDS; w++)
            if (B[a][b][w] & B[b][c][w] & B[c][a][w]) return false;
        return true;
    };

    ll best = 0;

    for (int p0 = 0; p0 < n; p0++) {
        vector<int> idx;
        for (int i = 0; i < n; i++)
            if (i != p0) idx.push_back(i);

        sort(idx.begin(), idx.end(), [&](int a, int b) {
            ll ax = pts[a].first - pts[p0].first, ay = pts[a].second - pts[p0].second;
            ll bx = pts[b].first - pts[p0].first, by = pts[b].second - pts[p0].second;
            int qa = (ay > 0 || (ay == 0 && ax > 0)) ? 0 : 1;
            int qb = (by > 0 || (by == 0 && bx > 0)) ? 0 : 1;
            if (qa != qb) return qa < qb;
            ll cr = ax * by - ay * bx;
            if (cr != 0) return cr > 0;
            return ax * ax + ay * ay < bx * bx + by * by;
        });

        int m = idx.size();
        vector<int> first_vtx(m * m, -1);
        auto FV = [&](int j, int k) -> int& { return first_vtx[j * m + k]; };
        vector<ll> dp(m * m, -1);
        auto DP = [&](int j, int k) -> ll& { return dp[j * m + k]; };

        // Base triangles
        for (int j = 0; j < m; j++)
            for (int k2 = j + 1; k2 < m; k2++) {
                ll cr = cross(pts[p0], pts[idx[j]], pts[idx[k2]]);
                if (cr <= 0) continue;
                if (!triangle_empty(p0, idx[j], idx[k2])) continue;
                DP(j, k2) = cr;
                FV(j, k2) = j;
                best = max(best, cr);
            }

        // Extend
        for (int k2 = 0; k2 < m; k2++)
            for (int j = 0; j < k2; j++) {
                if (DP(j, k2) < 0) continue;
                int f = FV(j, k2);
                for (int l = k2 + 1; l < m; l++) {
                    ll cr1 = cross(pts[idx[j]], pts[idx[k2]], pts[idx[l]]);
                    if (cr1 <= 0) continue;
                    ll cr2 = cross(pts[p0], pts[idx[k2]], pts[idx[l]]);
                    if (cr2 <= 0) continue;
                    if (!triangle_empty(p0, idx[k2], idx[l])) continue;
                    if (cross(pts[p0], pts[idx[f]], pts[idx[l]]) <= 0) continue;

                    ll new_area = DP(j, k2) + cr2;
                    if (new_area > DP(k2, l)) {
                        DP(k2, l) = new_area;
                        FV(k2, l) = f;
                    }
                    if (cross(pts[idx[k2]], pts[idx[l]], pts[p0]) > 0 &&
                        cross(pts[idx[l]], pts[p0], pts[idx[f]]) > 0)
                        best = max(best, new_area);
                }
            }
    }

    printf("%.1f\n", best / 2.0);
    return 0;
}

Python project_euler/problem_252/solution.py

"""
Project Euler Problem 252: Convex Holes

Find the largest area convex hole among 500 LCG-generated points.
A convex hole is a convex polygon whose vertices are from the point set
and whose interior contains no other points from the set.

Algorithm: For each anchor point p0, sort remaining points by polar angle,
then DP on (last two vertices) tracking first vertex for closure checks.
Empty-triangle tests use precomputed bitsets.

Answer: 104924.0
"""

import math
from functools import cmp_to_key


def cross(O, A, B):
    return (A[0] - O[0]) * (B[1] - O[1]) - (A[1] - O[1]) * (B[0] - O[0])


def solve():
    # Generate points
    s = 290797
    t = []
    for _ in range(1000):
        s = (s * s) % 50515093
        t.append(s % 2000 - 1000)
    pts = [(t[2 * k], t[2 * k + 1]) for k in range(500)]
    n = len(pts)

    WORDS = (n + 63) // 64

    # Precompute bitsets: B[i][j] = set of points strictly left of edge i->j
    B = [[[0] * WORDS for _ in range(n)] for _ in range(n)]
    for i in range(n):
        for j in range(n):
            if i == j:
                continue
            for p in range(n):
                if p == i or p == j:
                    continue
                if cross(pts[i], pts[j], pts[p]) > 0:
                    B[i][j][p // 64] |= 1 << (p % 64)

    def triangle_empty(a, b, c):
        for w in range(WORDS):
            if B[a][b][w] & B[b][c][w] & B[c][a][w]:
                return False
        return True

    best = 0

    for p0 in range(n):
        idx = [i for i in range(n) if i != p0]

        def cmp(a, b):
            ax = pts[a][0] - pts[p0][0]
            ay = pts[a][1] - pts[p0][1]
            bx = pts[b][0] - pts[p0][0]
            by_ = pts[b][1] - pts[p0][1]
            qa = 0 if (ay > 0 or (ay == 0 and ax > 0)) else 1
            qb = 0 if (by_ > 0 or (by_ == 0 and bx > 0)) else 1
            if qa != qb:
                return -1 if qa < qb else 1
            cr = ax * by_ - ay * bx
            if cr != 0:
                return 1 if cr > 0 else -1
            da = ax * ax + ay * ay
            db = bx * bx + by_ * by_
            if da != db:
                return -1 if da < db else 1
            return 0

        idx.sort(key=cmp_to_key(cmp))
        m = len(idx)
        dp = {}
        fv = {}

        # Base: triangles (p0, j, k)
        for j in range(m):
            for k in range(j + 1, m):
                cr = cross(pts[p0], pts[idx[j]], pts[idx[k]])
                if cr <= 0:
                    continue
                if not triangle_empty(p0, idx[j], idx[k]):
                    continue
                dp[(j, k)] = cr
                fv[(j, k)] = j
                best = max(best, cr)

        # Extend
        for k in range(m):
            for j in range(k):
                if (j, k) not in dp:
                    continue
                f = fv[(j, k)]
                area_jk = dp[(j, k)]
                for l in range(k + 1, m):
                    cr1 = cross(pts[idx[j]], pts[idx[k]], pts[idx[l]])
                    if cr1 <= 0:
                        continue
                    cr2 = cross(pts[p0], pts[idx[k]], pts[idx[l]])
                    if cr2 <= 0:
                        continue
                    if not triangle_empty(p0, idx[k], idx[l]):
                        continue
                    if cross(pts[p0], pts[idx[f]], pts[idx[l]]) <= 0:
                        continue
                    new_area = area_jk + cr2
                    if (k, l) not in dp or new_area > dp[(k, l)]:
                        dp[(k, l)] = new_area
                        fv[(k, l)] = f
                    if (cross(pts[idx[k]], pts[idx[l]], pts[p0]) > 0
                            and cross(pts[idx[l]], pts[p0], pts[idx[f]]) > 0):
                        best = max(best, new_area)

    print(f"{best / 2.0:.1f}")


if __name__ == "__main__":
    solve()