Problem 165: Intersections

Mathematical Development

Theorem 1 (Parametric intersection). Given two line segments $P_1P_2$ and $P_3P_4$ where $P_i = (x_i, y_i)$, define:

If $D \neq 0$, the lines containing the segments intersect at the point $P_1 + t(P_2 - P_1)$. This point is an interior point of both segments if and only if $0 < t < 1$ and $0 < u < 1$ (strictly).

Proof. Points on segment $P_1P_2$ are parameterized as $P_1 + t(P_2 - P_1)$ for $t \in [0,1]$, and points on segment $P_3P_4$ as $P_3 + u(P_4 - P_3)$ for $u \in [0,1]$. Setting these equal and solving the $2 \times 2$ linear system:

By Cramer’s rule, the determinant is $-D$ (with appropriate sign convention), giving the formulas above. The intersection is interior to both segments when $t$ and $u$ are strictly between 0 and 1, excluding endpoints. $\square$

Lemma 1 (Exact arithmetic sufficiency). Since all coordinates $t_n \in \{0, 1, \ldots, 499\}$, the numerators and denominator $D$ in the intersection formulas are integers bounded by $|D| \leq 2 \cdot 499^2 = 498{,}002$. The intersection point $(x, y)$ can be represented exactly as $\left(\frac{p}{q}, \frac{r}{q}\right)$ with $\gcd(|p|, |q|, |r|) = 1$, using 64-bit integers.

Proof. The coordinates are integers in $[0, 499]$. Each term in $D$ is a product of two differences of coordinates, each bounded by 499 in absolute value. So $|D| \leq 2 \cdot 499^2$. The numerator of the intersection point’s $x$-coordinate is:

where $t_{\text{num}}$ is the numerator of $t$. All quantities fit in 64-bit integers (products of at most three 499-bounded values: $\leq 499^3 < 2^{27}$). $\square$

Lemma 2 (Distinct point identification). Two intersection points $\left(\frac{p_1}{q_1}, \frac{r_1}{q_1}\right)$ and $\left(\frac{p_2}{q_2}, \frac{r_2}{q_2}\right)$ (in lowest terms with positive denominators) are identical if and only if $p_1 = p_2$, $q_1 = q_2$, and $r_1 = r_2$.

Proof. Rational numbers in lowest terms have a unique representation. By normalizing each coordinate pair $(p/q, r/q)$ to have $q > 0$ and $\gcd(|p|, |r|, q)$ appropriately reduced, equality of the triple $(p, q, r)$ is equivalent to geometric equality. $\square$

Editorial

The algorithm has two phases. First it generates the 5000 segments exactly as prescribed by the Blum Blum Shub recurrence. After that, the task is simply to examine every unordered pair of segments and determine whether their supporting lines meet at an interior point of both segments.

The important implementation choice is to keep every intersection in exact rational form. For a pair of non-parallel segments, the parameters $t$ and $u$ decide whether the crossing is strictly internal, and the resulting point is stored as a reduced triple representing $(x,y)$ over a common denominator. Inserting those normalized triples into a set automatically deduplicates repeated geometric intersections coming from different segment pairs.

Pseudocode

Generate the 5000 segments from the pseudorandom sequence.
Create an empty set of normalized intersection points.

For each unordered pair of segments:
    Compute the line determinant.
    If it is zero, the segments are parallel or collinear, so skip the pair.

    Compute the two intersection parameters for the segments.
    Keep the pair only when both parameters lie strictly between 0 and 1,
    which means the intersection is interior to both segments.

    Build the intersection point with integer numerators and denominator.
    Normalize the triple by making the denominator positive and dividing by the gcd.
    Insert the normalized point into the set.

Return the size of the set.

Complexity Analysis

• Time: $O(N^2)$ where $N = 5000$. The number of pairs is $\binom{5000}{2} = 12{,}497{,}500$. Each pair requires $O(1)$ arithmetic and a hash-set insertion. Total: $O(N^2) = O(2.5 \times 10^7)$.
• Space: $O(K)$ where $K$ is the number of distinct interior intersection points, for the hash set. In the worst case $K = O(N^2)$.

Answer

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

// Problem 165: Intersections
// Find distinct interior intersection points of 5000 line segments
// generated by the Blum Blum Shub generator.

typedef long long ll;
typedef __int128 lll;

struct Seg {
    ll x1, y1, x2, y2;
};

int main() {
    // Generate segments
    const int N = 5000;
    vector<Seg> segs(N);

    ll s = 290797;
    auto next_s = [&]() -> ll {
        s = (ll)((lll)s * s % 50515093LL);
        return s;
    };

    for (int i = 0; i < N; i++) {
        ll t1 = next_s() % 500;
        ll t2 = next_s() % 500;
        ll t3 = next_s() % 500;
        ll t4 = next_s() % 500;
        segs[i] = {t1, t2, t3, t4};
    }

    // For each pair, find interior intersection
    // Use exact rational arithmetic
    // Store points as (px/d, py/d) reduced to lowest terms with d > 0

    set<tuple<ll,ll,ll,ll>> points;

    auto gcd_func = [](ll a, ll b) -> ll {
        a = abs(a); b = abs(b);
        while (b) { a %= b; swap(a, b); }
        return a;
    };

    for (int i = 0; i < N; i++) {
        ll x1 = segs[i].x1, y1 = segs[i].y1;
        ll x2 = segs[i].x2, y2 = segs[i].y2;
        ll dx1 = x2 - x1, dy1 = y2 - y1;

        for (int j = i + 1; j < N; j++) {
            ll x3 = segs[j].x1, y3 = segs[j].y1;
            ll x4 = segs[j].x2, y4 = segs[j].y2;
            ll dx2 = x4 - x3, dy2 = y4 - y3;

            // D = dx1 * dy2 - dy1 * dx2
            ll D = dx1 * dy2 - dy1 * dx2;
            if (D == 0) continue; // parallel or collinear

            // t_num = (x3 - x1) * dy2 - (y3 - y1) * dx2
            // u_num = (x3 - x1) * dy1 - (y3 - y1) * dx1
            ll t_num = (x3 - x1) * dy2 - (y3 - y1) * dx2;
            ll u_num = (x3 - x1) * dy1 - (y3 - y1) * dx1;

            // Interior: 0 < t < 1 and 0 < u < 1
            // t = t_num / D, u = u_num / D
            // If D > 0: 0 < t_num < D and 0 < u_num < D
            // If D < 0: D < t_num < 0 and D < u_num < 0

            bool t_ok, u_ok;
            if (D > 0) {
                t_ok = (t_num > 0 && t_num < D);
                u_ok = (u_num > 0 && u_num < D);
            } else {
                t_ok = (t_num < 0 && t_num > D);
                u_ok = (u_num < 0 && u_num > D);
            }

            if (!t_ok || !u_ok) continue;

            // Compute intersection point: (x1 + t * dx1, y1 + t * dy1)
            // = (x1 + t_num * dx1 / D, y1 + t_num * dy1 / D)
            // = ((x1 * D + t_num * dx1) / D, (y1 * D + t_num * dy1) / D)

            lll px = (lll)x1 * D + (lll)t_num * dx1;
            lll py = (lll)y1 * D + (lll)t_num * dy1;
            lll dd = D;

            // Normalize: dd > 0
            if (dd < 0) { px = -px; py = -py; dd = -dd; }

            // Reduce to lowest terms
            ll px_l = (ll)px, py_l = (ll)py, dd_l = (ll)dd;
            ll g = gcd_func(gcd_func(abs(px_l), abs(py_l)), dd_l);
            px_l /= g; py_l /= g; dd_l /= g;

            points.insert({px_l, py_l, dd_l, 0});
        }
    }

    cout << points.size() << endl;
    return 0;
}

"""
Problem 165: Intersections

Count distinct interior intersection points of 5000 line segments
generated by the Blum Blum Shub pseudo-random number generator.
Uses exact rational arithmetic to identify distinct points.
"""

from math import gcd

def solve():
    N = 5000

    # Generate segments using BBS
    s = 290797
    coords = []
    for _ in range(4 * N):
        s = s * s % 50515093
        coords.append(s % 500)

    segments = []
    for i in range(N):
        x1, y1, x2, y2 = coords[4*i], coords[4*i+1], coords[4*i+2], coords[4*i+3]
        segments.append((x1, y1, x2, y2))

    points = set()

    for i in range(N):
        x1, y1, x2, y2 = segments[i]
        dx1, dy1 = x2 - x1, y2 - y1
        for j in range(i + 1, N):
            x3, y3, x4, y4 = segments[j]
            dx2, dy2 = x4 - x3, y4 - y3

            D = dx1 * dy2 - dy1 * dx2
            if D == 0:
                continue

            t_num = (x3 - x1) * dy2 - (y3 - y1) * dx2
            u_num = (x3 - x1) * dy1 - (y3 - y1) * dx1

            # Interior: 0 < t < 1 and 0 < u < 1
            if D > 0:
                if not (0 < t_num < D and 0 < u_num < D):
                    continue
            else:
                if not (D < t_num < 0 and D < u_num < 0):
                    continue

            # Intersection: (x1 + t_num/D * dx1, y1 + t_num/D * dy1)
            # = (x1*D + t_num*dx1, y1*D + t_num*dy1) / D
            px = x1 * D + t_num * dx1
            py = y1 * D + t_num * dy1
            dd = D

            if dd < 0:
                px, py, dd = -px, -py, -dd

            g = gcd(gcd(abs(px), abs(py)), dd)
            px //= g
            py //= g
            dd //= g

            points.add((px, py, dd))

    print(len(points))

solve()