FEM: Nicer Element Shape for Spherical Region

I'm trying to generate a mesh for later use in the Finite Element Method of the DSolve command. It is basically a parallelepiped with a spherical indentation. I'm trying to generate the mesh as follows:

Needs["NDSolveFEM"]
region = ImplicitRegion[!(Norm[{x, y, z}] < 1), {{x, -5, 5}, {y, -5, 5}, {z, 0, 5}}];
mesh = ToElementMesh[region];
Show[mesh["Wireframe"], PlotRange -> {{-2, 2}, {-2, 2}, {0, 1}}]


However, as you an see on the picture below this produces a really irregular and inaccurate (as I think) mesh for the spherical part. Is there a way to improve it and to generate more and homogeneously distributed elements on the spherical part?

The mesh seems to be fine and you can see that it is by doing:

region = ImplicitRegion[! (Norm[{x, y, z}] < 1), {{x, -5, 5}, {y, -5, 5}, {z, 0, 5}}];

m = DiscretizeRegion[region, {{-2, 2}, {-2, 2}, {0, 1}}]


To view as wireframe you can do:

Needs["NDSolveFEM"]

mesh = ToElementMesh[m] // Quiet;


Then:

Show[mesh["Wireframe"]]


If you want to keep your original region, you can vary MaxCellMeasure in DiscretizeRegion (of course, this will increase computation time and the number of mesh elements)

m2 = DiscretizeRegion[region, MaxCellMeasure -> 0.001]


Show[m2, PlotRange -> {{-2, 2}, {-2, 2}, {0, 1}}]


Edit

It seems there's a way to do what you want and hopefully keep performance hit to a minimum using MeshRefinementFunction:

m3 = DiscretizeRegion[region, MeshRefinementFunction -> Function[{vertices, vol},
Block[{x, y, z}, {x, y, z} = Mean[vertices];
If[-1 < x < 1 && -1 < y < 1 && 0 < z < 1, vol > 0.001, vol > 0.1]]]]


Then:

mesh3 = ToElementMesh[m3] // Quiet;

Show[mesh3["Wireframe"], PlotRange -> {{-2, 2}, {-2, 2}, {0, 1}}]


Looks reasonable.

• I'm not sure if it is a display problem or a discretization problem. You have chosen a smaller region, and the sphere part looks nice. But whenever my sphere is small in comparison to the big parallelepiped, e.g. sphere of radius 1 in a parallelepiped of dimensions (-5,5)x(-5,5)x(0,5), then it starts to look ugly as in my picture.
– Stan
Sep 2, 2014 at 13:07
• @Stan, please see my update. Sep 2, 2014 at 13:18
• yes this helps, but with the impact on the performance that you mentioned. I guess there is no easy way to generate a fine-grained mesh around the sphere and a coarse-grained one farther away? Something like to pre-generate a fine-grained mesh on the surface of the sphere separately, and then use it to generate the rest. Or to generate the parallelepiped and the sphere separately, and then to "subtract" the sphere from the block. Not sure if this makes sense.
– Stan
Sep 2, 2014 at 13:32
• @Stan, see my new edit. Hopefully this is what you want. Sep 2, 2014 at 13:47
• @Stan The Method -> "Continuation" option apparently helps with corners and cusps but it doesn't do anything in this case. Sep 2, 2014 at 13:54

Here are a few additions to @RunnyKine suggestions. If you are ever in doubt about the quality of a mesh (an ElementMesh to be exact) you can query the mesh.

Needs["NDSolveFEM"]
region = ImplicitRegion[! (Norm[{x, y, z}] < 1), {{x, -5, 5}, {y, -5,
5}, {z, 0, 5}}];
mesh = ToElementMesh[region];

Min[mesh["Quality"]]
0.004439742441262357


So the minimum mesh quality does not look too good. The overall mesh quality distribution seems sort of OK.

Histogram[mesh["Quality"]]


We can do better by increasing the sample points of the underlying RegionPlot (which is the only method available in 3D) with:

mesh = ToElementMesh[region,
"BoundaryMeshGenerator" -> {"RegionPlot", "SamplePoints" -> 25}];
Show[mesh["Wireframe"], PlotRange -> {{-2, 2}, {-2, 2}, {0, 1}}]


Looking at the quality again:

Histogram[mesh["Quality"]]
`

More information can be found in the documentation for ToElementMesh and in the mesh generation tutorial.

• I was hoping you would chime in. +1 Sep 2, 2014 at 14:13
• @RunnyKine, I was about not too :-) your answer is pretty good :-) Sep 2, 2014 at 14:18
• @user21: this is a very nice addition! Thanks!
– Stan
Sep 2, 2014 at 14:28
• I'm glad you did :) I like the analysis. Sep 2, 2014 at 14:32