luoguP3172 [CQOI2015]选数

题意

所求即为：
\(\sum\limits_{i_1=L}^{R}\sum\limits_{i_2=L}^{R}...\sum\limits_{i_k=L}^{R}[\gcd(i_1,i_2,...,i_k)=k]\)
套路地进行莫比乌斯反演：
\(\sum\limits_{i_1=\frac{L-1}{k}+1}^{\frac{R}{k}}\sum\limits_{i_2=\frac{L-1}{k}+1}^{\frac{R}{k}}...\sum\limits_{i_k=\frac{L-1}{k}+1}^{\frac{R}{k}}[\gcd(i_1,i_2,...,i_k)=1]\)
\(\sum\limits_{i_1=\frac{L-1}{k}+1}^{\frac{R}{k}}\sum\limits_{i_2=\frac{L-1}{k}+1}^{\frac{R}{k}}...\sum\limits_{i_k=\frac{L-1}{k}+1}^{\frac{R}{k}}\sum\limits_{x|\gcd(i_1,i_2,...,i_k)}\mu(x)\)
\(\sum\limits_{x=1}^{\frac{R}{k}}\mu(x)\sum\limits_{i_1=\frac{L-1}{k}+1}^{\frac{R}{k}}\sum\limits_{i_2=\frac{L-1}{k}+1}^{\frac{R}{k}}...\sum\limits_{i_k=\frac{L-1}{k}+1}^{\frac{R}{k}}[x|\gcd(i_1,i_2,...,i_k)]\)
\(\sum\limits_{x=1}^{\frac{R}{k}}\mu(x)\sum\limits_{i_1=\frac{L-1}{k*x}+1}^{\frac{R}{k*x}}\sum\limits_{i_2=\frac{L-1}{k*x}+1}^{\frac{R}{k*x}}...\sum\limits_{i_k=\frac{L-1}{k*x}+1}^{\frac{R}{k*x}}1\)
\(\sum\limits_{x=1}^{\frac{R}{k}}\mu(x)(\frac{R}{k*x}-\frac{L-1}{k*x})^n\)

杜教筛求\(\sum\limits_{x=1}^{\frac{R}{k}}\mu(x)\)就可以除法分块了

code:

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=1e5+10;
const int inf=1e9;
const ll mod=1000000007;
int n,K,L,R;
int mu[maxn],sum[maxn];
ll ans;
bool vis[maxn];
vector<int>prime;
unordered_map<int,int>mp;
inline ll power(ll x,ll k,ll mod)
{
	ll res=1;
	while(k)
	{
		if(k&1)res=res*x%mod;
		x=x*x%mod;k>>=1;
	}
	return res;
}
inline void pre_work(int n)
{
	vis[1]=1;mu[1]=1;
	for(int i=2;i<=n;i++)
	{
		if(!vis[i])prime.push_back(i),mu[i]=-1;
		for(unsigned int j=0;j<prime.size()&&i*prime[j]<=n;j++)
		{
			vis[i*prime[j]]=1;
			if(i%prime[j]==0)break;
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=1;i<=n;i++)sum[i]=sum[i-1]+mu[i];
}
inline int getsum(int x)
{
	if(x<=100000)return sum[x];
	if(mp.count(x))return mp[x];
	ll res=1;
	for(int l=2,r;l<=x;l=r+1)
	{
		r=x/(x/l);
		res=(res-(r-l+1)*getsum(x/l)%mod)%mod;
	}
	return mp[x]=(res%mod+mod)%mod;
}
int main()
{
	pre_work(100000);
	scanf("%d%d%d%d",&n,&K,&L,&R);
	L=(L-1)/K,R=R/K;
	for(int l=1,r;l<=R;l=r+1)
	{
		r=min(L/l?L/(L/l):inf,R/(R/l));
		ans=((ans+1ll*(getsum(r)-getsum(l-1))*power(R/l-L/l,n,mod)%mod)%mod+mod)%mod;
	}
	printf("%lld",ans);
	return 0;
}