0304 素波那契
* * * *
拉格朗日计划
* * * *
素波那契

给定正整数n,函数next_prime(n)的值是大于n的最小素数。

定义$a_1$为next_prime($10^14$),$a_n$为next_prime($a_{n-1}$)。

再定义$b_n$的定义为$f_{a_n}$,其中$f_n$是以$f_1=f_2=1$为初值的斐波那契数列。

求$\sum_{n=1}^{100000}b_n$模1234567891011的余数。

本题难度:



解答

为了减少码量,使用sympy的nextprime函数来生成$a_n$,$b_n$则直接用矩阵快速幂计算。

将斐波那契数列的递推公式写作 $$\begin{bmatrix}f_{n-1}\\f_n\end{bmatrix}=\begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix} \begin{bmatrix}f_{n-2} \\ f_{n-1}\end{bmatrix}$$ 记 $$A\begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix}$$ 用矩阵快速幂求出$A^{2^k}$缓存在列表B中,就可以把$a_n$写成二进制再计算矩阵乘法求出$b_n$,最终结果是$283988410192$。

from sympy import nextprime

target=10**14
bound=10**5
m=1234567891011

product=lambda a,b:[[sum(a[i][k]*b[k][j] for k in range(len(b)))%m for j in range(len(b[0]))] for i in range(len(a))]

B=[[[0,1],[1,1]]]
for i in range(len(bin(target+bound))+2):
    B.append(product(B[-1],B[-1]))

v=[[0],[1]]
r=0
p=0
x=target
for n in range(bound):
    x=nextprime(x)
    y=x-p
    i=0
    while y:
        if y%2==1:
            v=product(B[i],v)
        i+=1
        y>>=1
    r=(r+v[0][0])%m
    p=x
print(r)