Die milliardste Fibonacci-Zahl berechnen

Nach dem Titel fragst du dich vielleicht, warum man sich so etwas überhaupt antut. Ganz einfach: weil ich es kann. Manchmal muss man sich selbst herausfordern und etwas Ambitioniertes ausprobieren. Auch wenn es erst einmal simpel klingt, gibt es dabei viele Herausforderungen. Du kannst nicht einfach irgendeinen Algorithmus nehmen und ihn eine Milliarde Iterationen laufen lassen. Je nach Methode würde das Jahre dauern, vielleicht sogar Jahrzehnte. Ein weiteres Problem ist das Speichern der Zahl. Die milliardste Fibonacci-Zahl ist so absurd groß, dass sie menschlich kaum vorstellbar ist, und selbst Computer haben damit zu kämpfen. Also wollte ich sehen, ob ich sie berechnen kann. Aber nicht nur das: Es sollte auch schnell sein. Also legen wir los.

Zuerst reden wir kurz darüber, wie Fibonacci-Zahlen funktionieren. Falls du das Konzept noch nicht kennst: Es ist ziemlich einfach. Jede Fibonacci-Zahl entsteht als Summe der beiden vorherigen Fibonacci-Zahlen. Ab jetzt definieren wir $\text{Fib}(n)$ als die $n^{\text{th}}$ Fibonacci-Zahl. Mathematisch können wir das so schreiben:

\text{Fib}(n) = \text{Fib}(n-1) + \text{Fib}(n-2)

Die folgenden Gleichungen sind die Basisfälle der Fibonacci-Folge. Sie sind immer gleich:

\text{Fib}(0) = 0 \quad \text{and} \quad \text{Fib}(1) = 1

Die Rechnungen selbst sind also recht simpel. Sobald wir aber in den Milliardenbereich gehen, kann selbst so eine einfache Gleichung extrem teuer werden. Dazu kommt das kleine Problem der Computerdatentypen. Meistens sind sie auf eine bestimmte Größe begrenzt. Ein normaler Integer ist zum Beispiel nur 32 Bit lang. Das reicht gerade einmal für die Zahl 4 Milliarden. Da Fibonacci-Zahlen mit einer exponentiellen Rate wachsen, funktionieren unsere klassischen Datentypen für dieses Vorhaben leider nicht.

Es gibt mehrere Wege, dieses Problem anzugehen. Der erste wäre, eigene Datentypen mit beliebiger Größe zu implementieren. Das ist aber ziemlich viel Kopfschmerz, und ich bin faul. Wir könnten dafür auch eine kleine Library namens GMP nutzen (mehr dazu später). Für den Anfang gehen wir aber noch einfacher vor. Die perfekte Lösung ist hier Python. Python unterstützt Integer beliebiger Größe nativ. Der normale int-Typ reicht also für unser Programm. Ich nehme das außerdem als Gelegenheit zu zeigen, dass Python zwar langsamer sein kann als viele Sprachen, mit den richtigen Anpassungen aber für die meisten Vorhaben mehr als schnell genug ist.

Bevor wir starten: Das Ziel ist nicht nur, die milliardste Fibonacci-Zahl zu berechnen, sondern sie in relativ kurzer Zeit zu berechnen. Was genau „kurz“ heißt, besprechen wir später. Ich kann aber schon sagen, dass ich nicht jahrelang auf das Ergebnis warten will. Noch ein wichtiges Detail: Alle Benchmarks und Performance-Daten hier wurden auf meinem M2 Pro MacBook Pro mit 16 GB RAM ausgeführt. Damit können wir uns jetzt ein paar Ansätze anschauen, um die Fibonacci-Folge zu berechnen.

Rekursiv

Die erste und simpelste Lösung ist, den Code genau wie die Gleichung zu schreiben. Gleichzeitig ist das mit Abstand die schlechteste Lösung. Aber schauen wir uns zuerst den Code an:

pythondef fib(n):
    if n <= 1:
        return n
    return fib(n-1) + fib(n-2)

Es ist eine sehr einfache und schöne rekursive Funktion. So wie sie hier steht, ist sie in vielen Fällen aber nicht einmal eine echte Lösung. Für $n = 43$ braucht das Programm zum Beispiel schon über 40 Sekunden. Warum? Bei jedem Funktionsaufruf ruft die Funktion sich selbst zweimal auf, und jeder dieser Aufrufe tut das wieder. Das führt zu exponentiellem Wachstum, sodass wir insgesamt $2^n$ Operationen ausführen müssen. Mathematisch hat die Methode eine Laufzeitkomplexität von $O(2^n)$ . Bei einer Eingabe von 43 ruft sie sich also $8.796.093.022.208$ Mal auf. Du kannst dir vorstellen, wie $2^{1.000.000.000}$ aussehen würde. Da 43 ein kleines Stück von einer Milliarde entfernt ist, funktioniert diese Methode für uns nicht.

Performance der rekursiven Methode

Iterativ

Etwas effizienter und trotzdem relativ leicht zu verstehen. Der iterative Ansatz beseitigt einen der kritischsten Flaschenhälse des rekursiven Ansatzes, indem er $\text{fib}(0)$ und $\text{fib}(1)$ am Anfang definiert und dann einfach die letzten beiden berechneten Fibonacci-Zahlen speichert und aufsummiert. Hier ist der Code für diese Lösung:

pythondef fib(n):
    a, b = 0, 1  # Initialize the first two Fibonacci numbers
    for _ in range(n):
        a, b = b, a + b
    return a

Und ja, diesmal ist es tatsächlich eine Lösung. Sie ist zwar immer noch unglaublich langsam, würde die milliardste Fibonacci-Zahl aber irgendwann erreichen. Starke Betonung auf „irgendwann“. Im Vergleich zur rekursiven Lösung hat sie immerhin eine Laufzeitkomplexität von „nur“ $O(n)$ . Nicht großartig, aber eine massive Verbesserung gegenüber der rekursiven Methode. Sie bräuchte nur eine Milliarde Iterationen, um unser Ziel zu erreichen. Da ich aber auch nicht wirklich Zeit habe, auf eine Milliarde Iterationen zu warten, betrachte ich diese Methode trotzdem als Fehlschlag.

Performance der iterativen Methode

Matrix-Exponentiation

Jetzt wird es interessant. Statt das Problem nur auf hoher Ebene anzuschauen, gehen wir etwas tiefer. Eine Eigenschaft der Fibonacci-Folge, die uns sehr hilft, ist, dass sie über Matrix-Exponentiation berechnet werden kann. Das lässt sich so darstellen:

\begin{bmatrix}1 & 1 \\1 & 0\end{bmatrix}^n = \begin{bmatrix}F_{n+1} & F_n \\F_n & F_{n-1}\end{bmatrix}

Darauf basierend können wir die $n^{\text{th}}$ Fibonacci-Zahl berechnen, indem wir die Matrix in die $n^{\text{th}}$ Potenz nehmen. Das Gute daran: Computer sind schnell bei Matrixoperationen. Naiv würde Matrix-Exponentiation aber immer noch $O(n)$ Matrixmultiplikationen brauchen. Wenn man berücksichtigt, dass jede einzelne Matrixmultiplikation mehrere Additionen und Multiplikationen enthält, landen wir sogar bei mehr Gesamtoperationen als beim iterativen Ansatz. Mit ein paar cleveren Tricks kann man Matrix-Exponentiation allerdings in $O(\log n)$ berechnen.

pythondef fibonacci_matrix(n):
    F = [[1, 1], [1, 0]]  # Base matrix for Fibonacci sequence
    result = matrix_pow(F, n)
    return result[0][1]  # The nth Fibonacci number is in this position
@staticmethod
def matrix_mult(A, B):
    return [
        [A[0][0] * B[0][0] + A[0][1] * B[1][0], A[0][0] * B[0][1] +
            A[0][1] * B[1][1]],
        [A[1][0] * B[0][0] + A[1][1] * B[1][0], A[1][0] * B[0][1] +
            A[1][1] * B[1][1]]
    ]
def matrix_pow(self, mat, p):
    result = [[1, 0], [0, 1]]  # Identity matrix
    while p > 0:
        if p & 1:  # If p is odd
            result = matrix_mult(result, mat)
        mat = matrix_mult(mat, mat)  # Square the matrix
        p >>= 1  # Divide p by 2
    return result

Das nennt man Fast Exponentiation. Wenn wir das als Exponentiationsalgorithmus nutzen, bekommen wir die $n^{\text{th}}$ Fibonacci-Zahl mit einer Laufzeit von $O(\log n)$ . Das ist deutlich schneller als alle bisherigen Versuche und würde uns die milliardste Zahl tatsächlich in einer halbwegs vernünftigen Zeit liefern.

Performance der Matrix-Exponentiation

Aber es geht noch besser...

Verbesserte Matrix-Exponentiation

Diese Methode basiert ähnlich wie oben weiterhin auf Matrix-Exponentiation. Wir nutzen aber ein paar wichtige Verbesserungen, um diese Exponentiationen noch weiter zu beschleunigen. Zuerst wird die Matrix diesmal etwas anders definiert. Auch wenn die Form anders aussieht, funktioniert sie mathematisch gleich.

\begin{bmatrix}0 & 1 \\1 & 1\end{bmatrix}^n = \begin{bmatrix}F_{n-1} & F_n \\F_n & F_{n+1}\end{bmatrix}

pythondef mat_fib(n):
    a, b, c, d = 0, 1, 1, 1  # Initial matrix [[0, 1], [1, 1]]
    x, z = 0, 1  # Initial vector [0, 1]
    while n:
        n, r = divmod(n, 2)
        if r:  # If n is odd
            x, z = a*x + b*z, c*x + d*z  # Multiply current matrix with vector
        if n:  # If n is not zero
            # Square the current matrix
            a, b, c, d = a*a + b*c, b*(a + d), c*(a + d), c*b + d*d
    return x  # x holds the nth Fibonacci number

Der zentrale Unterschied dieser Methode ist, dass Operationen auf einzelnen Elementen statt auf der ganzen Matrix ausgeführt werden. Einfach gesagt reduziert das Rechenaufwand und Speicherverbrauch und sollte die Cache-Effizienz erhöhen. Obwohl die Laufzeitkomplexität weiterhin $O(\log n)$ ist, sinkt die tatsächliche Laufzeit, weil das Programm gleich viele Iterationen ausführt, aber pro Iteration weniger berechnet und besser mit dem Cache arbeitet.

Performance-Vergleich von normaler und verbesserter Matrixmethode

Doubling Fibonacci

Was wäre, wenn wir die Anzahl der Berechnungen noch weiter reduzieren könnten? Auch diese Methode startet mit ein paar mathematischen Definitionen. Zuerst definieren wir zwei Identitäten, aus denen wir unsere Folge berechnen können.

Double-Index-Identität für gerade Zahlen:

F(2k) = F(k) \times \left[ 2F(k+1) - F(k) \right]

Double-Index-Identität für ungerade Zahlen:

F(2k+1) = F(k+1)^2 + F(k)^2

Damit können wir ein Problem im Grunde in zwei kleinere Teilprobleme zerlegen. Wir halbieren den Index und dadurch in jeder Iteration auch die nötigen Schritte. Wenn wir zum Beispiel $\text{Fib}(100)$ berechnen wollen, könnten wir das so schreiben:

F(100) = F(50) \times \left[ 2F(51) - F(50) \right]

Das heißt, wir müssen keine Zahlen zwischen 100 und 51 berechnen. Auch diese Methode braucht insgesamt noch $O(\log n)$ Iterationen. Im Vergleich zu den beiden vorherigen Methoden braucht sie aber weniger Operationen pro Iteration, weil sie Matrixoperationen komplett vermeidet. Dadurch ist sie insgesamt effizienter.

Der Code dafür sieht so aus:

pythondef fibonacci_doubling(n):
    def _fibonacci_doubling(n):
        if n == 0:
            return (0, 1)
        else:
            a, b = _fibonacci_doubling(n >> 1)
            c = a * (2 * b - a)
            d = a * a + b * b
            if n & 1:
                return (d, c + d)
            else:
                return (c, d)
    return _fibonacci_doubling(n)[0]

Eine kleine Randnotiz: Wie einige vielleicht bemerkt haben, arbeitet diese Methode rekursiv. Im Gegensatz zu unserem ersten Ansatz ist sie aber extrem effizient, weil sie jede Zahl nur einmal berechnet und sehr viele überspringt.

Performance-Vergleich von Doubling und verbesserter Matrixmethode

Und das ist sie. Das ist die Methode, mit der wir zur milliardsten Fibonacci-Zahl kommen. Sie ist mit Abstand die schnellste Methode, die wir haben. Für meinen Geschmack ist sie aber immer noch nicht ganz schnell genug. Und ich habe noch ein Ass im Ärmel. Bis hierhin haben wir Pythons native Integer beliebiger Größe genutzt. Die sind, wie du dir denken kannst, langsam. Um das abzufedern, nutzen wir eine starke Library namens GMP. Auch mit GMP können wir mit Integern beliebiger Größe arbeiten. GMP wurde aber ursprünglich in purem C geschrieben und enthält sehr viele mathematische und algorithmische Optimierungen. Den Code auf diese Library umzuschreiben ist recht einfach. Es braucht nur ein paar zusätzliche Deklarationen in jeder Klasse, ansonsten bleibt der Code größtenteils gleich. Schauen wir, ob es schneller ist.

Performance-Vergleich von GMP und Nicht-GMP-Doubling

Und wow. GMP lässt Python komplett hinter sich. Bei sehr kleinen $n$ ist es wegen des zusätzlichen Overheads beim Erzeugen der MPZ-Objekte langsamer, aber bei größeren $n$ holt es das deutlich wieder rein. Jetzt können wir die milliardste Fibonacci-Zahl in nur ~7,38 Sekunden berechnen. Was ich damit mache? Noch keine Ahnung. Es war aber cool zu sehen, wie ein relativ simples Problem wie Fibonacci auf so viele verschiedene Arten gelöst werden kann, von einfacher bis komplexer Mathematik. Spannend war auch zu sehen, wie schnell Python sein kann, wenn man richtig damit arbeitet. Ja, eine Lösung komplett in C wäre wahrscheinlich noch schneller, aber sie hätte mich auch 10x mehr Zeit zum Programmieren gekostet. Nur zum Spaß habe ich außerdem jeden Ansatz auf GMP portiert. Unten findest du den Graphen mit allen Algorithmen im Vergleich auf logarithmischer Skala.

Performance-Vergleich aller GMP-Methoden (log)

Weitere Ideen

Auch wenn ich dieses Projekt als Erfolg sehe, habe ich ein paar Ideen, die für ein noch schnelleres Programm interessant sein könnten:

Multithreading für Doubling Fibonacci: Man könnte potenziell Multithreading nutzen, weil im Grunde zwei Probleme gleichzeitig laufen. Ist der Overhead es wert? Wird es überhaupt schneller? Fragen für einen möglichen zweiten Besuch...
GPUs nutzen: Könnte man das mit GPUs beschleunigen? Was ist mit Tensor Cores, die speziell für Matrixmultiplikationen gedacht sind? Wäre das schneller als auf der CPU? Könnte man es weiter parallelisieren?
Weitere mathematische Optimierungen: Gibt es noch schnellere Wege, die Fibonacci-Folge zu berechnen? Könnte man andere mathematische Identitäten nutzen, um die Anzahl der Berechnungen weiter zu reduzieren?

GitHub

Falls du dieses Projekt weiter erkunden, eigene Methoden zum Vergleichen implementieren oder sonst etwas damit machen willst, schau dir gern den gesamten hier gezeigten Code und mehr auf meiner GitHub-Seite an.

Skills Improved

Python +35

C +25