数论函数学习笔记

数论函数学习笔记

由多个博客学习组合而成

狄利克雷卷积

所谓的狄利克雷卷积，是定义在数论函数$({\mathbb{Z}}^{+}\to C的函数)$间的一种二元运算，可这样定义：

$(f\ast g)(n):={\sum }_{xy=n}f(x)g(y)$

也可以写成：

$(f\ast g)(n):={\sum }_{d|n}f(d)g(\frac{n}{d})$

为了方便讨论先定义一些常用的数论函数符号：

• 单位函数$\epsilon (n):=\{\begin{array}{l}1,whilen=1\\ 0,otherwise\end{array}$

• 幂函数$I{d}_k(n):=n^k$.当$k=1$时为恒等函数$Id(n)$，当$k=0$时为常数函数$1(n)$

• 除数函数${\sigma }_k(n):={\sum }_{d|n}{d}^k$，当$k=1$是为因数和函数$\sigma (n)$，当$k=0$时为因数个数函数${\sigma }_0(n)$

• 欧拉函数$\phi (n)$

值得注意的是，这些函数都是所谓的积性函数，即满足$f(1)=1$，且若a，b互质，则有$f(a)f(b)=f(ab)$.

其中前两者还是完全积性函数，它们不需要满足a，b互质的条件也符合等式

积性函数之间的狄利克雷卷积还有特殊的性质：

若$f,g$为积性函数，那么$f\ast g$也是积性函数

首先$(f\ast g)(1)=f(1)\ast g(1)=1$，设$a,b$互质则有：

$$\begin{gathered} (f\ast g)(a)\cdot (f\ast g)(b)=\sum _{d1|a}f(d1)g(\frac{a}{d1})\cdot \sum _{d2|b}f(d2)g(\frac{b}{d2}) \\ =\sum _{d1|a,d2|b}f(d1)f(d2)g(\frac{a}{d1})g(\frac{b}{d2}) \\ =\sum _{d1|a,d2|b}f(d1d2)g(\frac{ab}{d1d2}) \end{gathered}$$

因为$a,b$互质，所以就有了$(f\ast g)(a)\cdot (f\ast g)(b)=(f\ast g)(ab)$该等式成立

再考虑某些数论函数之间的关系：

除数函数与幂函数

首先根据定义我们有$(f\ast 1)(n)={\sum }_{d|n}f(d)\ast 1(\frac{n}{d})$

所以可以有$(I{d}_k\ast 1)(n)={\sum }_{d|n}I{d}_k(d)1(\frac{n}{d})={\sum }_{d|n}{d}^k={\sigma }_k(n)$

简写成$(I{d}_k\ast 1)(x)={\sigma }_k(x)$

欧拉函数与恒等函数

$(\phi \ast 1)(n)={\sum }_{d|n}\phi (d)$

当$n=p^m$时，有：

${\sum }_{d|n}\phi (d)=\phi (1)+{\sum }_{i=1}^m\phi (p^i)=1+{\sum }_{i=1}^m(p^i-p^{i-1})=p^m=n$

故$(\phi \ast 1)(p^m)=p^m$

现在令$n$为任意正整数，那么$n=\prod p^m$ 又因为$(\phi \ast 1)(n)$为积性函数，所以$(\phi \ast 1)(n)=\prod p^m=n$

所以有$(\phi \ast 1)=Id$

下面对于一些卷积的法则进行证明：

交换律：
$$\begin{gathered} (f\ast g)(n)=\sum _{xy=n}f(x)g(y) \\ =\sum _{xy=n}f(y)g(x) \\ =(g\ast f)(n) \end{gathered}$$

结合律：
$$\begin{gathered} (f\ast (g\ast h))(n)=\sum _{xy=n}f(x)\ast (g\ast h)(y) \\ =\sum _{xy=n}f(x)\ast \sum _{uv=y}g(u)\ast h(v) \\ =\sum _{xuv=n}f(x)g(u)h(v) \\ =\sum _{yv=n}h(v)\ast \sum _{xu=y}f(x)g(u) \\ =((f\ast g)\ast h)(n) \end{gathered}$$

分配律：
$$\begin{gathered} (f\ast (g+h))(n)=\sum _{xy=n}f(x)\ast (g+h)(y) \\ =\sum _{xy=n}f(x)g(y)+f(x)h(y) \\ =\sum _{xy=n}f(x)g(y)+\sum _{xy=n}f(x)h(y) \\ =(f\ast g)(n)+(f\ast h)(n) \end{gathered}$$

单位元的用法：
$$(f\ast \epsilon )(n)=\sum _{d|n}f(d)\epsilon (\frac{n}{d})=f(n)$$

逆元：

若$f\ast g=\epsilon$时我们把$g$称为$f$的狄利克雷逆元并写作$f^{-1}$
$$(f\ast f^{-1})=\epsilon$$

莫比乌斯反演

我们定义数论函数$1$的狄利克雷逆元为$\mu$
得到$\mu (x)=\{\begin{array}{l}1n=1\\ (-1{)}^{m}n=p_1{p}_2...p_m\\ 0otherwise\end{array}$

所以有$(f\ast 1)=g$可推出$f=(g\ast \mu )$

下面是几道经典题目的推导

P3455 [POI2007]ZAP-Queries

简要题意

求${\sum }_{x=1}^a{\sum }_{y=1}^b[gcd(x,y)=d]$

推导

$$\begin{gathered} \sum _{x=1}^a\sum _{y=1}^b[gcd(x,y)=d] \\ =\sum _{x=1}^{\lfloor \frac{a}{d}\rfloor }\sum _{y=1}^{\lfloor \frac{b}{d}\rfloor }[gcd(x,y)=1] \\ =\sum _{x=1}^{\lfloor \frac{a}{d}\rfloor }\sum _{y=1}^{\lfloor \frac{b}{d}\rfloor }\epsilon (gcd(x,y)) \\ =\sum _{x=1}^{\lfloor \frac{a}{d}\rfloor }\sum _{y=1}^{\lfloor \frac{b}{d}\rfloor }\sum _{k|gcd(x,y)}\mu (k) \\ =\sum _{k=1}^n\mu (k)\sum _{x=1}^{\lfloor \frac{a}{dk}\rfloor }\sum _{y=1}^{\lfloor \frac{b}{dk}\rfloor }1 \\ =\sum _{k=1}^n\mu (k)\lfloor \frac{a}{dk}\rfloor \lfloor \frac{b}{dk}\rfloor \end{gathered}$$

其中$n=max(a,b)$

推导到最后一步就可以使用整除分块加速

P2522 [HAOI2011]Problem b

简要题意

求${\sum }_{i=x}^a{\sum }_{j=y}^b[gcd(i,j)=d]$

推导

这道题需要想到容斥原理。

将这个式子写成
$$\sum _{i=1}^a\sum _{j=1}^b[gcd(i,j)=d]-\sum _{i=1}^a\sum _{j=1}^{y-1}[gcd(i,j)=d]-\sum _{i=1}^{x-1}\sum _{j=1}^b[gcd(i,j)=d]+\sum _{i=1}^{x-1}\sum _{j=1}^{y-1}[gcd(i,j)=d]$$
就跟上面那道一样了

P2257 YY的GCD

简要题意

求${\sum }_{x=1}^n{\sum }_{y=1}^m[gcd(x,y)=p]$其中$p$为任意质数

推导

首先这里筛出质数后称作$prime$数量为$tot$
$$\begin{gathered} \sum _{x=1}^n\sum _{y=1}^m[gcd(x,y)=p] \\ =\sum _{i=1}^{tot}\sum _{x=1}^n\sum _{y=1}^m[gcd(x,y)=prime[i]] \\ =\sum _{i=1}^{tot}\sum _{x=1}^{\lfloor \frac{n}{prime[i]}\rfloor }\sum _{y=1}^{\lfloor \frac{m}{prime[i]}\rfloor }[gcd(x,y)=1] \\ =\sum _{i=1}^{tot}\sum _{d=1}^{max(n,m)}\mu (d)\lfloor \frac{n}{prime[i]\ast d}\rfloor \lfloor \frac{m}{prime[i]\ast d}\rfloor \\ 这里T=prime[i]\ast d可以将式子转化为 \\ =\sum _{T=1}^{max(n,m)}\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{k|T,k\in prime}\mu (\frac{T}{k}) \end{gathered}$$
发现后面一个式子${\sum }_{k|T,k\in prime}\mu (\frac{T}{k})$是可以预处理出来的，这样就减少了复杂度

P1829 [国家集训队]Crash的数字表格 / JZPTAB

简要题意

求${\sum }_{i=1}^n{\sum }_{j=1}^{m}lcm(i,j)$

推导

$$\begin{gathered} \sum _{i=1}^n\sum _{j=1}^{m}lcm(i,j) \\ =\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{gcd(i,j)} \\ =\sum _{d=1}^{max(n,m)}\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{d}[gcd(i,j)=d] \\ =\sum _{d=1}^{max(n,m)}\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ijd[gcd(i,j)=1] \\ =\sum _{d=1}^{max(n,m)}d\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ij[gcd(i,j)=1] \\ 设sum(n,m)=\sum _{i=1}^n\sum _{j=1}^{m}ij[gcd(i,j)=1] \\ 化简这个式子 \\ \sum _{i=1}^n\sum _{j=1}^{m}ij[gcd(i,j)=1] \\ =\sum _{d=1}^{max(n,m)}\mu (d)d^2\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ij \\ 设f(n,m)=\sum _{i=1}^n\sum _{j=1}^{m}ij \\ 化简这个式子有 \\ \sum _{i=1}^n\sum _{j=1}^{m}ij \\ =\sum _{i=1}^{n}i\sum _{j=1}^{m}j \\ =\frac{(n+1)n}{2}\cdot \frac{(m+1)m}{2} \\ 所以sum(n,m)=\sum _{d=1}^{max(n,m)}\mu (d)\cdot d^2\cdot f(\frac{n}{d},\frac{m}{d}) \\ 原式=\sum _{d=1}^{max(n,m)}d\cdot sum(\frac{n}{d},\frac{m}{d}) \end{gathered}$$

这里的两个式子都可以用整除分块+前缀和的方式加速

代码

#include <bits/stdc++.h>

const int N = 1e7 + 5, mod = 20101009;

inline int add(int a, int b) { return (a += b) >= mod ? a - mod : a; }
inline int sub(int a, int b) { return (a -= b) < 0 ? a + mod : a; }
inline int mul(int a, int b) { return 1LL * a * b % mod; }

std::vector < int > prime;
int v[N], mu[N], sum[N];
void input(int n) {
    mu[1] = 1;
    for(int i = 2; i <= n; ++i) {
        if(!v[i]) prime.push_back(i), v[i] = i, mu[i] = mod - 1;
        for(auto x : prime) {
            if(x > n / i) break;
            v[i * x] = x;
            if(!(i % x)) {
                mu[x * i] = 0;
                break;
            }
            mu[i * x] = sub(mod, mu[i]);
        }
    }
    for(int i = 1; i <= n; ++i) sum[i] = add(sum[i - 1], mul(mu[i], mul(i, i)));
}

inline int sumij(int i, int j) {
    return mul((1LL * (i + 1) * i / 2) % mod, (1LL * (j + 1) * j / 2) % mod);
}

inline int f(int n, int m) {
    int re = 0;
    for(int l = 1, r; l <= m; l = r + 1) {
        r = std::min(n / (n / l), m / (m / l));
        int cnt1 = n / l, cnt2 = m / l;
        re = add(re, mul(sumij(cnt1, cnt2), sub(sum[r], sum[l - 1])));
    }
    return re;
}

int main() {
    std::ios::sync_with_stdio(false);
    std::cin.tie(nullptr);
    int n, m;
    std::cin >> n >> m;
    if(n < m) std::swap(n, m);
    input(n);
    int ans = 0;
    for(int d = 1; d <= n; ++d) {
        ans = add(ans, mul(d, f(n / d, m / d)));
    }
    std::cout << ans << '\n';
    return 0;
}

杜教筛

对于任意一个数论函数$g$，必满足：
$$\sum _{i=1}^n(f\ast g)(i)=\sum _{i=1}^n\sum _{d|i}g(d)f(\frac{i}{d})=\sum _{i=1}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$$

重要推论：
$$\begin{gathered} g(1)S(n)=\sum _{i=1}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )-\sum _{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor ) \\ =\sum _{i=1}^n(f\ast g)(i)-\sum _{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor ) \end{gathered}$$

如果可以构造一个数论函数$g$使得：

• 可以快速计算${\sum }_{i=1}^n(f\ast g)(i)$
• 可以快速计算$g$的前缀和，这样就能用数论分块求解${\sum }_{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$
  则我们可以取$m=\frac{2}{3}n$的复杂度求得$g(1)S(n)$

P3768 简单的数学题

$$\begin{gathered} \sum _{i=1}^n\sum _{j=1}^{n}ijgcd(i,j) \\ =\sum _{i=1}^n\sum _{j=1}^{n}ij\sum _{k|i,k|j}\phi (k) \\ =\sum _{k=1}^n\phi (k)\sum _{k|i}\sum _{k|i}ij \\ =\sum _{k=1}^n\phi (k)\times k^2\times (\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }i{)}^2 \\ =\sum _{k=1}^n\phi (k)\times k^2\times \sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }{i}^3 \end{gathered}$$

对$\phi (k)k^2$做杜教筛

把$k$看成$Id$函数有，$f=I{d}^2\phi$，卷一个$g=I{d}^2$可以得到$(f\ast g)(n)=n^3$

这里有一个结论就是$I{d}^k$卷一个$I{d}^k$可以看成是$(I{d}^k\ast I{d}^k)(n)=n^k(n-\phi (n))$

Powerful Number筛

Powerful Number（以下叫做PN）的定义：设$n={\prod }_{i=1}^m{p}_i^{e_i}$只有所有的$e_i>1$才能是PN

可以证明PN的数量级是$\Theta (\sqrt{n})$的，也就是可以暴力搜索所有小于$\sqrt{n}$的质数处理出PN

Powerful Number筛其实是杜教筛的升级版本，计算的也是$F(n)={\sum }_{i=1}^{n}f(i)$（$f(i)$是积性函数），但是依然不适用于所有场景。

PN筛的步骤：

• 找到一个积性函数$g(n)$在$n$为质数的时候$g(p)=f(p)$，并且$G(n)={\sum }_{i=1}^{n}g(i)$
• 找到一个快速计算$G$也就是$g$的前缀和的方法
• 令$f=g\ast h$，此时有$h$为积性函数且$h(1)=1$
• 最后$F(n)={\sum }_{d=1,disPN}^{n}h(d)G(\lfloor \frac{n}{d}\rfloor )$

下面是对于这条公式的证明（可跳过）：

首先对于$h$函数有一个特殊性质：

$$\begin{gathered} f(p)=g(p)h(1)+g(1)h(p)=g(p) \\ h(p)=0 \end{gathered}$$
因为$h$是一个积性函数，我们知道$h(notPN)=0$

于是
$$\begin{gathered} F(n)=\sum _{i=1}^n\sum _{d|i}h(d)\ast g(\lfloor \frac{i}{d}\rfloor ) \\ =\sum _{d=1}^{n}h(d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }g(i) \\ =\sum _{d=1}^{n}h(d)G(\lfloor \frac{n}{d}\rfloor ) \\ =\sum _{d=1,disPN}^{n}h(d)G(\lfloor \frac{n}{d}\rfloor ) \end{gathered}$$
得证

P5325【模板】Min_25 筛

给定$f(p^k)=p^k(p^k-1)$，求${\sum }_{i=1}^{n}f(i)$
$$\begin{gathered} 定义g(n)=Id(n)\phi (n) \\ 此时G(n)可以通过杜教筛快速得出 \\ G(n)=\sum _{i=1}^n{i}^2-\sum _{d=2}^{n}d\ast G(\lfloor \frac{n}{d}\rfloor ) \\ 那么我们知道h(p^k)可以枚举计算（法1） \\ 另外也有h(p^k)直接推出的公式 \\ f(p^k)=\sum _{i=0}^{k}g(p^{k-i})h(p^i) \\ {p}^k(p^k-1)=\sum _{i=0}^k{p}^{k-i}\phi (p^i)h(p^i) \\ {p}^k(p^k-1)=\sum _{i=0}^k{p}^{2k-2i-1}(p-1)h(p^i) \\ h(p^k)=p^k(p^k-1)-\sum _{i=0}^{k-1}{p}^{2k-2i-1}h(p^i) \\ h(p^k)-p^{2}h(p^{k-1})=p^k(p^k-1)-p^{k+1}(p^{k-1}-1)-p(p-1)h(p^{k-1}) \\ h(p^k)-ph(p^{k-1})=p^{k+1}-p^k \\ \frac{h(p^k)}{p^k}-\frac{h(p^{k-1})}{p^{k-1}}=p-1 \\ 有h(p)=0推出 \\ h(p^k)=(k-1)(p-1)p^k（法二） \\ 时间复杂度相似但是法二比较好写 \end{gathered}$$

void solve() {
    i64 n;
    std::cin >> n;
    int N1 = 10000000, N2 = sqrt(n);
    std::vector < int > prime, phi(N1 + 1), v(N1 + 1);
    std::vector < MInt > g(N1 + 1), s(N1 + 1);
    phi[1] = 1;
    for(int i = 2; i <= N1; ++i) {
        if(!v[i]) prime.emplace_back(i), v[i] = i, phi[i] = i - 1;
        for(auto x : prime) {
            if(x > N1 / i || x > v[i]) break;
            v[i * x] = x;
            phi[i * x] = phi[i] * (i % x == 0 ? x : x - 1);
        }
    }
    for(int i = 1; i <= N1; ++i) g[i] = 1LL * i * phi[i], s[i] = s[i - 1] + g[i];
    MInt inv6 = MInt(6).inverse(), inv2 = 500000004;
    std::unordered_map < i64 , MInt > mp;
    auto G = [&](auto self, i64 n) -> MInt {
        if(n <= N1) return s[n];
        if(mp.count(n)) return mp[n];
        MInt re = MInt(n) * MInt(n + 1) * MInt(2 * n + 1) * inv6;
        for(i64 l = 2, r; l <= n; l = r + 1) {
            r = std::min(n, n / (n / l));
            re -= (MInt(r) * MInt(r + 1) * inv2 - MInt(l - 1) * MInt(l) * inv2) * self(self, n / l);
        }
        return mp[n] = re;
    } ;
    // std::cout << g[2] << ' ' << g[3] << ' ' << g[6] << ' ' << g[9] << '\n';
    MInt ans = 0;
    auto dfs = [&](auto self, int id, i64 now, MInt re) -> void {
        i64 cnt = 1LL * prime[id] * prime[id];
        if(now > n / cnt) {
            // std::cout << now << ' ' << prime[id] << ' ' << n << '\n';
            ans += re * G(G, n / now);
            return ;
        }
        self(self, id + 1, now, re);
        MInt st = cnt * g[prime[id]];
        i64 las = n / now;
        while(cnt <= las) {
            MInt to = (MInt(cnt) * MInt(cnt - 1) - st);
            self(self, id + 1, now * cnt, re * to);
            if(cnt > las / prime[id]) break;
            st *= prime[id] * prime[id];
            st += to * g[prime[id]];
            cnt *= prime[id];
        }
    } ;
    dfs(dfs, 0, 1, 1);
    std::cout << ans << '\n';
}

void solve() {
    i64 n;
    std::cin >> n;
    int N1 = 10000000, N2 = sqrt(n);
    std::vector < int > prime, phi(N1 + 1), v(N1 + 1);
    std::vector < MInt > g(N1 + 1), s(N1 + 1);
    phi[1] = 1;
    for(int i = 2; i <= N1; ++i) {
        if(!v[i]) prime.emplace_back(i), v[i] = i, phi[i] = i - 1;
        for(auto x : prime) {
            if(x > N1 / i || x > v[i]) break;
            v[i * x] = x;
            phi[i * x] = phi[i] * (i % x == 0 ? x : x - 1);
        }
    }
    for(int i = 1; i <= N1; ++i) g[i] = 1LL * i * phi[i], s[i] = s[i - 1] + g[i];
    MInt inv6 = MInt(6).inverse(), inv2 = 500000004;
    std::unordered_map < i64 , MInt > mp;
    auto G = [&](auto self, i64 n) -> MInt {
        if(n <= N1) return s[n];
        if(mp.count(n)) return mp[n];
        MInt re = MInt(n) * MInt(n + 1) * MInt(2 * n + 1) * inv6;
        for(i64 l = 2, r; l <= n; l = r + 1) {
            r = std::min(n, n / (n / l));
            re -= (MInt(r) * MInt(r + 1) * inv2 - MInt(l - 1) * MInt(l) * inv2) * self(self, n / l);
        }
        return mp[n] = re;
    } ;
    // std::cout << g[2] << ' ' << g[3] << ' ' << g[6] << ' ' << g[9] << '\n';
    MInt ans = 0;
    auto dfs = [&](auto self, int id, i64 now, MInt re) -> void {
        i64 cnt = 1LL * prime[id] * prime[id];
        if(now > n / cnt) {
            // std::cout << now << ' ' << prime[id] << ' ' << n << '\n';
            ans += re * G(G, n / now);
            return ;
        }
        self(self, id + 1, now, re);
        MInt i = prime[id] - 1;
        i64 las = n / now;
        while(cnt <= las) {
            self(self, id + 1, now * cnt, re * i * MInt(cnt));
            if(cnt > las / prime[id]) break;
            i += MInt(prime[id] - 1);
            cnt *= prime[id];
        }
    } ;
    dfs(dfs, 0, 1, 1);
    std::cout << ans << '\n';
}

Min_25筛

使用要求

在$\Theta (\frac{n^{\frac{3}{4}}}{logn})$的时间复杂度下解决积性函数的前缀和问题。

要求$f(p)$是一个关于$p$的项数较少的多项式或可以快速求值，$f(p^c)$可以快速求值。

定义说明

$p$没有特殊说明时代表全体质数

$isprime(n)$表示$n$为质数时为$1$

$p_k$表示全体质数中第$k$小的质数，特别的$p_0=1$

l p f ( n ) : = [ 1 < n ] m i n { p : p ∣ n } + [ 1 = n ] lpf(n):=[1<n]min\{p:p|n\}+[1=n] lpf(n):=[1<n]min{p:p∣n}+[1=n]表示$n$的最小质因数。特别的，$lpf(1)=1$

$F_{prime}(n):={\sum }_{2\le p\le n}f(p)$

$F_k(n):={\sum }_{i=2}^n[p_k\le lpf(i)]f(i)$

使用方法

观察$F_k(n)$的定义，我们要求${\sum }_{i=1}^{n}f(i)=F_1(n)+f(1)=F_1(n)+1$

考虑如何求出$F_k(n)$
$$\begin{gathered} F_k(n)=\sum _{i=2}^n[p_k\le lpf(i)]f(i) \\ =\sum _{k\le i,p_i^2\le n}\sum _{1\le c,p_i^c\le n}f(p^c)([c>1]+F_{i+1}(\lfloor \frac{n}{p_i^c}\rfloor ))+\sum _{k\le i,p_i\le n}f(p_i) \\ =\sum _{k\le i,p_i^2\le n}(\sum _{1\le c,p_i^c\le n}f(p^c)([c>1]+F_{i+1}(\lfloor \frac{n}{p_i^c}\rfloor ))+F_{prime}(n)-F_{prime}(p_{k-1})) \\ =\sum _{k\le i,p_i^2\le n}(\sum _{1\le c,p_i^{c+1}\le n}(f(p^c)F_{i+1}(\lfloor \frac{n}{p_i^c}\rfloor )+f(p_i^{c+1}))+F_{prime}(n)-F_{prime}(p_{k-1})) \end{gathered}$$

最后一步推导的说明

对于满足$p_i^c\le n<p_i^{c+1}$的$c$，有$\lfloor \frac{n}{p_i^c}\rfloor <p_i<p_{i+1}$。

故$F_{i+1}(\lfloor \frac{n}{p_i^c}\rfloor )=0$

有$F_k(n)=0(p_k>n)$。

假设现在已经求出了所有的$F_{prime}(n)$，那么有两种方式可以求出所有的$F_k(n)$：

• 直接按照递推式计算
• 从大到小枚举$p$转移，仅当$p^2<n$时转移增加值不为零，故按照递推式后缀和优化即可。

计算$F_{prime}(n)$

一般情况下$f(p)$是一个关于$p$的低次多项式，可以表示成$f(p)=\sum a_i{p}^{c_i}$

那么对于每个$p^{c_i}$，其对$F_{prime}(n)$的贡献可以表示为$a_i{\sum }_{2\le p\le n}{p}^{c_i}$

分开考虑每个$p^{c_i}$的贡献，问题就转变为了：给定$n,s,g(p)=p^s$，对所有的$m=\lfloor \frac{n}{i}\rfloor$，求${\sum }_{p\le m}g(p)$。

设：
$$G_x(n)=\sum _{2\le i\le n}[isprime(i)or{p}_x<lpf(i)]g(i)$$
即埃氏筛第$k$轮晒完后剩下的数的$g$值之和。

对于一个合数$x$，必然有$lpf(x)\le \sqrt{x}$，则$F_{prime}=G_{\sqrt{n}}$故只需筛到$G_{\sqrt{n}}$即可

考虑$G$的边界值，$G_0(n)={\sum }_{i=2}^{n}g(i)$

对于转移，考虑埃氏筛的过程，分开讨论每部分的贡献，有：

• 对于$n<p_k^2$的部分，$G$值不变，即$G_x(n)=G_{x-1}(n)$
• 对于$p_k^2\le n$部分，被筛掉的数必有质因子$p_k$即$-g(p_k)G_{k-1}(\lfloor \frac{n}{p_k}\rfloor )$
• 对于第二个部分，会有一些$lpf(i)<p_k$的$i$被减去，加回来的贡献是$g(p_k)G_{k-1}(p_{k-1})$

则可以得到$G_k(n)=G_{k-1}(n)-[p_k^2\le n]g(p_k)(G_{k-1}(\lfloor \frac{n}{p_k}\rfloor )-G_{k-1}(p_{k-1}))$

总结

对于$F_k(n)$的计算，利用第一种方式写一个函数，这个实现难度较低，已经适用于大部分题目。

对于$F_{prime}的计算$，直接按递推式实现即可。

对于$p_k^2\le n$可以用线性筛预处理出$s_k:=F_{prime}(p_k)$来代替$F_k$递推式中的$F_{prime}(p_{k-1})$

相应的$G_{k-1}(p_{k-1})={\sum }_{i=1}^{k-1}g(p_i)$也可以预处理。

用埃氏筛求$f$的前缀和时，应当明确以下几点：

• 如何快速（一般是线性复杂度）筛出前$\sqrt{n}$个$f$值
• $f(p)$的多项式表示
• 如何快速求出$f(p^c)$

明确上述几点之后按顺序实现以下几部分即可：

1. 筛出$[1,\sqrt{n}]$内的质数与前$\sqrt{n}$个$f$值；
2. 对$f(p)$多项式表示中的每一项筛出对应的$G$，合并得到$F_{prime}$的所有$O(\sqrt{n})$有用点值
3. 按照$F_k$的递推式实现递归，求出$F_1(n)$

P5325【模板】Min_25 筛

const int N = 2e5;
i64 n, lim;

void solve() {
    std::cin >> n;
    lim = sqrt(n);
    std::vector < int > prime, v(N + 10);
    std::vector < MInt > sum1(N + 10), sum2(N + 10);
    for(int i = 2; i <= N; ++i) {
        if(!v[i]) {
            prime.emplace_back(i);
            v[i] = i;
        }
        for(auto x : prime) {
            if(x > N / i || x > v[i]) break;
            v[i * x] = x;
        }
    }
    for(int i = 0; i < prime.size(); ++i) {
        sum1[i + 1] = sum1[i] + prime[i];
        sum2[i + 1] = sum2[i] + MInt(1LL * prime[i] * prime[i]);
    }
    MInt inv2 = MInt(2).inverse(), inv6 = MInt(6).inverse();
    auto f1 = [&](i64 n) -> MInt {
        return MInt(n) * MInt(n + 1) * inv2;
    } ;
    auto f2 = [&](i64 n) -> MInt {
        return MInt(n) * MInt(n + 1) * MInt(2 * n + 1) * inv6;
    } ;
    std::vector < i64 > idpre(N + 10), idsuf(N + 10), w;
    std::vector < std::vector < MInt > > G(2, std::vector < MInt > (N + 10));
    for(i64 l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        w.emplace_back(n / l);
        G[0][w.size()] = f1(n / l) - 1, G[1][w.size()] = f2(n / l) - 1;
        if(n / l <= lim) idpre[n / l] = w.size();
        else idsuf[r] = w.size();
    }
    // std::cerr << '\n';
    auto id = [&](i64 x) -> i64 {
        if(x <= lim) return idpre[x];
        else return idsuf[n / x];
    } ;
    for(int i = 0; i < prime.size(); ++i) {
        int x = prime[i];
        for(int j = 1; j <= w.size() && 1LL * x * x <= w[j - 1]; ++j) {
            G[0][j] -= x * (G[0][id(w[j - 1] / x)] - sum1[i]);
            G[1][j] -= x * x * (G[1][id(w[j - 1] / x)] - sum2[i]);
        }
    }
    // std::cerr << '\n';
    // std::cerr << "123" << '\n';
    auto f = [&](MInt pc) -> MInt {
        return pc * (pc - 1);
    } ;
    auto F = [&](auto self, int k, i64 n) -> MInt {
        if(prime[k] > n) return 0;
        MInt ans = ((G[1][id(n)] - G[0][id(n)]) - (sum2[k] - sum1[k]));
        for(int i = k; i < prime.size() && prime[i] * prime[i] <= n; ++i) {
            for(i64 c = 1, sp = prime[i]; sp <= n / prime[i]; ++c, sp *= prime[i]) {
                ans += f(MInt(sp)) * self(self, i + 1, n / sp) + f(MInt(sp * prime[i]));
            }
        }
        return ans;
    } ;
    std::cout << F(F, 0, n) + 1 << '\n';
}

2024GDCPC 另一个计数问题

题解

注意到，对于所有 $x\le \lfloor \frac{n+1}{2}\rfloor$ 的结点 $x$，$x$ 与 $2x$ 间有一条边，而 $2x$ 与 $2$ 之间有一条边，因此所有 $2\sim \lfloor \frac{n+1}{2}\rfloor$ 范围内的结点连通。而对于 $x>\lfloor \frac{n+1}{2}\rfloor$ 的结点 $x$，$x$ 不与 $2\sim \lfloor \frac{n+1}{2}\rfloor$ 中的结点连通当且仅当 $x$ 为素数。因此整个图仅包含若干个单独的素数结点与一个大的连通块。

考虑计算不连通的点对之间的贡献。不难发现，只需要求出所有 $>\lfloor \frac{n+1}{2}\rfloor$ 的素数的和与平方和。使用 Min_25 筛或分块打表均可。

2024杭电多校第2场 1004 a*b problem

题意

求出$b=x_1^{t_1}{x}_2^{t_2}{x}_3^{t_3}...x_k^{t_k},1\le b\le n$在$gcd(x_1,x_2,...,x_k)=1$时，

$(x_1,x_2,...,x_k)$有序对的数量。

题解

设$f(n)$表示$n=x_1^{t_1}{x}_2^{t_2}{x}_3^{t_3}...x_k^{t_k}$时$(x_1,x_2,...,x_k)$有序对的数量。

即答案为${\sum }_{i=1}^{n}f(i)$。

如果在这里可以联想到Min_25筛的话，我们可以发现$f(p^c)$其实是与$p$无关的，即$p$的$0$次多项式。

此时我们知道$f(p)$就是$t_1,t_2,...,t_k$中$1$的数量，那么$F_{prime}(n)$就是$1$的数量与比$n$小质数个数的乘积

而$f(p^c)$我们知道每次放进$x_i^{t_i}$中，就是放一个$t_i$的倍数在其中。其实就是关于$c$的一个完全背包，但是这个完全背包有$gcd(x_1,x_2,...,x_k)=1$的限制，这个限制其实就是需要保证每个质数放的时候一定会至少有一个$0$的质量，魔改一下完全背包，加一个是否有放$0$就行。

剩下就变成板子了。

const int N = 2e5;
i64 n, lim;
int m;
MInt f[100001][34][2];

void solve() {
    std::cin >> n >> m;
    std::vector < int > t(m), p(34);
    lim = sqrt(n);
    for(auto &x : t) std::cin >> x, p[x]++;
    if(n == 1) {
        std::cout << "1\n";
        return ;
    }
    std::vector < int > prime, v(N + 10);
    for(int i = 2; i <= N; ++i) {
        if(!v[i]) {
            prime.emplace_back(i);
            v[i] = i;
        }
        for(auto x : prime) {
            if(x > N / i || x > v[i]) break;
            v[i * x] = x;
        }
    }
    MInt inv2 = MInt(2).inverse(), inv6 = MInt(6).inverse();
    std::vector < i64 > idpre(N + 10), idsuf(N + 10), w;
    std::vector < MInt > G(N + 10);
    for(i64 l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        w.emplace_back(n / l);
        G[w.size()] = MInt(n / l) - 1;
        if(n / l <= lim) idpre[n / l] = w.size();
        else idsuf[r] = w.size();
    }
    auto id = [&](i64 x) -> i64 {
        if(x <= lim) return idpre[x];
        else return idsuf[n / x];
    } ;
    for(int i = 0; i < prime.size(); ++i) {
        int x = prime[i];
        for(int j = 1; j <= w.size() && 1LL * x * x <= w[j - 1]; ++j) {
            G[j] -= G[id(w[j - 1] / x)] - i;
        }
    }
    std::vector < MInt > Fprime(w.size() + 1);
    for(int i = 1; i <= w.size(); ++i) {
        Fprime[i] = G[i] * p[1];
    }
    for(int i = 1; i <= m; ++i) {
        for(int j = 0; j <= 33; ++j) {
            f[i][j][0] = f[i][j][1] = 0;
        }
    }
    f[0][0][0] = 1;
    for(int i = 1; i <= m; ++i) {
        for(int j = 0; j <= 33; j += t[i - 1]) {
            if(j == 0) {
                for(int k = 0; k <= 33; ++k) {
                    f[i][k][1] += f[i - 1][k][0] + f[i - 1][k][1];
                }
            } else {
                for(int k = j; k <= 33; ++k) {
                    f[i][k][0] += f[i - 1][k - j][0];
                    f[i][k][1] += f[i - 1][k - j][1];
                }
            }
        }
    }
    auto dp = [&](int c) -> MInt {
        if(c > 33) return 0;
        return f[m][c][1];
    } ;
    auto F = [&](auto self, int k, i64 n) -> MInt {
        if(prime[k] > n) return 0;
        MInt ans = Fprime[id(n)] - k * p[1];
        for(int i = k; i < prime.size() && prime[i] * prime[i] <= n; ++i) {
            for(i64 c = 1, sp = prime[i]; sp <= n / prime[i]; ++c, sp *= prime[i]) {
                ans += dp(c) * self(self, i + 1, n / sp) + dp(c + 1);
            }
        }
        return ans;
    } ;
    std::cout << F(F, 0, n) + 1 << '\n';
}