1. 程式人生 > >LibreOJ 6053 簡單的函式(Min25篩)

LibreOJ 6053 簡單的函式(Min25篩)

題意如題……

又是一個典型的題。 

模板中需要替換的主要有幾個地方。一是在初始化g(x,j)的時候,需要知道把所有出了1以外的正整數當作質數去計算的和是什麼,以及1到第j個質數對應的f(i)的和。對於第一個,根據函式的定義,所有數字異或1再求和相當於是所有數字直接求和,那麼答案就是n*(n+1)/2-1。對於第二個,由於質數只是需要算到sqrt(N)的級別,所以可以考慮初始篩素數的時候求一下。還有種方法,可以看到出了2以外,所以質數都是奇數,異或1之後會減一,所以可以考慮同時算一個質數的個數,然後求的時候減去之前的質數個數即可。這裡運用了第二種方法。具體見程式碼:

#include<bits/stdc++.h>
#define mod 1000000007
#define LL long long
#define pb push_back
#define lb lower_bound
#define ub upper_bound
#define INF 0x3f3f3f3f
#define sf(x) scanf("%lld",&x)
#define sc(x,y,z) scanf("%d%d%d",&x,&y,&z)
#define clr(x,n) memset(x,0,sizeof(x[0])*(n+5))
#define file(x) freopen(#x".in","r",stdin),freopen(#x".out","w",stdout)
using namespace std;

const int N = 1e6 + 10;
const int inv = 5e8 +4;

LL p[N],w[N],g[N],h[N],ps[N],M,n;
int block,id[N],sz;
bool isp[N];

void init(int n)
{
    sz=0;
    for(int i=2;i<=n;i++)
    {
        if(!isp[i])p[++sz]=i,ps[sz]=(ps[sz-1]+i)%mod;
        for(int j=1;j<=sz&&p[j]*i<n;j++)
        {
            isp[i*p[j]]=1;
            if(i%p[j]==0) break;
        }
    }
}

void sieve_g(LL n)
{
    M=0;
    for(LL i=1,last;i<=n;i=last+1)
    {
        LL len=n/i; last=n/len;
        w[++M]=len; g[M]=(w[M]-1)%mod;
        h[M]=w[M]%mod*((w[M]+1)%mod)%mod;
        h[M]=(h[M]*inv%mod-1+mod)%mod;
        if(len<=block) id[len]=M;
    }
    for(int i=1;i<=sz;i++)
        for(int j=1;j<=M&&p[i]*p[i]<=w[j];j++)
        {
            int op=w[j]/p[i]<=block?id[w[j]/p[i]]:n/(w[j]/p[i]);
            (g[j]-=g[op]-i+1)%=mod; (h[j]-=p[i]*(h[op]-ps[i-1])%mod)%=mod;
        }
}

LL S(LL x,LL y)
{
    LL k,res=0;
    if (x<=1||p[y]>x) return 0;
    if (x>block) k=n/x; else k=id[x];
    res=(h[k]-ps[y-1]-g[k]+y-1)%mod;
    if (y==1) res+=2;
    for(int i=y;i<=sz&&p[i]*p[i]<=x;i++)
        for(LL e=1,s=p[i]*p[i];s<=x;s*=p[i],e++)
            (res+=(p[i]^e)*S(x*p[i]/s,i+1)%mod+(p[i]^(e+1))%mod)%=mod;
    return res;
}

int main()
{
    sf(n);
    block=ceil(sqrt(n));
    init(n<=1e6?1e4:1e6); sieve_g(n);
    printf("%lld\n",(S(n,1)+1+mod)%mod);
    return 0;
}