# How to get control the precision of contours in ContourPlot?

I want to get equal value contour of a function. But I want the curve to be precise.

A naive way would be using ContourPlot, then extract the contour data. To do this, we need a utility to extract contours like below

ClearAll[extractContourAssoc];
extractContourAssoc[contourPlot_] := Module[{},
tooltipData =
Cases[Normal@contourPlot, Tooltip[x__] -> {x}, Infinity];
contourValueList = tooltipData[[;; , -1]];
contourDataList =
Table[Cases[i, Line[pts__] -> pts, Infinity], {i, tooltipData}];


Then for example

func[x_, y_] := x^2 + y^2;
cplot = ContourPlot[func[x, y] == 1, {x, 0, 3}, {y, 0, 3}];


to extract the contour in cplot , we use

pts = First@First@extractContourAssoc@cplot;


However, the precision is limited, for example, it we substitute back the coordinate of curves.

func @@ # - 1 & /@ pts
(*{-0.0025047, -0.00131712, -0.00137604, -0.00225737, -0.00143394, \
-0.00116087, -0.00113379, -0.00280045, -0.00107967, -0.00161928, \
-0.00143129, -0.000948526, -0.00279327, -0.00140771, -0.000650118, \
-0.000687417, -0.000736281, -0.000528992, -0.00106009, -0.000248064, \
-0.00138826, -0.00156147, -0.000895443, -0.00156147, -0.000757219, \
-0.00138826, -0.000248064, -0.00106009, -0.00126354, -0.000736281, \
-0.00190364, -0.00140771, -0.00034149, -0.000606891, -0.000948526, \
-0.00143129, -0.00161928, -0.00107967, -0.00280045, -0.000860253, \
-0.00113379, -0.00116087, -0.00143394, -0.00225737, -0.00137604, \
-0.00131712, -0.000643397}*)


If the curve is precise enough, then the numbers should all be zero.

Is there a way to get more precise contours? Though there is MaxRecursion option, but it seems that there is no precision control in ContourPlot.

update

Sorry I found I made a mistake about extractContourAssoc. Line[pts__] -> pts should be Line[pts__] :> pts

ClearAll[extractContourAssoc];
extractContourAssoc[contourPlot_] := Module[{},
tooltipData =
Cases[Normal@contourPlot, Tooltip[x__] -> {x}, Infinity];
contourValueList = tooltipData[[;; , -1]];
contourDataList =
Table[Cases[i, Line[pts__] :> pts, Infinity], {i, tooltipData}];


And MaxRecursion actually works. For example, Set MaxRecursion->10 will have precision roughly about 5x10^-6. So the problem is solved.

• you can take your contourplot results and use as initial points for FindRoot Mar 4, 2018 at 13:53

Why not use an NDSolve based approach instead? For instance:

eqn = x[t]^2 + y[t]^2 == 1;
{xsol, ysol, diff} = Quiet @ NDSolveValue[
{D[eqn, t], x'[t]^2 + y'[t]^2 == 1, x[0]==1, y[0]==0},
{x, y, x[t]^2+y[t]^2-1},
{t, 0, 2 Pi}
];

ParametricPlot[{xsol[t], ysol[t]}, {t, 0, 2Pi}]


I included the difference you were interested in as an output of the NDSolveValue call:

Plot[diff, {t, 0, 2Pi}]


Finally, if you extract the Line object from the parametric plot, you will see that its errors are also less than $10^{-5}$ as well.

• Hi, @Carl Woll. It seems that your method needs to know at least one precise point at first, am I right? Feb 2, 2018 at 9:54

A typical way is to polish with FindRoot. Use the contour for one equation and the normal line for the second.

extractContourAssoc::noimpr = "No improvement possible in contour tagged ";
extractContourAssoc[contourPlot_, {f_, x_, y_}] :=
Module[{grad = D[f, {{x, y}}], contour},
tooltipData = Cases[Normal@contourPlot, Tooltip[t__] :> {t}, Infinity];
contourValueList = tooltipData[[;; , -1]];
contourDataList = Table[
Cases[
i,
Line[pts__] :> Block[{},
(* check the Tooltip tag to see how to use it *)
If[Head[i[[-1]]] === Equal,   (* tag is an equation with held variables *)
contour = i[[-1]] /. Equal -> Subtract // ReleaseHold,
If[NumericQ[i[[-1]]],        (* tag is numerical value *)
contour = f - i[[-1]],
Message[extractContourAssoc::noimpr, i[[-1]]];
contour = $Failed ]]; (* improve the position if possible *) If[contour =!=$Failed,
{x, y} /.
pts]
],
Infinity],
{i, tooltipData}];


Example:

func[x_, y_] := x^2 + y^2;
cplot = ContourPlot[func[x, y] == 1, {x, 0, 3}, {y, 0, 3}];


The one-argument form of extractContourAssoc is the OP's function:

pts = First@First@extractContourAssoc@cplot;  (* OP's extractContourAssoc *)
func @@ # - 1 & /@ pts // RealExponent // Histogram


The two-argument form is the one above that attempts to improve the position of the points with FindRoot:

pts = First@First@extractContourAssoc[cplot, {func[x, x], x, y}]; (* improved *)
func @@ # - 1 & /@ pts // RealExponent // Histogram[#, {1}] &