POJ 3070 - 快速矩阵幂求斐波纳契数列

这题并不复杂。

 

设$A=\begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}$

 

由题中公式：

$\begin{pmatrix}
f(n+1) & f(n)\\ 
f(n+1) & f(n-1)
\end{pmatrix} = {\begin{pmatrix}
1 & 1 \\ 
1 & 0
\end{pmatrix}}^{n}$

可知，若要求f(n)只要求矩阵A的n次方即可。设我们所需的矩阵为$Answer$.

 

对于此题，我们可以先将$Answer$矩阵设置为$E$。

再求出${A}^{{2}^{0}}$、${A}^{{2}^{1}}$、${A}^{{2}^{2}}$ ... ${A}^{{2}^{m}}$ (${2}^{m}\leq n <{2}^{m+1}$)

其中，后一个矩阵为前一个矩阵的平方。

把他们储存起来。

对上述矩阵从后到前遍历。

当遍历到第i项时，若${2}^{i} \leq n$，则将$Answer$矩阵与此矩阵项相乘，积作为新的$Answer$矩阵。然后，将$n$减去${2}^{i}$，再接着遍历下一项，直至$n = 0$。

遍历结束后的$Answer$矩阵即为我们所需的矩阵。

 

 1 #include <cstddef>
 2 #include <cstdio>
 3 #include <cstring>
 4 
 5 struct matrix {
 6     unsigned m[2][2];
 7 };
 8 
 9 #define multiply(a,b,r) (((r)[0][0]=(a)[0][0]*(b)[0][0]+(a)[0][1]*(b)[1][0]),((r)[0][1]=(a)[0][0]*(b)[0][1]+(a)[0][1]*(b)[1][1]),((r)[1][0]=(a)[1][0]*(b)[0][0]+(a)[1][1]*(b)[1][0]),((r)[1][1]=(a)[1][0]*(b)[0][1]+(a)[1][1]*(b)[1][1]))
10 
11 int fibo_mod_by_10000(unsigned int n) {
12     if (n == 0)
13         return 0; 14 unsigned int mask = 0u, m = 0u; 15 16 while ((mask & n) != n) { 17 mask <<= 1u; 18 mask += 1u; 19 ++m; 20  } 21 22 matrix * ms = new matrix[m + 1u]; 23 ms[1u].m[0][0] = 1u; 24 ms[1u].m[0][1] = 1u; 25 ms[1u].m[1][0] = 1u; 26 ms[1u].m[1][1] = 0u; 27 28 for (unsigned int i = 1u; i < m; ++i) { 29 multiply(ms[i].m, ms[i].m, ms[i + 1].m); 30 ms[i + 1].m[0][0] %= 10000u; 31 ms[i + 1].m[0][1] %= 10000u; 32 ms[i + 1].m[1][0] %= 10000u; 33 ms[i + 1].m[1][1] %= 10000u; 34  } 35 36  matrix result, tmp; 37 memcpy(&result, &(ms[m]), sizeof(matrix)); 38 n -= (1u << (m - 1u)); 39 40 while (n != 1u && n != 0u) { 41 while ((1u << (m - 1u)) > n) 42 --m; 43 memcpy(&tmp, &result, sizeof(matrix)); 44  multiply(tmp.m, ms[m].m, result.m); 45 result.m[0][0] %= 10000u; 46 result.m[0][1] %= 10000u; 47 result.m[1][0] %= 10000u; 48 result.m[1][1] %= 10000u; 49 n -= (1u << (m - 1u)); 50  } 51  unsigned r; 52 delete[] ms; 53 if (n == 1u) 54 return result.m[0][0]; 55 else 56 return result.m[0][1]; 57 } 58 59 int main() 60 { 61 int i; 62 while ((scanf("%d", &i)), (i != -1)) 63 printf("%d\n", fibo_mod_by_10000(i)); 64 return 0; 65 }

 

转载于:https://www.cnblogs.com/Hoyy/p/6705664.html