P6222-「P6156 简单题」加强版【莫比乌斯反演】

正题

题目链接:https://www.luogu.com.cn/problem/P6222

---

题目大意

给出(k)，(T)组询问给出(n)求

[sum_{i=1}^nsum_{j=1}^n(i+j)^k imes gcd(i,j) imes mu(gcd(i,j))^2 ]

---

解题思路

开始忘记了(k)次幂能线性筛后面全推错了，既然可以线性筛(k)次幂就把(gcd(i,j))提到前面来。

[sum_{d=1}^nmu(d)^2dsum_{i=1}^nsum_{j=1}^n[gcd(i,j)=d](i+j)^k ]

这里有一个比较巧妙的思路就是，可以把(mu(d)^2d)提前莫反之后一起套进莫反里。
设有

[sum_{d|n}g(d)=mu(n)^2nRightarrow g(n)=sum_{d|n}mu(frac{n}{d})mu(d)^2d ]

然后莫反原式就有一个比较简单的式子了

[sum_{d=1}^ng(d)d^ksum_{i=1}^{lfloorfrac{n}{d} floor}sum_{j=1}^{lfloorfrac{n}{d} floor}(i+j)^k ]

设(S(n)=sum_{i=1}^nsum_{j=1}^n(i+j)^k)，这个东西是可以线性预处理的。

先预线性筛+前缀和处理出(s(n)=sum_{i=1}^ni^k)，就有(S(n)=left(sum_{i=n+1}^{2n}s(i) ight)-left(sum_{i=1}^ns(i) ight))，再用一个前缀和就好了。

然后就有式子

[sum_{d=1}^ng(d)d^kS(lfloorfrac{n}{d} floor) ]

这个知道怎么筛(g(d)d^k)之后好像就可以(O(Tsqrt n))做了？发现(g)里面有个(n)还是很难搞。
然后可以一顿操作把下取整搞掉，设(S'(x)=S(x)-S(x-1))，那么就有

[sum_{d=1}^ng(d)d^ksum_{i=1}^{lfloorfrac{n}{d} floor}S'(i) ]

然后把(i)提出来就有

[sum_{i=1}^nsum_{d|i}g(d)d^kS'(frac{n}{d}) ]

然后就会发现这是两个函数狄利克雷卷积的前缀和？
(f(n)=g(n)n^k=n^ksum_{d|n}mu(frac{n}{d})mu(d)^2d)是积性函数，所以我们可以用埃氏筛的方法来快速搞这个东西。

首先我们需要知道对于一个质数(p)，(f(p^e))如何快速计算。
这里是有结论的

• 若(e=0)则(f(p^0)=f(1)=1)
• 若(e=1)则(f(p^1)=p imes (1 imes mu(1)^2mu(p)+p imes mu(1)mu(p)^2)=p(p-1))
• 若(e=2)则(f(p^2)=p^2 imes (mu(p^2)mu(1)^2+mu(p)mu(p)^2p+mu(1)mu(p^2)^2p^2)=-p^2)
• 若(egeq 3)那么若(mu(d) eq0)那么一定有(mu(frac{n}{d})=0)，所以(f(p^e)=0)

然后就可以开始做了，怎么快速计算狄利克雷卷积？埃氏筛给过我们方法，对于一个积性函数(f)，我们可以拆成若干个(f_p)满足

[f_p(x)=[p^e=x]f(x) (pin Pri) ]

（(Pri)是质数集）

然后有(f=prod_{pin Pri}f_p)（乘法表示狄利克雷积）

所以把所有的(f_p)乘到(S')里去就好了，这个是枚举所有质数的倍数来搞的，时间复杂度(O(nlog log n))

所以总共的时间复杂度就是(O(T+nlog log n))

---

code

#include<cstdio>
#include<cstring>
#include<algorithm>
#define uit unsigned int
using namespace std;
const int N=2e7+10;
int t,n,k,cnt;
uit pri[N],pw[N],s[N];
bool v[N];
uit power(uit x,int b){
    uit ans=1;
    while(b){
        if(b&1)ans=ans*x;
        x=x*x;b>>=1;
    }
    return ans;
}
void init(){
    s[1]=1;
    for(int i=2;i<=n*2;i++){
        if(!v[i])pri[++cnt]=i,s[i]=power(i,k);
        for(int j=1;j<=cnt&&i*pri[j]<=n*2;j++){
            v[i*pri[j]]=1;s[i*pri[j]]=s[i]*s[pri[j]];
            if(i%pri[j]==0)break;
        }
    }
    for(int i=1;i<=n*2;i++)
        pw[i]=s[i],s[i]=s[i-1]+s[i];
    for(int i=1;i<=n;i++)
        s[i]=s[i*2]+s[i*2-1]-s[i]*2;
    for(int j=1;j<=cnt&&pri[j]<=n;j++)
        for(int x=pri[j],i=n/x;i>=1;i--){
            s[i*x]+=s[i]*(x-1)*pw[x];
            if(i%x==0)s[i*x]-=s[i/x]*x*pw[x]*pw[x];
        }
    for(int i=1;i<=n;i++)s[i]+=s[i-1];
    return;
}
int main()
{
    scanf("%d%d%d",&t,&n,&k);
    init();
    while(t--){
        scanf("%d",&n);
        printf("%u
",s[n]);
    }
    return 0;
}