【luogu3768】简单的数学题【杜教筛】【欧拉函数】
题目传送门
题意
求 ,,为质数。
题解
之前一直想着用函数搞,结果死都只能搞到。后来,我fa♂现了bzoj3518点组计数,这道题用了欧拉函数的一个性质,好像叫做欧拉反演?
于是新日暮里的大门打开了。
于是开始推式子。
我们令
则原式
于是我们只需要能够快速求出就可以分块求和计算了。
这不就是杜教筛干的活吗?
因为
所以
因此原式
把带人
得
然后就递归计算即可。细节详见代码。
这道题就解决了!
代码
#include<cstdio>
typedef long long ll;
const ll N=5000005,M=1000017;
ll mod,n,ans,p[N],phi[N];
bool vis[N];
struct hashmap{
int cnt,head[M],nxt[M*10];
ll to[M*10],v[M*10];
bool count(ll f){
ll x=f%M;
for(int i=head[x];i;i=nxt[i]){
if(to[i]==f){
return true;
}
}
return false;
}
ll get(ll f){
ll x=f%M;
for(int i=head[x];i;i=nxt[i]){
if(to[i]==f){
return v[i];
}
}
return 0;
}
void set(ll f,ll val){
ll x=f%M;
to[++cnt]=f;
nxt[cnt]=head[x];
v[cnt]=val;
head[x]=cnt;
}
}mp;
inline ll sqr(ll x){
x%=mod;
return x*x%mod;
}
ll calc(ll x){
if(!x){
return 0;
}
x%=mod;
return sqr((1+x)*x/2)%mod;
}
ll get(ll n){
n%=mod;
if(n%3!=1){
return n*(n+1)/6%mod*(2*n+1)%mod;
}else if(n&1){
return (n+1)*(2*n+1)/6%mod*n%mod;
}else{
return n*(2*n+1)/6%mod*(n+1)%mod;
}
}
ll solve(ll n){
if(n<=5000000){
return phi[n];
}
if(mp.count(n)){
return mp.get(n);
}
ll res=calc(n);
for(ll i=2,last;i<=n;i=last+1){
last=n/(n/i);
res-=solve(n/i)*(get(last)-get(i-1))%mod;
res%=mod;
}
mp.set(n,res);
return res;
}
int main(){
scanf("%lld%lld",&mod,&n);
phi[1]=1;
for(int i=2;i<=5000000;i++){
if(!vis[i]){
p[++p[0]]=i;
phi[i]=i-1;
}
for(int j=1;j<=p[0]&&i*p[j]<=5000000;j++){
vis[i*p[j]]=true;
if(i%p[j]){
phi[i*p[j]]=phi[i]*(p[j]-1);
}else{
phi[i*p[j]]=phi[i]*p[j];
}
}
}
for(int i=1;i<=5000000;i++){
phi[i]=phi[i]*i%mod*i%mod;
phi[i]+=phi[i-1];
phi[i]%=mod;
}
for(ll i=1,last;i<=n;i=last+1){
last=n/(n/i);
ans+=(solve(last)-solve(i-1))*calc(n/i)%mod;
ans%=mod;
}
printf("%lld\n",(ans+mod)%mod);
return 0;
}
最后附上蒟蒻博主龙hou飞yan凤wu舞chi的公式草稿一张2333 鼠标写字真是太方便啦