Skip to main content

memoization - How to memoize a function with Options?


Consider the following modified fibonacci function:


tempHist = $History
$History = 0

ClearAll[fibon]
Options[fibon] = {k -> 1}
fibon[0] = 0;
fibon[1] = 1;
fibon[n_, OptionsPattern[]]:= fibon[n-1] + OptionValue[k] * fibon[n-2]

fibon[30] // Timing (* {4.111894, 832040} *)
fibon[30, k -> 2] // Timing (* {4.161676, 1149851} *)

The modification involves a coefficent that multiplies the second term of the recursive definition. Can the idiom f[x_]:=f[x]= ... be used to implement memoization in this case?



I tried with no luck the following naive approaches for c-aching fibon


ClearAll[cfibon]
Options[cfibon] = {k -> 1}
cfibon[0] = 0;
cfibon[1] = 1;
cfibon[n_, OptionsPattern[]]:= cfibon[n, OptionsPattern[]] = cfibon[n-1] + OptionValue[k] * cfibon[n-2]

cfibon[30] // Timing (* {2.1999999997746045*^-5, 832040} *)
cfibon[30, k -> 2] // Timing (* {2.899999999961267*^-5, 832040} *)


Something is being memoized, but this is not the wanted behavior (note that the two calls return the same value 832040, the second call is meant return 1149851). Here is another try:


ClearAll[cfibon]
Options[cfibon] = {k -> 1}
cfibon[0] = 0;
cfibon[1] = 1;
cfibon[n_, OptionsPattern[]]:= cfibon[n, k -> OptionValue[k]] = cfibon[n-1] + OptionValue[k] * cfibon[n-2]

cfibon[30] // Timing
(* {10.664260999999998, 832040} *)


cfibon[30, k -> 2] // Timing
(* {10.677505, 1149851} *)

Now, the results are as expected, by the computation is slower than the non-cached version! So...


How to memoize a function with Options?



Answer



The case at hand


Here is one possibility:


ClearAll[fibon]
Options[fibon] = {k -> 1}

fibon[0, OptionsPattern[]] = 0;
fibon[1, OptionsPattern[]] = 1;

fibon[n_, opts : OptionsPattern[]] /; ! OrderedQ[{opts}] :=
fibon[n, Sequence @@ Sort[{opts}]];
fibon[n_, opts : OptionsPattern[]] :=
fibon[n, opts] =
fibon[n - 1, opts] + OptionValue[k]*fibon[n - 2, opts]

So, for example:



fibon[30]

(* 832040 *)

fibon[30,k->2]

(* 357913941 *)

Note that I added a rule which sorts the options, and I explicitly assumed the two sets of options the same when they are only different by the order of options. This assumption could be improved by accounting for some options not present explicitly but being the same due to defaults.


General meta-programming solution



The above version has the following limitations:



  • Default option values not accounted for

  • Has to explicitly pass the options to inner calls to self

  • It is not general enough (which is ok for a specific problem, but not for a general approach). Although one can see that some of the code does not really depend on the specifics of the problem, the above form does not allow one to easily factor that out.


Here I will address these limitations of the previous version, and show how we can achieve a fairly complete automation of this task for a wide class of problems by using meta-programming and code generation.


Here is a code-generating function that would generate the full memoizing code from simpler version not involving options:


ClearAll[optionMemo];
SetAttributes[optionMemo, HoldAll];

optionMemo[sym_Symbol, {defs___}, defOptions : {___?OptionQ}] :=
Module[{},
ClearAll[sym];
Options[sym] = defOptions;
sym[args : Except[_?OptionQ] ..., opts : OptionsPattern[]] :=
With[{
newopts =
Sort @ DeleteDuplicates[
{opts}~Join~Options[sym],
First@#1 == First@#2 &]

},
sym[args, Sequence @@ newopts] /; newopts =!= {opts}
];
sym[args : Except[_?OptionQ] ..., opts : OptionsPattern[]]/;!OrderedQ[{opts}] :=
sym[args, Sequence @@ Sort[{opts}]];

sym[args : Except[_?OptionQ] ..., opts : OptionsPattern[]] :=
Block[{sym, options = Options[sym]},
Options[sym] = options;
defs;

sym[args]
]
];

It takes the name of the symbol, the list of code statements which produce the relevant definitions, and the list of options for this symbol. It works by generating the boilerplate code, but also it relies on a rather interesting technique where a function Block-s itself and creates some definitions locally inside Block. What we buy by this construct here:



  • We don't have to pass options to inner recursive calls

  • We can still globally call OptionVaulue as before, with a single argument, and this will work.

  • The memoization only happens during the computation, to speed it up, but all the memoized values are cleared by Block at the end, automatically;



So, here is how we use it: first generate the boilerplate code


optionMemo[
fibon,
{
fibon[0] = 0,
fibon[1] = 1,
fibon[m_] := fibon[m] = fibon[m - 1] + OptionValue[k]*fibon[m - 2]
},
{k -> 1, i -> 2}
]


and now use it as before:


fibon[30]

(* 832040 *)

fibon[30,k->2]

(* 357913941 *)


The optionMemo construct should be general enough to handle many similar problems.


Comments

Popular posts from this blog

plotting - Plot 4D data with color as 4th dimension

I have a list of 4D data (x position, y position, amplitude, wavelength). I want to plot x, y, and amplitude on a 3D plot and have the color of the points correspond to the wavelength. I have seen many examples using functions to define color but my wavelength cannot be expressed by an analytic function. Is there a simple way to do this? Answer Here a another possible way to visualize 4D data: data = Flatten[Table[{x, y, x^2 + y^2, Sin[x - y]}, {x, -Pi, Pi,Pi/10}, {y,-Pi,Pi, Pi/10}], 1]; You can use the function Point along with VertexColors . Now the points are places using the first three elements and the color is determined by the fourth. In this case I used Hue, but you can use whatever you prefer. Graphics3D[ Point[data[[All, 1 ;; 3]], VertexColors -> Hue /@ data[[All, 4]]], Axes -> True, BoxRatios -> {1, 1, 1/GoldenRatio}]

plotting - Filling between two spheres in SphericalPlot3D

Manipulate[ SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, Mesh -> None, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], {n, 0, 1}] I cant' seem to be able to make a filling between two spheres. I've already tried the obvious Filling -> {1 -> {2}} but Mathematica doesn't seem to like that option. Is there any easy way around this or ... Answer There is no built-in filling in SphericalPlot3D . One option is to use ParametricPlot3D to draw the surfaces between the two shells: Manipulate[ Show[SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], ParametricPlot3D[{ r {Sin[t] Cos[1.5 Pi], Sin[t] Sin[1.5 Pi], Cos[t]}, r {Sin[t] Cos[0 Pi], Sin[t] Sin[0 Pi], Cos[t]}}, {r, 1, 2 - n}, {t, 0, Pi}, PlotStyle -> Yellow, Mesh -> {2, 15}]], {n, 0, 1}]

plotting - Mathematica: 3D plot based on combined 2D graphs

I have several sigmoidal fits to 3 different datasets, with mean fit predictions plus the 95% confidence limits (not symmetrical around the mean) and the actual data. I would now like to show these different 2D plots projected in 3D as in but then using proper perspective. In the link here they give some solutions to combine the plots using isometric perspective, but I would like to use proper 3 point perspective. Any thoughts? Also any way to show the mean points per time point for each series plus or minus the standard error on the mean would be cool too, either using points+vertical bars, or using spheres plus tubes. Below are some test data and the fit function I am using. Note that I am working on a logit(proportion) scale and that the final vertical scale is Log10(percentage). (* some test data *) data = Table[Null, {i, 4}]; data[[1]] = {{1, -5.8}, {2, -5.4}, {3, -0.8}, {4, -0.2}, {5, 4.6}, {1, -6.4}, {2, -5.6}, {3, -0.7}, {4, 0.04}, {5, 1.0}, {1, -6.8}, {2, -4.7}, {3, -1....