Skip to content

Commit

Permalink
mesh: allow for marking facets with Q2 for hex for now
Browse files Browse the repository at this point in the history
  • Loading branch information
ksagiyam committed Oct 18, 2023
1 parent 4b202ef commit 4b5593f
Show file tree
Hide file tree
Showing 2 changed files with 24 additions and 2 deletions.
12 changes: 10 additions & 2 deletions firedrake/mesh.py
Original file line number Diff line number Diff line change
Expand Up @@ -835,6 +835,8 @@ def mark_entities(self, tf, label_value, label_name=None):
must be "DP" or "DQ" (degree 0) to mark cell entities and
"P" (degree 1) in 1D or "HDiv Trace" (degree 0) in 2D or 3D
to mark facet entities.
Can use "Q" (degree 2) functions for 3D hex meshes until
we support "HDiv Trace" elements on hex.
:arg lable_value: The value used in the label.
:arg label_name: The name of the label to store entity selections.
Expand Down Expand Up @@ -1278,7 +1280,8 @@ def mark_entities(self, tf, label_value, label_name=None):
height = 0
label_name = label_name or dmcommon.CELL_SETS_LABEL
elif (elem.family() == "HDiv Trace" and elem.degree() == 0 and self.cell_dimension() > 1) or \
(elem.family() == "Lagrange" and elem.degree() == 1 and self.cell_dimension() == 1):
(elem.family() == "Lagrange" and elem.degree() == 1 and self.cell_dimension() == 1) or \
(elem.family() == "Q" and elem.degree() == 2 and self.ufl_cell().cellname() == "hexahedron"):
# facets
height = 1
label_name = label_name or dmcommon.FACE_SETS_LABEL
Expand Down Expand Up @@ -2363,6 +2366,8 @@ def mark_entities(self, f, label_value, label_name=None):
must be "DP" or "DQ" (degree 0) to mark cell entities and
"P" (degree 1) in 1D or "HDiv Trace" (degree 0) in 2D or 3D
to mark facet entities.
Can use "Q" (degree 2) functions for 3D hex meshes until
we support "HDiv Trace" elements on hex.
:arg lable_value: The value used in the label.
:arg label_name: The name of the label to store entity selections.
Expand Down Expand Up @@ -3980,6 +3985,8 @@ def RelabeledMesh(mesh, indicator_functions, subdomain_ids, **kwargs):
"DP"/"DQ" (degree 0) functions to mark cell entities and
"P" (degree 1) functions in 1D or "HDiv Trace" (degree 0) functions
in 2D or 3D to mark facet entities.
Can use "Q" (degree 2) functions for 3D hex meshes until
we support "HDiv Trace" elements on hex.
:arg subdomain_ids: list of subdomain ids associated with
the indicator functions in indicator_functions; thus,
must have the same length as indicator_functions.
Expand Down Expand Up @@ -4032,7 +4039,8 @@ def RelabeledMesh(mesh, indicator_functions, subdomain_ids, **kwargs):
height = 0
dmlabel_name = dmcommon.CELL_SETS_LABEL
elif (elem.family() == "HDiv Trace" and elem.degree() == 0 and mesh.topological_dimension() > 1) or \
(elem.family() == "Lagrange" and elem.degree() == 1 and mesh.topological_dimension() == 1):
(elem.family() == "Lagrange" and elem.degree() == 1 and mesh.topological_dimension() == 1) or \
(elem.family() == "Q" and elem.degree() == 2 and mesh.topology.ufl_cell().cellname() == "hexahedron"):
# facets
height = 1
dmlabel_name = dmcommon.FACE_SETS_LABEL
Expand Down
14 changes: 14 additions & 0 deletions tests/regression/test_mark_entities.py
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,20 @@ def test_mark_entities_overlapping_facet_subdomains():
assert abs(v - 1.0) < 1.e-10


def test_mark_entities_mesh_mark_entities_3d_hex():
label_name = "test_label"
label_value = 999
mesh = UnitCubeMesh(2, 2, 2, hexahedral=True)
x, y, z = SpatialCoordinate(mesh)
V = FunctionSpace(mesh, "Q", 2)
f = Function(V).interpolate(conditional(And(x > .6, And(y > .4, y < .6)), 1, 0))
mesh.mark_entities(f, label_value, label_name=label_name)
plex = mesh.topology.topology_dm
label = plex.getLabel(label_name)
assert label.getStratumIS(label_value).getSize() == 2
assert all(label.getStratumIS(label_value).getIndices() == [51, 57])


def test_mark_entities_mesh_mark_entities_2d():
# +---+---+
# | \ | \ | mark \
Expand Down

0 comments on commit 4b5593f

Please sign in to comment.