Search code examples
recursionwolfram-mathematicadynamic-programmingexponential

Dynamic Programming in Mathematica: how to automatically localize and / or clear memoized function's definitions


In Mathematica 8.0, suppose I have some constants:


a:=7
b:=9
c:=13
d:=.002
e:=2
f:=1

and I want to use them to evaluate some interlinked functions



g[0,k_]:=0
g[t_,0]:=e
g[t_,k_]:=g[t-1,k]*a+h[t-1,k-1]*b

h[0,k_]:=0
h[t_,0]:=f
h[t_,k_]:=h[t-1,k]*c+g[t-1,k-1]*d

But this is really slow and in need of dynamic programming, or else you get an exponential slowdown:


g[0, k_] := 0
g[t_, 0] := e
g[t_, k_] := g[t, k] = g[t - 1, k]*a + h[t - 1, k - 1]*b

h[0, k_] := 0
h[t_, 0] := f
h[t_, k_] := h[t, k] = h[t - 1, k]*c + g[t - 1, k - 1]*d

Now it's really fast, but if we ever want to change the constants(say, to use this in a Manipulate function) we have to Clear g and h each time. If we had complex interdependencies, it might be really annoying to clear them all out every single time we wanted a value from g and h.

Is there an easy way to run g and h in a Module or Block or similar, so that I can get a fresh result back each time it's evaluated without the exponential slowdown? Or even a fast way to build up a table of results for both g and h in a nice way? As said, I want to be able to compute g and h in a Manipulate function.


Solution

  • Here is one way, using Block as you requested:

    ClearAll[defWrap];
    SetAttributes[defWrap, HoldFirst];
    defWrap[fcall_] :=
      Block[{g, h},
         (* Same defintions with memoization as you had, but within Block*)
    
         g[0, k_] := 0;
         g[t_, 0] := e;
         g[t_, k_] := g[t, k] = g[t - 1, k]*a + h[t - 1, k - 1]*b;   
         h[0, k_] := 0;
         h[t_, 0] := f;
         h[t_, k_] := h[t, k] = h[t - 1, k]*c + g[t - 1, k - 1]*d;
    
         (* Our function call, but within a dynamic scope of Block *)
         fcall];
    

    We will use this to give definitions to f and h as

    ClearAll[g, h];
    g[tt_, kk_] := defWrap[g[tt, kk]];
    h[tt_, kk_] := defWrap[h[tt, kk]];
    

    We call now:

    In[1246]:= g[20,10]//Timing
    Out[1246]= {0.,253809.}
    
    In[1247]:= h[20,10]//Timing
    Out[1247]= {6.50868*10^-15,126904.}
    

    There are no global memoized definitions left after each call - Block takes care to destroy them just before the execution exits Block. In particular, here I will change the parameters and call them again:

    In[1271]:= 
    a:=1
    b:=2
    c:=3
    d:=.01
    e:=4
    f:=5
    
    In[1279]:= g[20,10]//Timing
    Out[1279]= {0.015,0.808192}
    
    In[1280]:= h[20,10]//Timing
    Out[1280]= {0.,1.01024}
    

    An alternative to this scheme would be to explicitly pass all parameters like a,b,c,d,e,f to functions, extending their formal parameter lists (signatures), but this has a disadvantage that the older memoized definitions corresponding to different past parameter values would not be automatically cleared. Another problem with that approach is that the resulting code will be more fragile, w.r.t changes in the number of parameters, etc.

    EDIT

    However, if you want to build a table of results, this could be somewhat faster since you do it once and for all, and in this case you do want to keep all memoized definitions. So, here is the code:

    ClearAll[g, h];
    g[0, k_, _] := 0;
    g[t_, 0, {a_, b_, c_, d_, e_, f_}] := e;
    g[t_, k_, {a_, b_, c_, d_, e_, f_}] := 
         g[t, k, {a, b, c, d, e, f}] = 
            g[t - 1, k, {a, b, c, d, e, f}]*a +  h[t - 1, k - 1, {a, b, c, d, e, f}]*b;
    
    h[0, k_, _] := 0;
    h[t_, 0, {a_, b_, c_, d_, e_, f_}] := f;
    h[t_, k_, {a_, b_, c_, d_, e_, f_}] := 
         h[t, k, {a, b, c, d, e, f}] = 
            h[t - 1, k, {a, b, c, d, e, f}]*c +  g[t - 1, k - 1, {a, b, c, d, e, f}]*d;
    

    You call it, passing the parameters explicitly:

    In[1317]:= g[20,10,{a,b,c,d,e,f}]//Timing
    Out[1317]= {0.,253809.}
    

    (I was using the original parameters). You can check that the memoized definitions remain in the global rule base, in this method. Next time you call a function with exact same parameters, it will fetch the memoized definition rather than recompute. Apart from the problems with this approach that I outlined above, you should also watch for the memory usage, since nothing gets cleared.