zoukankan      html  css  js  c++  java
  • 【51nod1847】 奇怪的数学题

    就当我是 A 了此题吧...

    首先预备知识有点多(因为题目要处理的东西都挺毒瘤):

    1. 杜教筛运用(当然你可以用其他筛?)
    2. 第二类斯特林数相关定理
    3. 下降阶乘幂相关定理
    4. min25 筛运用

    好了可以关掉本题解了

    咳咳上面我除了杜教筛都不会,于是熬了好久才 做掉的

    题意

    我们首先考虑题目让我们求什么:

    [ANS=sum_{i=1}^nsum_{j=1}^n sgcd(i,j) ]

    上面的 (sgcd(i,j)) 表示 (gcd(i,j)) 的次大约数

    emmm?有问题么?虽然题目让我们求 两个数的次大公约数 ,但是这个次大公约数必然能被他们的 gcd 整除,那么这个次大公约数也就是 gcd 的次大约数了

    noteskey

    Part 1

    我们考虑用 f(i) 表示一个数的次大约数:

    [ANS=sum_{i=1}^nsum_{j=1}^n f(gcd(i,j)) ]

    然后这里有 gcd 啊,立马开始化式子:

    [ANS=sum_{d=1}^nf(d)sum_{i=1}^{n/d}sum_{j=1}^{n/d} [gcd(i,j)=1] ]

    这时候我们就有两个选择了,要么反演变成 μ ,要么转成 φ ,个人感觉可能 φ 方便点(毕竟大家都写的 φ ),但是 μ 也没有差到哪里去

    求 φ

    [ANS=sum_{d=1}^nf(d)sum_{i=1}^{n/d}sum_{j=1}^{n/d} sum_[gcd(i,j)=1] ]

    [ANS=sum_{d=1}^nf(d)~ (2·sum_{i=1}^{n/d}φ(i)-1) ]

    然后我们发现这里可以数论分块,上杜教筛!

    求 μ

    [ANS=sum_{d=1}^nf(d)sum_{i=1}^{n/d}sum_{j=1}^{n/d} sum_[gcd(i,j)=1] ]

    [ANS=sum_{d=1}^nf(d) sum_{k=1}^{n/d}μ(k) sum_{i=1}^{n/(kd)}sum_{j=1}^{n/(kd)} ]

    [ANS=sum_{d=1}^nf(d) sum_{k=1}^{n/d}μ(k) {n over kd}^{~2} ]

    然后我们发现还是可以数论分块,上杜教筛!


    然后这篇题解到此结束,完结撒花

    有什么不对的样子...emmm 哦! 前面的 f 还没说怎么求啊!

    Part 2

    我们重新考虑 f 是个啥:

    如果说 x 是个质数,那么 f(x) 显然是 1 ,然后我们考虑一个合数的 f 值其实就是其本身除去最小值因子

    总的来说,一个数 x 的 f 值就是 x 除去其最小质因子(当然 1 是有特例,算 0 的)

    那么按照套路我们肯定是要求 f 的前缀和的(不然怎么数论分块啊!)

    所以我们肯定也是要求出类似 (sum_{i=1}^n i^k) 这样的东西(我们考虑除去一个最小值因子之后的贡献时肯定要用到的嘛)

    然后我们看看这玩意儿怎么求先

    k 次幂和的快速计算

    可能看到这玩意儿你会想到拉格朗日插值...

    但是很遗憾这道题的模数很不友好(就是那个自然溢出的玩意儿),我们只能另寻他路

    所以首先我们要用的就是第二类斯特林数的性质了:

    [n^k = sum_{i=1}^k left{ ~ ^k_i ight}i! (~ ^n_i) ]

    那么前缀和就是:

    [sum_{i=1}^n i^k = sum _{j=1}^n sum_{i=1}^k left{ ~ ^k_i ight}i! (~ ^j_i ) ]

    [sum_{i=1}^n i^k =sum_{i=1}^k left{ ~ ^k_i ight} sum _{j=1}^n i! ( ~^j_i) ]

    然后看到后面的东西其实是一个下降阶乘幂

    [sum_{i=1}^n i^k =sum_{i=1}^k left{ ~ ^k_i ight} sum _{j=1}^n j^underline{i} ]

    分析这个下降阶乘幂的前缀和:

    [sum _{j=1}^n j^underline{i}= {(n+1)^underline{i+1}over i+1} ]

    至于证明...要用积分(我们看到这玩意儿挺像多项式求积分的形式的,就是感觉底数 n+1 有点不对劲,不过毕竟这不是简单的多项式而是下降阶乘幂,所以就默认它存在吧)

    那么原来的式子带进去:

    [sum_{i=1}^n i^k =sum_{i=1}^k left{ ~ ^k_i ight} {(n+1)^underline{i+1}over i+1} ]

    这样的话我们就可以预处理斯特林数然后 (O(k)) 跑出 (i^k) 的前缀和了 (注意这里的 k 上限 50 ,所以不会出事)

    那么我们考虑质数 合数贡献分开来算,就是,就是说我们只要让合数的 f 值前缀和加上质数个数(每个 f 贡献为 1 ,相当于求质数个数)就是要求的 f 前缀和了

    质数个数可能还比较好做,但是还要考虑合数的次大约数前缀和怎么求:

    我们考虑之前说的一个数除去它的最小质因子,这有点像欧拉筛哈,但是欧拉筛是线性的,对于这样的问题无能为力

    于是我们考虑 min25 筛中其实也有每个数被最小质因子筛去的这一性质,于是我们欢乐地套上 min25 筛就过了此题

    code

    这里用的是第一种化简(也就是求 φ 的),但其实这题关键在于求出 f 的前缀和所以求 φ 还是 μ 无所谓了啦

    //by Judge
    #include<bits/stdc++.h>
    #define fp(i,a,b) for(register int i=(a),I=(b)+1;i<I;++i)
    #define ui unsigned int
    using namespace std;
    const int M=1e6+3;
    typedef ui arr[M<<1];
    ui n,m,mm,K,cnt,cntval,sq; ui ans,S[60][60];
    arr v,p,loc1,loc2,g,rec,phi,pcnt,pw,G; map<ui,ui> mp;
    inline ui qpow(ui x,int p){ ui s=1;
    	for(;p;p>>=1,x=x*x)
    		if(p&1) s=s*x; return s;
    }
    inline ui Pos(ui x){ //判断返回两种 id 中哪种的 val 
    	return x<=sq?loc1[x]:loc2[n/x];
    }
    inline ui f(ui x){  // 求 f 前缀和 
    	return G[Pos(x)]+pcnt[Pos(x)];
    }
    inline void prep(int n){ v[1]=phi[1]=1; //预处理各种变量 
    	for(int i=2;i<=n;++i){ if(!v[i]) p[++cnt]=i,phi[i]=i-1;
    		for(int j=1;j<=cnt&&i*p[j]<=n;++j){ v[i*p[j]]=1;
    			if(!(i%p[j])){phi[i*p[j]]=phi[i]*p[j];break;}
    			phi[i*p[j]]=phi[i]*(p[j]-1);
    		} phi[i]+=phi[i-1]; //线性筛 phi  
    	} S[0][0]=1;
    	fp(i,1,K){ S[i][1]=1;
    		fp(j,2,i) S[i][j]=S[i-1][j]*(ui)j+S[i-1][j-1]; //计算第二类斯特林数 
    	}
    	fp(i,1,cnt) pw[i]=qpow(p[i],K); //计算 p[i] 的 K 次幂 
    }
    inline void init_loc(ui n){ cntval=0; //预处理出要处理的 n/i 的 val 
    	for(ui l=1,r;l<=n;l=r+1)
    		r=n/(n/l),rec[++cntval]=n/l;
    	reverse(rec+1,rec+1+cntval);
    	fp(i,1,cntval)
    		if(rec[i]<=sq) loc1[rec[i]]=i;
    		else loc2[n/rec[i]]=i;
    }
    inline ui sum(ui l,ui r){ //计算 l~ r 的和 
    	ui tmp1=l+r,tmp2=r-l+1;
    	if(tmp1&1) tmp2>>=1;
    	else tmp1>>=1;
    	return tmp1*tmp2;
    }
    inline ui get_phi(ui n){ if(n<=mm) return phi[n]; //杜教筛求 phi 前缀和 
    	if(mp[n]) return mp[n]; ui ans=sum(1,n);
    	for(ui l=2,r;l<=n;l=r+1) r=n/(n/l),
    		ans-=get_phi(n/l)*(r-l+1); return mp[n]=ans;
    }
    inline ui calc_sumk(long long n){  //计算 n^K 的前缀和 
    	ui ans=0,tmpval,tmp;
    	fp(i,1,K){ // K 次幂 ,循环 K 次 
    		tmpval=1,tmp=(n-i+1)%(i+1);
    		for(long long j=n-i+1;j<=n+1;++j,++tmp){
    			if(tmp>=i+1) tmp-=i+1;
    			if(!tmp) tmpval*=(ui)j/(i+1);
    			else tmpval*=(ui)j; //计算下降幂 
    		}
    		ans+=S[K][i]*tmpval;
    	} return ans;
    }
    inline void get_g(){ // g 存 rec[i] 范围内 质数 ^ k 次幂 的前缀和 , pcnt 存 rec[i] 内质数个数(也就是 所有质数 f 的前缀和) 
    	fp(i,2,cntval) g[i]=calc_sumk(rec[i])-1,pcnt[i]=rec[i]-1;
    	fp(j,1,m) for(int i=cntval;i&&1ll*p[j]*p[j]<=rec[i];--i)
    		g[i]-=pw[j]*(g[Pos(rec[i]/p[j])]-g[Pos(p[j]-1)]),  // g 当前为最小质因子大于 p[j] 或者为质数 的数的 k 次幂前缀和 
    		pcnt[i]-=pcnt[Pos(rec[i]/p[j])]-pcnt[Pos(p[j]-1)],
    		G[i]+=g[Pos(rec[i]/p[j])]-g[Pos(p[j]-1)]; // G 中记录合数的 f 前缀和,因为是从大到小筛所以这里的 g 表示的是 g[rec[i]/p[j]][j-1] 
    }
    int main(){ cin>>n>>K;
    	prep(mm=pow(n,2/3.0)),sq=sqrt(n);
    	m=upper_bound(p+1,p+1+cnt,sq)-1-p;
    	init_loc(n),get_g(),ans=0; ui pre=0,now;
    	for(ui l=2,r;l<=n;l=r+1,pre=now) r=n/(n/l),
    		now=f(l),ans+=(get_phi(n/l)*2-1)*(now-pre);
    	return !printf("%u
    ",ans);
    }
    
  • 相关阅读:
    Atitit.播放系统规划新版本 v4 q18 and 最近版本回顾
    Atitit.播放系统规划新版本 v4 q18 and 最近版本回顾
    atitit.极光消息推送服务器端开发实现推送  jpush v3. 总结o7p
    atitit.极光消息推送服务器端开发实现推送  jpush v3. 总结o7p
    Atitit.文件搜索工具 attilax 总结
    Atitit.文件搜索工具 attilax 总结
    Atitit.软件命名空间  包的命名统计 及命名表(2000个名称) 方案java package
    Atitit.软件命名空间  包的命名统计 及命名表(2000个名称) 方案java package
    Atitit..状态机与词法分析  通用分词器 分词引擎的设计与实现 attilax总结
    Atitit..状态机与词法分析  通用分词器 分词引擎的设计与实现 attilax总结
  • 原文地址:https://www.cnblogs.com/Judge/p/10620755.html
Copyright © 2011-2022 走看看