Memoization in Mathematica

I wrote this short tutorial on memoization years ago. Since it turned out to be somewhat popular, and my old website is going away, I am reproducing it here.

Memoization is a very well known programming pattern in Mathematica. Memoization is an optimisation technique: it means making a function “remember” its last return value for a given input argument.

For example, let us look at the recursively defined Fibonacci function:

fibonacci[0] = 0;
fibonacci[1] = 1;
fibonacci[n_] := fibonacci[n-1] + fibonacci[n-2]

In[4]:= fibonacci[30] // Timing
Out[4]= {5.391, 832040}

It is easy to see that the execution time of this function depends exponentially on the argument, so above a threshold it gets very-very slow. Now let us change the function definition slightly: in Mathematica a function may redefine itself!

In[5]:= fibonacci[n_] := fibonacci[n] = fibonacci[n-1] + fibonacci[n-2]

In[6]:= fibonacci[30] // Timing
Out[6]= {2.50234*10^-16, 832040}

This new Fibonacci function is a lot faster since it needs to compute its return value only once for each input argument. On the downside, it never “forgets”, so this technique is not very practical with functions that act on real numbers (as the memory may fill up quickly).

The following function is an attempt to remedy this problem. It caches only up to 200 values. It is intended to be used with functions that take a single numerical quantity as their argument.

SetAttributes[memo, HoldAll]
SetAttributes[memoStore, HoldFirst]
SetAttributes[memoVals, HoldFirst]

In[4]:= memoVals[_]={};

memoStore[f_, x_] := With[{vals = memoVals[f]},
    If[Length[vals] > 200,
        f /: memoStore[f, First[vals]] =. ;
        memoVals[f] ^= Append[Rest[memoVals[f]], x]
        memoVals[f] ^= Append[memoVals[f], x]
    f /: memoStore[f, x] = f[x]

In[6]:= memo[f_Symbol][x_?NumericQ] := memoStore[f, x]

Note that memo[] associates the cached results with the function on which it is acting (the argument f). Therefore, when this function is Cleared, the cached results are removed too. If memoization has been used on a function, care must be taken to always Clear[] it before its definition is changed (otherwise the old cached return values will shadow the new ones). Example usage:

fun[a_?NumericQ] := First@Module[{f, x},
    f /. NDSolve[{f''[x] == -a f[x], f[0]==1, f'[0]==0}, f, {x,0,10}]]

In[9]:= Plot3D[fun[a][x], {a, 1, 2}, {x, 0, 10}] // Timing
Out[9]= {10.015, <removed plot>}

In[10]:= Plot3D[memo[fun][a][x], {a, 1, 2}, {x, 0, 10}] // Timing
Out[10]= {0.453, <removed plot>}

Comments !