104. 首尾均为全数字的斐波那契数(Pandigital Fibonacci ends)

斐波那契数列可以通过以下递归关系来定义:

$$ F_n=F_{n-1}+F_{n-2},where\ F_1=1\ and\ F_2=1. $$
可以发现\(F_{541}\)总共包含\(113\)位数字,而且是第一个结尾九个数字构成全数字的斐波那契数(全数字即包含一至九九个数字,但不一定是按照这个顺序),而\(F_{2749}\)包含\(575\)位数且是第一个开关九个数字是全数字的斐波那契数。

假设\(F_k\)是第一个结尾与开头九个数字均为全数字的斐波那契数,求\(k\)的值。

分析:这道题有一个显而易见的暴力解法,即逐一求出斐波那契数,把它们转换成字符串,分别截取前九位和后九位数字,判断它们是否是全数字,如果均为全数字则为题目所求。但是因为涉及到的斐波那契数都很大,多达几百上千位,这个暴力解法的效率较低,很难在一分钟时间内得到答案,我们需要一个更聪明的办法,技巧在于我们并不需要知道整个斐波那契数是多少,只需要知道前九位和后九位就可以了,从而可以避免大数运算提高计算效率。下面我们分别考虑前九位和后九位数字:

一、斐波那契数的后九位

我们知道要得到一个数字的后九位,只需要用这个数字对\(10^9\)取余即可,所以要求斐波那契数后九位,我们可对期递推关系式两侧分别对\(10^9\)取余,有:

$$ F_n \% 10^9=(F_{n-1}\%10^9+F_{n-2}\%10^9)\%10^9 $$

我们可把\(F_n\%10^9\)看作是一个新的数列,假设为\(MF_n\),则有:

$$ MF_n=(MF_{n-1}+MF_{n-2})\%10^9 $$

这个新数列中的每一项均为对应的斐波那契数对\(10^9\)取余的结果,因此我们只需要计算这个数列,就可以知道每一个斐波那契数的后九位,然后就可以判断它是否为全数字了。由于我们每次都只是对一个不超过九位数的数字进行运算,所以运算的效率要更高。具体编程时,我们仍然可以使用python语言中变量交换的方式来更新每一项。

二、斐波那契数的前九位

取一个数的前几位的思路也很简单,我们只需要将小数点向前移,移动恰当的位置后取整即可。举例来说,我们要取12345数字的前三位,我们需要把小数点往前移两位,小数点需要移动的位数等于数字总的位数减去需要取的前几位数字,比如12345总的位数是五位,需要取前三位,则小数点需要移动的位数等于五减去三等于二。小数向前移动实际上就是用数字除以十的乘方,比如需要向前移两位就除以十的二次方。

我们当然可以算出每个斐波那契数然后用上面的方式取出它前九位,但同样会涉及大数运算,效率很低。在第二十五题中,我们提到了斐波那契数的近似值,即当\(n\)足够大时,有:

$$ F_n\approx \frac{\varphi^n}{\sqrt{5}},where\ \varphi=\frac{1+\sqrt{5}}{2} $$

这里的近似误差随着\(n\)的增大而迅速减小,当\(n\ge4\)时,误差在10%以内;当\(n\ge8\)时,误差缩小到1%以内。对于这道题目涉及到的数量级,误差基本可以忽略不计。则对\(F_n\)取以十为底的对数,我们有:

$$ d=log(F_n)=log(\frac{\varphi^n}{\sqrt{5}})=nlog(\varphi)-\frac{1}{2}log(5) $$

\(F_n\)的位数\(w=int(d)+1\),则要求它的前九位,有:

$$ fnd=10^d/10^{w-9}=10^d/10^{int(d)+1-9}=10^{d-int(d)+8} $$

然后我们对\(fnd\)取整,即为斐波那契数的前九位。

三、结论

因此,在实际计算斐波那契数的前九位和后九位时,我们并不需要计算出整个斐波那契数。对于后九位,使用上面介绍的取余运算,不断迭代计算各个斐波那契数的后九位,并检验其是否为全数字。如果后九位是全数字,再把其对应的\(k\)值带入到上面的公式中计算斐波那契数的前九位,如果也为全数字,则为题目所求。代码如下:

# time cost = 456 ms ± 4.23 ms

from math import log10

def fib_top_digits(k):
    phi = (1+5**0.5)/2
    d = k*log10(phi) - 0.5*log10(5)
    res = int(pow(10,d-int(d)+8))
    return res

def main(N=10**9):
    is_pandigital = lambda n:set(str(n)) == set('123456789')
    a,b,k = 1,1,2
    while True:
        a,b = b,(a+b)%N
        k = k + 1
        if is_pandigital(b) and is_pandigital(fib_top_digits(k)):
            return k