My program can't run that fast even with memoization

The other answer is the correct way to calculate the fibonacci sequence, indeed, but you should also know why your memoization wasn't working. To be specific:

fibonacci_cache = {}

This line being inside the function means you were emptying your cache every time fibonacci_memo was called.

You shouldn't be computing the Fibonacci sequence, not even by dynamic programming. Since the Fibonacci sequence satisfies a linear recurrence relation with constant coefficients and constant order, then so will be the sequence of their sums.

Definitely don't cache all the values. That will give you an unnecessary consumption of memory. When the recurrences have constant order, you only need to remember as many previous terms as the order of the recurrence.

Further more, there is a way to turn recurrences of constant order into systems recurrences of order one. The solution of the latter is given by a power of a matrix. This gives a faster algorithm, for large values of n. Each step will be more expensive, though. So, the best method would use a combination of the two, choosing the first method for small values of n and the latter for large inputs.

O(n) using the recurrence for the sum

Denote S_n=F_0+F_1+...+F_n the sum of the first Fibonacci numbers F_0,F_1,...,F_n.

Observe that

• S_{n+1}-S_n=F_{n+1}
• S_{n+2}-S_{n+1}=F_{n+2}
• S_{n+3}-S_{n+2}=F_{n+3}

Since F_{n+3}=F_{n+2}+F_{n+1} we get that S_{n+3}-S_{n+2}=S_{n+2}-S_n. So

S_{n+3}=2S_{n+2}-S_n

with the initial conditions S_0=F_0=1, S_1=F_0+F_1=1+1=2, and S_2=S_1+F_2=2+2=4.

One thing that you can do is compute S_n bottom up, remembering the values of only the previous three terms at each step. You don't need to remember all of the values of S_k, from k=0 to k=n. This gives you an O(n) algorithm with O(1) amount of memory.

O(ln(n)) by matrix exponentiation

You can also get an O(ln(n)) algorithm in the following way:

Call X_n to be the column vector with components S_{n+2},S_{n+1},S_{n}

So, the recurrence above gives the recurrence

X_{n+1}=AX_n

where A is the matrix

[
 [2,0,-1],
 [1,0,0],
 [0,1,0],
]

Therefore, X_n=A^nX_0. We have X_0. To multiply by A^n we can do exponentiation by squaring.

For the sake of completeness here are implementations of the general ideas described in @NotDijkstra's answer plus my humble optimizations including the "closed form" solution implemented in integer arithmetic.

We can see that the "smart" methods are not only an order of magnitude faster but also seem to scale better compatible with the fact (thanks @NotDijkstra) that Python big ints use better than naive multiplication.

import numpy as np
import operator as op
from simple_benchmark import BenchmarkBuilder, MultiArgument
B = BenchmarkBuilder()

def pow(b,e,mul=op.mul,unit=1):
    if e == 0:
        return unit
    res = b
    for bit in bin(e)[3:]:
        res = mul(res,res)
        if bit=="1":
            res = mul(res,b)
    return res

def mul_fib(a,b):
    return (a[0]*b[0]+5*a[1]*b[1])>>1 , (a[0]*b[1]+a[1]*b[0])>>1

def fib_closed(n):
    return pow((1,1),n+1,mul_fib)[1]

def fib_mat(n):
    return pow(np.array([[1,1],[1,0]],'O'),n,op.matmul)[0,0]

def fib_sequential(n):
    t1,t2 = 1,1
    for i in range(n-1):
        t1,t2 = t2,t1+t2
    return t2

def sum_fib_direct(n):
    t1,t2,res = 1,1,1
    for i in range(n):
        t1,t2,res = t2,t1+t2,res+t2
    return res
    
def sum_fib(n,method="closed"):
    if method == "direct":
        return sum_fib_direct(n)
    return globals()[f"fib_{method}"](n+2)-1

methods = "closed mat sequential direct".split()

def f(method):
    def f(n):
        return sum_fib(n,method)
    f.__name__ = method
    return f

for method in methods:
    B.add_function(method)(f(method))

B.add_arguments('N')(lambda:(2*(1<<k,) for k in range(23)))

r = B.run()
r.plot()

import matplotlib.pylab as P
P.savefig(fib.png)

I am not sure how you are taking anything near a second. Here is the memoized version without fanciness:

class fibs(object):
    def __init__(self):
        self.thefibs = {0:0, 1:1}

    def __call__(self, n):
        if n not in self.thefibs:
            self.thefibs[n] = self(n-1)+self(n-2)
        return self.thefibs[n]

dog = fibs()

sum([dog(i) for i in range(40) if dog(i) < 4000000])