莫比乌斯反演总结

---

此总结写于18年10月左右, 忘记放在CSDN上共享了, 现在补一下.
同时提供markdown的源文件供下载修改成自己的模板, 也可在其中查看博客里面渲染失败的公式
链接: https://pan.baidu.com/s/1QjdjJa2ek-7NRBRIB4uDjw 提取码: 8t4c

---

莫比乌斯反演

目录

• 约定

• 基本概念

• 杜教筛

• 常用性质

• 小技巧

• 线性筛各种积性函数

• $\gcd$ 相关题目

• $lcm$ 相关题目

• 约数个数函数相关题目

• 约数和函数相关题目

• 欧拉函数相关题目

• 莫比乌斯函数相关题目

• 反演的巧妙应用

• 综合应用

• min_25筛及其应用

---

约定

1.任何一个数可以被表示成: $P=p_1^{q_1}\cdot p_2^{q_2}\cdot p_3^{q_3}\cdots p_t^{q_t}$

并约定: $q_1+q_2+q_3+...+q_t=tsum$.

2.$d(x)$表示$x$的约数个数,即约数个数函数.

3.$\sigma (x)$表示$x$的约数和,即约数和函数.

4.$w(x)$表示$x$的不同质因子个数.

5.未特别注明且 $n$ 不一定能整除 $k$ 的 $\frac{n}{k}$ 均表示 $\lfloor \frac{n}{k}\rfloor$ .

6.未注明范围的题均是 $10^5$ ~ $10^{12}$ 都可以做的题.

7.题目中未特别注明的题,都根据范围考虑取模或不取模.

---

基本概念

莫比乌斯函数

莫比乌斯函数是一个由容斥系数所构成的函数. $\mu (n)$ 的定义如下:

1. 当 $n=1$ 时, $\mu (n)=1$
2. 当 $n={\prod }_{i=1}^t{p}_i$ 且 $p_i$ 为不同的素数时, $\mu (n)=(-1{)}^t$ . 即 $n$ 分解质因子后没有幂次大于等于平方的质因子,此时函数值根据分解的个数决定.
3. 当 $n$ 含有任何质因子的幂次大于等于 $2$ 时, $\mu (n)=0$

莫比乌斯反演

• 定理 : $F(n)$ 和 $f(n)$ 是定义在非负整数集合上的两个函数,并且满足条件 :

  ​ $$F(n)=\sum _{d|n}f(d)$$

  那么存在一个结论 :

  ​ $$f(n)=\sum _{d|n}\mu (d)\cdot F(\frac{n}{d})$$

• 更常用的形式: 如果 $F(n)$ 和 $f(n)$ 满足条件:

  ​ $F(n)={\sum }_{n|d}f(d)$

​ 那么存在一个结论:

​ $f(n)={\sum }_{n|d}\mu (\frac{d}{n})\cdot F(d)$

---

杜教筛

作用: 以低于线性的时间复杂度来计算积性函数的前缀和(常为 $O(n^{\frac{2}{3}})$).

前提 : 对于一个积性函数 $f(i)$ ,如果知道 ${\sum }_{d|n}f(d)=g(n)$ ,且 $g(n)$ 非常容易求出来,那么就可以用杜教筛快速求出 ${\sum }_{i=1}^{n}f(i)$ .

预处理 : 通常预处理出 ${\sum }_{i=1}^{n}f(i)$ 的前 $n^{\frac{2}{3}}$ 项 .

原理 : 假设现在要求 ${\sum }_{i=1}^n\phi (i)$ , 我们知道 ${\sum }_{d|n}\phi (d)=n$ .

将 ${\sum }_{d|n}\phi (d)=n$ 拆分移项可以得到 : $\phi (n)=n-{\sum }_{d|n,d<n}\phi (d)$

将式子左右两端同时取和得 : ${\sum }_{i=1}^n\phi (i)={\sum }_{i=1}^{n}i-{\sum }_{i=1}^n{\sum }_{d|i,d<i}\phi (d)$ .

考虑将右边的式子更换枚举项由 $i$ 变成 $\frac{i}{d}$ 得: ${\sum }_{i=1}^n\phi (i)={\sum }_{i=1}^{n}i-{\sum }_{i=2}^n{\sum }_{d=1}^{\frac{n}{i}}\phi (d)$

此时容易发现可以可以分块+记忆化处理.

常用式子总结:

${\sum }_{i=1}^n\mu (i)=1-{\sum }_{i=2}^n{\sum }_{d=1}^{\frac{n}{i}}\mu (d)$

${\sum }_{i=1}^{n}i\cdot \mu (i)=1-{\sum }_{i=2}^{n}i\cdot {\sum }_{d=1}^{\frac{n}{i}}d\cdot \mu (d)$

${\sum }_{i=1}^n\phi (i)={\sum }_{i=1}^{n}i-{\sum }_{i=2}^n{\sum }_{d=1}^{\frac{n}{i}}\phi (d)$

${\sum }_{i=1}^{n}i\cdot \phi (i)={\sum }_{i=1}^n{i}^2-{\sum }_{i=2}^{n}i\cdot {\sum }_{d=1}^{\frac{n}{i}}d\cdot \phi (d)$

${\sum }_{i=1}^n{i}^2\cdot \phi (i)={\sum }_{i=1}^n{i}^3-{\sum }_{i=2}^n{i}^2\cdot {\sum }_{d=1}^{\frac{n}{i}}{d}^2\cdot \phi (d)$

---

常用性质

1.设 $g(n)={\sum }_{i=1}^n\lceil \frac{n}{i}\rceil$ , 那么 $g(n)=g(n-1)+d(n-1)+1$ , $d(i)$为约数个数函数.

2.$d(i\cdot j)={\sum }_{x|i}{\sum }_{y|j}[\gcd (x,y)==1]$.

3.$\sigma (i\cdot j)={\sum }_{x|i}{\sum }_{y|j}(x\cdot \frac{j}{y})[\gcd (x,y)==1]$ .

4.设 $f(n)$ 为 $1$ ~ $n$ 中无平方因子的数的数量.则 $f(n)={\sum }_{i=1}^n\mu (i^2)={\sum }_{i=1}^{\sqrt{n}}\mu (i)\cdot \frac{n}{i^2}$

5.${\sum }_{i=1}^{n}i\cdot \lfloor \frac{n}{i}\rfloor ={\sum }_{i=1}^n{\sum }_{j=1}^{\frac{n}{i}}j={\sum }_{i=1}^n\sigma (i)$ .

6.${\sum }_{d|n}\mu (d)=[n==1]$

7.${\sum }_{d|n}\frac{\mu (d)}{d}=\frac{\phi (n)}{n}$

8.${\sum }_{d|n}\phi (d)=n$

9.${\sum }_{i=1}^{n}i\cdot [\gcd (i,n)==1]=\frac{n\cdot \phi (n)+[n==1]}{2}$

10.$2^{w(n)}={\sum }_{d|n}{\mu }^2(d)$

11.$\phi (i^k)=i^{k-1}\cdot \phi (i)$

12.${\sum }_{i=1}^n{\sum }_{j=1}^n[\gcd (i,j)==1]=2\cdot {\sum }_{i=1}^n\phi (i)-1$

---

小技巧

1.对于要求${\sum }_{i=1}^n\frac{f[i]}{x}$ , ${\prod }_{i=1}^n\frac{f[i]}{x}$的情况(f为给定的数组,x为给定的值),如果$f[i]$的范围较小, 可以把$f$数组里面的值装入桶里,然后做一次前缀和,再根据x值来分块计算.从而将式子转化成: ${\sum }_{i=1}^{max}i\cdot (sum[x\cdot i+x-1]-sum[x\cdot i-1])$ , ${\prod }_{i=1}^{max}{i}^{(sum[x\cdot i+x-1]-sum[x\cdot i-1])}$ .

2.对于一个数组 $a_1,a_2,a_3\cdots a_n$ 求每一个元素的逆元, 可以先做一次前缀积得到数组 $S$.然后我们可以快速求出 $S_n$ ~ $S_1$ 的逆元,那么 $a_n$ 的逆元就等于 $S_n$ 的逆元乘以 $S_{n-1}$ .

下面以斐波那契数列为例求逆元:

void init()
{
    fib[0]=0;fib[1]=1;
    multi_fib[0]=multi_fib[1]=1;
    invfib[0]=invfib[1]=1;
    for(int i=2;i<=maxn;i++){
        fib[i]=(fib[i-1]+fib[i-2])%mod;
        multi_fib[i]=multi_fib[i-1]*fib[i]%mod;//预处理fib的前缀积
    }
    invfib[maxn]=qpow(multi_fib[maxn],mod-2);
    for(int i=maxn;i>=2;i--){
        invfib[i-1]=invfib[i]*fib[i]%mod;//先给invfib[i-1]赋值为前i-1项前缀积的逆元
        invfib[i]=invfib[i]*multi_fib[i-1]%mod;//更新invfib为fib[i]的逆元.
    }
}

3.有些函数本身就可以分块的处理出来的,在范围比较大的时候,可以用杜教筛的思想预处理出前 $n^{\frac{2}{3}}$ 项,然后分块的处理后 $n^{\frac{1}{3}}$ 项.

---

线性筛各种积性函数

关键: 线性筛积性函数 $f(n)$,只需要知道 $f(p)$ 的值和 $f(p^k)$ 即可做到线性筛.

1.莫比乌斯函数 $\mu (i)$:

​ $f(p)=-1$ , $f(p^k)=0(k>1)$.

void make()
{
    mu[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]) prime[++prime[0]]=i,mu[i]=-1;
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            prime[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            else mu[i*prime[j]]=-mu[i];
        }
    }
}

2.欧拉函数 $\phi (i)$ :

​ $f(p)=p-1,f(p^k)=f(p^{k-1})\cdot p$

void make()
{
    phi[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]) prime[++prime[0]]=i,phi[i]=i-1;
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            prime[i*prime[j]]=1;
            if(i%prime[j]==0) { phi[i*prime[j]]=phi[i]*prime[j]; break;}
            else phi[i*prime[j]]=phi[i]*phi[prime[j]];
        }
    }
}

3.约数个数函数 $d(i):$

​ $f(p)=2,f(p^k)=k+1$

void make()
{
    d[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i;
            d[i]=2; pow_cnt[i]=1;//pow_cnt记录最小质因子的幂次.
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            int to=i*prime[j]; prime[to]=1;
            if(i%prime[j]==0){
                pow_cnt[to]=pow_cnt[i]+1;
                d[to]=d[i]/(pow_cnt[i]+1)*(pow_cnt[to]+1);
                break;
            }
            else{
                pow_cnt[to]=1; d[to]=d[i]*d[prime[j]];
            }
        }
    }
}

4.约数和函数 $\sigma (i):$

​ $f(p)=p+1,f(p^k)=p^0+p^1+\cdots +p^k$

void make()
{
    d_sum[1]=1;//表示约数和函数
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i;  d_sum[i]=i+1;//素数的话约数和就是i+1.
            prime_pow_count[i]=i+1;//素数对应的最小质因子幂和就是i+1.
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            int to=i*prime[j]; prime[to]=1;
            if(i%prime[j]==0){//如果不是互素的话,说明当前的i中至少含有一个prime[j].
                prime_pow_count[to]=prime_pow_count[i]*prime[j]+1;
                //先除去之前的,再重新乘上去即可.
                d_sum[to]=d_sum[i]/prime_pow_count[i]*prime_pow_count[to];
                break; 
            }
            else {
                prime_pow_count[to]=prime[j]+1;
                d_sum[to]=d_sum[i]*d_sum[prime[j]];
                //由于约数和函数是积性函数,所以互质的情况下等于直接相乘.
            }
        }
    }
}

5.$\sigma (i^2)$ :

​ $f(p)=p^2+p+1,f(p^k)=p^0+p^1+p^2+\cdots +p^k+\cdots +p^{2\cdot k}$

void make()
{
    d_sum_pow_2[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i; d_sum_pow_2[i]=i*i+i+1;//素数的话约数和是i^2+i+1.
            prime_pow_2_count[i]=d_sum_pow_2[i];//素数对应的最小质因子幂和就是i^2+i+1.
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            int to=i*prime[j]; prime[to]=1;
            if(i%prime[j]==0){//如果不是互素的话,说明当前的i中至少含有一个prime[j].
            prime_pow_2_count[to]=prime_pow_2_count[i]*prime[j]*prime[j]+prime[j]+1;
             //先除去之前的,在重新乘上去即可.
          d_sum_pow_2[to]=d_sum_pow_2[i]/prime_pow_2_count[i]*prime_pow_2_count[to];
             break;
            }
            else {
                prime_pow_2_count[to]=prime[j]*prime[j]+prime[j]+1;
                d_sum_pow_2[to]=d_sum_pow_2[i]*d_sum_pow_2[prime[j]];
                //由于约数和函数是积性函数,所以互质的情况下等于直接相乘.
            }
        }
    }
}

6.$f(k)={\sum }_{d|k}d\ast \mu (d)$ :

​ $f(p)=1-p$ , $f(p^k)=f(p)$ .

void make()
{
    f[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]) prime[++prime[0]]=i,f[i]=1-i;
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            prime[i*prime[j]]=1;
            if(i%prime[j]==0){ f[i*prime[j]]=f[i]; break; } 
            else f[i*prime[j]]=f[i]*f[prime[j]];
        }
    }
}

7.$f(k)={\sum }_{d|k}{d}^x\cdot \mu (\frac{k}{d})$ :

​ $f(p)=p^x-1$ , $f(p^k)=f(p^{k-1})\cdot p^k$
​ 具体推的过程可以先取 $f(p)$ 观察出结果,再取 $f(p^2)$ 和 $f(p)$ 比较,从而得出结论

void make()
{
    f[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i;
            g[i]=qpow(i,x);//记录p^x
            f[i]=g[i]-1;
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            prime[i*prime[j]]=1;
            if(i%prime[j]==0){ f[i*prime[j]]=(f[i]*g[prime[j]])%mod; break; }
            else f[i*prime[j]]=(f[i]*f[prime[j]])%mod;
        }
    }
}

8.$f(k)={\sum }_{d|k}\phi (d)\cdot \mu (\frac{k}{d})$ :

​ $f(p)=p-2$ , $f(p^k)=p^k+p^{k-2}-2\cdot p^{k-1}$

void make()
{
    phimu[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i;
            low[i]=i; phimu[i]=i-2;//low数组代表i的最小质因子乘积.
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            int to=i*prime[j]; prime[to]=1;
            if(i%prime[j]==0){
                low[to]=low[i]*prime[j];
                if(low[i]==i){
                    //如果当前这个合数是一个质因子的幂的形式,单独讨论.(也可以直接在下面讨论)
                    if(i==prime[j]) phimu[to]=prime[j]*prime[j]+1-2*prime[j];
                    //单独讨论为2的情况.
                    else phimu[to]=phimu[i]*prime[j];
                }
                else phimu[to]=phimu[i/low[i]]*phimu[low[i]*prime[j]];
                //将最小质因子部分分割出去.
                break;
            }
            else{
                low[to]=prime[j];
                phimu[to]=phimu[i]*phimu[prime[j]];//互质就直接相乘.
            }
        }
    }
}

9.$x|i^k,x=\prod {p_i}^{q_i}$ ,令 $f_k(x)=\prod p_i^{\lceil \frac{qi}{k}\rceil }$ .

​ 比如求 ${\sum }_{i=1}^A[d|i^k]$ ,那么答案就等于 $\lfloor \frac{A}{f_k(d)}\rfloor$ .

​ 可以发现 $k=1$ 的时候, $f_1(x)=x$ ,下面给出 $f_2(x),f_3(x)$ 的线性筛代码:

void make()
{
    f2[1]=f3[1]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]){
            prime[++prime[0]]=i;
            low[i]=i;//low数组代表i的最小质因子乘积p^k.
            pow_cnt[i]=1;//pow_cnt记录i的最小质因子的幂.
            f2[i]=i; f3[i]=i;
        }
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            int to=i*prime[j];  prime[to]=1;
            if(i%prime[j]==0){
                low[to]=low[i]*prime[j];  pow_cnt[to]=pow_cnt[i]+1;
                //如果当前这个合数是一个质因子的幂的形式,单独讨论.(也可以直接在下面讨论)
                if(low[i]==i){
                    f2[to]=f2[i];  f3[to]=f3[i];
                    if(pow_cnt[to]%2==1) f2[to]*=prime[j];//如果恰好多1说明要加1
                    if(pow_cnt[to]%3==1) f3[to]*=prime[j];
                }
                else{
                    f2[to]=f2[i/low[i]]*f2[low[i]*prime[j]];//将最小质因子部分分割出去.
                    f3[to]=f3[i/low[i]]*f3[low[i]*prime[j]];
                }
                break;
            }
            else{
                low[to]=prime[j];  pow_cnt[to]=1;
                f2[to]=f2[i]*f2[prime[j]];//互质就直接相乘.
                f3[to]=f3[i]*f3[prime[j]];
            }
        }
    }
}

---

$\gcd$ 相关题目

---

题意:

给定 $n$ , $m$ , $d$ 求 ${\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==d]$ .

题解:

根据反演常用套路,设 :

$f(d)$ 为 ${\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==d]$ . $F(n)$ 为 ${\sum }_{i=1}^n{\sum }_{j=1}^{m}d|\gcd (i,j)$ .

则: ans = $f(d)$

​ = ${\sum }_{d|k}\mu (\frac{k}{d})\cdot F(k)$

可以观察出: 对答案有贡献的 $k$ 均是 $d$ 的倍数,而 $d$ 的倍数是有上限的,因此考虑更换枚举项为 $\frac{k}{d}$ .

​ = ${\sum }_{k=1}^{min(\frac{n}{d},\frac{m}{d})}\mu (k)\cdot F(k\cdot d)$

由于$F(k)=\frac{n}{k}\cdot \frac{m}{k}$. 考虑后面分块处理即可解决问题.

​ — $P3455$

---

题意:

给定 $T$ , $n$ , $m$ , $T$ 次询问, 输出 ${\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==Prime]$

$T\le 10^4$ , $N,M\le 10^6$ .

题解:

根据反演常用套路,设 :

$f(d)$ 为 ${\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==d]$ . $F(n)$ 为 ${\sum }_{i=1}^n{\sum }_{j=1}^{m}d|\gcd (i,j)$ .

设: $S$ 为质数集.

​ ans = ${\sum }_{p\in S}{\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==p]$

​ = ${\sum }_{p\in S}f(p)$

​ = ${\sum }_{p\in S}{\sum }_{p|d}\mu (\frac{d}{p})\cdot F(d)$ (一般把分块的部分和 $\mu$ 函数分开,所以需要更换枚举项)

​ = ${\sum }_{p\in S}{\sum }_{d=1}^{min(\frac{n}{p},\frac{m}{p})}\mu (d)\cdot F(d\cdot p)$ (由枚举 $d$ 变成 $\frac{d}{p}$ )

由于此时枚举 $p$ 的部分在最外面无法预处理,所以更换枚举项,改成枚举 $d\ast p$.

​ = ${\sum }_{k=1}^{min(n,m)}F(k)\cdot {\sum }_{p|k,p\in S}\mu (\frac{k}{p})$

此时容易发现后面一部分可以预处理处理.前一部分可以分块.

​ — $P2257$

---

题意:

求: ${\sum }_{i=1}^n{\sum }_{j=1}^n\gcd (i,j)$ , $n\le 10^{10}$

题解:

对于两个$\sum$ 前后都是 $n$ 的情况,通常将后面的 $\sum$ 上限换成 $i$ ,然后减去多余的部分.

​ ans = $(2\cdot {\sum }_{i=1}^n{\sum }_{j=1}^i\gcd (i,j))-{\sum }_{i=1}^{n}i$

设 ans${}^{{}^{\prime }}$ = ${\sum }_{i=1}^n{\sum }_{j=1}^i\gcd (i,j)$

​ ans${}^{{}^{\prime }}$ = ${\sum }_{d=1}^{n}d\cdot {\sum }_{i=1}^{\frac{n}{d}}{\sum }_{j=1}^i[\gcd (i,j)==1]$

此时容易发现 ${\sum }_{j=1}^i[\gcd (i,j)==1]$ 是 $i$ 以内与 $i$ 互质的数的个数 ,因此原式再次化简为:

​ = ${\sum }_{d=1}^{n}d\cdot {\sum }_{i=1}^{\frac{n}{d}}\phi (i)$

​ = ${\sum }_{i=1}^n\phi (i)\cdot {\sum }_{d=1}^{\frac{n}{i}}d$

此时可用杜教筛处理出${\sum }_{i=1}^n\phi (i)$ ,即可解决该题.

​ — $51Nod1237$

---

题意:

${\sum }_{i=1}^n{\sum }_{j=1}^m[\gcd (i,j)==k][tsum\le P]$ (tsum为k的质数幂和,P为题目给定)

$1\le n,m\le 10^5$,$0\le P\le 10^5$. Q次询问, $Q\le 10^5$.( 1的tsum为0 ).

题解:

用f(d)表示满足$d=\gcd (i,j),且1\le i\le n,1\le j\le m$的对数.

用F(d)表示满足$d|\gcd (i,j),且1\le i\le n,1\le j\le m$的对数.

可知: $F(d)=\lfloor \frac{n}{d}\rfloor$ $\cdot \lfloor \frac{m}{d}\rfloor$ , $f(x)={\sum }_{x|d}\mu (\frac{d}{x})\cdot F(d)$.

设k为满足$tsum\le P$的数. 则枚举每一个k既是答案:

​ ans = ${\sum }_{k}f(k)$

​ = ${\sum }_k{\sum }_{k|d}\mu (\frac{d}{k})\cdot F(d)$

​ = ${\sum }_{d=1}^{min(n,m)}F(d)\cdot {\sum }_{k|d}\mu (\frac{d}{k})$ (更换枚举项)

此时发现后面的可以线性筛出来后做前缀和即可.由于tsum最多不超过20.所以线性筛的时候记录

一下即可.

const int maxn = 5e5+7;
int prime[maxn+5],mu[maxn+5],pow_cnt[maxn+5];//pow_cnt记录i的质因子幂和.(出现+1即可).
long long sum[21][maxn+5];
void init()
{
    make_mu();
    for(int i=1;i<=maxn;i++)
        for(int j=i;j<=maxn;j+=i)
            sum[pow_cnt[i]][j]+=mu[j/i];//质因子个数为pow_cnt[i],能贡献的所有j.
    for(int i=1;i<=20;i++)
        for(int j=1;j<=maxn;j++)
            sum[i][j]+=sum[i-1][j];
    for(int i=0;i<=20;i++)
        for(int j=1;j<=maxn;j++)
            sum[i][j]+=sum[i][j-1];
}
int main()
{
    scanf("%d",&t);
    while(t--){
        scanf("%d%d%d",&n,&m,&p);
        if(p>20){printf("%lld\n",1LL*n*m);continue;}
        long long ans=0;
        for(int l=1,r;l<=min(n,m);l=r+1){
            r=min(n/(n/l),m/(m/l));
            ans+=1LL*(n/l)*(m/l)*(sum[p][r]-sum[p][l-1]);
        }
        printf("%lld\n",ans);
    }
}

​ — $hdu4746$

---

$lcm$ 相关题目

---

题意:

求 ${\sum }_{i=1}^n{\sum }_{j=1}^{m}LCM(i,j)$ , $n,m\le 10^7$ .

题解:

由于 $LCM(i,j)=\frac{i\cdot j}{\gcd (i,j)}$ .

​ ans = ${\sum }_{i=1}^n{\sum }_{j=1}^m\frac{i\cdot j}{\gcd (i,j)}$

​ = ${\sum }_{d=1}^{min(n,m)}{\sum }_{i=1}^{\frac{n}{d}}{\sum }_{j=1}^{\frac{m}{d}}i\cdot j\cdot d\cdot [\gcd (i,j)==1]$ (考虑枚举 $\gcd (i,j)$)

​ = ${\sum }_{d=1}^{min(n,m)}{\sum }_{i=1}^{\frac{n}{d}}{\sum }_{j=1}^{\frac{m}{d}}i\cdot j\cdot d\cdot {\sum }_{k|\gcd (i,j)}\mu (k)$ (考虑将$[\gcd (i,j)==1]$ 更换为 $\mu$)

​ = ${\sum }_{d=1}^{\min (n,m)}d\cdot {\sum }_{k=1}^{\min (\frac{n}{d},\frac{m}{d})}{k}^2\cdot \mu (k)\cdot {\sum }_{i=1}^{\frac{n}{k\cdot d}}i\cdot {\sum }_{j=1}^{\frac{m}{k\cdot d}}j$ (考虑枚举 k)

此时枚举 $d$ 即可过题. 但还可以通过更换枚举项将复杂度再次降低.

由于 ${\sum }_{i=1}^{x}i$ 可以 $O(1)$ 的计算出,所以我们用 $f(x)$ 来代表 ${\sum }_{i=1}^{x}i$ .

​ = ${\sum }_{d=1}^{min(n,m)}f(\frac{n}{d})\cdot f(\frac{m}{d})\cdot d\cdot {\sum }_{k|d}k\cdot \mu (k)$

考虑预处理出 $d\cdot {\sum }_{k|d}k\cdot \mu (k)$ .设$g(d)={\sum }_{k|d}k\cdot \mu (k)$ .由于 $g(d)$ 是积性函数,所以可以线性筛出.

最后做次前缀和即可解决这道题.

— $P1829$

---

题意:

求: ${\sum }_{i=1}^n{\sum }_{j=1}^{n}LCM(i,j)$ , $n\le 10^{10}$

题解:

对于两个$\sum$ 前后都是 $n$ 的情况,通常将后面的 $\sum$ 上限换成 $i$ ,然后减去多余的部分.

​ ans = $(2\cdot {\sum }_{i=1}^n{\sum }_{j=1}^i\frac{i\cdot j}{\gcd (i,j)})-{\sum }_{i=1}^{n}i$

设 ans${}^{{}^{\prime }}$ = ${\sum }_{i=1}^n{\sum }_{j=1}^i\frac{i\cdot j}{\gcd (i,j)}$

​ ans${}^{{}^{\prime }}$ = ${\sum }_{d=1}^{n}d\cdot {\sum }_{i=1}^{\frac{n}{d}}{\sum }_{j=1}^{i}i\cdot j\cdot [\gcd (i,j)==1]$

此时容易发现 ${\sum }_{j=1}^{i}j\cdot [\gcd (i,j)==1]$ 是 $i$ 以内与 $i$ 互质的数之和 ,因此原式再次化简为:

​ = ${\sum }_{d=1}^{n}d\cdot {\sum }_{i=1}^{\frac{n}{d}}\frac{i^2\cdot \phi (i)}{2}$

​ = ${\sum }_{i=1}^n\frac{i^2\cdot \phi (i)}{2}{\sum }_{d=1}^{\frac{n}{i}}d$

此时可用杜教筛处理出${\sum }_{i=1}^n\frac{i^2\cdot \phi (i)}{2}$ ,即可解决该题.

​ — $51Nod1238$

---

题意:

求满足 $a\le lcm(x,y)\le b$ 且 $x\le y$ 的二元组 $(x,y)$ 的数量. $1\le a,b\le 10^{11}$

题解:

容易发现可以答案可以转化为求满足 $1\le lcm(x,y)\le b$ 的数量减去 $1\le lcm(x,y)\le (a-1)$ 的数量.

问题转化为求满足 $1\le lcm(x,y)\le n$ 的数量,我们先不考虑 $x\le y$ ,算出这种情况下的答案再去掉重复的即可.

因此要求的式子变成了: ${\sum }_{i=1}^n{\sum }_{j=1}^n[lcm(i,j)\le n]$

​ = ${\sum }_{d=1}^n{\sum }_{i=1}^{\frac{n}{d}}{\sum }_{j=1}^{\frac{n}{d}}[i\cdot j\cdot d\le n]\cdot [\gcd (i,j)==1]$

​ = ${\sum }_{d=1}^n{\sum }_{k=1}^{\frac{n}{d}}\mu (k)\cdot {\sum }_{i=1}^{\frac{n}{k\cdot d}}{\sum }_{j=1}^{\frac{n}{k\cdot d}}[i\cdot j\cdot d\cdot k^2\le n]$

此时更换 $k$ 和 $d$ 的位置,同时将 $k^2$ 移到右边可得:

​ = ${\sum }_{k=1}^{\sqrt{n}}{\sum }_{d=1}^{\frac{n}{k^2}}{\sum }_{i=1}^{\frac{n}{k^2\cdot d}}{\sum }_{j=1}^{\frac{n}{k^2\cdot d}}[i\cdot j\cdot d\le \frac{n}{k^2}]$

观察发现后面一部分求的是形如 $[a\cdot b\cdot c\le n]$ 的东西.

因此令 $a\le b\le c$ ,那么 $a$ 只需要枚举到 $n^{\frac{1}{3}}$ ,然后 $b$ 枚举到 $\sqrt{\lfloor \frac{n}{a}\rfloor }$ , $c$ 的范围可以直接算出.

然后统计三个数相同,两个数相同,三个数都不同的情况即可解决该题.

long long cal(long long n)//计算1~n中LCM(i,j)<=n的对数.
{
    if(!n) return 0;
    long long res=0;
    for(int k=1;1LL*k*k<=n;k++){
        if(mu[k]==0) continue;
        long long x=n/(1LL*k*k);
        long long now=0;//计算a*b*c<=x的对数.
        for(int i=1;1LL*i*i*i<=x;i++){//枚举第一位.
            for(int j=i+1;1LL*j*j*i<=x;j++)//枚举第二位
                now+=1LL*(x/(1LL*i*j)-j)*6+3;//a<b<c时累计num*6,a<b=c时累计3.
            now+=1LL*(x/(1LL*i*i)-i)*3+1;//a=b<c时累计num*3,a=b=c时累计1.
        }
        if(mu[k]==1) res+=now;
        else res-=now;
    }
    return (res+n)>>1;//程序计算的是无顺序的,答案要求有顺序,所以(ans+n)/2;
}

​ — $51Nod1222$

---

约数个数函数相关题目

---

题意:

求 ${\sum }_{i=1}^n{\sum }_{j=1}^{m}d(i\cdot j)$ $n,m\le 10^5$

题解:

​ $ans={\sum }_{i=1}^n{\sum }_{j=1}^m{\sum }_{x|i}{\sum }_{y|j}[\gcd (x,y)==1]$

​ $={\sum }_{x=1}^n{\sum }_{y=1}^m\frac{n}{x}\cdot \frac{m}{y}\cdot [\gcd (x,y)==1]$ (由上一步更换枚举项得到)

​ $={\sum }_{x=1}^n{\sum }_{y=1}^m\frac{n}{x}\cdot \frac{m}{y}\cdot {\sum }_{k|\gcd (x,y)}\mu (k)$

​ $={\sum }_{k=1}^{min(n,m)}\mu (k)\cdot {\sum }_{x=1}^{\frac{n}{k}}\frac{n}{k\cdot x}\cdot {\sum }_{y=1}^{\frac{m}{k}}\frac{m}{k\cdot y}$

可以发现后一部分可以分块处理,因此问题得到解决.

​ — $BZOJ3994$

---

题意:

求 ${\sum }_{i=1}^n\sigma (i)$ , $n<2^{63}$ , $10^5$ 组数据 , $15s$ . (神仙题)
题目给的是 $\sigma$ , 但求的是约数个数.

题解:

#include <bits/stdc++.h>
using namespace std;
const int maxn = 1e6+7;
typedef long long ll;
typedef __int128 lll;
int t;
ll n;
struct Node{
    ll x,y;
    Node(ll _x=0,ll _y=0):x(_x),y(_y){}
    inline Node operator + (const Node &B){return Node(x+B.x,y+B.y);}
}S[maxn];
int top=0;
inline bool inner(ll x,ll y){return n<x*y;}
inline bool steep(ll x,Node v){return (lll)n*v.x<=(lll)x*x*v.y; }
lll solve()
{
    int i,crn = cbrt(n);//获取n的立方根.
    ll srn = sqrt(n);//获取n的平方根.
    ll x = n/srn;
    ll y = srn+1;
    lll res=0;struct Node L,R,M;
    S[++top]=(Node(1,0));S[++top]=(Node(1,1));
    while(true){
        for(L=S[top--];inner(x+L.x,y-L.y);x+=L.x,y-=L.y)
            res+=x*L.y+(L.y+1)*(L.x-1)/2;
        if(y<=crn) break;
        for(R=S[top];!inner(x+R.x,y-R.y);R=S[--top]) L=R;
        while(M=(L+R),1){
            if(inner(x+M.x,y-M.y)) S[++top]=R=M;
            else{
                if(steep(x+M.x,R))  break;
                L=M;
            }
        }
    }
    for(i=1;i<y;i++) res+=n/i;
    return res*2-srn*srn;
}
inline void print(lll x)
{
    static char buf[36];
    if(!x) {putchar(48); return ;}int i=0;
    for(;x;buf[++i]=x%10|48,x/=10);
    for(;i;--i) putchar(buf[i]);
}
int main()
{
    scanf("%d",&t);
    while(t--){
        scanf("%lld",&n);
        print(solve());
        putchar(10);
    }
    return 0;
}

​ — $SPOJ$ $DIVCNT1$

---

题意:

求 ${\sum }_{i=1}^{n}d(i^2)$ $n\le 10^{12}$

题解:

容易发现我们对于 $d(i^2)$ 没办法直接处理,所以对于 $d(n^2)$ ,我们考虑那些 $n^2$ 有的因子,而 $n$ 没有的因子,可以知道,对于 $n=p_1^{q_1}\cdot {p_2}^{q_2}\cdots {p_t}^{q_t}$ ,每一个 $n$ 有的约数对应质因子指数上分别加上 $0$ 或 $q_i$ ,可以构成一个 $n^2$ 的因子.

即: $d(n^2)={\sum }_{d|n}{2}^{w(d)}$ ( $w(d)$表示 $d$ 的不同质因子个数 ) .

​ $ans={\sum }_{i=1}^n{\sum }_{d|n}{2}^{w(d)}$

考虑 $2^{w(d)}$ , 本质上是 $d$ 的所有约数中无平方因数的个数,他们有一个特点, $\mu$ 要么是$-1$要么是$1$ .

​ $={\sum }_{i=1}^n{\sum }_{d|i}{\sum }_{k|d}{\mu }^2(k)$

​ $={\sum }_{i=1}^n{\sum }_{k|i}{\mu }^2(k)\cdot d(\frac{i}{k})$

​ $={\sum }_{k=1}^n{\mu }^2(k)\cdot {\sum }_{i=1}^{\frac{n}{k}}d(i)$

发现 ${\sum }_{k=1}^n{\mu }^2(k)$ 等于 $1$ ~ $n$ 中无平方因子数的个数,等于 ${\sum }_{i=1}^{\sqrt{n}}\mu (i)\cdot \frac{n}{i^2}$

因此分块即可解决此题.

— $SPOJ$ $DIVCNT2$

---

约数和函数相关题目

---

题意:

求 ${\sum }_{i=1}^n{\sum }_{j=1}^n\sigma (i\cdot j)$ $n\le 10^9$

题解:

​ $ans={\sum }_{i=1}^n{\sum }_{j=1}^n{\sum }_{x|i}{\sum }_{y|j}x\cdot \frac{j}{y}\cdot [\gcd (x,y)==1]$

​ $={\sum }_{x=1}^n{\sum }_{y=1}^{n}x\cdot {\sum }_{j=1}^{\frac{n}{y}}j\cdot [\gcd (x,y)==1]$

​ $={\sum }_{x=1}^n{\sum }_{y=1}^{n}x\cdot {\sum }_{j=1}^{\frac{n}{y}}j\cdot {\sum }_{k|\gcd (x,y)}\mu (k)$

​ $={\sum }_{k=1}^{n}k\cdot \mu (k)\cdot {\sum }_{x=1}^{\frac{n}{k}}x\cdot {\sum }_{y=1}^{\frac{n}{k}}{\sum }_{j=1}^{\frac{n}{k\cdot y}}j$

由于 ${\sum }_{i=1}^{n}i\cdot \lfloor \frac{n}{i}\rfloor ={\sum }_{i=1}^n{\sum }_{j=1}^{\frac{n}{i}}j={\sum }_{i=1}^n\sigma (i)$

所以: $={\sum }_{k=1}^{n}k\cdot \mu (k)\cdot ({\sum }_{i=1}^{\frac{n}{k}}\sigma (i){)}^2$

因此分块加杜教筛即可解决该题. (后面的部分可以用杜教筛的思想预处理前 $n^{\frac{2}{3}}$ 项,然后分块出后 $n^{\frac{1}{3}}$ 项)

​ — $51Nod1220$

---

题意:

求 ${\sum }_{i=1}^n{\sum }_{j=1}^{n}max(i,j)\cdot \sigma (i\cdot j)$ $n\le 10^6$ , $5\cdot 10^4$ 组数据.

题解:

易得: $ans=2\cdot {\sum }_{i=1}^n{\sum }_{j=1}^{i}i\cdot \sigma (i\cdot j)-{\sum }_{i=1}^{n}i\cdot \sigma (i^2)$

可以发现后面一部分是可以线性筛出来的.所以式子变成:

​ $an{s}^{{}^{\prime }}={\sum }_{i=1}^n{\sum }_{j=1}^{i}i\cdot \sigma (i\cdot j)$

​ $={\sum }_{i=1}^n{\sum }_{j=1}^{i}i\cdot {\sum }_{x|i}{\sum }_{y|j}x\cdot \frac{j}{y}\cdot [\gcd (x,y)==1]$

​ $={\sum }_{x=1}^n{x}^2\cdot {\sum }_{i=1}^{\frac{x}{i}}i\cdot {\sum }_{y=1}^{i\cdot x}{\sum }_{j=1}^{\frac{i\cdot x}{y}}j\cdot [\gcd (x,y)==1]$

​ $={\sum }_{x=1}^n{x}^2\cdot {\sum }_{i=1}^{\frac{x}{i}}i\cdot {\sum }_{y=1}^{i\cdot x}\sigma (y)\cdot [\gcd (x,y)==1]$

​ $={\sum }_{k=1}^n{k}^2\cdot \mu (k)\cdot {\sum }_{x=1}^{\frac{n}{k}}{x}^2\cdot {\sum }_{i=1}^{\frac{n}{k\cdot x}}i\cdot {\sum }_{y=1}^{i\cdot x}\sigma (y)$

​ $={\sum }_{k=1}^n{k}^2\cdot \mu (k){\sum }_{i=1}^{\frac{n}{k}}i\cdot \sigma (i)\cdot {\sum }_{j=1}^i\sigma (j)$ (由上一步更换枚举项可得)

此时可以分块做,但由于数据组数过多,令 $T=i\cdot k$ 可以再次变换:

​ $={\sum }_{T=1}^n{\sum }_{k|T}T\cdot k\cdot \mu (k)\cdot \sigma (\frac{T}{k})\cdot {\sum }_{i=1}^{\frac{T}{k}}\sigma (i)$

此时可以 $n\cdot \ln n$ 预处理出解然后 $O(1)$ 的输出答案.

​ — $51Nod1584$

---

欧拉函数相关题目

---

题意:

给定 $n,m$ ,求 ${\sum }_{i=1}^n{\sum }_{j=1}^m\phi (i\cdot j)$ , $n\le 10^5$ , $m\le 10^9$ .

---

题解:

首先观察到 $n$ 的范围较小,因此考虑令 $S(n,m)={\sum }_{i=1}^m\phi (n\cdot i)$ , 则 $ans $ = ${\sum }_{i=1}^{n}S(i,m)$

现在考虑化简 ${\sum }_{i=1}^m\phi (n\cdot i)$ . 我们发现没有直接的公式可以化简 $\phi (n\cdot i)$ , 因此考虑从 $\phi (x)$ 的定义考虑.

由于 $\phi (x)$ $=$ $p_1^{q_1}\cdot p_2^{q_2}\cdot p_3^{q_3}\cdots p_t^{q_t}\cdot (1-\frac{1}{p_1})\cdot (1-\frac{1}{p_2})\cdot (1-\frac{1}{p_3})\cdots (1-\frac{1}{p_t})$

​ $=$ $p_1^{q_1-1}\cdot p_2^{q_2-1}\cdot p_3^{q_3-1}\cdots p_t^{q_t-1}\cdot (p_1-1)\cdot (p_2-1)\cdot (p_3-1)\cdots (p_t-1)$

所以我们可以把 $x$ 表示成 $y\cdot w$ ,分别对应上面的前一部分和后一部分.因此,原式可以化简成:

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (w\cdot i)$

此时再考虑将 $w$ 和 $i$ 的公因子都放到 $w$ 上:

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (\frac{w}{\gcd (w,i)}\cdot i\cdot \gcd (w,i))$

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (\frac{w}{\gcd (w,i)})\cdot \phi (i\cdot \gcd (w,i))$

由于 $\gcd (w,i)$ 有的因子 $i$ 也有,所以可以考虑将 $gcd(w,i)$ 提取出来

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (\frac{w}{\gcd (w,i)})\cdot \phi (i)\cdot \gcd (w,i)$

此时考虑利用欧拉函数的性质 ${\sum }_{d|n}\phi (d)=n$ 将 $gcd(w,i)$ 转化成 $\phi ()$ 的形式

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (\frac{w}{\gcd (w,i)})\cdot \phi (i)\cdot {\sum }_{k|\gcd (w,i)}\phi (k)$

此时观察发现 $\phi (\frac{w}{\gcd (w,i)})$ 和 $\phi (k)$ 之间有一定关系,所以可以再次转化:

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (i)\cdot {\sum }_{k|\gcd (w,i)}\phi (\frac{w}{k})$

​ $=$ $y\cdot {\sum }_{i=1}^m\phi (i)\cdot {\sum }_{k|w,k|i}\phi (\frac{w}{k})$

此时可以考虑更换枚举项的位置,把 $k$ 的位置提前,同时:

​ $=$ $y\cdot {\sum }_{k|w}\phi (\frac{w}{k})\cdot {\sum }_{i=1}^{\frac{m}{k}}\phi (k\cdot i)$

​ $=$ $y\cdot {\sum }_{k|w}\phi (\frac{w}{k})\cdot S(k,\frac{m}{k})$

此时发现后面一部分是一个递归的形式,且当 $n=1$ 时 $S(n,m)={\sum }_{i=1}^m\phi (i)$ 可以杜教筛算出来.

因此记忆化后即可解决该题.

​ — $BZOJ3512$

---

莫比乌斯函数相关题目

---

题意:

求 ${\sum }_{i=1}^m\mu (n\cdot i)$ , $n\le 10^{12}$ , $m\le 2\cdot 10^9$

题解:

设 $S(n,m)={\sum }_{i=1}^m\mu (n\cdot i)$ .

则 ans = $\mu (n)\cdot {\sum }_{i=1}^m\mu (i)\cdot [\gcd (n,i)==1]$

​ = $\mu (n)\cdot {\sum }_{i=1}^m\mu (i)\cdot {\sum }_{k|\gcd (n,i)}\mu (k)$

​ = $\mu (n)\cdot {\sum }_{k|n}\mu (k)\cdot {\sum }_{i=1}^{\frac{m}{k}}\mu (k\cdot i)$

​ = $\mu (n)\cdot {\sum }_{k|n}\mu (k)\cdot S(k,\frac{m}{k})$

发现化简得到一个递归的形式,当 $n==1$ 时为出口,用杜教筛求出 ${\sum }_{i=1}^m\mu (i)$ 即可.

小技巧: 只在最外面分解一次 $n$ 的质因子,然后在求 $S(n,m)$ 时,对于一个 $n$ ,暴力判断哪些原始 $n$ 的质因子可以整除现在的 $n $ ,然后二进制枚举 $k$ 即可.

​ — $18$ 年徐州网络赛 $easymath$

---

反演的巧妙应用

---

题意:

给定长度为n的数组a, 求${\sum }_{i=1}^n{\sum }_{j=1}^n\gcd (a[i],a[j])\cdot (\gcd (a[i],a[j])-1)$ .

$1\le n\le 10^4$,$1\le a[i]\le 10^4$.

题解:

考虑直接枚举$\gcd (a[i],a[j])==d$的对数统计. 用max表示a[i]的最大值.

用f(d) 表示$\gcd (a[i],a[j])==d$ 的对数, 用F(d)表示$d|\gcd (a[i],a[j])$ 的对数.

F(d)可以较为轻松的预处理出来.

容易得到式子: ans = ${\sum }_{d=1}^{max}(d^2-d)\cdot f(d)$

​ = ${\sum }_{d=1}^{max}(d^2-d)\cdot {\sum }_{d|x}\mu (\frac{x}{d})\cdot F(x)$

由于此题式子化简后系数是$d^2-d$ ,且范围在$10^4$,所以可以直接预处理后面一部分来计算.

当计算的是$gcd(a[i],a[j])$时,还可以继续化简.

— $hdu5212$

---

题意:

求: ${\prod }_{i=1}^n{\prod }_{j=1}^{m}fib[\gcd (i,j)]$ $n,m\le 10^6$

题解:

根据反演基本套路,可得: ans = ${\prod }_{d=1}^{min(n,m)}fib[d{]}^{f(d)}$

​ ans = ${\prod }_{d=1}^{min(n,m)}fib[d{]}^{{\sum }_{d|k}}\mu (\frac{k}{d})\cdot F(k)$

​ ans = ${\prod }_{k=1}^{min(n,m)}{\prod }_{d|k}fib[d{]}^{\mu (\frac{k}{d})\cdot F(k)}$

​ ans = ${\prod }_{k=1}^{min(n,m)}({\prod }_{d|k}fib[d{]}^{\mu (\frac{k}{d})}{)}^{F(k)}$

此时可以发现 ${\prod }_{d|k}fib[d{]}^{\mu (\frac{k}{d})}$ 这一部分可以预处理,然后做一个前缀积,求次逆元,就可以和 $F(k)$ 一起分块处理了.

​ — $BZOJ4816$

---

题意:

有一个长为 $l$ 的数组 $a$ ,初始全是 $0$ ,现在有 $q$ 个操作,第一种是 $1$ $n$ $d$ $v$ ,表示给下标 $x$ 满足 $\gcd (x,n)=d$ 的 $a_x$ 增加 $v$ ,第二种是 $2$ $x$ ,查询 ${\sum }_{i=1}^x{a}_i$ . $1\le n,d,v\le 2\cdot 10^5$ .

题解:

首先列出第一种操作的式子: $a_x+=v\cdot [\gcd (x,n)==d]$

尝试化简: $v\cdot [\gcd (x,n)==d]$ = $v\cdot [\gcd (\frac{x}{d},\frac{n}{d})==1]$ = $v\cdot {\sum }_{k|\frac{x}{d},k|\frac{n}{d}}\mu (k)$

此时观察可以发现,对于 $k|\frac{n}{d}$ ,能贡献的 $x$ 一定是 $k\cdot d$ 的倍数.但是我们又不好直接去计算每一个 $x$ .

因此考虑构造一个函数 $f(d)$ : $a_x={\sum }_{d|x}f(d)$ . 此时对于每一个给定的 $n,d$ ,我们只需要把对应的 $k\cdot d$ 给加上 $v\cdot \mu (k)$ 即可. 此时第二种操作的式子可以表示为: ${\sum }_{i=1}^x{a}_i={\sum }_{i=1}^x{\sum }_{d|i}f(d)={\sum }_{d=1}^{x}f(d)\cdot (\frac{x}{d})$ ,可以分块做.

体会:

如果一个题中,每次会更新一个值的倍数,我们可以通过反演构造来转化问题.

​ — $hdu4947$

---

综合应用

---

题意:

求 : ${\sum }_{i=1}^A{\sum }_{j=1}^B{\sum }_{k=1}^C\phi (\gcd (i,j^2,k^3))$ $A,B,C\le 10^7$

题解:

设 $f(d)={\sum }_{i=1}^A{\sum }_{j=1}^B{\sum }_{k=1}^C[\gcd (i,j^2,k^3)==d]$

$F(d)={\sum }_{i=1}^A{\sum }_{j=1}^B{\sum }_{k=1}^{C}d|\gcd (i,j^2,k^3)$

令 ans = ${\sum }_{d=1}^{min(A,B^2,C^3)}\phi (d)\cdot f(d)$

​ ans = ${\sum }_{d=1}^{min(A,B^2,C^3)}\phi (d)\cdot {\sum }_{d|k}\mu (\frac{k}{d})\cdot F(k)$

​ ans = ${\sum }_{k=1}^{min(A,B^2,C^3)}F(k)\cdot {\sum }_{d|k}\phi (d)\cdot \mu (\frac{k}{d})$

可以发现,后面一部分是积性函数,可以线性筛出来.所以现在问题转化为求 $F(k)$ .

容易发现 $F(d)=\lfloor \frac{A}{f_1(d)}\rfloor \lfloor \frac{B}{f_2(d)}\rfloor \lfloor \frac{C}{f_3(d)}\rfloor$ (由于 $F(d)$ 有上下,所以外层枚举 $k$ 的时候应取 $min(A,B,C)$.)

当 $x|i^k,x=\prod {p_i}^{q_i}$ ,令 $f_k(x)=\prod p_i^{\lceil \frac{qi}{k}\rceil }$

可以发现 $f_k(x)$ 也是积性函数,可以线性筛出来.因此该题得到解决.

​ — $hdu6428$

---

题意:

给定一个包含n个正整数的集合,有两种计算集合权值的办法:

1. 计算集合的所有排列的所有区间的gcd之和.
2. 对于从集合中选出k个数( $k\in [1,n]$ )的所有方案,选出的数的gcd⋅k之和.

分别计算两种方式的权值并取模,$1\le n\le 10^5$,元素num范围$1\le num\le 10^5$.

题解:

用cnt[i]代表集合中是i及i的倍数的数的个数. max表示集合中最大的元素.

f(k,x)表示从集合中选取k个元素且这k个元素的gcd为x的种数.

F(k,x)表示从集合中选取k个元素且这k个元素的gcd为x或x的倍数的种数.

显而易见: F(k,x) = $C_{cnt[x]}^k$.

第一种情况下的权值:

​ ans = ${\sum }_{k=1}^n{\sum }_{x=1}^{max}f(k,x)\cdot k!\cdot (n-k+1)!\cdot x$

​ = ${\sum }_{k=1}^{n}k!\cdot (n-k+1)!\cdot {\sum }_{x=1}^{max}x\cdot {\sum }_{x|d}\mu (\frac{d}{x})\cdot F(k,d)$ (反演常用套路.)

​ = ${\sum }_{k=1}^{n}k!\cdot (n-k+1)!\cdot {\sum }_{d=1}^{max}F(k,d)\cdot {\sum }_{x|d}\mu (\frac{d}{x})\cdot x$ (更换枚举项)

​ = ${\sum }_{k=1}^{n}k!\cdot (n-k+1)!\cdot {\sum }_{d=1}^{max}F(k,d)\cdot \phi (d)$

​ = ${\sum }_{d=1}^{max}\phi (d)\cdot {\sum }_{k=1}^{n}k!\cdot (n-k+1)\cdot F(k,d)$ (前后调换位置)

​ = ${\sum }_{d=1}^{max}\phi (d)\cdot cnt[d]!\cdot {\sum }_{k=1}^{cnt[d]}\frac{(n-k+1)!}{(cnt[d]-k)!}$ (可有F(k,d)的公式展开后化简得)

由手动推导或者公式:${\sum }_{k=1}^x\frac{(n-k+1)!}{(x-k)!}=\frac{(n+1)!}{(x-1)!\cdot (n-x+2)}$ 可将原式再次化简:

​ =${\sum }_{d=1}^{max}\phi (d)\cdot cnt[d]!\cdot \frac{(n+1)!}{(cnt[d]-1)!\cdot (n-cnt[d]+2)}$

第二种情况下的权值:

​ ans = ${\sum }_{k=1}^n{\sum }_{x=1}^{max}f(k,x)\cdot k\cdot x$

​ = ${\sum }_{k=1}^{n}k\cdot {\sum }_{d=1}^{max}\phi (d)\cdot F(k,d)$ (化简过程同上)

​ = ${\sum }_{d=1}^{max}\phi (d)\cdot {\sum }_{k=1}^{cnt[d]}k\cdot C_{cnt[d]}^k$ (k的范围根据其上限调整了)

​ = ${\sum }_{d=1}^{max}\phi (d)\cdot cnt[d]\cdot {\sum }_{k=1}^{cnt[d]}{C}_{cnt[d]-1}^k$ (将k带入组合数里面然后化简一下)

​ = ${\sum }_{d=1}^{max}\phi (d)\cdot cnt[d]\cdot 2^{cnt[d]-1}$ (再次化简)

体会:

给定集合选数gcd的题一定要先列出式子,然后再一步一步化简.

#include <bits/stdc++.h>
using namespace std;
const int maxn = 1e5+7;
const int mod = 258280327;
long long fac[maxn*2+5],inv[maxn*2+5];//预处理阶乘和阶乘的逆元.
int prime[maxn+5],phi[maxn+5],vis[maxn+5],cnt[maxn+5];//vis记录每个数字出现了多少次.cnt记录每个数字的倍数出现了多少次.
long long pow_2[maxn+5];
long long qpow(long long x,long long y)
{
    long long res=1;
    while(y){
        if(y&1) res=(res*x)%mod;
        x=(x*x)%mod;
        y>>=1;
    }
    return res;
}
void init()
{
    phi[1]=1;fac[0]=1;pow_2[0]=1;
    for(int i=2;i<=maxn;i++){
        if(!prime[i]) prime[++prime[0]]=i,phi[i]=i-1;
        for(int j=1;j<=prime[0]&&prime[j]<=maxn/i;j++){
            prime[i*prime[j]]=1;
            if(i%prime[j]==0) {phi[i*prime[j]]=phi[i]*prime[j];break;}
            else phi[i*prime[j]]=phi[i]*phi[prime[j]];
        }
    }
    for(int i=1;i<=2*maxn;i++) fac[i]=fac[i-1]*i%mod;
    inv[2*maxn]=qpow(fac[2*maxn],mod-2);
    for(int i=2*maxn-1;i>=0;i--) inv[i]=(1LL*inv[i+1]*(i+1))%mod;
    for(int i=1;i<=maxn;i++) pow_2[i]=1LL*pow_2[i-1]*2%mod;
}
long long C(int n,int m){if(n<0||m>n) return 0;return fac[n]*inv[m]%mod*inv[n-m]%mod;}
int main()
{
    int n,x,mx;init();
    while(~scanf("%d",&n)){
        memset(vis,0,sizeof(vis));
        memset(cnt,0,sizeof(cnt));mx=0;
        for(int i=1;i<=n;i++) scanf("%d",&x),vis[x]++,mx=max(mx,x);
        for(int i=1;i<=mx;i++)//与数组和gcd相关的题大部分都要做这一步.
            for(int j=i;j<=mx;j+=i)
                cnt[i]+=vis[j];//记录为i的倍数的数有多少个.
        long long ans1=0,ans2=0;
        for(int i=1;i<=mx;i++){
            if(!cnt[i]) continue;
            ans1=(ans1+1LL*phi[i]*fac[cnt[i]]%mod*fac[n-cnt[i]+1]%mod*C(n+1,cnt[i]-1)%mod)%mod;
            ans2=(ans2+1LL*cnt[i]*pow_2[cnt[i]-1]%mod*phi[i]%mod)%mod;
        }
        if(ans1>ans2) printf("Mr. Zstu ");
        else if(ans1==ans2) printf("Equal ");
        else printf("Mr. Hdu ");
        printf("%lld\n",max(ans1,ans2));
    }
    return 0;
}

​ — $hdu5321$

---

题意:

定义 $rad(n)$ 表示 $n$ 的因子中最大的无平方因子数. $f(n)=rad(n)\cdot \phi (\frac{n}{rad(n)})$

$g(n)={\sum }_{d|n}f(d)$ , 求 ${\sum }_{i=1}^{n}g(i)$ , 其中 $n\le 10^{12}$ .

题解:

首先,我们只考虑如何化简求 $g(n)$ , 观察发现,对于 $f(n)$ 这个函数来说,貌似没办法直接化简,因此,我们可以尝试拆分 $n$ 为素因子,对每一个素因子考虑贡献来计算答案,此时式子可以表示为 :​

​ $g(n)={\prod }_{i=1}^t{\sum }_{j=0}^{q_i}f({p_i}^j)$

​ $=$ ${\prod }_{i=1}^t(1+{\sum }_{j=1}^{q_i}f({p_i}^j))$

​ $=$ ${\prod }_{i=1}^t(1+{\sum }_{j=1}^{q_i}(p_i\cdot \phi (p_i^{j-1})))$

​ $=$ ${\prod }_{i=1}^t(1+p_i^{q_i})$

由 $g(n)$ 的式子化简可以发现,对于 $g(n)$ 有贡献的值要么不取某一个素因子,要么取完某一个素因子.

所以: $g(n)$ $=$ ${\sum }_{d|n}d\cdot [\gcd (d,\frac{n}{d})==1]$

​ $=$ ${\sum }_{i=1}^n{\sum }_{j=1}^{n}i\cdot [\gcd (i,j)==1]\cdot [i\cdot j==n]$

此时发现: ${\sum }_{i=1}^{n}g(i)={\sum }_{i=1}^n{\sum }_{j=1}^{n}i\cdot [\gcd (i,j)==1]\cdot [i\cdot j\le n]$

​ $=$ ${\sum }_{i=1}^n{\sum }_{j=1}^n{\sum }_{k|\gcd (i,j)}\mu (k)\cdot i\cdot [i\cdot j\le n]$

​ $=$ ${\sum }_{k=1}^{n}k\cdot \mu (k)\cdot {\sum }_{i=1}^{\frac{n}{k}}i\cdot {\sum }_{j=1}^{\frac{n}{k}}[i\cdot j\cdot k^2\le n]$

容易发现将前面的枚举规定严格后,后面的部分就可以直接去掉.

​ $=$ ${\sum }_{k=1}^{\sqrt{n}}k\cdot \mu (k)\cdot {\sum }_{i=1}^{\frac{n}{k^2}}i\cdot \lfloor \frac{n}{k^2\cdot i}\rfloor$

注意到最后一部分是 ${\sum }_{i=1}^{n}i\cdot \lfloor \frac{n}{i}\rfloor$ 的形式,即为 ${\sum }_{i=1}^n\sigma (i)$ ,因此最后式子化为:

​ $=$ ${\sum }_{k=1}^{\sqrt{n}}k\cdot \mu (k)\cdot {\sum }_{i=1}^{\frac{n}{k^2}}\sigma (i)$

— $ZOJ3881$

---

---

min_25筛

概念

---

问题:

对于积性函数 $f(x)$ ,求 ${\sum }_{i=1}^{n}f(i),n\le 10^{10}$

要求:

$f(x)$ 满足 $f(p)$ 可以用多项式表示,且可以快速求出 $f(p^k)$ .(p指质数)

预处理:

设 $P$ 为质数集合, $P_i$ 表示第 $i$ 个质数.

先不考虑怎么求答案,考虑求这个式子: ${\sum }_{i=2}^{n}f(i)\cdot [i\in P]$

设 $g(n,j)={\sum }_{i=2}^{n}f(i)\cdot [i\in P$ 或 $i$ 的最小质因子大于 $P_j]$

考虑 $g(n,j)$ 的转移方程.

当 $P_j^2>n$ 时, $P_j$ 不会产生贡献,有 $g(n,j)=g(n,j-1)$

当 $P_j^2\le n$ 时, 要减去最小质因子为 $P_j$ 的贡献,即减去 $f(P_j)\cdot (g(\frac{n}{P_j},j-1)-g(P_j-1,j-1))$

显然可以知道: 我们要求的 ${\sum }_{i=1}^{n}f(i)\cdot [i\in P]=g(n,|P|)$ ,由于无法直接求 $g(n,j)$ ,所以初始把所有的数都当作质数来求 $g(n,j)$ ,然后只取其中的质数,即可得到答案.因此在非递归写法中,初始化每一个 $g$ 数组为把 $2$ ~ $i$ 中每一个数都当作质数的情况下的值,(假设 $f(p)=p^k$ ,那么 $g(n,0)={\sum }_{i=2}^n{i}^k$ ) 然后再转移即可.注意 $g(P_j-1,j-1)$ 是可以直接求的,等于${\sum }_{i=1}^{j-1}f(P_i)$

计算答案:

当做完上一步之后,考虑求真正的前缀和.设: $S(n,j)={\sum }_{i=2}^{n}f(i)\cdot [i\in P$ 或 $i$ 的最小质因子大于等于 $P_j]$

此时要求的答案就是 $S(n,1)$ .把 $S(n,j)$ 分为质数,合数,1来分别计算,则:

​ $S(n,j)=g(n,|P|)-{\sum }_{i=1}^{j-1}f(P_i)+{\sum }_{k\ge j}{\sum }_e(f(P_k^e)\cdot S(\frac{n}{P_k^e},k+1)+f(P_k^{e+1}))$

预处理 $g$ 后,直接递归计算 $S(n,1)$ 即可.

复杂度: $O(\frac{n^{\frac{3}{4}}}{\log n})$

---

题目

---

题意:

求 ${\sum }_{i=1}^{n}d(i^k)$ $(mod$ $2^{64})$ $n,k\le 10^{10}$

题解:

考虑用min_25筛解决此题:

$f(1)=1,f(p)=k+1,f(p^c)=k\cdot c+1$ .

本来 $g(n,0)={\sum }_{i=2}^{n}k+1$ ,考虑在计算的过程中先不乘以 $k+1$ ,再最后计算 $S$ 的时候再乘上去,所以初始化 $g(n,0)={\sum }_{i=2}^n=n-1$ . 即 $g(n,|P|)\cdot (k+1)={\sum }_{i=2}^{n}f(i)\cdot [i\in P]$

#include <bits/stdc++.h>
using namespace std;
typedef unsigned long long ull;
const ull maxn = 233333;//maxn代表的是分块的块数,设置为sqrt(n).
ull n,k,cntp,m,sn,id1[maxn+5],id2[maxn+5],w[maxn+5];
ull prime[maxn+5],g[maxn+5];
bool vis[maxn+5];
ull S(ull x,ull y)
{
    if(x<prime[y]||x<=1) return 0;
    ull id=(x<=sn)?id1[x]:id2[n/x];//找到x对应的块.
    ull res = g[id]*(k+1)-(y-1)*(k+1);//加上前一部分.
    for(ull i=y;i<=cntp&&prime[i]*prime[i]<=x;i++){
        ull t=prime[i],t2=prime[i]*prime[i];
        for(ull j=1;t2<=x;j++,t=t2,t2*=prime[i])
            res+=((k*j+1)*S(x/t,i+1)+(k*(j+1)+1));
    }
    return res;
}
void solve()
{
    scanf("%lld%lld",&n,&k); sn=sqrt(n);
    cntp=m=0;//cntp代表质数的数量,m代表块的数量.
    for(ull i=2;i<=sn;i++){
        if(!vis[i]){
            prime[++cntp]=i;
        }
        for(ull j=1;j<=cntp&&i*prime[j]<=sn;j++){
            vis[i*prime[j]]=1;
            if(i%prime[j]==0) break;
        }
    }
    //g数组预处理第一维所有可能取值,然后第二层转移.
    for(ull l=1,r;l<=n;l=r+1){
        r=n/(n/l);//分块得右边界
        w[++m]=n/l;//记录每一块的id和值.
        g[m]=(w[m]-1);//初始化的上界为w[m],同时还要减掉1.
        if(w[m]<=sn) id1[w[m]]=m;//给两种情况下的块分别取id.
        else id2[n/w[m]]=m;
    }
    for(ull i=1;i<=cntp;i++){
        //枚举每一个prime[i],考虑其能贡献的所有w[j],然后更新一层.
        for(ull j=1;j<=m&&prime[i]*prime[i]<=w[j];j++){
            ull d=w[j]/prime[i];
            ull id=(d<=sn)?id1[d]:id2[n/d];//id表示n/pj对应的块.
            g[j]-=(g[id]-i+1);//减去转移方程对应的部分.
        }
    }
    printf("%llu\n",S(n,1)+1);//1要特殊处理.
}
int main()
{
    int t;scanf("%d",&t);
    while(t--) solve();
    return 0;
}

​ — $SPOJ$ $DIVCNTK$

---

题意:

定义一个函数 $f(x)$ , 满足:

1.$f(1)=1$

2.$f(p^c)=p⨁c$ ( $p$ 为质数, $⨁$ 表示异或)

3.$f(a\cdot b)=f(a)\cdot f(b)$ ( $a$ 和 $b$ 互质)

求 ${\sum }_{i=1}^{n}f(i)$ $mod(10^9+7)$

题解:

考虑用 min_25 筛 解决此题.

对应质数 $p$ ,有 $f(p)=p$ $xor$ $1$ . 除了$f(2)=3$ , 其余 $f(p)=p-1$ .

那么转化成积性函数 $g(x)=x,h(x)=1$ ,求出 $g(n,|P|),h(n,|P|)$

则 ${\sum }_{i=1}^{n}f(i)\cdot [i\in P]=g(n,|P|)-h(n,|P|)$ .注意遇到$f(2)$ 时要额外 $+2$ .

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod = 1e9+7;
const int maxn = 2e5+1e4+7;//maxn代表的是分块的块数,设置为sqrt(n).
const ll inv2 = 500000004;
int prime[maxn+5];
bool vis[maxn+5];
ll n,m,sn,cntp,id1[maxn+5],id2[maxn+5],w[maxn+5];
ll g[maxn+5],h[maxn+5],sump[maxn+5];
ll S(ll x,ll y)
{
    if(x<prime[y]||x<=1) return 0;
    ll id = (x<=sn)?id1[x]:id2[n/x];
    ll res = ((g[id]-h[id])-(sump[y-1]-(y-1)))%mod;
    if(y==1) res+=2;//有1的情况需要额外+2.
    for(ll i=y;i<=cntp&&prime[i]*prime[i]<=x;i++){
        ll t=prime[i],t2=prime[i]*prime[i];
        for(ll j=1;t2<=x;j++,t=t2,t2*=prime[i])
            res=(res+(prime[i]^j)*S(x/t,i+1)%mod+(prime[i]^(j+1))%mod)%mod;
    }
    return res;
}
int main()
{
    scanf("%lld",&n); sn=sqrt(n);
    cntp=m=0;//cntp代表质数的数量,m代表块的数量.
    for(ll i=2;i<=sn;i++){
        if(!vis[i]){
            prime[++cntp]=i;
            sump[cntp]=(sump[cntp-1]+i)%mod;//记录质数前缀和.
        }
        for(ll j=1;j<=cntp&&i*prime[j]<=maxn;j++){
            vis[i*prime[j]]=1;
            if(i%prime[j]==0) break;
        }
    }
    //g数组预处理第一维所有可能取值,然后第二层转移.
    for(ll l=1,r;l<=n;l=r+1){
        r=n/(n/l);//分块得右边界
        w[++m]=n/l;//记录每一块的id和值.
        g[m]=(w[m]%mod*((w[m]+1)%mod)%mod*inv2%mod+mod-1)%mod;//g[i]=i的和.
        h[m]=(w[m]-1)%mod;//h[i]=1的和.
        if(w[m]<=sn) id1[w[m]]=m;
        else id2[n/w[m]]=m;
    }
    for(ll i=1;i<=cntp;i++){
        for(ll j=1;j<=m&&prime[i]*prime[i]<=w[j];j++){
            ll d = w[j]/prime[i];
            ll id = (d<=sn)?id1[d]:id2[n/d];
            g[j]=((g[j]-prime[i]*(g[id]-sump[i-1])%mod)%mod+mod)%mod;
            h[j]=((h[j]-(h[id]-i+1))%mod+mod)%mod;
            //分别更新g[j],h[j].
        }
    }
    printf("%lld\n",(S(n,1)+mod+1+mod)%mod);
    return 0;
}

​ — $LOJ6053$

---

题意:

一个合数的真因数指不包含本身的因数,现在求 $l$ ~ $r$ 之间的合数的最大真因数之和. $l,r\le 5\cdot 10^9$

题解:

考虑在min_25筛中,初始化 $g(n,0)={\sum }_{i=2}^{n}i$ ,可以发现在转移方程中:

KaTeX parse error: Expected '}', got 'EOF' at end of input: …,j-1) & \text{P_j^2>nKaTeX parse error: Expected 'EOF', got '}' at position 1: }̲ \\ g(n,j-1)-f(…P_j^2\le{n}KaTeX parse error: Expected 'EOF', got '}' at position 1: }̲ \end{cases}

$g(\frac{n}{P_j},j-1)-g(P_j-1,j-1)$ 对应的就是满足最小质因子是 $P_j$ 的情况下的和.因此直接统计这一部分即可.

注意: 此题是只求合数的.直接套用时需要注意.

#include <bits/stdc++.h>
using namespace std;
const int maxn = 2e5+7;
typedef long long ll;
ll sn,cntp,m,w[maxn+5],g[maxn+5],id1[maxn+5],id2[maxn+5];
ll prime[maxn+5],psum[maxn+5];
bool vis[maxn+5];
ll solve(ll n)
{
    if(n<=2) return 0; sn=sqrt(n);
    cntp=m=0;
    memset(vis,false,sizeof(vis));
    for(ll i=2;i<=sn;i++){
        if(!vis[i]) prime[++cntp]=i,psum[cntp]=psum[cntp-1]+i;
        for(int j=1;j<=cntp&&i*prime[j]<=sn;j++){
            vis[i*prime[j]]=1;
            if(i%prime[j]==0) break;
        }
    }
    for(ll l=1,r;l<=n;l=r+1){
        r=n/(n/l); w[++m]=n/l;
        g[m]=(w[m]+2)*(w[m]-1)/2;//与i*(i+1)/2-1等效.
        if(w[m]<=sn) id1[w[m]]=m;
        else id2[n/w[m]]=m;
    }
    ll ans=0;
    for(ll i=1;i<=cntp;i++){
        for(ll j=1;j<=m&&prime[i]*prime[i]<=w[j];j++){
            ll d=w[j]/prime[i];
            ll id=(d<=sn)?id1[d]:id2[n/d];
            g[j]-=prime[i]*(g[id]-psum[i-1]);
            if(j==1) ans+=g[id]-psum[i-1];//只有j=1的情况才对答案有贡献.
        }
    }
    return ans;
}
int main()
{
    ll l,r;
    while(~scanf("%lld%lld",&l,&r))
        printf("%lld\n",solve(r)-solve(l-1));
    return 0;
}

​ — 最大真因数 ( $oj$ 挂掉了)

---

题意:

求 ${\sum }_{i=1}^N{\sum }_{j=1}^{N}sgcd(i,j{)}^k$ $mod$ $2^{32}$.其中 $sgcd(i,j)$ 表示 $i,j$ 的公约数中第二大的,如果 $\gcd (i,j)==1$ ,那么 $sgcd(i,j)=0$ . $N\le 10^9,k\le 50$

题解:

考虑 $f(n)=\frac{n}{minp(n)}$ . 反演一波将答案转化为 ${\sum }_{d=2}^{n}f(d{)}^k\cdot (2\cdot {\sum }_{i=1}^{\frac{n}{d}}\phi (i)-1)$

容易发现如果可以求出 $f(d{)}^k$ 的前缀和,那么后面一部分杜教筛出来,即可分块处理得到答案.

考虑如何求 $f(d{)}^k$ 的前缀和 ,可以发现 $f(d)$ 等价于 $d$ 的最大真因数(素数的部分单独讨论). 因此在求最大真因数的基础上加一个求自然数幂和即可解决问题.

#include<cmath>
#include<cstdio>
#define ll long long
#define ul unsigned int
#define fo(i, x, y) for(int i = x; i <= y; i ++)
using namespace std;
const int N = 1e6 + 5;
int n, k, sqr, m;
int bz[N]; ul p[N], phi[N];
int w[N], i1[N], i2[N];
ul s[51][51], g[N], sk[N], sp[N], sum[N], h[N], c[N];
ul ksm(ul x, ul y) {
    ul s = 1;
    for(; y; y /= 2, x = x * x)
        if(y & 1) s = s * x;
    return s;
}
void sieve(int n) {
    fo(i, 2, n) bz[i] = 0; p[0] = 0;
    fo(i, 2, n) {
        if(!bz[i]){//sk代表p^k,sp代表sk的前缀和.
            p[++ p[0]] = i; phi[i] = i - 1; sk[p[0]] = ksm(i, k);
            sp[p[0]] = sp[p[0] - 1] + sk[p[0]]; c[i] = sk[p[0]];
        }//c代表i^k.
        for(int j = 1; i * p[j] <= n; j ++) {
            int k = i * p[j]; bz[k] = 1;
            c[k] = c[i] * c[p[j]];
            if(i % p[j] == 0) { phi[k] = phi[i] * p[j]; break;}
            phi[k] = phi[i] * phi[p[j]];
        }
    }
    phi[1] = c[1] = 1;//做一次前缀和.
    fo(i, 2, n) phi[i] += phi[i - 1], c[i] += c[i - 1];
}
ul zm(int n) {//计算1~n的自然幂和.
    if(n <= 1e6) return c[n];
    ul S = 0;
    fo(j, 1, k) {
        ul p = 1;
        fo(i, n - j + 1, n + 1)
            p *= (i % (j + 1) == 0 ? i / (j + 1) : i);
        S += s[k][j] * p;
    }
    return S;
}
int bx[N]; ul bs[N];
ul dg(int x) {//杜教筛出phi的前缀和.
    if(x <= 1e6) return phi[x];
    int k = x <= sqr ? i1[x] : i2[n / x];
    if(bx[k]) return bs[k];
    bx[k] = 1; bs[k] = (ul) x * (x + 1) / 2;
    for(int i = 2, j; i <= x; i = j + 1) {
        j = x / (x / i);
        bs[k] -= dg(x / i) * (j - i + 1);
    }
    return bs[k];
}
ul ans;
int main() {
    scanf("%d %d", &n, &k); sqr = sqrt(n);
    sieve(N - 5);
    s[0][0] = 1; fo(i, 1, k) fo(j, 1, i) s[i][j] = s[i - 1][j - 1] + s[i - 1][j] * j;
    for(int i = 1, j; i <= n; i = j + 1) {
        j = n / (n / i); w[++ m] = n / i;
        if(w[m] <= sqr) i1[w[m]] = m; else i2[n / w[m]] = m;
        h[m] = w[m] - 1;//初始化质数的.
        g[m] = zm(w[m]) - 1;//初始化合数的.
    }
    p[0] = 0; while(p[p[0] + 1] <= sqr) p[0] ++;
    fo(j, 1, p[0]) for(int i = 1; i <= m && p[j] * p[j] <= w[i]; i ++) {
        int k = (w[i] / p[j] <= sqr) ? i1[w[i] / p[j]] : i2[n / (w[i] / p[j])];
        g[i] -= sk[j] * (g[k] - sp[j - 1]); h[i] -= h[k] - j + 1;
        sum[i] += g[k] - sp[j - 1];
    }
    fo(i, 1, m) sum[i] += h[i];//统计上质数的贡献.
    fo(i, 1, m) ans += (sum[i] - sum[i + 1]) * (dg(n / w[i]) * 2 - 1);
    printf("%u", ans);
}

​ — $51Nod1847$

---