• [洛谷P4213]【模板】杜教筛(Sum)(杜教筛讲解)


    题面

    https://www.luogu.com.cn/problem/P4213

    题解

    前置知识

    杜教筛

    对于一些数论函数,杜教筛可以在低于线性的时间复杂度内求出其前缀和。

    适用范围:对于数论函数f,g,h,如果g、h的前缀和能够(O(1))求解,并且有

    [h=f{ imes}g$$(其中${ imes}$表示狄利克雷卷积) 那么,我们可以在低于线性的时间复杂度内求出f的前缀和。原理如下:(以下sf,sg,sh分别表示f,g,h的前缀和) $$h(n)={sumlimits_{d1d2=n}}f(d1)g(d2)]

    [{sumlimits_{i=1}^{n}}h(i)={sumlimits_{d1d2{leq}n}}f(d1)g(d2) ]

    [sh(i)={sumlimits_{d2=1}^{n}}g(d2){sumlimits_{d1=1}^{{lfloor}{frac{n}{d1}}{ floor}}}f(d1) ]

    [={sumlimits_{d=1}^{n}}g(d)sf({lfloor}{frac{n}{d}}{ floor}) ]

    将d=1一项提出,得

    [g(1)sf(n)=sh(n)-{sumlimits_{i=2}^{n}g(i)sf({lfloor}{frac{n}{i}}{ floor})} ]

    为了方便,不妨设(g(1)=1)。那么有

    [sf(n)=sh(n)-{sumlimits_{i=2}^{n}}g(i)sf({lfloor}{frac{n}{i}}{ floor}) ]

    即可递归计算了。使用数论分块可以将计算(sf(n))的时间降到(O({sqrt{n}}))(不算递归下去的时间);并且,递归计算的所有子问题sf(n'),一定可以写成(n'={lfloor}{frac{n}{t}}{ floor})的形式,其中t是某个正整数。所以根据数论分块,记忆化后所有需要计算sf的下标只有({lfloor}{frac{n}{i}}{ floor}(1{leq}i{leq}n))去重后的那些数,计算的总时间是它们的平方根之和,也就是

    [{sumlimits_{i=1}^{{sqrt{n}}}}{sqrt{i}}+{sumlimits_{i=1}^{sqrt{n}}}{sqrt{frac{n}{i}}} ]

    两边分别做积分近似,可以得到这样的时间复杂度为(O(n^{frac{3}{4}}))

    如果函数f是积性函数,那么我们还可以进行进一步的优化:设置阈值M,将1~M的f值线性筛出来,如果计算sf(n)时(n{leq}M)就直接调用这些值,否则再进行记忆化搜索。

    考虑M的最优值。如果(M<{sqrt{n}}),总时间为

    [O(M)+{sumlimits_{i=M+1}^{sqrt{n}}}+{sumlimits_{i=1}^{sqrt{n}}}{sqrt{frac{n}{i}}} ]

    发现第三项还是(O(n^{frac{3}{4}})),所以这样不行。

    如果(M{geq}{sqrt{n}}),总时间为

    [O(M)+{sumlimits_{i=1}^{frac{n}{M}}}{sqrt{frac{n}{i}}} ]

    积分近似得第二项是(O({frac{n}{sqrt{m}}})),所以置(M=n^{frac{2}{3}})的时候,总时间复杂度最小,为(O(n^{frac{2}{3}}))

    代码

    #include<bits/stdc++.h>
    
    using namespace std;
    
    #define N 2147483647
    #define M 5000000
    #define T 1300  
    #define ll long long
    #define rg register
    
    inline ll read(){
    	ll s = 0,ww = 1;
    	char ch = getchar();
    	while(ch < '0' || ch > '9'){if(ch == '-')ww = -1;ch = getchar();}
    	while('0' <= ch && ch <= '9'){s = 10 * s + ch - '0';ch = getchar();}
    	return s * ww;
    }
    
    inline void write(ll x){
    	if(x < 0)putchar('-'),x = -x;
    	if(x > 9)write(x / 10);
    	putchar('0' + x % 10);
    }
    
    ll pri[M/5+5];
    ll pn;
    bool isp[M+5];
    ll smu[M+5],sphi[M+5];
    
    inline void Eular(){
    	smu[1] = sphi[1] = 1;
    	for(rg ll i = 2;i <= M;i++)isp[i] = 1;
    	for(rg ll i = 2;i <= M;i++){
    		if(isp[i])pri[++pn] = i,smu[i] = -1,sphi[i] = i - 1;
    		for(rg ll j = 1;i * pri[j] <= M;j++){
    			isp[i*pri[j]] = 0;
    			if(i % pri[j])smu[i*pri[j]] = -smu[i],sphi[i*pri[j]] = sphi[i] * sphi[pri[j]];
    			else{
    				smu[i*pri[j]] = 0;
    				sphi[i*pri[j]] = sphi[i] * pri[j];
    				break;
    			}
    		}
    	}
    	for(rg ll i = 2;i <= M;i++)smu[i] += smu[i-1],sphi[i] += sphi[i-1]; 
    }
    
    ll n;
    
    unordered_map<ll,ll>Hmu,Hphi;
    
    inline ll sumphi(ll x){
    	if(x <= M)return sphi[x]; //小于M的部分,预处理出来,使用数组存储 
    	if(Hphi.count(x))return Hphi[x]; //大于M的部分,记忆化搜索,使用Hash表存储 
    	ll ans = 1ll * (ll)x * ((ll)x + 1) >> 1;
    	ll L,R = 1;
    	while(R < x){
    		L = R + 1,R = x / (x / L);
    		ans -= 1ll * (R - L + 1) * sumphi(x / L);
    	}
    	return Hphi[x] = ans;
    }
    
    inline ll summu(ll x){
    	if(x <= M)return smu[x]; //同上
    	if(Hmu.count(x))return Hmu[x]; //同上
    	ll ans = 1;
    	ll L,R = 1;
    	while(R < x){
    		L = R + 1,R = x / (x / L);
    		ans -= 1ll * (R - L + 1) * summu(x / L);
    	}
    	return Hmu[x] = ans;
    }
    
    int main(){
    	Eular();
    	ll c = read();
    	while(c--){
    		n = read();
    		write(sumphi(n)),putchar(' '),write(summu(n)),putchar('
    ');
    	}	
    	return 0;
    }
    
    
  • 相关阅读:
    2017 ICPC beijing E
    1629 B君的圆锥
    1298 圆与三角形
    通过String获取字符数组
    Java中的代码点与代码单元
    数据库事务隔离级别
    oracle修改密码、添加用户及授权
    Python起航
    软件测试常见概念
    TestNG--@Factory
  • 原文地址:https://www.cnblogs.com/xh092113/p/12291710.html
Copyright © 2020-2023  润新知