automatic processing of numerical results in `Plot`

First I want to solve an equation F(x,y)=0F(x,y)=0F(x,y)=0 for yyy by supplying a value of xx. (suppose obtaining the analytic form of y(x)y(x) is too difficult) Then I want to plot root yy (numerically calculated) as a function of xx by using the following:

Plot[y /. FindRoot[…/.{x->x0},{y, 0.2}],{x0, 0, 1}]

and I got something like the following

I omit … here since it is terribly complicated. Update: F(x,y) is of the form


and mm and nn can be as high as 19, which basically makes Solve impractical.

The result is satisfying except a small number of points near 1.0. Setting another initial starting value of yy in FindRoot might be an option but it is very tedious and often I cannot find a value of yy that fits the whole range of xx.

My question is: suppose I stick with the initial value of yy, is there a way to just eliminate that anomaly point after I Plot the numerical results? Or is there a better way to deal with this kind of numerical problem in general?




For some functions you can use the previously found root to start looking for the next. But I can’t know how your function behaves since it isn’t supplied.
– C. E.
Oct 23 ’13 at 7:58



@Anon,sounds like a good idea, could you write some code on it? For the function, it is rational function very similar to a Pade approximant. But it is too complicated to be readable so I didn’t put it in the question.
– wdg
Oct 23 ’13 at 8:13



@wdg You could try to use the form FindRoot[…,{x,xstart,xmin,xmax}] and restrict the search to [0,1]. Maybe this helps.
– halirutan
Oct 23 ’13 at 8:16



@halirutan, it does help!
– wdg
Oct 23 ’13 at 8:22


3 Answers


Using NDSolve to create an interpolation often works well. Hard to tell if it will work with your function.

With[{f = x^2 y + y^3 – 1/3},
sol = NDSolveValue[{Dt[f == 0, x] /. y -> y[x],
y[0] == (y /. FindRoot[f /. x -> 0, {y, 0.6}])}, y, {x, 0, 1}]]
(* InterpolatingFunction[{{0., 1.}}, <>] *)

Plot[sol[x], {x, 0, 1}]



A very nice trick!
– wdg
Oct 23 ’13 at 14:09

In case it helps to determine the new initial search value from the old, which is sometimes the case, it can be done like this:

f[x_, y_] := ArcTan[x + y]
x0 = 0;
Plot[(x0 = FindRoot[f[x, y], {x, x0}][[1, 2]]), {y, -100, 100}]



This method works well if the value of y is changing slowly from step to step, for example when creating a table of values for ListPlot. However, Plot will sample the function at y values which hop around all over the place as it goes through the recursive refinement process, so I don’t think this method is reliable.
– Simon Woods
Oct 23 ’13 at 9:00



@SimonWoods You are absolutely right! I stored the y-values in this example using Sow and the largest differences were very large which is not something I had anticipated. It’s still better than using a static initial value, I think, but as you say it’s better to use Table and ListPlot to do it this way.
– C. E.
Oct 23 ’13 at 12:31

Like this:?

Here I give F(x,y)= Abs@x + Abs@(x + y)-1

Plot[y /. Solve[Abs@x + Abs@(x + y) == 1, y, Reals] //
Evaluate, {x, -2, 2}, AspectRatio -> Automatic]

Another methods:

Using the ContourPlot command

ContourPlot[Abs[x] + Abs[x + y] == 1, {x, -2, 2}, {y, -2, 2},
AxesOrigin -> {0, 0}, Frame -> False, Axes -> True]



ContourPlot is a good idea. In fact I often use ContourPlot to verify numerical results obtained through other means.
– wdg
Oct 23 ’13 at 13:21