4816: [Sdoi2017]数字表格

4816: [Sdoi2017]数字表格

Time Limit: 50 Sec Memory Limit: 128 MB
Submit: 537 Solved: 236
[Submit][Status][Discuss]
Description

Doris刚刚学习了fibonacci数列。用f[i]表示数列的第i项，那么
f[0]=0
f[1]=1
f[n]=f[n-1]+f[n-2],n>=2
Doris用老师的超级计算机生成了一个n×m的表格，第i行第j列的格子中的数是f[gcd(i,j)]，其中gcd(i,j)表示i,
j的最大公约数。Doris的表格中共有n×m个数，她想知道这些数的乘积是多少。答案对10^9+7取模。
Input

有多组测试数据。

第一个一个数T，表示数据组数。
接下来T行，每行两个数n,m
T<=1000,1<=n,m<=10^6
Output

输出T行，第i行的数是第i组数据的结果
Sample Input

3

2 3

4 5

6 7

Sample Output

1

6

960
HINT

Source

鸣谢infinityedge上传

[Submit][Status][Discuss]

枚举每个公约数对答案的贡献，大力推一波式子
记$N=min(n,m)$，上一个莫比乌斯反演
$Ans = \prod_{i=1}^n\prod_{j=1}^mf[gcd(i,j)]=\prod_{d=1}^Nf[d]^{\sum_{t=1}^{\lfloor \frac Nd \rfloor}\mu(t) \lfloor \frac n{dt} \rfloor \lfloor \frac m{dt} \rfloor}$
$\lfloor \frac Nd \rfloor$不同的取值只有$\sqrt N$种
当$d$确定时，$\lfloor \frac n{dt} \rfloor，\lfloor \frac m{dt} \rfloor$的不同取值分别为$\sqrt {\lfloor \frac n{dt} \rfloor},\sqrt {\lfloor \frac m{dt} \rfloor}$种
两个位置都直接暴力枚举不同的方案，复杂度为
$\sum_{i=1}^{\sqrt n}O(\sqrt i) + \sum_{i=1}^{\sqrt n}O(\sqrt {\lfloor \frac ni \rfloor}) \approx O(n^{\frac 34})$
可以预处理比较小的值，加上用hash记忆化搜索，卡卡就过了

#include<iostream>
#include<cstdio>
#include<cstring>
#include<vector>
#include<queue>
#include<algorithm>
#include<cmath>
#include<stack>
//#include<ctime>
#define min(a,b) ((a) < (b) ? (a) : (b))
using namespace std;

const int G = 1005;
const int maxn = 1E6 + 10;
typedef long long LL;
const LL mo = 1000000007;
const int Mod = 4999987;

struct data{
    int N,M,Ans; data(){}
    data(int N,int M,int Ans): N(N),M(M),Ans(Ans){}
};

int T,n,m,tot,Ans,pri[maxn],f[maxn],p[maxn],Inv[maxn],c[maxn*4],g[G][G],mi[maxn];
bool not_pri[maxn]; LL mu[maxn];

vector <data> v[Mod];

inline int Mul(const LL &x,const LL &y) {return x * y % mo;}
inline int Add(const int &x,const int &y) {return x + y < mo ? x + y : x + y - mo;}
inline int Dec(const int &x,const int &y) {return x - y >= 0 ? x - y : x - y + mo;}

inline int ksm(int x,int y)
{
    int ret = 1;
    for (; y; y >>= 1LL)
    {
        if (y & 1LL) ret = Mul(ret,x);
        x = Mul(x,x);
    }
    return ret;
}

inline int Get(int N,int M)
{
    LL ret = 0;
    for (int i = 1,last; i <= N; i = last + 1)
    {
        int x = N / i,y = M / i; LL A;
        last = min(N / x,M / y); A = 1LL * x * y;
        ret += (mu[last] - mu[i - 1]) * A;
    }
    return ret % (mo - 1);
}

inline int Search(int N,int M)
{
    int pos = mi[N] + mi[M]; if (pos >= Mod) pos -= Mod;
    for (int i = 0; i < v[pos].size(); i++)
        if (v[pos][i].N == N && v[pos][i].M == M) return v[pos][i].Ans;
    return -1;
}

inline void Insert(int N,int M,int Ans)
{
    int pos = mi[N] + mi[M]; if (pos >= Mod) pos -= Mod;
    v[pos].push_back(data(N,M,Ans));
}

inline LL Calc(int N,int M)
{
    if (N < G && M < G) return g[N][M];
    LL ret = Search(N,M);
    if (ret != -1) return ret; ret = 0;
    for (int i = 1,last; i <= N; i = last + 1)
    {
        int x = N / i,y = M / i; LL A;
        last = min(N / x,M / y); A = 1LL * x * y;
        ret += (mu[last] - mu[i - 1]) * A;
    }
    ret %= (mo - 1); Insert(N,M,ret); return ret;
}

inline int getint()
{
    char ch = getchar(); int ret = 0;
    while (ch < '0' || '9' < ch) ch = getchar();
    while ('0' <= ch && ch <= '9')
        ret = ret * 10 + ch - '0',ch = getchar();
    return ret;
}

int main()
{
    #ifdef DMC
        freopen("DMC.txt","r",stdin);
        //freopen("1.out","w",stdout);
    #endif

    mu[1] = f[1] = p[0] = p[1] = mi[1] = 1;
    for (int i = 2; i < maxn; i++)
    {
        mi[i] = 1LL * i * i % Mod;
        f[i] = Add(f[i - 1],f[i - 2]);
        if (!not_pri[i]) pri[++tot] = i,mu[i] = -1;
        for (int j = 1; j <= tot; j++)
        {
            int Nex = i * pri[j];
            if (Nex >= maxn) break;
            not_pri[Nex] = 1;
            if (i % pri[j] == 0) break;
            mu[Nex] = -mu[i];
        }
    }
    for (int i = 2; i < maxn; i++)
        mu[i] += mu[i - 1],p[i] = Mul(p[i - 1],f[i]);
    Inv[maxn - 1] = ksm(p[maxn - 1],mo - 2);
    for (int i = maxn - 2; i >= 0; i--)
        Inv[i] = Mul(Inv[i + 1],f[i + 1]);
    for (int N = 1; N < G; N++)
    {
        int *gg = g[N];
        for (int M = N; M < G; M++)
            gg[M] = Get(N,M);
    }

    T = getint();
    while (T--)
    {
        n = getint(); m = getint();
        if (n > m) swap(n,m); Ans = 1;
        for (int i = 1,last; i <= n; i = last + 1)
        {
            int x = n / i,y = m / i,A;
            last = min(n / x,m / y);
            A = Mul(p[last],Inv[i - 1]);
            Ans = Mul(Ans,ksm(A,Calc(x,y)));
        }
        printf("%d\n",Ans);
    }
    return 0;
}