1. 程式人生 > >LOJ 6053 簡單的函式 - Min_25篩

LOJ 6053 簡單的函式 - Min_25篩

不知道為啥腦子一抽就來溫習了一波Min_25篩
這題F( p )=p-1(除了p=2,這個特判一下即可),最後把2加回來。

#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
#define Rep(i,v) rep(i,0,(int)v.size()-1)
#define lint long long
#define mod 1000000007
#define ull unsigned lint
#define db long double
#define pb push_back
#define mp make_pair
#define fir first #define sec second #define gc getchar() #define debug(x) cerr<<#x<<"="<<x #define sp <<" " #define ln <<endl using namespace std; typedef pair<int,int> pii; typedef set<int>::iterator sit; #define P(x) (x>=mod?x-=mod:0) #define Q(x) (x<0?x+=mod:0)
inline int fast_pow(int x,int k,int ans=1) { for(;k;k>>=1,x=(lint)x*x%mod) (k&1)?ans=(lint)ans*x%mod:0;return ans; } namespace min_25{ const int MXS=100010;lint n; int s,A[MXS],B[MXS],pri[MXS],cnt; int a[MXS],b[MXS],mip[MXS],np[MXS]; inline int MySqrt(lint n) { int x=(int)sqrt(n); while
((lint)x*x<n) x++; while((lint)x*x>n) x--; return x; } inline int init() { memset(np,0,sizeof(int)*(s+1)); memset(A,0,sizeof(int)*(s+1)); memset(B,0,sizeof(int)*(s+1)); return cnt=0; } inline int prelude() { rep(i,2,s) if(!np[i]) rep(j,i,s/i) np[i*j]=1; rep(i,2,s) if(!np[i]) pri[++cnt]=i;return pri[++cnt]=s+1; } inline int getg(int xs) { // rep(i,1,s) debug(i)sp,debug(a[i])sp,debug(b[i])ln;cerr ln; rep(p,2,s) if(!np[p]) { lint r=(lint)p*p;int ed1=s/p,ed2=(int)min(n/r,(lint)s); rep(i,1,ed1) b[i]-=(lint)mip[p]*(b[i*p]-a[p-1])%mod,P(b[i]),Q(b[i]); rep(i,ed1+1,ed2) b[i]-=(lint)mip[p]*(a[n/i/p]-a[p-1])%mod,P(b[i]),Q(b[i]); for(int i=s;i>=r;i--) a[i]-=(lint)mip[p]*(a[i/p]-a[p-1])%mod,P(a[i]),Q(a[i]); // debug(p)ln;rep(i,1,s) debug(i)sp,debug(a[i])sp,debug(b[i])ln;cerr ln; } rep(i,1,s) A[i]=(A[i]+(lint)a[i]*xs)%mod,B[i]=(B[i]+(lint)b[i]*xs)%mod;return 0; } inline int F(int p,int c) { return p^c; } inline int G(lint x) { return x<=s?A[x]:B[n/x]; } int S(lint x,int y,int curd=0) { // debug(x)sp,debug(y)sp,debug(pri[y])ln; if(pri[y]>x) return 0; int res=G(x)-G(pri[y]-1);Q(res); rep(i,y,cnt&&(lint)pri[i]*pri[i]<=x) { lint v=pri[i]; for(int j=1;v<=x/pri[i];v*=pri[i],j++) res=(res+(lint)F(pri[i],j)*S(x/v,i+1,curd+1)+F(pri[i],j+1))%mod; // debug(pri[i])sp,debug(j)sp,debug(j+1)ln, // debug(qwq)sp,debug(F(pri[i],j))sp,debug(F(pri[i],j+1))ln ln; } // debug(x)sp,debug(pri[y])sp,debug(res)ln ln; return res; } inline int solve(lint _n) { s=MySqrt(n=_n),init(),prelude(); rep(i,2,s) if(!np[i]) mip[i]=1; rep(i,1,s) a[i]=i,b[i]=(n/i)%mod;getg(mod-1); // rep(i,1,s) debug(i)sp,debug(A[i])sp,debug(B[i])ln; rep(i,2,s) if(!np[i]) mip[i]=i;int inv2=fast_pow(2,mod-2),x; rep(i,1,s) a[i]=(i*(i+1ll)/2)%mod,x=(n/i)%mod,b[i]=x*(x+1ll)%mod*inv2%mod; getg(1); // rep(i,1,s) debug(i)sp,debug(A[i])sp,debug(B[i])ln;cerr ln; int Ans=S(n,1)+1+(n>=2?2:0);return P(Ans),Q(Ans),Ans; } } int main() { lint n;cin>>n;return printf("%d\n",min_25::solve(n)); }