[SDOI2017] 数字表格

标签：https   open   mat   getchar   running   space   bug   return   其他   

Portal

还是一样, 我们应用枚举GCD的套路:
\[ Ans = \prod_{i}\prod_{j} f(gcd(i, j)) \= \prod_{d} f(d)^ {\sum_{i}\sum_{j} [(i, j) = d]}\= \prod_{d} f(d)^ {\sum_{s} \mu(s)\frac{n}{ds}\frac{m}{ds}} \\]
之后我们继续应用合并元的套路:
\[ Ans = \prod_{d} \prod_{i | d} f(i) ^ {\mu(\frac{d}{i})\frac{n}{d}\frac{m}{d}} \\]
之后考虑整除分块, 因为要对\(a^{\frac{n}{d}\frac{m}{d}}\)整除分块,所以必须把所有其他项视为整体\(g(n)\)
\[ = \prod_{d} f(i) ^ {\sum_{i | d}\mu(\frac{d}{i})\frac{n}{d}\frac{m}{d}} \= \prod_{d} (\prod_{i | d} f(i) ^ {\mu(\frac{d}{i})}) ^ {\frac{n}{d}\frac{m}{d}} \Let ~ g(n) = \prod_{d | n} f(d) ^ {\mu(\frac{n}{d})} \Ans = \prod_{d = 1}^{min(n, m)} g(d) ^ {\frac{n}{d} \frac{m}{d}}\\]
发现\(n\)的范围很小,于是我们可以枚举约数/倍数暴力计算. 这里枚举倍数的复杂度较小,所以我们采用枚举倍数.

然后就做完了,总复杂度:\[O(T\sqrt n ~log_{2}Mod + n ~ln~n + n log_2n)\]

Code

#include<bits/stdc++.h>
using namespace std;
#define rep(i, a, b) for(int i = (a), i##_end_ = (b); i <= i##_end_; ++i)
#define drep(i, a, b) for(int i = (a), i##_end_ = (b); i >= i##_end_; --i)
#define clar(a, b) memset((a), (b), sizeof(a))
#define debug(...) fprintf(stderr, __VA_ARGS__)
typedef long long LL;
typedef long double LD;
int read() {
    char ch = getchar();
    int x = 0, flag = 1;
    for (;!isdigit(ch); ch = getchar()) if (ch == ‘-‘) flag *= -1;
    for (;isdigit(ch); ch = getchar()) x = x * 10 + ch - 48;
    return x * flag;
}
void write(int x) {
    if (x < 0) putchar(‘-‘), x = -x;
    if (x >= 10) write(x / 10);
    putchar(x % 10 + 48);
}

const int Maxn = 1e6 + 9, Mod = 1000000007;
LL fpm(LL base, LL tims) {
    LL r = 1;
    for (; tims; tims >>= 1) {
        if (tims & 1) r = 1ll * r * base % Mod;
        base = 1ll * base * base % Mod;
    }
    return r;
}

static int prime[Maxn], tot, mu[Maxn];
static bool isnprime[Maxn];
void linearSieve() {
    mu[1] = 1;
    rep (i, 2, Maxn - 1) {
        if (!isnprime[i]) prime[++tot] = i, mu[i] = -1;
        for (int k, j = 1; j <= tot && (k = prime[j] * i) < Maxn; ++j) {
            isnprime[k] = 1;
            if (i % prime[j] == 0) {
                mu[k] = 0;
                break;
            } else mu[k] = -mu[i];
        }
    }
}

static int Fib[Maxn], invFib[Maxn];
static int F[Maxn], prodF[Maxn], invProdF[Maxn];
void init() {
    linearSieve();

    Fib[0] = 0, Fib[1] = 1;
    rep (i, 2, Maxn - 1) Fib[i] = (Fib[i - 1] + Fib[i - 2]) % Mod;
    rep (i, 1, Maxn - 1) invFib[i] = fpm(Fib[i], Mod - 2), F[i] = 1;

    rep (i, 1, Maxn - 1)
        for (int j = i; j < Maxn; j += i) {
            if (mu[j / i] == -1) F[j] = 1ll * F[j] * invFib[i] % Mod;
            if (mu[j / i] == 1) F[j] = 1ll * F[j] * Fib[i] % Mod;
        }

    prodF[0] = invProdF[0] = 1;
    rep (i, 1, Maxn - 1) {
        prodF[i] = prodF[i - 1] * 1ll * F[i] % Mod;
        invProdF[i] = fpm(prodF[i], Mod - 2);
    }
}

void solve() {
    int T = read();
    while (T--) {
        int n = read(), m = read();
        int Limit = min(n, m); LL ans = 1;
        for (int l = 1, r; l <= Limit; l = r + 1) {
            r = min(n / (n / l), m / (m / l));
            ans = ans * 1ll * fpm(1ll * prodF[r] * invProdF[l - 1] % Mod, 1ll * (n / l) * (m / l)) % Mod;
        }
        printf("%lld\n", ans);
    }
}

int main() {
    freopen("loj2000.in", "r", stdin);
    freopen("loj2000.out", "w", stdout);

    init();
    solve();

#ifdef Qrsikno
    debug("\nRunning time: %.3lf(s)\n", clock() * 1.0 / CLOCKS_PER_SEC);
#endif
    return 0;
}

[SDOI2017] 数字表格

标签：https   open   mat   getchar   running   space   bug   return   其他   

原文地址：https://www.cnblogs.com/qrsikno/p/10201651.html