Mathematica Stack Exchange is a question and answer site for users of Mathematica. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

Let's say I have the NDSolve example for documentation involving splitting 2nd order into set of 1st order ODEs:

NDSolve[{x'[t] == y[t], y'[t] == -x[t], x[0] == 1, y[0] == 0}, {x, 
y}, {t, 0, 10}]

I need to solve my second order system like this for various reasons.

My ODE is more complicated, and depends on some other parameters call them $(a,b)$, and so when I solved it as a second order system I used memoization in the form below, so that if a solution for a given parameter set had already been computed, Mathematica wouldn't redo it all again. I would like to do this, or something like it for the first order system and where NDSolve is giving the two variables $x[a,b],y[a,b]$.

 x[a,b] = x /.NDSolve[{eq[a,b],x[0]== init[a,b],x'[0]==dinit[a,b],x,{r, 0, 10}][[1]];
share|improve this question
Have you, by any chance, seen the NDSolve`StateData tutorial? – J. M. Sep 30 '12 at 1:38
up vote 5 down vote accepted

When I see this right, you have one single function call (NDSolve) depending on two parameters which results in a vector (the interp. functions). Here is an equivalent minimal example function for this call

f[a_, b_] := {b, a}

The underlying issue is now, that you always get both answers from f (x and y) even if you need only x for a specific choice of a and b. You can solve this by memorizing both, x and y, no matter which function you call.

x[a_, b_] := {x[a, b], y[a, b]} = f[a, b];
y[a_, b_] := {x[a, b], y[a, b]} = f[a, b];

A quit check shows you that you memorize the values for x and y even if you compute only x

Table[x[a, b], {a, 0, 2}, {b, 0, 2}];

Lets try this with your NDSolve example. I only wrapped Module around to localize x and y.

ClearAll[x, y, f];
f[a_, b_] := Module[{x, y}, {x, y} /. 
   First@NDSolve[{x'[t] == y[t], y'[t] == -x[t], x[0] == a, y[0] == b}, 
    {x, y}, {t, 0, 10}]];

x[a_, b_] := {x[a, b], y[a, b]} = f[a, b];
y[a_, b_] := {x[a, b], y[a, b]} = f[a, b];

Table[x[a, b], {a, 0, 2}, {b, 0, 2}];


Finally, the repetition of the exact same function definition for x and y looks kind of like we are not experts. Lets improve this. The only thing which is different is the head of the left side in the :=. So we write the line one time and substitute x and y resp.

SetDelayed[#[a_, b_], {x[a, b], y[a, b]} = f[a, b]] & /@ {x, y}
share|improve this answer
one small problem with this is on the initial call, if I do say: x[1,2][3.0] instead of getting a numeric answer back I get two interpolating functions. The second call works fine, and gives the result desired. This is in contrast to standard memoization. I think it can be solved by just changing to x[a_, b_] := ({x[a, b], y[a, b]} = f[a, b]; x[a,b]) etc – fpghost Sep 30 '12 at 12:14
@fpghost Yes, I haven't addressed that. I only worked out the problem with the memoization. I'll have a look at it later. – halirutan Sep 30 '12 at 14:08

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.