bzoj3930 [CQOI2015]选数

题目链接：bzoj3930
题目大意：
从区间[L,H]（L和H为整数）中选取N个整数，总共有(H-L+1)^N种方案。小z会告诉你一个整数K，你需要回答他选出来的这N个数的最大公约数恰好为K的选取方案有多少个。由于方案数较大，你只需要输出其除以1000000007的余数即可。

注：观察样例解释，这N个数是可以一样的！！！我TM想了一晚上！！！

题解：
莫比乌斯反演+杜教筛
设选的这N个数分别为$a_1、a_2、a_3、···、a_n$
那么即求

$$\sum_{a_1=L}^{H}\sum_{a_2=L}^{H}···\sum_{a_n=L}^{H}[gcd(a_1,a_2,···,a_n)=K]$$
套路变形
$$\sum_{a_1=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}\sum_{a_2=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}···\sum_{a_n=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}[gcd(a_1,a_2,···,a_n)=1]$$
反演一下即 $$\sum_{a_1=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}\sum_{a_2=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}···\sum_{a_n=\lfloor\frac{L}{K}\rfloor}^{\lfloor\frac{H}{K}\rfloor}\sum_{t|(a_1,a_2,···,a_n)}\mu(t)$$
交换枚举倍数和约数即
$$\sum_{t=1}^{\lfloor\frac{H}{K}\rfloor}\sum_{a_1=\lfloor\frac{L}{Kt}\rfloor}^{\lfloor\frac{H}{Kt}\rfloor}\sum_{a_2=\lfloor\frac{L}{Kt}\rfloor}^{\lfloor\frac{H}{Kt}\rfloor}···\sum_{a_n=\lfloor\frac{L}{Kt}\rfloor}^{\lfloor\frac{H}{Kt}\rfloor}\mu(t)$$
即 $$\sum_{t=1}^{\lfloor\frac{H}{K}\rfloor}\mu(t)\times (\lfloor\frac{H}{Kt}\rfloor-\lfloor\frac{L-1}{Kt}\rfloor)^n$$

后面那一坨表示的是$a_i$有多少种取值。所以应该从一开始的有$R-L+1$种取值开始看(不然取整什么的会有误差？)，就是$\lfloor\frac{H}{Kt}\rfloor-\lfloor\frac{L-1}{Kt}\rfloor$。
好吧。。这里有点迷。。主要是试了其它打法不对在这里-1才对= =
看别人的好像弄了$f(i)F(i)$之类的东西。。不过最后化到的也是一样的(?)

剩下的就分块就好了，后面的快速幂求，关于$\mu$的前缀和，，额嗯直接上杜教筛就好了
然后看了一下status..woc4ms??我跑了几百ms啊(为什么居然分块杜教筛什么也那么快？)。再翻翻别人blog。。组合数学？容斥？算了不说了，去看看究竟是smg了。。

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<map>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long LL;
#define N 1001000

map<int,LL> M;
const LL mod=1000000007;
bool ispri[N];int cnt,lim,pri[N/4];LL mu[N];
int mymin(int x,int y){return (x<y)?x:y;}
void pre()
{
    cnt=0;mu[1]=1;
    for (int i=2;i<=lim;i++)
    {
        if (!ispri[i]){pri[++cnt]=i;mu[i]=-1;}
        for (int j=1;j<=cnt && i*pri[j]<=lim;j++)
        {
            ispri[i*pri[j]]=true;
            if (i%pri[j]==0)
            {
                mu[i*pri[j]]=0;
                break;
            }mu[i*pri[j]]=-mu[i];
        }
    }
    for (int i=2;i<=lim;i++) mu[i]+=mu[i-1];
}
LL qpow(LL x,int t)
{
    LL ret=1;
    while (t)
    {
        if (t&1) ret=(ret*x)%mod;
        x=(x*x)%mod;t>>=1;
    }return ret;
}
LL solve(int n)
{
    if (n<=lim) return mu[n];
    if (M.count(n)) return M[n];
    int i,r;LL ans=1LL;
    for (int i=2;i<=n;i=r+1)
    {
        r=mymin(n,n/(n/i));
        ans-=(LL)(r-i+1)*solve((LL)(n/i))%mod;
        ans%=mod;if (ans<0) ans+=mod;
    }M[n]=ans;
    return ans;
}
int main()
{
    //freopen("a.in","r",stdin);
    //freopen("a.out","w",stdout);
    int n,K,L,R,l,r,nt;LL ans=0;
    scanf("%d%d%d%d",&n,&K,&L,&R);
    lim=1e6;l=(L-1)/K;r=R/K;pre();
    for (int i=1;i<=r;i=nt+1)
    {
        if (l/i) nt=mymin(l/(l/i),r/(r/i));
        else nt=r/(r/i);
        ans+=(solve(nt)-solve(i-1)+mod)%mod*qpow((LL)(r/i-l/i),n)%mod;
        ans%=mod;if (ans<0) ans+=mod;
    }
    printf("%lld\n",ans);
    return 0;
}