洛谷很早以前就写过了，今天交到bzoj发现TLE了。
检查了一下发现自己复杂度是错的。
题目传送门：洛谷P3704。
题意简述：
求 (prod_{i=1}^{N}prod_{j=1}^{M}F_{gcd(i,j)}mod mod) ，其中 (F_{i}) 是斐波那契数列的第 (i) 项， (mod=10^9+7) 。
(T) 组数据。
题解：
喜闻乐见的推式子时间。
不失一般性，假设 (Nle M) 。
[egin{aligned}&prod_{i=1}^{N}prod_{j=1}^{M}F_{gcd(i,j)} \=&prod_{k=1}^{N}{F_{k}}^{left(sum_{i=1}^{N};sum_{j=1}^{M};left[gcd(i,j)=k ight] ight)}end{aligned}]
右上角的指数部分是老套路了。
[egin{align*}&= sum_{i=1}^{N}sum_{j=1}^{M}left[gcd(i,j)=k ight]\&= sum_{i=1}^{leftlfloorfrac{N}{k} ight floor}sum_{j=1}^{leftlfloorfrac{M}{k} ight floor}left[gcd(i,j)=1 ight]\&= sum_{d=1}^{leftlfloorfrac{N}{k} ight floor}mu(d)leftlfloorfrac{N}{kd} ight floorleftlfloorfrac{M}{kd} ight floorend{align*}]
所以
[egin{align*} &= prod_{k=1}^{N}{F_{k}}^{left(sum_{d=1}^{leftlfloorfrac{N}{k} ight floor}mu(d)leftlfloorfrac{N}{kd} ight floorleftlfloorfrac{M}{kd} ight floor ight)}\ &= prod_{T=1}^{N}left(prod_{k|T}{F_{k}}^{mu(frac{T}{k})} ight)^{leftlfloorfrac{N}{T} ight floorleftlfloorfrac{M}{T} ight floor} end{align*}]
令 (f(n)=prod_{d|n}{F_{d}}^{mu(frac{n}{d})}) 。
则
[=prod_{T=1}^{N}{f(T)}^{leftlfloorfrac{N}{T} ight floorleftlfloorfrac{M}{T} ight floor}]
外层数论分块求出。内层的 (f(T)) 直接暴力预处理，每个数直接乘到它的倍数中，复杂度 (Theta(nlog n))。
注意实现的时候的时间复杂度，我因为实现多了快速幂的一个 (log) 被卡了。
正确的时间复杂度应该是 (Theta(N(log N+log mod)+Tsqrt{N}log mod)) 。
1 #include<cstdio> 2 #include<algorithm> 3 using namespace std; 4 5 #define mod 1000000007 6 #define LL long long 7 8 int Pow(int b, LL e) { 9 if (e < 0) e += mod - 1; 10 int a = 1; 11 for (; e; b = (LL)b * b % mod, e >>= 1) 12 if (e & 1) a = (LL)a * b % mod; 13 return a; 14 } 15 16 bool ip[1000005]; 17 int p[80005], pc; 18 int mu[1000005]; 19 int f[1000005], fr[1000005]; 20 21 void init() { 22 23 ip[1] = 1; 24 mu[1] = 1; 25 26 for (int i = 2; i <= 1000000; ++i) { 27 if (!ip[i]) { 28 p[++pc] = i; 29 mu[i] = -1; 30 } 31 for (int j = 1; j <= pc && (LL)p[j] * i <= 1000000; ++j) { 32 register int k = p[j] * i; 33 ip[k] = 1; 34 if (i % p[j]) mu[k] = -mu[i]; 35 else break; 36 } 37 } 38 39 for (int i = 1; i <= 1000000; ++i) 40 f[i] = 1, fr[i] = 1; 41 42 int A = 1, B = 0; 43 for (int i = 1; i <= 1000000; ++i) { 44 B = (A + B) % mod; 45 A = (B - A + mod) % mod; 46 int G[3] = {Pow(B, -1), 1, B}; 47 for (int j = i, k = 1; j <= 1000000; j += i, ++k) { 48 f[j] = (LL)f[j] * G[mu[k] + 1] % mod, 49 fr[j] = (LL)fr[j] * G[1 - mu[k]] % mod; 50 } 51 } 52 53 f[0] = fr[0] = 1; 54 for (int i = 1; i <= 1000000; ++i) 55 f[i] = (LL)f[i - 1] * f[i] % mod, 56 fr[i] = (LL)fr[i - 1] * fr[i] % mod; 57 } 58 59 int main() { 60 init(); 61 int T; 62 scanf("%d", &T); 63 while (T--) { 64 int N, M; 65 scanf("%d%d", &N, &M); 66 if (N > M) swap(N, M); 67 int A = 1; 68 for (int i = 1, j; i <= N; i = j + 1) { 69 j = min(N / (N / i), M / (M / i)); 70 A = (LL)A * Pow((LL)f[j] * fr[i - 1] % mod, (LL)(N / i) * (M / i)) % mod; 71 } 72 printf("%d ", A); 73 } 74 return 0; 75 }