1. 程式人生 > >Loj 6053(EES篩法)

Loj 6053(EES篩法)

problem

思路

滿足EES篩法的要求 考慮f(p)f(p)p1p^1是啥,當p>2p>2時,f(p)=p1f(p)=p-1,而f(2)=3f(2)=3 所以預處理出p1,p0p^1,p^0的字首和即可 並注意,當i>=2i>=2時,這個字首和少計算了2(即p=2的貢獻3-1),所以加上一下

在dfs中維護時,沒什麼特別的,先計數第一類貢獻(當前數乘以更大的素數的貢獻),再計數第二類貢獻(最大的質因子有更大的指數)

注意別忘了取模

//phi(1^2)=1
//phi(p^2)=p^2-p
//phi((p^e)^2)=phi((p^{e-1})^2)*p*p; #include<bits/stdc++.h> using namespace std; typedef long long ll; const int mod=1e9+7; const int ni6=166666668; const int ni2=500000004; ll n,M; vector<int> pre[2],hou[2],primes; inline int add(const int x, const int v) { return x + v >= mod ? x + v -
mod : x + v; } inline int dec(const int x, const int v) { return x - v < 0 ? x - v + mod : x - v; } //這裡res是n/列舉的數 int dfs(ll res, int last, ll f){ //最大質因子是prime[last-1] 但將1放在外面值顯然一樣 int t=dec((res > M ? hou[1][n/res] : pre[1][res]),pre[1][primes[last]-1]); int ret= t*f%mod;//這裡需修改 for
(int i=last;i<(int) primes.size();++i){ int p = primes[i]; if((ll)p*p > res) break; int e=1; for(ll q=p,nres=res,nf=f*(p^e)%mod;q*p<=res;q*=p){//nf需修改 ret = add(ret,dfs(nres/=p,i+1,nf));//列舉更大的數 e++; nf = f* (p^e)%mod; ret =add(ret,nf);//指數大於1時,記上貢獻 } } return ret; } inline int ff(ll x){ x%=mod; return x*(x+1)%mod*ni2%mod; } int solve(ll n){ M=sqrt(n); for(int i=0;i<2;++i){ pre[i].clear();pre[i].resize(M+1); hou[i].clear();hou[i].resize(M+1); } primes.clear();primes.reserve(M+1); for(int i=1;i<=M;++i){ pre[0][i]=i-1; hou[0][i]=(n/i-1)%mod; pre[1][i]=dec(ff(i),1); hou[1][i]=dec(ff(n/i),1); } for(int p=2;p<=M;++p){ if(pre[0][p]==pre[0][p-1]) continue; primes.push_back(p); const ll q=(ll)p*p,m=n/p; const int pnt0=pre[0][p-1],pnt1=pre[1][p-1]; const int mid=M/p; const int End=min((ll)M,n/q); for(int i=1;i<=mid;++i){ hou[0][i]=dec(hou[0][i],dec(hou[0][i*p],pnt0)); hou[1][i]=dec(hou[1][i],dec(hou[1][i*p],pnt1)*(ll)p%mod); //hou[2][i]=dec(hou[2][i],dec(hou[2][i*p],pnt2)*q%mod); } for(int i=mid+1;i<=End;++i){ hou[0][i]=dec(hou[0][i],dec(pre[0][m/i],pnt0)); hou[1][i]=dec(hou[1][i],dec(pre[1][m/i],pnt1)*(ll)p%mod); //hou[2][i]=dec(hou[2][i],dec(pre[2][m/i],pnt2)*q%mod); } for(int i=M;i>=q;--i){ pre[0][i]=dec(pre[0][i],dec(pre[0][i/p],pnt0)); pre[1][i]=dec(pre[1][i],dec(pre[1][i/p],pnt1)*(ll)p%mod); //pre[2][i]=dec(pre[2][i],dec(pre[2][i/p],pnt2)*q%mod); } } primes.push_back(M+1); for (int i = 1; i <= M;++i) { pre[1][i] = i>=2 ? add(2,dec(pre[1][i], pre[0][i])) : (dec(pre[1][i], pre[0][i]));//p-1 hou[1][i] = n>=2*i ? add(2,dec(hou[1][i], hou[0][i])) :(dec(hou[1][i], hou[0][i])); } return n>1 ? add(dfs(n,0,1),1) : 1;//別忘了加上1的貢獻 } int main() { //freopen("in.txt","r",stdin); ios::sync_with_stdio(false); cin>>n; cout<<solve(n)<<endl; return 0; }