U390993 斐波那契数列Ⅰ（加强版）- Official题解

# U390643 斐波那契数列Ⅰ（加强版）- 题解

[> Problem Link <](https://www.luogu.com.cn/problem/U390993)

## 设计说明

样例点是按照以下梯度设计的：
| 算法 | 期望得分 |
| :----------: | :----------: |
| 压位高精度 | $30$ 分 |
| FFT+矩阵快速幂 | $100$ 分 |

## 算法

### 压位高精度 | $30$ 分

[> 详见弱化版题解 <](https://www.luogu.com.cn/article/rxo1572a)

### FFT+矩阵快速幂 | $100$ 分

写出递推矩阵  
FFT优化乘法计算即可

## 工程细节

1. 矩阵简化
>斐波那契数列有递推矩阵:  
$
\begin{pmatrix}
  f_{n+1} \\
  f_n \\
\end{pmatrix}
=
\begin{pmatrix}
  f_1 \\
  f_0 \\
\end{pmatrix}
\begin{pmatrix}
  1&1\\
  1&0\\
\end{pmatrix}^n
$

>而这个递推矩阵一定可以写成：  
$
\begin{pmatrix}
  a+b&b\\
  b&a\\
\end{pmatrix}
$  
的形式。

>所以矩阵运算为:  
$
\begin{pmatrix}
  a+b&b\\
  b&a\\
\end{pmatrix}
\begin{pmatrix}
  a^{\prime}+b^{\prime}&b^{\prime}\\
  b^{\prime}&a^{\prime}\\
\end{pmatrix}
=
\begin{pmatrix}
  aa^{\prime}+ab^{\prime}+a^{\prime}b+2bb^{\prime}&ab^{\prime}+a^{\prime}b+bb^{\prime}\\
  ab^{\prime}+a^{\prime}b+bb^{\prime}&aa^{\prime}+bb^{\prime}\\
\end{pmatrix}
$  
仅存储 $a,b$ 即可,定义他们之间的乘法运算为：  
$
\begin{pmatrix}
  a&b
\end{pmatrix}
\begin{pmatrix}
  b^{\prime}&a^{\prime}
\end{pmatrix}
=
\begin{pmatrix}
  ab^{\prime}+a^{\prime}b+bb^{\prime}&aa^{\prime}+bb^{\prime}
\end{pmatrix}
$  
所以 
$f_n=
\text{后一个元素}\{\begin{pmatrix}
  0&1
\end{pmatrix}^n\}
$

2. 共轭复数优化 FFT
>略

时间复杂度：$\Theta(n\log n)$  
证明：略

## 代码实现

这是std： 
```cpp
#include <bits/stdc++.h>
#define PI 3.141592653589793238462643383279
#define MAXN 1048576
typedef int Int;
typedef double F;
typedef std::complex<F> C;
C I(0, 1);
void FFT(C*a, Int l, Int op = 1) {
	Int n = 1 << l, *rev = new Int[n];
	rev[0] = 0;
	for (Int i = 0; i != n; ++i) {
		rev[i] = (rev[i >> 1] >> 1) | (i & 1) << (l - 1);
		if (i < rev[i])std::swap(a[i], a[rev[i]]);
	}
	delete rev;
	for (Int i = 1; i != n; i <<= 1) {
		C wi(cos(PI / i), sin(PI / i)*op);
		for (Int j = 0; j != n; j += (i << 1)) {
			C w(1, 0);
			for (Int k = 0; k != i; ++k, w *= wi) {
				C x = a[j + k], y = w * a[i + j + k];
				a[j + k] = x + y, a[i + j + k] = x - y;
			}
		}
	}
	if (op == -1)
		for (Int i = 0; i != n; ++i)a[i] /= n;
}
void move(C*d, C*a, C*b, Int n) {
	d[n] = d[0];
	for (Int i = 0; i != n; ++i) {
		a[i] = (d[i] + std::conj(d[n - i])) * F(0.5);
		b[i] = (d[i] - std::conj(d[n - i])) * F(-0.5) * I;
	}
}
C D[MAXN], A1[MAXN], B1[MAXN], A2[MAXN], B2[MAXN];
inline Int R(F x) {
	return Int(x + 0.5);
}
class matrix {
	public:
		Int size = 1;
//		C v[MAXN];
		C *v;
		matrix() {
			v = new C[MAXN];
		}
		~matrix() {
			delete v;
		}
		/*这上面几行必须要加，不能直接C v[MAXN];
		不然编译器会优化，导致编译文件高达33MB*/
		void E() {
			v[0] = {1, 0};
		}
		void operator*=(const matrix &x) {
			Int l = 0, n = 1;
			while (n < size + x.size + 32)++l, n <<= 1;
			memcpy(D, v, n * sizeof(C));
			FFT(D, l), move(D, A1, B1, n);
			memcpy(D, x.v, n * sizeof(C));
			FFT(D, l), move(D, A2, B2, n);
			for (Int i = 0; i != n; ++i) {
				C B1B2 = B1[i] * B2[i];
				v[i] = A1[i] * A2[i] + B1B2 +
				       I * (A1[i] * B2[i] + A2[i] * B1[i] + B1B2);
			}
			FFT(v, l, -1);
			Int s = 0;
			for (Int i = 0; i != n; ++i) {
				v[i + 1] += C{F(R(v[i].real()) / 10), F(R(v[i].imag()) / 10)};
				v[i] = C{F(R(v[i].real()) % 10), F(R(v[i].imag()) % 10)};
				if (Int(v[i].imag()))s = i + 1;
			}
			size = s;
		}
} Ans, A;
void pow(Int x) {
	for (; x; x >>= 1, A *= A)
		if (x & 1)Ans *= A;
}
int main() {
	Int n;
	std::cin >> n;
	Ans.E();
	A.v[0] = {0, 1};
	pow(n);
	for (Int i = Ans.size - 1; i != -1; --i)
		putchar('0' + Int(Ans.v[i].imag()));
	return 0;
}
```