bzoj 4816: 洛谷 P3704: [SDOI2017]数字表格
洛谷很早以前就写过了，今天交到bzoj发现TLE了。
检查了一下发现自己复杂度是错的。
题目传送门：洛谷P3704。
题意简述：
求 \(\prod_{i=1}^{N}\prod_{j=1}^{M}F_{\gcd(i,j)}\bmod mod\) ，其中 \(F_{i}\) 是斐波那契数列的第 \(i\) 项， \(mod=10^9+7\) 。
\(T\) 组数据。
题解：
喜闻乐见的推式子时间。
不失一般性，假设 \(N\le M\) 。
\[\begin{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\right]\right)}\end{aligned}\]
右上角的指数部分是老套路了。
\[\begin{align*}&= \sum_{i=1}^{N}\sum_{j=1}^{M}\left[\gcd(i,j)=k\right]\\&= \sum_{i=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{M}{k}\right\rfloor}\left[\gcd(i,j)=1\right]\\&= \sum_{d=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\mu(d)\left\lfloor\frac{N}{kd}\right\rfloor\left\lfloor\frac{M}{kd}\right\rfloor\end{align*}\]
所以
\[\begin{align*} &= \prod_{k=1}^{N}{F_{k}}^{\left(\sum_{d=1}^{\left\lfloor\frac{N}{k}\right\rfloor}\mu(d)\left\lfloor\frac{N}{kd}\right\rfloor\left\lfloor\frac{M}{kd}\right\rfloor\right)}\\ &= \prod_{T=1}^{N}\left(\prod_{k|T}{F_{k}}^{\mu(\frac{T}{k})}\right)^{\left\lfloor\frac{N}{T}\right\rfloor\left\lfloor\frac{M}{T}\right\rfloor} \end{align*}\]
令 \(f(n)=\prod_{d|n}{F_{d}}^{\mu(\frac{n}{d})}\) 。
则
\[=\prod_{T=1}^{N}{f(T)}^{\left\lfloor\frac{N}{T}\right\rfloor\left\lfloor\frac{M}{T}\right\rfloor}\]
外层数论分块求出。内层的 \(f(T)\) 直接暴力预处理，每个数直接乘到它的倍数中，复杂度 \(\Theta(n\log n)\)。
注意实现的时候的时间复杂度，我因为实现多了快速幂的一个 \(\log\) 被卡了。
正确的时间复杂度应该是 \(\Theta(N(\log N+\log mod)+T\sqrt{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\n", A); 73 } 74 return 0; 75 }
