编辑:这是对我上一个答案的更正
首先根据所描述的条件写出差分方程:
f[0] = a
f[1] = a + b
f[2] = f[1] + f[0]
= 2a + b = a + b + a
f[3] = f[2] + f[1] = f[1] + f[0] + f[1]
= 3a + 2b
= a + b + a + a + b
f[4] = f[3] + f[2]
= 3a + 2b + 2a + b = 5a + 3b
f[5] = f[4] + f[3]
= 5a + 3b + 3a + 2b = 8a + 5b
f[6] = f[5] + f[4]
= 8a + 5b + 5a + 3b = 13a + 8b
...
f[n] = f[n-1] + f[n-2]
如果我们分开a和b,我们实际上可以简化这个问题:
f_a[n] = a*(f[n-1] + f[n-2]) with f[0] = 1 and f[1] = 1
f_b[n] = b*(f[n-1] + f[n-2]) with f[0] = 0 and f[1] = 1
现在,如果我们计算差分方程的解,我们应该得到以下假设 s=sqrt(5) 那 n \in N(是一个自然数):
w1a = ((1+s)/2)ˆ{n+1}
w2a = ((1-s)/2)ˆ{n+1}
w1b = ((1+s)/2)ˆn
w2b = ((1-s)/2)ˆn
f_a[n] = (1/s) * [w1a - w2a] * a
f_b[n] = (1/s) * [w1b - w2b] * b
简化:
l = (1+s)/2
g = (1-s)/2
f[n] = f_a[n] + f_b[n]
= (1/s) * [lˆn(al+b) - gˆn(ag+b)]
您可以在此处找到有关如何求解微分方程的更多信息:https://www.cl.cam.ac.uk/teaching/2003/Probability/prob07.pdf
您可以在 Python 函数中实现这些等式以获得该函数的任何值。
from math import sqrt
def f(n, a, b):
s = sqrt(5)
l = (1+s)/2
g = (1-s)/2
fn = (1/s) * ((a*l + b) * (l**n) - (a*g + b) * (g**n))
return int(round(fn, 0))
迭代搜索索引
如果您应用对数函数,您现在可以找到求解该方程的特定 f(n) 的 n(请参阅下面的部分)。但是,如果时间复杂度对您来说不是问题,并且考虑到 f[n] 对于 n 呈指数增长(这意味着在达到或超过 100k 之前您不需要进行太多搜索),您也可以简单地找到 n通过执行以下搜索为给定的a 和b 提供f[n]:
def search_index(a, b, value):
n = 0
while(True):
fn = f(n, a, b)
if fn == value:
return n
elif fn > value:
return -1
else:
n += 1
def brute_search(range_a, range_b, value):
for a in range(range_a + 1):
for b in range(range_b + 1):
if (a == 0) and (b == 0):
a = 1
res = search_index(a, b, value)
if res != -1:
return a, b, res
return -1
brute_search(1000, 1000, 100000)
>>> (80, 565, 12) # a = 80, b = 565 and n = 12
通过这种(相当糟糕的)方法,我们发现对于a=80 和b=565,n=12 将返回f_n = 100k。如果您想为a 和b 的值范围找到所有可能的解决方案,您可以通过以下方式修改brute_search:
def brute_search_many_solutions(range_a, range_b, value):
solutions = []
for a in range(range_a + 1):
for b in range(range_b + 1):
if (a == 0) and (b == 0):
a = 1
res = search_index(a, b, value)
if res != -1:
solutions.append((a, b, res))
return solutions
解析解
将先前的差分方程f_n 转换为现在n 是a、b 和f_n 的函数,我们得到:
n \aprox log((f_n * s) / (a * l + b)) / log(l)
此结果是一个近似值,可以指导您的搜索。您可以通过以下方式使用它:
def find_n(a, b, value):
s = sqrt(5)
l = (1+s)/2
g = (1-s)/2
return int(round(log(value * s / (a * l + b)) / log(l), 0))
def search(a, b, value):
n = find_n(a, b, value)
sol = f(n, a, b)
if sol == value:
return(a, b, n)
elif sol > value:
for i in range(n-1, 0, -1):
sol = f(i, a, b)
if sol == value:
return(a, b, i)
elif sol < value:
return(-1, 'no solution exits for a={} and b={}'.format(a, b))
else: # this should probably never be reached as find_n should
# provide an upper bound. But I still need to prove it
i = n
while(sol < value):
i += 1
sol = f(i, a, b)
if sol == value:
return(a, b, i)
elif sol > value:
return(-1, 'no solution exits for a={} and b={}'.format(a, b))
search(80, 565, 100000)
>>> (80, 565, 12) # a = 80, b = 565 and n = 12
注意:我很想在这里使用 LaTeX 的数学符号,但不幸的是我没有找到一种简单的方法来做到这一点......很可能,这个问题更适合 Stack Exchange,而不是 Stack Overflow。