Plotting the components of a function that returns a list in different colors without redundant evaluations of the function



I have a function f which takes a number as input, and returns a list of numbers (the length of the list is constant). f is hard to calculate (each evaluation takes a long time).

I want to plot the different components of f in different colors.

If I use this command:

Plot[f[x], {x, -2, 2}]

all the lines are drawn in the same color.
If I use this command:

Plot[{f[x][[1]], f[x][[2]], f[x][[3]]}, {x, -2, 2}]

(assuming the list has three components) the lines are drawn in different colors, but the function is called three times the necessary amount.

Note that this is a numeric function, it cannot be evaluated with a symbolic argument (i.e. the function definition begins with f[x_Real]:=), so there is no use in using Evaluate like in this question.


Posted 2012-06-05T11:31:36.493

Reputation: 1 451

5The following code suggests that even Plot[f[t],...] evaluates f multiple times. f[x_Real] := (i++; {x, x + 1, x + 2}) and i = 0; Plot[f[t], {t, 0, 1}, PlotPoints -> 10, MaxRecursion -> 0], then Print[i]. The result is 31, not 10-ish. So, solutions just using Plot may not work as well. – Yu-Sung Chang – 2012-06-05T12:47:03.710

@Yu-SungChang, sorry, I didn't see your comment. I wrote an answer which basically states the same. I give you +1 for your comment ;-) – halirutan – 2012-06-05T13:04:24.720

Relevant SO question: "How to select the “best” new point when sampling a near-parabolic function?"

– Alexey Popkov – 2012-06-07T08:33:23.260

Another related SO question: "Telling Plot to style vector-valued black-box functions in Mathematica"

– Joe – 2012-11-15T08:29:52.220



At risk of stating the obvious, if you are willing to give up the adaptive sampling, exclusions, etc. of Plot you could use ListLinePlot:

f[x_?NumericQ] := x + Mod[x, {1, 2, 3}]

ListLinePlot[Transpose@Table[f[x], {x, 0, 10, 0.01}], PlotStyle -> Thick] 

Mathematica graphics

Better I think is to restyle the Graphics data produced by plot, as Heike did for Plotting piecewise function with distinct colors in each section and which I refactored in my answer. Applied here:

Module[{i = 1},
  Plot[f[x], {x, 0, 10}, PlotStyle -> Thick] /.
    x_Line :> {ColorData[1][i++], x}

enter image description here

Even nicer is Simon Woods' method which styles the plot while it is created, posted in answer to:

Also useful and very interesting is the solution by wxffles in the follow-up question:


Posted 2012-06-05T11:31:36.493

Reputation: 259 163

+1 for the Module[{i=0}...ColorData[1][i++] alone! I'd have spent hours fiddling with MapIndexed's... – Aisamu – 2014-04-28T12:13:17.733

2(+1) Because MMA uses adaptive sampling for each component of the function separately, memoizing doesn't usually save a lot of evaluation: the arguments used to plot the first component only slightly overlap the arguments used to plot the second component, etc. (Sow the values of x with an EvaluationMonitor and then ListPlot them afterwards to see this.) Thus, this appears to be the only really effective workaround when good plotting arguments (values of $x$) are known in advance. – whuber – 2012-06-05T14:21:11.987

1@whuber: My thoughts exactly. It seems that there is no way to make MMA's adaptive algorithm take into account all the values of a certain component that were already calculated as a result of evaluating f for the other components, and that's a shame. If evaluating f three (or more) times the necessary amount is more costly than the evaluations that the adaptive algorithm saves us (which is almost always the case for functions that take reasonably long to evaluate), this answer seems to be the best workaround. Accepted. – Joe – 2012-06-06T08:32:18.077

@Joe thanks for the Accept. Sorry I didn't have something better for you. – Mr.Wizard – 2012-06-06T08:35:48.897


Please consider the following simple example

i = 0;
f[x_?NumericQ] := {x^2, 2 x^2, x^2 + 1};
Plot[f[x], {x, 0, 1}, EvaluationMonitor :> i++];

(* Out[30]= 471 *)

and now the same code for only the first component of f

i = 0;
Plot[x^2, {x, 0, 1}, EvaluationMonitor :> i++];

(* Out[33]= 157 *)

If you now recall, that $3\cdot 157 = 471$ you see that even in your Plot[f[x],...] call the function is called more often than one would naively expect. Therefore, I claim you can skip hard thinking about this issue.


Posted 2012-06-05T11:31:36.493

Reputation: 109 574

You are correct. However, if I still want to plot all the components of f with a number of evaluations more close to 157, I need to use the approach described in @Mr.Wizard's answer.

– Joe – 2012-06-06T08:38:07.990


You could memoize the function:

f[x_] := f[x] = x + Mod[x, {1, 2, 3}]

Plot[{f[x][[1]], f[x][[2]], f[x][[3]]}, {x, 0, 10}, PlotStyle -> Thick]

Mathematica graphics

If you don't want to use memory globally for this or memoize for all calls to f you could make a copy:

  g[x_] := g[x] = f[x];
  Plot[{g[x][[1]], g[x][[2]], g[x][[3]]}, {x, 0, 10}, PlotStyle -> Thick]


Posted 2012-06-05T11:31:36.493

Reputation: 259 163

This is a nice workaround in case there is no way to address the problem directly, and it has the slight disadvantage of memory use. I would like to wait and see if someone comes up with a more direct approach. – Joe – 2012-06-05T11:55:56.427

@Joe please do wait; I'm hoping someone has a perfect solution. – Mr.Wizard – 2012-06-05T12:00:12.647

From @whuber's comment, I learned that memoizing doesn't actually save function evaluations.

– Joe – 2012-06-06T08:42:56.377

1@Joe it can save evaluations, but not as many as reliably as I had thought. f[x_Real] := (i++; x + {Sin[x], Sqrt[x], Mod[x, 3]}) then:

i = 0; Plot[{f[x][[1]], f[x][[2]], f[x][[3]]}, {x, 0, 10}, PlotStyle -> Thick]; i Compare: i = 0; Module[{g}, g[x_] := g[x] = f[x]; Plot[{g[x][[1]], g[x][[2]], g[x][[3]]}, {x, 0, 10}, PlotStyle -> Thick]]; i – Mr.Wizard – 2012-06-06T09:19:44.570

You're right @Mr.Wizard, in your example memoizing saves evaluations, but it does not save two thirds of the evaluations (I got 715 instead of 953). Of course we can't expect exactly two thirds, since there may be some overhead, by I think that even asymptotically we won't get two thirds. – Joe – 2012-06-06T09:52:10.737


If your functions are continuous, you may use something like (naive example follows):

f[x_Real] := {2 x, 3 x};
Plot[f[x], {x, -2, 2}] /. {Hue[x__], Line[y__], Line[z__]} -> {Hue[x],
    Line[y], Hue[0, 1, 1], Line[z]}

enter image description here

Dr. belisarius

Posted 2012-06-05T11:31:36.493

Reputation: 112 848

Doesn't work well on my example, does it? (This was one of the first things I tried.) – Mr.Wizard – 2012-06-05T12:23:01.040

I also tried to get really clever and use Plot[{##&[], Sequence @@ f[x], ##&[]} . . . but that's no good at all. Maybe you can find a way to make it work. – Mr.Wizard – 2012-06-05T12:24:37.357

@Mr.Wizard Works only for continuous functions, as I said ... – Dr. belisarius – 2012-06-05T12:25:59.823