• Luogu 3704 [SDOI2017]数字表格


    列一下式子:

          $prod_{i = 1}^{n}prod_{j = 1}^{m}fib_{gcd(i, j)}$

    很套路的变成这样:

          $prod_{d = 1}^{min(n, m)}fib_{d}^{sum_{i = 1}^{left lfloor frac{n}{d} ight floor}sum_{j = 1}^{left lfloor frac{m}{d} ight floor}[gcd(i, j) == 1]}$

    右上角的那个东西:

          $sum_{i = 1}^{left lfloor frac{n}{d} ight floor}sum_{j = 1}^{left lfloor frac{m}{d} ight floor}[gcd(i, j) == 1]$

    太熟悉了。

          $sum_{d = 1}^{min(n, m)}sum_{t = 1}^{min(left lfloor frac{n}{d} ight floor,left lfloor frac{m}{d} ight floor)}mu (t) * left lfloor frac{n}{td} ight floor * left lfloor frac{m}{td} ight floor$

    代回去之后按照套路枚举$T = dt$:

          $prod_{T = 1}^{min(n, m)}sum_{d | T}fib^{left lfloor frac{n}{T} ight floorleft lfloor frac{m}{T} ight floor mu (frac{T}{d})}_{d}$

    这样子的话我们记$h(i) = sum_{d | T}fib^{mu (frac{T}{d})}_{d}$

    原式就变为:

          $prod_{T = 1}^{min(n, m)}h(T)^{left lfloor frac{n}{T} ight floorleft lfloor frac{m}{T} ight floor}$。

    发现外面已经可以整除分块了。

    然而这个$h(i)$怎么办,这玩意...不能线性筛的呀。

    喂喂,不能线性筛就暴力算吧,暴力...似乎并不慢啊,其实是一个$O(nlogn)$。

    注意到$mu (i) == -1$的时候其实是乘上一个逆元。

    时间复杂度$O(MaxNlogMaxN + T sqrt{n} )$。

    复杂度写的并不严格。

    Code:

    #include <cstdio>
    #include <cstring>
    using namespace std;
    typedef long long ll;
    
    const int N = 1e6 + 5;
    const ll P = 1e9 + 7;
    
    int testCase, pCnt = 0, pri[N];
    ll mu[N], fib[N], h[N];
    bool np[N];
    
    template <typename T>
    inline void read(T &X) {
        X = 0; char ch = 0; T op = 1;
        for(; ch > '9' || ch < '0'; ch = getchar())
            if(ch == '-') op = -1;
        for(; ch >= '0' && ch <= '9'; ch = getchar())
            X = (X << 3) + (X << 1) + ch - 48;
        X *= op;
    }
    
    inline ll pow(ll x, ll y) {
        ll res = 1LL;
        for(; y > 0; y >>= 1) {
            if(y & 1) res = res * x % P;
            x = x * x % P;
        }
        return res;
    }
    
    ll inv(ll x) {
        return pow(x, P - 2);
    }
    
    void sieve() {
        fib[1] = 1LL;
        for(int i = 2; i < N; i++) 
            fib[i] = (fib[i - 1] + fib[i - 2]) % P;
        
    /*    for(int i = 1; i <= 20; i++)
            printf("%lld ", fib[i]);
        printf("
    ");   */
        
        mu[1] = 1;
        for(int i = 2; i < N; i++) {
            if(!np[i]) pri[++pCnt] = i, mu[i] = -1;
            for(int j = 1; j <= pCnt && i * pri[j] < N; j++) {
                np[i * pri[j]] = 1;
                if(i % pri[j] == 0) {
                    mu[i *  pri[j]] = 0;
                    break;
                }
                mu[i * pri[j]] = -mu[i];
            }
        }
        
    /*    for(int i = 1; i <= 20; i++)
            printf("%d ", mu[i]);
        printf("
    ");   */
        
        for(int i = 0; i < N; i++) h[i] = 1LL;
        for(int i = 1; i < N; i++) {
            if(!mu[i]) continue;
            for(int j = i; j < N; j += i) 
                h[j] = h[j] * ((mu[i] == 1) ? fib[j / i] : inv(fib[j / i])) % P;
        }
        
    /*    for(int i = 1; i <= 20; i++)
            printf("%lld ", h[i]);
        printf("
    ");   */
        
        for(int i = 1; i < N; i++)
            h[i] = 1LL * h[i] * h[i - 1] % P;
    }
    
    inline int min(int x, int y) {
        return x > y ? y : x;
    }
    
    int main() {
        sieve();
        for(read(testCase); testCase--; ) {
            int n, m; read(n), read(m);
            int rep = min(n, m); ll ans = 1LL;
            for(int l = 1, r; l <= rep; l = r + 1) {
                r = min((n / (n / l)), (m / (m / l)));
                ans = ans * pow(h[r] * inv(h[l - 1]) % P, 1LL * (n / l) * (m / l) % (P - 1)) % P;
            }
            printf("%lld
    ", ans);
        }
        return 0;
    }
    View Code
  • 相关阅读:
    [POJ2752]Seek the Name, Seek the Fame
    [HDOJ1247]Hat’s Words
    [POJ2001]Shortest Prefixes
    不用代码,10分钟打造属于自己的第一款小程序
    不用代码,10分钟打造属于自己的第一款小程序
    不用代码,10分钟打造属于自己的第一款小程序
    移动端iPhone系列适配问题的一些坑
    移动端iPhone系列适配问题的一些坑
    【前端统计图】echarts多条折线图和横柱状图实现
    【前端统计图】echarts多条折线图和横柱状图实现
  • 原文地址:https://www.cnblogs.com/CzxingcHen/p/9633747.html
Copyright © 2020-2023  润新知