题解 P3172 【[CQOI2015]选数】

· · 题解

很显然，题目让我们求的是(以下R=H看得顺眼点)

ans=\sum_{a_1=L}^R...\sum_{a_n=L}^R[gcd(a_1,...,a_n)==k] =\sum_{a_1=\lceil \frac L k \rceil}^{\lfloor \frac R k \rfloor}...\sum_{a_n=\lceil \frac L k \rceil}^{\lfloor \frac R k \rfloor}[gcd(a_1,...,a_n)==1]

接下来我们令l=\lceil \frac L k \rceil,r=\lfloor \frac R k \rfloor

=\sum_{a_1=l}^r...\sum_{a_n=l}^r[gcd(a_1,...,a_n)==1]

把最后一项莫比乌斯反演一下

=\sum_{a_1=l}^r...\sum_{a_n=l}^r \sum_{d|gcd(a_1,...,a_n)}\mu(d)

把枚举d提到最前面

=\sum_{d=1}^r\mu(d)\sum_{a_1=l}^r[d|a_1]...\sum_{a_n=l}^r[d|a_n]

辣么很显然，\sum_{a=l}^r[d|a]=[l,r]里d的倍数的个数

=\sum_{d=1}^r\mu(d)(\lfloor \frac r x \rfloor-\lfloor\frac {l-1} x \rfloor)^n

在整除分块一波就OK了

关于线性筛，它死了

这里我们用杜教筛，先上公式：

g(1)S(n)=\sum_{i=1}^{n}h(i)-\sum_{d=2}^{n}g(d)\cdot S(\lfloor\frac{n}{d}\rfloor)

这里我们去g=I,h=\mu *I=id那么式子变成

\to S(n)=1-\sum_{d=2}^{n}S(\lfloor\frac{n}{d}\rfloor)

用整除分块就可以做了

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int maxn=5e6+10;
bool isprime[maxn];
int miu[maxn],prime[maxn],tot,Sum[maxn];
const long long mod=1e9+7;
void init(){
    memset(isprime,0xff,sizeof isprime);
    tot=0;miu[1]=1;
    for(int i=2;i<maxn;i++){
        if(isprime[i])
            prime[++tot]=i,miu[i]=-1;
        for(int j=1;j<=tot&&i*prime[j]<maxn;j++){
            isprime[i*prime[j]]=false;
            if(i%prime[j])
                miu[i*prime[j]]=-miu[i];
            else{
                miu[i*prime[j]]=0;
                break;
            }
        }
    }
    Sum[0]=0;
    for(int i=1;i<maxn;i++)
        Sum[i]=(Sum[i-1]+miu[i])%mod;
}
map<int,int>S;
inline int calc(int x){
    if(x<maxn)return Sum[x];
    if(S[x])return S[x];
    int ans=1;
    for(int l=2,r;l<=x;l=r+1){
        r=x/(x/l);
        ans-=calc(x/l)*(r-l+1)%mod;
        ans%=mod;
    }
    return S[x]=ans%mod;
}
template <class T>inline void read(T &x){
    x=0;char c=getchar();
    for(;!isdigit(c);c=getchar());
    for(;isdigit(c);c=getchar())x=x*10+c-'0';
    return;
}
template <class T>void write(T x){
    if(x<0)putchar('-'),x=-x;
    if(x>=10)write(x/10);
    putchar(x%10+'0');
}
long long qpow(long long a,long long b){
    long long ret=1;
    while(b){
        if(b&1)ret=(ret*a)%mod;
        a=(a*a)%mod,b>>=1;
    }
    return ret;
}
void solve(long long N,long long K,long long L,long long R){
    L=(L-1)/K,R=R/K;
    long long ans=0;
    for(int l=1,r;l<=R;l=r+1){
        r=min(R/(R/l),(L/l)?L/(L/l):0x3f3f3f3f3f3f3f3f);
        ans=(ans+1ll*(calc(r)-calc(l-1))*qpow((R/l)-(L/l),N))%mod;
        //printf("%lld %lld %lld\n",l,r,ans);
    }
    write((ans+mod)%mod);putchar('\n');
}

signed main(){
    init();
    long long N,K,L,R;
    read(N);read(K);read(L);read(R);
    solve(N,K,L,R);
}