• 题解-Sakuya's task


    题面

    Sakuya's task

    [left(sum_{i=1}^nsum_{j=1}^n varphi(gcd(i,j)) ight)mod 10^9+7 ]

    数据范围:(1le nle 10^{10})


    蒟蒻语

    考场爆零真开森。

    本来以为要卷 (1*1),没想到真要卷 (1*1),只不过要一个一个卷……

    考场上还以为要洲阁 ( t Min\_25)yayayaya.png


    正解

    先莫反操作一发:

    [egin{split} &sum_{i=1}^nsum_{j=1}^n varphi(gcd(i,j))\ =&sum_{d=1}^n varphi(d)sum_{i=1}^{lfloorfrac{n}{d} floor}sum_{j=1}^{lfloorfrac{n}{d} floor}epsilon(gcd(i,j))\ =&sum_{d=1}^nvarphi(d)sum_{k=1}^{lfloorfrac{n}{d} floor}mu(k)lfloorfrac{n}{dk} floor^2\ =&sum_{T=1}^nlfloorfrac{n}{T} floor^2sum_{d|T}varphi(d)mu(frac{T}{d})\ end{split} ]

    整除分块左边,右边杜教。

    第一次杜教:(f_1=varphi)(g_1=1)(f_1*g_1=id)

    第二次杜教:(f_2=varphi*mu)(g_2=1)(f_2*g_2=varphi=f_1)

    (f_2) 会多次调用 (f_1),但是内部调用的函数 (x) 集相等,所以可以一起求:

    //Dusieve
    bool vis[iN+1];
    int duphi[iN+1],dupm[iN+1];
    int Phi(ll x){return x<=N?phi[x]:duphi[n/x];}
    int Pm(ll x){return x<=N?pm[x]:dupm[n/x];}
    void Dusieve(ll x){
    	if(x<=N||vis[n/x]) return;
    	vis[n/x]=true;
    	for(ll l=2,r;l<=x;l=r+1){
    		r=x/(x/l),Dusieve(x/l);
    		(duphi[n/x]-=(ll)(r-l+1)*Phi(x/l)%mod)%=mod;
    		(dupm[n/x]-=(ll)(r-l+1)*Pm(x/l)%mod)%=mod;
    	}
    	(duphi[n/x]+=(ll)x%mod*(x%mod+1)/2%mod)%=mod;
    	(dupm[n/x]+=duphi[n/x])%=mod;
    	(duphi[n/x]+=mod)%=mod,(dupm[n/x]+=mod)%=mod;
    }
    

    还有个问题:怎么线性筛 (varphi*mu)

    其实可以狄利克雷前缀和一下,但是这里有个更妙的方法:

    (mu)(varphi) 为积性,(varphi*mu) 必为积性。

    根据 (mu) 函数的性质与找规律可得:

    [(varphi*mu)(1)=1\ (varphi*mu)(p)=p-2\ (varphi*mu)(p^2)=p(varphi*mu)(p)+(varphi*mu)(1)\ (varphi*mu)(p^3)=p(varphi*mu)(p^2)\ ]

    然后根据积性函数性质,就可以线性筛了。

    时间复杂度 (Theta(n^{frac{2}{3}}))


    代码

    取模坑死蒟蒻,细节会有注释。

    #include <bits/stdc++.h>
    using namespace std;
    
    //Start
    typedef long long ll;
    typedef double db;
    #define mp(a,b) make_pair(a,b)
    #define x first
    #define y second
    #define be(a) a.begin()
    #define en(a) a.end()
    #define sz(a) int((a).size())
    #define pb(a) push_back(a)
    const int inf=0x3f3f3f3f;
    const ll INF=0x3f3f3f3f3f3f3f3f;
    
    //Data
    const int mod=1e9+7;
    ll n; int ans;
    
    //Sieve
    const int N=1e7,iN=1e3;
    bool np[N+1];
    int phi[N+1],pm[N+1],cnt,p[N];
    void Sieve(){
    	np[1]=true,phi[1]=pm[1]=1;
    	for(int i=2;i<=N;i++){
    		if(!np[i]) p[cnt++]=i,phi[i]=i-1,pm[i]=i-2;
    		for(int j=0;j<cnt&&i*p[j]<=N;j++){
    			np[i*p[j]]=1;
    			if(i%p[j]==0){
    				phi[i*p[j]]=(ll)phi[i]*p[j]%mod;
    				if((i/p[j])%p[j]==0) pm[i*p[j]]=(ll)pm[i]*p[j]%mod;
    				else pm[i*p[j]]=((ll)pm[i]*p[j]+pm[i/p[j]])%mod;
    				break;
    			}
    			phi[i*p[j]]=(ll)phi[i]*phi[p[j]]%mod;
    			pm[i*p[j]]=(ll)pm[i]*pm[p[j]]%mod;
    		}
    	}
    	for(int i=2;i<=N;i++)
    		(phi[i]+=phi[i-1])%=mod,(pm[i]+=pm[i-1])%=mod;
    }
    
    //Dusieve
    bool vis[iN+1];
    int duphi[iN+1],dupm[iN+1];
    int Phi(ll x){return x<=N?phi[x]:duphi[n/x];}
    int Pm(ll x){return x<=N?pm[x]:dupm[n/x];}
    void Dusieve(ll x){
    	if(x<=N||vis[n/x]) return;
    	vis[n/x]=true;
    	for(ll l=2,r;l<=x;l=r+1){
    		r=x/(x/l),Dusieve(x/l);
    		(duphi[n/x]-=(ll)(r-l+1)*Phi(x/l)%mod)%=mod;
    		(dupm[n/x]-=(ll)(r-l+1)*Pm(x/l)%mod)%=mod;
    	}
    	(duphi[n/x]+=(ll)x%mod*(x%mod+1)/2%mod)%=mod;
    	(dupm[n/x]+=duphi[n/x])%=mod;
    	(duphi[n/x]+=mod)%=mod,(dupm[n/x]+=mod)%=mod;
    }
    
    //Main
    int main(){
    	ios::sync_with_stdio(0);
    	cin.tie(0),cout.tie(0);
    	cin>>n;
    	Sieve(),Dusieve(n);
    	// cout<<Pm(n)<<'
    ';
    	for(ll l=1,r;l<=n;l=r+1){
    		r=n/(n/l);
    		(ans+=(ll)(n/l%mod)*(n/l%mod)%mod*(Pm(r)-Pm(l-1)+mod)%mod)%=mod;
    		/*
    		杜教筛是在开始整除分块前开始的,但是为什么这里可以直接Pm调用呢?
    		蒟蒻的回答:因为杜教筛内部处理了所有n的整除分块的答案。
    		*/
    	}
    	cout<<ans<<'
    ';
    	return 0;
    }
    

    祝大家学习愉快!

  • 相关阅读:
    软件补丁问题(网络流24题)
    飞行员配对方案问题(网络流24题)
    [NOIp普及组2011]瑞士轮
    如何在Linux上安装QQ
    [洛谷P2420] 让我们异或吧
    高斯消元
    [SCOI2014]方伯伯的OJ
    [USACO1.3]虫洞
    KMP算法讲解
    [洛谷P1382] 楼房
  • 原文地址:https://www.cnblogs.com/George1123/p/13307244.html
Copyright © 2020-2023  润新知