From 8c8c89cecfba035d602c8f6c69511a7a8485a94b Mon Sep 17 00:00:00 2001
From: Thomas Gibson <gibsonthomas1120@hotmail.com>
Date: Wed, 24 Mar 2021 15:25:40 -0500
Subject: [PATCH] Reunify discretization types

---
 doc/discretization.rst                        |   4 -
 examples/plot-connectivity.py                 |   4 +-
 examples/simple-dg.py                         |   4 +-
 examples/to_firedrake.py                      |   4 +-
 meshmode/array_context.py                     |   4 +-
 meshmode/discretization/__init__.py           | 141 +++++++++++--
 meshmode/discretization/connection/direct.py  |   4 +-
 meshmode/discretization/connection/modal.py   |  76 +++----
 .../discretization/connection/refinement.py   |   4 +-
 meshmode/discretization/modal.py              |  61 ------
 meshmode/discretization/nodal.py              | 185 ------------------
 meshmode/dof_array.py                         |   6 +-
 meshmode/interop/firedrake/connection.py      |  18 +-
 meshmode/interop/nodal_dg.py                  |   6 +-
 meshmode/mesh/visualization.py                |   4 +-
 test/test_array.py                            |   6 +-
 test/test_chained.py                          |   4 +-
 test/test_firedrake_interop.py                |   8 +-
 test/test_mesh.py                             |   4 +-
 test/test_meshmode.py                         |  51 +++--
 test/test_modal.py                            |  19 +-
 test/test_partition.py                        |   8 +-
 test/test_refinement.py                       |   4 +-
 23 files changed, 248 insertions(+), 381 deletions(-)
 delete mode 100644 meshmode/discretization/modal.py
 delete mode 100644 meshmode/discretization/nodal.py

diff --git a/doc/discretization.rst b/doc/discretization.rst
index 77073db5..9389cd85 100644
--- a/doc/discretization.rst
+++ b/doc/discretization.rst
@@ -6,10 +6,6 @@ Abstract interface
 
 .. automodule:: meshmode.discretization
 
-.. automodule:: meshmode.discretization.nodal
-
-.. automodule:: meshmode.discretization.modal
-
 Element Groups for Composite Polynomial Discretization
 ------------------------------------------------------
 
diff --git a/examples/plot-connectivity.py b/examples/plot-connectivity.py
index a96fc6f7..0faaa9b3 100644
--- a/examples/plot-connectivity.py
+++ b/examples/plot-connectivity.py
@@ -18,11 +18,11 @@ def main():
     mesh = generate_icosahedron(1, order=order)
     #mesh = generate_torus(3, 1, order=order)
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.poly_element import \
             PolynomialWarpAndBlendGroupFactory
 
-    discr = NodalDiscretization(
+    discr = Discretization(
             actx, mesh, PolynomialWarpAndBlendGroupFactory(order))
 
     from meshmode.discretization.visualization import make_visualizer
diff --git a/examples/simple-dg.py b/examples/simple-dg.py
index cc04ef3f..e1d63ef6 100644
--- a/examples/simple-dg.py
+++ b/examples/simple-dg.py
@@ -59,11 +59,11 @@ class DGDiscretization:
     def __init__(self, actx, mesh, order):
         self.order = order
 
-        from meshmode.discretization.nodal import NodalDiscretization
+        from meshmode.discretization import Discretization
         from meshmode.discretization.poly_element import \
                 PolynomialWarpAndBlendGroupFactory
         self.group_factory = PolynomialWarpAndBlendGroupFactory(order=order)
-        self.volume_discr = NodalDiscretization(actx, mesh, self.group_factory)
+        self.volume_discr = Discretization(actx, mesh, self.group_factory)
 
         assert self.volume_discr.dim == 2
 
diff --git a/examples/to_firedrake.py b/examples/to_firedrake.py
index 0a961c98..395211ca 100644
--- a/examples/to_firedrake.py
+++ b/examples/to_firedrake.py
@@ -61,11 +61,11 @@ def main():
 
     order = 3
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.poly_element import \
         InterpolatoryQuadratureSimplexGroupFactory
     group_factory = InterpolatoryQuadratureSimplexGroupFactory(order=order)
-    discr = NodalDiscretization(actx, mesh, group_factory)
+    discr = Discretization(actx, mesh, group_factory)
 
     # Get our solution: we will use
     # Real(e^z) = Real(e^{x+iy})
diff --git a/meshmode/array_context.py b/meshmode/array_context.py
index 20bfedad..29c0042e 100644
--- a/meshmode/array_context.py
+++ b/meshmode/array_context.py
@@ -176,8 +176,8 @@ class _BaseFakeNumpyLinalgNamespace:
 
 
 class ArrayContext:
-    """An interface that allows a subclass of
-    :class:`meshmode.discretization.DiscretizationBase` to create and interact
+    """An interface that allows a
+    :class:`~meshmode.discretization.Discretization` to create and interact
     with arrays of degrees of freedom without fully specifying their types.
 
     .. automethod:: empty
diff --git a/meshmode/discretization/__init__.py b/meshmode/discretization/__init__.py
index 0bdf5151..ebd5a5b6 100644
--- a/meshmode/discretization/__init__.py
+++ b/meshmode/discretization/__init__.py
@@ -26,8 +26,9 @@ THE SOFTWARE.
 import numpy as np
 
 from abc import ABCMeta, abstractproperty, abstractmethod
-from pytools import memoize_method
-from meshmode.array_context import ArrayContext
+from pytools import memoize_in, memoize_method
+from pytools.obj_array import make_obj_array
+from meshmode.array_context import ArrayContext, make_loopy_program
 
 # underscored because it shouldn't be imported from here.
 from meshmode.dof_array import DOFArray as _DOFArray
@@ -37,7 +38,7 @@ __doc__ = """
 .. autoclass:: ModalElementGroupBase
 .. autoclass:: NodalElementGroupBase
 .. autoclass:: InterpolatoryElementGroupBase
-.. autoclass:: DiscretizationBase
+.. autoclass:: Discretization
 """
 
 
@@ -48,7 +49,7 @@ class NoninterpolatoryElementGroupError(TypeError):
 
 
 class ElementGroupBase(object, metaclass=ABCMeta):
-    """Container for data of any subclass of :class:`DiscretizationBase`
+    """Container for data of any :class:`Discretization`
     corresponding to one :class:`meshmode.mesh.MeshElementGroup`.
 
     .. attribute :: mesh_el_group
@@ -140,7 +141,7 @@ class ElementGroupBase(object, metaclass=ABCMeta):
 # {{{ Nodal element group base
 
 class NodalElementGroupBase(ElementGroupBase):
-    """Container for :class:`meshmode.discretization.nodal.NodalDiscretization`
+    """Container for :class:`meshmode.discretization.Discretization`
     data corresponding to one :class:`meshmode.mesh.MeshElementGroup`.
 
     .. autoattribute:: nunit_dofs
@@ -231,7 +232,7 @@ class InterpolatoryElementGroupBase(NodalElementGroupBase):
 # {{{ modal element group base
 
 class ModalElementGroupBase(ElementGroupBase):
-    """Container for :class:`meshmode.discretization.modal.ModalDiscretization`
+    """Container for :class:`meshmode.discretization.Discretization`
     data corresponding to one :class:`meshmode.mesh.MeshElementGroup`. This
     is a subclass of :class:`ElementGroupBase`, equipped with a function space
     and a hierarchical basis that is orthonormal with respect to the
@@ -286,10 +287,10 @@ class ModalElementGroupBase(ElementGroupBase):
 # }}}
 
 
-# {{{ discretization base
+# {{{ discretization
 
-class DiscretizationBase:
-    """An unstructured composite discretization base class.
+class Discretization:
+    """An unstructured composite discretization class.
 
     .. attribute:: real_dtype
     .. attribute:: complex_dtype
@@ -304,6 +305,9 @@ class DiscretizationBase:
     .. automethod:: zeros
     .. automethod:: empty_like
     .. automethod:: zeros_like
+    .. automethod:: num_reference_derivative
+    .. automethod:: quad_weights
+    .. automethod:: coordinates
     """
 
     def __init__(self, actx: ArrayContext, mesh, group_factory,
@@ -407,11 +411,122 @@ class DiscretizationBase:
     def zeros_like(self, array: _DOFArray):
         return self.zeros(array.array_context, dtype=array.entry_dtype)
 
-# }}}
+    def num_reference_derivative(self, ref_axes, vec):
+        actx = vec.array_context
+        ref_axes = list(ref_axes)
+
+        if any([isinstance(grp, ModalElementGroupBase) for grp in self.groups]):
+            raise NotImplementedError("Differentiation matrices not currently "
+                                      "implemented for modal element groups yet.")
+
+        @memoize_in(actx, (Discretization, "reference_derivative_prg"))
+        def prg():
+            return make_loopy_program(
+                """{[iel,idof,j]:
+                    0<=iel<nelements and
+                    0<=idof,j<nunit_dofs}""",
+                "result[iel,idof] = sum(j, diff_mat[idof, j] * vec[iel, j])",
+                name="diff")
+
+        def get_mat(grp):
+            mat = None
+            for ref_axis in ref_axes:
+                next_mat = grp.diff_matrices()[ref_axis]
+                if mat is None:
+                    mat = next_mat
+                else:
+                    mat = np.dot(next_mat, mat)
+
+            return mat
+
+        return _DOFArray(actx, tuple(
+                actx.call_loopy(
+                    prg(), diff_mat=actx.from_numpy(get_mat(grp)), vec=vec[grp.index]
+                    )["result"]
+                for grp in self.groups))
+
+    @memoize_method
+    def quad_weights(self):
+        """:returns: A :class:`~meshmode.dof_array.DOFArray` with quadrature weights.
+        """
+        actx = self._setup_actx
+
+        if any([isinstance(grp, ModalElementGroupBase) for grp in self.groups]):
+            raise TypeError("Cannot compute quadrature weights for "
+                            "modal element groups")
+
+        @memoize_in(actx, (Discretization, "quad_weights_prg"))
+        def prg():
+            return make_loopy_program(
+                "{[iel,idof]: 0<=iel<nelements and 0<=idof<nunit_dofs}",
+                "result[iel,idof] = weights[idof]",
+                name="quad_weights")
+
+        return _DOFArray(None, tuple(
+                actx.freeze(
+                    actx.call_loopy(
+                        prg(),
+                        weights=actx.from_numpy(grp.weights),
+                        nelements=grp.nelements,
+                        )["result"])
+                for grp in self.groups))
 
+    @memoize_method
+    def coordinates(self):
+        r"""
+        :returns: object array of shape ``(ambient_dim,)`` containing
+            :class:`~meshmode.dof_array.DOFArray`\ s of global coordinates.
+        """
 
-# For backwards compatibility, we need to be sure we export the
-# alias for NodalDiscretization
-from meshmode.discretization.nodal import Discretization  # noqa: F401
+        actx = self._setup_actx
+
+        if any([isinstance(grp, ModalElementGroupBase) for grp in self.groups]):
+            raise TypeError("Modal element groups do not have `unit_nodes`")
+
+        @memoize_in(actx, (Discretization, "coordinates_prg"))
+        def prg():
+            return make_loopy_program(
+                """{[iel,idof,j]:
+                    0<=iel<nelements and
+                    0<=idof<ndiscr_nodes and
+                    0<=j<nmesh_nodes}""",
+                """
+                    result[iel, idof] = \
+                        sum(j, resampling_mat[idof, j] * nodes[iel, j])
+                    """,
+                name="coordinates")
+
+        def resample_mesh_nodes(grp, iaxis):
+            # TODO: would be nice to have the mesh use an array context already
+            nodes = actx.from_numpy(grp.mesh_el_group.nodes[iaxis])
+
+            grp_unit_nodes = grp.unit_nodes.reshape(-1)
+            meg_unit_nodes = grp.mesh_el_group.unit_nodes.reshape(-1)
+
+            tol = 10 * np.finfo(grp_unit_nodes.dtype).eps
+            if (grp_unit_nodes.shape == meg_unit_nodes.shape
+                    and np.linalg.norm(grp_unit_nodes - meg_unit_nodes) < tol):
+                return nodes
+
+            return actx.call_loopy(
+                    prg(),
+                    resampling_mat=actx.from_numpy(grp.from_mesh_interp_matrix()),
+                    nodes=nodes,
+                    )["result"]
+
+        return make_obj_array([
+            _DOFArray(None, tuple([
+                actx.freeze(resample_mesh_nodes(grp, iaxis)) for grp in self.groups
+                ]))
+            for iaxis in range(self.ambient_dim)])
+
+    def nodes(self):
+        from warnings import warn
+        warn("`Discretization.nodes` will be dropped in version 2022.x. "
+             "use `Discretization.coordinates` instead.",
+             DeprecationWarning, stacklevel=2)
+        return self.coordinates()
+
+# }}}
 
 # vim: fdm=marker
diff --git a/meshmode/discretization/connection/direct.py b/meshmode/discretization/connection/direct.py
index f94064ec..613a95ad 100644
--- a/meshmode/discretization/connection/direct.py
+++ b/meshmode/discretization/connection/direct.py
@@ -110,8 +110,8 @@ class DiscretizationConnectionElementGroup:
 # {{{ connection classes
 
 class DiscretizationConnection:
-    """Abstract interface for transporting a DOF vector from one subclass
-    of :class:`meshmode.discretization.DiscretizationBase` to another.
+    """Abstract interface for transporting a DOF vector from one
+    :class:`meshmode.discretization.Discretization` to another.
     Possible applications include:
 
     *   upsampling/downsampling on the same mesh
diff --git a/meshmode/discretization/connection/modal.py b/meshmode/discretization/connection/modal.py
index 3994bf88..e1f5bd47 100644
--- a/meshmode/discretization/connection/modal.py
+++ b/meshmode/discretization/connection/modal.py
@@ -29,11 +29,12 @@ import modepy as mp
 
 from meshmode.array_context import make_loopy_program
 from meshmode.dof_array import DOFArray
-from meshmode.discretization import InterpolatoryElementGroupBase
+from meshmode.discretization import (
+    NodalElementGroupBase,
+    ModalElementGroupBase,
+    InterpolatoryElementGroupBase)
 from meshmode.discretization.poly_element import QuadratureSimplexElementGroup
 from meshmode.discretization.connection.direct import DiscretizationConnection
-from meshmode.discretization.modal import ModalDiscretization
-from meshmode.discretization.nodal import NodalDiscretization
 
 from pytools import memoize_in, keyed_memoize_in
 from pytools.obj_array import obj_array_vectorized_n_args
@@ -51,13 +52,13 @@ class NodalToModalDiscretizationConnection(DiscretizationConnection):
 
     .. attribute:: from_discr
 
-        an instance of
-            :class:`meshmode.discretization.nodal.NodalDiscretization`
+        an instance of :class:`meshmode.discretization.Discretization` containing
+        :class:`meshmode.discretization.NodalElementGroupBase` element groups
 
     .. attribute:: to_discr
 
-        an instance of
-            :class:`meshmode.discretization.modal.ModalDiscretization`
+        an instance of :class:`meshmode.discretization.Discretization` containing
+        :class:`meshmode.discretization.ModalElementGroupBase` element groups
 
     .. attribute:: groups
 
@@ -71,24 +72,23 @@ class NodalToModalDiscretizationConnection(DiscretizationConnection):
 
     def __init__(self, from_discr, to_discr, allow_approximate_quad=False):
         """
-        :arg from_discr: an instance of
-            :class:`meshmode.discretization.nodal.NodalDiscretization`
-        :arg to_discr: an instance of
-            :class:`meshmode.discretization.modal.ModalDiscretization`
-        :arg allow_approximate_quad: an optional :class:`bool` flag indicating
-            whether to proceed with numerically approximating (via quadrature)
-            modal coefficients, despite using an insufficient (not high enough
-            approximation order for exact integration) quadrature method.
-            The default value is *False*.
+        :arg from_discr: a :class:`meshmode.discretization.Discretization`
+            containing :class:`meshmode.discretization.NodalElementGroupBase`
+            element groups.
+        :arg to_discr: a :class:`meshmode.discretization.Discretization`
+            containing :class:`meshmode.discretization.ModalElementGroupBase`
+            element groups.
         """
 
-        if not isinstance(from_discr, NodalDiscretization):
-            raise TypeError("from_discr must be a NodalDiscretization "
-                            "object to use this connection.")
+        if not all([isinstance(grp, NodalElementGroupBase)
+                    for grp in from_discr.groups]):
+            raise ValueError("`from_discr` must be defined on nodal "
+                             "element groups to use this connection.")
 
-        if not isinstance(to_discr, ModalDiscretization):
-            raise TypeError("to_discr must be a ModalDiscretization "
-                            "object to use this connection.")
+        if not all([isinstance(grp, ModalElementGroupBase)
+                    for grp in to_discr.groups]):
+            raise ValueError("`to_discr` must be defined on modal "
+                             "element groups to use this connection.")
 
         if to_discr.mesh != from_discr.mesh:
             raise ValueError("Both `from_discr` and `to_discr` must be on "
@@ -264,13 +264,13 @@ class ModalToNodalDiscretizationConnection(DiscretizationConnection):
 
     .. attribute:: from_discr
 
-        an instance of
-            :class:`meshmode.discretization.modal.ModalDiscretization`
+        an instance of :class:`meshmode.discretization.Discretization` containing
+        :class:`meshmode.discretization.ModalElementGroupBase` element groups
 
     .. attribute:: to_discr
 
-        an instance of
-            :class:`meshmode.discretization.nodal.NodalDiscretization`
+        an instance of :class:`meshmode.discretization.Discretization` containing
+        :class:`meshmode.discretization.NodalElementGroupBase` element groups
 
     .. attribute:: groups
 
@@ -284,19 +284,23 @@ class ModalToNodalDiscretizationConnection(DiscretizationConnection):
 
     def __init__(self, from_discr, to_discr):
         """
-        :arg from_discr: an instance of
-            :class:`meshmode.discretization.modal.ModalDiscretization`
-        :arg to_discr: an instance of
-            :class:`meshmode.discretization.nodal.NodalDiscretization`
+        :arg from_discr: a :class:`meshmode.discretization.Discretization`
+            containing :class:`meshmode.discretization.ModalElementGroupBase`
+            element groups.
+        :arg to_discr: a :class:`meshmode.discretization.Discretization`
+            containing :class:`meshmode.discretization.NodalElementGroupBase`
+            element groups.
         """
 
-        if not isinstance(from_discr, ModalDiscretization):
-            raise TypeError("from_discr must be a ModalDiscretization "
-                            "object to use this connection.")
+        if not all([isinstance(grp, ModalElementGroupBase)
+                    for grp in from_discr.groups]):
+            raise ValueError("`from_discr` must be defined on modal "
+                             "element groups to use this connection.")
 
-        if not isinstance(to_discr, NodalDiscretization):
-            raise TypeError("to_discr must be a NodalDiscretization "
-                            "object to use this connection.")
+        if not all([isinstance(grp, NodalElementGroupBase)
+                    for grp in to_discr.groups]):
+            raise ValueError("`to_discr` must be defined on nodal "
+                             "element groups to use this connection.")
 
         if to_discr.mesh != from_discr.mesh:
             raise ValueError("Both `from_discr` and `to_discr` must be on "
diff --git a/meshmode/discretization/connection/refinement.py b/meshmode/discretization/connection/refinement.py
index a42d59ba..f044873f 100644
--- a/meshmode/discretization/connection/refinement.py
+++ b/meshmode/discretization/connection/refinement.py
@@ -118,8 +118,8 @@ def make_refinement_connection(actx, refiner, coarse_discr, group_factory):
         :class:`meshmode.mesh.refinement.Refiner`
 
     :arg coarse_discr: An instance of
-        :class:`meshmode.discretization.nodal.NodalDiscretization`
-        associated with the mesh given to the refiner
+        :class:`meshmode.discretization.Discretization` associated
+        with the mesh given to the refiner
 
     :arg group_factory: An instance of
         :class:`meshmode.discretization.poly_element.ElementGroupFactory`. Used
diff --git a/meshmode/discretization/modal.py b/meshmode/discretization/modal.py
deleted file mode 100644
index 3cc27408..00000000
--- a/meshmode/discretization/modal.py
+++ /dev/null
@@ -1,61 +0,0 @@
-__copyright__ = """
-Copyright (C) 2021 University of Illinois Board of Trustees
-"""
-
-__license__ = """
-Permission is hereby granted, free of charge, to any person obtaining a copy
-of this software and associated documentation files (the "Software"), to deal
-in the Software without restriction, including without limitation the rights
-to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
-copies of the Software, and to permit persons to whom the Software is
-furnished to do so, subject to the following conditions:
-
-The above copyright notice and this permission notice shall be included in
-all copies or substantial portions of the Software.
-
-THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
-THE SOFTWARE.
-"""
-
-
-import numpy as np
-
-from meshmode.array_context import ArrayContext
-from meshmode.discretization import DiscretizationBase
-from meshmode.discretization.poly_element import (
-        ModalSimplexGroupFactory,
-        ModalTensorProductGroupFactory)
-
-__doc__ = """
-.. autoclass:: ModalDiscretization
-"""
-
-
-class ModalDiscretization(DiscretizationBase):
-    """An unstructured composite modal discretization."""
-
-    def __init__(self, actx: ArrayContext, mesh, group_factory,
-            real_dtype=np.float64):
-        """
-        :param actx: A :class:`ArrayContext` used to perform computation needed
-            during initial set-up of the mesh.
-        :param mesh: A :class:`meshmode.mesh.Mesh` over which the discretization is
-            built.
-        :param group_factory: An :class:`ElementGroupFactory`. Note that the
-            element groups must be subclasses of :class:`ModalElementGroupBase`.
-        :param real_dtype: The :mod:`numpy` data type used for representing real
-            data, either :class:`numpy.float32` or :class:`numpy.float64`.
-        """
-
-        # Some sanity checking so we don't do something horrendous later on
-        if not isinstance(group_factory, (ModalSimplexGroupFactory,
-                                          ModalTensorProductGroupFactory)):
-            raise ValueError("You must use a modal group factory when using "
-                             "a ModalDiscretization.")
-
-        super().__init__(actx, mesh, group_factory, real_dtype)
diff --git a/meshmode/discretization/nodal.py b/meshmode/discretization/nodal.py
deleted file mode 100644
index c519127f..00000000
--- a/meshmode/discretization/nodal.py
+++ /dev/null
@@ -1,185 +0,0 @@
-__copyright__ = """
-Copyright (C) 2013-2021 Andreas Kloeckner
-Copyright (C) 2021 University of Illinois Board of Trustees
-"""
-
-__license__ = """
-Permission is hereby granted, free of charge, to any person obtaining a copy
-of this software and associated documentation files (the "Software"), to deal
-in the Software without restriction, including without limitation the rights
-to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
-copies of the Software, and to permit persons to whom the Software is
-furnished to do so, subject to the following conditions:
-
-The above copyright notice and this permission notice shall be included in
-all copies or substantial portions of the Software.
-
-THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
-THE SOFTWARE.
-"""
-
-
-import numpy as np
-
-from pytools import memoize_in, memoize_method
-from pytools.obj_array import make_obj_array
-from meshmode.array_context import ArrayContext, make_loopy_program
-
-# underscored because it shouldn't be imported from here.
-from meshmode.dof_array import DOFArray as _DOFArray
-
-from meshmode.discretization import DiscretizationBase
-from meshmode.discretization.poly_element import (
-        ModalSimplexGroupFactory,
-        ModalTensorProductGroupFactory)
-
-__doc__ = """
-.. autoclass:: NodalDiscretization
-"""
-
-
-class NodalDiscretization(DiscretizationBase):
-    """An unstructured composite nodal discretization.
-
-    .. attribute:: ndofs
-
-    .. automethod:: nodes
-    .. automethod:: num_reference_derivative
-    .. automethod:: quad_weights
-    """
-
-    def __init__(self, actx: ArrayContext, mesh, group_factory,
-            real_dtype=np.float64):
-        """
-        :param actx: A :class:`ArrayContext` used to perform computation needed
-            during initial set-up of the mesh.
-        :param mesh: A :class:`meshmode.mesh.Mesh` over which the discretization is
-            built.
-        :param group_factory: An :class:`ElementGroupFactory`. Note that the
-            element groups must be subclasses of :class:`NodalElementGroupBase`.
-        :param real_dtype: The :mod:`numpy` data type used for representing real
-            data, either :class:`numpy.float32` or :class:`numpy.float64`.
-        """
-
-        # Some sanity checking so we don't do something horrendous later on
-        if isinstance(group_factory, (ModalSimplexGroupFactory,
-                                      ModalTensorProductGroupFactory)):
-            raise ValueError("Cannot use a modal element group factory with "
-                             "a NodalDiscretization.")
-
-        super().__init__(actx, mesh, group_factory, real_dtype)
-
-    def num_reference_derivative(self, ref_axes, vec):
-        actx = vec.array_context
-        ref_axes = list(ref_axes)
-
-        @memoize_in(actx, (NodalDiscretization, "reference_derivative_prg"))
-        def prg():
-            return make_loopy_program(
-                """{[iel,idof,j]:
-                    0<=iel<nelements and
-                    0<=idof,j<nunit_dofs}""",
-                "result[iel,idof] = sum(j, diff_mat[idof, j] * vec[iel, j])",
-                name="diff")
-
-        def get_mat(grp):
-            mat = None
-            for ref_axis in ref_axes:
-                next_mat = grp.diff_matrices()[ref_axis]
-                if mat is None:
-                    mat = next_mat
-                else:
-                    mat = np.dot(next_mat, mat)
-
-            return mat
-
-        return _DOFArray(actx, tuple(
-                actx.call_loopy(
-                    prg(), diff_mat=actx.from_numpy(get_mat(grp)), vec=vec[grp.index]
-                    )["result"]
-                for grp in self.groups))
-
-    @memoize_method
-    def quad_weights(self):
-        """:returns: A :class:`~meshmode.dof_array.DOFArray` with quadrature weights.
-        """
-        actx = self._setup_actx
-
-        @memoize_in(actx, (NodalDiscretization, "quad_weights_prg"))
-        def prg():
-            return make_loopy_program(
-                "{[iel,idof]: 0<=iel<nelements and 0<=idof<nunit_dofs}",
-                "result[iel,idof] = weights[idof]",
-                name="quad_weights")
-
-        return _DOFArray(None, tuple(
-                actx.freeze(
-                    actx.call_loopy(
-                        prg(),
-                        weights=actx.from_numpy(grp.weights),
-                        nelements=grp.nelements,
-                        )["result"])
-                for grp in self.groups))
-
-    @memoize_method
-    def nodes(self):
-        r"""
-        :returns: object array of shape ``(ambient_dim,)`` containing
-            :class:`~meshmode.dof_array.DOFArray`\ s of node coordinates.
-        """
-
-        actx = self._setup_actx
-
-        @memoize_in(actx, (NodalDiscretization, "nodes_prg"))
-        def prg():
-            return make_loopy_program(
-                """{[iel,idof,j]:
-                    0<=iel<nelements and
-                    0<=idof<ndiscr_nodes and
-                    0<=j<nmesh_nodes}""",
-                """
-                    result[iel, idof] = \
-                        sum(j, resampling_mat[idof, j] * nodes[iel, j])
-                    """,
-                name="nodes")
-
-        def resample_mesh_nodes(grp, iaxis):
-            # TODO: would be nice to have the mesh use an array context already
-            nodes = actx.from_numpy(grp.mesh_el_group.nodes[iaxis])
-
-            grp_unit_nodes = grp.unit_nodes.reshape(-1)
-            meg_unit_nodes = grp.mesh_el_group.unit_nodes.reshape(-1)
-
-            tol = 10 * np.finfo(grp_unit_nodes.dtype).eps
-            if (grp_unit_nodes.shape == meg_unit_nodes.shape
-                    and np.linalg.norm(grp_unit_nodes - meg_unit_nodes) < tol):
-                return nodes
-
-            return actx.call_loopy(
-                    prg(),
-                    resampling_mat=actx.from_numpy(grp.from_mesh_interp_matrix()),
-                    nodes=nodes,
-                    )["result"]
-
-        return make_obj_array([
-            _DOFArray(None, tuple([
-                actx.freeze(resample_mesh_nodes(grp, iaxis)) for grp in self.groups
-                ]))
-            for iaxis in range(self.ambient_dim)])
-
-
-class Discretization(NodalDiscretization):
-
-    def __init__(self, *args, **kwargs):
-
-        from warnings import warn
-        warn("`Discretization` will be dropped in favor of "
-             "`NodalDiscretization` in version 2022.x.",
-             DeprecationWarning, stacklevel=2)
-
-        super().__init__(*args, **kwargs)
diff --git a/meshmode/dof_array.py b/meshmode/dof_array.py
index 3dacc378..d865cec6 100644
--- a/meshmode/dof_array.py
+++ b/meshmode/dof_array.py
@@ -57,8 +57,8 @@ __doc__ = """
 # {{{ DOFArray
 
 class DOFArray:
-    r"""This array type holds degree-of-freedom arrays for use with subclasses of
-    :class:`~meshmode.discretization.DiscretizationBase`,
+    r"""This array type holds degree-of-freedom arrays for use with
+    :class:`~meshmode.discretization.Discretization`,
     with one entry in the :class:`DOFArray` for each
     :class:`~meshmode.discretization.ElementGroupBase`.
     The arrays contained within a :class:`DOFArray`
@@ -67,7 +67,7 @@ class DOFArray:
     :attr:`~meshmode.discretization.ElementGroupBase.nelements`
     of the associated group.
     ``ndofs_per_element`` is typically, but not necessarily, the same as
-    :attr:`~meshmode.discretization.NodalElementGroupBase.nunit_dofs`
+    :attr:`~meshmode.discretization.ElementGroupBase.nunit_dofs`
     of the associated group. The entries in this array are further arrays managed by
     :attr:`array_context`.
 
diff --git a/meshmode/interop/firedrake/connection.py b/meshmode/interop/firedrake/connection.py
index 2836f844..bbdf8369 100644
--- a/meshmode/interop/firedrake/connection.py
+++ b/meshmode/interop/firedrake/connection.py
@@ -42,8 +42,8 @@ from meshmode.discretization.poly_element import (
     PolynomialWarpAndBlendGroupFactory,
     PolynomialRecursiveNodesGroupFactory,
     ElementGroupFactory)
-from meshmode.discretization import InterpolatoryElementGroupBase
-from meshmode.discretization.nodal import NodalDiscretization
+from meshmode.discretization import (
+    Discretization, InterpolatoryElementGroupBase)
 
 from pytools import memoize_method
 
@@ -96,7 +96,7 @@ class FiredrakeConnection:
 
     .. attribute:: discr
 
-        A :class:`meshmode.discretization.nodal.NodalDiscretization`.
+        A :class:`meshmode.discretization.Discretization`.
 
     .. attribute:: group_nr
 
@@ -135,7 +135,7 @@ class FiredrakeConnection:
     """
     def __init__(self, discr, fdrake_fspace, mm2fd_node_mapping, group_nr=None):
         """
-        :param discr: A :class:`meshmode.discretization.nodal.NodalDiscretization`
+        :param discr: A :class:`meshmode.discretization.Discretization`
         :param fdrake_fspace: A
             :class:`firedrake.functionspaceimpl.WithGeometry`.
             Must use ufl family ``"Discontinuous Lagrange"``.
@@ -157,9 +157,9 @@ class FiredrakeConnection:
             if *group_nr* is *None* when *discr* has more than one group.
         """
         # {{{ Validate input
-        if not isinstance(discr, NodalDiscretization):
+        if not isinstance(discr, Discretization):
             raise TypeError("'discr' must be of type "
-                            "meshmode.discretization.nodal.NodalDiscretization, "
+                            "meshmode.discretization.Discretization, "
                             "not '%s'`." % type(discr))
         from firedrake.functionspaceimpl import WithGeometry
         if not isinstance(fdrake_fspace, WithGeometry):
@@ -639,7 +639,7 @@ PolynomialWarpAndBlendGroupFactory` is used.
           to :class:`~meshmode.mesh.BTAG_ALL`).
 
         If not *None*, creates a
-        :class:`~meshmode.discretization.nodal.NodalDiscretization` on a submesh
+        :class:`~meshmode.discretization.Discretization` on a submesh
         of ``fdrake_fspace.mesh()`` created from the cells with at least
         one vertex on a facet marked with the marker
         *restrict_to_boundary*.
@@ -726,7 +726,7 @@ PolynomialWarpAndBlendGroupFactory` is used.
     # Create to_discr
     mm_mesh, orient = import_firedrake_mesh(fdrake_fspace.mesh(),
                                             cells_to_use=cells_to_use)
-    to_discr = NodalDiscretization(actx, mm_mesh, grp_factory)
+    to_discr = Discretization(actx, mm_mesh, grp_factory)
 
     # get firedrake unit nodes and map onto meshmode reference element
     group = to_discr.groups[0]
@@ -771,7 +771,7 @@ def build_connection_to_firedrake(discr, group_nr=None, comm=None):
     space and allow for conversion back and forth
     by resampling at the nodes.
 
-    :param discr: A :class:`~meshmode.discretization.nodal.NodalDiscretization`
+    :param discr: A :class:`~meshmode.discretization.Discretization`
         to intialize the connection with
     :param group_nr: The group number of the discretization to convert.
         If *None* there must be only one group. The selected group
diff --git a/meshmode/interop/nodal_dg.py b/meshmode/interop/nodal_dg.py
index 6d739b15..0e438d9d 100644
--- a/meshmode/interop/nodal_dg.py
+++ b/meshmode/interop/nodal_dg.py
@@ -102,7 +102,7 @@ class NodalDGContext:
 
         self.octave.eval(f"StartUp{mesh.dim}D;")
 
-    def get_discr(self, actx) -> meshmode.discretization.nodal.NodalDiscretization:
+    def get_discr(self, actx) -> meshmode.discretization.Discretization:
         """Get a discretization with nodes exactly matching the ones used
         by the nodal-DG code.
 
@@ -144,10 +144,10 @@ class NodalDGContext:
 
         mesh = Mesh(vertices=vertices, groups=[egroup], is_conforming=True)
 
-        from meshmode.discretization.nodal import NodalDiscretization
+        from meshmode.discretization import Discretization
         from meshmode.discretization.poly_element import (
                 PolynomialGivenNodesGroupFactory)
-        return NodalDiscretization(actx, mesh,
+        return Discretization(actx, mesh,
                 PolynomialGivenNodesGroupFactory(order, unit_nodes))
 
     def push_dof_array(self, name, ary: meshmode.dof_array.DOFArray):
diff --git a/meshmode/mesh/visualization.py b/meshmode/mesh/visualization.py
index 20bb70b5..a6048626 100644
--- a/meshmode/mesh/visualization.py
+++ b/meshmode/mesh/visualization.py
@@ -309,8 +309,8 @@ def vtk_visualize_mesh(actx, mesh, filename, vtk_high_order=True):
 
     from meshmode.discretization.poly_element import \
             PolynomialWarpAndBlendGroupFactory
-    from meshmode.discretization.nodal import NodalDiscretization as NDiscr
-    discr = NDiscr(actx, mesh, PolynomialWarpAndBlendGroupFactory(order))
+    from meshmode.discretization import Discretization
+    discr = Discretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(order))
 
     from meshmode.discretization.visualization import make_visualizer
     vis = make_visualizer(actx, discr, order, force_equidistant=vtk_high_order)
diff --git a/test/test_array.py b/test/test_array.py
index b6c157d1..2f089a05 100644
--- a/test/test_array.py
+++ b/test/test_array.py
@@ -27,7 +27,7 @@ from meshmode.array_context import (  # noqa
         pytest_generate_tests_for_pyopencl_array_context
         as pytest_generate_tests)
 
-from meshmode.discretization.nodal import NodalDiscretization
+from meshmode.discretization import Discretization
 from meshmode.discretization.poly_element import PolynomialWarpAndBlendGroupFactory
 from meshmode.dof_array import flatten, unflatten, DOFArray
 
@@ -44,7 +44,7 @@ def test_array_context_np_workalike(actx_factory):
     mesh = generate_regular_rect_mesh(
             a=(-0.5,)*2, b=(0.5,)*2, n=(8,)*2, order=3)
 
-    discr = NodalDiscretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(3))
+    discr = Discretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(3))
 
     for sym_name, n_args in [
             ("sin", 1),
@@ -88,7 +88,7 @@ def test_dof_array_arithmetic_same_as_numpy(actx_factory):
     mesh = generate_regular_rect_mesh(
             a=(-0.5,)*2, b=(0.5,)*2, n=(3,)*2, order=1)
 
-    discr = NodalDiscretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(3))
+    discr = Discretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(3))
 
     def get_real(ary):
         return ary.real
diff --git a/test/test_chained.py b/test/test_chained.py
index dea34243..70296c6f 100644
--- a/test/test_chained.py
+++ b/test/test_chained.py
@@ -70,10 +70,10 @@ def create_discretization(actx, ndim,
         raise ValueError(f"unsupported dimension: {ndim}")
 
     # create discretization
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.poly_element import \
             InterpolatoryQuadratureSimplexGroupFactory
-    discr = NodalDiscretization(actx, mesh,
+    discr = Discretization(actx, mesh,
             InterpolatoryQuadratureSimplexGroupFactory(order))
 
     return discr
diff --git a/test/test_firedrake_interop.py b/test/test_firedrake_interop.py
index c4313783..045978c8 100644
--- a/test/test_firedrake_interop.py
+++ b/test/test_firedrake_interop.py
@@ -29,7 +29,7 @@ from pyopencl.tools import (  # noqa
 
 from meshmode.array_context import PyOpenCLArrayContext
 
-from meshmode.discretization.nodal import NodalDiscretization
+from meshmode.discretization import Discretization
 from meshmode.discretization.poly_element import (
     InterpolatoryQuadratureSimplexGroupFactory)
 
@@ -192,7 +192,7 @@ def test_to_fd_consistency(ctx_factory, mm_mesh, fspace_degree):
     actx = PyOpenCLArrayContext(queue)
 
     factory = InterpolatoryQuadratureSimplexGroupFactory(fspace_degree)
-    discr = NodalDiscretization(actx, mm_mesh, factory)
+    discr = Discretization(actx, mm_mesh, factory)
     fdrake_connection = build_connection_to_firedrake(discr)
     fdrake_fspace = fdrake_connection.firedrake_fspace()
     # Check consistency
@@ -538,7 +538,7 @@ def test_to_fd_transfer(ctx_factory, fspace_degree, mesh_name, mesh_pars, dim):
 
         # Make discr and connect it to firedrake
         factory = InterpolatoryQuadratureSimplexGroupFactory(fspace_degree)
-        discr = NodalDiscretization(actx, mm_mesh, factory)
+        discr = Discretization(actx, mm_mesh, factory)
 
         fdrake_connection = build_connection_to_firedrake(discr)
         fdrake_fspace = fdrake_connection.firedrake_fspace()
@@ -658,7 +658,7 @@ def test_to_fd_idempotency(ctx_factory, mm_mesh, fspace_degree):
 
     # Make a function space and a function with unique values at each node
     factory = InterpolatoryQuadratureSimplexGroupFactory(fspace_degree)
-    discr = NodalDiscretization(actx, mm_mesh, factory)
+    discr = Discretization(actx, mm_mesh, factory)
     fdrake_connection = build_connection_to_firedrake(discr)
     fdrake_mesh = fdrake_connection.firedrake_fspace().mesh()
     dtype = fdrake_mesh.coordinates.dat.data.dtype
diff --git a/test/test_mesh.py b/test/test_mesh.py
index 28c9b026..1bf76851 100644
--- a/test/test_mesh.py
+++ b/test/test_mesh.py
@@ -42,10 +42,10 @@ def test_nonequal_rect_mesh_generation(actx_factory, visualize=False):
     mesh = mgen.generate_regular_rect_mesh(
             a=(0, 0)*2, b=(5, 3), n=(10, 6,), order=3)
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.poly_element import \
             PolynomialWarpAndBlendGroupFactory as GroupFactory
-    discr = NodalDiscretization(actx, mesh, GroupFactory(3))
+    discr = Discretization(actx, mesh, GroupFactory(3))
 
     if visualize:
         from meshmode.discretization.visualization import make_visualizer
diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index ca7cc0e1..1a50edfc 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -112,8 +112,8 @@ def test_parallel_vtk_file(actx_factory, dim):
     else:
         raise ValueError("unknown dimensionality")
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    discr = NodalDiscretization(actx, mesh,
+    from meshmode.discretization import Discretization
+    discr = Discretization(actx, mesh,
             InterpolatoryQuadratureSimplexGroupFactory(target_order))
 
     from meshmode.discretization.visualization import make_visualizer
@@ -193,8 +193,8 @@ def test_visualizers(actx_factory, dim, group_cls):
     else:
         group_factory = LegendreGaussLobattoTensorProductGroupFactory
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    discr = NodalDiscretization(actx, mesh, group_factory(target_order))
+    from meshmode.discretization import Discretization
+    discr = Discretization(actx, mesh, group_factory(target_order))
 
     nodes = thaw(actx, discr.nodes())
     f = actx.np.sqrt(sum(nodes**2)) + 1j*nodes[0]
@@ -419,7 +419,7 @@ def test_boundary_interpolation(actx_factory, group_factory, boundary_tag,
     else:
         group_cls = SimplexElementGroup
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.connection import (
             make_face_restriction, check_connection)
 
@@ -468,7 +468,7 @@ def test_boundary_interpolation(actx_factory, group_factory, boundary_tag,
 
         # }}}
 
-        vol_discr = NodalDiscretization(actx, mesh, group_factory(order))
+        vol_discr = Discretization(actx, mesh, group_factory(order))
         print("h=%s -> %d elements" % (
                 h, sum(mgrp.nelements for mgrp in mesh.groups)))
 
@@ -533,7 +533,7 @@ def test_all_faces_interpolation(actx_factory, group_factory,
     else:
         group_cls = SimplexElementGroup
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.connection import (
             make_face_restriction, make_face_to_all_faces_embedding,
             check_connection)
@@ -574,7 +574,7 @@ def test_all_faces_interpolation(actx_factory, group_factory,
 
         # }}}
 
-        vol_discr = NodalDiscretization(actx, mesh, group_factory(order))
+        vol_discr = Discretization(actx, mesh, group_factory(order))
         print("h=%s -> %d elements" % (
                 h, sum(mgrp.nelements for mgrp in mesh.groups)))
 
@@ -656,7 +656,7 @@ def test_opposite_face_interpolation(actx_factory, group_factory,
     else:
         group_cls = SimplexElementGroup
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.connection import (
             make_face_restriction, make_opposite_face_connection,
             check_connection)
@@ -705,7 +705,7 @@ def test_opposite_face_interpolation(actx_factory, group_factory,
 
         # }}}
 
-        vol_discr = NodalDiscretization(actx, mesh, group_factory(order))
+        vol_discr = Discretization(actx, mesh, group_factory(order))
         print("h=%s -> %d elements" % (
                 h, sum(mgrp.nelements for mgrp in mesh.groups)))
 
@@ -865,9 +865,9 @@ def test_merge_and_map(actx_factory, group_cls, visualize=False):
     mesh4 = mesh3.copy()
 
     if visualize:
-        from meshmode.discretization.nodal import NodalDiscretization
+        from meshmode.discretization import Discretization
         actx = actx_factory()
-        discr = NodalDiscretization(actx, mesh4, discr_grp_factory)
+        discr = Discretization(actx, mesh4, discr_grp_factory)
 
         from meshmode.discretization.visualization import make_visualizer
         vis = make_visualizer(actx, discr, 3, element_shrink_factor=0.8)
@@ -985,10 +985,10 @@ def test_sanity_qhull_nd(actx_factory, dim, order):
     mesh = from_vertices_and_simplices(dtri.points.T, dtri.simplices,
             fix_orientation=True)
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    low_discr = NodalDiscretization(actx, mesh,
+    from meshmode.discretization import Discretization
+    low_discr = Discretization(actx, mesh,
             PolynomialEquidistantSimplexGroupFactory(order))
-    high_discr = NodalDiscretization(actx, mesh,
+    high_discr = Discretization(actx, mesh,
             PolynomialEquidistantSimplexGroupFactory(order+1))
 
     from meshmode.discretization.connection import make_same_mesh_connection
@@ -1151,10 +1151,10 @@ def test_box_mesh(actx_factory, visualize=False):
     mesh = mgen.generate_box_mesh(3*(np.linspace(0, 1, 5),))
 
     if visualize:
-        from meshmode.discretization.nodal import NodalDiscretization
+        from meshmode.discretization import Discretization
 
         actx = actx_factory()
-        discr = NodalDiscretization(actx, mesh,
+        discr = Discretization(actx, mesh,
                 PolynomialWarpAndBlendGroupFactory(7))
 
         from meshmode.discretization.visualization import make_visualizer
@@ -1404,8 +1404,8 @@ def test_vtk_overwrite(actx_factory):
 
     mesh = mgen.generate_torus(10.0, 2.0, order=target_order)
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    discr = NodalDiscretization(
+    from meshmode.discretization import Discretization
+    discr = Discretization(
             actx, mesh,
             InterpolatoryQuadratureSimplexGroupFactory(target_order))
 
@@ -1478,8 +1478,8 @@ def test_mesh_without_vertices(actx_factory):
     mesh = refine_uniformly(mesh, 1)
 
     # make sure the world doesn't end
-    from meshmode.discretization.nodal import NodalDiscretization
-    discr = NodalDiscretization(actx, mesh,
+    from meshmode.discretization import Discretization
+    discr = Discretization(actx, mesh,
             InterpolatoryQuadratureSimplexGroupFactory(4))
     thaw(actx, discr.nodes())
 
@@ -1605,9 +1605,8 @@ def test_mesh_multiple_groups(actx_factory, ambient_dim, visualize=False):
         import matplotlib.pyplot as plt
         plt.savefig("test_mesh_multiple_groups_2d_elements.png", dpi=300)
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    discr = NodalDiscretization(actx, mesh,
-            PolynomialWarpAndBlendGroupFactory(order))
+    from meshmode.discretization import Discretization
+    discr = Discretization(actx, mesh, PolynomialWarpAndBlendGroupFactory(order))
 
     if visualize:
         group_id = discr.empty(actx, dtype=np.int32)
@@ -1697,9 +1696,9 @@ def test_mesh_with_interior_unit_nodes(actx_factory, ambient_dim):
     assert mesh.facial_adjacency_groups
     assert mesh.nodal_adjacency
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.poly_element import QuadratureSimplexGroupFactory
-    discr = NodalDiscretization(actx, mesh,
+    discr = Discretization(actx, mesh,
             QuadratureSimplexGroupFactory(order))
 
     from meshmode.discretization.connection import make_face_restriction
diff --git a/test/test_modal.py b/test/test_modal.py
index 8def0ea5..b6ed033b 100644
--- a/test/test_modal.py
+++ b/test/test_modal.py
@@ -49,8 +49,7 @@ from meshmode.discretization.poly_element import (
     QuadratureSimplexGroupFactory
     )
 
-from meshmode.discretization.nodal import NodalDiscretization
-from meshmode.discretization.modal import ModalDiscretization
+from meshmode.discretization import Discretization
 from meshmode.discretization.connection.modal import (
     NodalToModalDiscretizationConnection,
     ModalToNodalDiscretizationConnection
@@ -89,8 +88,8 @@ def test_inverse_modal_connections(actx_factory, nodal_group_factory):
         a=(0, 0)*2, b=(5, 3), n=(10, 6,), order=order, group_cls=group_cls)
 
     # Make discretizations
-    nodal_disc = NodalDiscretization(actx, mesh, nodal_group_factory(order))
-    modal_disc = ModalDiscretization(actx, mesh, modal_group_factory(order))
+    nodal_disc = Discretization(actx, mesh, nodal_group_factory(order))
+    modal_disc = Discretization(actx, mesh, modal_group_factory(order))
 
     # Make connections
     nodal_to_modal_conn = NodalToModalDiscretizationConnection(
@@ -130,8 +129,8 @@ def test_modal_coefficients_by_projection(actx_factory, quad_group_factory):
         a=(0, 0)*2, b=(5, 3), n=(10, 6,), order=order, group_cls=group_cls)
 
     # Make discretizations
-    nodal_disc = NodalDiscretization(actx, mesh, quad_group_factory(order))
-    modal_disc = ModalDiscretization(actx, mesh, modal_group_factory(m_order))
+    nodal_disc = Discretization(actx, mesh, quad_group_factory(order))
+    modal_disc = Discretization(actx, mesh, modal_group_factory(m_order))
 
     # Make connections one using quadrature projection
     nodal_to_modal_conn_quad = NodalToModalDiscretizationConnection(
@@ -190,8 +189,8 @@ def test_quadrature_based_modal_connection_reverse(actx_factory, quad_group_fact
         a=(0, 0)*2, b=(5, 3), n=(10, 6,), order=order, group_cls=group_cls)
 
     # Make discretizations
-    nodal_disc = NodalDiscretization(actx, mesh, quad_group_factory(order))
-    modal_disc = ModalDiscretization(actx, mesh, modal_group_factory(m_order))
+    nodal_disc = Discretization(actx, mesh, quad_group_factory(order))
+    modal_disc = Discretization(actx, mesh, modal_group_factory(m_order))
 
     # Make connections one using quadrature projection
     nodal_to_modal_conn_quad = NodalToModalDiscretizationConnection(
@@ -257,8 +256,8 @@ def test_modal_truncation(actx_factory, nodal_group_factory,
         h = 1/mesh_par
 
         # Make discretizations
-        nodal_disc = NodalDiscretization(actx, mesh, nodal_group_factory(order))
-        modal_disc = ModalDiscretization(actx, mesh, modal_group_factory(order))
+        nodal_disc = Discretization(actx, mesh, nodal_group_factory(order))
+        modal_disc = Discretization(actx, mesh, modal_group_factory(order))
 
         # Make connections (nodal -> modal)
         nodal_to_modal_conn = NodalToModalDiscretizationConnection(
diff --git a/test/test_partition.py b/test/test_partition.py
index 4232349e..77db5abe 100644
--- a/test/test_partition.py
+++ b/test/test_partition.py
@@ -103,8 +103,8 @@ def test_partition_interpolation(actx_factory, dim, mesh_pars,
             for i_remote_part in neighbors:
                 connected_parts.add((i_local_part, i_remote_part))
 
-        from meshmode.discretization.nodal import NodalDiscretization
-        vol_discrs = [NodalDiscretization(actx, part_meshes[i], group_factory(order))
+        from meshmode.discretization import Discretization
+        vol_discrs = [Discretization(actx, part_meshes[i], group_factory(order))
                         for i in range(num_parts)]
 
         from meshmode.mesh import BTAG_PARTITION
@@ -348,8 +348,8 @@ def _test_mpi_boundary_swap(dim, order, num_groups):
     queue = cl.CommandQueue(cl_ctx)
     actx = PyOpenCLArrayContext(queue)
 
-    from meshmode.discretization.nodal import NodalDiscretization
-    vol_discr = NodalDiscretization(actx, local_mesh, group_factory)
+    from meshmode.discretization import Discretization
+    vol_discr = Discretization(actx, local_mesh, group_factory)
 
     from meshmode.distributed import get_connected_partitions
     connected_parts = get_connected_partitions(local_mesh)
diff --git a/test/test_refinement.py b/test/test_refinement.py
index 42ebf9b9..ca748e4f 100644
--- a/test/test_refinement.py
+++ b/test/test_refinement.py
@@ -202,7 +202,7 @@ def test_refinement_connection(
     # discretization order
     order = 5
 
-    from meshmode.discretization.nodal import NodalDiscretization
+    from meshmode.discretization import Discretization
     from meshmode.discretization.connection import (
             make_refinement_connection, check_connection)
 
@@ -249,7 +249,7 @@ def test_refinement_connection(
 
             return result
 
-        discr = NodalDiscretization(actx, mesh, group_factory(order))
+        discr = Discretization(actx, mesh, group_factory(order))
 
         refiner = refiner_cls(mesh)
         flags = refine_flags(mesh)
-- 
GitLab