luogu3172 [CQOI2015]选数 莫比乌斯反演+杜教筛

标签：数论   http   题目   math   for   res   pow   main   cpp   

link

题目大意：有N个数，每个数都在区间[L,H]之间，请求出所有数的gcd恰好为K的方案数

推式子

首先可以把[L,H]之间的数字gcd恰好为K转化为[(L-1)/K+1,H/K]之间数字gcd恰好为1

然后就可以反演了

下面手误把所有的H都打成了R

\(\sum_{i_1=L}^R\sum_{i_2=L}^R\dots\sum_{i_N=L}^R[\gcd(i_1,i_2,\dots,i_N)=1]\)

\(\sum_{i_1=L}^R\sum_{i_2=L}^R\dots\sum_{i_N=L}^R\sum_{d|i_1,d|i_2,\dots,d|i_N}\mu(d)\)

\(\sum_{d=1}^R\mu(d)\left(\frac{R}d-\frac{L-1}d\right)^N\)

显然可以打数论分块

但是这道题H范围达到10的9次方，果断杜教筛

注意由于我们的d是枚举到R的（因为右边式子的关系不是相乘，而是相减，所以大于L的项也有效）而>=L的数整除L-1的值为0，所以会导致除法爆炸，所以需要特判，详细请见代码

#include <cstdio>
#include <utility>
#include <map>
using namespace std;

int p = 1000000007, fuck = 2000000;
int n, k, l, h;
bool vis[2000010];
int prime[2000000], tot, mu[2000010];
map<int, int> memory;

int qpow(int x, int y)
{
    int res = 1;
    for (x %= p; y > 0; y >>= 1, x = x * (long long)x % p) if (y & 1) res = res * (long long)x % p;
    return res;
}

int chumu(int x)
{
    if (x <= fuck) return mu[x];
    if (memory.count(x)) return memory[x];
    int res = 1;
    for (int i = 2, j; i <= x; i = j + 1)
    {
        j = x / (x / i);
        res -= (j - i + 1) * chumu(x / i);
    }
    return memory[x] = res;
}

int main()
{
    scanf("%d%d%d%d", &n, &k, &l, &h);
    h /= k, l = (l - 1) / k;
    mu[1] = 1;
    for (int i = 2; i <= fuck; i++)
    {
        if (vis[i] == false) prime[++tot] = i, mu[i] = -1;
        for (int j = 1; j <= tot && i * prime[j] <= fuck; j++)
        {
            vis[i * prime[j]] = true;
            if (i % prime[j] == 0) break;
            mu[i * prime[j]] = -mu[i];
        }
        mu[i] += mu[i - 1];
    }
    int res = 0;
    for (int i = 1, j; i <= h; i = j + 1)
    {
        j = min(h / (h / i), l / i == 0 ? h : l / (l / i));
        res += qpow((h / i - l / i), n) * (long long)(chumu(j) - chumu(i - 1)) % p;
        res %= p;
        if (res < 0) res += p;
    }
    printf("%d\n", res);
    return 0;
}

luogu3172 [CQOI2015]选数 莫比乌斯反演+杜教筛

标签：数论   http   题目   math   for   res   pow   main   cpp   

原文地址：https://www.cnblogs.com/oier/p/10323200.html