每日算法 — 使用java实现莫比乌斯反演：数论分块与积性函数线性筛优化

可莉爸爸 AI算法

引言

莫比乌斯反演是数论中一项极为优美的技术，它建立了两个数论函数之间的对偶关系，广泛应用于组合计数、容斥原理和算法竞赛中。本文将用Java实现莫比乌斯反演及其核心计算框架，重点讲解数论分块与积性函数线性筛两大优化技术。

理论基础

狄利克雷卷积与莫比乌斯函数

对于两个数论函数 $f$ 和 $g$，其狄利克雷卷积定义为：
$$(f * g)(n) = \sum_{d|n} f(d) \cdot g\left(\frac{n}{d}\right)$$

莫比乌斯函数 $\mu(n)$ 定义为：
– $\mu(1) = 1$
– 若 $n$ 含有平方因子，则 $\mu(n) = 0$
– 若 $n$ 是 $k$ 个不同素数的乘积，则 $\mu(n) = (-1)^k$

核心定理：若 $F(n) = \sum_{d|n} f(d)$，则 $f(n) = \sum_{d|n} \mu(d) \cdot F\left(\frac{n}{d}\right)$

线性筛法预处理

线性筛（欧拉筛）可在 $O(n)$ 时间内同时预处理出每个数的最小质因子、莫比乌斯函数值和欧拉函数值。

import java.util.*;

/**
 * 数论工具类：线性筛预处理 + 莫比乌斯反演
 */
public class NumberTheory {

    private int maxN;
    private boolean[] isComposite;
    private int[] primes;        // 素数列表
    private int[] mu;            // 莫比乌斯函数
    private int[] phi;           // 欧拉函数
    private int[] minPrime;      // 最小质因子
    private int primeCount;

    /**
     * 线性筛预处理
     * 时间复杂度: O(n)
     * 同时计算：素数表、莫比乌斯函数、欧拉函数、最小质因子
     */
    public NumberTheory(int n) {
        this.maxN = n;
        isComposite = new boolean[n + 1];
        primes = new int[n + 1];
        mu = new int[n + 1];
        phi = new int[n + 1];
        minPrime = new int[n + 1];
        primeCount = 0;

        mu[1] = 1;
        phi[1] = 1;

        for (int i = 2; i <= n; i++) {
            if (!isComposite[i]) {
                primes[primeCount++] = i;
                mu[i] = -1;        // 素数: μ(p) = -1
                phi[i] = i - 1;   // φ(p) = p - 1
                minPrime[i] = i;
            }

            for (int j = 0; j < primeCount && i * primes[j] <= n; j++) {
                int p = primes[j];
                isComposite[i * p] = true;
                minPrime[i * p] = p;

                if (i % p == 0) {
                    // p 是 i 的最小质因子
                    mu[i * p] = 0;  // 含平方因子
                    phi[i * p] = phi[i] * p; // φ(i*p) = φ(i) * p
                    break;
                } else {
                    mu[i * p] = -mu[i]; // 增加一个不同质因子，符号翻转
                    phi[i * p] = phi[i] * (p - 1); // φ(i*p) = φ(i) * φ(p)
                }
            }
        }
    }

    public int getMu(int n) { return mu[n]; }
    public int getPhi(int n) { return phi[n]; }
    public int[] getPrimes() { return Arrays.copyOf(primes, primeCount); }

    /**
     * 打印前n个数的莫比乌斯函数值
     */
    public void printMuTable(int n) {
        System.out.println("n\tμ(n)\tφ(n)");
        for (int i = 1; i <= n; i++) {
            System.out.printf("%d\t%d\t%d%n", i, mu[i], phi[i]);
        }
    }

莫比乌斯反演核心实现

    /**
     * 莫比乌斯反演：由F(n)求f(n)
     * 前提: F(n) = sum_{d|n} f(d)
     * 则:   f(n) = sum_{d|n} μ(d) * F(n/d)
     * @param F 已知函数F的值数组（F[i] = F(i)）
     * @return f数组
     */
    public long[] mobiusInversion(long[] F) {
        int n = F.length - 1;
        long[] f = new long[n + 1];

        for (int i = 1; i <= n; i++) {
            long sum = 0;
            // 枚举i的所有因子d
            for (int d = 1; d * d <= i; d++) {
                if (i % d == 0) {
                    sum += mu[d] * F[i / d];
                    if (d * d != i) {
                        sum += mu[i / d] * F[d];
                    }
                }
            }
            f[i] = sum;
        }

        return f;
    }

    /**
     * 莫比乌斯反演（倍数形式）
     * 前提: F(n) = sum_{n|d} f(d)
     * 则:   f(n) = sum_{n|d} μ(d/n) * F(d)
     * @param F 已知函数F的值数组
     * @return f数组
     */
    public long[] mobiusInversionMultiple(long[] F) {
        int n = F.length - 1;
        long[] f = new long[n + 1];
        System.arraycopy(F, 0, f, 0, n + 1);

        // 倒序进行：从大到小消去贡献
        for (int i = n; i >= 1; i--) {
            for (int j = 2 * i; j <= n; j += i) {
                f[i] -= f[j];
            }
        }

        return f;
    }

数论分块优化

数论分块（也称为整除分块）是处理 $\sum_{i=1}^{n} \lfloor \frac{n}{i} \rfloor$ 类求和的核心优化，可将 $O(n)$ 优化到 $O(\sqrt{n})$。

    /**
     * 数论分块核心技巧
     * 关键观察：对于固定的n，floor(n/i) 的值只有O(sqrt(n))种不同的取值
     * 当 i 在 [l, r] 范围内时，floor(n/i) = floor(n/l) 保持不变
     * 其中 r = floor(n / floor(n/l))
     */

    /**
     * 计算 sum_{i=1}^{n} floor(n / i)
     * 普通做法: O(n)
     * 数论分块: O(sqrt(n))
     */
    public long sumFloorDivision(long n) {
        long sum = 0;
        long l = 1;
        while (l <= n) {
            long q = n / l;           // 当前块内的商
            long r = n / q;           // 当前块的右端点
            // 块 [l, r] 内 floor(n/i) 都等于 q
            sum += q * (r - l + 1);
            l = r + 1;
        }
        return sum;
    }

    /**
     * 莫比乌斯反演 + 数论分块
     * 计算 sum_{i=1}^{n} f(floor(n/i))，其中f可通过莫比乌斯反演得到
     * 
     * 经典应用：计算 sum_{i=1}^{n} sum_{j=1}^{m} gcd(i,j)
     * 转化为：sum_{d=1}^{min(n,m)} φ(d) * floor(n/d) * floor(m/d)
     */
    public long sumGcdPairs(int n, int m) {
        long sum = 0;
        int limit = Math.min(n, m);
        int l = 1;

        while (l <= limit) {
            int q1 = n / l;
            int q2 = m / l;
            // 找到最大的r使得 floor(n/r)=q1 且 floor(m/r)=q2
            int r = Math.min(n / q1, m / q2);
            r = Math.min(r, limit);

            // 计算 φ 的前缀和 [l, r]
            long phiPrefixSum = prefixPhi(r) - prefixPhi(l - 1);
            sum += phiPrefixSum * q1 * q2;

            l = r + 1;
        }

        return sum;
    }

    // 欧拉函数前缀和（可预处理）
    private long[] phiPrefix;

    public void precomputePhiPrefix() {
        phiPrefix = new long[maxN + 1];
        for (int i = 1; i <= maxN; i++) {
            phiPrefix[i] = phiPrefix[i - 1] + phi[i];
        }
    }

    private long prefixPhi(int n) {
        if (phiPrefix == null) precomputePhiPrefix();
        return phiPrefix[n];
    }

完整应用：互质对计数

    /**
     * 计算1 <= i <= n, 1 <= j <= m 中满足 gcd(i,j) = 1 的数对个数
     * 利用莫比乌斯反演：
     * 令 F(d) = 满足 d|gcd(i,j) 的数对个数 = floor(n/d) * floor(m/d)
     * 令 f(1) = 满足 gcd(i,j) = 1 的数对个数
     * 由 F(d) = sum_{d|k} f(k)，反演得：
     * f(1) = sum_{d=1}^{min(n,m)} μ(d) * floor(n/d) * floor(m/d)
     */
    public long countCoprimePairs(int n, int m) {
        long sum = 0;
        int limit = Math.min(n, m);
        int l = 1;

        while (l <= limit) {
            int q1 = n / l;
            int q2 = m / l;
            int r = Math.min(n / q1, m / q2);
            r = Math.min(r, limit);

            // μ的前缀和
            long muSum = prefixMu(r) - prefixMu(l - 1);
            sum += muSum * q1 * q2;

            l = r + 1;
        }

        return sum;
    }

    private long[] muPrefix;

    public void precomputeMuPrefix() {
        muPrefix = new long[maxN + 1];
        for (int i = 1; i <= maxN; i++) {
            muPrefix[i] = muPrefix[i - 1] + mu[i];
        }
    }

    private long prefixMu(int n) {
        if (muPrefix == null) precomputeMuPrefix();
        return muPrefix[n];
    }

杜教筛：亚线性时间求前缀和

当 $n$ 达到 $10^{11}$ 级别时，线性筛预处理不可行。杜教筛利用递归+哈希表可在 $O(n^{2/3})$ 时间内求出积性函数前缀和。

import java.util.HashMap;
import java.util.Map;

/**
 * 杜教筛实现：计算莫比乌斯函数和欧拉函数的前缀和
 */
public class DuJiaoSieve {

    private NumberTheory nt;
    private long n;
    private Map<Long, Long> muSumCache;   // S_μ(n) 缓存
    private Map<Long, Long> phiSumCache;  // S_φ(n) 缓存
    private int precomputeLimit;          // 线性筛预处理的上限

    public DuJiaoSieve(long n) {
        this.n = n;
        // 预处理前 n^(2/3) 个数
        this.precomputeLimit = (int) Math.pow(n, 2.0 / 3.0);
        this.nt = new NumberTheory(precomputeLimit);
        nt.precomputeMuPrefix();
        nt.precomputePhiPrefix();

        muSumCache = new HashMap<>();
        phiSumCache = new HashMap<>();

        // 基础情况
        muSumCache.put(0L, 0L);
        muSumCache.put(1L, 1L);
        phiSumCache.put(0L, 0L);
        phiSumCache.put(1L, 1L);
    }

    /**
     * 杜教筛计算莫比乌斯函数前缀和 S_μ(n) = sum_{i=1}^{n} μ(i)
     * 利用公式: sum_{i=1}^{n} sum_{d|i} μ(d) = 1 = sum_{d=1}^{n} μ(d) * floor(n/d)
     * 即: S_μ(n) = 1 - sum_{i=2}^{n} S_μ(floor(n/i)) * (某个系数)
     * 
     * 更简洁的推导：
     * sum_{i=1}^{n} [i==1] = sum_{i=1}^{n} sum_{d|i} μ(d)
     * 1 = sum_{d=1}^{n} μ(d) * floor(n/d)
     * S_μ(n) = 1 - sum_{i=2}^{n} S_μ(floor(n/i)) (通过整除分块)
     */
    public long muPrefixSum(long x) {
        if (x <= precomputeLimit) {
            return nt.prefixMu((int) x);
        }
        if (muSumCache.containsKey(x)) {
            return muSumCache.get(x);
        }

        long result = 1; // sum_{i=1}^{n} [i==1] 的实质就是1
        long l = 2;
        while (l <= x) {
            long q = x / l;
            long r = x / q;
            result -= (r - l + 1) * muPrefixSum(q);
            l = r + 1;
        }

        muSumCache.put(x, result);
        return result;
    }

    /**
     * 杜教筛计算欧拉函数前缀和 S_φ(n) = sum_{i=1}^{n} φ(i)
     * 利用公式: sum_{d|n} φ(d) = n
     * sum_{i=1}^{n} i = sum_{i=1}^{n} sum_{d|i} φ(d) = sum_{d=1}^{n} φ(d) * floor(n/d)
     * n(n+1)/2 = sum_{d=1}^{n} φ(d) * floor(n/d)
     * S_φ(n) = n(n+1)/2 - sum_{i=2}^{n} S_φ(floor(n/i))
     */
    public long phiPrefixSum(long x) {
        if (x <= precomputeLimit) {
            return nt.prefixPhi((int) x);
        }
        if (phiSumCache.containsKey(x)) {
            return phiSumCache.get(x);
        }

        long result = x * (x + 1) / 2;
        long l = 2;
        while (l <= x) {
            long q = x / l;
            long r = x / q;
            result -= (r - l + 1) * phiPrefixSum(q);
            l = r + 1;
        }

        phiSumCache.put(x, result);
        return result;
    }
}

完整测试程序

/**
 * 莫比乌斯反演完整测试
 */
public class MobiusDemo {

    public static void main(String[] args) {
        System.out.println("=== 线性筛预处理 ===\n");
        NumberTheory nt = new NumberTheory(100);
        nt.printMuTable(30);

        System.out.println("\n=== 莫比乌斯反演测试 ===");
        // 构造 F(n) = sum_{d|n} d（即因子和函数 σ(n)）
        long[] F = new long[101];
        for (int i = 1; i <= 100; i++) {
            for (int j = i; j <= 100; j += i) {
                F[j] += i;
            }
        }
        // 反演求 f(n) = n（因为 σ = id * 1，反演得 id = σ * μ）
        long[] f = nt.mobiusInversion(F);
        System.out.println("验证: 反演后 f(n) 应等于 n");
        for (int i = 1; i <= 20; i++) {
            System.out.printf("f(%d) = %d (期望: %d)%n", i, f[i], i);
        }

        System.out.println("\n=== 数论分块测试 ===");
        long sumDiv = nt.sumFloorDivision(100);
        System.out.println("sum_{i=1}^{100} floor(100/i) = " + sumDiv);
        // 手动验证: 100+50+33+25+20+16+14+12+11+10+9+8+7+7+6+6+5+5+5+5+4+4+4+4+4+3+3+3+3+3+3+3+3+2+2+2+2+2+2+2+2+2+2+2+2+2+2+2+2+2

        System.out.println("\n=== 互质对计数 ===");
        nt.precomputeMuPrefix();
        long coprime = nt.countCoprimePairs(100, 100);
        System.out.println("1<=i,j<=100 中 gcd(i,j)=1 的数对个数: " + coprime);

        System.out.println("\n=== GCD之和 ===");
        nt.precomputePhiPrefix();
        long gcdSum = nt.sumGcdPairs(100, 100);
        System.out.println("sum_{i=1}^{100} sum_{j=1}^{100} gcd(i,j) = " + gcdSum);

        System.out.println("\n=== 杜教筛大数测试 ===");
        DuJiaoSieve djs = new DuJiaoSieve(1_000_000_000L);
        long muSum = djs.muPrefixSum(1_000_000_000L);
        long phiSum = djs.phiPrefixSum(1_000_000_000L);
        System.out.println("S_μ(10^9) = " + muSum);
        System.out.println("S_φ(10^9) = " + phiSum);
    }
}

复杂度分析

算法 预处理 单次查询 空间
线性筛 $O(n)$ $O(1)$ $O(n)$
数论分块 无需 $O(\sqrt{n})$ $O(1)$
杜教筛 $O(n^{2/3})$ $O(n^{2/3})$ $O(n^{2/3})$
  • 莫比乌斯反演（因子枚举）：$O(\sqrt{n})$ 计算单个值
  • 莫比乌斯反演（倍数形式）：$O(n \log n)$ 计算全部值
  • 互质对计数 + 数论分块：$O(\sqrt{n})$ 单次查询
  • 杜教筛：$O(n^{2/3})$ 处理 $10^{11}$ 级别数据

总结

本文通过Java实现了莫比乌斯反演的完整计算框架：

  1. 线性筛：$O(n)$ 同时预处理素数、莫比乌斯函数、欧拉函数
  2. 莫比乌斯反演：建立 $F$ 与 $f$ 之间的双向转换
  3. 数论分块：将 $O(n)$ 求和优化到 $O(\sqrt{n})$
  4. 杜教筛：$O(n^{2/3})$ 处理超大规模数据的前缀和

读者可以进一步扩展：实现Min_25筛处理更复杂的积性函数、应用莫比乌斯反演解决GCD卷积问题、或构建交互式数论计算工具。

Tagged with: