# RegionPlot with FindMaximum

I don't understand why each time I'm launching this computation Mathematica buggs and I get a message from Windows that the Kernel stopped working...

I'm trying to get a map of where the value of Eig is below 0 (or above).

Eig comes from a matrix, I'm writing it as an appendix if you're interested :

Eig = {(-x + x^2 - y - alpha x y - gamma x y + alpha x^2 y +
gamma x^2 y + y^2 + alpha x y^2 +
gamma x y^2 - \[Sqrt]((x - x^2 + y + alpha x y + gamma x y -
alpha x^2 y - gamma x^2 y - y^2 - alpha x y^2 -
gamma x y^2)^2 -
4 (x y - x^2 y + alpha x^2 y - 2 alpha x^3 y + alpha x^4 y -
x y^2 + gamma x y^2 - alpha x^2 y^2 - 2 beta x^2 y^2 -
beta^2 x^2 y^2 - gamma x^2 y^2 + alpha gamma x^2 y^2 +
alpha x^3 y^2 + 2 beta x^3 y^2 + 2 beta^2 x^3 y^2 -
2 alpha gamma x^3 y^2 - beta^2 x^4 y^2 +
alpha gamma x^4 y^2 - 2 gamma x y^3 + 2 beta x^2 y^3 +
2 beta^2 x^2 y^3 + gamma x^2 y^3 - 2 alpha gamma x^2 y^3 -
2 beta^2 x^3 y^3 + 2 alpha gamma x^3 y^3 + gamma x y^4 -
beta^2 x^2 y^4 + alpha gamma x^2 y^4))), (-x + x^2 - y -
alpha x y - gamma x y + alpha x^2 y + gamma x^2 y + y^2 +
alpha x y^2 +
gamma x y^2 + \[Sqrt]((x - x^2 + y + alpha x y + gamma x y -
alpha x^2 y - gamma x^2 y - y^2 - alpha x y^2 -
gamma x y^2)^2 -
4 (x y - x^2 y + alpha x^2 y - 2 alpha x^3 y + alpha x^4 y -
x y^2 + gamma x y^2 - alpha x^2 y^2 - 2 beta x^2 y^2 -
beta^2 x^2 y^2 - gamma x^2 y^2 + alpha gamma x^2 y^2 +
alpha x^3 y^2 + 2 beta x^3 y^2 + 2 beta^2 x^3 y^2 -
2 alpha gamma x^3 y^2 - beta^2 x^4 y^2 +
alpha gamma x^4 y^2 - 2 gamma x y^3 + 2 beta x^2 y^3 +
2 beta^2 x^2 y^3 + gamma x^2 y^3 - 2 alpha gamma x^2 y^3 -
2 beta^2 x^3 y^3 + 2 alpha gamma x^3 y^3 + gamma x y^4 -
beta^2 x^2 y^4 + alpha gamma x^2 y^4)))}
Eig[[2]]
gamma = 0
f2[b0_?NumericQ, b1_?NumericQ] :=
Block[{alpha = b0, beta = b1},
FindMaximum[{Re[Eig[[2]]],
0.001 < x < 0.99 && 0.01 < y < 0.99 && x + y < 0.99}, {x,
y}][[1]]];

RegionPlot[f2[b0, b1] <= 0, {b0, 0, 3}, {b1, 0, 3}]



The corresponding densityplot is working :

DensityPlot[f2[b0, b1], {b0, 0, 3}, {b1, 0, 3},
PlotLegends -> Automatic, PlotRange -> All]


But as one can see, it's not so smooth due to the FindMaximum that is not super efficient, so I assume RegionPlot has some difficulties to find smooth borders for the domains.

Do you have an idea how to solve that please ?

And for example would it be possible to smooth the FindMaximum function to see if it would give better results if it was smooth ?

Appendix : the matrix :

Eig = Eigenvalues[{{1/(1 - x - y) + alpha + 1/x,
1/(1 - x - y) + beta}, {1/(1 - x - y) + beta,
1/y + 1/(1 - x - y) +gamma}}]

--->
{(1/(2 x y (-1 + x + y)))(-x + x^2 - y - alpha x y - gamma x y +
alpha x^2 y + gamma x^2 y + y^2 + alpha x y^2 +
gamma x y^2 - \[Sqrt]((x - x^2 + y + alpha x y + gamma x y -
alpha x^2 y - gamma x^2 y - y^2 - alpha x y^2 -
gamma x y^2)^2 -
4 (x y - x^2 y + alpha x^2 y - 2 alpha x^3 y + alpha x^4 y -
x y^2 + gamma x y^2 - alpha x^2 y^2 - 2 beta x^2 y^2 -
beta^2 x^2 y^2 - gamma x^2 y^2 + alpha gamma x^2 y^2 +
alpha x^3 y^2 + 2 beta x^3 y^2 + 2 beta^2 x^3 y^2 -
2 alpha gamma x^3 y^2 - beta^2 x^4 y^2 +
alpha gamma x^4 y^2 - 2 gamma x y^3 + 2 beta x^2 y^3 +
2 beta^2 x^2 y^3 + gamma x^2 y^3 - 2 alpha gamma x^2 y^3 -
2 beta^2 x^3 y^3 + 2 alpha gamma x^3 y^3 + gamma x y^4 -
beta^2 x^2 y^4 + alpha gamma x^2 y^4))), (1/(
2 x y (-1 + x + y)))(-x + x^2 - y - alpha x y - gamma x y +
alpha x^2 y + gamma x^2 y + y^2 + alpha x y^2 +
gamma x y^2 + \[Sqrt]((x - x^2 + y + alpha x y + gamma x y -
alpha x^2 y - gamma x^2 y - y^2 - alpha x y^2 -
gamma x y^2)^2 -
4 (x y - x^2 y + alpha x^2 y - 2 alpha x^3 y + alpha x^4 y -
x y^2 + gamma x y^2 - alpha x^2 y^2 - 2 beta x^2 y^2 -
beta^2 x^2 y^2 - gamma x^2 y^2 + alpha gamma x^2 y^2 +
alpha x^3 y^2 + 2 beta x^3 y^2 + 2 beta^2 x^3 y^2 -
2 alpha gamma x^3 y^2 - beta^2 x^4 y^2 +
alpha gamma x^4 y^2 - 2 gamma x y^3 + 2 beta x^2 y^3 +
2 beta^2 x^2 y^3 + gamma x^2 y^3 - 2 alpha gamma x^2 y^3 -
2 beta^2 x^3 y^3 + 2 alpha gamma x^3 y^3 + gamma x y^4 -
beta^2 x^2 y^4 + alpha gamma x^2 y^4)))}


and I multiplicated by $$xy(1-x-y)$$ in order to simplify the problem.

EDIT 1 : It looks like the RegionFunction option doesn't work either for the DensityPlot :

DensityPlot[f2[b0, b1], {b0, 0, 3}, {b1, 0, 3},
PlotLegends -> Automatic, PlotRange -> All,
RegionFunction -> Function[{b0, b1}, f2[b0, b1] > 0.001]]


EDIT 2 : close to a solution but still weird : : So I found a way to have something that looks loke what I want. But It seems that at some point the result is not exact, since the diagram is a bit chaotic...

DensityPlot[
If[f2[b0, b1] > 0.0001, 1, 0], {b0, 0, 2.5}, {b1, -2.5, 2.5},
PlotLegends -> Automatic, PlotRange -> All]


• I assume Eig comes from an Eigenvalues call. If so, could you provide the matrix? – Carl Woll May 16 '19 at 15:11
• You're almost right :) – J.A May 16 '19 at 15:12