Plot of function-generated lists: Why does this workaround not work?

Question

My problem is a classic one: I want to plot a function-generated list, and there's no way to pre-evaluate (because the actual function includes NMinimize and thus cannot be calculated on symbolic values). I use the following replacement for the actual function:

f[x_?NumericQ]:={x,x^2}

Now doing

Plot[f[x], {x, 0, 1}, Evaluated->True]

will still generate only one colour, as evaluation with a symbolic value will not generate a list. So I thought I could bypass this by just saving the result in a temporary variable and creating a wrapper function to read that out again (after all, unlike Mathematica, I do know the shape of the generated list):

Module[{fwrap1,fwrap2,tmp},
  fwrap1[x_?NumericQ] := (tmp = f[x])[[1]];
  fwrap2[x_?NumericQ] := tmp[[2]];
  Plot[{fwrap1[x], fwrap2[x]}, {x, 0, 1}]]

I expected this to work as follows:

• First, fwrap1 is called, which obtains the complete list from f, stores it in tmp and returns the first item.
• Then, fwrap2 is called, which just reads the second value which fwrap1 stored in tmp.

However while this indeed gives two colours as expected, the second function is replaced with a constant value (the one at maximal x).

By having fwrap2 increment a counter, I can verify that it is called 53 times, and by recording the minimal and maximal x, I can verify that it is indeed called on the complete interval (well, at least at both end points). Moreover, adding lines of the form

Print[{fwrap1[0.5],fwrap2[0.5]}];

inside the module (with different values each time) gives the expected output.

So why does this code not work?

Answer 1

I'm turning my previous comment into an answer. In your case, you assumed that, for each value of $x$, fwrap1 is called, then fwrap2.

However, from the behavior you observe I suspect that Plot calculates all values for fwrap1 first, repeatedly rewriting the value of tmp and discarding the intermediate values you wanted to save; only afterwards does it move on to fwrap2.

This would seem to make sense with the fact that only the value of tmp corresponding to one of the boundaries of the $x$ range is retained.

Answer 2

f[x_?NumericQ] := {x, x^2}

These three expressions all produce the same output. Quiet is used in each case to suppress the message "Part::partw: "Part 2 of f[x] does not exist."

Module[{fwrap1, fwrap2, tmp},
  fwrap1[x_] := (tmp = f[x])[[1]];
  fwrap2[x_] := tmp[[2]];
  Plot[{fwrap1[x], fwrap2[x]}, {x, 0, 1},
   Evaluated -> True,
   PlotLegends -> {x, x^2}]] //
 Quiet

Module[{fwrap1, fwrap2, tmp},
  fwrap1[x_] = (tmp = f[x])[[1]];
  fwrap2[x_] = tmp[[2]];
  Plot[{fwrap1[x], fwrap2[x]}, {x, 0, 1},
   PlotLegends -> {x, x^2}]] //
 Quiet

Module[{tmp},
  Plot[
   {(tmp = f[x])[[1]], tmp[[2]]},
   {x, 0, 1},
   Evaluated -> True,
   PlotLegends -> {x, x^2}]] //
 Quiet

EDIT:

Looking at the number of times at f[x] is evaluated

f[x_?NumericQ] := (n++; {x, x^2})

n = 0;
Module[{fwrap1, fwrap2, tmp}, fwrap1[x_] := (tmp = f[x])[[1]];
  fwrap2[x_] := tmp[[2]];
  Plot[{fwrap1[x], fwrap2[x]}, {x, 0, 1}, Evaluated -> True, 
   PlotLegends -> {x, x^2}]] // Quiet; n

(*  159  *)

n = 0;
Module[{fwrap1, fwrap2, tmp}, fwrap1[x_] = (tmp = f[x])[[1]];
   fwrap2[x_] = tmp[[2]];
   Plot[{fwrap1[x], fwrap2[x]}, {x, 0, 1}, PlotLegends -> {x, x^2}]] // 
  Quiet;
n

(*  159  *)

n = 0;
Module[{tmp}, 
   Plot[{(tmp = f[x])[[1]], tmp[[2]]}, {x, 0, 1}, Evaluated -> True, 
    PlotLegends -> {x, x^2}]] // Quiet;
n

(*  159  *)

n = 0; Plot[{f[x][[1]], f[x][[2]]}, {x, 0, 1}];
n

(*  238  *)