Plotting boundary of piecewise function in ContourPlot

Question

I am trying to do a ContourPlot of the following function:

$$\psi=\begin{cases}\left(\sqrt{x^{2}+y^{2}}-\frac{0.25^{2}}{\sqrt{x^{2}+y^{2}}}\right)\frac{y}{\sqrt{x^{2}+y^{2}}} & x^{2}+y^{2} \geq 0.25^{2} \\ 0 & \text{otherwise}\end{cases}$$

However, when I plot this in Mathematica with the following code, I get a white line where the boundary should be (cf. image below):

ContourPlot[Piecewise[{{(Sqrt[x^2 + y^2] - 0.25^2/Sqrt[x^2 + y^2]) (y/Sqrt[x^2 + y^2]), x^2 + y^2 > 0.25^2}}, 0], {x, -1, 1}, {y, -1, 1}]

How can I fix this to get a contour boundary at $x^{2}+y^{2}=0.25^{2}$ instead of the white "undefined"-looking boundary I have at the moment?

Answer 1

ContourPlot[
 Piecewise[{{(Sqrt[x^2 + y^2] - 0.25^2/Sqrt[x^2 + y^2]) (y/
       Sqrt[x^2 + y^2]), x^2 + y^2 >= 0.25^2}}, 0], {x, -1, 
  1}, {y, -1, 1}, Exclusions -> None, PlotPoints -> 200]

Answer 2

ContourPlot works best on nonconstant, continuous functions, so one ought to expect to have to do some extra work when the function is discontinuous or, as in this case, constant on a region. In this case, the default behavior makes a choice that rather forces the user to make a choice about how to represent the function accurately. The particular issue is that the default choice of contours will include the contour f == 0. The problem is that the points where f == 0 does not form a curve but is the union of the horizontal line y == 0 and the disk about the origin. The disk in fact is colored the same as the region just above y == 0, which should indicate to that the function is actually positive on the disk, which is incorrect.

There are two workarounds that occur to me: (1) Avoid the contour f == 0 or (2) color the disk black.

For (1), specify the contours explicitly.

ContourPlot[
 Piecewise[{
  {(Sqrt[x^2 + y^2] - 0.25^2/Sqrt[x^2 + y^2]) (y/Sqrt[x^2 + y^2]), x^2 + y^2 > 0.25^2}},
  0],
 {x, -1, 1}, {y, -1, 1},
 Contours -> Range[-0.9, 0.9, 0.2], Exclusions -> None]

The bad part of this is that the locally constant behavior of the function is hidden. One can address this to some extent by increasing the number of contours or including two contours close together on either side of zero with something like

Contours -> Flatten[Range[-1, 1, 0.1] /. zero_ /; zero == 0 :> 10^-4 {-1, 1}]

Increasing PlotPoints or MaxRecursion will smooth the contours near the boundary of the disk. To my mind the contours near f == 0, being irregular, are a bit of a distraction.

For (2), add a black disk to the plot. Since the radius of the disk is 1/4 and the length of the domain in the x direction is 2, starting with a number of plot points that subdivides {-1, 1} into a multiple of 8 subintervals will produce a contour y == 0 that had an endpoint that matches the boundary of the disk.

ContourPlot[
 Piecewise[{
  {(Sqrt[x^2 + y^2] - 0.25^2/Sqrt[x^2 + y^2]) (y/Sqrt[x^2 + y^2]), x^2 + y^2 > 0.25^2}},
  0],
 {x, -1, 1}, {y, -1, 1},
 PlotPoints -> 2*8 + 1, Exclusions -> None,
 Epilog -> {Black, Tooltip[Disk[{0, 0}, 1/4], 0]}]