Etydhv

I am trying to plot a function in a region below a level curve of the function and within a cell. I have been doing this by calculating an ElementMesh using ImplicitRegion and ToElementMesh, but the result has holes.

Here is the cell (it's just a square),

cell = Parallelogram[-0.5`, -0.5`, 1.`, 0.`, 0.`, 1.`];
Graphics[Transparent, EdgeForm[Thick], cell]

and the function,

f[kx_, ky_, n_] := 
 Sort[Eigenvalues[(-1. + kx)^2 + (-1. + ky)^2, -0.23, 0., -0.23, 
 0.12, 0., 0., 0., 
 0., -0.23, (-1. + kx)^2 + (0. + ky)^2, -0.23, 0.12, -0.23, 
 0.12, 0., 0., 0., 0., -0.23, (-1. + kx)^2 + (1. + ky)^2, 0., 
 0.12, -0.23, 0., 0., 0., -0.23, 0.12, 
 0., (0. + kx)^2 + (-1. + ky)^2, -0.23, 0., -0.23, 0.12, 
 0., 0.12, -0.23, 
 0.12, -0.23, (0. + kx)^2 + (0. + ky)^2, -0.23, 0.12, -0.23, 
 0.12, 0., 0.12, -0.23, 0., -0.23, (0. + kx)^2 + (1. + ky)^2, 
 0., 0.12, -0.23, 0., 0., 0., -0.23, 0.12, 
 0., (1. + kx)^2 + (-1. + ky)^2, -0.23, 0., 0., 0., 0., 
 0.12, -0.23, 
 0.12, -0.23, (1. + kx)^2 + (0. + ky)^2, -0.23, 0., 0., 0., 
 0., 0.12, -0.23, 0., -0.23, (1. + kx)^2 + (1. + ky)^2]][[
 n]];
Plot3D[f[x, y, 4], x, y [Element] cell, PlotPoints -> 50]

and what the region should look like,

isovalue = 1.29897233417072;
ContourPlot[f[x, y, 4], x, y [Element] cell, 
 Contours -> isovalue, ColorFunction -> GrayLevel, 
 PlotPoints -> 100]

This is what I have tried

reg = ToElementMesh[
 ImplicitRegion[
 f[x, y, 4] < isovalue && x, y [Element] cell, x, y], 
 "MaxBoundaryCellMeasure" -> 0.01, MeshQualityGoal -> 1, 
 PerformanceGoal -> "Quality", MaxCellMeasure -> 0.01, 
 "BoundaryMeshGenerator" -> "Continuation"];
RegionPlot[reg]

The region is no more accurate when I decrease MaxCellMeasure or MaxBoundaryCellMeasure. I also tried the solution suggested here.

I hope I interpreted your question correctly that you want a more accurate ElementMesh representation of the region.

First we create a high quality Graphics of the region of interest.

isovalue = 1.29897233417072;
(* Add some margins to plot range to get connected region. *)
tolerance = 0.05;
plot = ContourPlot[
 f[x, y, 4],
 x, y ∈ Cuboid[-0.5, -0.5 - tolerance, 0.5, 0.5 + tolerance],
 Contours -> isovalue,
 ColorFunction -> GrayLevel,
 (* We need high quality plot for ImageMesh later. *)
 PlotPoints -> 200,
 Frame -> None
]

Create MeshRegion from Graphics object.

mreg = ImageMesh[ColorNegate[plot]]

And convert it to ElementMesh.

Needs["NDSolve`FEM`"]
mesh = ToElementMesh[mreg,"MeshOrder"->1]
(* ElementMesh[7., 353., 7., 353., TriangleElement["<" 1057 ">"]] *)

mesh["Wireframe"]

Another approach is:

reg = ToElementMesh[
 ImplicitRegion[
 f[x, y, 4] < isovalue && x, y [Element] cell, x, y], 
 "MaxBoundaryCellMeasure" -> 0.01, MeshQualityGoal -> 1, 
 PerformanceGoal -> "Quality", MaxCellMeasure -> 0.01, 
 "BoundaryMeshGenerator" -> "RegionPlot", "SamplePoints" -> 41];

reg["Wireframe"]

One thing to be a bit careful about is the question if the holes intersect the boundary. From the mesh it does not look like it but the math might say it.