[航海协会]求和

求和

题目概述

题解

既然是数学题，那我们就先来化化式子。
显然，看到$\gcd (i,j)$，那我们不妨枚举这个最大公因数是多少，再看看有多少对数它们的$\gcd$是这个。
有,
$$Ans=\sum _{d=1}^n(\sum _{i=1}^n\sum _{j=1}^n[gcd(i,j)=d])(\sum _{i=1}^K{f}_i(d))$$对于$f_i(x)$这一项，显然是当所有质因子的次数都不超过$i$的时候为$(-1{)}^{d(x)}$，其它时候都为$0$。
而这有显然是一个积性函数，我们可以很容易地用筛子筛出来。
所以我们主要考虑的是前面这个${\sum }_{i=1}^n{\sum }_{j=1}^n[gcd(i,j)=d]$怎么算。
这主要有两种方法，一种是大家都很熟悉的莫比乌斯反演，直接枚举这个，进行容斥，${\sum }_{d|t}\mu (\frac{t}{d})\lfloor \frac{n}{t}{\rfloor }^2={\sum }_{i=1}^{\lfloor \frac{n}{d}\rfloor }\mu (i)\lfloor \frac{\lfloor \frac{n}{d}\rfloor }{i}\rfloor$。
我们定义$h(x)={\sum }_{i=1}^x\mu (i)\lfloor \frac{x}{i}{\rfloor }^2$，则$Ans={\sum }_{d=1}^{n}h(\lfloor \frac{n}{d}\rfloor )({\sum }_{i=1}^K{f}_i(d))$。
显然，前面的$h(x)$是可以数论分块$O\left(\sqrt{x}\right)$计算，后面有可以整除分块按$\frac{n}{d}$计算。
后面的$f_i(d)$可以$\min 25$筛，按$\lfloor \frac{n}{d}\rfloor$位置计算前缀和，然后我们就得到了一个$O\left(n^{1-ϵ}\right)$的亚线性算法。
然后吗？你就$T$了。

好的，这也就是说我们得换一种方法计算我们的$h(\lfloor \frac{n}{d}\rfloor )$，同样的思路，还是计算有多少对在$[1,\frac{n}{d}]$以内的数对互质，这不是可以欧拉函数吗？
显然，$h(x)=2({\sum }_{i=1}^x\varphi (x))-1$，之间枚举数对内较大一个数是那个不就行了吗，注意减去重复的$(1,1)$。
这$\varphi (x)$不也是积性函数吗？我们这个前缀和可以与我们的$f_i(x)$一起计算了。
$\varphi (x)$的计算就是经典的$\min 25$筛问题了，先计算质数和与质数个数，减去后就得到质数处上的位置，再把其他质因子加回去即可。
至于$f_i(x)$的计算，是与$\mu (x)$的计算类似的你，你只需要在意把质数加回去时每个质数的次项不超过$i$即可。

时间复杂度$O\left(\frac{K{n}^{\frac{3}{4}}}{\ln n}\right)$。
注意$f$数组维度的顺序，这会严重地影响常数，就跟矩阵乘法一样。

源码

 #pragma GCC optimize(2)
 #pragma GCC optimize(3)
 #pragma GCC optimize("Ofast")
 #include<bits/stdc++.h>
 using namespace std;
 typedef long long LL;
 typedef pair<int,int> pii;
 typedef unsigned int uint;
 #define MAXN 200005
 #define pb push_back
 #define mkpr make_pair
 #define fir first
 #define sec second
 template<typename _T>
 void read(_T &x){
    _T f=1;x=0;char s=getchar();
    while(s<'0'||s>'9'){if(s=='-')f=-1;s=getchar();}
    while('0'<=s&&s<='9'){x=(x<<3)+(x<<1)+(s^48);s=getchar();}
    x*=f;
 }
 template<typename _T>
 _T Fabs(_T x){return x<0?-x:x;}
 int add(int x,int y,int p){return x+y<p?x+y:x+y-p;}
 void Add(int &x,int y,int p){x=add(x,y,p);}
 int qkpow(int a,int s,int p){int t=1;while(s){if(s&1)t=1ll*a*t%p;s>>=1;}return t;}
 int K,prime[MAXN],cntp,idx;
 LL n,a[MAXN],id1[MAXN],id2[MAXN];
 uint F[MAXN],H[MAXN][42],ans,G[MAXN],P[MAXN],pre[MAXN];
 bool oula[MAXN];
 void init(){
    int n1=sqrt(n);
    for(int i=2;i<=n1;i++){
        if(!oula[i])prime[++cntp]=i,pre[cntp]=pre[cntp-1]+i;
        for(int j=1;j<=cntp&&1ll*i*prime[j]<=n1;j++){
            oula[i*prime[j]]=1;
            if(i%prime[j]==0)break;
        }
    }
 }
 inline int getId(LL x){return x<=n/x?id1[x]:id2[n/x];}
 int main(){
    //freopen("sum.in","r",stdin);
    //freopen("sum.out","w",stdout);
    read(n);read(K);init();
    for(LL l=1,r;l<=n;l=r+1){
        r=n/(n/l);++idx;
        if(l<n/l)id2[l]=idx;else id1[n/l]=idx;
        a[idx]=n/r;F[idx]=1-a[idx];P[idx]=a[idx]-1;
        G[idx]=1ll*(a[idx]+1)*a[idx]/2LL-1;
    }
    for(int i=1;i<=cntp;i++){
        for(int k=1;k<=idx&&a[k]>=1ll*prime[i]*prime[i];k++)
            F[k]-=F[getId(a[k]/prime[i])]+(i-1),
            P[k]-=P[getId(a[k]/prime[i])]-(i-1),
            G[k]-=prime[i]*(G[getId(a[k]/prime[i])]-pre[i-1]);
    }
    for(int j=1;j<=idx;j++)for(int i=1;i<=K;i++)H[j][i]=F[j];
    for(int i=1;i<=idx;i++)F[i]=G[i]-P[i];
    for(int i=cntp,up=1;i>0;i--){
        while(a[up+1]>=1ll*prime[i]*prime[i])up++;
        for(int k=1;k<=up;k++){
            LL now=prime[i],nowp=1ll*prime[i]*now,nw=prime[i]-1;
            for(int j=1;nowp<=a[k];j++,now=nowp,nowp*=prime[i]){
                int t=getId(a[k]/now);F[k]+=nowp-now;
                F[k]+=1ll*nw*(F[t]-pre[i]+i);nw=nowp-now;
            }
        }
        for(int k=1;k<=up;k++){
            LL now=prime[i],nowp=1ll*prime[i]*now;uint *A=H[k];
            for(int l=1;nowp<=a[k];l++,now=nowp,nowp*=prime[i]){
                int t=getId(a[k]/now);uint *B=H[t];
                if(l+1&1)for(int j=l+1;j<=K;j++)A[j]--;
                else for(int j=l+1;j<=K;j++)A[j]++;
                if(l&1)for(int j=l;j<=K;j++)A[j]-=B[j]+i;
                else for(int j=l;j<=K;j++)A[j]+=B[j]+i;
            }
        }
    }
    for(int j=1;j<=idx;j++)for(int i=K;i>0;i--)H[j][i]-=H[j][i-1];
    ans=(F[1]+F[1]+1)*K;
    for(LL l=2,r;l<=n;l=r+1){
        r=n/(n/l);int z=getId(n/l);uint tmp=F[z]+F[z]+1,sum=0;
        int x=getId(r),y=getId(l-1);uint *A=H[x],*B=H[y];
        for(int j=1;j<=K;j++)sum+=(A[j]-B[j])*(K+1-j);
        ans+=sum*tmp;
    }
    printf("%u\n",ans%(((uint)1)<<30));
    return 0;
 }

谢谢！！！