BZOJ 3930 CQOI2015 選數 莫比烏斯反演

題目見 http://pan.baidu.com/s/1o6zajc2

此外不知道H-L<=10^5這個條件是幹嘛的。。。。

#include <map>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define M 10001000
#define INF 0x3f3f3f3f
#define MOD 1000000007
using namespace std;
int mu[M],prime[1001001],tot;
bool not_prime[M];
map<int,long long> mu_sum;
long long n,d,l,r;
void Linear_Shaker()
{
	int i,j;
	mu[1]=1;
	for(i=2;i<=10000000;i++)
	{
		if(!not_prime[i])
		{
			mu[i]=-1;
			prime[++tot]=i;
		}
		for(j=1;prime[j]*i<=10000000;j++)
		{
			not_prime[prime[j]*i]=true;
			if(i%prime[j]==0)
			{
				mu[prime[j]*i]=0;
				break;
			}
			mu[prime[j]*i]=-mu[i];
		}
	}
	for(i=1;i<=10000000;i++)
		mu[i]+=mu[i-1];
}
long long Mu_Sum(int x)
{
	if(x<=10000000)
		return mu[x];
	if(mu_sum.find(x)!=mu_sum.end())
		return mu_sum[x];
	long long i,last,re=1;
	for(i=1;i<=x;i=last+1)
	{
		last=x/(x/i);
		if(x/i-1)
			re-=(Mu_Sum(last)-Mu_Sum(i-1))*(x/i-1);
	}
	return mu_sum[x]=re;
}
long long Quick_Power(long long x,int y)
{
	long long re=1;
	while(y)
	{
		if(y&1) (re*=x)%=MOD;
		(x*=x)%=MOD; y>>=1;
	}
	return re;
}
long long Solve()
{
	long long i,last,re=0;
	for(i=1;i<=r;i=last+1)
	{
		last=min(r/(r/i),l/i?(l/(l/i)):INF);
		re+=(Mu_Sum(last)-Mu_Sum(i-1))*Quick_Power(r/i-l/i,n);
		re%=MOD;
	}
	return (re%MOD+MOD)%MOD;
}
int main()
{
	cin>>n>>d>>l>>r;
	l=(l-1)/d;r=r/d;
	Linear_Shaker();
	cout<<Solve()<<endl;
}