狄利克雷卷积到筛法

狄利克雷卷积到筛法

前置符号

[P]={1,if the proposition P is true0,if the proposition P is false[P] = \begin{cases} 1,&\text{if the proposition } P \text{ is true}\\ 0,&\text{if the proposition } P \text{ is false} \end{cases}

集合的表示:
S={x:P}S = \{x:P\}
表示满足条件 PPxx 构成的集合

数论函数与狄利克雷卷积

指一类函数:其定义域为 N\mathbb N^*,值域为某个数集。

两个数论函数的加法:
(f+g)(n)=f(n)+g(n)(f+g)(n) = f(n) + g(n)
逐项相加即可。

数乘:
(cf)(n)=c(f(n))(cf)(n) = c(f(n))
两个数论函数的狄利克雷卷积(用 * 表示):
(fg)(n)=dnf(d)×g(nd)=d1×d2=nf(d1)×f(d2)\begin{aligned} (f*g)(n) &= \sum_{d\mid n}f(d)\times g\left(\frac n d\right)\\ &= \sum_{d_1 \times d_2 = n}f(d_1)\times f(d_2) \end{aligned}
注意狄利克雷卷积其实就是关于乘法的卷积(d1d2=nd_1d_2=n)。

其具有如下性质:

  1. 交换律:fg=gff*g = g*f,乘法满足交换律,得证。

  2. 结合律:(fg)h=f(gh)(f*g)*h=f*(g*h)
    (ij)k=n(f(i)g(j))h(k)=i(jk)=nf(i)(g(j)h(k))\sum_{(ij)k = n}(f(i)g(j))h(k)=\sum_{i(jk) = n}f(i)(g(j)h(k))

  3. 对加法的分配律:(f+g)h=fh+gh(f+g)*h = f*h + g*h
    (f(n)+g(n))h(n)=dn(f(d)+d(d))h(nd)=dnf(d)h(nd)+dng(d)h(nd)=f(n)h(n)+g(n)h(n)\begin{aligned} &(f(n) + g(n)) * h(n)\\ =&\sum_{d\mid n}(f(d) + d(d))h\left(\frac nd\right)\\ =&\sum_{d\mid n}f(d)h\left(\frac nd\right) + \sum_{d\in n}g(d)h\left(\frac nd\right)\\ =&f(n)*h(n) + g(n)*h(n) \end{aligned}

  4. 与数乘的结合律:(cf)g=c(fg)(cf)*g = c(f*g)

  5. 单位元ϵf=f\epsilon * f = fϵ(n)=[n=1]\epsilon(n) = [n = 1]。即这个单位元函数当 n=1n = 1 时返回 11,当 n1n\not=1 时返回 00
    fϵ=dnf(d)ϵ(nd)\begin{aligned} f * \epsilon &= \sum_{d\mid n} f(d)\epsilon\left(\frac nd\right) \end{aligned}
    不难验证其正确性:只有当 d=nd = nϵ(n/d)\epsilon(n/d) 才等于 11,只保留了 f(n)f(n) 这一项。

  6. 狄利克雷逆元:f1f=ϵf^{-1}*f=\epsilon。考虑如何求这个逆元:
    ϵ=f1fϵ=dnf1(d)f(nd)ϵ(n)=f1(n)f(1)+dnd<nf1(d)f(nd)f1(n)=1f(1)(ϵ(n)dnd<nf1(d)f(nd))\begin{aligned} \epsilon &= f^{-1}*f\\ \epsilon &= \sum_{d\mid n}f^{-1}(d)f\left(\frac nd\right)\\ \epsilon(n) &= f^{-1}(n)f(1) + \sum_{d\mid n\land d < n} f^{-1}(d)f\left(\frac nd\right)\\ f^{-1}(n) &=\frac1{f(1)}\left(\epsilon(n) – \sum_{d\mid n\land d < n} f^{-1}(d)f\left(\frac nd\right)\right) \end{aligned}
    于是这个逆元就被求出来了。不难发现只要 ff 满足 f(1)0f(1)\not=0,其就会具有逆元

狄利克雷卷积就是研究函数的乘法结构

容斥原理到莫比乌斯反演

首先我们有算术基本定理
n=k=1ω(n)pkakn=\prod_{k = 1}^{\omega(n)} p_k^{a_k}
注意到此时 ω(n)\omega(n) 描述了一个数 nn 的约数种类数

定义欧米茄函数
Ω(n)=k=1ω(n)ak\Omega(n) = \sum_{k = 1}^{\omega(n)} a_k
它描述了一个数约数的总个数。

有了如上定义,我们就可以引入莫比乌斯函数了:
μ(n)=11(n)=[ω(n)=Ω(n)]×(1)ω(n)\mu(n) = 1^{-1}(n) = [\omega(n) =\Omega(n)]\times (-1)^{\omega(n)}
其中 11 的意义为常函数 f(n)=1f(n) = 1

分析一下这个神奇的定义:

  • 首先,此处定义莫比乌斯函数为 11 的逆元,注意到其具有如下性质
    dnμ(d)={1,n=10,n1\sum_{d\mid n}\mu(d) = \begin{cases} 1,&n = 1\\ 0,&n\not= 1 \end{cases}
    证明:设 n=k=1ωpkn' = \displaystyle\prod_{k=1}^{\omega}p_k,不难发现
    dnμ(d)=dnμ(d)=k=0ω(ωk)(1)k=(1+(1))ω\sum_{d\mid n}\mu(d) = \sum_{d\mid n'}\mu(d) = \sum_{k = 0}^\omega \binom{\omega}{k}(-1)^k = (1 + (-1))^\omega
    根据二项式定理,发现 ω=0\omega = 0n=1n = 1 时式子为 11,否则为 00

    所以我们知道 μ1=ϵ\mu*1 = \epsilon

  • 其次,分析后面的那坨式子:

    • ω(n)=Ω(n)\omega(n) = \Omega(n) 如果不成立,说明整个函数值就为 00 了。即 nn 具有平方因子
    • 如果 n=1n = 1,则显然 μ(1)=1\mu(1) = 1
    • 如果 ω(n)=Ω(n)\omega(n) = \Omega(n) 成立,那么函数的值即取决于函数具有奇数个还是偶数个质因子,前者 1-1,后者 11

很重要的一点:
idkt=ϵ\mathrm{id}^k * t = \epsilon
其中 idk(n)=nk\mathrm{id}^k(n) = n^kt(n)=μ(n)nkt(n) = \mu(n)n^k

简证:
(idkt)(n)=dndkμ(nd)(nd)k=nkdnμ(nd)=nkϵ(n)=ϵ(n)\begin{aligned} &(\mathrm{id}^k*t)(n)\\ =&\sum_{d\mid n}d^k\mu\left(\frac n d\right)\left(\frac nd\right)^k\\ =&n^k\sum_{d\mid n}\mu\left(\frac n d\right)\\ =&n^k\epsilon(n)\\ =&\epsilon(n) \end{aligned}
还有:ϕ=μid\phi = \mu * \mathrm{id},简证:
(μid)(n)=dnμ(d)(nd)=dnμ(nd)d=d=1n[dn]μ(nd)d\begin{aligned} &(\mu*\mathrm{id})(n)\\ =&\sum_{d\mid n}\mu(d)\left(\frac nd\right)\\ =&\sum_{d\mid n}\mu\left(\frac nd\right)d\\ =&\sum_{d = 1}^n[d\mid n]\mu\left(\frac nd\right)d\\ \end{aligned}

再来看看容斥原理:
i=1nSi=T{1,2,,n}(1)T1jTSj\left|\bigcup_{i=1}^nS_i\right|=\sum_{T\subseteq\{1,2,\cdots,n\}}(-1)^{|T| – 1}\left|\bigcap_{j\in T}S_j\right|
其描述的就是多个集合的并内的元素个数应该如何计算。

证明:考虑一个任意一个元素 xxkk 个集合 SiS_i 中,则在选取的时候:

  • 选取集合个数为 11 时,这个元素被加了 kk
  • 选取集合个数为 22 时,这个元素被减了 (k2)\displaystyle\binom k 2 次,因为一共是存在 (k2)\displaystyle\binom k 2 种选法使得选出的两个集合内包含 xx
  • 选取集合个数为 33 时,这个元素又被加了 (k3)\displaystyle\binom k 3
  • 选取集合个数为 44 时,这个元素被减了 (k4)\displaystyle\binom k 4
  • ……

总的来看,元素 xx 对答案产生的贡献即为
i=1k(ki)×(1)i+1\sum_{i = 1}^k \binom k i\times(-1)^{i + 1}
由二项式定理知其结果为 11,所以每个元素都正正好好被统计了一次,命题得证。

回到刚才的唯一分解式
n=k=1ω(n)pkakn=\prod_{k = 1}^{\omega(n)} p_k^{a_k}
我们将容斥原理换一种方式表达:
(1)0×k=1ωAi+k=1ω((1)k×1j1<j2<<jkω(p=1kAjp))=(-1)^0\times \bigcup_{k = 1}^{\omega}A_i + \sum_{k = 1}^\omega\left((-1)^k\times\sum_{1\le j_1< j_2 <\cdots <j_k\le \omega}\left(\bigcap_{p = 1}^kA_{j_p}\right) \right)= \varnothing
定义
f(a)=kωpk[aAk]f(a) = \prod_{k}^\omega p_k^{[a\in A_k]}

Bnd={a:df(a)}B_{\frac nd} = \{a:d\mid f(a)\}

其中 dnd\mid n

理解一下上面的式子:AkA_k 表示好多好多个集合,其中 AkA_k 被染色成了唯一分解式中的 pkp_k,然后 Ak1Ak2A_{k_1}\cap A_{k_2} 被染成了 pk1×pk2p_{k_1}\times p{k_2},依次类推。不难发现 f(a)f(a) 即代表了一个元素所在区域的颜色(如下图)

BndB_{\frac nd} 表示的就是图中的一个划分。举个例子:n=pk1ak1pk2ak2pk3ak3n = p_{k_1}^{a_{k_1}}p_{k_2}^{a_{k_2}}p_{k_3}^{a_{k_3}}(对应上面的图),当 d=pk1pk2d = p_{k_1}p_{k_2} 时,满足 df(a)d \mid f(a) 的也只有上图绿色的两块区域了,所以 Bnd=Ak1Ak2B_{\frac n d} = A_{k_1}\cap A_{k_2}

有了上述定义,我们有下面的式子:
dn(μ(d)×Bnd)=\sum_{d\mid n}(\mu(d)\times B_{\frac nd}) = \varnothing
发现其具有的就是容斥原理的形式。

进一步地,我们定义
g(nd)={a:d=f(a)}\begin{aligned} g\left(\frac n d\right) &= \{a:d = f(a)\} \end{aligned}
所以有
(1g)(nd)={a:df(a)}(1*g)\left(\frac nd\right) = \{a:d\mid f(a)\}

=g=μ1g=dn(μ(d)×(1g)(nd))\varnothing = g = \mu * 1 * g = \sum_{d\mid n}\left(\mu(d)\times(1*g)\left(\frac nd\right)\right)

联系上述式子尽量理解其含义

最终,莫比乌斯反演的表达式是如下的:
g(n)=dnf(d)    f(n)=dn(μ(d)×g(nd))=dn(g(d)×μ(nd))g=1f    f=μg\begin{aligned} g(n) = \sum_{d\mid n}f(d) & \iff f(n) = \sum_{d\mid n}\left(\mu(d)\times g\left(\frac nd\right) \right)=\sum_{d\mid n}\left(g(d)\times \mu\left(\frac nd\right) \right)\\ g = 1 * f&\iff f = \mu * g \end{aligned}
莫比乌斯反演的本质就是容斥原理

例题:洛谷P3455 [POI2007]ZAP-Queries

题意:求
i=1aj=1b[gcd(i,j)=d]\sum_{i = 1}^a\sum_{j = 1}^b[\gcd(i,j) = d]
解:

莫比乌斯反演时,一般要将含 gcd\gcd 的式子化为 gcd(i,j)=1\gcd(i,j) = 1 的形式才好反演。因此考虑设 a=a/da' = a/db=b/db' = b/d,原式化为
i=1aj=1b[gcd(i,j)=1]\sum_{i = 1}^{a'}\sum_{j = 1}^{b'}[\gcd(i,j) = 1]
然后就可以愉快的反演了。我们知道
μ1=ϵ\mu * 1 = \epsilon
所以
ϵ(n)=dnμ(d)\epsilon(n) = \sum_{d\mid n}\mu(d)
稍微变一下:
ϵ(gcd(i,j))=dgcd(i,j)μ(d)\epsilon(\gcd(i,j)) = \sum_{d\mid \gcd(i,j)}\mu(d)
所以原式化为
i=1aj=1bdgcd(i,j)μ(d)\sum_{i = 1}^{a'}\sum_{j = 1}^{b'}\sum_{d\mid \gcd(i,j)}\mu(d)
考虑写成枚举 dd 的形式:
i=1aj=1bdgcd(i,j)μ(d)=i=1aj=1bd=1min(a,b)[dgcd(i,j)]μ(d)=d=1min(a,b)μ(d)i=1aj=1b[dgcd(i,j)]\begin{aligned} &\sum_{i = 1}^{a'}\sum_{j = 1}^{b'}\sum_{d\mid \gcd(i,j)}\mu(d)\\ =&\sum_{i = 1}^{a'}\sum_{j = 1}^{b'}\sum_{d = 1}^{\min(a',b')}[d\mid\gcd(i,j)]\mu(d)\\ =&\sum_{d=1}^{\min(a',b')}\mu(d)\sum_{i = 1}^{a'}\sum_{j = 1}^{b'}[d\mid\gcd(i,j)] \end{aligned}
考虑到若要使 dgcd(i,j)d\mid \gcd(i,j)iijj 必须同为 dd 的倍数,根据乘法原理,原式可以化为
d=1min(a,b)μ(d)adbd\sum_{d=1}^{\min(a',b')}\mu(d)\left\lfloor\frac{a'}{d}\right\rfloor\left\lfloor\frac{b'}{d}\right\rfloor
预处理出 μ\mu 的前缀和,再套一个整除分块即可。

#include <cstdio>
#include <cctype>
#define FOR(i, a, b) for (int i = a; i <= b; ++i)

const int maxn = 5e4 + 7;

inline int read()
{
    char c = getchar();
    int s = 0;
    while (!isdigit(c))
        c = getchar();
    while (isdigit(c))
        s = 10 * s + c - '0', c = getchar();
    return s;
}

int mu[maxn], prime[maxn], tot, vis[maxn];

void getmu()
{
    mu[1] = 1;
    for (int i = 2; i <= maxn - 5; ++i)
    {
        if (!vis[i])
        {
            mu[i] = -1;
            prime[++tot] = i;
        }
        for (int j = 1; j <= tot && i * prime[j] <= maxn - 5; ++j)
        {
            vis[i * prime[j]] = 1;
            if (i % prime[j] == 0)
            {
                mu[i * prime[j]] = 0;
                break;
            }
            else
                mu[i * prime[j]] = -mu[i];
        }
    }
    FOR(i, 1, maxn - 5)
        mu[i] += mu[i - 1];
    return;
}

inline int min(int a, int b) {return a < b ? a : b;}

int main()
{
    getmu();
    int T = read();
    while (T--)
    {
        int a = read(), b = read(), d = read();
        a /= d, b /= d;
        int top = min(a, b), ans = 0;
        for (int l = 1, r = 0; l <= top && r <= top; l = r + 1)
        {
            r = min(top, min(a / (a / l), b / (b / l)));
            ans += (a / l) * (b / l) * (mu[r] - mu[l - 1]);
        }
        printf("%d\n", ans);
    }
    return 0;
}

洛谷 P2522 [HAOI2011]Problem b 与此题类似,只不过需要套一个容斥。

例题:洛谷 P2257 YY的GCD

题意:求
i=1nj=1m[gcd(i,j)P]\sum_{i = 1}^n\sum_{j = 1}^m[\gcd(i,j) \in P]
其中 PP 为质数集。

思路:

不妨枚举这个质数,如下:
kPi=1nj=1m[gcd(i,j)=k]\sum_{k\in P}\sum_{i = 1}^n\sum_{j = 1}^m[\gcd(i,j) = k]
然后就是像上题一样化简:
kPi=1nj=1m[gcd(i,j)=k]=kPi=1nj=1m[gcd(i,j)=1],n=nk,m=mk=kPi=1nj=1mdgcd(i,j)μ(d)=kPd=1min(n,m)μ(d)nkdmkd\begin{aligned} &\sum_{k\in P}\sum_{i = 1}^n\sum_{j = 1}^m[\gcd(i,j) = k]\\ =&\sum_{k\in P}\sum_{i = 1}^{n'}\sum_{j = 1}^{m'}[\gcd(i,j) = 1],\quad n' = \left\lfloor\frac nk \right\rfloor,m' = \left\lfloor\frac mk \right\rfloor\\ =&\sum_{k\in P}\sum_{i = 1}^{n'}\sum_{j = 1}^{m'}\sum_{d\mid \gcd(i,j)}\mu(d)\\ =&\sum_{k\in P}\sum_{d = 1}^{\min(n',m')}\mu(d)\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor \end{aligned}
此时这个式子看上去已经差不多了,但是实际运行的时候还是会 T 掉,考虑设 T=kdT = kd,有
kPd=1min(n,m)μ(d)nkdmkd=kPd=1min(n,m)μ(Tk)nTmT\begin{aligned} &\sum_{k\in P}\sum_{d = 1}^{\min(n',m')}\mu(d)\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor\\ =&\sum_{k\in P}\sum_{d= 1}^{\min(n',m')}\mu\left(\frac Tk\right)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\ \end{aligned}
换成枚举 TT,提前:
kPd=1min(n,m)μ(Tk)nTmT=T=1nnTmTkTkPμ(Tk)\begin{aligned} &\sum_{k\in P}\sum_{d= 1}^{\min(n',m')}\mu\left(\frac Tk\right)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\ =&\sum_{T = 1}^n\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\sum_{k\mid T\land k\in P}\mu\left(\frac Tk\right) \end{aligned}
kTkPμ(Tk)\displaystyle\sum_{k\mid T\land k\in P}\mu\left(\frac Tk\right) 是可以预处理的,预处理然后整除分块就可以了。实际实现的时候记得开 long long 并卡卡常(不要处处 long long

#include <cstdio>
#include <cctype>
#define FOR(i, a, b) for (register int i = a; i <= b; ++i)

const int maxn = 1e7 + 7;

inline int read()
{
    char c = getchar();
    int s = 0;
    while (!isdigit(c))
        c = getchar();
    while (isdigit(c))
        s = 10 * s + c - '0', c = getchar();
    return s;
}

inline int min(int a, int b) {return a < b ? a : b;}

typedef long long ll;

int mu[maxn], p[maxn], vis[maxn], tot;
ll sum[maxn];

void getmu()
{
    mu[1] = 1;
    FOR(i, 2, maxn - 5)
    {
        if (!vis[i])
        {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && i * p[j] <= maxn - 5; ++j)
        {
            vis[i * p[j]] = 1;
            if (i % p[j] == 0)
            {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
    FOR(j, 1, tot)
        for (int i = 1; i * p[j] <= maxn - 5; ++i)
            sum[i * p[j]] += mu[i];
    FOR(i, 1, maxn - 5)
        sum[i] += sum[i - 1];
    return;
}

signed main()
{
    getmu();
    int kase = read();
    while (kase--)
    {
        int n = read(), m = read();
        int top = min(n, m);
        ll ans = 0;
        for (int l = 1, r = 0; l <= top; l = r + 1)
        {
            r = min(top, min(n / (n / l), m / (m / l)));
            ans += (ll)(n / l) * (ll)(m / l) * (sum[r] - sum[l - 1]);
        }
        printf("%lld\n", ans);
    }
    return 0;
}

积性函数

定义:
gcd(a,b)=1    f(a×b)=f(a)×f(b)\gcd(a, b) = 1 \implies f(a\times b) = f(a)\times f(b)
由定义知 f(1)=1f(1) = 1

一些魔术:
gcd(a,b)=1    dabf(d)=d1ad2b(f(d1)×f(d2))=(d1af(d1))×(d2bf(d2))\gcd(a,b) = 1\implies\sum_{d\mid ab}f(d) = \sum_{d_1 \mid a \land d_2\mid b}(f(d_1)\times f(d_2)) = \left(\sum_{d_1\mid a}f(d_1)\right)\times\left(\sum_{d_2\mid b}f(d_2)\right)
上面的式子应该比较好理解,利用这个性质我们可以将一个 nn 拆分开来求和再乘在一起。通过上式可以推出下式:
gcd(a,b)=1    (fg)(ab)=dab(f(d)×g(abd))=d1ad2b(f(d1)g(ad1)f(d2)g(bd2))=d1a(f(d1)g(ad1))×d2b(f(d2)g(bd2))=(fg)(a)×(fg)(b)\begin{aligned} \gcd(a,b) = 1\implies (f*g)(ab) &= \sum_{d\mid ab}\left(f(d)\times g\left(\frac{ab}{d}\right)\right)\\ &=\sum_{d_1\mid a\land d_2\mid b}\left(f(d_1)g\left(\frac{a}{d_1}\right)f(d_2)g\left(\frac{b}{d_2}\right) \right)\\ &=\sum_{d_1\mid a}\left(f(d_1)g\left(\frac{a}{d_1}\right)\right)\times\sum_{d_2\mid b}\left(f(d_2)g\left(\frac{b}{d_2}\right) \right)\\ &= (f*g)(a)\times(f*g)(b) \end{aligned}
其中 f,gf,g 为积性函数。

完全积性函数:
a,bN,f(a×b)=f(a)×f(b)\forall a,b\in\mathbb N^*,\:f(a\times b) = f(a)\times f(b)
例子:1(n)=11(n) = 1idk(n)=nk\mathrm{id}^k(n) = n^kid(n)=n\mathrm{id}(n) = n

不难发现其具有如下性质
f(pn)=(f(p))nf(p^n) = \left(f(p)\right)^n
显然。
(f×g)(f×h)=dn(f(d)×g(d)×f(nd)×h(nd))=f×dn(g(d)×h(nd))=f×(gh)\begin{aligned} (f\times g)*(f\times h) &= \sum_{d\mid n}\left(f(d)\times g(d)\times f\left(\frac n d \right)\times h\left(\frac nd \right) \right)\\ &=f\times\sum_{d\mid n}\left(g(d)\times h\left(\frac nd \right) \right)\\ &=f\times(g*h) \end{aligned}

除数函数:
τ(n)=dn1σ(n)=dnd\begin{aligned} \tau(n) &= \sum_{d\mid n}1\\ \sigma(n) &= \sum_{d\mid n}d \end{aligned}
前面的 τ(n)\tau(n) 描述了 nn 的约数个数,σ(n)\sigma(n) 则描述的是 nn 的约数和。

一般的,
σm(n)=dndm=1nm\sigma_m(n) = \sum_{d\mid n}d^m = 1 * n^m
不难发现,τ(n)=σ0(n)\tau(n) = \sigma_0(n)σ(n)=σ1(n)\sigma(n) = \sigma_1(n),且 σm(n)\sigma_m(n) 为积性函数,验证如下:

假定 gcd(a,b)=1\gcd(a,b) = 1
σm(a)×σm(b)=(1am)×(1bm)=d1ad1m×d2bd2m=d1ad2b(d1d2)m=σm(ab)\begin{aligned} \sigma_m(a)\times \sigma_m(b) &= (1*a^m)\times(1*b^m)\\ &=\sum_{d_1\mid a}d_1^m\times \sum_{d_2\mid b}d_2^m\\ &=\sum_{d_1\mid a\land d_2\mid b}(d_1d_2)^m\\ &=\sigma_m(ab) \end{aligned}
欧拉函数:
ϕ(n)=1kngcd(k,n)=11\phi(n) = \sum_{\substack{1\le k\le n\\\gcd(k, n) = 1}}1
更一般地:
ϕm(n)=1kjngcd(k1,k2,,km,n)=11=μidm\phi_m(n) = \sum_{\substack{1\le k_j\le n\\ \gcd(k_1,k_2,\cdots,k_m,n) = 1}}1 = \mu *\mathrm{id}^m
其意义为从 1,2,,n1,2,\cdots,n 中选出 mm 个不相同且两两互质的数的方案数。

下面证明 ϕm=μidm\phi_m = \mu * \mathrm{id}^m

f(i)f(i) 为任意完全积性函数
1kjngcd(k1,k2,,km,n)=1(f1(k1)f2(k2)fm(km))=1kjn([gcd(k1,k2,,km,n)=1]1jmfj(kj))=1kjndgcd(k1,k2,km,n)(μ(d)1jmfj(kj))=1kjndgcd(k1,k2,km,n)(μ(d)1jmkjd=kjfj(kjd))=1kjndgcd(k1,k2,km,n)(μ(d)1jmkjd=kjfj(d)1jmkjd=kjfj(kj))\begin{aligned} \sum_{\substack{1\le k_j\le n\\ \gcd(k_1,k_2,\cdots,k_m,n) = 1}}(f_1(k_1)f_2(k_2)\cdots f_m(k_m))&=\sum_{1\le k_j\le n}\left([\gcd(k_1, k_2, \cdots, k_m, n) = 1]\prod_{1\le j\le m}f_j(k_j) \right)\\ &=\sum_{\substack{1\le k_j\le n\\d\mid \gcd(k_1,k_2\cdots,k_m, n)}}\left(\mu(d)\prod_{1\le j\le m}f_j(k_j) \right)\\ &=\sum_{\substack{1\le k_j\le n\\d\mid \gcd(k_1,k_2\cdots,k_m, n)}}\left(\mu(d)\prod_{\substack{1\le j\le m\\k_j'd = k_j}}f_j(k_j'd) \right)\\ &=\sum_{\substack{1\le k_j\le n\\d\mid \gcd(k_1,k_2\cdots,k_m, n)}}\left(\mu(d)\prod_{\substack{1\le j\le m\\k_j'd = k_j}}f_j(d)\prod_{\substack{1\le j\le m\\k_j'd = k_j}}f_j(k_j') \right)\\ \end{aligned}
定义 sj(i)=k=1ifj(k)s_j(i) = \displaystyle\sum_{k = 1}^if_j(k),继续化简:
1kjngcd(k1,k2,,km,n)=1(f1(k1)f2(k2)fm(km))=1kjndgcd(k1,k2,km,n)(μ(d)1jmkjd=kjfj(d)1jmkjd=kjfj(kj))=dn(μ(d)1jmfj(d)1jmsj(nd))=(μ(n)×1jmfj(n))1jmsj(n)\begin{aligned} \sum_{\substack{1\le k_j\le n\\ \gcd(k_1,k_2,\cdots,k_m,n) = 1}}(f_1(k_1)f_2(k_2)\cdots f_m(k_m))&=\sum_{\substack{1\le k_j\le n\\d\mid \gcd(k_1,k_2\cdots,k_m, n)}}\left(\mu(d)\prod_{\substack{1\le j\le m\\k_j'd = k_j}}f_j(d)\prod_{\substack{1\le j\le m\\k_j'd = k_j}}f_j(k_j') \right)\\ &=\sum_{d\mid n}\left(\mu(d)\prod_{1\le j\le m}f_j(d)\prod_{1\le j\le m}s_j\left(\frac nd\right)\right)\\ &=\left(\mu(n)\times\prod_{1\le j\le m}f_j(n)\right) * \prod_{1\le j\le m}s_j(n) \end{aligned}
若此时限定 fj(n)=1f_j(n) = 1,则此时 1jmsj(n)=nm\displaystyle \prod_{1\le j\le m}s_j(n) = n^m
1kjngcd(k1,k2,,km,n)=1(f1(k1)f2(k2)fm(km))=1kjngcd(k1,k2,,km,n)=11=(μ(n)×1jmfj(n))1jmsj(n)=μidm\begin{aligned} \sum_{\substack{1\le k_j\le n\\ \gcd(k_1,k_2,\cdots,k_m,n) = 1}}(f_1(k_1)f_2(k_2)\cdots f_m(k_m))&=\sum_{\substack{1\le k_j\le n\\ \gcd(k_1,k_2,\cdots,k_m,n) = 1}}1\\ &=\left(\mu(n)\times\prod_{1\le j\le m}f_j(n)\right) * \prod_{1\le j\le m}s_j(n)\\ &=\mu*\mathrm{id}^m \end{aligned}
例:给定 aj(1jm)a_j\:(1\le j\le m)
ParseError: KaTeX parse error: Undefined control sequence: \ at position 258: …(k_j = k_j'd)\\\̲ ̲=&\sum_{{1\le k…

杜教筛

杜教筛可以求解积性函数的前缀和。

设现在要求积性函数 ff 的前缀和,设 S(n)=i=1nf(i)S(n) = \displaystyle\sum_{i=1}^nf(i)

再随便找一个积性函数 gg,考虑他们的狄利克雷卷积的前缀和
i=1n(fg)(i)=i=1ndif(i)g(id)=i=1nj=1nig(i)f(j)=i=1n(g(i)i=1nif(j))=i=1ng(i)S(ni)\begin{aligned} &\sum_{i =1}^n(f*g)(i)\\ =&\sum_{i=1}^n\sum_{d\mid i}f(i)g\left(\frac id\right)\\ =&\sum_{i = 1}^n\sum_{j = 1}^{\left\lfloor\frac ni\right\rfloor}g(i)f(j) \\ =&\sum_{i = 1}^n\left(g(i)\sum_{i =1}^{\left\lfloor\frac ni\right\rfloor}f(j)\right)\\ =&\sum_{i=1}^ng(i)S\left(\left\lfloor\frac ni\right\rfloor\right) \end{aligned}

现在我们需要求 S(n)S(n),所以把它提出来,移项:
g(1)S(n)=i=1n(fg)(i)i=2ng(i)S(ni)g(1)S(n) = \sum_{i = 1}^n(f*g)(i) – \sum_{i=2}^ng(i)S\left(\left\lfloor\frac ni\right\rfloor\right)
假如我们可以快速求解 i=1n(fg)(i)\displaystyle\sum_{i = 1}^n(f*g)(i) 并使用数论分块,就可以较快地求得 g(1)S(n)g(1)S(n)

伪代码:

int getsumf(int n)
{
    int ans = fgs(n);// f * g 的前缀和
    for (int l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        ans -= (g_sum[r] - g_sum[l - 1]) * getsumf(n / l);
    }
    return ans;
}

可见其是一个递归调用的过程。其复杂度为 O(n34)O(n^{\frac 34}),证明如下:

设求出 S(n)S(n) 的时间复杂度为 T(n)T(n),要求出 S(n)S(n) 需要求出 n\sqrt{n}S(ni)\displaystyle S\left(\left\lfloor\frac ni\right\rfloor\right) 的值,故
T(n)=i=1nO(i)+O(ni)=O(n34)T(n) = \sum_{i = 1}^{\sqrt n}O(\sqrt i) + O\left(\sqrt\frac ni\right) = O\left(n^{\frac34}\right)
进一步地,我们可以先筛出前 mm 个答案,再使用杜教筛,优化后的复杂度为
T(n)=i=1nmni=O(nm)T(n) = \sum_{i = 1}^{\left\lfloor\frac nm\right\rfloor}\sqrt\frac ni = O\left(\frac{n}{\sqrt m}\right)
m=n23m = n^\frac23 时,T(n)=n23T(n) = n^\frac23

例题:洛谷 P4213 【模板】杜教筛(Sum)

题意:求
ParseError: KaTeX parse error: \cr valid only within a tabular/array environment
由狄利克雷卷积我们知道 μ1=ϵ\mu*1 = \epsilon,故 ϵ(n)=dnμ(d)\epsilon(n) = \displaystyle\sum_{d\mid n}\mu(d)

所以
S1(n)1(1)=i=1nϵ(i)i=2nS1(ni)S1(n)=1i=2nS1(ni)\begin{aligned} S_1(n)1(1) &= \sum_{i=1}^n\epsilon(i) – \sum_{i=2}^nS_1\left(\left\lfloor\frac ni\right\rfloor\right)\\ S_1(n)&=1-\sum_{i=2}^nS_1\left(\left\lfloor\frac ni\right\rfloor\right) \end{aligned}
由之前,我们也知道 ϕ1=id\phi * 1 = \mathrm {id}

所以
S2(n)1(1)=i=1nid(i)i=2n1(i)S2(ni)S2(n)=(n+1)n2i=2nS2(ni)\begin{aligned} S_2(n)1(1) &= \sum_{i=1}^n\mathrm{id}(i) – \sum_{i=2}^n1(i)S_2\left(\left\lfloor\frac ni\right\rfloor\right)\\ S_2(n)&=\frac{(n + 1)n}{2} – \sum_{i=2}^nS_2\left(\left\lfloor\frac ni\right\rfloor\right) \end{aligned}

#include <cstdio>
#include <cctype>
#include <unordered_map>
#define FOR(i, a, b) for (int i = a; i <= b; ++i)

inline int read()
{
    char c = getchar();
    int s = 0;
    bool x = 0;
    while (!isdigit(c))
        x = x | (c == '-'), c = getchar();
    while (isdigit(c))
        s = 10 * s + c - '0', c = getchar();
    return x?-s:s;
}

typedef unsigned long long ll;

const int maxn = 1e7 + 15;
std::unordered_map<int, ll> Mu, Phi;
ll mu[maxn], phi[maxn];
int vis[maxn], p[maxn], tot;

void init()
{
    mu[1] = 1, phi[1] = 1;
    FOR(i, 2, maxn - 5)
    {
        if (!vis[i])
            p[++tot] = i, mu[i] = -1, phi[i] = i - 1;
        for (int j = 1; j <= tot && i * p[j] <= maxn - 5; ++j)
        {
            vis[i * p[j]] = 1;
            if (i % p[j] == 0)
            {
                phi[i * p[j]] = phi[i] * p[j];
                break;
            }
            mu[i * p[j]] = -mu[i];
            phi[i * p[j]] = phi[i] * (p[j] - 1);
        }
    }
    FOR(i, 2, maxn - 5)
        phi[i] += phi[i - 1], mu[i] += mu[i - 1];
    return;
}

ll getSumPhi(ll n)
{
    if (n <= maxn)
        return phi[n];
    if (Phi[n])
        return Phi[n];
    ll ans = (ll)(n + 1ll) * (ll)n / 2ll;
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        ans -= (r - l + 1LL) * getSumPhi(n / l);
    }
    return Phi[n] = ans;
}

ll getSumMu(ll n)
{
    if (n <= maxn)
        return mu[n];
    if (Mu[n])
        return Mu[n];
    ll ans = 1;
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        ans -= (r - l + 1LL) * getSumMu(n / l);
    }
    return Mu[n] = ans;
}

int main()
{
    init();
    int T = read();
    while (T--)
    {
        int n = read();
        printf("%lld %lld\n", getSumPhi(n), getSumMu(n));
    }
    return 0;
}
最后修改日期:2021年1月28日

作者

留言

撰写回覆或留言

发布留言必须填写的电子邮件地址不会公开。