Mathematica 中的动态规划：如何自动定位和/或清除记忆函数的定义

Question

在 Mathematica 8.0 中，假设我有一些常量：

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

我想用它们来评估一些相互关联的功能

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

但这真的很慢，需要动态编程，否则你会得到指数级的减速：

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

现在它真的很快，但是如果我们想要更改常量（例如，在 Manipulate 函数中使用它），我们必须每次都使用 Clear g 和 h。如果我们有复杂的相互依赖关系，每次我们想要来自g 和h 的值时都将它们全部清除可能真的很烦人。

是否有一种简单的方法可以在 Module 或 Block 或类似名称中运行 g 和 h，以便在每次评估时都能得到新的结果而不会出现指数级减速？或者甚至是一种快速的方法来以一种很好的方式为g 和h 建立一个结果表？如前所述，我希望能够在 Manipulate 函数中计算 g 和 h。

Answer 1

这是一种方法，按照您的要求使用Block：

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];

我们将使用它来定义f 和h 为

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

我们现在打电话：

In[1246]:= g[20,10]//Timing
Out[1246]= {0.,253809.}

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

每次调用后都没有留下全局记忆定义 - Block 在执行退出之前小心销毁它们 Block。特别是，这里我将参数改一下，再次调用：

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}

此方案的替代方案是将所有参数（如 a,b,c,d,e,f）显式传递给函数，扩展它们的形式参数列表（签名），但这有一个缺点，即对应于不同过去参数值的旧记忆定义不会自动清除。这种方法的另一个问题是生成的代码将更加脆弱，w.r.t 参数数量的变化等。

编辑

但是，如果您想构建一个结果表，这可能会更快一些，因为您只需一次性完成，并且在这种情况下，您确实希望保留所有已记忆的定义。所以，这里是代码：

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;

你调用它，显式传递参数：

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

（我使用的是原始参数）。您可以在此方法中检查记忆定义是否保留在全局规则库中。下次您调用具有完全相同参数的函数时，它将获取记忆定义而不是重新计算。除了我上面概述的这种方法的问题之外，您还应该注意内存使用情况，因为没有任何东西被清除。

Answer 2

使用辅助符号记忆

可以修改问题中展示的记忆技术，以便在需要清除缓存时无需重新建立g 和h 的定义。想法是将记忆值存储在辅助符号上，而不是直接存储在 g 和 h 上：

g[0,k_] = 0;
g[t_,0] = e;
g[t_,k_] := memo[g, t, k] /. _memo :> (memo[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_] := memo[h, t, k] /. _memo :> (memo[h, t, k] = h[t-1,k]*c+g[t-1,k-1]*d)

这些定义与g 和h 的原始记忆版本基本相同，只是引入了一个新符号memo。有了这些定义，就可以简单地使用Clear@memo 清除缓存——无需重新定义g 和h。更好的是，可以通过将memo 放在Block 中来本地化缓存，因此：

Block[{memo, a = 7, b = 9, c = 13, d = 0.002, e = 2, f = 1}
, Table[g[t, k], {t, 0, 100}, {k, 0, 100}]
]

退出块时缓存被丢弃。

使用建议

进行记忆

原始和修订的记忆技术需要在函数 g 和 h 内进行侵入性更改。有时，在事后引入记忆是很方便的。一种方法是使用advising 的技术——一种类似于OO 编程中的子类化的函数式编程。 particular implementation of function advice 经常出现在 StackOverflow 的页面中。然而，该技术也是侵入性的。让我们考虑另一种向g 和h 添加建议而不改变它们的全局定义的技术。

诀窍是在Block 中临时重新定义g 和h。重新定义将首先检查缓存中的结果，如果失败，则从块外部调用原始定义。让我们回到 g 和 h 的原始定义，它们完全不知道记忆：

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

此技术的基本架构如下所示：

Module[{gg, hh}
, copyDownValues[g, gg]
; copyDownValues[h, hh]
; Block[{g, h}
  , m:g[a___] := m = gg[a]
  ; m:h[a___] := m = hh[a]
  ; (* ... do something with g and h ... *)
  ]
]

引入了临时符号gg 和hh 来保存g 和h 的原始定义。然后g 和h 在本地重新绑定到新的缓存定义，这些定义根据需要委托给原始定义。这是copyDownValues的定义：

ClearAll@copyDownValues
copyDownValues[from_Symbol, to_Symbol] :=
  DownValues[to] =
    Replace[
      DownValues[from]
    , (Verbatim[HoldPattern][from[a___]] :> d_) :> (HoldPattern[to[a]] :> d)
    , {1}
    ]

为了保持这篇文章的简短（呃），这个“复制”功能只关注向下值。一般建议工具还需要考虑上值、子值、符号属性等。

这种通用模式很容易实现自动化，即使很乏味。下面的宏函数memoize 做到了这一点，几乎没有评论：

ClearAll@memoize
SetAttributes[memoize, HoldRest]
memoize[symbols:{_Symbol..}, body_] :=
  Module[{pairs, copy, define, cdv, sd, s, m, a}
  , pairs = Rule[#, Unique[#, Temporary]]& /@ symbols
  ; copy = pairs /. (f_ -> t_) :> cdv[f, t]
  ; define = pairs /. (f_ -> t_) :> (m: f[a___]) ~sd~ (m ~s~ t[a])
  ; With[{ temps = pairs[[All, 2]]
         , setup1 = Sequence @@ copy
         , setup2 = Sequence @@ define }
    , Hold[Module[temps, setup1; Block[symbols, setup2; body]]] /.
        { cdv -> copyDownValues, s -> Set, sd -> SetDelayed }
    ] // ReleaseHold
  ]

经过一番折腾，我们现在可以从外部对g 和h 的非缓存版本进行记忆：

memoize[{g, h}
, Block[{a = 7, b = 9, c = 13, d = .002, e = 2, f = 1}
  , Table[g[t, k], {t, 0, 100}, {k, 0, 100}]
  ]
]

综上所述，我们现在可以创建一个响应式Manipulate 块：

Manipulate[
  memoize[{g, h}
  , Table[g[t, k], {t, 0, tMax}, {k, 0, kMax}] //
      ListPlot3D[#, InterpolationOrder -> 0, PlotRange -> All, Mesh -> None] &
  ]
, {{tMax, 10}, 5, 25}
, {{kMax, 10}, 5, 25}
, {{a, 7}, 0, 20}
, {{b, 9}, 0, 20}
, {{c, 13}, 0, 20}
, {{d, 0.002}, 0, 20}
, {{e, 2}, 0, 20}
, {{f, 1}, 0, 20}
, LocalizeVariables -> False
, TrackedSymbols -> All
]

LocalizeVariables 和 TrackedSymbols 选项是 g 和 h 对全局符号 a 到 f 的依赖关系的产物。