[CQOI2015]选数

丛雨 · 2021-01-09 17:49:09 · 题解

传送门

首先解读题意，我们得到下式：

令l=\lfloor\frac{H}{K}\rfloor,r=\lfloor\frac{L-1}{K}\rfloor

接着进行莫比乌斯反演

把d放在外层枚举

然后我们会发现上式可以用整除分块

但问题就来了，r(H)的范围是1e9

有的题解中使用了杜教筛，但事实上完全不需要。

注意到这个条件，同时我们又可以知道\lim\limits_{d\rightarrow\infty}(\lfloor\frac{r}{d}\rfloor-\lfloor\frac{l-1}{d}\rfloor)=0

所以当d足够大时可以不用计算\sum\mu(d)(并非单调)

于是我们设计出了以下算法：

先整除分块

• \lfloor\frac{r}{d}\rfloor\not=\lfloor\frac{l-1}{d}\rfloor
  • d\leq X$($X$为我们设置的闸量)时，通过预处理的$\sum\mu$得到值,预处理时间复杂度为$O(X)
  • d > X$枚举$i\in[l,r]$,通过筛出的$prime$得到$\mu(i)$,单次时间复杂度为$O(\frac{\sqrt{i}}{\ln i})$,当$X=10^7$时,期望出现次数$=(H-L+1)\int_{L}^{H}\frac{1}{x}dx=2

期望时间复杂度O(X+\sqrt H\log_2N)

#include<bits/stdc++.h>
using namespace std;
# define ll long long
# define read read1<ll>()
# define Type template<typename T>
Type T read1(){
    T t=0;
    char k;
    bool vis=0;
    do (k=getchar())=='-'&&(vis=1);while('0'>k||k>'9');
    while('0'<=k&&k<='9')t=(t<<3)+(t<<1)+(k^'0'),k=getchar();
    return vis?-t:t;
}
# define fre(k) freopen(k".in","r",stdin);freopen(k".out","w",stdout)
# define N 10000000
bool vis[10000005];
int mu[10000005],prime[1000005],s;
void init(){
    mu[1]=1;
    for(int i=2;i<=N;++i){
        if(!vis[i]){
            mu[i]=-1;
            prime[++prime[0]]=i;
        }
        for(int j=1;j<=prime[0]&&prime[j]*i<=N;++j){
            vis[prime[j]*i]=1;
            if(!(i%prime[j])){
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=2;i<=N;++i)mu[i]+=mu[i-1];
}
int Mu(int x){
    int t=1;
    for(int i=1;prime[i]*prime[i]<=x;++i)
        if(!(x%prime[i])){
            x/=prime[i];t=-t;
            if(!(x%prime[i]))return 0; 
        }
    if(x!=1)t=-t;
    return t;
}
# define mod 1000000007
ll qkpow(int n,int k){
    if(!k)return 1;
    ll t=qkpow(n,k>>1);
    t=t*t%mod;
    if(k&1)t=t*n%mod;
    return t;
}
int main(){
    init();
    int n=read,k=read,L=(read-1)/k,R=read/k;
    int ans=0;
    for(int l=1,r=0;l<=R;l=r+1){
        r=min(R/(R/l),l>L?1000000000:L/(L/l));
        if(R/l==L/l)continue;
        ll t=qkpow(R/l-L/l,n),v=0;
        if(r<=N)v=mu[r]-mu[l-1];
        else for(int i=l;i<=r;++i)v+=Mu(i);
        ans=(ans+t*v)%mod;
    }
    cout<<(ans+mod)%mod;
    return 0;
}