acm-(数论、莫比乌斯反演、杜教筛)51nod 1227 平均最小公倍数

最新推荐文章于 2024-09-12 19:02:20 发布

&*^*&

最新推荐文章于 2024-09-12 19:02:20 发布

阅读量137

点赞数

文章标签： acm竞赛算法

本文链接：https://blog.csdn.net/qq_30361651/article/details/109184626

版权

数论同时被 3 个专栏收录

11 篇文章 0 订阅

订阅专栏

莫比乌斯反演

9 篇文章 0 订阅

订阅专栏

杜教筛

3 篇文章 0 订阅

订阅专栏

传送门
设函数 $F(n)=\sum_{k=1}^n\frac 1k\sum_{i=1}^klcm(i,k)$ ，题目就是让求 $F (b) - F (a - 1)$ ，于是对 $F (n)$ 先套路反演一波，得到 $F(n)=\frac{n+\sum_{i=1}^nf(i)}2$ ( $f(n)=\sum_{d\mid n}d\varphi(d)$ )，设 $S(n)=\sum_{i=1}^nf(i)$ ，看上去就像杜教筛的样子，由于 $f=(id\cdot \varphi)*I$ ，故我们可以设 $g=f*id=((\varphi*id)*id)*I=\sigma_2$ ，又由于 $\sum_{i=1}^ng(i)=\sum_{i=1}^ni^2\lfloor \frac ni\rfloor$ ，显然可以分块 $O(\sqrt n)$ 求解，不过我们可以预处理出 $g$ 函数的 $n^{\frac 23}$ 的前缀和，这样能加快速度，然后在此基础上套一个杜教筛式子，我们有 $S(n)=\sum_{i=1}^ng(i)-\sum_{i=2}^niS(\lfloor\frac ni\rfloor)$ ，由于计算 $\sum_{i=1}^ng(i)$ 比较快，总复杂度应该还是 $O(n^{\frac 23})$ ，注意 $f$ 函数也要预处理一下前缀和，否则杜教筛复杂度会退化。

int nn,N,limt,prim[maxn],tot,s[maxn],e[maxn],f[maxn],g[maxm],dp[maxm],ss[maxn],
	dps[maxm],es[maxn],fs[maxn];
bool flag[maxn];
int id(int x){
	return x<=limt?x:nn/x+limt;
}
void init(int n){
	s[1]=1,ss[1]=1;
	FOR(i,2,n+1){
		if(!flag[i])prim[tot++]=i,e[i]=1ll*i*(i-1)%mod,s[i]=e[i]+1,f[i]=1,
					es[i]=1ll*i*i%mod,ss[i]=es[i]+1,fs[i]=1;
		for(register int j=0;j<tot && prim[j]*i<=n;++j){
			flag[i*prim[j]]=1;
			if(i%prim[j]==0){
				e[i*prim[j]]=1ll*e[i]*sqr(prim[j])%mod;
				f[i*prim[j]]=f[i];
				s[i*prim[j]]=(s[i]+1ll*e[i*prim[j]]*f[i]%mod)%mod;
				es[i*prim[j]]=1ll*es[i]*sqr(prim[j])%mod;
				fs[i*prim[j]]=fs[i];
				ss[i*prim[j]]=(ss[i]+1ll*es[i*prim[j]]*fs[i]%mod)%mod;
				break;
			}
			e[i*prim[j]]=1ll*prim[j]*(prim[j]-1)%mod;
			f[i*prim[j]]=s[i];
			s[i*prim[j]]=1ll*s[i]*(e[i*prim[j]]+1)%mod;
			es[i*prim[j]]=sqr(prim[j]);
			fs[i*prim[j]]=ss[i];
			ss[i*prim[j]]=1ll*ss[i]*(es[i*prim[j]]+1)%mod;
		}
	}
	FOR(i,1,n+1)add(s[i],s[i-1]),add(ss[i],ss[i-1]);
}
int SS(int n){
	if(n<=N)return ss[n];
	register int idd=id(n),i=1,j,ans=0,pre=0,cur=0;
	if(dps[idd])return dps[idd];
	while(i<=n){
		j=n/(n/i);
		add(ans,1ll*((cur=sm2(j))-pre+mod)*(n/i)%mod);
		pre=cur;
		i=j+1;
	}
	return dps[idd]=ans;
}
int S(int n){
	if(n<=N)return s[n];
	register int idd=id(n),ans=SS(n),i=2,j,cur=0,pre=1;
	if(dp[idd])return dp[idd];
	while(i<=n){
		j=n/(n/i);
		dec(ans,1ll*((cur=sm1(j))-pre+mod)*S(n/i)%mod);
		i=j+1;
		pre=cur;
	}
	return dp[idd]=ans;
}
int cal(int n){
	if(n<=0)return 0;
	nn=n;
	limt=sqrt(nn);
	register int i =1,j,ans=0,c=nn/limt+limt;
	FOR(i,1,c+1)dps[i]=dp[i]=0;
	return 1ll*(S(n)+n)*qpow(2,mod-2,mod)%mod;
}
int main(){
	int a,b;
	rd(&a,&b);
	N=pow(b,2.0/3.0);
	init(N=max(N,10000));
	wrn((cal(b)-cal(a-1)+mod)%mod);
}