# 3D FEM - How to define Point-Markers for an ElementMesh object created with ElementMeshRegionProduct

I want to create a 3D grid based on a 2D mesh, on the surface of which I want to assign point-marker for certain regions. The mesh created in the sample code is a fluidic Y-structure, at whose inputs (marker 4,3) and outputs (marker 2) I want to define corresponding Dirichlet boundary conditions.

As soon as I create the 3D mesh from the 2D mesh, all markers assigned in the 2D mesh are set to 0. After creating the 3D mesh, how can I add markers for the point markers at the inputs and outputs afterwards?

Thank you very much.

Needs["NDSolveFEM"]

coord = {{0., -25.}, {200., -25.}, {200., 25.}, {0., 25.}, {-75., 100.}, {-75., 75.}, {0., 0.}, {-75., -75.}, {-75., -100.}};
iID = {{1, 2}, {2, 3}, {3, 4}, {4, 5}, {5, 6}, {6, 7}, {7, 8}, {8, 9}, {9, 1}};
marker = {1, 2, 1, 1, 3, 1, 1, 4, 1};

(*2D-Boundary-Mesh*)
bmesh2D = ToBoundaryMesh["Coordinates" -> coord, "BoundaryElements" -> {LineElement[iID, marker]}]
bmesh2D["Wireframe"["MeshElement" -> "BoundaryElements", "MeshElementMarkerStyle" -> Blue]]

mesh2D = ToElementMesh[bmesh2D];
mesh2D["PointElementMarkerUnion"]
Show[mesh2D["Wireframe"], mesh2D["Wireframe"["MeshElement" -> "PointElements", "MeshElementMarkerStyle" -> Blue]]]

(*create 3D-Mesh*)
mesh = ElementMeshRegionProduct[mesh2D, ToElementMesh[Line[{{0}, {50}}], "MaxCellMeasure" -> 5.0]]
mesh["PointElementMarkerUnion"]
mesh["Wireframe"]

bcs = {DirichletCondition[{u[x, y, z] == 0, v[x, y, z] == 0.,
w[x, y, z] == 0.}, ElementMarker == 1],
DirichletCondition[p[x, y, z] == pIn1, ElementMarker == 4],
DirichletCondition[p[x, y, z] == pIn2, ElementMarker == 3],
DirichletCondition[p[x, y, z] == pOut, ElementMarker == 2]}

• I am not sure this is possible to do automatically; what outcome would you have expected? Nov 2, 2022 at 7:50

Thanks to the hint from user21, I found a solution to my question using ElementMeshResetPointElementMarker.

Needs["NDSolveFEM"]

(*2D-Mesh*)
coord = {{0., -25.}, {200., -25.}, {200., 25.}, {0., 25.}, {-75., 100.}, {-75., 75.}, {0., 0.}, {-75., -75.}, {-75., -100.}};
iID = {{1, 2}, {2, 3}, {3, 4}, {4, 5}, {5, 6}, {6, 7}, {7, 8}, {8, 9}, {9, 1}};
marker = {1, 2, 1, 1, 3, 1, 1, 4, 1};

(*2D-Boundary-Mesh*)
bmesh2D = ToBoundaryMesh["Coordinates" -> coord, "BoundaryElements" -> {LineElement[iID, marker]}]
mesh2D = ToElementMesh[bmesh2D];
Show[mesh2D["Wireframe"], mesh2D["Wireframe"["MeshElement" -> "PointElements", "MeshElementMarkerStyle" -> Blue]]]

(*3D-Mesh*)
mesh = ElementMeshRegionProduct[mesh2D, ToElementMesh[Line[{{0},{50}}], "MaxCellMeasure" -> 5.0]]
mesh["Wireframe"]

(*PointMarkerFunction*)
xMin = -75.0; xMax = 200.;
pMF = Compile[{{coords, _Real, 2}},
Block[{x = #[[1]], y = #[[2]], z = #[[3]], epsilon}, epsilon = 0.01;
Which[((xMax - x) <= epsilon) && (0. < z < 50.), 2,
((x - xMin) <= epsilon) && ((75. < y < 100) && (0. < z < 50.)), 3,
((x - xMin) <= epsilon) && ((-100. < y < -75.) && (0. < z < 50.)), 4,
True, 1]] & /@ coords];

peInci = Flatten[ElementIncidents[mesh["PointElements"]]];
peCoord = mesh["Coordinates"][[#]] & /@ peInci;
peMarker = pMF[peCoord];

(*Mesh-Reset*)
m = ElementMeshResetPointElementMarker[mesh, peCoord, peMarker]
m["Wireframe"]
m["PointElementMarkerUnion"]

(*Plot Point-Marker*)
pIDs = m["PointElementMarkerUnion"];
Manipulate[
Show[m["Wireframe"[
"MeshElementStyle" ->
Directive[Opacity[0.2], FaceForm[LightBlue], EdgeForm[]]]],
m["Edgeframe"],
m["Wireframe"[ElementMarker == pIDs[[id]],
"MeshElement" -> "PointElements",
"MeshElementStyle" -> Directive[Red]]]], {{id, 1,
"ElementMarker ID"}, 1, Length[pIDs], 1, Appearance -> "Open"},
SaveDefinitions -> True]


Unfortunately, that is not possible directly right now. However, you can do it manually. For that you would regenerate the mesh with ToElementMesh and specify a BoundaryMarkerFunction and a PointMarkerFunction. So, the None needs to be replaced with a function that does what you want. (See the documentation of ToElementMesh in the options section to see how those are written)

m2 = ToElementMesh["Coordinates" -> mesh["Coordinates"],
"MeshElements" -> mesh["MeshElements"],
"BoundaryMarkerFunction" -> None, "PointMarkerFunction" -> None]


Update: The direct setting of boundary and point markers will possible in version 13.3 (not in 13.2).

• I created a PointMarker and BoundaryMarker function as I believe is described in the documentation you mentioned and then recreated the mesh again - unfortunately without success: ´xMin = -75.0; xMax = 200.; boundaryMarkerFunction = Compile[{{boundaryElementCoords, _Real, 3}, {boundaryElementPointMarkres, _Integer, 2}}, Which[Union[#] == {2}, 2, Union[#] == {3}, 3, Union[#] == {4}, 4, True, 1 ] & /@ boundaryElementPointMarkres];´ Nov 3, 2022 at 12:29
• pointMarkerFunction = Compile[{{coords, _Real, 2}, {pMarker, _Integer, 1}}, Block[{x = #[[1]], y = #[[2]], z = #[[3]], epsilon}, epsilon = 0.1; Which[(xMax - x) <= epsilon, 2, ((x - xMin) <= epsilon) &&(75. <= y <= 100), 3,((x - xMin) <= epsilon) && (-100. <= y <= -75.), 4, True, 1]] & /@ coords]; m2 = ToElementMesh["Coordinates" -> mesh["Coordinates"], "MeshElements" -> mesh["MeshElements"], "PointMarkerFunction" -> pointMarkerFunction, "BoundaryMarkerFunction" -> boundaryMarkerFunction] Nov 3, 2022 at 12:36
• @eragim, I think you need to describe the condition with the coordinates; since there are no boundary markers or point markers. Nov 3, 2022 at 15:23