Essential Number Theory and Linear Algebra Algorithms
- Number Theory Fundamentals
Extended Euclidean Algorithm and Linear Diophantine Equations
For integers a and b, the equation ax + by = d has integer solutions if and only if the greatest common divisor gcd(a, b) divides d. This is known as Bézout's Identity. The extended Euclidean algorithm alows us to find x and y such that ax + by = gcd(a, b).
The logic relies on the observation that ax + by = bx' + (a mod b)y'. Substituting a mod b = a - floor(a/b) * b, we can derive the recursive relationship for the coefficients.
long long extended_gcd(long long a, long long b, long long& coeff_a, long long& coeff_b) {
if (b == 0) {
// Base case: a * 1 + 0 * 0 = a
coeff_a = 1;
coeff_b = 0;
return a;
}
long long x1, y1;
long long gcd_val = extended_gcd(b, a % b, x1, y1);
// Update coefficients based on the recursive step
coeff_a = y1;
coeff_b = x1 - (a / b) * y1;
return gcd_val;
}
Modular Arithmetic and Multiplicative Inverses
If m divides (a - b), we say a ≡ b (mod m). Modular arithmetic preserves properties under addition, subtraction, and multiplication. However, division is distinct; it requires the concept of a multiplicative inverse.
The multiplicative inverse of x modulo p, denoted x-1, satisfies x * x-1 ≡ 1 (mod p). To solve ax ≡ 1 (mod b), we use the extended Euclidean algorithm to find a solution to ax + by = 1.
#include <iostream>
using namespace std;
long long extended_gcd(long long a, long long b, long long& x, long long& y) {
if (b == 0) { x = 1; y = 0; return a; }
long long g = extended_gcd(b, a % b, x, y);
long long t = x; x = y; y = t - a / b * y;
return g;
}
int main() {
long long a, mod;
if (!(cin >> a >> mod)) return 0;
long long x, y;
extended_gcd(a, mod, x, y);
// Normalize the result to the smallest positive integer
long long inverse = (x % mod + mod) % mod;
cout << inverse << endl;
return 0;
}
Prime Number Sieves
The Sieve of Eratosthenes is a efficient method to find all primes up to n with O(n log log n) complexity. The Linear Sieve (or Euler's Sieve) improves upon this by ensuring each composite number is marked only once by its smallest prime factor, resulting in O(n) complexity.
#include <vector>
#include <iostream>
using namespace std;
const int MAX_N = 10005;
vector<int> primes;
bool is_composite[MAX_N];
void linear_sieve(int limit) {
for (int i = 2; i <= limit; ++i) {
if (!is_composite[i]) {
primes.push_back(i);
}
for (int prime : primes) {
long long multiple = (long long)i * prime;
if (multiple > limit) break;
is_composite[multiple] = true;
if (i % prime == 0) break; // Key: ensures smallest prime factor marks the number
}
}
}
int main() {
int n = 100;
linear_sieve(n);
for (int p : primes) cout << p << " ";
return 0;
}
Factorial Prime Factorization
To decompose n! into prime factors, we calculate the exponent of each prime p. The exponent is given by the sum: floor(n/p) + floor(n/p2) + .... This represents the count of multiples of p, p2, etc., in 1 to n.
- Linear Algebra and Exponentiation
Fast Exponentiation (Binary Exponentiation)
This technique computes ab mod m in O(log b) time by decomposing the exponent b into its binary representation.
long long fast_pow(long long base, long long exp, long long mod) {
long long result = 1;
base %= mod;
while (exp > 0) {
if (exp & 1) result = (result * base) % mod;
base = (base * base) % mod;
exp >>= 1;
}
return result;
}
Gaussian Elimination
Gaussian elimination solves systems of linear equations Ax = B. It transforms the augmented matrix [A|B] into row-echelon form using elementary row operations (swapping rows, multiplying a row by a scalar, adding a multiple of one row to another).
Application: Finding a Sphere's Center
Given n+1 points on an n-dimensional sphere, the center can be found by solving linear equations derived from equating the squared distances from the center to any two points. This transforms quadratic equations into linear ones.
Application: Equipment Selection
Given equipment with attribute vectors and costs, we want to maximize the number of items purchased while minimizing cost, ensuring selected attribute vectors are linearly independent. We use Gaussian elimination with a greedy strategy: when selecting a pivot row, we choose the cheapest available option.
#include <iostream>
#include <vector>
#include <cmath>
using namespace std;
const double EPS = 1e-9;
const int MAX_DIM = 510;
int n, m;
double attributes[MAX_DIM][MAX_DIM]; // Matrix A
int costs[MAX_DIM]; // Cost vector
int solve_equipment() {
int current_rank = 0;
int total_cost = 0;
for (int col = 0; col < m && current_rank < n; ++col) {
int selected_row = -1;
// Greedy: Find the row with non-zero coefficient in this column and minimal cost
for (int row = current_rank; row < n; ++row) {
if (fabs(attributes[row][col]) > EPS &&
(selected_row == -1 || costs[row] < costs[selected_row])) {
selected_row = row;
}
}
if (selected_row == -1) continue; // No valid vector for this column
// Swap current row with selected row
swap(attributes[current_rank], attributes[selected_row]);
swap(costs[current_rank], costs[selected_row]);
total_cost += costs[current_rank];
// Eliminate this variable from other rows
for (int row = 0; row < n; ++row) {
if (row != current_rank && fabs(attributes[row][col]) > EPS) {
double factor = attributes[row][col] / attributes[current_rank][col];
for (int k = col; k < m; ++k) {
attributes[row][k] -= attributes[current_rank][k] * factor;
}
}
}
current_rank++;
}
cout << current_rank << " " << total_cost << endl;
return 0;
}
- Permutations, Combinations, and Dynamic Programming
Theoretical Foundations
Key formulas for permutations P(n, m) and combinations C(n, m) are standard. An important property for modular arithmetic is Lucas' Theorem, which allows computing C(n, m) modulo a prime p by decomposing n and m in to base p digits.
Combinatorial Grouping Problem
Problem: Group n distinguishable people. Person i must be in a group of size at most ai. Groups are indistinguishable, and order within a group does not matter. Find the number of valid grouping schemes modulo 998244353.
Analysis: We use Dynamic Programming. Let dp[i][j] be the number of ways to form groups of sizes strictly greater than i using j remaining people. We process group sizes from large to small.
#include <iostream>
#include <vector>
using namespace std;
const long long MOD = 998244353;
const int MAX_N = 2010;
int n;
int limits[MAX_N]; // a[i]
int freq[MAX_N]; // frequency of limit = i
long long dp[MAX_N][MAX_N];
long long comb[MAX_N][MAX_N];
long long fact[MAX_N], inv_fact[MAX_N];
long long fast_pow(long long base, long long exp) {
long long res = 1;
while (exp) {
if (exp & 1) res = res * base % MOD;
base = base * base % MOD;
exp >>= 1;
}
return res;
}
void precompute() {
// Combinations
for (int i = 0; i <= n; ++i) {
comb[i][0] = 1;
for (int j = 1; j <= i; ++j) {
comb[i][j] = (comb[i-1][j] + comb[i-1][j-1]) % MOD;
}
}
// Factorials and Inverses
fact[0] = 1;
for (int i = 1; i <= n; ++i) fact[i] = fact[i-1] * i % MOD;
inv_fact[n] = fast_pow(fact[n], MOD - 2);
for (int i = n - 1; i >= 0; --i) inv_fact[i] = inv_fact[i + 1] * (i + 1) % MOD;
}
int main() {
ios::sync_with_stdio(false);
cin >> n;
for (int i = 1; i <= n; ++i) {
int val; cin >> val;
freq[val]++;
}
precompute();
// Base case: group size > n means 0 people left
dp[n + 1][0] = 1;
for (int size = n; size >= 1; --size) {
long long inv_pow_size = 1;
// Iterate number of groups of current 'size'
for (int k = 0; size * k <= n; ++k) {
// Term: (size*k)! / ((size!)^k * k!)
long long ways_to_partition = fact[size * k] * inv_pow_size % MOD * inv_fact[k] % MOD;
for (int rem = 0; rem + size * k <= n; ++rem) {
// Transition: choose size*k people from remaining, partition them, add to previous state
long long ways_choose = comb[rem + size * k][size * k];
dp[size][rem] = (dp[size][rem] +
dp[size + 1][rem + size * k - freq[size]] * ways_choose % MOD *
ways_to_partition % MOD) % MOD;
}
inv_pow_size = inv_pow_size * inv_fact[size] % MOD;
}
}
cout << dp[1][0] << endl;
return 0;
}
