数学与数论基础

#算法 · #数学 · #GCD · #素数 · #快速幂 · #组合数学

算法题目中涉及数学的地方远比想象中多。GCD/LCM、素数、快速幂、组合数学——这些不是"纯数学"，而是编程中的高频工具。

---

1. 最大公约数与最小公倍数

1.1 欧几里得算法（GCD）

// 辗转相除法：O(log min(a,b))
func gcd(a, b int) int {
    for b != 0 {
        a, b = b, a%b
    }
    return a
}

// 递归版本
func gcdRec(a, b int) int {
    if b == 0 {
        return a
    }
    return gcdRec(b, a%b)
}

证明思路：gcd(a, b) = gcd(b, a % b)。因为 a = qb + r，任何能整除 a 和 b 的数必能整除 r。

1.2 LCM（最小公倍数）

func lcm(a, b int) int {
    return a / gcd(a, b) * b   // 先除后乘防溢出
}

1.3 扩展欧几里得算法

求 ax + by = gcd(a, b) 的整数解（x, y），用于解模线性方程和求逆元。

func exgcd(a, b int) (gcd, x, y int) {
    if b == 0 {
        return a, 1, 0
    }
    gcd, x1, y1 := exgcd(b, a%b)
    x = y1
    y = x1 - (a/b)*y1
    return
}

---

2. 素数

2.1 素数判定 O(√n)

func isPrime(n int) bool {
    if n <= 1 {
        return false
    }
    if n <= 3 {
        return true
    }
    if n%2 == 0 || n%3 == 0 {
        return false
    }
    // 6k ± 1 优化
    for i := 5; i*i <= n; i += 6 {
        if n%i == 0 || n%(i+2) == 0 {
            return false
        }
    }
    return true
}

2.2 埃氏筛 O(n log log n)

func eratosthenes(n int) []bool {
    isPrime := make([]bool, n+1)
    for i := 2; i <= n; i++ {
        isPrime[i] = true
    }
    for i := 2; i*i <= n; i++ {
        if isPrime[i] {
            for j := i * i; j <= n; j += i {
                isPrime[j] = false
            }
        }
    }
    return isPrime
}

2.3 欧拉线性筛 O(n)

func eulerSieve(n int) []int {
    primes := []int{}
    isComposite := make([]bool, n+1)
    for i := 2; i <= n; i++ {
        if !isComposite[i] {
            primes = append(primes, i)
        }
        for _, p := range primes {
            if i*p > n {
                break
            }
            isComposite[i*p] = true
            if i%p == 0 {   // 保证每个合数只被其最小质因子筛去
                break
            }
        }
    }
    return primes
}

---

3. 快速幂与模运算

3.1 快速幂（二分幂）

// x^n mod m: O(log n)
func powMod(x, n, m int64) int64 {
    result := int64(1)
    x %= m
    for n > 0 {
        if n&1 == 1 {
            result = result * x % m
        }
        x = x * x % m
        n >>= 1
    }
    return result
}

3.2 模逆元

当 gcd(a, m) = 1 时，存在 x 使 a*x ≡ 1 (mod m)。

// 费马小定理：当 m 为素数时，a^(m-2) ≡ a^(-1) (mod m)
func modInverse(a, m int64) int64 {
    return powMod(a, m-2, m)
}

// 通用方法：扩展欧几里得
func modInverseEx(a, m int) int {
    _, x, _ := exgcd(a, m)
    return (x%m + m) % m
}

3.3 大数阶乘模素数

// n! mod p, 使用威尔逊定理和分块
// 对于 p ≤ 1e7: 预处理阶乘
// 对于大 p: Lucas 定理

---

4. 组合数学

4.1 组合数公式

C(n, k) = n! / (k! * (n-k)!)
C(n, k) = C(n-1, k-1) + C(n-1, k)   // 杨辉三角递推
C(n, k) = C(n, n-k)                   // 对称性

4.2 组合数计算

// 方法1: 杨辉三角 — O(n²), 适合 n ≤ 5000
func preprocessPascal(n int) [][]int {
    C := make([][]int, n+1)
    for i := range C {
        C[i] = make([]int, n+1)
        C[i][0], C[i][i] = 1, 1
        for j := 1; j < i; j++ {
            C[i][j] = C[i-1][j-1] + C[i-1][j]
        }
    }
    return C
}

// 方法2: 阶乘 + 逆元 — O(n + log MOD), 适合 MOD 为素数
var fact, invFact []int64
const MOD = 1_000_000_007

func preprocess(n int) {
    fact = make([]int64, n+1)
    invFact = make([]int64, n+1)
    fact[0] = 1
    for i := 1; i <= n; i++ {
        fact[i] = fact[i-1] * int64(i) % MOD
    }
    invFact[n] = powMod(fact[n], MOD-2, MOD)
    for i := n - 1; i >= 0; i-- {
        invFact[i] = invFact[i+1] * int64(i+1) % MOD
    }
}

func nCr(n, r int) int64 {
    if r < 0 || r > n {
        return 0
    }
    return fact[n] * invFact[r] % MOD * invFact[n-r] % MOD
}

4.3 卡特兰数

C₀ = 1
Cₙ₊₁ = Cₙ * 2*(2n+1) / (n+2)

C₁ = 1, C₂ = 2, C₃ = 5, C₄ = 14, C₅ = 42, ...

应用：n 对括号的合法序列数 = Cₙ
      n 个节点的不同 BST 数 = Cₙ
      不穿过对角线的网格路径数 = Cₙ

func catalan(n int) int64 {
    // C(n) = C(2n, n) / (n+1)
    return nCr(2*n, n) * powMod(int64(n+1), MOD-2, MOD) % MOD
}

---

4.4 中国剩余定理（CRT）

中国剩余定理解决"同余方程组"问题——在密码学（RSA 加速）和算法竞赛（大数模运算拆分）中十分重要。

问题定义

求 x 满足:
  x ≡ a₁ (mod m₁)
  x ≡ a₂ (mod m₂)
  ...
  x ≡ aₖ (mod mₖ)

其中 m₁, m₂, ..., mₖ 两两互质。

解法

令 M = m₁ × m₂ × ... × mₖ
令 Mᵢ = M / mᵢ
令 tᵢ = Mᵢ⁻¹ (mod mᵢ)  （Mᵢ 模 mᵢ 的逆元）

解:  x = Σ aᵢ × Mᵢ × tᵢ   (mod M)

验证: x mod mⱼ = aⱼ，因为:
  - 当 i ≠ j: Mᵢ 包含 mⱼ 作为因子 → aᵢ×Mᵢ×tᵢ mod mⱼ = 0
  - 当 i = j: Mⱼ×tⱼ ≡ 1 (mod mⱼ) → aⱼ×Mⱼ×tⱼ ≡ aⱼ (mod mⱼ)

// 中国剩余定理 — O(k log M)
func crt(a, m []int64) int64 {
    // 计算 M = m₁ × m₂ × ... × mₖ
    M := int64(1)
    for _, mi := range m {
        M *= mi
    }

    var x int64
    for i := range a {
        Mi := M / m[i]
        // 求 Mi 模 m[i] 的逆元
        ti := modInverseEx(int(Mi), int(m[i]))
        x = (x + a[i]*Mi*int64(ti)) % M
    }
    return (x + M) % M
}

// 示例：
// x ≡ 2 (mod 3), x ≡ 3 (mod 5), x ≡ 2 (mod 7)
// 解: x = 23

应用

场景	说明
RSA 解密加速	用 CRT 将 cᵈ mod n 拆为 cᵈ mod p 和 cᵈ mod q，提速 ~4x
大数组合数	将 MOD 拆为素数幂乘积，用 CRT 合并结果
哈希碰撞	有时用 CRT 合并多个小哈希的冲突解决

---

4.5 欧拉函数与约数

欧拉函数 φ(n)

定义：φ(n) = [1, n] 中与 n 互质的数的个数。

性质 1（乘积性）: 若 gcd(a, b) = 1，则 φ(ab) = φ(a) × φ(b)
性质 2（素数幂）: φ(pᵏ) = pᵏ - pᵏ⁻¹ = pᵏ⁻¹(p - 1)

计算: φ(n) = n × Π(1 - 1/pᵢ)，其中 pᵢ 是 n 的不同质因子

// 单个欧拉函数 O(√n)
func phi(n int64) int64 {
    result := n
    for i := int64(2); i*i <= n; i++ {
        if n%i == 0 {
            for n%i == 0 {
                n /= i
            }
            result = result / i * (i - 1)  // result *= (1 - 1/i)
        }
    }
    if n > 1 {
        result = result / n * (n - 1)
    }
    return result
}

// 线性筛求 [1, n] 的 φ 值 — O(n)
func phiRange(n int) []int {
    phi := make([]int, n+1)
    primes := []int{}
    isComposite := make([]bool, n+1)
    phi[1] = 1

    for i := 2; i <= n; i++ {
        if !isComposite[i] {
            primes = append(primes, i)
            phi[i] = i - 1  // 素数 p: φ(p) = p-1
        }
        for _, p := range primes {
            if i*p > n {
                break
            }
            isComposite[i*p] = true
            if i%p == 0 {
                phi[i*p] = phi[i] * p       // p 是 i 的因子
                break
            }
            phi[i*p] = phi[i] * (p - 1)     // p 与 i 互质
        }
    }
    return phi
}

约数和与约数个数

// 1 ~ n 每个数的约数个数和 — 调和级数
// Σ(d(n)) ≈ n log n + (2γ-1)n  （γ ≈ 0.577，欧拉常数）

// 批量求约数
func getAllDivisors(n int) [][]int {
    divs := make([][]int, n+1)
    for i := 1; i <= n; i++ {
        for j := i; j <= n; j += i {
            divs[j] = append(divs[j], i)
        }
    }
    return divs
}
// 复杂度：O(n log n)，因为调和级数 Σ n/i ≈ n log n

---

4.6 矩阵快速幂

矩阵快速幂是普通快速幂的推广——将"数的幂"替换为"方阵的幂"，底数运算符从乘法变为矩阵乘法。这使 O(n) 的 Fibonacci 数列降到 O(log n)。

// 2x2 矩阵乘法
func matMul(a, b [2][2]int64, mod int64) [2][2]int64 {
    return [2][2]int64{
        {
            (a[0][0]*b[0][0] + a[0][1]*b[1][0]) % mod,
            (a[0][0]*b[0][1] + a[0][1]*b[1][1]) % mod,
        },
        {
            (a[1][0]*b[0][0] + a[1][1]*b[1][0]) % mod,
            (a[1][0]*b[0][1] + a[1][1]*b[1][1]) % mod,
        },
    }
}

// 矩阵快速幂 O(k³ log n)，k 为矩阵宽度
func matPow(base [2][2]int64, n int64, mod int64) [2][2]int64 {
    result := [2][2]int64{{1, 0}, {0, 1}} // 单位矩阵
    b := base
    for n > 0 {
        if n&1 == 1 {
            result = matMul(result, b, mod)
        }
        b = matMul(b, b, mod)
        n >>= 1
    }
    return result
}

// Fibonacci 数列 O(log n)
// [F(n+1), F(n)]^T = [[1,1],[1,0]]^n × [F(1), F(0)]^T
// [Fib(n+1)]   = [[1,1]]ⁿ  ×  [1]
// [ Fib(n) ]     [[1,0]]      [0]
func fib(n int64, mod int64) int64 {
    if n <= 1 {
        return n
    }
    base := [2][2]int64{{1, 1}, {1, 0}}
    result := matPow(base, n-1, mod)
    return result[0][0]  // 即 F(n)
}

矩阵快速幂的推广——线性递推

任意常系数 k 阶线性递推:
  a(n) = c₁×a(n-1) + c₂×a(n-2) + ... + cₖ×a(n-k)

都可以用 k×k 矩阵快速幂在 O(k³ log n) 求解。

例如：a(n) = 2a(n-1) + a(n-3)
  [a(n)  ]   [[2,0,1]]ⁿ⁻³   [a(3)]
  [a(n-1)] = [[1,0,0]]    × [a(2)]
  [a(n-2)]   [[0,1,0]]      [a(1)]

---

5. 概率与期望

算法竞赛中的概率期望问题通常不是让你算数值，而是让你推导期望的线性性质来简化计算。

5.1 期望的线性性质

E[X + Y] = E[X] + E[Y]，即使 X 和 Y 不独立也成立！

这是算法中最重要的期望性质，可将复杂问题拆解。

// 示例：掷 n 次骰子，期望有几个 6？
// E[每个骰子出现6] = 1/6，n个 → E[总数] = n/6
// 不需要枚举所有 6ⁿ 种情况！

5.2 贡献法 — 期望的核心技巧

问题：给定一个排列，求所有连续子数组的最大值之和的期望。

暴力：O(n²) 枚举所有子数组
贡献法：计算每个元素 a[i] 作为最大值出现在哪些子数组中。

步骤：
1. 找到左边第一个大于 a[i] 的位置 L（单调栈）
2. 找到右边第一个大于 a[i] 的位置 R
3. 以 a[i] 为最大值的子数组数 = (i-L) × (R-i)
4. 贡献 = a[i] × (i-L) × (R-i)

// 贡献法模板 — O(n)
func sumOfSubarrayMaximums(arr []int) int64 {
    n := len(arr)
    left := make([]int, n)   // 左边界（第一个 > arr[i] 的位置）
    right := make([]int, n)  // 右边界（第一个 > arr[i] 的位置）

    // 单调栈求左右边界
    stack := []int{}
    for i := 0; i < n; i++ {
        for len(stack) > 0 && arr[stack[len(stack)-1]] <= arr[i] {
            stack = stack[:len(stack)-1]
        }
        if len(stack) == 0 {
            left[i] = -1
        } else {
            left[i] = stack[len(stack)-1]
        }
        stack = append(stack, i)
    }
    // ... right 边界类似 ...

    var total int64
    for i := 0; i < n; i++ {
        count := int64(i-left[i]) * int64(right[i]-i)
        total += int64(arr[i]) * count
    }
    return total
}

5.3 马尔可夫链与状态转移

问题：在一个 n 个节点的图上随机游走，第一次到达目标节点的期望步数？

解法：设 E[u] = 从 u 出发到达目标的期望步数

转移方程：
  E[target] = 0
  E[u] = 1 + Σ (E[v] / deg(u))   对于每条边 (u, v)

这是马尔可夫链期望，可用高斯消元 O(n³) 解线性方程组。

5.4 常见概率问题速查

问题类型	核心技巧	复杂度
随机排列期望	期望线性性，每个元素独立贡献	O(n)
几何分布	E = 1/p（首次成功期望）	O(1)
抛硬币	状态机 + DP	O(n)
随机游走	马尔可夫链 + 高斯消元	O(n³)
蓄水池抽样	数学归纳法证明等概率	O(n)
Monte Carlo	随机采样近似计算	与精度相关

---

6. 其他常见数学工具

6.1 整数平方根

func sqrt(n int) int {
    x := n
    for x*x > n {
        x = (x + n/x) / 2
    }
    return x
}

6.2 蓄水池抽样

// 从数据流中等概率选 k 个元素
func reservoirSample(stream []int, k int) []int {
    res := make([]int, k)
    copy(res, stream[:k])
    for i := k; i < len(stream); i++ {
        j := rand.Intn(i + 1)
        if j < k {
            res[j] = stream[i]
        }
    }
    return res
}

6.3 洗牌算法（Fisher-Yates）

func shuffle(arr []int) {
    for i := len(arr) - 1; i > 0; i-- {
        j := rand.Intn(i + 1)
        arr[i], arr[j] = arr[j], arr[i]
    }
}

---

7. 常用数学常数与公式

名称	公式 / 值	用途
等差数列和	n(a₁+aₙ)/2	连续序列求和
等比数列和	a₁(1-rⁿ)/(1-r)	倍增场景
平方和	n(n+1)(2n+1)/6	DP 推导
异或性质	1⊕2⊕3⊕...⊕n 有 O(1) 公式	异或问题
容斥原理	|A∪B| = |A|+|B|-|A∩B|	计数去重

---

参考

• CP-Algorithms
• CLRS（《算法导论》）第 31 章：数论算法
• OEIS（整数数列在线百科）