4

Everyone knows that the Fibonacci sequence goes

F[0]=1, F[1]=1, F[2]=2, F[3]=3, F[4]=5, F[5]=8,

with F[n] = F[n-1]+F[n-2].

Now, how do you compute a number in the Fibonacci sequence when taken modulo 1000000007 = 10^9+7?

Needs to run as efficiently as possible, and in Python language :)

For example F[10**15] should take less than a second or so.

I know matrix exponentiation works, but how do you correct Matrix Exponentiation to reflect MODULO? (Another example, see https://www.nayuki.io/page/fast-fibonacci-algorithms)

Nayuki
  • 17,911
  • 6
  • 53
  • 80
William Yi
  • 69
  • 1
  • 4
  • Precompute the whole thing and store it in a table? :-) – Platinum Azure Oct 18 '14 at 16:31
  • 1
    This is not a duplicate. Even an answer proving how the modulus might be insignificant and then falling back to the usual suspects adds unique value to this question. Also added +1 to cancel downvote. – Platinum Azure Oct 18 '14 at 16:34
  • 5
    @PlatinumAzure your vote is not meant to compensate for the opinion of downvoters. Upvote only if you found a question/answer useful, not just to make justice. I didn't downvote, FWIW. – Stefano Sanfilippo Oct 18 '14 at 16:37
  • Another method is Doubling Method for Fibonacci, but idk how to fix that to MOD – William Yi Oct 18 '14 at 16:48
  • Python has exact "big" integer math, so integer overflow is not a problem. What is wrong with the simplest, most obvious approach? – Paul Oct 24 '14 at 05:06
  • The only way this is at all interesting is if you want `F[10**15]%M` where `F[0]=1,F[1]=1,F[2]=2,...` is the Fibonacci sequence. If you only want `F[0],...,F[k] s.t. F[k]<10**15`, that is trivial. – Paul Oct 24 '14 at 05:36
  • related: http://stackoverflow.com/questions/1525521/nth-fibonacci-number-in-sublinear-time – Paul Oct 24 '14 at 05:39
  • 1
    The output and recurrence relationship is unclear to me. Do you want `output = F[10**15] mod 1000000007 ` or do you want to apply modulo to the recurrence relationship so that, say, `G[n] = (G[n-1]+G[n-2]) mod 1000000007` where I have moved from F to G because G is apparently only piecewise related to Fibonacci. Also, you might ask a different version of this question as a math question on http://math.stackexchange.com – Paul Oct 24 '14 at 05:48
  • 1
    another modulo 1000000007 question: http://stackoverflow.com/questions/4874688/what-can-i-do-to-improve-my-fibonacci-number-generator – Paul Oct 25 '14 at 00:20

1 Answers1

2

The tricks needed:

1) Use the closed form of Fibonacci numbers, this is much faster than recursion. http://mathworld.wolfram.com/FibonacciNumber.html (formula 6)

2) modulo essentially factors over multiplication and addition, and sort of factors over division (you have to compute the multiplicative inverse in the mod space with the extended euclidean algorithm first) so you can basically just mod as you go. https://en.wikipedia.org/wiki/Modulo_operation#Equivalencies https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Modular_integers

code:

def rootiply(a1,b1,a2,b2,c):
    ''' multipy a1+b1*sqrt(c) and a2+b2*sqrt(c)... return a,b'''
    return a1*a2 + b1*b2*c, a1*b2 + a2*b1

def rootipower(a,b,c,n):
    ''' raise a + b * sqrt(c) to the nth power... returns the new a,b and c of the result in the same format'''
    ar,br = 1,0
    while n != 0:
        if n%2:
            ar,br = rootiply(ar,br,a,b,c)
        a,b = rootiply(a,b,a,b,c)
        n /= 2
    return ar,br

def rootipowermod(a,b,c,k,n):
    ''' compute root powers, but modding as we go'''
    ar,br = 1,0
    while k != 0:
        if k%2:
            ar,br = rootiply(ar,br,a,b,c) 
            ar,br = ar%n,br%n
        a,b = rootiply(a,b,a,b,c)
        a,b = a%n, b%n
        k /= 2
    return ar,br

def fib(k):
    ''' the kth fibonacci number'''
    a1,b1 = rootipower(1,1,5,k)
    a2,b2 = rootipower(1,-1,5,k)
    a = a1-a2
    b = b1-b2
    a,b = rootiply(0,1,a,b,5)
    # b should be 0!
    assert b == 0
    return a/2**k/5

def powermod(a,k,n):
    ''' raise a**k, modding as we go by n'''
    r = 1
    while k!=0:
        if k%2:
            r = (a*r)%n
        a = (a**2)%n
        k/=2
    return r

def mod_inv(a,n):
    ''' compute the multiplicative inverse of a, mod n'''
    t,newt,r,newr = 0,1,n,a
    while newr != 0:
        quotient = r / newr
        t, newt = newt, t - quotient * newt
        r, newr = newr, r - quotient * newr
    if r > 1: return "a is not invertible"
    if t < 0: t = t + n
    return t

def fibmod(k,n):
    ''' compute the kth fibonacci number mod n, modding as we go for efficiency'''
    a1,b1 = rootipowermod(1,1,5,k,n)
    a2,b2 = rootipowermod(1,-1,5,k,n)
    a = a1-a2
    b = b1-b2
    a,b = rootiply(0,1,a,b,5)
    a,b = a%n,b%n
    assert b == 0
    return (a*mod_inv(5,n)*mod_inv(powermod(2,k,n),n))%n

if __name__ == "__main__":
    assert rootipower(1,2,3,3) == (37,30) # 1+2sqrt(3) **3 => 13 + 4sqrt(3) => 39 + 30sqrt(3)
    assert fib(10)==55
    #print fib(10**15)%(10**9+7) # takes forever because the integers involved are REALLY REALLY REALLY BIG
    print fibmod(10**15,10**9+7) # much faster because we never deal with integers bigger than 10**9+7
Him
  • 5,257
  • 3
  • 26
  • 83