# Finer mesh for selected subregion of a solid 3D cylinder

I am trying to create a finer mesh in a subregion of a solid cylinder in order to improve the resolution there. I don't want to have too many elements for the entire cylinder because it slows down the computation. From other questions here and the tutorial on Element Mesh Creation https://reference.wolfram.com/language/FEMDocumentation/tutorial/ElementMeshCreation.html, I see that this can be done for 2D shapes and composite shapes (different materials).

1. I'd like a finer mesh along the axis of the cylinder, with a smaller radius (e.g. 1/3 of the radius).

I suppose I can break it in two separate cylinders--the inner solid cylinder and the outer "pipe" and then combine them with RegionUnion. But, then I am not sure how to define the boundary surfaces for the boundary conditions (outer cylindrical surface and the whole end faces). I hope there is a more elegant solution.

1. How to create finer mesh just for the end faces? The boundary surfaces are defined with ElementMarker. I tried defining separate boundary mesh for each face and then turning it ToElementMesh, but it didn't work. Below is the code for the whole solid cylinder.
OpenCascadeShapeType[cyl];
cyl = OpenCascadeShape[c1 = Cylinder[{{0, 0, 0}, {0.01, 0, 0}}, 0.003]];
mesh = ToElementMesh[bmesh]
Show[mesh["Wireframe"], Graphics3D[c1]]


When I visualize the solution with SliceContourPlot3D, the resolution is insufficient to capture small features near the axis.

• You'll probably need to find an appropriate MeshRefinementFunction. Perhaps start by searching the site with that keyword. Jul 17, 2023 at 0:17

Version 1

One way to do it could be to generate a 2D disk mesh with an appropriate refinement and the extrude that to 3D.

(m2 = ToElementMesh[Disk[{0, 0}, 0.003],
MeshRefinementFunction ->
Function[{vertices, area},
area > 0.00000005 (0.001 +
500 Norm[Mean[vertices]])]])["Wireframe"]


mesh = ElementMeshRegionProduct[m2,
ToElementMesh[Line[{{0}, {0.01}}]]]

mesh["Wireframe"]


Here is a view inside:

mesh["Wireframe"["MeshElement" -> "MeshElements",
"ElementMeshDirective" ->
Directive[EdgeForm[Black], FaceForm[LightGray]],
PlotRange -> {{-0.003, 0}, All, {0, 0.005}}]]


Version 2

A variant of this approach is to use a graded mesh for the element mesh region product, like so:

m1 = ToGradedMesh[{Line[{{0}, {0.01}}], <|
"Alignment" -> "BothEnds"|>}];
MeshRegion[m1]


mesh = ElementMeshRegionProduct[m2, m1];

mesh["Wireframe"["MeshElement" -> "MeshElements",
"ElementMeshDirective" ->
Directive[EdgeForm[Black], FaceForm[LightGray]],
PlotRange -> {{-0.003, 0}, All, {0, 0.005}}]]


Note, how the mesh is becoming denser towards the bottom (and top not shown) Here is a view of the boundaries:

mesh["Wireframe"[
"MeshElementStyle" -> {FaceForm[Green], FaceForm[Red],
FaceForm[Blue]}, PlotRange -> {{-0.003, 0}, All, {0, 0.01}}]]


Version 3

Here is a version that has a refined boundary but a coarse interior:

Needs["NDSolveFEM"]
c1 = Cylinder[{{0, 0, 0}, {0.01, 0, 0}}, 0.003];
mesh = ToElementMesh[c1, "MaxBoundaryCellMeasure" -> 0.1,
"MaxCellMeasure" -> Infinity];
Show[mesh["Wireframe"], Graphics3D[c1]]


Cross section:

mesh["Wireframe"["MeshElement" -> "MeshElements",
"ElementMeshDirective" ->
Directive[EdgeForm[Black], FaceForm[LightGray]],
PlotRange -> {{0, 0.005}, {0, 0.003}, {0, 0.003}}]]


Version 4

Yet another way to do it would be by using OpenCascade and create an inner cylinder that is refined to a different level.

Needs["NDSolveFEM"]
Needs["OpenCascadeLink"]
c1 = Cylinder[{{0, 0, 0}, {0.01, 0, 0}}, 0.003];
c2 = Cylinder[{{0, 0, 0}, {0.01, 0, 0}}, 0.003/2];

faces = Flatten[OpenCascadeShapeFaces /@ {s1, s2}];

bmesh["Wireframe"]


Create to material markers p1 and p2:

p1 = {0.01/2, 0, 0};
p2 = {0.01/2, 0, 0.005/2};
Show[
Graphics3D[{{Red, PointSize[0.02], Point[p1]}, {Green,
PointSize[0.02], Point[p2]}}],
bmesh["Wireframe"],
Graphics3D[{Opacity[0.2], c1}]]


mesh = ToElementMesh[bmesh,
"RegionMarker" -> {{p1, 1, 10^-12}, {p2, 2, 10^-10}}];
Show[mesh["Wireframe"], Graphics3D[c1]]


Version 5

Here is a version with a MeshRefinementFunction. By testing I found that the mesh changes when a "MaxCellMeasure" around 10^-11 was given. So I create a function that re-scales around this value from 0 to the boundary at 0.003.

Needs["NDSolveFEM"]
rescale = Rescale[Sqrt[y^2 + z^2], {0, 0.003}, {10^-14, 10^-11}]


The Sqrt[y^2 + z^2] is just the radial distance from the x-axis. This I then injected into the mesh refinement function:

c1 = Cylinder[{{0, 0, 0}, {0.01, 0, 0}}, 0.003];
cf = With[{r = rescale},
Compile[{{coordinates, _Real, 2}, {vol, _Real, 0}},
Block[{com, x, y, z, rvol},
com = Mean[coordinates];
{x, y, z} = com;
rvol = r;
If[vol > rvol, True, False]
]
]
];
mesh = ToElementMesh[c1, MeshRefinementFunction -> cf]
Show[mesh["Wireframe"], Graphics3D[c1]]


• how would you define the boundaries--cylindrical surface, top-end face, and bottom-end face? I tried mesh["Wireframe"[ "MeshElementStyle" -> {FaceForm[Green], FaceForm[Red], FaceForm[Blue]}, PlotRange -> {{-0.003, 0}, All, {0, 0.005}}]]` but it gives a series of cross-sections along the length (instead of top and bottom). In addition, I cannot assign names to the boundaries in order to use them in the boundary conditions for solving the PDE. Jul 17, 2023 at 22:08
• @user93372, I hope i have addressed the first part of your question, I do not understand what you mean with 'assign names' though. Jul 18, 2023 at 18:04