[BZOJ4816][Sdoi2017]数字表格

[BZOJ4816][Sdoi2017]数字表格

试题描述

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取模。

输入

有多组测试数据。

第一个一个数T，表示数据组数。

接下来T行，每行两个数n,m

T<=1000,1<=n,m<=10^6

输出

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

输入示例

3
2 3
4 5
6 7

输出示例

1
6
960

数据规模及约定

见“输入”

题解

要求这么个东西

考虑枚举 gcd(i, j)

于是分块套分块，好像能卡过。。。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <algorithm>
using namespace std;

int read() {
	int x = 0, f = 1; char c = getchar();
	while(!isdigit(c)){ if(c == '-') f = -1; c = getchar(); }
	while(isdigit(c)){ x = x * 10 + c - '0'; c = getchar(); }
	return x * f;
}

#define maxn 1000010
#define maxt 2828
#define MOD 1000000007
#define LL long long

void gcd(LL a, LL b, LL& x, LL& y) {
	if(!b){ x = 1; y = 0; return ; }
	gcd(b, a % b, y, x); y -= a / b * x;
	return ;
}
int inv(LL a) {
	LL x, y;
	gcd(a, MOD, x, y);
	return (x % MOD + MOD) % MOD;
}

int prime[maxn], cp, mu[maxn], smu[maxn], f[maxn], tf[maxn], invt[maxn];
bool vis[maxn];
void init() {
	mu[1] = 1; smu[1] = 1;
	for(int i = 2; i < maxn; i++) {
		if(!vis[i]) prime[++cp] = i, mu[i] = -1;
		for(int j = 1; i * prime[j] < maxn && j <= cp; j++) {
			vis[i*prime[j]] = 1;
			if(i % prime[j] == 0){ mu[i*prime[j]] = 0; break; }
			mu[i*prime[j]] = -mu[i];
		}
		smu[i] = smu[i-1] + mu[i];
	}
	int f1 = 0, f2 = 1; tf[1] = invt[1] = 1; invt[0] = 1;
	for(int i = 2; i < maxn; i++) {
		swap(f1, f2); f2 += f1; if(f2 >= MOD) f2 -= MOD;
		tf[i] = (LL)tf[i-1] * f2 % MOD;
		invt[i] = inv(tf[i]);
	}
	return ;
}

LL get[maxt][maxt];
LL g(int n, int m) {
	if(n < maxt && m < maxt) {
		LL& ans = get[n][m];
		if(ans) return ans;
		for(int i = 1, lst; i <= n; i = lst + 1) {
			lst = min(n / (n / i), m / (m / i));
			ans += (LL)(n / i) * (m / i) * (smu[lst] - smu[i-1]);
		}
		return ans;
	}
	LL ans = 0;
	for(int i = 1, lst; i <= n; i = lst + 1) {
		lst = min(n / (n / i), m / (m / i));
		ans += (LL)(n / i) * (m / i) * (smu[lst] - smu[i-1]);
	}
	return ans;
}

int Pow(int a, LL b) {
	int t = a, ans = 1;
	while(b) {
		if(b & 1) ans = (LL)ans * t % MOD;
		t = (LL)t * t % MOD; b >>= 1;
	}
	return ans;
}

int main() {
	init();
	
	int T = read();
	while(T--) {
		int n = read(), m = read(), ans = 1;
		if(n > m) swap(n, m);
		for(int i = 1, lst; i <= n; i = lst + 1) {
			lst = min(n / (n / i), m / (m / i));
			ans = (LL)ans * Pow((LL)tf[lst] * invt[i-1] % MOD, g(n / i, m / i)) % MOD;
		}
		printf("%d\n", ans);
	}
	
	return 0;
}