# Cissoid of Diocles not smooth plotted by ContourPlot

Table[ContourPlot[y^2 (10-x)==x^3,{x,0,.1},{y,-.1,.1},PlotPoints->5],{k,1,150,10}] How to get the smooth graph?

When PlotPoints is larger than 500, performance is too slow, however still not good.  • increase your plot points? btw. for what do you need k? Jun 27 '13 at 9:50
• @Stefan PlotPoints -> 500 two slow, however still not good ---! x=0, y should be 0. Jun 27 '13 at 9:54
• i've used 50. still too slow? Jun 27 '13 at 9:57
• @Stefan but when using PlotPoints->50 output is not the right graph. Performance too slow. Jun 27 '13 at 10:00
• Try an odd number of PlotPoints. You might also try MaxRecursion->6 or so. Also, I wouldn't generally expect ContourPlot, which samples points on a 2D grid, to perform as well as Plot, which samples points on a 1D grid. Jun 27 '13 at 10:18

First, I would not generally expect ContourPlot, which samples points on a 2D grid, to perform as well as Plot, which samples points on a 1D grid. That said, you can improve the plot by ensuring that you sample the function at the origin by simply using an odd number of PlotPoints:

ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
PlotPoints -> 25, MaxRecursion -> 6] To make it clear what's going on here, let's zoom into your image and show the sample points.

{pic, {pts}} = Reap[ContourPlot[
y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
PlotPoints -> 50,EvaluationMonitor :> Sow[{x, y}],
ContourStyle -> {Thick, Blue}]];
Show[{pic, Graphics[{Point[pts], PointSize[Large], Red, Point[{0, 0}]}]},
PlotRange -> -0.05] The origin is the large red dot and it is simply missed. Now, here's the same picture with PlotPoints set to 25 and MaxRecursion left at it's default value of 2. • Nice answer, but what's a 1D grid? Jun 27 '13 at 10:46
• @cormullion Is grid not an $n$-dimensional concept for arbitrary, positive, integer $n$? As in GridGraph[{5}] or Grid[{{1,2,3,4,5}}]. Jun 27 '13 at 10:51

ContourPlot generally works on the principle of the Intermediate Value Theorem, and works best when it can find points for which the function being plotted lies on either side of the contour level (which is hard, for example, for equations of the form f^2 == 0). For convenience I will talk about contours of the form f == 0 and speak about the function changing sign instead of taking values above and below the contour level.

ContourPlot uses PlotPoints to divide the domain initially into rectangular cells. ContourPlot decides which cells contain part of the contour and which contain connected components; in the adaptive refinement step (up to MaxRecursion steps), it also decides which might contain more of the contour, and which are not worth pursuing. In a symmetric graph on a symmetric domain, whether there are an even or an odd number of PlotPoint can make a remarkable difference. So can offsetting the domain to make it slightly asymmetrical.

It might be surprising that one can produce an excellent graph with just PlotPoints -> 2:

ContourPlot[y^2 (10 - x) == x^3, {x, 0, 0.1}, {y, -0.1, 0.1},
PlotPoints -> 2, MaxRecursion -> #] & /@ {2, 4} // GraphicsRow This works for an even number of PlotPoints up to 14 and then fails for 16 or higher, as in the OP. I will return to this later. First I wish illustrate another point made above. The connectedness of the contour seems (naturally) related to the values of the function on the sample grid. If a function does not change sign, then the cell contains a disconnected segment of the contour or no segment. This can change as sample points are added under adaptive refinement. This can be illustrated by manually adding extra sample points (as discussed in this question. Here is the OP's 50-point plot and the same plot with an extra sample point added to each cell straddling the x-axis:

ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
PlotPoints -> {{50}, #}, MaxRecursion -> 0, ImageSize -> 150,
BaseStyle -> Tiny, PlotRange -> {{0, 0.04}, Automatic},
PlotLabel -> Row[{Length@#, " extra points"}]] & /@
{{}, Table[{x0, 0}, {x0, 0.05/49, 0.1, 0.1/49}]} // GraphicsRow Adaptive refinement slowly connects the components generated by the extra sample points to the curve (MaxRecursion is set to 2, 4, 6):

Table[ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
PlotPoints -> {{50}, Table[{x0, 0}, {x0, 0.05/49, 0.1, 0.1/49}]},
MaxRecursion -> mr, PlotLabel -> mr, ImageSize -> 150,
BaseStyle -> Tiny, PlotRange -> {{0, 0.04}, Automatic}],
{mr, 2, 6, 2}] // GraphicsRow I will try to explain what is happening below.

For an even number of PlotPoints up to 14, a very large proportion of the plot domain is empty (all x less than 0.08 or rectangle at least 80% of the area). There must be some threshold that triggers ContourPlot to subdivide the whole domain. In cases of 2, 4,..., 14 points, the first subdivision effectively adds the points corresponding to PlotPoints being set to 3, 5,..., 15 points resp. In the OP's case, a low number of plot points happens to work well because of the symmetry of the plot.

The case of PlotPoints -> 16

The function adaptivePlotGrid[plotpoints, xMax, maxRecursion] is a tool for exploring the adaptive refinement procedure (code below). It shows the orginal sample grid (thin brown lines), the new sample points (red) and the previously generated sample points (light gray) for a succession of refinements.

Here we can see the adaptive refinement working its way back to the origin for PlotPoints -> 16. When an original cell is subdivided at y == 0, the (approximate) contour in that cell is incorporated into the contour plot. We can see that the origin is first incorporated into the plot at step MaxRecursion -> 13; subsequent refinements result in better approximations to the contour.

adaptivePlotGrid[16, 0.09, 15] The case of PlotPoints -> 50

The OP's case is similar to the case of PlotPoints -> 50, but it has many more cells to work through to get back to the origin. The origin is still not reached when the limit of MaxRecursion -> 15 is reached.

Aside: An issue with adding extra sample points is discussed in the linked question. The cells with an extra sample point in them (apparently) will not be subdivided until an adjacent cell is subdivided and the shared side is split. When this happens, the components of the contour on the adjacent cells are joined (assuming the changes in sign of the function at the sample points indicate they join, as in the OP's case).

The case of PlotPoints -> 2

The case of PlotPoints -> 2 may interesting, too. Here we see the initial graph is blank. The first refinement adds the origin (as other points, on and off the x-axis, not shown due to the PlotRange). Then the adaptive refinement is on its way to plotting a more and more accurate contour.

adaptivePlotGrid[2, 0.04, 3] adaptivePlotGrid[pp_, maxX_: 0.04, maxMR_: 5] :=
Module[{pts, mr1 = 0, mr2 = maxMR, dmr = 1},
pts[mr1 - dmr] = {};
Partition[
Table[
{pic, {pts[mr]}} =
Reap[ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
PlotPoints -> pp, MaxRecursion -> mr,
EvaluationMonitor :> Sow[{x, y}], PlotLabel -> mr,
PlotRange -> {{-0.002, maxX}, {-0.015, 0.015}},
ImageSize -> 150, BaseStyle -> Tiny]];
Show[{Graphics[{LightBrown,
Line[Table[{{x0, -0.1}, {x0, 0.1}}, {x0, 0., 0.1, 0.1/(pp - 1)}]],
Line[Table[{{0., y0}, {0.1, y0}}, {y0, -0.1, 0.1, 0.2/(pp - 1)}]],
PointSize[Tiny], LightGray, Point[pts[mr - dmr]],
Red, Point[Complement[pts[mr], pts[mr - dmr]]]}],
pic,
Graphics[{PointSize[Tiny],
Red, Point[Cases[Complement[pts[mr], pts[mr - dmr]], {_, 0.}]]}]},
Options[pic]],
{mr, mr1, mr2, dmr}],
4, 4, 1, {}] // GraphicsGrid
]

• Note: There seems to be a bug in V10 using PlotLabel with BaseStyle -> Tiny. (Reported to WRI.) Jan 2 '15 at 19:32

I knew about MaxBend for Plot and ParametricPlot, but I've never thought of applying it to ContourPlot until I saw this answer today by Mr.Wizard.

Setting MaxBend to a very low value will force ContourPlot to refine the curve as many times as allowed by MaxRecursion.

ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
Method -> {MaxBend -> 0.01}] It is already an improvement over the default plot. Increase MaxRecursion from 2 to 3 and the plot is quite good.

ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
Method -> {MaxBend -> 0.01}, MaxRecursion -> 3] The advice about using an appropriate number (odd, in this case) for PlotPoints, which has been pointed out by Mark McClure, is excellent. But perhaps it is worth pointing out that MaxBend can overcome the deficiency inherent in a small, even number of PlotPoints:

ContourPlot[y^2 (10 - x) == x^3, {x, 0, .1}, {y, -.1, .1},
Method -> {MaxBend -> 0.01}, PlotPoints -> 10] Yet another method in V10, "BoundaryMeshGenerator" -> "Continuation" in the FEM mesh generation functionality. It's quite simple if one is satisfied with a RegionPlot instead of a ContourPlot.

Needs["NDSolveFEM"]

mesh = ToElementMesh[
ImplicitRegion[y^2 (10 - x) <= x^3, {{x, 0, 1/10}, {y, -1/10, 1/10}}],
"BoundaryMeshGenerator" -> "Continuation"];
RegionPlot[mesh] Or the curve can be extracted from the mesh elements and the spurious boundary at the edge of the plot range can be deleted.

Show[
DeleteCases[
Normal@ mesh["Wireframe"[
"MeshElement" -> "BoundaryElements",
"MeshElementStyle" ->
Directive[RGBColor[0.368417, 0.506779, 0.709798],
AbsoluteThickness[1.6]]]],
Line[{{x1_ /; Abs[x1 - 0.1] < 10^-8, _Real},
{x2_ /; Abs[x2 - 0.1] < 10^-8, _Real}}],
Infinity
],
Frame -> True, PlotRange -> {{0, 1/10}, {-1/10, 1/10}},
PlotRangePadding -> Scaled[.02], AspectRatio -> 1]
` 