【bzoj4816】[Sdoi2017]数字表格 莫比乌斯反演

题目描述

Doris刚刚学习了fibonacci数列。用f[i]表示数列的第i项，那么
f[0]=0
f[1]=1
f
=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<=1000,1<=n,m<=10^6
输出

输出T行，第i行的数是第i组数据的结果
样例输入

3

2 3

4 5

6 7

样例输出

1

6

960

题解

莫比乌斯反演

$\prod\limits_{i=1}^n\prod\limits_{j=1}^mf(\gcd(i,j))\\=\prod\limits_{d=1}^{min(n,m)}(\prod\limits_{i=1}^n\prod\limits_{j=1}^m[\gcd(i,j)=d]·f(d))\\=\prod\limits_{d=1}^{min(n,m)}f(d)^{\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=d]}\\=\prod\limits_{d=1}^{min(n,m)}f(d)^{\sum\limits_{i=1}^{\lfloor\frac nd\rfloor}\sum\limits_{j=1}^{\lfloor \frac md\rfloor}[\gcd(i,j)=1]}\\=\prod\limits_{d=1}^{min(n,m)}f(d)^{\sum\limits_{i=1}^{\lfloor\frac nd\rfloor}\sum\limits_{j=1}^{\lfloor \frac md\rfloor}\sum\limits_{k|i\&k\&j}\mu(d)}\\=\prod\limits_{d=1}^{min(n,m)}f(d)^{\sum\limits_{k=1}^{min(\lfloor\frac nd\rfloor,\lfloor\frac md\rfloor)}\mu(k)\lfloor\frac n{dk}\rfloor\lfloor\frac m{dk}\rfloor}$

到了这一步我们可以选择分块套分块，不过显然时间复杂度不足以应对多组询问。

继续令$D=dk$，可以得到：

$\prod\limits_{d=1}^{min(n,m)}f(d)^{\sum\limits_{k=1}^{min(\lfloor\frac nd\rfloor,\lfloor\frac md\rfloor)}\mu(k)\lfloor\frac n{dk}\rfloor\lfloor\frac m{dk}\rfloor}\\=\prod\limits_{D=1}^{min(n,m)}(\prod\limits_{d|D}f(d)^{\mu(\frac Dd)})^{\lfloor\frac nD\rfloor\lfloor\frac mD\rfloor}\\=\prod\limits_{D=1}^{min(n,m)}(t(D))^{\lfloor\frac nD\rfloor\lfloor\frac mD\rfloor}\\(t(D)=\prod\limits_{d|D}f(d)^{\mu(\frac Dd)})$

于是线性筛出$\mu$，递推出f，预处理出f的逆元，进而使用$O(n\ln n)$的时间预处理出t数组。这里有一个小技巧：先枚举$\frac Dd$，当它的$\mu$值等于0时不作任何处理。亲测可以有效减少时间。

然后预处理出t数组之后分块处理询问即可。

时间复杂度为O(跑得过)$O(n\ln n+T\sqrt n)$，实际上跑了37s。

#include <cstdio>
#include <cstring>
#include <algorithm>
#define N 1000010
using namespace std;
typedef long long ll;
const int k = 1000000 , mod = 1000000007;
int mu
, prime
, tot;
ll f
, inv
, t
, sum
;
bool np
;
ll pow(ll x , ll y)
{
ll ans = 1;
while(y)
{
if(y & 1) ans = ans * x % mod;
x = x * x % mod , y >>= 1;
}
return ans;
}
int main()
{
int i , j , T , n , m;
ll ans;
sum[0] = t[1] = f[1] = inv[1] = mu[1] = 1;
for(i = 2 ; i <= k ; i ++ )
{
t[i] = 1 , f[i] = (f[i - 1] + f[i - 2]) % mod , inv[i] = pow(f[i] , mod - 2) % mod;
if(!np[i]) mu[i] = -1 , prime[++tot] = i;
for(j = 1 ; j <= tot && i * prime[j] <= k ; j ++ )
{
np[i * prime[j]] = 1;
if(i % prime[j] == 0)
{
mu[i * prime[j]] = 0;
break;
}
else mu[i * prime[j]] = -mu[i];
}
}
for(j = 1 ; j <= k ; j ++ )
{
if(!mu[j]) continue;
for(i = 1 ; i * j <= k ; i ++ )
{
if(mu[j] == 1) t[i * j] = t[i * j] * f[i] % mod;
else t[i * j] = t[i * j] * inv[i] % mod;
}
}
for(i = 1 ; i <= k ; i ++ ) sum[i] = sum[i - 1] * t[i] % mod;
scanf("%d" , &T);
while(T -- )
{
scanf("%d%d" , &n , &m) , ans = 1;
for(i = 1 ; i <= n && i <= m ; i = j + 1)
j = min(n / (n / i) , m / (m / i)) , ans = ans * pow(sum[j] * pow(sum[i - 1] , mod - 2) % mod , (ll)(n / i) * (m / i)) % mod;
printf("%lld\n" , ans);
}
return 0;
}