从[L,H]中可重复地选出N个数，问其gcd等于K的方案(mod 10^9+7)，1<=N,K<=109，1<=L<=H<=109，H-L<=10^5。
解
法一：直接
转换为Mobius问题，不难有
[ans=sum_{i_1=L}^Hsum_{i_2=L}^H...sum_{i_N=L}^H(gcd(i_1,i_2,...,i_N)==k)
]
于是设
[f(k)=sum_{i_1=L}^Hsum_{i_2=L}^H...sum_{i_n=L}^H(gcd(i_1,i_2,...,i_n)==k)
]
[F(k)=sum_{i_1=L}^Hsum_{i_2=L}^H...sum_{i_N=L}^H(k|gcd(i_1,i_2,...,i_N))=([frac{H}{k}]-[frac{L-1}{k}])^N
]
由Mobius反演定理代入有
[ans=sum_{k|d}([frac{H}{d}]-[frac{L-1}{d}])^Nmu(d/k)=
]
[sum_{d=1}^{[H/k]}([frac{H}{dk}]-[frac{L-1}{dk}])^Nmu(d)
]
显然左式可以整除分块，而至于后式只需要快速维护其前缀和，注意到N很大，于是用杜教筛优化即可。
参考代码：
#include <iostream>
#include <cstdio>
#define il inline
#define ri register
#define ll long long
#define control 6000000
#define yyb 1000000007
using namespace std;
template<class f1,class f2>
struct chain{
chain*next;f1 x;f2 y;
};
template<class f1,class f2,class f3,f3 key>
struct unordered_map{
chain<f1,f2>*head[key],*p;
il void insert(f1 x,f2 y){
p=new chain<f1,f2>,p->x=x,p->y=y;
x%=key,p->next=head[x],head[x]=p;
}
il chain<f1,f2>* find(f1 x){
for(p=head[x%key];p!=NULL;p=p->next)
if(p->x==x)return p;
return NULL;
}
};
bool check[control+1];
int prime[600000],pt,mb[control+1];
il void prepare(int);
il int djs(int),pow(int,int),min(int,int);
int main(){
int n,k,l,h,i,j,ans(0);
scanf("%d%d%d%d",&n,&k,&l,&h);
h/=k,(--l)/=k,prepare(control);
for(i=1;i<=h;i=j+1){
j=h/(h/i);if(l/i)j=min(j,l/(l/i));
(ans+=(ll)(djs(j)-djs(i-1))*pow(h/i-l/i,n)%yyb)%=yyb;
}printf("%d",(ans+yyb)%yyb);
return 0;
}
il int min(int a,int b){
return a<b?a:b;
}
il int pow(int x,int y){
int ans(1);while(y){
if(y&1)ans=(ll)ans*x%yyb;
x=(ll)x*x%yyb,y>>=1;
}return ans;
}
unordered_map<int,int,int,999983>m;
il int djs(int n){
if(n<=control)return mb[n];
if(m.find(n)!=NULL)return m.find(n)->y;int i,j,ans(1);
for(i=2;i<=n;i=j+1)j=n/(n/i),ans-=(j-i+1)*djs(n/j);
return m.insert(n,ans),ans;
}
il void prepare(int n){
int i,j;check[1]|=mb[1]|=true;
for(i=2;i<=n;++i){
if(!check[i])prime[++pt]=i,mb[i]=-1;
for(j=1;j<=pt&&prime[j]<=n/i;++j){
check[i*prime[j]]|=true;
if(!(i%prime[j]))break;
mb[i*prime[j]]=~mb[i]+1;
}
}for(i=1;i<=n;++i)mb[i]+=mb[i-1];
}
法二：优化
思路一：
[gcd(a,b)leq|a-b|(a
eq b)
]
证明：
(|x-y|=|k_1gcd(x,y)-k_2gcd(x,y)|=gcd(x,y)|k_1-k_2|)
故得证。
注意到(H-Lleq10^5),于是猜测gcd范围应在里面，于是我们有结论(gcd(a,b)leq|a-b|(a eq b)),所以只要对a,b相同的情况单独讨论，而显然这种方案有且仅有一种，只要看其是否落在([L,H])内即可，于是接下来就可以直接对法一进行优化了,但是注意该结论一定是不存在全相同对数。
思路二：
而注意到问题比较简单，而有可以使用Mobius反演，于是我们考虑容斥，于是令(L=[L-1]/k+1,H/=k),接下来问题变成只要对互质对数考虑即可，所以设(D(n)=([frac{H}{n}]-[frac{L-1}{n}])^N-[frac{H}{n}]+[frac{L-1}{n}])，表示有公约数n的互质的不全相同的对数，根据容斥原理，我们有
[ans=sum_{d=1}^{|H-L|}mu(d)D(d)
]
进行整除分块即可，注意后面全相同对数的情况，因为这也是一种方案。
参考代码：
#include <iostream>
#include <cstdio>
#define il inline
#define ri register
#define ll long long
#define yyb 1000000007
using namespace std;
bool check[100001];
int prime[2500],pt,mb[100001];
il void prepare(int);
il int pow(int,int),min(int,int);
int main(){
int n,k,l,h,l_,i,j,cjx,li,ans(0);
scanf("%d%d%d%d",&n,&k,&l,&h);
++(--l/=k),h/=k,l_=l-1;
li=h-l,prepare(li);
for(i=1;i<=li;i=cjx+1){
j=h/i-l_/i,cjx=h/(h/i);
if(l_/i)cjx=min(cjx,l_/(l_/i));
(ans+=(ll)(mb[cjx]-mb[i-1])*(pow(j,n)-j)%yyb)%=yyb;
}if(l==1)++ans;printf("%d",(ans+yyb)%yyb);
return 0;
}
il int min(int a,int b){
return a<b?a:b;
}
il int pow(int x,int y){
int ans(1);while(y){
if(y&1)ans=(ll)ans*x%yyb;
x=(ll)x*x%yyb,y>>=1;
}return ans;
}
il void prepare(int n){
int i,j;check[1]|=mb[1]|=true;
for(i=2;i<=n;++i){
if(!check[i])prime[++pt]=i,mb[i]=-1;
for(j=1;j<=pt&&prime[j]<=n/i;++j){
check[i*prime[j]]|=true;
if(!(i%prime[j]))break;
mb[i*prime[j]]=~mb[i]+1;
}
}for(i=1;i<=n;++i)mb[i]+=mb[i-1];
}