【学术篇】51nod 1238 最小公倍数之和

这是一道杜教筛的入(du)门(liu)题目...


题目大意

ni=1nj=1lcm(i,j)ni=1nj=1lcm(i,j)

一看就是辣鸡反演一类的题目, 那就化式子呗..

ni=1nj=1lcm(i,j)=ni=1nj=1ijgcd(i,j)=ni=1nj=1nk=1ijk[gcd(i,j)=k]=nk=1nki=1nkj=1ijk[gcd(i,j)=1]=nk=1knki=1inkj=1j[gcd(i,j)=1]=nk=1knki=1inkj=1jnkd=1[d|i][d|j]μ(d)=nk=1knkd=1μ(d)nki=1[d|i]inkj=1[d|j]ji=dq,j=dq,=nk=1knkd=1μ(d)nkp=1dpnkq=1dq=nk=1knkd=1μ(d)d2nkp=1pnkq=1qni=1i=n(n+1)2,=nk=1knkd=1μ(d)d2(nk(nk+1)2)2t=kd,=nt=1(nt(nt+1)2)2d|tμ(d)tdg(x)=d|xμ(d)xd,=nt=1(nt(nt+1)2)2g(t)

然后前面这一堆我们可以分块O(1)搞出来, 所以问题的关键就变成了g(x)的前缀和.
因为题目中n1010, 所以这里要用O(n23)的杜教筛.
但是杜教筛的框架我们是知道的,

f(x),sf(x)g(x)h(x),使(fg)(x)=h(x),h(x)sh(x),sf(x)=sh(x)nd=2sf(nd)g(d)g(1)

所以我们的任务又变成了找到合适的g(x)h(x).
在某篇歪果仁的blog上, 他说:

With a little luck, we can notice that Id(l) = g * Id2(l),

  • 注: id(x)=x(以下称为n)id2(x)=x2(以下称为n2)

But 这运气也太好了点吧? 我怎么就找不到这种函数呢?!
不过我们还是要考虑一下的... (这波化式子的时候注意乘号和卷积号)

g(x)=d|xμ(d)xd=d|xd2μ(d)xd=(n2μn)(x)(gn2)(x)=((n2μn)n2)(x),,.x2,( )=((n2μn2)n)(x)=(n2(μ1)n)(x)=(n2ϵn)(x)=n(x)

所以我们就证明出来了(gn2)(x)=n(x)
这样就可以直接杜教筛了= =
然后就是注意细节了OvO
比如非常坑人的1010(1e9+7)会爆long long!!
所以要开unsigned...
就这样吧= =
代码:

#include <map>
#include <cstring>
#include <iostream>
using namespace std;
typedef unsigned long long LL;
#define ri register int
#define rll register unsigned long long 
const int p=1e9+7;
int qpow(int a,int b,int s=1){
	for(;b;b>>=1,a=1LL*a*a%p)
		if(b&1) s=1LL*s*a%p;
	return s;
}
const int inv2=qpow(2,p-2),inv4=qpow(4,p-2),inv6=qpow(6,p-2);
const int N=4800000,_N=N+10;
const int SQ=233333;
LL n;
map<LL,int> mp;
map<LL,int>::iterator it;
int g[_N],pr[N/2],tot;
bool notp[_N];
void shai(){
	g[1]=1; ri i,j; rll k;
	for(i=2;i<=N;++i){
		if(!notp[i])
			pr[++tot]=i,g[i]=(-1LL*i*(i-1)%p+p)%p;
		for(j=1;j<=tot&&(k=1LL*i*pr[j])<=N;++j){
			notp[k]=1;
			if(i%pr[j]==0){g[k]=1LL*g[i]*pr[j]%p;break;}
			g[k]=1LL*g[i]*g[pr[j]]%p;
		}
	}
	for(i=2;i<=N;++i) g[i]=(g[i]+g[i-1])%p;
}
inline int calcsqr(LL x){x%=p;return x*(x+1)%p*(2*x+1)%p*inv6%p;}
inline int calcsum(LL x){x%=p;return x*x%p*(x+1)%p*(x+1)%p*inv4%p;}
int calc(LL x){
	if(x<=N) return g[x];
	it=mp.find(x);
	if(it!=mp.end()) return it->second;
	int ans=1LL*x%p*(x+1)%p*inv2%p;
	LL last=0;
	for(rll i=2;i<=x;i=last+1){
		last=x/(x/i);
		ans=(ans-1LL*calc(x/i)*(calcsqr(last)-calcsqr(i-1)+p)%p+p)%p;
	}
	return mp[x]=(ans%p+p)%p;
}
int main(){
	shai();cin>>n;
	LL last=0; int ans=0;
	for(rll i=1;i<=n;i=last+1){
		last=n/(n/i);
		ans=(ans+1LL*calcsum(n/i)*(calc(last)-calc(i-1)+p)%p)%p;
	}
	cout<<ans;
}

这破题我连推带写加总结写了整整一天!!!
窝还是太弱了..

posted @   Enzymii  阅读(256)  评论(0编辑  收藏  举报
编辑推荐:
· 浏览器原生「磁吸」效果!Anchor Positioning 锚点定位神器解析
· 没有源码,如何修改代码逻辑?
· 一个奇形怪状的面试题:Bean中的CHM要不要加volatile?
· [.NET]调用本地 Deepseek 模型
· 一个费力不讨好的项目,让我损失了近一半的绩效!
阅读排行:
· 微软正式发布.NET 10 Preview 1:开启下一代开发框架新篇章
· 没有源码,如何修改代码逻辑?
· PowerShell开发游戏 · 打蜜蜂
· 在鹅厂做java开发是什么体验
· WPF到Web的无缝过渡:英雄联盟客户端的OpenSilver迁移实战
点击右上角即可分享
微信分享提示