From 3d16443e10490ea43193b40941a5beb478ad34c3 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Fri, 8 May 2020 16:02:02 -0500
Subject: [PATCH 01/16] Small changes to get SBP-DG cases running with Grudge

---
 meshmode/discretization/__init__.py | 3 +--
 meshmode/mesh/generation.py         | 9 +++++----
 2 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/meshmode/discretization/__init__.py b/meshmode/discretization/__init__.py
index 6a4fcda7..70a93798 100644
--- a/meshmode/discretization/__init__.py
+++ b/meshmode/discretization/__init__.py
@@ -356,9 +356,8 @@ class Discretization(object):
                     "stride:auto,stride:auto,stride:auto")
             return knl
 
-        result = self.empty(dtype=self.real_dtype, extra_dims=(self.ambient_dim,))
-
         with cl.CommandQueue(self.cl_context) as queue:
+            result = self.empty(queue=queue, dtype=self.real_dtype, extra_dims=(self.ambient_dim,))
             for grp in self.groups:
                 if grp.nelements == 0:
                     continue
diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 3540dadd..de5d9581 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -658,7 +658,7 @@ def generate_urchin(order, m, n, est_rel_interp_tolerance, min_rad=0.2):
 # {{{ generate_box_mesh
 
 def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
-        group_factory=None):
+        group_factory=None, boundary_tags=[]):
     """Create a semi-structured mesh.
 
     :param axis_coords: a tuple with a number of entries corresponding
@@ -776,14 +776,15 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 
     from meshmode.mesh import Mesh
     return Mesh(vertices, [grp],
-            is_conforming=True)
+            is_conforming=True, boundary_tags=boundary_tags)
 
 # }}}
 
 
 # {{{ generate_regular_rect_mesh
 
-def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1):
+def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
+                               boundary_tags=[]):
     """Create a semi-structured rectangular mesh.
 
     :param a: the lower left hand point of the rectangle
@@ -797,7 +798,7 @@ def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1):
     axis_coords = [np.linspace(a_i, b_i, n_i)
             for a_i, b_i, n_i in zip(a, b, n)]
 
-    return generate_box_mesh(axis_coords, order=order)
+    return generate_box_mesh(axis_coords, order=order, boundary_tags=boundary_tags)
 
 # }}}
 
-- 
GitLab


From fc8cbe68b841c7c05673816adad421cd79170701 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Fri, 8 May 2020 16:38:44 -0500
Subject: [PATCH 02/16] Placate Flake8

---
 meshmode/discretization/__init__.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/meshmode/discretization/__init__.py b/meshmode/discretization/__init__.py
index 70a93798..9340a982 100644
--- a/meshmode/discretization/__init__.py
+++ b/meshmode/discretization/__init__.py
@@ -357,7 +357,8 @@ class Discretization(object):
             return knl
 
         with cl.CommandQueue(self.cl_context) as queue:
-            result = self.empty(queue=queue, dtype=self.real_dtype, extra_dims=(self.ambient_dim,))
+            result = self.empty(queue=queue, dtype=self.real_dtype,
+                                extra_dims=(self.ambient_dim,))
             for grp in self.groups:
                 if grp.nelements == 0:
                     continue
-- 
GitLab


From 5d563ef7974c92d1060d22e4d86d1df6488e9c03 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Tue, 12 May 2020 09:49:19 -0500
Subject: [PATCH 03/16] Reverts incorrect change to nodes result array

---
 meshmode/discretization/__init__.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/meshmode/discretization/__init__.py b/meshmode/discretization/__init__.py
index 9340a982..6a4fcda7 100644
--- a/meshmode/discretization/__init__.py
+++ b/meshmode/discretization/__init__.py
@@ -356,9 +356,9 @@ class Discretization(object):
                     "stride:auto,stride:auto,stride:auto")
             return knl
 
+        result = self.empty(dtype=self.real_dtype, extra_dims=(self.ambient_dim,))
+
         with cl.CommandQueue(self.cl_context) as queue:
-            result = self.empty(queue=queue, dtype=self.real_dtype,
-                                extra_dims=(self.ambient_dim,))
             for grp in self.groups:
                 if grp.nelements == 0:
                     continue
-- 
GitLab


From ce9852dd24e5ef61566af08aa48fe09b2acac1d0 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Tue, 12 May 2020 16:42:23 -0500
Subject: [PATCH 04/16] Beefs up boundary changes for box mesh generation

---
 meshmode/mesh/__init__.py   |  3 +-
 meshmode/mesh/generation.py | 56 +++++++++++++++++++++++++++++++++++--
 2 files changed, 55 insertions(+), 4 deletions(-)

diff --git a/meshmode/mesh/__init__.py b/meshmode/mesh/__init__.py
index 6d6191aa..16983a06 100644
--- a/meshmode/mesh/__init__.py
+++ b/meshmode/mesh/__init__.py
@@ -1141,7 +1141,8 @@ def _compute_facial_adjacency_from_vertices(groups, boundary_tags,
                     tag_mask = 0
                     for tag in tags:
                         tag_mask |= boundary_tag_bit(tag)
-                    fagrp.neighbors[idx] = -(-(fagrp.neighbors[idx]) | tag_mask)
+                    fagrp.neighbors[idx] = -(tag_mask | boundary_tag_bit(
+                                             BTAG_REALLY_ALL))
 
         else:
             raise RuntimeError("unexpected number of adjacent faces")
diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index de5d9581..04b5288f 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -658,7 +658,7 @@ def generate_urchin(order, m, n, est_rel_interp_tolerance, min_rad=0.2):
 # {{{ generate_box_mesh
 
 def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
-        group_factory=None, boundary_tags=[]):
+        group_factory=None, boundary_dict={}):
     """Create a semi-structured mesh.
 
     :param axis_coords: a tuple with a number of entries corresponding
@@ -666,6 +666,10 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         specifying the coordinates to be used along that axis.
     :param group_factory: One of :class:`meshmode.mesh.SimplexElementGroup`
         or :class:`meshmode.mesh.TensorProductElementGroup`.
+    :param boundary_dict: an optional dictionary for boundary configuration.
+        The keys correspond to custom boundary tags, with the values giving
+        a list of the faces on which they should be applied in terms of cardinal
+        directions.
 
     .. versionchanged:: 2017.1
 
@@ -774,8 +778,53 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
             vertices.reshape(dim, -1), el_vertices, order,
             group_factory=group_factory)
 
+    # compute facial adjacency for Mesh if there is tag information
+    facial_adjacency_groups = None
+    face_vertex_indices_to_tags = {}
+    boundary_tags = list(boundary_dict.keys())
+    if boundary_tags:
+        for tag_idx, tag in enumerate(boundary_tags):
+            # Need to map the correct face vertices to the boundary tags
+            for face in boundary_dict[tag]:
+                if face == "east":
+                    face_id = 1
+                    dim_crit = 0
+                    node_crit = axis_coords[0][0]
+                elif face == "west":
+                    face_id = 1
+                    dim_crit = 0
+                    node_crit = axis_coords[0][-1]
+                elif face == "north":
+                    face_id = 0
+                    dim_crit = 1
+                    node_crit = axis_coords[1][-1]
+                elif face == "south":
+                    face_id = 0
+                    dim_crit = 1
+                    node_crit = axis_coords[1][0]
+                for ielem in range(0, grp.nelements):
+                    if grp.nodes[dim_crit][ielem][0] == node_crit:
+                        fvi = grp.vertex_indices[ielem,
+                                                 grp.face_vertex_indices()[
+                                                     face_id]]
+                        if frozenset(fvi) not in face_vertex_indices_to_tags:
+                            face_vertex_indices_to_tags[frozenset(fvi)] = []
+                        face_vertex_indices_to_tags[frozenset(fvi)].append(tag)
+        # Need to add BTAG ALL and BTAG_REALLY_ALL to list of tags
+        from meshmode.mesh import BTAG_ALL, BTAG_REALLY_ALL
+        boundary_tags.append(BTAG_ALL)
+        boundary_tags.append(BTAG_REALLY_ALL)
+        from meshmode.mesh import _compute_facial_adjacency_from_vertices
+        facial_adjacency_groups = _compute_facial_adjacency_from_vertices(
+                [grp], boundary_tags, np.int32, np.int8,
+                face_vertex_indices_to_tags)
+    else:
+        facial_adjacency_groups = _compute_facial_adjacency_from_vertices(
+                [grp], boundary_tags, np.int32, np.int8)
+
     from meshmode.mesh import Mesh
     return Mesh(vertices, [grp],
+            facial_adjacency_groups=facial_adjacency_groups,
             is_conforming=True, boundary_tags=boundary_tags)
 
 # }}}
@@ -784,7 +833,7 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 # {{{ generate_regular_rect_mesh
 
 def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
-                               boundary_tags=[]):
+                               boundary_dict={}):
     """Create a semi-structured rectangular mesh.
 
     :param a: the lower left hand point of the rectangle
@@ -798,7 +847,8 @@ def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
     axis_coords = [np.linspace(a_i, b_i, n_i)
             for a_i, b_i, n_i in zip(a, b, n)]
 
-    return generate_box_mesh(axis_coords, order=order, boundary_tags=boundary_tags)
+    return generate_box_mesh(axis_coords, order=order,
+                             boundary_dict=boundary_dict)
 
 # }}}
 
-- 
GitLab


From fbca9d3033decf10e5c9ac17927335cc3e14b096 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Wed, 13 May 2020 11:47:27 -0500
Subject: [PATCH 05/16] Get Meshmode tests to pass by reinstating BTAG_ALL
 everywhere

---
 meshmode/mesh/__init__.py   | 3 +--
 meshmode/mesh/generation.py | 5 +----
 2 files changed, 2 insertions(+), 6 deletions(-)

diff --git a/meshmode/mesh/__init__.py b/meshmode/mesh/__init__.py
index 16983a06..6d6191aa 100644
--- a/meshmode/mesh/__init__.py
+++ b/meshmode/mesh/__init__.py
@@ -1141,8 +1141,7 @@ def _compute_facial_adjacency_from_vertices(groups, boundary_tags,
                     tag_mask = 0
                     for tag in tags:
                         tag_mask |= boundary_tag_bit(tag)
-                    fagrp.neighbors[idx] = -(tag_mask | boundary_tag_bit(
-                                             BTAG_REALLY_ALL))
+                    fagrp.neighbors[idx] = -(-(fagrp.neighbors[idx]) | tag_mask)
 
         else:
             raise RuntimeError("unexpected number of adjacent faces")
diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 04b5288f..2376cb44 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -782,6 +782,7 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     facial_adjacency_groups = None
     face_vertex_indices_to_tags = {}
     boundary_tags = list(boundary_dict.keys())
+    from meshmode.mesh import _compute_facial_adjacency_from_vertices
     if boundary_tags:
         for tag_idx, tag in enumerate(boundary_tags):
             # Need to map the correct face vertices to the boundary tags
@@ -814,13 +815,9 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         from meshmode.mesh import BTAG_ALL, BTAG_REALLY_ALL
         boundary_tags.append(BTAG_ALL)
         boundary_tags.append(BTAG_REALLY_ALL)
-        from meshmode.mesh import _compute_facial_adjacency_from_vertices
         facial_adjacency_groups = _compute_facial_adjacency_from_vertices(
                 [grp], boundary_tags, np.int32, np.int8,
                 face_vertex_indices_to_tags)
-    else:
-        facial_adjacency_groups = _compute_facial_adjacency_from_vertices(
-                [grp], boundary_tags, np.int32, np.int8)
 
     from meshmode.mesh import Mesh
     return Mesh(vertices, [grp],
-- 
GitLab


From 12c3cb67500c7799e1e8b9382bd97681cf6ed492 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Thu, 14 May 2020 10:05:46 -0500
Subject: [PATCH 06/16] Changes naming convention(s), adds check for
 dimensionality

---
 meshmode/mesh/generation.py | 52 ++++++++++++++++++++++++++++---------
 1 file changed, 40 insertions(+), 12 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 2376cb44..6fc9c8f5 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -658,7 +658,7 @@ def generate_urchin(order, m, n, est_rel_interp_tolerance, min_rad=0.2):
 # {{{ generate_box_mesh
 
 def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
-        group_factory=None, boundary_dict={}):
+        group_factory=None, face_to_boundary_tag={}):
     """Create a semi-structured mesh.
 
     :param axis_coords: a tuple with a number of entries corresponding
@@ -666,10 +666,10 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         specifying the coordinates to be used along that axis.
     :param group_factory: One of :class:`meshmode.mesh.SimplexElementGroup`
         or :class:`meshmode.mesh.TensorProductElementGroup`.
-    :param boundary_dict: an optional dictionary for boundary configuration.
+    :param face_to_boundary_tag: an optional dictionary for boundary configuration.
         The keys correspond to custom boundary tags, with the values giving
-        a list of the faces on which they should be applied in terms of cardinal
-        directions.
+        a list of the faces on which they should be applied in terms of coordinate
+        directions (+x, -x, +y, -y, +z, -z).
 
     .. versionchanged:: 2017.1
 
@@ -781,28 +781,52 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     # compute facial adjacency for Mesh if there is tag information
     facial_adjacency_groups = None
     face_vertex_indices_to_tags = {}
-    boundary_tags = list(boundary_dict.keys())
+    boundary_tags = list(face_to_boundary_tag.keys())
     from meshmode.mesh import _compute_facial_adjacency_from_vertices
     if boundary_tags:
         for tag_idx, tag in enumerate(boundary_tags):
             # Need to map the correct face vertices to the boundary tags
-            for face in boundary_dict[tag]:
-                if face == "east":
+            for face in face_to_boundary_tag[tag]:
+                if face == "-x":
                     face_id = 1
                     dim_crit = 0
                     node_crit = axis_coords[0][0]
-                elif face == "west":
+                elif face == "+x":
                     face_id = 1
                     dim_crit = 0
                     node_crit = axis_coords[0][-1]
-                elif face == "north":
+                elif face == "-y":
+                    # Check to make sure we are at least 2D.
+                    if dim < 2:
+                        raise ValueError("Attempting to implement y- "
+                                "boundary condition on 1D box mesh.")
                     face_id = 0
                     dim_crit = 1
                     node_crit = axis_coords[1][-1]
-                elif face == "south":
+                elif face == "+y":
+                    # Check to make sure we are at least 2D.
+                    if dim < 2:
+                        raise ValueError("Attempting to implement y- "
+                                "boundary condition on 1D box mesh.")
                     face_id = 0
                     dim_crit = 1
                     node_crit = axis_coords[1][0]
+                elif face == "-z":
+                    # Check to make sure we are at least 3D.
+                    if dim < 3:
+                        raise ValueError("Attempting to implement z- "
+                                "boundary condition on 2D box mesh.")
+                    face_id = 0
+                    dim_crit = 2
+                    node_crit = axis_coords[2][-1]
+                elif face == "+z":
+                    # Check to make sure we are at least 3D.
+                    if dim < 3:
+                        raise ValueError("Attempting to implement z- "
+                                "boundary condition on 2D box mesh.")
+                    face_id = 0
+                    dim_crit = 2
+                    node_crit = axis_coords[2][0]
                 for ielem in range(0, grp.nelements):
                     if grp.nodes[dim_crit][ielem][0] == node_crit:
                         fvi = grp.vertex_indices[ielem,
@@ -830,13 +854,17 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 # {{{ generate_regular_rect_mesh
 
 def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
-                               boundary_dict={}):
+                               face_to_boundary_tag={}):
     """Create a semi-structured rectangular mesh.
 
     :param a: the lower left hand point of the rectangle
     :param b: the upper right hand point of the rectangle
     :param n: a tuple of integers indicating the total number of points
       on [a,b].
+    :param face_to_boundary_tag: an optional dictionary for boundary configuration.
+        The keys correspond to custom boundary tags, with the values giving
+        a list of the faces on which they should be applied in terms of coordinate
+        directions (+x, -x, +y, -y, +z, -z).
     """
     if min(n) < 2:
         raise ValueError("need at least two points in each direction")
@@ -845,7 +873,7 @@ def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
             for a_i, b_i, n_i in zip(a, b, n)]
 
     return generate_box_mesh(axis_coords, order=order,
-                             boundary_dict=boundary_dict)
+                             face_to_boundary_tag=face_to_boundary_tag)
 
 # }}}
 
-- 
GitLab


From 96494e51b80c471ef4ad7376af1fed1cd32a7a6f Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Fri, 15 May 2020 13:12:06 -0500
Subject: [PATCH 07/16] Cuts down new code by using axis list

---
 meshmode/mesh/generation.py | 57 +++++++++++--------------------------
 1 file changed, 17 insertions(+), 40 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 6fc9c8f5..335751e9 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -782,51 +782,28 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     facial_adjacency_groups = None
     face_vertex_indices_to_tags = {}
     boundary_tags = list(face_to_boundary_tag.keys())
+    axes = ["x", "y", "z"]
+    face_ids = [1, 0, 0]
     from meshmode.mesh import _compute_facial_adjacency_from_vertices
     if boundary_tags:
         for tag_idx, tag in enumerate(boundary_tags):
             # Need to map the correct face vertices to the boundary tags
             for face in face_to_boundary_tag[tag]:
-                if face == "-x":
-                    face_id = 1
-                    dim_crit = 0
-                    node_crit = axis_coords[0][0]
-                elif face == "+x":
-                    face_id = 1
-                    dim_crit = 0
-                    node_crit = axis_coords[0][-1]
-                elif face == "-y":
-                    # Check to make sure we are at least 2D.
-                    if dim < 2:
-                        raise ValueError("Attempting to implement y- "
-                                "boundary condition on 1D box mesh.")
-                    face_id = 0
-                    dim_crit = 1
-                    node_crit = axis_coords[1][-1]
-                elif face == "+y":
-                    # Check to make sure we are at least 2D.
-                    if dim < 2:
-                        raise ValueError("Attempting to implement y- "
-                                "boundary condition on 1D box mesh.")
-                    face_id = 0
-                    dim_crit = 1
-                    node_crit = axis_coords[1][0]
-                elif face == "-z":
-                    # Check to make sure we are at least 3D.
-                    if dim < 3:
-                        raise ValueError("Attempting to implement z- "
-                                "boundary condition on 2D box mesh.")
-                    face_id = 0
-                    dim_crit = 2
-                    node_crit = axis_coords[2][-1]
-                elif face == "+z":
-                    # Check to make sure we are at least 3D.
-                    if dim < 3:
-                        raise ValueError("Attempting to implement z- "
-                                "boundary condition on 2D box mesh.")
-                    face_id = 0
-                    dim_crit = 2
-                    node_crit = axis_coords[2][0]
+                for i_ax, axis in enumerate(axes):
+                    if face == "-" + axis:
+                        if dim < i_ax + 1:
+                            raise ValueError("Boundary condition dimension "
+                            "mismatch")
+                        face_id = face_ids[i_ax]
+                        dim_crit = i_ax
+                        node_crit = axis_coords[i_ax][0]
+                    elif face == "+" + axis:
+                        if dim < i_ax + 1:
+                            raise ValueError("Boundary condition dimension "
+                            "mismatch")
+                        face_id = face_ids[i_ax]
+                        dim_crit = i_ax
+                        node_crit = axis_coords[i_ax][-1]
                 for ielem in range(0, grp.nelements):
                     if grp.nodes[dim_crit][ielem][0] == node_crit:
                         fvi = grp.vertex_indices[ielem,
-- 
GitLab


From d1647c9fac3556dbab53e5223c795fc3245d681a Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Fri, 15 May 2020 13:16:20 -0500
Subject: [PATCH 08/16] Adds example of new face_to_boundary_tag input to
 docstring

---
 meshmode/mesh/generation.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 335751e9..d5153ec6 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -669,7 +669,9 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     :param face_to_boundary_tag: an optional dictionary for boundary configuration.
         The keys correspond to custom boundary tags, with the values giving
         a list of the faces on which they should be applied in terms of coordinate
-        directions (+x, -x, +y, -y, +z, -z).
+        directions (+x, -x, +y, -y, +z, -z). One example of this would be:
+
+        face_to_boundary_tag={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
 
     .. versionchanged:: 2017.1
 
-- 
GitLab


From c9835aacd63732a99df6c71c2a615a811e0e267a Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Fri, 15 May 2020 13:33:59 -0500
Subject: [PATCH 09/16] Adds test for 2D box with custom boundary tags

---
 test/test_meshmode.py | 50 +++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 50 insertions(+)

diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index 915b89e8..697e0bbd 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -123,6 +123,56 @@ def test_boundary_tags():
 # }}}
 
 
+# {{{ test custom boundary tags on box mesh
+
+def test_box_boundary_tags():
+    from meshmode.mesh.generation import generate_regular_rect_mesh
+    from meshmode.mesh import is_boundary_tag_empty
+    from meshmode.mesh import check_bc_coverage
+    # Test this capability with a 2D mesh.
+    mesh = generate_regular_rect_mesh(a=(0, -1), b=(1, 1),
+                                      n=(20, 20), order=3,
+                                      face_to_boundary_tag={
+                                          "btag_test_1": ["+x", "-y"],
+                                          "btag_test_2": ["+y", "-x"]})
+
+    # ensure boundary is covered
+    assert not is_boundary_tag_empty(mesh, "btag_test_1")
+    assert not is_boundary_tag_empty(mesh, "btag_test_2")
+    check_bc_coverage(mesh, ['btag_test_1', 'btag_test_2'])
+
+    # correct answers
+    num_on_bdy = 38
+
+    # check how many elements are marked on each boundary
+    num_marked_bdy_1 = 0
+    num_marked_bdy_2 = 0
+    btag_1_bit = mesh.boundary_tag_bit('btag_test_1')
+    btag_2_bit = mesh.boundary_tag_bit('btag_test_2')
+    for igrp in range(len(mesh.groups)):
+        bdry_fagrp = mesh.facial_adjacency_groups[igrp].get(None, None)
+
+        if bdry_fagrp is None:
+            continue
+
+        for i, nbrs in enumerate(bdry_fagrp.neighbors):
+            if (-nbrs) & btag_1_bit:
+                num_marked_bdy_1 += 1
+            if (-nbrs) & btag_2_bit:
+                num_marked_bdy_2 += 1
+
+    # raise errors if wrong number of elements marked
+    if num_marked_bdy_1 != num_on_bdy:
+        raise ValueError("%i marked on custom boundary 1, should be %i" %
+                         (num_marked_bdy_1, num_on_bdy))
+
+    if num_marked_bdy_2 != num_on_bdy:
+        raise ValueError("%i marked on custom boundary 2, should be %i" %
+                         (num_marked_bdy_2, num_on_bdy))
+
+# }}}
+
+
 # {{{ convergence of boundary interpolation
 
 @pytest.mark.parametrize("group_factory", [
-- 
GitLab


From 589969fe854e230c13a1680b3bd444aa27f16071 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Mon, 18 May 2020 15:08:14 -0500
Subject: [PATCH 10/16] Makes a number of changes to the box tag interface and
 test

---
 meshmode/mesh/generation.py | 39 ++++++++++++++++++++++++-------------
 test/test_meshmode.py       | 24 +++++++++++------------
 2 files changed, 37 insertions(+), 26 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 3955fbdd..71c8dfd2 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -673,7 +673,7 @@ def generate_urchin(order, m, n, est_rel_interp_tolerance, min_rad=0.2):
 # {{{ generate_box_mesh
 
 def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
-        group_factory=None, face_to_boundary_tag={}):
+        group_factory=None, boundary_tag_to_face={}):
     """Create a semi-structured mesh.
 
     :param axis_coords: a tuple with a number of entries corresponding
@@ -681,12 +681,12 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         specifying the coordinates to be used along that axis.
     :param group_factory: One of :class:`meshmode.mesh.SimplexElementGroup`
         or :class:`meshmode.mesh.TensorProductElementGroup`.
-    :param face_to_boundary_tag: an optional dictionary for boundary configuration.
+    :param boundary_tag_to_face: an optional dictionary for boundary configuration.
         The keys correspond to custom boundary tags, with the values giving
         a list of the faces on which they should be applied in terms of coordinate
-        directions (+x, -x, +y, -y, +z, -z). One example of this would be:
+        directions (+x, -x, +y, -y, +z, -z, +w, -w). One example of this would be:
 
-        face_to_boundary_tag={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
+        boundary_tag_to_face={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
 
     .. versionchanged:: 2017.1
 
@@ -798,26 +798,39 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     # compute facial adjacency for Mesh if there is tag information
     facial_adjacency_groups = None
     face_vertex_indices_to_tags = {}
-    boundary_tags = list(face_to_boundary_tag.keys())
-    axes = ["x", "y", "z"]
-    face_ids = [1, 0, 0]
+    boundary_tags = list(boundary_tag_to_face.keys())
+    axes = ["x", "y", "z", "w"]
+    face_ids = [1, 0, 0, 0]
     from meshmode.mesh import _compute_facial_adjacency_from_vertices
     if boundary_tags:
         for tag_idx, tag in enumerate(boundary_tags):
             # Need to map the correct face vertices to the boundary tags
-            for face in face_to_boundary_tag[tag]:
+            for face in boundary_tag_to_face[tag]:
+                # Check if face to boundary input is formatted properly
+                if face.startswith("-"):
+                    pass
+                elif face.startswith("+"):
+                    pass
+                else:
+                    raise ValueError("Face ", face, " given for boundary tag ",
+                            tag, " is not formatted properly. Use +x, -x, +y, "
+                            "+z, -z, +w, or -w")
                 for i_ax, axis in enumerate(axes):
                     if face == "-" + axis:
                         if dim < i_ax + 1:
                             raise ValueError("Boundary condition dimension "
-                            "mismatch")
+                            "mismatch: facial dimension ", face, " for tag ",
+                            tag, " is higher than the dimension of the "
+                            "problem")
                         face_id = face_ids[i_ax]
                         dim_crit = i_ax
                         node_crit = axis_coords[i_ax][0]
                     elif face == "+" + axis:
                         if dim < i_ax + 1:
                             raise ValueError("Boundary condition dimension "
-                            "mismatch")
+                            "mismatch: facial dimension ", face, " for tag ",
+                            tag, " is higher than the dimension of the "
+                            "problem")
                         face_id = face_ids[i_ax]
                         dim_crit = i_ax
                         node_crit = axis_coords[i_ax][-1]
@@ -848,14 +861,14 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 # {{{ generate_regular_rect_mesh
 
 def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
-                               face_to_boundary_tag={}):
+                               boundary_tag_to_face={}):
     """Create a semi-structured rectangular mesh.
 
     :param a: the lower left hand point of the rectangle
     :param b: the upper right hand point of the rectangle
     :param n: a tuple of integers indicating the total number of points
       on [a,b].
-    :param face_to_boundary_tag: an optional dictionary for boundary configuration.
+    :param boundary_tag_to_face: an optional dictionary for boundary configuration.
         The keys correspond to custom boundary tags, with the values giving
         a list of the faces on which they should be applied in terms of coordinate
         directions (+x, -x, +y, -y, +z, -z).
@@ -867,7 +880,7 @@ def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
             for a_i, b_i, n_i in zip(a, b, n)]
 
     return generate_box_mesh(axis_coords, order=order,
-                             face_to_boundary_tag=face_to_boundary_tag)
+                             boundary_tag_to_face=boundary_tag_to_face)
 
 # }}}
 
diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index fe6a1ea2..609e426b 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -124,26 +124,25 @@ def test_boundary_tags():
 
 
 # {{{ test custom boundary tags on box mesh
-
-def test_box_boundary_tags():
+@pytest.mark.parametrize(("dim", "nelem"), [
+    (2, 20),
+    ])
+def test_box_boundary_tags(dim, nelem):
     from meshmode.mesh.generation import generate_regular_rect_mesh
     from meshmode.mesh import is_boundary_tag_empty
     from meshmode.mesh import check_bc_coverage
-    # Test this capability with a 2D mesh.
     mesh = generate_regular_rect_mesh(a=(0, -1), b=(1, 1),
-                                      n=(20, 20), order=3,
-                                      face_to_boundary_tag={
-                                          "btag_test_1": ["+x", "-y"],
-                                          "btag_test_2": ["+y", "-x"]})
+                                      n=(nelem, nelem), order=3,
+                                      boundary_tag_to_face={
+                                      "btag_test_1": ["+x", "-y"],
+                                      "btag_test_2": ["+y", "-x"]})
+    # correct answer
+    num_on_bdy = dim*(nelem-1)
 
-    # ensure boundary is covered
     assert not is_boundary_tag_empty(mesh, "btag_test_1")
     assert not is_boundary_tag_empty(mesh, "btag_test_2")
     check_bc_coverage(mesh, ['btag_test_1', 'btag_test_2'])
 
-    # correct answers
-    num_on_bdy = 38
-
     # check how many elements are marked on each boundary
     num_marked_bdy_1 = 0
     num_marked_bdy_2 = 0
@@ -160,16 +159,15 @@ def test_box_boundary_tags():
                 num_marked_bdy_1 += 1
             if (-nbrs) & btag_2_bit:
                 num_marked_bdy_2 += 1
-
     # raise errors if wrong number of elements marked
     if num_marked_bdy_1 != num_on_bdy:
         raise ValueError("%i marked on custom boundary 1, should be %i" %
                          (num_marked_bdy_1, num_on_bdy))
-
     if num_marked_bdy_2 != num_on_bdy:
         raise ValueError("%i marked on custom boundary 2, should be %i" %
                          (num_marked_bdy_2, num_on_bdy))
 
+
 # }}}
 
 
-- 
GitLab


From 25a465c4ed06f9896ed8fc8779b9451d9420a104 Mon Sep 17 00:00:00 2001
From: Cory Mikida <cmikida2@illinois.edu>
Date: Wed, 27 May 2020 15:16:14 -0500
Subject: [PATCH 11/16] Adds 3D meshmode test, further ensures catching of bad
 face tag inputs

---
 meshmode/mesh/generation.py | 49 ++++++++++++++++++++++++++++---------
 test/test_meshmode.py       | 24 +++++++++++++-----
 2 files changed, 55 insertions(+), 18 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 71c8dfd2..5dd3b22e 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -728,6 +728,9 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 
     el_vertices = []
 
+    vertex_indices_left = {}
+    vertex_indices_right = {}
+
     if dim == 1:
         for i in range(shape[0]-1):
             # a--b
@@ -738,6 +741,11 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
             el_vertices.append((a, b,))
 
     elif dim == 2:
+        # For custom boundary tagging.
+        vertex_indices_left["x"] = vertex_indices[0, :].flatten()
+        vertex_indices_right["x"] = vertex_indices[-1, :].flatten()
+        vertex_indices_left["y"] = vertex_indices[:, 0].flatten()
+        vertex_indices_right["y"] = vertex_indices[:, -1].flatten()
         for i in range(shape[0]-1):
             for j in range(shape[1]-1):
 
@@ -757,6 +765,13 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
                     el_vertices.append((d, c, b))
 
     elif dim == 3:
+        # For custom boundary tagging.
+        vertex_indices_left["x"] = vertex_indices[0, :, :].flatten()
+        vertex_indices_right["x"] = vertex_indices[-1, :, :].flatten()
+        vertex_indices_left["y"] = vertex_indices[:, 0, :].flatten()
+        vertex_indices_right["y"] = vertex_indices[:, -1, :].flatten()
+        vertex_indices_left["z"] = vertex_indices[:, :, 0].flatten()
+        vertex_indices_right["z"] = vertex_indices[:, :, -1].flatten()
         for i in range(shape[0]-1):
             for j in range(shape[1]-1):
                 for k in range(shape[2]-1):
@@ -800,7 +815,6 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
     face_vertex_indices_to_tags = {}
     boundary_tags = list(boundary_tag_to_face.keys())
     axes = ["x", "y", "z", "w"]
-    face_ids = [1, 0, 0, 0]
     from meshmode.mesh import _compute_facial_adjacency_from_vertices
     if boundary_tags:
         for tag_idx, tag in enumerate(boundary_tags):
@@ -815,6 +829,18 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
                     raise ValueError("Face ", face, " given for boundary tag ",
                             tag, " is not formatted properly. Use +x, -x, +y, "
                             "+z, -z, +w, or -w")
+                if face.endswith("x"):
+                    pass
+                elif face.endswith("y"):
+                    pass
+                elif face.endswith("z"):
+                    pass
+                elif face.endswith("w"):
+                    pass
+                else:
+                    raise ValueError("Face ", face, " given for boundary tag ",
+                            tag, " is not formatted properly. Use +x, -x, +y, "
+                            "+z, -z, +w, or -w")
                 for i_ax, axis in enumerate(axes):
                     if face == "-" + axis:
                         if dim < i_ax + 1:
@@ -822,26 +848,25 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
                             "mismatch: facial dimension ", face, " for tag ",
                             tag, " is higher than the dimension of the "
                             "problem")
-                        face_id = face_ids[i_ax]
-                        dim_crit = i_ax
-                        node_crit = axis_coords[i_ax][0]
+                        vert_crit = vertex_indices_left[axis]
                     elif face == "+" + axis:
                         if dim < i_ax + 1:
                             raise ValueError("Boundary condition dimension "
                             "mismatch: facial dimension ", face, " for tag ",
                             tag, " is higher than the dimension of the "
                             "problem")
-                        face_id = face_ids[i_ax]
-                        dim_crit = i_ax
-                        node_crit = axis_coords[i_ax][-1]
+                        vert_crit = vertex_indices_right[axis]
                 for ielem in range(0, grp.nelements):
-                    if grp.nodes[dim_crit][ielem][0] == node_crit:
+                    # Loop through potential boundary-facing faces.
+                    for i in range(0, dim+1):
+                        # Element vertices:
                         fvi = grp.vertex_indices[ielem,
                                                  grp.face_vertex_indices()[
-                                                     face_id]]
-                        if frozenset(fvi) not in face_vertex_indices_to_tags:
-                            face_vertex_indices_to_tags[frozenset(fvi)] = []
-                        face_vertex_indices_to_tags[frozenset(fvi)].append(tag)
+                                                     i]]
+                        if all(ind in vert_crit for ind in fvi):
+                            if frozenset(fvi) not in face_vertex_indices_to_tags:
+                                face_vertex_indices_to_tags[frozenset(fvi)] = []
+                            face_vertex_indices_to_tags[frozenset(fvi)].append(tag)
         # Need to add BTAG ALL and BTAG_REALLY_ALL to list of tags
         from meshmode.mesh import BTAG_ALL, BTAG_REALLY_ALL
         boundary_tags.append(BTAG_ALL)
diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index 609e426b..a0248712 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -126,18 +126,29 @@ def test_boundary_tags():
 # {{{ test custom boundary tags on box mesh
 @pytest.mark.parametrize(("dim", "nelem"), [
     (2, 20),
+    (3, 20),
     ])
 def test_box_boundary_tags(dim, nelem):
     from meshmode.mesh.generation import generate_regular_rect_mesh
     from meshmode.mesh import is_boundary_tag_empty
     from meshmode.mesh import check_bc_coverage
-    mesh = generate_regular_rect_mesh(a=(0, -1), b=(1, 1),
-                                      n=(nelem, nelem), order=3,
-                                      boundary_tag_to_face={
-                                      "btag_test_1": ["+x", "-y"],
-                                      "btag_test_2": ["+y", "-x"]})
+    if dim == 2:
+        a = (0, -1)
+        b = (1, 1)
+        n = (nelem, nelem)
+        btag_to_face = {"btag_test_1": ["+x", "-y"],
+                        "btag_test_2": ["+y", "-x"]}
+    elif dim == 3:
+        a = (0, -1, -1)
+        b = (1, 1, 1)
+        n = (nelem, nelem, nelem)
+        btag_to_face = {"btag_test_1": ["+x", "-y", "-z"],
+                        "btag_test_2": ["+y", "-x", "+z"]}
+    mesh = generate_regular_rect_mesh(a=a, b=b,
+                                      n=n, order=3,
+                                      boundary_tag_to_face=btag_to_face)
     # correct answer
-    num_on_bdy = dim*(nelem-1)
+    num_on_bdy = dim*(dim-1)*(nelem-1)**(dim-1)
 
     assert not is_boundary_tag_empty(mesh, "btag_test_1")
     assert not is_boundary_tag_empty(mesh, "btag_test_2")
@@ -159,6 +170,7 @@ def test_box_boundary_tags(dim, nelem):
                 num_marked_bdy_1 += 1
             if (-nbrs) & btag_2_bit:
                 num_marked_bdy_2 += 1
+
     # raise errors if wrong number of elements marked
     if num_marked_bdy_1 != num_on_bdy:
         raise ValueError("%i marked on custom boundary 1, should be %i" %
-- 
GitLab


From 7bf6643ff7bb86deb17176d3008cb91c8dac06b6 Mon Sep 17 00:00:00 2001
From: "[6~" <inform@tiker.net>
Date: Thu, 28 May 2020 22:54:37 -0500
Subject: [PATCH 12/16] test_box_boundary_tags: Add 1D test

---
 test/test_meshmode.py | 18 ++++++++++++++----
 1 file changed, 14 insertions(+), 4 deletions(-)

diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index 03863784..e1b829eb 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -125,6 +125,7 @@ def test_boundary_tags():
 
 # {{{ test custom boundary tags on box mesh
 @pytest.mark.parametrize(("dim", "nelem"), [
+    (1, 20),
     (2, 20),
     (3, 20),
     ])
@@ -132,7 +133,13 @@ def test_box_boundary_tags(dim, nelem):
     from meshmode.mesh.generation import generate_regular_rect_mesh
     from meshmode.mesh import is_boundary_tag_empty
     from meshmode.mesh import check_bc_coverage
-    if dim == 2:
+    if dim == 1:
+        a = (0,)
+        b = (1,)
+        n = (nelem,)
+        btag_to_face = {"btag_test_1": ["+x"],
+                        "btag_test_2": ["-x"]}
+    elif dim == 2:
         a = (0, -1)
         b = (1, 1)
         n = (nelem, nelem)
@@ -148,7 +155,10 @@ def test_box_boundary_tags(dim, nelem):
                                       n=n, order=3,
                                       boundary_tag_to_face=btag_to_face)
     # correct answer
-    num_on_bdy = dim*(dim-1)*(nelem-1)**(dim-1)
+    if dim == 1:
+        num_on_bdy = 1
+    else:
+        num_on_bdy = dim*(dim-1)*(nelem-1)**(dim-1)
 
     assert not is_boundary_tag_empty(mesh, "btag_test_1")
     assert not is_boundary_tag_empty(mesh, "btag_test_2")
@@ -157,8 +167,8 @@ def test_box_boundary_tags(dim, nelem):
     # check how many elements are marked on each boundary
     num_marked_bdy_1 = 0
     num_marked_bdy_2 = 0
-    btag_1_bit = mesh.boundary_tag_bit('btag_test_1')
-    btag_2_bit = mesh.boundary_tag_bit('btag_test_2')
+    btag_1_bit = mesh.boundary_tag_bit("btag_test_1")
+    btag_2_bit = mesh.boundary_tag_bit("btag_test_2")
     for igrp in range(len(mesh.groups)):
         bdry_fagrp = mesh.facial_adjacency_groups[igrp].get(None, None)
 
-- 
GitLab


From ce77b6c262e6f444469e8baf3516ee73a866f4c8 Mon Sep 17 00:00:00 2001
From: "[6~" <inform@tiker.net>
Date: Thu, 28 May 2020 22:55:21 -0500
Subject: [PATCH 13/16] generate_box_mesh+BCs: Simplify logic, don't hardcode
 face count

---
 meshmode/mesh/generation.py | 139 ++++++++++++++++--------------------
 1 file changed, 62 insertions(+), 77 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 5dd3b22e..7fe8d1a4 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -673,7 +673,7 @@ def generate_urchin(order, m, n, est_rel_interp_tolerance, min_rad=0.2):
 # {{{ generate_box_mesh
 
 def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
-        group_factory=None, boundary_tag_to_face={}):
+        group_factory=None, boundary_tag_to_face=None):
     """Create a semi-structured mesh.
 
     :param axis_coords: a tuple with a number of entries corresponding
@@ -681,18 +681,26 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         specifying the coordinates to be used along that axis.
     :param group_factory: One of :class:`meshmode.mesh.SimplexElementGroup`
         or :class:`meshmode.mesh.TensorProductElementGroup`.
-    :param boundary_tag_to_face: an optional dictionary for boundary configuration.
+    :param boundary_tag_to_face: an optional dictionary for tagging boundaries.
         The keys correspond to custom boundary tags, with the values giving
         a list of the faces on which they should be applied in terms of coordinate
-        directions (+x, -x, +y, -y, +z, -z, +w, -w). One example of this would be:
+        directions (``+x``, ``-x``, ``+y``, ``-y``, ``+z``, ``-z``, ``+w``, ``-w``).
+        For example::
 
-        boundary_tag_to_face={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
+            boundary_tag_to_face={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
 
     .. versionchanged:: 2017.1
 
         *group_factory* parameter added.
+
+    .. versionchanged:: 2020.1
+
+        *boundary_tag_to_face* parameter added.
     """
 
+    if boundary_tag_to_face is None:
+        boundary_tag_to_face = {}
+
     for iaxis, axc in enumerate(axis_coords):
         if len(axc) < 2:
             raise ValueError("need at least two points along axis %d"
@@ -728,9 +736,6 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 
     el_vertices = []
 
-    vertex_indices_left = {}
-    vertex_indices_right = {}
-
     if dim == 1:
         for i in range(shape[0]-1):
             # a--b
@@ -741,11 +746,6 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
             el_vertices.append((a, b,))
 
     elif dim == 2:
-        # For custom boundary tagging.
-        vertex_indices_left["x"] = vertex_indices[0, :].flatten()
-        vertex_indices_right["x"] = vertex_indices[-1, :].flatten()
-        vertex_indices_left["y"] = vertex_indices[:, 0].flatten()
-        vertex_indices_right["y"] = vertex_indices[:, -1].flatten()
         for i in range(shape[0]-1):
             for j in range(shape[1]-1):
 
@@ -765,13 +765,6 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
                     el_vertices.append((d, c, b))
 
     elif dim == 3:
-        # For custom boundary tagging.
-        vertex_indices_left["x"] = vertex_indices[0, :, :].flatten()
-        vertex_indices_right["x"] = vertex_indices[-1, :, :].flatten()
-        vertex_indices_left["y"] = vertex_indices[:, 0, :].flatten()
-        vertex_indices_right["y"] = vertex_indices[:, -1, :].flatten()
-        vertex_indices_left["z"] = vertex_indices[:, :, 0].flatten()
-        vertex_indices_right["z"] = vertex_indices[:, :, -1].flatten()
         for i in range(shape[0]-1):
             for j in range(shape[1]-1):
                 for k in range(shape[2]-1):
@@ -801,8 +794,7 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
                         el_vertices.append((a011, a111, a101, a110))
 
     else:
-        raise NotImplementedError("box meshes of dimension %d"
-                % dim)
+        raise NotImplementedError("box meshes of dimension %d" % dim)
 
     el_vertices = np.array(el_vertices, dtype=np.int32)
 
@@ -810,70 +802,63 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
             vertices.reshape(dim, -1), el_vertices, order,
             group_factory=group_factory)
 
-    # compute facial adjacency for Mesh if there is tag information
+    # {{{ compute facial adjacency for mesh if there is tag information
+
     facial_adjacency_groups = None
     face_vertex_indices_to_tags = {}
     boundary_tags = list(boundary_tag_to_face.keys())
     axes = ["x", "y", "z", "w"]
-    from meshmode.mesh import _compute_facial_adjacency_from_vertices
+
     if boundary_tags:
-        for tag_idx, tag in enumerate(boundary_tags):
-            # Need to map the correct face vertices to the boundary tags
-            for face in boundary_tag_to_face[tag]:
-                # Check if face to boundary input is formatted properly
-                if face.startswith("-"):
-                    pass
-                elif face.startswith("+"):
-                    pass
-                else:
-                    raise ValueError("Face ", face, " given for boundary tag ",
-                            tag, " is not formatted properly. Use +x, -x, +y, "
-                            "+z, -z, +w, or -w")
-                if face.endswith("x"):
-                    pass
-                elif face.endswith("y"):
-                    pass
-                elif face.endswith("z"):
-                    pass
-                elif face.endswith("w"):
-                    pass
-                else:
-                    raise ValueError("Face ", face, " given for boundary tag ",
-                            tag, " is not formatted properly. Use +x, -x, +y, "
-                            "+z, -z, +w, or -w")
-                for i_ax, axis in enumerate(axes):
-                    if face == "-" + axis:
-                        if dim < i_ax + 1:
-                            raise ValueError("Boundary condition dimension "
-                            "mismatch: facial dimension ", face, " for tag ",
-                            tag, " is higher than the dimension of the "
-                            "problem")
-                        vert_crit = vertex_indices_left[axis]
-                    elif face == "+" + axis:
-                        if dim < i_ax + 1:
-                            raise ValueError("Boundary condition dimension "
-                            "mismatch: facial dimension ", face, " for tag ",
-                            tag, " is higher than the dimension of the "
-                            "problem")
-                        vert_crit = vertex_indices_right[axis]
-                for ielem in range(0, grp.nelements):
-                    # Loop through potential boundary-facing faces.
-                    for i in range(0, dim+1):
-                        # Element vertices:
-                        fvi = grp.vertex_indices[ielem,
-                                                 grp.face_vertex_indices()[
-                                                     i]]
-                        if all(ind in vert_crit for ind in fvi):
-                            if frozenset(fvi) not in face_vertex_indices_to_tags:
-                                face_vertex_indices_to_tags[frozenset(fvi)] = []
-                            face_vertex_indices_to_tags[frozenset(fvi)].append(tag)
-        # Need to add BTAG ALL and BTAG_REALLY_ALL to list of tags
-        from meshmode.mesh import BTAG_ALL, BTAG_REALLY_ALL
-        boundary_tags.append(BTAG_ALL)
-        boundary_tags.append(BTAG_REALLY_ALL)
+        from pytools import indices_in_shape
+        vert_index_to_tuple = {
+                vertex_indices[itup]: itup
+                for itup in indices_in_shape(shape)}
+
+    for tag_idx, tag in enumerate(boundary_tags):
+        # Need to map the correct face vertices to the boundary tags
+        for face in boundary_tag_to_face[tag]:
+            if len(face) != 2:
+                raise ValueError("face identifier '%s' does not "
+                        "consist of exactly two characters" % face)
+
+            side, axis = face
+            try:
+                axis = axes.index(axis)
+            except ValueError:
+                raise ValueError("unrecognized axis in face identifier '%s'" % face)
+            if axis >= dim:
+                raise ValueError("axis in face identifier '%s' does not exist in %dD"
+                        % (face, dim))
+
+            if side == "-":
+                vert_crit = 0
+            elif side == "+":
+                vert_crit = shape[axis] - 1
+            else:
+                raise ValueError("first character of face identifier '%s' is not"
+                        "'+' or '-'" % face)
+
+            for ielem in range(0, grp.nelements):
+                for ref_fvi in grp.face_vertex_indices():
+                    fvi = grp.vertex_indices[ielem, ref_fvi]
+                    fvi_tuples = [vert_index_to_tuple[i] for i in fvi]
+
+                    if all(fvi_tuple[axis] == vert_crit for fvi_tuple in fvi_tuples):
+                        key = frozenset(fvi)
+                        face_vertex_indices_to_tags.setdefault(key, []).append(tag)
+
+    if boundary_tags:
+        from meshmode.mesh import (
+                _compute_facial_adjacency_from_vertices, BTAG_ALL, BTAG_REALLY_ALL)
+        boundary_tags.extend([BTAG_ALL, BTAG_REALLY_ALL])
         facial_adjacency_groups = _compute_facial_adjacency_from_vertices(
                 [grp], boundary_tags, np.int32, np.int8,
                 face_vertex_indices_to_tags)
+    else:
+        facial_adjacency_groups = None
+
+    # }}}
 
     from meshmode.mesh import Mesh
     return Mesh(vertices, [grp],
-- 
GitLab


From 4c9c25daac6e3b724c71597cdb80a86b8ef6aa77 Mon Sep 17 00:00:00 2001
From: "[6~" <inform@tiker.net>
Date: Thu, 28 May 2020 23:00:59 -0500
Subject: [PATCH 14/16] generate_regular_rect_mesh: Allow passing group_factory

---
 meshmode/mesh/generation.py | 13 +++++++------
 1 file changed, 7 insertions(+), 6 deletions(-)

diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py
index 7fe8d1a4..ca5ddfb4 100644
--- a/meshmode/mesh/generation.py
+++ b/meshmode/mesh/generation.py
@@ -685,6 +685,7 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
         The keys correspond to custom boundary tags, with the values giving
         a list of the faces on which they should be applied in terms of coordinate
         directions (``+x``, ``-x``, ``+y``, ``-y``, ``+z``, ``-z``, ``+w``, ``-w``).
+
         For example::
 
             boundary_tag_to_face={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
@@ -871,17 +872,16 @@ def generate_box_mesh(axis_coords, order=1, coord_dtype=np.float64,
 # {{{ generate_regular_rect_mesh
 
 def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
-                               boundary_tag_to_face={}):
+                               boundary_tag_to_face=None,
+                               group_factory=None):
     """Create a semi-structured rectangular mesh.
 
     :param a: the lower left hand point of the rectangle
     :param b: the upper right hand point of the rectangle
     :param n: a tuple of integers indicating the total number of points
       on [a,b].
-    :param boundary_tag_to_face: an optional dictionary for boundary configuration.
-        The keys correspond to custom boundary tags, with the values giving
-        a list of the faces on which they should be applied in terms of coordinate
-        directions (+x, -x, +y, -y, +z, -z).
+    :param boundary_tag_to_face: an optional dictionary for tagging boundaries.
+        See :func:`generate_box_mesh`.
     """
     if min(n) < 2:
         raise ValueError("need at least two points in each direction")
@@ -890,7 +890,8 @@ def generate_regular_rect_mesh(a=(0, 0), b=(1, 1), n=(5, 5), order=1,
             for a_i, b_i, n_i in zip(a, b, n)]
 
     return generate_box_mesh(axis_coords, order=order,
-                             boundary_tag_to_face=boundary_tag_to_face)
+                             boundary_tag_to_face=boundary_tag_to_face,
+                             group_factory=group_factory)
 
 # }}}
 
-- 
GitLab


From c4a07b960074cef22f8bab5f2dfaaee371b41fb2 Mon Sep 17 00:00:00 2001
From: "[6~" <inform@tiker.net>
Date: Thu, 28 May 2020 23:01:33 -0500
Subject: [PATCH 15/16] test_box_boundary_tags: Prepare for tensor product
 elements

---
 test/test_meshmode.py | 15 ++++++++++++---
 1 file changed, 12 insertions(+), 3 deletions(-)

diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index e1b829eb..b408175d 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -33,6 +33,7 @@ from pyopencl.tools import (  # noqa
         pytest_generate_tests_for_pyopencl
         as pytest_generate_tests)
 
+from meshmode.mesh import SimplexElementGroup, TensorProductElementGroup
 from meshmode.discretization.poly_element import (
         InterpolatoryQuadratureSimplexGroupFactory,
         PolynomialWarpAndBlendGroupFactory,
@@ -124,12 +125,19 @@ def test_boundary_tags():
 
 
 # {{{ test custom boundary tags on box mesh
+
 @pytest.mark.parametrize(("dim", "nelem"), [
     (1, 20),
     (2, 20),
-    (3, 20),
+    (3, 10),
+    ])
+@pytest.mark.parametrize("group_factory", [
+    SimplexElementGroup,
+
+    # FIXME: Not implemented: TPE.face_vertex_indices
+    # TensorProductElementGroup
     ])
-def test_box_boundary_tags(dim, nelem):
+def test_box_boundary_tags(dim, nelem, group_factory):
     from meshmode.mesh.generation import generate_regular_rect_mesh
     from meshmode.mesh import is_boundary_tag_empty
     from meshmode.mesh import check_bc_coverage
@@ -153,7 +161,8 @@ def test_box_boundary_tags(dim, nelem):
                         "btag_test_2": ["+y", "-x", "+z"]}
     mesh = generate_regular_rect_mesh(a=a, b=b,
                                       n=n, order=3,
-                                      boundary_tag_to_face=btag_to_face)
+                                      boundary_tag_to_face=btag_to_face,
+                                      group_factory=group_factory)
     # correct answer
     if dim == 1:
         num_on_bdy = 1
-- 
GitLab


From 6513f19ffd26a8ad2692af2382ca65e85ad2ac94 Mon Sep 17 00:00:00 2001
From: "[6~" <inform@tiker.net>
Date: Fri, 29 May 2020 00:10:45 -0500
Subject: [PATCH 16/16] Placate flake8

---
 test/test_meshmode.py | 5 +----
 1 file changed, 1 insertion(+), 4 deletions(-)

diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index b408175d..ded91fe2 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -39,7 +39,7 @@ from meshmode.discretization.poly_element import (
         PolynomialWarpAndBlendGroupFactory,
         PolynomialEquidistantSimplexGroupFactory,
         )
-from meshmode.mesh import BTAG_ALL
+from meshmode.mesh import Mesh, BTAG_ALL
 from meshmode.discretization.connection import \
         FACE_RESTR_ALL, FACE_RESTR_INTERIOR
 import meshmode.mesh.generation as mgen
@@ -606,7 +606,6 @@ def test_3d_orientation(ctx_factory, what, mesh_gen_func, visualize=False):
 def test_merge_and_map(ctx_factory, visualize=False):
     from meshmode.mesh.io import generate_gmsh, FileSource
     from meshmode.mesh.generation import generate_box_mesh
-    from meshmode.mesh import TensorProductElementGroup
     from meshmode.discretization.poly_element import (
             PolynomialWarpAndBlendGroupFactory,
             LegendreGaussLobattoTensorProductGroupFactory)
@@ -1088,7 +1087,6 @@ def no_test_quad_mesh_3d():
 
 def test_quad_single_element():
     from meshmode.mesh.generation import make_group_from_vertices
-    from meshmode.mesh import Mesh, TensorProductElementGroup
 
     vertices = np.array([
                 [0.91, 1.10],
@@ -1116,7 +1114,6 @@ def test_quad_single_element():
 
 def test_quad_multi_element():
     from meshmode.mesh.generation import generate_box_mesh
-    from meshmode.mesh import TensorProductElementGroup
     mesh = generate_box_mesh(
             (
                 np.linspace(3, 8, 4),
-- 
GitLab