[SDOI2017]数字表格

[SDOI2017]数字表格

标签: 数学 莫比乌斯


题目描述

Doris刚刚学习了fibonacci数列。用f[i]f[i]f[i]表示数列的第iii项,那么
$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$取模。

输入输出格式

输入格式:

有多组测试数据。
第一个一个数T,表示数据组数。
接下来T行,每行两个数n,m

输出格式:

输出T行,第i行的数是第i组数据的结果

输入输出样例

输入样例#1:
3
2 3
4 5
6 7

输出样例#1:
1
6
960

说明

对10%的数据,$1≤n,m≤100$
对30%的数据,$1≤n,m≤1000$
另外存在30%的数据,$T≤3$
对100%的数据,$T≤1000,1≤n,m≤10^6$
时间限制:5s
内存限制:128MB

题解

做了一些题目后感觉莫比乌斯题目套路都比较明显了~~
一般都是关于$gcd$这个神奇的东西
一句话描述题目:
$$Ans=\prod_{i=1}^{n}\prod_{j=1}^{m}Fi(gcd(i,j))$$
其中$Fi$表示斐波那契数列
那么按照套路,我们设:
$$f(x)=\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)==x]$$
$$F(x)=\sum_{x|d}f(d)=\left \lfloor \frac{n}{x} \right \rfloor\left \lfloor \frac{m}{x} \right \rfloor$$
那么我们枚举$gcd$
$$Ans=\prod_{x=1}^{min(n,m)}Fi(x)^{f(x)}$$
指数部分很简单,我们直接反演即可
$$Ans=\prod_{x=1}^{min(n,m)}Fi(x)^{\sum_{x|d}u(\frac{d}{x})\frac{n}{d}\frac{m}{d}}$$
更改枚举项
$$Ans=\prod_{x=1}^{min(n,m)}Fi(x)^{\sum_{d=1}^{min(n,m)}u(d)\frac{n}{dx}\frac{m}{dx}}$$
化简到这里还是不能直接做,我们知道$Fi$和$u$都可以预处理,因此我们要把这两个”撮合”到一起
外部更改枚举项,我们枚举$dx$(即x的倍数)
$$Ans=\prod_{T=1}^{min(n,m)}(\prod_{k|T}Fi(k)^{u(\frac{T}{k})})^{\frac{n}{T}\frac{m}{T}}$$
我们设$g(T)=\prod_{k|T}Fi(k)^{u(\frac{T}{k})}$
这个$g(T)$是可以预处理的
那么
$$Ans=\prod_{T=1}^{min(n,m)}g(T)^{\frac{n}{T}\frac{m}{T}}$$
最后用整除分块即可解决本题了

代码

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#define N 1000010
#define ll long long
#define debug system("pause")
using namespace std;
const ll p = 1e9 + 7;

int t, n, m, cnt = 0;
int pri[N], u[N];
ll f[N], g[N], sum[N];
bool vis[N];

ll pow_mul(ll a, ll b)
{
    ll ans = 1;
    while(b)
    {
        if(b & 1) ans = ans * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return ans;
}

ll calc(ll a, ll b)
{
    if(b == 0) return 1;
    if(b == 1) return a;
    if(b == -1) return pow_mul(a, p - 2);
}

void get()
{
    u[1] = 1; int t = N - 10;
    for(int i = 2; i <= t; i ++)
    {
        if(!vis[i]) u[i] = -1, pri[++ cnt] = i;
        for(int j = 1; j <= cnt && i * pri[j] <= t; j ++)
        {
            vis[i * pri[j]] = true;
            if(i % pri[j]) u[i * pri[j]] = -u[i];
            else break;
        }
    }
    g[0] = 0; g[1] = 1;
    for(int i = 2; i <= t; i ++)
        g[i] = (g[i - 1] + g[i - 2]) % p;
    for(int i = 0; i <= t; i ++) 
        f[i] = sum[i] = 1;
    for(int i = 1; i <= t; i ++)
        for(int j = 1; i * j <= t; j ++)
            f[i * j] = (f[i * j] * calc(g[i], u[j])) % p;
    for(int i = 1; i <= t; i ++) sum[i] = sum[i - 1] * f[i] % p;       
}

int main()
{
    get();
    scanf("%d", &t);
    while(t --)
    {
        scanf("%d%d", &n, &m);
        if(n > m) swap(n, m);
        ll ans = 1;
        for(int l = 1, r; l <= n; l = r + 1)
        {
            r = min(m / (m / l), n / (n / l));
            ll tmp = (sum[r] * pow_mul(sum[l - 1], p - 2) % p);
            ans = ans * pow_mul(tmp, 1ll * (n / l) * (m / l) % (p - 1)) % p;
        }
        printf("%lld\n", ans % p);
    }
    return 0;
} 

心得

巧妙的更改枚举项,把能预处理的合并在一起即可优化时间复杂度


 上一篇
莫比乌斯反演总结 莫比乌斯反演总结
莫比乌斯反演总结标签: 总结 1.序做了一段时间的莫比乌斯反演题目,决定还是总结一下莫比乌斯反演做题方法(套路)吧。便于以后的复习与理解 2.莫比乌斯函数莫比乌斯函数本质是一个容斥函数,它的定义如下:当$d=1$时$u(d)=1$当$d=
2019-03-02
下一篇 
HNOI2009 图的同构记数 HNOI2009 图的同构记数
HNOI2009 图的同构记数标签: 数学 置换 polya定理 题目描述小雪在了解到以上情况后,自认为直接挑战终极难题还有不少困难,于是决定先从简单的问题做起,具体来说,他对图同构问题产生了浓厚的兴趣。A图与B图被认为是同构的是指:A图
2019-02-23