# What is my problem with MeshFunction?

I want plot the max value in a sine plot. We can use the following code

Plot[Sin[x], {x, 0, 50}, Mesh -> {{.99}},
MeshFunctions -> {#2 &},
MeshStyle -> {PointSize[Large], Red}]


Why can't Mesh -> 1 be used to plot the red point?

But in a mathematics way, I want to use this code

Plot[Sin[x], {x, 0, 50},
Mesh -> {{1}},
MeshFunctions -> {Boole[(Cos[#] == 0) && (-Sin[#] < 0)] &},
MeshStyle -> {PointSize[Large], Red}]


You can see the MeshFunction doesn't work? Why?

Then I try some other test. For example,I want to emphasize the point above 0.5 to draw on Red. I use this code

Plot[Sin[x], {x, 0, 50},
Mesh -> {{1}},
MeshFunctions -> {Boole[Greater[#2, 0.5]] &}]


It doesn't work again. So I guess the equation PrimePi[z] == 2 will give 2 <= z < 3. To see whether this region can plot in MeshFunction, I tried the following:

f[x_, y_] := (x^2 + 3 y^2)*E^(1 - x^2 - y^2)

Plot3D[f[x, y], {x, -2, 2}, {y, -2, 2},
Mesh -> {{1}},
MeshFunctions -> {PrimePi[#3] &}]


You can see the plot is weird.

So what is my problem with MeshFunction? Can I use MeshFunction to plot the maximum value in a plot?

• The first problem is related to a similar problem in ContourPlot: 23363, 32734. – Michael E2 Aug 25 '14 at 12:21
• Would Plot[Sin[x], {x, 0, 50}, Mesh -> {Pi/2 + 2 Pi Range[0, 50/(2 Pi)]}, MeshStyle -> {PointSize[Large], Red}] be sufficient? Or do you seek a particular method of solution? – Michael E2 Aug 25 '14 at 13:34

Another alternative is to use ConditionalExpression using the second-order condition for a local maximum as the second argument:

f = Sin;
Plot[f[x], {x, 0, 20 Pi},
Mesh -> {{0}},
MeshFunctions -> {ConditionalExpression[f'[#], f''[#] < 0] &},
MeshStyle -> {PointSize[Large], Red}]


f = Sin[#] - 1/2 Cos[Pi #] &;
Plot[f[x], {x, 0, 10 Pi},
Mesh -> {{0}},
MeshFunctions -> {ConditionalExpression[f'[#], f''[#] < 0] &},
MeshStyle -> {PointSize[Large], Red}]


You can add additional constraints to the second argument of ConditionalExpression

e.g. f[#]>0:

 Plot[f[x], {x, 0, 20 Pi},
Mesh -> {{0}},
MeshFunctions -> {ConditionalExpression[f'[#],f''[#] < 0 && f[#] > 0] &},
MeshStyle -> {PointSize[Large], Red}]


Update ... or, essentially any constraint:

cond = ((# - 10 Pi)/(2 Pi))^2 + (Pi #2)^2 &;
rplt = RegionPlot[{#, ! #}, {x, 0, 20  Pi}, {y, -2, 2},
PlotLegends -> "Expressions"] &@(4 < cond[x, y] < 16);
meshF = ConditionalExpression[f'[#], f''[#] < 0 && ff@(4 < cond[#, f[#]] < 16)] &;
Row[Table[Plot[f[x], {x, 0, 20 Pi},
Mesh -> {{0}}, ImageSize -> 400,
MeshFunctions -> {meshF},
MeshStyle -> {PointSize[Large], Red},
Epilog -> {Opacity[.6], rplt[[1, 1]]}],
{ff, {Identity, Not}}],
rplt[[2, 1, 1]]]


• wow..smart! Well done! – partida Aug 26 '14 at 13:55

MeshFunctions, according to the documentation, "should normally be chosen to be continuous monotonic functions." Failing that, the mesh functions should be transverse to the mesh levels (i.e., cross them, not have a local extremum); in this case, however, one might have trouble with sampling missing a small region where the mesh function very briefly crosses and recrosses the mesh value.

The way the mesh algorithm works is to find two points where the function values bracket the mesh value, and then refine by subdividing the interval. If the mesh function does not cross the mesh value, it is (probably) impossible for the mesh point to be found.

In short, there's no way to do what the OP wants to in Plot in the way the OP wants to do it. But there are other ways. In the Plot3D, it is worse in that PrimePi is locally constant. So the equation PrimePi[f[x, y]] == 2 has a region of solutions and cannot be used to define a mesh line.

One solution

Solve for the maxima and use the default # & as the mesh function. This mesh function is increasing. One can either list the x coordinates in a simple case like Sin[x] or use Solve or NSolve to solve for them.

Plot[Sin[x], {x, 0, 50}, Mesh -> {Pi/2 + 2 Pi Range[0, 50/(2 Pi)]},
MeshStyle -> {PointSize[Large], Red}]

Plot[Sin[x], {x, 0, 50},
Mesh -> {x /. NSolve[D[Sin[x], x] == 0 && Sin[x] > 0 && 0 < x < 50, x]},
MeshStyle -> {PointSize[Large], Red}]


Other possibilities

Other solutions tend to be more complicated. One can try to create a mesh function that is continuous, monotonic and takes on, say, integer values exactly at maxima. For instance,

MeshFunctions -> (#/(2 Pi) - 1/4 + (#2 - 1)/(2 Pi) &),
Mesh -> {Range[0, 50/(2 Pi)]}


An alternative in this case is to apply the double angle formula and factor out of the derivative a factor that vanishes exactly at maxima.

doubleangle = Cos[t_] :> Cos[t/2]^2 - Sin[t/2]^2;
Cos[x] /. doubleangle // Factor
(* (Cos[x/2] - Sin[x/2]) (Cos[x/2] + Sin[x/2]) *)

Plot[Sin[x], {x, 0, 50},
Mesh -> {{0}},
MeshFunctions -> (Cos[#/2] - Sin[#/2] &),
MeshStyle -> {PointSize[Large], Red}]


But they're certainly more work and rely on the simplicity of Sin[x] to greater extent than the first solutions.

• I don't understand the meaning of the OP.can you explain for it? – partida Aug 27 '14 at 1:42
• @user15961 "OP" means the Original Poster -- that is, you. What I was referring to was trying to use the mesh functions depending on #2 &, Boole, or Sin[#1] & to detect extrema, which is what I thought you were trying to do. It seems reasonable at first, but as you noticed, it wasn't working, – Michael E2 Aug 27 '14 at 1:49

You have a numerical accuracy problem more than anything else. The mesh points are not computed all that accurately, and you have to allow for it in your code.

For example,

Plot[Sin[x], {x, 0, N[8 Pi]},
Mesh -> {{1}},
MeshFunctions -> (Boole[Chop[Cos[#1], .005] == 0 && #2 > 0] &),
MeshStyle -> {PointSize[Large], Red}]


works fine and gives

### Edit

This method is sensitive to the size of plot domain. To get good results form the OP's full domain of {x, 0, 50}, one must loosen the chop tolerance to about 0.01.

Plot[Sin[x], {x, 0, 50},
Mesh -> {{1}},
MeshFunctions -> (Boole[Chop[Cos[#1], .01] == 0 && #2 > 0] &),
MeshStyle -> {PointSize[Large], Red}]


• Note that it gives some extra points, but it looks fine: Cases[Normal@Plot[Sin[x], {x, 0, N[8 Pi]}, Mesh -> {{1}}, MeshFunctions -> (Boole[Chop[Cos[#1], .005] == 0 && #2 > 0] &), MeshStyle -> {PointSize[Large], Red}], _Point, Infinity]. – Michael E2 Aug 25 '14 at 14:45
• @MichaelE2. Do you think the multiplicity invalidates my answer? – m_goldberg Aug 25 '14 at 15:19
• No, I think it is valid because with plotting the goal should be for the graphics to look right, which is so. It is, in some sense, less efficient than ideal. I'm working a more expansive explanation of what's happening with the mesh functions. – Michael E2 Aug 25 '14 at 15:27
• It occurs to me that the success of this method depends on number of PlotPoints being sufficiently high. For the OP's domain {x, 0, 50}, I need around PlotPoints -> 80. – Michael E2 Aug 25 '14 at 16:35
• @MichaelE2. It is easier and more efficient to loosen the chop tolerance. Boole[Chop[Cos[#1], .01] == 0 && #2 > 0] & works fine of the domain {x, 0, N[15 Pi]} without increasing PlotPoints. – m_goldberg Aug 25 '14 at 16:40

I think Michael E2 is exactly true. "should normally be chosen to be continuous monotonic functions." I have tried draw like this. Because it is not continuous monotonic function it is drawed lines along with the polygon meshs. I think MeshFunction seems improbable for your purpose.

f[x_, y_] := (x^2 + 3 y^2)*E^(1 - x^2 - y^2)

f1 = Plot3D[f[x, y], {x, -2, 2}, {y, -2, 2}, Mesh -> {{1}},
MeshFunctions -> {PrimePi[#3] &},
MeshStyle -> {Thickness[0.006], Red}];
f2 = Graphics3D[Replace[f1[[1]],
EdgeForm[___] :> EdgeForm[{Opacity[0.3], Blue}], 4],
BoxRatios -> {1, 1, 0.4}];
Grid[{{f1, f2}}]


g1 = Plot3D[f[x, y], {x, -2, 2}, {y, -2, 2}, Mesh -> {{2}},
MeshFunctions -> {#3 &},
MeshStyle -> {Thickness[0.006], Red}];
g2 = Graphics3D[Replace[g2[[1]],
EdgeForm[___] :> EdgeForm[{Opacity[0.3], Blue}], 4],
BoxRatios -> {1, 1, 0.4}];
Grid[{{g1, g2}}]


This is a verbose epilog to the nice answers given:

fun = Sin[x];
lim = 4 Pi;

max = Round @ First @ FindMaximum[fun, {x, 0}];
min = Round @ First @ FindMinimum[fun, {x, 0}];

xp =
FindInstance[(fun == max || fun == min) && 0 <= x < lim, x, Reals, 15]//Values//Flatten;


yp = Table[fun, {x, xp}]


plo =
Plot[fun, {x, 0, lim},
GridLines -> {xp, {min, 0, max}},
GridLinesStyle -> Dashed,
Frame -> True,
FrameTicks -> {{{min, 0, max}, None}, {xp, None}},
Epilog -> {PointSize @ 0.02, Point @ Transpose[{xp, yp}]}];

nlp =
With[{d = D[fun, x]},
NumberLinePlot[{d > 0, d < 0, d == 0}, {x, 0, lim},
PlotStyle -> {Green, Red, Black},
Spacings -> 0,
PlotLegends -> {"Increasing", "Decreasing", "Stationary"}]];

Show[plo, nlp, ImageSize -> 500]


• Perhaps using Round on the maximum and minimum only makes sense if you know beforehand that that they are integers? (If you already know them, why not use their values?) Also, the OP wants only the maxima marked. – Michael E2 Aug 26 '14 at 21:36
• @MichaelE2 Even worse: I haven't answered the question which is about MeshFunctions. Once in a while I misuse SE as a personal cloud to see "my natural solution" in the context of other solutions. If this is a misuse please let me know, and I'll delete. – eldo Aug 26 '14 at 22:35
• I suppose you might point out in the text of your answer that another way to get the desired points is through Epilog instead of MeshFunctions. One sometimes finds such alternative approaches on the site. – Michael E2 Aug 26 '14 at 23:14