bzoj 3930: [CQOI2015]选数 莫比乌斯反演+杜教筛

题意

我们知道，从区间[L,H]（L和H为整数）中选取N个整数，总共有(H-L+1)^N种方案。小z很好奇这样选出的数的最大公约数的规律，他决定对每种方案选出的N个整数都求一次最大公约数，以便进一步研究。然而他很快发现工作量太大了，于是向你寻求帮助。你的任务很简单，小z会告诉你一个整数K，你需要回答他最大公约数刚好为K的选取方案有多少个。由于方案数较大，你只需要输出其除以1000000007的余数即可。
1≤N,K≤10^9，1≤L≤H≤10^9，H-L≤10^5

分析

显然按照一般套路就是设$f(d)$表示$gcd==d$的方案数，$g(d)$表示$d|gcd$的方案数。
容易得到$g(d)=(\lfloor\frac{H}{d}\rfloor-\lfloor\frac{L-1}{d}\rfloor)^n$还有

$$f(d)=\sum_{d|k}g(k)\mu(\frac{k}{d})=\sum_{i=1}^{\lfloor\frac{H}{d}\rfloor}\mu(i)g(di)$$
$$f(d)=\sum_{i=1}^{\lfloor\frac{H}{d}\rfloor}\mu(i)(\lfloor\frac{H}{di}\rfloor-\lfloor\frac{L-1}{di}\rfloor)$$
看到下取整就知道可以分块啦。
处理 $\mu$的前缀和的话，小的预处理，大的杜教筛即可。

据说有一种超水的dp做法，但是懒得学啦。

代码

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<map>
using namespace std;

typedef long long LL;

const int N=10000005;
const int inf=0x7fffffff;
const int MOD=1000000007;

int tot,n,k,l,r,prime[N/10],mu[N];
bool not_prime[N];
map<int,int> w;

void prework(int n)
{
    mu[1]=1;
    for (int i=2;i<=n;i++)
    {
        if (!not_prime[i]) prime[++tot]=i,mu[i]=-1;
        for (int j=1;j<=tot&&i*prime[j]<=n;j++)
        {
            not_prime[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
    for (int i=1;i<=n;i++) mu[i]+=mu[i-1];
}

int ksm(int x,int y)
{
    int ans=1;
    while (y)
    {
        if (y&1) ans=(LL)ans*x%MOD;
        x=(LL)x*x%MOD;y>>=1;
    }
    return ans;
}

int get_mu(int n)
{
    if (n<=10000000) return mu[n];
    if (w[n]) return w[n];
    int ans=1;
    for (int i=2,last;i<=n;i=last+1)
    {
        last=n/(n/i);
        ans=(ans-(LL)(last-i+1)*get_mu(n/i)%MOD)%MOD;;
    }
    return w[n]=ans;
}

int solve()
{
    int ans=0;l--;
    for (int i=1,last;i<=r;i=last+1)
    {
        last=min(l/i?l/(l/i):inf,r/(r/i));
        ans=(ans+(LL)(get_mu(last)-get_mu(i-1))*ksm(r/(k*i)-l/(k*i),n)%MOD)%MOD;
    }
    return (ans+MOD)%MOD;
}

int main()
{
    scanf("%d%d%d%d",&n,&k,&l,&r);
    prework(10000000);
    printf("%d",solve());
}