Why this matters when you build with AI
An AI tool can write a number theory solution that reads well and passes the sample. The usual failures are specific: a 32-bit sum that wraps, a modulus that is not prime when the code divides by it, a sieve sized for the sample and not for the limit, or a loop over every pair when the statement allows 2000 values. None of these show up on a short example. A problem set at the right difficulty is a way to practise catching them before a contest does.
This set uses only the ideas from the lessons in this module. Each problem names the lesson it draws on in its hint. Work on the problems first, then compare with the reference approach. A solution that you cannot test is not finished. For every problem, write a brute force for small inputs and check that your fast solution agrees with it on random small cases.
The problems are original. We wrote every statement, sample and limit ourselves. We did not take problems from any judge or contest. The reference programs are in C++17.
How to work the set under timed conditions
Give yourself 150 minutes for all fourteen problems. Work through them in order, and move on after 15 minutes on one problem, then come back to it. For each problem do these steps in order.
- Read the statement and write the limits on paper. Note the largest value that can appear.
- Write a brute force for small inputs first. It must be simple enough to trust.
- Write the fast solution. Keep 64-bit integers for sums and products.
- Run both on the sample, then on at least 100 random small inputs.
- Time your fast solution at the largest limit the statement allows.
- Compare the measured time with the complexity target in the statement.
The complexity target is the time a good solution should meet at the largest limit. A solution that is correct but too slow at the largest limit is not finished.
Here is a worked trace for problem 1 on its sample, before the problems start. The array is 12, 18, 35, 49. The gcd of all four is found step by step: gcd(12, 18) = 6, gcd(6, 35) = 1, gcd(1, 49) = 1. The pairs with gcd 1 are (12, 35), (12, 49), (18, 35) and (18, 49), so the count is 4. The pairs (12, 18) with gcd 6 and (35, 49) with gcd 7 are not counted.
The problems
Problem 1: gcd and coprime pairs
Statement. You get N numbers. Print their gcd. Then print how many pairs of positions i < j have gcd(a_i, a_j) = 1.
Input. The first line holds N. The second line holds the N numbers.
Output. Two lines. The first is the gcd of all numbers. The second is the number of coprime pairs.
Limits. 1 ≤ N ≤ 2000 and 1 ≤ a_i ≤ 10^9.
Example.
input
4
12 18 35 49
output
1
4Hint. Start the gcd fold at the first number, or at 0. Lesson: gcd and lcm.
Target. O(N² log A) with A the largest value. For N = 2000 that is about 2 million gcd calls.
Show the answer
Fold gcd over the array, starting from 0, since gcd(0, x) = x. Then test every pair with a gcd call. The fold must not start at 1. The value gcd(1, x) is always 1, so a fold that starts at 1 prints 1 for every input.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref1(istream& in, ostream& out) {
int n;
in >> n;
vector<ll> a(n);
for (ll& x : a) in >> x;
ll g = 0;
for (ll x : a) g = std::gcd(g, x);
ll cnt = 0;
for (int i = 0; i < n; i++)
for (int j = i + 1; j < n; j++)
if (std::gcd(a[i], a[j]) == 1) cnt++;
out << g << "\n" << cnt << "\n";
}
int main() {
ref1(cin, cout);
return 0;
}Checked: the program prints 1 and 4 on the sample. A subtraction-based gcd brute force agreed on 300 random inputs.
Problem 2: solutions of a linear equation
Statement. For each query you get a, b and c. Count the pairs of non-negative integers (x, y) with a times x plus b times y equal to c.
Input. The first line holds Q. Each of the next Q lines holds a, b and c.
Output. One count per line. The count can be as large as 10^18 + 1, so it needs a 64-bit integer.
Limits. 1 ≤ Q ≤ 10^4, 1 ≤ a, b ≤ 10^9, 0 ≤ c ≤ 10^18.
Example.
input
3
3 5 15
2 3 12
4 6 7
output
2
3
0Hint. Use extended Euclid. The equation has solutions only when gcd(a, b) divides c. Lesson: gcd, lcm and extended Euclid.
Target. O(log C) per query.
Show the answer
Let g = gcd(a, b). If g does not divide c, the answer is 0. Otherwise divide a, b and c by g. Find the smallest non-negative x0 with a1 times x0 ≡ c1 (mod b1), using the inverse of a1 modulo b1 from extended Euclid. The other solutions of the congruence are x0 + k times b1. Each needs y ≥ 0, which means x ≤ c1 / a1. Count the values from x0 to that bound in steps of b1.
The product c1 times the inverse can reach 10^18 times 10^9, which overflows 64 bits. The reference uses __int128 for that product. The extended Euclid coefficients can be negative, so reduce them modulo b1 before use.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ll ext_gcd(ll a, ll b, ll& x, ll& y) { // returns g = gcd(a, b) and a*x + b*y = g
if (b == 0) {
x = 1;
y = 0;
return a;
}
ll x1, y1;
ll g = ext_gcd(b, a % b, x1, y1);
x = y1;
y = x1 - (a / b) * y1;
return g;
}
void ref2(istream& in, ostream& out) {
int q;
in >> q;
while (q--) {
ll a, b, c;
in >> a >> b >> c;
ll x, y;
ll g = ext_gcd(a, b, x, y);
if (c % g != 0) {
out << 0 << "\n";
continue;
}
ll a1 = a / g, b1 = b / g, c1 = c / g;
ll inv = 0;
if (b1 > 1) {
ext_gcd(a1, b1, x, y); // a1 * x + b1 * y = 1
inv = ((x % b1) + b1) % b1;
}
ll x0 = (ll)((__int128)(c1 % b1) * inv % b1); // smallest x >= 0 with a1 * x = c1 (mod b1)
ll xmax = c1 / a1; // largest x with a1 * x <= c1
if (x0 > xmax) out << 0 << "\n";
else out << (xmax - x0) / b1 + 1 << "\n";
}
}
int main() {
ref2(cin, cout);
return 0;
}Checked: the sample output is 2, 3 and 0. A brute force over x and y agreed on 300 random small cases.
Problem 3: primes in a range
Statement. For each query, print how many primes p satisfy L ≤ p ≤ R.
Input. The first line holds Q. Each of the next Q lines holds L and R.
Output. One count per line.
Limits. 1 ≤ Q ≤ 10^5, 1 ≤ L ≤ R ≤ 10^7.
Example.
input
2
10 20
1 100
output
4
25Hint. Sieve once up to the largest R, then store prefix counts. Lesson: sieves.
Target. O(R log log R) for the sieve, then O(1) per query.
Show the answer
Sieve up to the largest R in the input. Store pi[x], the number of primes up to x, as a prefix array. Each query is pi[R] - pi[L - 1]. The entry pi[0] must be 0, which covers queries with L = 1.
The prefix array is an int for each value up to 10^7, which is about 40 MB. A sieve of char flags adds about 10 MB. If the memory limit is tight, store the prefix counts in a smaller type or sieve only the odd numbers.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref3(istream& in, ostream& out) {
int q;
in >> q;
vector<pair<int, int>> qs(q);
int mx = 2;
for (auto& [l, rr] : qs) {
in >> l >> rr;
mx = max(mx, rr);
}
vector<char> comp(mx + 1, 0);
comp[0] = comp[1] = 1;
for (ll i = 2; i * i <= mx; i++)
if (!comp[i])
for (ll j = i * i; j <= mx; j += i) comp[j] = 1;
vector<int> pi(mx + 1, 0);
for (int x = 1; x <= mx; x++) pi[x] = pi[x - 1] + (comp[x] ? 0 : 1);
for (auto [l, rr] : qs) out << pi[rr] - pi[l - 1] << "\n";
}
int main() {
ref3(cin, cout);
return 0;
}Checked: the sample prints 4 and 25. A trial division brute force agreed on 300 random small inputs.
Problem 4: distinct prime factors of each value
Statement. You get N numbers. For each one, print how many different primes divide it. The value 1 has no prime factors, so its count is 0.
Input. The first line holds N. The second line holds the N numbers.
Output. The N counts on one line, separated by spaces.
Limits. 1 ≤ N ≤ 10^6 and 1 ≤ x ≤ 10^6.
Example.
input
5
1 12 30 97 1000
output
0 2 3 1 2Hint. Build the smallest prime factor table once. Then divide out each prime completely before you count the next one. Lesson: factorisation and divisors.
Target. O(M log log M) for the table with M the largest value, then O(log x) per value.
Show the answer
Build spf, the smallest prime factor of each value up to the largest input. For each x, repeatedly take p = spf[x], count one prime, and divide out every copy of p. The count is of distinct primes, so the inner loop removes all copies before the count moves on. For x = 1 the loop never runs, which gives 0.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref4(istream& in, ostream& out) {
int n;
in >> n;
vector<int> a(n);
int mx = 2;
for (int& x : a) {
in >> x;
mx = max(mx, x);
}
vector<int> spf(mx + 1, 0);
for (int i = 2; i <= mx; i++)
if (spf[i] == 0)
for (int j = i; j <= mx; j += i)
if (spf[j] == 0) spf[j] = i;
for (int i = 0; i < n; i++) {
int x = a[i], cnt = 0;
while (x > 1) {
int p = spf[x];
cnt++;
while (x % p == 0) x /= p;
}
out << cnt << (i + 1 < n ? ' ' : '\n');
}
}
int main() {
ref4(cin, cout);
return 0;
}Checked: the sample prints 0 2 3 1 2. A trial division brute force agreed on 300 random small inputs.
Problem 5: the value with the most divisors
Statement. Find the value x in the range from L to R that has the most divisors. If several values tie, print the smallest. Also print its number of divisors.
Input. One line with L and R.
Output. One line with x and its divisor count.
Limits. 1 ≤ L ≤ R ≤ 10^6.
Example.
input
10 20
output
12 6Hint. Count divisors for all values up to R with a loop over multiples. Then scan from L with a strict comparison. Lesson: factorisation and divisors, and harmonic-sum tricks.
Target. O(R log R).
Show the answer
Count the divisors of every value up to R with a loop over multiples, as in the harmonic-sum lesson. Then scan from L to R. Replace the best value only when the count is strictly larger. A >= comparison would keep the largest value among ties, which is the wrong answer here.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref5(istream& in, ostream& out) {
int L, R;
in >> L >> R;
vector<int> d(R + 1, 0);
for (int i = 1; i <= R; i++)
for (int j = i; j <= R; j += i) d[j]++;
int best = L;
for (int x = L + 1; x <= R; x++)
if (d[x] > d[best]) best = x; // strict, so the smallest x wins a tie
out << best << " " << d[best] << "\n";
}
int main() {
ref5(cin, cout);
return 0;
}Checked: the sample prints 12 and 6. On the range 1 to 10^6 the program printed 720720 with 240 divisors. A brute force that counts divisors by trying every number up to x agreed on 300 random small ranges.
Problem 6: a power tower modulo m
Statement. For each query you get a, b, c and m. Print a raised to the power b raised to the power c, taken modulo m. The exponent b^c is far too large to compute, so use the reduction rules for exponents.
Input. The first line holds T. Each of the next T lines holds a, b, c and m.
Output. One value per line.
Limits. 1 ≤ T ≤ 10^4 and 1 ≤ a, b, c, m ≤ 10^9. The result is the value of a to the power b^c, reduced modulo m.
Example.
input
3
2 3 2 10
3 2 3 1000
5 3 4 97
output
2
561
19Hint. Let phi be Euler's totient of m. Reduce the exponent modulo phi, but add phi back when the true exponent is at least phi. Lesson: Euler's totient, and fermat and Euler in practice.
Target. O(T times sqrt(m)) for the totient by trial division, plus O(log) for each power.
Show the answer
Compute phi(m) by trial division. Let E = b^c. The rule is that a^E ≡ a^(E mod phi + phi) (mod m) whenever E ≥ phi. This holds because phi(m) is at least the largest exponent of any prime in m. If E < phi, use E itself. The check "is b^c at least phi" stops after about 30 multiplications, because b ≥ 2 in that branch, and b = 1 is handled first.
The rule must include the extra phi. For a = 2, m = 8 and E = 5, phi(8) = 4 and E mod 4 = 1. The power 2^1 mod 8 is 2, but 2^5 mod 8 is 0. Adding phi gives the exponent 5, which is correct.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ll pw(ll b, ll e, ll m) { // b^e mod m, m >= 1
ll r = 1 % m;
b %= m;
while (e > 0) {
if (e & 1) r = r * b % m;
b = b * b % m;
e >>= 1;
}
return r;
}
static ll phi_of(ll m) {
ll r = m;
for (ll p = 2; p * p <= m; p++)
if (m % p == 0) {
while (m % p == 0) m /= p;
r -= r / p;
}
if (m > 1) r -= r / m;
return r;
}
static bool pow_at_least(ll b, ll c, ll t) {
if (t <= 1) return true;
if (b == 1) return false;
ll cur = 1;
for (ll i = 0; i < c; i++) {
cur *= b;
if (cur >= t) return true;
}
return false;
}
void ref6(istream& in, ostream& out) {
int T;
in >> T;
while (T--) {
ll a, b, c, m;
in >> a >> b >> c >> m;
ll ph = phi_of(m);
ll e = pw(b, c, ph); // b^c mod phi(m)
if (pow_at_least(b, c, ph)) e += ph; // exponent reduction rule
out << pw(a, e, m) << "\n";
}
}
int main() {
ref6(cin, cout);
return 0;
}Checked: the sample prints 2, 561 and 19. A brute force that computes b^c exactly and then multiplies a by itself agreed on 300 random small cases. The trial division step is the slow part in the worst case, about 3 x 10^8 modulo operations for 10^4 queries with m near 10^9, which is near one second in C++.
Problem 7: a geometric sum modulo m
Statement. For each query, print the sum of a to the power i, for i from 0 to n, taken modulo m. The modulus m need not be prime, so you cannot divide by a minus 1.
Input. The first line holds T. Each of the next T lines holds a, n and m.
Output. One value per line.
Limits. 1 ≤ T ≤ 10^4, 0 ≤ a ≤ 10^9, 0 ≤ n ≤ 10^18, 1 ≤ m ≤ 10^9.
Example.
input
2
3 10 1000
2 5 7
output
573
0Hint. Use divide and conquer on the number of terms. For k terms, keep the sum and the power a to the k. Lesson: modular arithmetic and fast power.
Target. O(log n) per query, with no division.
Show the answer
Let S(k) be the sum of a^i for i from 0 to k - 1, and P(k) be a^k, both modulo m. Then S(2h) = S(h) + P(h) times S(h) and P(2h) = P(h) squared. For an odd count, S(k) = S(k - 1) + P(k - 1) and P(k) = P(k - 1) times a. The answer is S(n + 1), because the sum includes i = n. The recursion has depth about 2 log2(n), so it is safe for n up to 10^18.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static pair<ll, ll> geo(ll a, ll k, ll m) {
if (k == 0) return {0, 1 % m};
if (k % 2 == 1) {
auto [s, p] = geo(a, k - 1, m);
return {(s + p) % m, p * (a % m) % m};
}
auto [s, p] = geo(a, k / 2, m);
return {(s + p * s) % m, p * p % m};
}
void ref7(istream& in, ostream& out) {
int T;
in >> T;
while (T--) {
ll a, n, m;
in >> a >> n >> m;
out << geo(a, n + 1, m).first << "\n";
}
}
int main() {
ref7(cin, cout);
return 0;
}Checked: the sample prints 573 and 0. For the first query, the sum 1 + 3 + ... + 3^10 is 88573, which ends in 573. A loop brute force agreed on 300 random small cases, including a up to 10^9.
Problem 8: multiplicative order modulo a prime
Statement. For each query, print the smallest k ≥ 1 such that a to the power k is 1 modulo p.
Input. The first line holds T. Each of the next T lines holds p and a.
Output. One value per line.
Limits. 1 ≤ T ≤ 1000. p is an odd prime with p ≤ 10^9, and 1 ≤ a < p.
Example.
input
3
3 5
7 3
1009 2
output
2
6
504Hint. The order divides p - 1. Start from p - 1, and for each prime factor q of p - 1, divide out q while the power stays 1. Lesson: fermat and Euler in practice.
Target. O(sqrt(p)) per query for the factoring by trial division, plus O(log p) for each prime factor of p - 1.
Show the answer
By Fermat's little theorem, the order divides p - 1. Start with ord = p - 1. For each distinct prime q of p - 1, while q divides ord and a to the power (ord / q) is 1, replace ord by ord / q. This greedy step reaches the exact order, because removing a prime factor is allowed only when the power still equals 1.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ll pw(ll b, ll e, ll m) { // b^e mod m, m >= 1
ll r = 1 % m;
b %= m;
while (e > 0) {
if (e & 1) r = r * b % m;
b = b * b % m;
e >>= 1;
}
return r;
}
static ll order_mod_prime(ll p, ll a) {
ll ord = p - 1, t = p - 1;
vector<ll> qs;
for (ll q = 2; q * q <= t; q++)
if (t % q == 0) {
qs.push_back(q);
while (t % q == 0) t /= q;
}
if (t > 1) qs.push_back(t);
for (ll q : qs)
while (ord % q == 0 && pw(a, ord / q, p) == 1) ord /= q;
return ord;
}
void ref8(istream& in, ostream& out) {
int T;
in >> T;
while (T--) {
ll p, a;
in >> p >> a;
out << order_mod_prime(p, a) << "\n";
}
}
int main() {
ref8(cin, cout);
return 0;
}Checked: the sample prints 2, 6 and 504. A loop that multiplies a until it reaches 1 agreed on 300 random small primes.
Problem 9: sums divided by a value modulo a prime
Statement. You get an array A of N values, each between 0 and p - 1. For each query, take the sum of A from position l to position r, and divide it by x modulo p.
Input. The first line holds N, Q and p. The second line holds the N values. Each of the next Q lines holds l, r and x.
Output. One value per query.
Limits. 1 ≤ N, Q ≤ 10^5, p is an odd prime with p ≤ 10^9 + 7, 0 ≤ A_i < p, and 1 ≤ x < p.
Example.
input
4 2 7
3 1 2 4
1 2 3
2 4 5
output
6
0Hint. Keep prefix sums modulo p. For the division, use the Fermat inverse. Lesson: modular inverse and division.
Target. O(N + Q log p).
Show the answer
Keep prefix sums modulo p. The range sum is (pre[r] - pre[l - 1] + p) mod p, where the extra p keeps the difference non-negative. Divide by x by multiplying with the inverse x^(p - 2) modulo p, which is valid because p is prime and x is not 0 modulo p.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ll pw(ll b, ll e, ll m) { // b^e mod m, m >= 1
ll r = 1 % m;
b %= m;
while (e > 0) {
if (e & 1) r = r * b % m;
b = b * b % m;
e >>= 1;
}
return r;
}
void ref9(istream& in, ostream& out) {
int n, q;
ll p;
in >> n >> q >> p;
vector<ll> pre(n + 1, 0);
for (int i = 1; i <= n; i++) {
ll v;
in >> v;
pre[i] = (pre[i - 1] + v) % p;
}
while (q--) {
int l, r;
ll x;
in >> l >> r >> x;
ll s = (pre[r] - pre[l - 1] + p) % p;
out << s * pw(x, p - 2, p) % p << "\n";
}
}
int main() {
ref9(cin, cout);
return 0;
}Checked: the sample prints 6 and 0. For the first query the sum is 4, and 4 divided by 3 modulo 7 is 6, because 3 times 6 = 18 = 4 (mod 7). A brute force that searches for y with y times x equal to the sum agreed on 300 random small cases.
Problem 10: combining congruences that may share factors
Statement. You get K congruences of the form x ≡ r_i (mod m_i). The moduli can share factors. Print the smallest non-negative x that satisfies all of them, or -1 if none exists.
Input. The first line holds K. Each of the next K lines holds r_i and m_i.
Output. One value.
Limits. 1 ≤ K ≤ 6, 0 ≤ r_i < m_i ≤ 1000. The lcm of the moduli is at most 10^18.
Example.
input
3
2 6
8 10
3 7
output
38Hint. Merge the congruences one at a time. Two congruences are compatible only when their residues agree modulo the gcd of the moduli. Lesson: Chinese remainder theorem.
Target. O(K log m).
Show the answer
Keep the current solution as x ≡ r (mod M). To add a congruence x ≡ ri (mod mi), let g = gcd(M, mi). If ri - r is not a multiple of g, the system has no solution. Otherwise solve M times t ≡ ri - r (mod mi) by dividing everything by g and using the inverse of M / g modulo mi / g. The new solution is r + M times t, with modulus M / g times mi.
The difference ri - r can be negative, and C++ keeps the sign of %. Normalise it with ((d % g) + g) % g before the divisibility check, and normalise t before using it. A negative t gives a result outside the range from 0 to the lcm.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ll ext_gcd(ll a, ll b, ll& x, ll& y) { // returns g = gcd(a, b) and a*x + b*y = g
if (b == 0) {
x = 1;
y = 0;
return a;
}
ll x1, y1;
ll g = ext_gcd(b, a % b, x1, y1);
x = y1;
y = x1 - (a / b) * y1;
return g;
}
void ref10(istream& in, ostream& out) {
int k;
in >> k;
ll M = 1, r = 0; // solution so far: x = r (mod M)
bool ok = true;
for (int i = 0; i < k; i++) {
ll ri, mi;
in >> ri >> mi;
if (!ok) continue;
ll g = std::gcd(M, mi);
ll diff = ri - r;
if (((diff % g) + g) % g != 0) { // the two congruences disagree mod g
ok = false;
continue;
}
ll m2 = mi / g;
ll t = 0;
if (m2 > 1) {
ll x, y;
ext_gcd((M / g) % m2, m2, x, y); // (M/g) * x + m2 * y = 1
ll inv = ((x % m2) + m2) % m2;
ll d = ((diff / g) % m2 + m2) % m2;
t = (ll)((__int128)d * inv % m2);
}
ll newM = M / g * mi;
r = (ll)(((__int128)r + (__int128)M * t) % newM);
M = newM;
}
out << (ok ? r : -1) << "\n";
}
int main() {
ref10(cin, cout);
return 0;
}Checked: the sample prints 38, which is 8 modulo 30 and satisfies 38 mod 6 = 2, 38 mod 10 = 8 and 38 mod 7 = 3. The congruences 2 mod 6, 4 mod 10 and 3 mod 15 have no common solution, and the program prints -1. A brute force that checks x from 0 to the lcm agreed on 300 random small cases.
Problem 11: sum of totients
Statement. Print the sum of phi(i) for i from 1 to N, where phi(i) is Euler's totient.
Input. One line with N.
Output. One integer.
Limits. 1 ≤ N ≤ 10^7.
Example.
input
10
output
32Hint. Sieve the totients up to N with the prime-marking loop. The total is about 3 N squared divided by pi squared, so it needs 64 bits. Lesson: Euler's totient.
Target. O(N log log N).
Show the answer
Start phi[i] = i. For each prime p, which is an index with phi[p] still equal to p, subtract phi[m] / p from every multiple m of p. The sum uses 64-bit integers, because the total for N = 10^7 is about 3 x 10^13.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref11(istream& in, ostream& out) {
int n;
in >> n;
vector<int> phi(n + 1);
for (int i = 0; i <= n; i++) phi[i] = i;
for (int p = 2; p <= n; p++)
if (phi[p] == p) // p is prime
for (int m = p; m <= n; m += p) phi[m] -= phi[m] / p;
ll s = 0;
for (int i = 1; i <= n; i++) s += phi[i];
out << s << "\n";
}
int main() {
ref11(cin, cout);
return 0;
}Checked: the sample prints 32. Our run for N = 10^7 printed 30396356427242 in about 0.37 seconds on our machine. A brute force that counts coprime pairs for each i agreed on 300 random values up to 300.
Problem 12: coprime pairs in a grid
Statement. For each query with n and m, count the pairs (i, j) with 1 ≤ i ≤ n, 1 ≤ j ≤ m and gcd(i, j) = 1.
Input. The first line holds Q. Each of the next Q lines holds n and m.
Output. One count per line.
Limits. 1 ≤ Q ≤ 1000, 1 ≤ n, m ≤ 10^7.
Example.
input
2
3 4
10 10
output
9
63Hint. The count equals the sum of mu(d) times floor(n / d) times floor(m / d) over d. Group d into blocks where both quotients stay the same. Lesson: multiplicative functions and Mobius inversion, and counting by divisors.
Target. O(M) for the Mobius table with M the largest limit, plus O(sqrt(max(n, m))) per query.
Show the answer
Build the Mobius function mu up to the largest limit with a linear sieve, and store its prefix sums. For each query, walk the blocks of d where floor(n / d) and floor(m / d) are both constant. For a block [l, r], the quotients are a = n / l and b = m / l, and r = min(n / a, m / b). Add (prefix[r] - prefix[l - 1]) times a times b. Stop at min(n, m), because beyond that one quotient is 0.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref12(istream& in, ostream& out) {
int q;
in >> q;
vector<pair<ll, ll>> qs(q);
ll mx = 1;
for (auto& [n, m] : qs) {
in >> n >> m;
mx = max({mx, n, m});
}
vector<int> mu(mx + 1, 0), primes;
vector<char> comp(mx + 1, 0);
mu[1] = 1;
for (ll i = 2; i <= mx; i++) { // linear sieve for the Mobius function
if (!comp[i]) {
primes.push_back((int)i);
mu[i] = -1;
}
for (int p : primes) {
if (i * p > mx) break;
comp[i * p] = 1;
if (i % p == 0) {
mu[i * p] = 0;
break;
}
mu[i * p] = -mu[i];
}
}
vector<ll> Mu(mx + 1, 0);
for (ll i = 1; i <= mx; i++) Mu[i] = Mu[i - 1] + mu[i];
for (auto [n, m] : qs) {
ll lim = min(n, m), ans = 0;
for (ll l = 1, r; l <= lim; l = r + 1) {
ll a = n / l, b = m / l;
r = min(n / a, m / b);
ans += (Mu[r] - Mu[l - 1]) * a * b;
}
out << ans << "\n";
}
}
int main() {
ref12(cin, cout);
return 0;
}Checked: the sample prints 9 and 63. A brute force over all pairs agreed on 300 random small queries. Our run with 1000 random queries up to 10^7 took about 0.3 seconds, including the sieve.
Problem 13: pairs with a common factor
Statement. You get N numbers. Count the pairs of positions i < j where gcd(a_i, a_j) is greater than 1.
Input. The first line holds N. The second line holds the N numbers.
Output. One integer.
Limits. 1 ≤ N ≤ 2 x 10^5 and 1 ≤ a_i ≤ 10^6.
Example.
input
5
6 10 15 7 9
output
5Hint. Count the pairs that are coprime with Mobius inversion over the divisor buckets, then subtract from all pairs. Lesson: multiplicative functions and Mobius inversion, and counting by divisors.
Target. O(N + M log M) with M the largest value.
Show the answer
Let c_d be the number of values divisible by d. The number of pairs i < j with gcd exactly 1 is the sum over d of mu(d) times c_d times (c_d - 1) / 2. This is Mobius inversion applied to "both values are divisible by d". The answer is all pairs, N times (N - 1) / 2, minus that count. Equal values count as pairs. Two copies of the value 6 have gcd 6, so they are not coprime.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
void ref13(istream& in, ostream& out) {
int n;
in >> n;
vector<int> a(n);
int M = 1;
for (int& x : a) {
in >> x;
M = max(M, x);
}
vector<ll> freq(M + 1, 0);
for (int x : a) freq[x]++;
vector<int> mu(M + 1, 0), primes;
vector<char> comp(M + 1, 0);
mu[1] = 1;
for (int i = 2; i <= M; i++) {
if (!comp[i]) {
primes.push_back(i);
mu[i] = -1;
}
for (int p : primes) {
if ((ll)i * p > M) break;
comp[i * p] = 1;
if (i % p == 0) {
mu[i * p] = 0;
break;
}
mu[i * p] = -mu[i];
}
}
ll coprime = 0;
for (int d = 1; d <= M; d++) {
if (mu[d] == 0) continue;
ll c = 0; // divisor bucket: how many a_i are divisible by d
for (int m = d; m <= M; m += d) c += freq[m];
coprime += mu[d] * (c * (c - 1) / 2);
}
ll total = (ll)n * (n - 1) / 2;
out << total - coprime << "\n";
}
int main() {
ref13(cin, cout);
return 0;
}Checked: the sample prints 5. The pairs with a common factor are (6, 10), (6, 15), (6, 9), (10, 15) and (9, 15). A brute force over pairs agreed on 300 random small arrays.
Problem 14: number of divisors of a large number
Statement. For each query, print the number of divisors of n.
Input. The first line holds Q. Each of the next Q lines holds n.
Output. One count per line.
Limits. 1 ≤ Q ≤ 100 and 1 ≤ n ≤ 10^18.
Example.
input
3
1
999999937
1000000000000000000
output
1
2
361Hint. Remove the small primes by trial division. For the rest, test primality with Miller-Rabin, and split composites with Pollard rho. Multiply the exponents plus one over the prime factors. Lesson: Miller-Rabin and Pollard rho.
Target. O(Q times n^(1/4)) expected steps for the splitting, with a small constant.
Show the answer
Trial-divide by primes below 1000. If the rest is greater than 1, factor it. The primality test is Miller-Rabin with the bases 2 to 37. These are the first twelve primes, a set that is widely used as deterministic for inputs below 2^64. We did not derive that bound ourselves. The splitting uses Pollard rho with Floyd cycle finding. If the result is n itself, try the next constant c. The products in the modular steps use unsigned __int128, because a product of two 64-bit numbers needs 128 bits.
Once the prime factors are sorted, count each prime's exponent e and multiply the answer by e + 1.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
static ull mulmod(ull a, ull b, ull m) { return (unsigned __int128)a * b % m; }
static ull powmod(ull a, ull e, ull m) {
ull r = 1 % m;
a %= m;
while (e) {
if (e & 1) r = mulmod(r, a, m);
a = mulmod(a, a, m);
e >>= 1;
}
return r;
}
static bool is_prime(ull n) {
if (n < 2) return false;
static const ull bases[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
for (ull p : bases)
if (n % p == 0) return n == p;
ull d = n - 1;
int s = 0;
while (d % 2 == 0) {
d /= 2;
s++;
}
for (ull a : bases) {
ull x = powmod(a, d, n);
if (x == 1 || x == n - 1) continue;
bool composite = true;
for (int r = 1; r < s; r++) {
x = mulmod(x, x, n);
if (x == n - 1) {
composite = false;
break;
}
}
if (composite) return false;
}
return true;
}
static ull rho(ull n) { // returns a non-trivial divisor of an odd composite n
if (n % 2 == 0) return 2;
for (ull c = 1;; c++) {
ull x = 2, y = 2, d = 1;
auto f = [&](ull v) { return (mulmod(v, v, n) + c) % n; };
while (d == 1) {
x = f(x);
y = f(f(y));
d = std::gcd(x > y ? x - y : y - x, n);
}
if (d != n) return d;
}
}
static void factor_big(ull n, vector<ull>& out) {
if (n == 1) return;
if (is_prime(n)) {
out.push_back(n);
return;
}
ull d = rho(n);
factor_big(d, out);
factor_big(n / d, out);
}
void ref14(istream& in, ostream& out) {
int q;
in >> q;
while (q--) {
ull n;
in >> n;
vector<ull> pf;
for (ull p = 2; p < 1000 && p * p <= n; p++)
while (n % p == 0) {
pf.push_back(p);
n /= p;
}
if (n > 1) factor_big(n, pf);
sort(pf.begin(), pf.end());
ull cnt = 1;
size_t i = 0;
while (i < pf.size()) {
size_t j = i;
while (j < pf.size() && pf[j] == pf[i]) j++;
cnt *= (j - i + 1);
i = j;
}
out << cnt << "\n";
}
}
int main() {
ref14(cin, cout);
return 0;
}Checked: the sample prints 1, 2 and 361. For 999999937 times 999999929, a product of two primes near 10^9, the program printed 4. Both factors were checked with trial division in Python. A trial division brute force agreed on 300 random inputs up to 10^9, which covers composites with two factors above 1000. Our run with 100 random numbers up to 10^18 took about 8 milliseconds.
Common mistakes and the symptoms they cause
- A gcd fold that starts at 1. The value gcd(1, x) is 1 for every x, so the output is 1 for every input. Start the fold at 0.
- A 32-bit sum. Sums for N = 10^7 in problem 11 and n = 10^12 in the lesson are above 2 x 10^9. The symptom is a wrong value or a negative number.
- A sieve sized for the sample. Arrays built from the sample's largest value overflow on the real tests. The symptom is a crash or results that differ between runs.
- Dividing by a non-prime modulus. Problem 7 has composite m. Using the formula with a minus 1 in the denominator is invalid, and the symptom is a wrong answer for some m.
- Forgetting the extra phi in problem 6. The symptom shows only for some inputs, such as a = 2 and m = 8 with a large exponent.
- Negative remainders. In C++,
%keeps the sign of the left operand. Problem 10 can produce a negative t if it is not normalised. The symptom is -1 or a value outside the range. - Empty containers in generators. While we built the random tests, a static vector was emptied on each call, so
r() % sizedivided by zero. The program stopped with exit code 136 and printed nothing. Check container sizes before you take a remainder. - Testing only the sample. Each sample here is small enough that most bugs stay hidden. Use random tests and the largest limits.
Checklist for every problem
- Write the limits and the largest value next to the statement.
- Write a brute force first and check it on the sample.
- Use 64-bit integers for any sum or product that can exceed 2 x 10^9.
- Check which modulus is prime before you divide.
- Keep the answer in range after every negative remainder.
- Time the fast solution at the largest allowed input.
- Compare with the complexity target and note any gap.
What we checked
Every reference program was compiled with g++ -std=c++17 -O2 -Wall -Wextra using g++ 15.2.0, with no warnings. Each one was compared with its brute force on 300 random small inputs. The generators produced inputs within the statements' limits at small size. The samples in each statement are the outputs of the reference programs, and we checked the small ones by hand. The timings are from one machine, and your numbers will differ. We did not run the solutions on any judge.
The Miller-Rabin base set is a claim we did not verify ourselves. We tested the code on random inputs up to 10^9 and on one product of two primes near 10^9. The claim that these bases are deterministic below 2^64 is one we found in a Python bug-tracker comment. We did not check it against the original paper, so check it before you rely on it.
Summary
The problems cover the main tools of the module in rising order. Problems 1 to 5 use gcd, sieves and divisor counts. Problems 6 to 10 use modular powers, orders, inverses and the Chinese remainder theorem. Problems 11 to 14 use totients, Mobius sums, divisor buckets and factoring of large numbers. For every problem, the reference solution has a complexity target, and the random tests give the evidence that it is correct.