题意:
给你
n
,
k
,
l
,
r
n,k,l,r
n,k,l,r,让你从
[
l
,
r
]
[l,r]
[l,r]选一个数,选
n
n
n次,总方案数是
(
r
−
l
+
1
)
n
(r-l+1)^n
(r−l+1)n,问选出的
n
n
n个数的gcd恰好是
k
k
k的方案数。
n
,
k
,
l
,
r
<
=
1
e
9
,
r
−
l
<
=
1
e
5
n,k,l,r<=1e9,r-l<=1e5
n,k,l,r<=1e9,r−l<=1e5
题解:
考虑反演。
我们设
f
(
x
)
f(x)
f(x)表示gcd是
x
x
x的倍数的方案数,设
g
(
x
)
g(x)
g(x)表示gcd是
x
x
x的方案数。那么我们有
f
(
n
)
=
∑
n
∣
d
r
g
(
d
)
=
(
⌊
r
k
⌋
−
⌊
l
−
1
k
⌋
)
n
f(n)=\sum_{n|d}^rg(d)=(\lfloor\frac{r}{k}\rfloor-\lfloor\frac{l-1}{k}\rfloor)^n
f(n)=n∣d∑rg(d)=(⌊kr⌋−⌊kl−1⌋)n 其中
⌊
r
k
⌋
−
⌊
l
−
1
k
⌋
\lfloor\frac{r}{k}\rfloor-\lfloor\frac{l-1}{k}\rfloor
⌊kr⌋−⌊kl−1⌋表示在
[
l
,
r
]
[l,r]
[l,r]之间有多少个数是
k
k
k的倍数,只有从这些数中选最后才看
g
c
d
gcd
gcd是
k
k
k或
k
k
k的倍数,而一次有这么多种情况,一共要选
n
n
n次,于是就是
n
n
n次方。然后接下来开始莫比乌斯反演。
g
(
n
)
=
∑
n
∣
d
f
(
d
)
μ
(
⌊
d
n
⌋
)
g(n)=\sum_{n|d}f(d)\mu(\lfloor\frac{d}{n}\rfloor)
g(n)=n∣d∑f(d)μ(⌊nd⌋)我们更换枚举对象,改为枚举
n
n
n的倍数,那么有
g
(
n
)
=
∑
d
=
1
⌊
r
n
⌋
f
(
d
n
)
μ
(
d
)
g(n)=\sum_{d=1}^{\lfloor\frac{r}{n}\rfloor}f(dn)\mu(d)
g(n)=d=1∑⌊nr⌋f(dn)μ(d)我们要求的是
g
(
k
)
g(k)
g(k)
g
(
k
)
=
∑
d
=
1
⌊
r
k
⌋
μ
(
d
)
(
⌊
r
d
k
⌋
−
⌊
l
−
1
d
k
⌋
)
n
g(k)=\sum_{d=1}^{\lfloor\frac{r}{k}\rfloor}\mu(d)(\lfloor\frac{r}{dk}\rfloor-\lfloor\frac{l-1}{dk}\rfloor)^n
g(k)=d=1∑⌊kr⌋μ(d)(⌊dkr⌋−⌊dkl−1⌋)n
对于这个式子,对后半部分整除分块+快速幂,同时用杜教筛处理
μ
\mu
μ的前缀和,然后就可以了。
代码:
#include <bits/stdc++.h>
using namespace std;
const long long mod=1e9+7;
long long n,k,l,r;
int vis[1000010],p[1000010],cnt,mu[1000010];
long long sum[1000010],ans;
map<long long,long long> mp;
inline long long ksm(long long x,long long y)
{
long long res=1;
while(y)
{
if(y&1)
res=res*x%mod;
x=x*x%mod;
y>>=1;
}
return res;
}
inline long long calc(long long x)
{
if(x<=1000000)
return sum[x];
if(mp[x])
return mp[x];
long long res=1,j;
for(long long i=2;i<=x;i=j+1)
{
j=x/(x/i);
res=(res-(j-i+1)*calc(x/i)%mod+mod)%mod;
}
mp[x]=res;
return res;
}
int main()
{
scanf("%lld%lld%lld%lld",&n,&k,&l,&r);
mu[1]=1;
sum[1]=1;
for(int i=2;i<=1000000;++i)
{
if(!vis[i])
{
p[++cnt]=i;
mu[i]=-1;
}
for(int j=1;j<=cnt&&p[j]*i<=1000000;++j)
{
vis[i*p[j]]=1;
if(i%p[j]==0)
{
mu[i*p[j]]=0;
break;
}
mu[i*p[j]]=-mu[i];
}
sum[i]=sum[i-1]+mu[i];
}
long long j,ji=r/k,jj=(l-1)/k;
for(long long i=1;i<=ji;i=j+1)
{
if(jj/i)
j=min(jj/(jj/i),ji/(ji/i));
else
j=ji/(ji/i);
ans=(ans+(calc(j)-calc(i-1)+mod)%mod*ksm(ji/i-jj/i,n)%mod)%mod;
}
printf("%lld\n",ans);
return 0;
}