PROBLEM LINK:
Practice
Contest: Division 1
Contest: Division 2
Contest: Division 3
Contest: Division 4
Author: ziad_el_gafy
Tester & Editorialist: iceknight1093
DIFFICULTY:
2990
PREREQUISITES:
Elementary number theory, Segmented sieve of Eratosthenes
PROBLEM:
Given a simple polygon with N vertices and perimeter at most 10^7, count the number of lattice points (x, y) on its perimeter such that (x+y) is prime.
The coordinates of the points can be as large as 10^{14}.
EXPLANATION:
The perimeter of the polygon is at most 10^7, so there are certainly at most 10^7 lattice points lying on it.
Let’s first try to find all such lattice points.
The polygon is simple and we’re given its vertices in order, so the perimeter can be broken up into a bunch of line segments.
It suffices to find lattice points lying on each such segment.
So, consider the segment from (x_1, y_1) to (x_2, y_2).
Then, if d = \gcd(x_2 - x_1, y_2 - y_1), there are exactly g+1 lattice points lying on this segment (including the endpoints); in fact their coordinates will be of the form
for 0 \leq k \leq d.
A quick explanation of why this is the case can be found here.
This allows us to easily find all the lattice points lying on the perimeter in \mathcal{O}(N + P) (where P is the perimeter).
Now, we need to test each one for primality.
A brute-force primality check for an integer \alpha would take \mathcal{O}(\sqrt \alpha) time, which is obviously too slow here since we need to run it several times.
If the maximum value of \alpha were small enough (say, \leq 10^7), we could use a sieve of Eratosthenes to precompute all primes upto the limit and then answer queries in \mathcal{O}(1).
However, here \alpha can be as large as 2\cdot 10^{14}, so this has no chance of working.
The important observation here comes once again from the fact that the perimeter is small: although the values themselves might be large, the range of values that we need to check primality for is quite small!
In particular, let m_x and M_x be the smallest and largest x-coordinates of some lattice points, respectively. SImilarly define m_y and M_y.
Then, the smallest integer we care about is m_x + m_y, and the largest is M_x + M_y.
However, since the perimeter is \leq 10^7, we have M_x - m_x \leq 10^7 and M_y - m_y \leq 10^7, and so (M_x + M_y) - (m_x + m_y) \leq 2\cdot 10^7.
This range being small allows us to use a modification of the sieve of Eratosthenes known as the segmented sieve, to compute primality of everything in this range.
The idea to precompute primes in range [L, R] is as follows:
- First, use a normal sieve to compute all the primes upto \sqrt{R}.
- Then, for each prime p, mark all of its multiples in range [L, R] as not prime.
- Finally, everything in range that wasn’t marked is a prime.
This works because every non-prime in this range will have a prime factor that’s \leq \sqrt{R}.
The time complexity of this is \mathcal{O}((R-L+1 + \sqrt{R})\log\log R).
Notice that this immediately solves the problem: the segmented sieve precomputes all primes we care about, then primality of each lattice point can be checked in \mathcal{O}(1) time.
TIME COMPLEXITY
\mathcal{O}(N + P + (R-L+1 + \sqrt{R})\log\log R) per test case, where P is the perimeter of the polygon and L and R are its minimum/maximum possible coordinate sums, bounded by 2\cdot 10^{14}.
CODE:
Author's code (C++)
#pragma GCC optimize("Ofast")
#pragma GCC optimize("O3")
#include <bits/stdc++.h>
using namespace std;
#define el '\n'
#define F first
#define S second
typedef long long ll;
typedef long double ld;
bool multipleTestCases = 0, sublime = 1;
const int N = 15000000 + 5;
bool isP[N], isPrime[N];
vector<int> primes;
ll n, mnX = 1e15, mnY = 1e15, mxX, mxY, ans, lowerLimit, upperLimit;
vector<pair<ll, ll>> v;
void sieve() {
isPrime[0] = isPrime[1] = 1;
for (int i = 2; i < N; i++) {
isP[i] = isPrime[i] = 1;;
}
for (int i = 2; i * i < N; i++) {
if (isP[i]) {
for (int j = i * i; j < N; j += i) {
isP[j] = 0;
}
}
}
for (int i = 2; i < N; i++) {
if (isP[i]) {
primes.push_back(i);
}
}
}
void processPrimes() {
lowerLimit = mnX + mnY, upperLimit = mxX + mxY;
if (upperLimit >= N) {
for (auto &p : primes) {
ll md = lowerLimit % p;
ll start = (md ? lowerLimit + (p - md) : lowerLimit);
if (start == p) {
start += p;
}
for (ll i = start; i <= upperLimit; i += p) {
isPrime[i - lowerLimit] = 0;
}
}
}
else {
for (int i = lowerLimit; i <= upperLimit; i++) {
isPrime[i - lowerLimit] = isP[i];
}
}
}
void doWork() {
sieve();
cin >> n;
v.resize(n);
for (auto &i : v) {
cin >> i.F;
}
for (auto &i : v) {
cin >> i.S;
mnX = min(mnX, i.F);
mxX = max(mxX, i.F);
mnY = min(mnY, i.S);
mxY = max(mxY, i.S);
}
processPrimes();
for (int i = 0; i < n; i++) {
pair<ll, ll> cur = v[i], nxt = v[(i + 1) % n];
ll dx = nxt.F - cur.F, dy = nxt.S - cur.S, g = __gcd(abs(dx), abs(dy));
ll xUnit = dx / g, yUnit = dy / g;
ll x = cur.F, y = cur.S;
while (x != nxt.F or y != nxt.S) {
ans += isPrime[x + y - lowerLimit];
x += xUnit;
y += yUnit;
}
}
cout << ans << el;
}
signed main() {
#ifdef ONLINE_JUDGE
ios_base::sync_with_stdio(0), cin.tie(0), cout.tie(0);
#else
if (sublime) {
freopen("input.txt", "r", stdin);
freopen("output.txt", "w", stdout);
}
#endif
int tests = 1;
if (multipleTestCases) {
cin >> tests;
}
for (int tc = 1; tc <= tests; tc++) {
doWork();
}
}
Editorialist's code (Python)
mx = int(1.5e7)
isprime = [1]*mx
isprime_seg = [1]*mx
n = int(input())
x = list(map(int, input().split()))
y = list(map(int, input().split()))
L, R = min(x) + min(y), max(x) + max(y)
for i in range(2, mx):
if isprime[i] == 0: continue
for j in range(2*i, mx, i): isprime[j] = 0
lo = i * (max(2, (L+i-1)//i))
while lo <= R:
isprime_seg[lo - L] = 0
lo += i
import math
ans = 0
for i in range(n):
x1, y1 = x[i], y[i]
x2, y2 = x[i-1], y[i-1]
g = math.gcd(x2 - x1, y2 - y1)
for k in range(g):
X, Y = x1 + (x2 - x1)//g * k, y1 + (y2 - y1)//g * k
ans += isprime_seg[X+Y-L]
print(ans)