diff --git a/examples/advection/surface.py b/examples/advection/surface.py
index df179a7615f1e643f8e25bdde656cfed4edfc3e7..d528f5f3eee29366d9c963c3dbd5fe410e1819fa 100644
--- a/examples/advection/surface.py
+++ b/examples/advection/surface.py
@@ -136,24 +136,28 @@ def main(ctx_factory, dim=2, order=4, product_tag=None, visualize=False):
else:
raise ValueError("unsupported dimension")
- quad_tag_to_group_factory = {}
+ discr_tag_to_group_factory = {}
if product_tag == "none":
product_tag = None
+ else:
+ product_tag = dof_desc.DISCR_TAG_QUAD
from meshmode.discretization.poly_element import \
PolynomialWarpAndBlendGroupFactory, \
QuadratureSimplexGroupFactory
- quad_tag_to_group_factory[dof_desc.QTAG_NONE] = \
- PolynomialWarpAndBlendGroupFactory(order)
+ discr_tag_to_group_factory[dof_desc.DISCR_TAG_BASE] = \
+ PolynomialWarpAndBlendGroupFactory(order)
if product_tag:
- quad_tag_to_group_factory[product_tag] = \
- QuadratureSimplexGroupFactory(order=4*order)
+ discr_tag_to_group_factory[product_tag] = \
+ QuadratureSimplexGroupFactory(order=4*order)
from grudge import DiscretizationCollection
- discr = DiscretizationCollection(actx, mesh,
- quad_tag_to_group_factory=quad_tag_to_group_factory)
+ discr = DiscretizationCollection(
+ actx, mesh,
+ discr_tag_to_group_factory=discr_tag_to_group_factory
+ )
volume_discr = discr.discr_from_dd(dof_desc.DD_VOLUME)
logger.info("ndofs: %d", volume_discr.ndofs)
diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py
index b99e21ad5056ed466524f9c9eafeb4f5ab0c756e..5843b99318cb65de77e53e9a78669a2c001dee2c 100644
--- a/examples/advection/var-velocity.py
+++ b/examples/advection/var-velocity.py
@@ -141,15 +141,17 @@ def main(ctx_factory, dim=2, order=4, product_tag=None, visualize=False):
QuadratureSimplexGroupFactory
if product_tag:
- quad_tag_to_group_factory = {
- product_tag: QuadratureSimplexGroupFactory(order=4*order)
- }
+ discr_tag_to_group_factory = {
+ product_tag: QuadratureSimplexGroupFactory(order=4*order)
+ }
else:
- quad_tag_to_group_factory = {}
+ discr_tag_to_group_factory = {}
from grudge import DiscretizationCollection
- discr = DiscretizationCollection(actx, mesh, order=order,
- quad_tag_to_group_factory=quad_tag_to_group_factory)
+ discr = DiscretizationCollection(
+ actx, mesh, order=order,
+ discr_tag_to_group_factory=discr_tag_to_group_factory
+ )
# }}}
diff --git a/examples/dagrt-fusion.py b/examples/dagrt-fusion.py
index cb51964f6442192ed30ea953fe1dc1be2720bdab..bd399c3a38542280f49cb66b00d6979541085a31 100755
--- a/examples/dagrt-fusion.py
+++ b/examples/dagrt-fusion.py
@@ -523,7 +523,7 @@ def test_stepper_equivalence(ctx_factory, order=4):
fused_steps = fused_stepper.run(ic, t_start, dt, t_end)
- for t_ref, (u_ref, v_ref) in stepper.run(ic, t_start, dt, t_end):
+ for t_ref, (u_ref, _v_ref) in stepper.run(ic, t_start, dt, t_end):
step += 1
logger.debug("step %d/%d", step, nsteps)
t, (u, v) = next(fused_steps)
@@ -620,7 +620,7 @@ class ExecutionMapperWithMemOpCounting(ExecutionMapperWrapper):
expr.op, expr.field, profile_data)
else:
- assert False, ("unknown operator: %s" % expr.op)
+ raise AssertionError("unknown operator: %s" % expr.op)
else:
logger.debug("assignment not profiled: %s", expr)
@@ -811,7 +811,7 @@ def test_stepper_mem_ops(ctx_factory, use_fusion):
step = 0
nsteps = int(np.ceil((t_end + 1e-9) / dt))
- for (_, _, profile_data) in stepper.run(
+ for (_, _, profile_data) in stepper.run( # noqa: B007
ic, t_start, dt, t_end, return_profile_data=True):
step += 1
logger.info("step %d/%d", step, nsteps)
@@ -984,7 +984,7 @@ def test_stepper_timing(ctx_factory, use_fusion):
import time
t = time.time()
nsteps = int(np.ceil((t_end + 1e-9) / dt))
- for (_, _, profile_data) in stepper.run(
+ for (_, _, profile_data) in stepper.run( # noqa: B007
ic, t_start, dt, t_end, return_profile_data=True):
step += 1
tn = time.time()
@@ -1146,7 +1146,7 @@ def mem_ops_results(actx, dims):
result = {}
- for (_, _, profile_data) in stepper.run(
+ for (_, _, profile_data) in stepper.run( # noqa: B007
ic, t_start, dt, t_end, return_profile_data=True):
pass
@@ -1164,7 +1164,7 @@ def mem_ops_results(actx, dims):
result["nonfused_bytes_read_by_scalar_assignments"] \
+ result["nonfused_bytes_written_by_scalar_assignments"]
- for (_, _, profile_data) in fused_stepper.run(
+ for (_, _, profile_data) in fused_stepper.run( # noqa: B007
ic, t_start, dt, t_end, return_profile_data=True):
pass
diff --git a/examples/wave/wave-op-var-velocity.py b/examples/wave/wave-op-var-velocity.py
index 799ee2ed371c3172f4a16a65dbe6bc22a1277989..5ac1233d34e59f6fc7073307c4adb2ff3bb2c494 100644
--- a/examples/wave/wave-op-var-velocity.py
+++ b/examples/wave/wave-op-var-velocity.py
@@ -33,7 +33,7 @@ from meshmode.dof_array import thaw
from meshmode.mesh import BTAG_ALL, BTAG_NONE # noqa
from grudge.discretization import DiscretizationCollection
-from grudge.dof_desc import QTAG_NONE, DOFDesc
+from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD, DOFDesc
import grudge.op as op
from grudge.shortcuts import make_visualizer
from grudge.symbolic.primitives import TracePair
@@ -43,7 +43,7 @@ from grudge.symbolic.primitives import TracePair
def wave_flux(dcoll, c, w_tpair):
dd = w_tpair.dd
- dd_quad = dd.with_qtag("vel_prod")
+ dd_quad = dd.with_discr_tag(DISCR_TAG_QUAD)
u = w_tpair[0]
v = w_tpair[1:]
@@ -78,13 +78,13 @@ def wave_operator(dcoll, c, w):
dir_bval = flat_obj_array(dir_u, dir_v)
dir_bc = flat_obj_array(-dir_u, dir_v)
- dd_quad = DOFDesc("vol", "vel_prod")
+ dd_quad = DOFDesc("vol", DISCR_TAG_QUAD)
c_quad = op.project(dcoll, "vol", dd_quad, c)
w_quad = op.project(dcoll, "vol", dd_quad, w)
u_quad = w_quad[0]
v_quad = w_quad[1:]
- dd_allfaces_quad = DOFDesc("all_faces", "vel_prod")
+ dd_allfaces_quad = DOFDesc("all_faces", DISCR_TAG_QUAD)
return (
op.inverse_mass(dcoll,
@@ -161,11 +161,13 @@ def main():
from meshmode.discretization.poly_element import \
QuadratureSimplexGroupFactory, \
PolynomialWarpAndBlendGroupFactory
- dcoll = DiscretizationCollection(actx, mesh,
- quad_tag_to_group_factory={
- QTAG_NONE: PolynomialWarpAndBlendGroupFactory(order),
- "vel_prod": QuadratureSimplexGroupFactory(3*order),
- })
+ dcoll = DiscretizationCollection(
+ actx, mesh,
+ discr_tag_to_group_factory={
+ DISCR_TAG_BASE: PolynomialWarpAndBlendGroupFactory(order),
+ DISCR_TAG_QUAD: QuadratureSimplexGroupFactory(3*order),
+ }
+ )
# bounded above by 1
c = 0.2 + 0.8*bump(actx, dcoll, center=np.zeros(3), width=0.5)
diff --git a/grudge/discretization.py b/grudge/discretization.py
index a46981c1b1f8e813f3159af5215d99f75216dd7e..2b86c994ec8fa58728f1bb0ef4df255f62549e20 100644
--- a/grudge/discretization.py
+++ b/grudge/discretization.py
@@ -21,14 +21,18 @@ THE SOFTWARE.
"""
from pytools import memoize_method
-from grudge.dof_desc import \
- QTAG_NONE, QTAG_MODAL, DTAG_BOUNDARY, DOFDesc, as_dofdesc
+from grudge.dof_desc import (
+ DISCR_TAG_BASE, DISCR_TAG_MODAL,
+ DTAG_BOUNDARY, DOFDesc, as_dofdesc
+)
import numpy as np # noqa: F401
from meshmode.array_context import ArrayContext
from meshmode.discretization.connection import \
FACE_RESTR_INTERIOR, FACE_RESTR_ALL, make_face_restriction
from meshmode.mesh import BTAG_PARTITION
+from warnings import warn
+
__doc__ = """
.. autoclass:: DiscretizationCollection
@@ -37,6 +41,8 @@ __doc__ = """
class DiscretizationCollection:
"""
+ .. automethod :: __init__
+
.. automethod :: discr_from_dd
.. automethod :: connection_from_dds
@@ -49,52 +55,76 @@ class DiscretizationCollection:
"""
def __init__(self, array_context, mesh, order=None,
- quad_tag_to_group_factory=None, mpi_communicator=None):
+ discr_tag_to_group_factory=None, mpi_communicator=None,
+ # FIXME: `quad_tag_to_group_factory` is deprecated
+ quad_tag_to_group_factory=None):
"""
- :param quad_tag_to_group_factory: A mapping from quadrature tags (typically
- strings--but may be any hashable/comparable object) to a
+ :param discr_tag_to_group_factory: A mapping from discretization tags
+ (typically one of: :class:`grudge.dof_desc.DISCR_TAG_BASE`,
+ :class:`grudge.dof_desc.DISCR_TAG_MODAL`, or
+ :class:`grudge.dof_desc.DISCR_TAG_QUAD`) to a
:class:`~meshmode.discretization.poly_element.ElementGroupFactory`
- indicating with which quadrature discretization the operations are
+ indicating with which type of discretization the operations are
to be carried out, or *None* to indicate that operations with this
- quadrature tag should be carried out with the standard volume
+ discretization tag should be carried out with the standard volume
discretization.
"""
+ if (quad_tag_to_group_factory is not None
+ and discr_tag_to_group_factory is not None):
+ raise ValueError(
+ "Both `quad_tag_to_group_factory` and `discr_tag_to_group_factory` "
+ "are specified. Use `discr_tag_to_group_factory` instead."
+ )
+
+ # FIXME: `quad_tag_to_group_factory` is deprecated
+ if (quad_tag_to_group_factory is not None
+ and discr_tag_to_group_factory is None):
+ warn("`quad_tag_to_group_factory` is a deprecated kwarg and will "
+ "be dropped in version 2022.x. Use `discr_tag_to_group_factory` "
+ "instead.",
+ DeprecationWarning, stacklevel=2)
+ discr_tag_to_group_factory = quad_tag_to_group_factory
+
self._setup_actx = array_context
from meshmode.discretization.poly_element import \
PolynomialWarpAndBlendGroupFactory
- if quad_tag_to_group_factory is None:
+ if discr_tag_to_group_factory is None:
if order is None:
- raise TypeError("one of 'order' and "
- "'quad_tag_to_group_factory' must be given")
+ raise TypeError(
+ "one of 'order' and 'discr_tag_to_group_factory' must be given"
+ )
- quad_tag_to_group_factory = {
- QTAG_NONE: PolynomialWarpAndBlendGroupFactory(order=order)}
+ discr_tag_to_group_factory = {
+ DISCR_TAG_BASE: PolynomialWarpAndBlendGroupFactory(order=order)}
else:
if order is not None:
- quad_tag_to_group_factory = quad_tag_to_group_factory.copy()
- if QTAG_NONE in quad_tag_to_group_factory:
- raise ValueError("if 'order' is given, "
- "'quad_tag_to_group_factory' must not have a "
- "key of QTAG_NONE")
-
- quad_tag_to_group_factory[QTAG_NONE] = \
+ discr_tag_to_group_factory = discr_tag_to_group_factory.copy()
+ if DISCR_TAG_BASE in discr_tag_to_group_factory:
+ raise ValueError(
+ "if 'order' is given, 'discr_tag_to_group_factory' must "
+ "not have a key of DISCR_TAG_BASE"
+ )
+
+ discr_tag_to_group_factory[DISCR_TAG_BASE] = \
PolynomialWarpAndBlendGroupFactory(order=order)
# Modal discr should always comes from the base discretization
- quad_tag_to_group_factory[QTAG_MODAL] = \
+ discr_tag_to_group_factory[DISCR_TAG_MODAL] = \
_generate_modal_group_factory(
- quad_tag_to_group_factory[QTAG_NONE]
+ discr_tag_to_group_factory[DISCR_TAG_BASE]
)
- self.quad_tag_to_group_factory = quad_tag_to_group_factory
+ self.discr_tag_to_group_factory = discr_tag_to_group_factory
from meshmode.discretization import Discretization
- self._volume_discr = Discretization(array_context, mesh,
- self.group_factory_for_quadrature_tag(QTAG_NONE))
+ self._volume_discr = Discretization(
+ array_context, mesh,
+ self.group_factory_for_discretization_tag(DISCR_TAG_BASE)
+ )
# {{{ management of discretization-scoped common subexpressions
@@ -112,6 +142,16 @@ class DiscretizationCollection:
self.mpi_communicator = mpi_communicator
+ @property
+ def quad_tag_to_group_factory(self):
+ warn("`DiscretizationCollection.quad_tag_to_group_factory` "
+ "is deprecated and will go away in 2022. Use "
+ "`DiscretizationCollection.discr_tag_to_group_factory` "
+ "instead.",
+ DeprecationWarning, stacklevel=2)
+
+ return self.discr_tag_to_group_factory
+
def get_management_rank_index(self):
return 0
@@ -123,7 +163,7 @@ class DiscretizationCollection:
== self.get_management_rank_index()
def _set_up_distributed_communication(self, mpi_communicator, array_context):
- from_dd = DOFDesc("vol", QTAG_NONE)
+ from_dd = DOFDesc("vol", DISCR_TAG_BASE)
boundary_connections = {}
@@ -135,13 +175,14 @@ class DiscretizationCollection:
raise RuntimeError("must supply an MPI communicator when using a "
"distributed mesh")
- grp_factory = self.group_factory_for_quadrature_tag(QTAG_NONE)
+ grp_factory = \
+ self.group_factory_for_discretization_tag(DISCR_TAG_BASE)
local_boundary_connections = {}
for i_remote_part in connected_parts:
local_boundary_connections[i_remote_part] = self.connection_from_dds(
from_dd, DOFDesc(BTAG_PARTITION(i_remote_part),
- QTAG_NONE))
+ DISCR_TAG_BASE))
from meshmode.distributed import MPIBoundaryCommSetupHelper
with MPIBoundaryCommSetupHelper(mpi_communicator, array_context,
@@ -156,9 +197,12 @@ class DiscretizationCollection:
return boundary_connections
def get_distributed_boundary_swap_connection(self, dd):
- if dd.quadrature_tag != QTAG_NONE:
+ if dd.discretization_tag is not DISCR_TAG_BASE:
# FIXME
- raise NotImplementedError("Distributed communication with quadrature")
+ raise NotImplementedError(
+ "Distributed communication with discretization tag "
+ f"{dd.discretization_tag} is not implemented."
+ )
assert isinstance(dd.domain_tag, DTAG_BOUNDARY)
assert isinstance(dd.domain_tag.tag, BTAG_PARTITION)
@@ -169,26 +213,27 @@ class DiscretizationCollection:
def discr_from_dd(self, dd):
dd = as_dofdesc(dd)
- qtag = dd.quadrature_tag
+ discr_tag = dd.discretization_tag
- if qtag is QTAG_MODAL:
+ if discr_tag is DISCR_TAG_MODAL:
return self._modal_discr(dd.domain_tag)
if dd.is_volume():
- if qtag is not QTAG_NONE:
- return self._quad_volume_discr(qtag)
+ if discr_tag is not DISCR_TAG_BASE:
+ return self._discr_tag_volume_discr(discr_tag)
return self._volume_discr
- if qtag is not QTAG_NONE:
+ if discr_tag is not DISCR_TAG_BASE:
no_quad_discr = self.discr_from_dd(DOFDesc(dd.domain_tag))
from meshmode.discretization import Discretization
return Discretization(
- self._setup_actx,
- no_quad_discr.mesh,
- self.group_factory_for_quadrature_tag(qtag))
+ self._setup_actx,
+ no_quad_discr.mesh,
+ self.group_factory_for_discretization_tag(discr_tag)
+ )
- assert qtag is QTAG_NONE
+ assert discr_tag is DISCR_TAG_BASE
if dd.domain_tag is FACE_RESTR_ALL:
return self._all_faces_volume_connection().to_discr
@@ -204,24 +249,26 @@ class DiscretizationCollection:
from_dd = as_dofdesc(from_dd)
to_dd = as_dofdesc(to_dd)
- to_qtag = to_dd.quadrature_tag
- from_qtag = from_dd.quadrature_tag
+ to_discr_tag = to_dd.discretization_tag
+ from_discr_tag = from_dd.discretization_tag
# {{{ mapping between modal and nodal representations
- if (from_qtag is QTAG_MODAL and to_qtag is not QTAG_MODAL):
+ if (from_discr_tag is DISCR_TAG_MODAL
+ and to_discr_tag is not DISCR_TAG_MODAL):
return self._modal_to_nodal_connection(to_dd)
- if (to_qtag is QTAG_MODAL and from_qtag is not QTAG_MODAL):
+ if (to_discr_tag is DISCR_TAG_MODAL
+ and from_discr_tag is not DISCR_TAG_MODAL):
return self._nodal_to_modal_connection(from_dd)
# }}}
- assert (to_qtag is not QTAG_MODAL and from_qtag is not QTAG_MODAL)
+ assert (to_discr_tag is not DISCR_TAG_MODAL
+ and from_discr_tag is not DISCR_TAG_MODAL)
- if (
- not from_dd.is_volume()
- and from_qtag == to_qtag
+ if (not from_dd.is_volume()
+ and from_discr_tag == to_discr_tag
and to_dd.domain_tag is FACE_RESTR_ALL):
faces_conn = self.connection_from_dds(
DOFDesc("vol"),
@@ -235,11 +282,11 @@ class DiscretizationCollection:
faces_conn, self.discr_from_dd(to_dd),
self.discr_from_dd(from_dd))
- # {{{ simplify domain + qtag change into chained
+ # {{{ simplify domain + discr_tag change into chained
if (from_dd.domain_tag != to_dd.domain_tag
- and from_qtag is QTAG_NONE
- and to_qtag is not QTAG_NONE):
+ and from_discr_tag is DISCR_TAG_BASE
+ and to_discr_tag is not DISCR_TAG_BASE):
from meshmode.discretization.connection import \
ChainedDiscretizationConnection
@@ -261,10 +308,10 @@ class DiscretizationCollection:
# {{{ generic to-quad
- # QTAG_MODAL is handled above
+ # DISCR_TAG_MODAL is handled above
if (from_dd.domain_tag == to_dd.domain_tag
- and from_qtag is QTAG_NONE
- and to_qtag is not QTAG_NONE):
+ and from_discr_tag is DISCR_TAG_BASE
+ and to_discr_tag is not DISCR_TAG_BASE):
from meshmode.discretization.connection.same_mesh import \
make_same_mesh_connection
@@ -276,11 +323,11 @@ class DiscretizationCollection:
# }}}
- if from_qtag is not QTAG_NONE:
+ if from_discr_tag is not DISCR_TAG_BASE:
raise ValueError("cannot interpolate *from* a "
"(non-interpolatory) quadrature grid")
- assert to_qtag is QTAG_NONE
+ assert to_discr_tag is DISCR_TAG_BASE
if from_dd.is_volume():
if to_dd.domain_tag is FACE_RESTR_ALL:
@@ -288,12 +335,12 @@ class DiscretizationCollection:
if to_dd.domain_tag is FACE_RESTR_INTERIOR:
return self._interior_faces_connection()
elif to_dd.is_boundary_or_partition_interface():
- assert from_qtag is QTAG_NONE
+ assert from_discr_tag is DISCR_TAG_BASE
return self._boundary_connection(to_dd.domain_tag.tag)
elif to_dd.is_volume():
from meshmode.discretization.connection import \
make_same_mesh_connection
- to_discr = self._quad_volume_discr(to_qtag)
+ to_discr = self._discr_tag_volume_discr(to_discr_tag)
from_discr = self._volume_discr
return make_same_mesh_connection(self._setup_actx, to_discr,
from_discr)
@@ -304,22 +351,32 @@ class DiscretizationCollection:
else:
raise ValueError("cannot interpolate from: " + str(from_dd))
- def group_factory_for_quadrature_tag(self, quadrature_tag):
+ def group_factory_for_quadrature_tag(self, discretization_tag):
+ warn("`DiscretizationCollection.group_factory_for_quadrature_tag` "
+ "is deprecated and will go away in 2022. Use "
+ "`DiscretizationCollection.group_factory_for_discretization_tag` "
+ "instead.",
+ DeprecationWarning, stacklevel=2)
+
+ return self.group_factory_for_discretization_tag(discretization_tag)
+
+ def group_factory_for_discretization_tag(self, discretization_tag):
"""
OK to override in user code to control mode/node choice.
"""
+ if discretization_tag is None:
+ discretization_tag = DISCR_TAG_BASE
- if quadrature_tag is None:
- quadrature_tag = QTAG_NONE
-
- return self.quad_tag_to_group_factory[quadrature_tag]
+ return self.discr_tag_to_group_factory[discretization_tag]
@memoize_method
- def _quad_volume_discr(self, quadrature_tag):
+ def _discr_tag_volume_discr(self, discretization_tag):
from meshmode.discretization import Discretization
- return Discretization(self._setup_actx, self._volume_discr.mesh,
- self.group_factory_for_quadrature_tag(quadrature_tag))
+ return Discretization(
+ self._setup_actx, self._volume_discr.mesh,
+ self.group_factory_for_discretization_tag(discretization_tag)
+ )
# {{{ modal to nodal connections
@@ -327,10 +384,10 @@ class DiscretizationCollection:
def _modal_discr(self, domain_tag):
from meshmode.discretization import Discretization
- discr_base = self.discr_from_dd(DOFDesc(domain_tag, QTAG_NONE))
+ discr_base = self.discr_from_dd(DOFDesc(domain_tag, DISCR_TAG_BASE))
return Discretization(
self._setup_actx, discr_base.mesh,
- self.group_factory_for_quadrature_tag(QTAG_MODAL)
+ self.group_factory_for_discretization_tag(DISCR_TAG_MODAL)
)
@memoize_method
@@ -370,10 +427,11 @@ class DiscretizationCollection:
@memoize_method
def _boundary_connection(self, boundary_tag):
return make_face_restriction(
- self._setup_actx,
- self._volume_discr,
- self.group_factory_for_quadrature_tag(QTAG_NONE),
- boundary_tag=boundary_tag)
+ self._setup_actx,
+ self._volume_discr,
+ self.group_factory_for_discretization_tag(DISCR_TAG_BASE),
+ boundary_tag=boundary_tag
+ )
# }}}
@@ -382,15 +440,16 @@ class DiscretizationCollection:
@memoize_method
def _interior_faces_connection(self):
return make_face_restriction(
- self._setup_actx,
- self._volume_discr,
- self.group_factory_for_quadrature_tag(QTAG_NONE),
- FACE_RESTR_INTERIOR,
-
- # FIXME: This will need to change as soon as we support
- # pyramids or other elements with non-identical face
- # types.
- per_face_groups=False)
+ self._setup_actx,
+ self._volume_discr,
+ self.group_factory_for_discretization_tag(DISCR_TAG_BASE),
+ FACE_RESTR_INTERIOR,
+
+ # FIXME: This will need to change as soon as we support
+ # pyramids or other elements with non-identical face
+ # types.
+ per_face_groups=False
+ )
@memoize_method
def opposite_face_connection(self):
@@ -408,15 +467,16 @@ class DiscretizationCollection:
@memoize_method
def _all_faces_volume_connection(self):
return make_face_restriction(
- self._setup_actx,
- self._volume_discr,
- self.group_factory_for_quadrature_tag(QTAG_NONE),
- FACE_RESTR_ALL,
-
- # FIXME: This will need to change as soon as we support
- # pyramids or other elements with non-identical face
- # types.
- per_face_groups=False)
+ self._setup_actx,
+ self._volume_discr,
+ self.group_factory_for_discretization_tag(DISCR_TAG_BASE),
+ FACE_RESTR_ALL,
+
+ # FIXME: This will need to change as soon as we support
+ # pyramids or other elements with non-identical face
+ # types.
+ per_face_groups=False
+ )
# }}}
@@ -451,7 +511,6 @@ class DiscretizationCollection:
@property
def order(self):
- from warnings import warn
warn("DiscretizationCollection.order is deprecated, "
"consider using the orders of element groups instead. "
"'order' will go away in 2021.",
@@ -463,7 +522,6 @@ class DiscretizationCollection:
class DGDiscretizationWithBoundaries(DiscretizationCollection):
def __init__(self, *args, **kwargs):
- from warnings import warn
warn("DGDiscretizationWithBoundaries is deprecated and will go away "
"in 2022. Use DiscretizationCollection instead.",
DeprecationWarning, stacklevel=2)
diff --git a/grudge/dof_desc.py b/grudge/dof_desc.py
index 0d83c37b91275a421c879d0315eba864f7ce5432..621f245e72ecc384e913bf81a78b20b5fc2f1345 100644
--- a/grudge/dof_desc.py
+++ b/grudge/dof_desc.py
@@ -29,14 +29,18 @@ from meshmode.discretization.connection import \
FACE_RESTR_INTERIOR, FACE_RESTR_ALL
from meshmode.mesh import \
BTAG_PARTITION, BTAG_ALL, BTAG_REALLY_ALL, BTAG_NONE
+from warnings import warn
+import sys
__doc__ = """
.. autoclass:: DTAG_SCALAR
.. autoclass:: DTAG_VOLUME_ALL
.. autoclass:: DTAG_BOUNDARY
-.. autoclass:: QTAG_NONE
-.. autoclass:: QTAG_MODAL
+
+.. autoclass:: DISCR_TAG_BASE
+.. autoclass:: DISCR_TAG_QUAD
+.. autoclass:: DISCR_TAG_MODAL
.. autoclass:: DOFDesc
.. autofunction:: as_dofdesc
@@ -50,19 +54,19 @@ __doc__ = """
# {{{ DOF description
class DTAG_SCALAR: # noqa: N801
- """A tag denoting scalar values."""
+ """A domain tag denoting scalar values."""
class DTAG_VOLUME_ALL: # noqa: N801
"""
- A tag denoting degrees of freedom defined
+ A domain tag denoting values defined
in all cell volumes.
"""
class DTAG_BOUNDARY: # noqa: N801
- """A tag describing the meaning of degrees of freedom
- on element boundaries which are adjacent to elements
+ """A domain tag describing the values on element
+ boundaries which are adjacent to elements
of another :class:`~meshmode.mesh.Mesh`.
.. attribute:: tag
@@ -96,21 +100,45 @@ class DTAG_BOUNDARY: # noqa: N801
return "<{}({})>".format(type(self).__name__, repr(self.tag))
-class QTAG_NONE: # noqa: N801
- """A quadrature tag indicating the use of a
+class DISCR_TAG_BASE: # noqa: N801
+ """A discretization tag indicating the use of a
basic discretization grid. This tag is used
- to distinguish the base discretization (`QTAG_NONE`)
- from quadrature (e.g. overintegration) grids.
+ to distinguish the base discretization from quadrature
+ (e.g. overintegration) or modal (:class:`DISCR_TAG_MODAL`)
+ discretizations.
"""
-class QTAG_MODAL: # noqa: N801
- """A quadrature tag indicating the use of a
+class DISCR_TAG_QUAD: # noqa: N801
+ """A discretization tag indicating the use of a
+ quadrature discretization grid. This tag is used
+ to distinguish the quadrature discretization
+ (e.g. overintegration) from modal (:class:`DISCR_TAG_MODAL`)
+ or base (:class:`DISCR_TAG_BASE`) discretizations.
+
+ For working with multiple quadrature grids, it is
+ recommended to create appropriate subclasses of
+ :class:`DISCR_TAG_QUAD` and define appropriate
+ :class:`DOFDesc` objects corresponding to each
+ subclass. For example:
+
+ .. code-block:: python
+
+ class CustomQuadTag(DISCR_TAG_QUAD):
+ "A custom quadrature discretization tag."
+
+ dd = DOFDesc(DTAG_VOLUME_ALL, CustomQuadTag)
+
+ """
+
+
+class DISCR_TAG_MODAL: # noqa: N801
+ """A discretization tag indicating the use of a
basic discretization grid with modal degrees of
freedom. This tag is used to distinguish the
- modal discretization (`QTAG_MODAL`) from
- the base (nodal) discretization (e.g. `QTAG_NONE`)
- or discretizations on quadrature grids.
+ modal discretization from the base (nodal)
+ discretization (e.g. :class:`DISCR_TAG_BASE`) or
+ discretizations on quadrature grids (:class:`DISCR_TAG_QUAD`).
"""
@@ -118,7 +146,7 @@ class DOFDesc:
"""Describes the meaning of degrees of freedom.
.. attribute:: domain_tag
- .. attribute:: quadrature_tag
+ .. attribute:: discretization_tag
.. automethod:: __init__
@@ -130,7 +158,7 @@ class DOFDesc:
.. automethod:: uses_quadrature
- .. automethod:: with_qtag
+ .. automethod:: with_discr_tag
.. automethod:: with_dtag
.. automethod:: __eq__
@@ -138,7 +166,9 @@ class DOFDesc:
.. automethod:: __hash__
"""
- def __init__(self, domain_tag, quadrature_tag=None):
+ def __init__(self, domain_tag, discretization_tag=None,
+ # FIXME: `quadrature_tag` is deprecated
+ quadrature_tag=None):
"""
:arg domain_tag: One of the following:
:class:`DTAG_SCALAR` (or the string ``"scalar"``),
@@ -154,11 +184,11 @@ class DOFDesc:
:class:`~meshmode.mesh.BTAG_PARTITION`,
or *None* to indicate that the geometry is not yet known.
- :arg quadrature_tag:
- *None* to indicate that the quadrature grid is not known, or
- :class:`QTAG_NONE` to indicate the use of the basic discretization
- grid, or a string to indicate the use of the thus-tagged quadratue
- grid.
+ :arg discretization_tag:
+ *None* or :class:`DISCR_TAG_BASE` to indicate the use of the basic
+ discretization grid, :class:`DISCR_TAG_MODAL` to indicate a
+ modal discretization, or :class:`DISCR_TAG_QUAD` to indicate
+ the use of a quadrature grid.
"""
if domain_tag is None:
@@ -180,14 +210,41 @@ class DOFDesc:
else:
raise ValueError("domain tag not understood: %s" % domain_tag)
- if domain_tag is DTAG_SCALAR and quadrature_tag is not None:
- raise ValueError("cannot have nontrivial quadrature tag on scalar")
+ if (quadrature_tag is not None and discretization_tag is not None):
+ raise ValueError(
+ "Both `quadrature_tag` and `discretization_tag` are specified. "
+ "Use `discretization_tag` instead."
+ )
+
+ # FIXME: `quadrature_tag` is deprecated
+ if (quadrature_tag is not None and discretization_tag is None):
+ warn("`quadrature_tag` is a deprecated kwarg and will be dropped "
+ "in version 2022.x. Use `discretization_tag` instead.",
+ DeprecationWarning, stacklevel=2)
+ discretization_tag = quadrature_tag
- if quadrature_tag is None:
- quadrature_tag = QTAG_NONE
+ if domain_tag is DTAG_SCALAR and discretization_tag is not None:
+ raise ValueError("cannot have nontrivial discretization tag on scalar")
+
+ if discretization_tag is None:
+ discretization_tag = DISCR_TAG_BASE
+
+ # FIXME: String tags are deprecated
+ if isinstance(discretization_tag, str):
+ warn("Support for string values of `discretization_tag` will "
+ "be dropped in version 2022.x. Use one of the `DISCR_TAG_` "
+ "tags instead.",
+ DeprecationWarning, stacklevel=2)
self.domain_tag = domain_tag
- self.quadrature_tag = quadrature_tag
+ self.discretization_tag = discretization_tag
+
+ @property
+ def quadrature_tag(self):
+ warn("`DOFDesc.quadrature_tag` is deprecated and will be dropped "
+ "in version 2022.x. Use `DOFDesc.discretization_tag` instead.",
+ DeprecationWarning, stacklevel=2)
+ return self.discretization_tag
def is_scalar(self):
return self.domain_tag is DTAG_SCALAR
@@ -208,29 +265,49 @@ class DOFDesc:
FACE_RESTR_INTERIOR])
def uses_quadrature(self):
- if self.quadrature_tag is None:
- return False
- if self.quadrature_tag is QTAG_NONE:
+ # FIXME: String tags are deprecated
+ # Check for string first, otherwise
+ # `issubclass` will raise an exception whenever
+ # its first argument is not a class.
+ # This can go away once support for strings is dropped
+ # completely.
+ if isinstance(self.discretization_tag, str):
+ # All strings are interpreted as quadrature-related tags
+ return True
+ elif issubclass(self.discretization_tag, DISCR_TAG_QUAD):
+ return True
+ elif issubclass(self.discretization_tag,
+ (DISCR_TAG_BASE, DISCR_TAG_MODAL)):
return False
+ else:
+ raise ValueError(
+ f"Unsure how to interpret tag: {self.discretization_tag}"
+ )
- return True
+ def with_qtag(self, discr_tag):
+ warn("`DOFDesc.with_qtag` is deprecated and will be dropped "
+ "in version 2022.x. Use `DOFDesc.with_discr_tag` instead.",
+ DeprecationWarning, stacklevel=2)
+ return self.with_discr_tag(discr_tag)
- def with_qtag(self, qtag):
- return type(self)(domain_tag=self.domain_tag, quadrature_tag=qtag)
+ def with_discr_tag(self, discr_tag):
+ return type(self)(domain_tag=self.domain_tag,
+ discretization_tag=discr_tag)
def with_dtag(self, dtag):
- return type(self)(domain_tag=dtag, quadrature_tag=self.quadrature_tag)
+ return type(self)(domain_tag=dtag,
+ discretization_tag=self.discretization_tag)
def __eq__(self, other):
return (type(self) == type(other)
and self.domain_tag == other.domain_tag
- and self.quadrature_tag == other.quadrature_tag)
+ and self.discretization_tag == other.discretization_tag)
def __ne__(self, other):
return not self.__eq__(other)
def __hash__(self):
- return hash((type(self), self.domain_tag, self.quadrature_tag))
+ return hash((type(self), self.domain_tag, self.discretization_tag))
def __repr__(self):
def fmt(s):
@@ -241,20 +318,44 @@ class DOFDesc:
return "DOFDesc({}, {})".format(
fmt(self.domain_tag),
- fmt(self.quadrature_tag))
+ fmt(self.discretization_tag))
DD_SCALAR = DOFDesc(DTAG_SCALAR, None)
DD_VOLUME = DOFDesc(DTAG_VOLUME_ALL, None)
-DD_VOLUME_MODAL = DOFDesc(DTAG_VOLUME_ALL, QTAG_MODAL)
+DD_VOLUME_MODAL = DOFDesc(DTAG_VOLUME_ALL, DISCR_TAG_MODAL)
def as_dofdesc(dd):
if isinstance(dd, DOFDesc):
return dd
- return DOFDesc(dd, quadrature_tag=None)
+ return DOFDesc(dd, discretization_tag=None)
+
+# }}}
+
+
+# {{{ Deprecated tags
+
+_deprecated_name_to_new_name = {"QTAG_NONE": "DISCR_TAG_BASE",
+ "QTAG_MODAL": "DISCR_TAG_MODAL"}
+
+
+def __getattr__(name):
+ if name in _deprecated_name_to_new_name:
+ warn(f"'{name}' is deprecated and will be dropped "
+ f"in version 2022.x. Use '{_deprecated_name_to_new_name[name]}' "
+ "instead.",
+ DeprecationWarning, stacklevel=2)
+ return globals()[_deprecated_name_to_new_name[name]]
+
+ raise AttributeError(f"module {__name__} has no attribute {name}")
+
+
+if sys.version_info < (3, 7):
+ for name in _deprecated_name_to_new_name:
+ globals()[name] = globals()[_deprecated_name_to_new_name[name]]
# }}}
diff --git a/grudge/execution.py b/grudge/execution.py
index a0624479210f84b72ccbec58b9372e6282025211..4b86fcb9e359c828a99f8f6d94a4714a999cd705 100644
--- a/grudge/execution.py
+++ b/grudge/execution.py
@@ -628,7 +628,7 @@ class BoundOperator:
else:
pass
- for key, val in context.items():
+ for val in context.values():
look_for_array_contexts(val)
if array_contexts:
@@ -719,7 +719,7 @@ def process_sym_operator(dcoll, sym_operator, post_bind_mapper=None, dumper=None
dumper("before-qcheck", sym_operator)
sym_operator = mappers.QuadratureCheckerAndRemover(
- dcoll.quad_tag_to_group_factory)(sym_operator)
+ dcoll.discr_tag_to_group_factory)(sym_operator)
# Work around https://github.com/numpy/numpy/issues/9438
#
diff --git a/grudge/models/advection.py b/grudge/models/advection.py
index 55ad4fc70add026ee43c25b8a2e6993183985e1f..980425e5be60079d46bd2067936cdad83f164ec7 100644
--- a/grudge/models/advection.py
+++ b/grudge/models/advection.py
@@ -200,11 +200,12 @@ def v_dot_n_tpair(velocity, dd=None):
dd = DOFDesc(FACE_RESTR_INTERIOR)
ambient_dim = len(velocity)
- normal = sym.normal(dd.with_qtag(None), ambient_dim, dim=ambient_dim - 2)
+ normal = sym.normal(dd.with_discr_tag(None),
+ ambient_dim, dim=ambient_dim - 2)
return sym.int_tpair(velocity.dot(normal),
- qtag=dd.quadrature_tag,
- from_dd=dd.with_qtag(None))
+ qtag=dd.discretization_tag,
+ from_dd=dd.with_discr_tag(None))
def surface_advection_weak_flux(flux_type, u, velocity):
@@ -231,7 +232,7 @@ class SurfaceAdvectionOperator(AdvectionOperatorBase):
def flux(self, u):
from grudge.dof_desc import DD_VOLUME
- surf_v = sym.project(DD_VOLUME, u.dd.with_qtag(None))(self.v)
+ surf_v = sym.project(DD_VOLUME, u.dd.with_discr_tag(None))(self.v)
return surface_advection_weak_flux(self.flux_type, u, surf_v)
def sym_operator(self):
diff --git a/grudge/models/em.py b/grudge/models/em.py
index 244d324e0749ef48d19a67af9ab0be5c3156aa0c..3f4bc13d42a017fcee36ba7b05237c512ca11451 100644
--- a/grudge/models/em.py
+++ b/grudge/models/em.py
@@ -327,7 +327,9 @@ class MaxwellOperator(HyperbolicOperator):
1 / sym.sqrt(self.epsilon * self.mu)
)
- def max_eigenvalue(self, t, fields=None, discr=None, context={}):
+ def max_eigenvalue(self, t, fields=None, discr=None, context=None):
+ if context is None:
+ context = {}
if self.fixed_material:
return self.max_eigenvalue_expr()
else:
diff --git a/grudge/models/poisson.py b/grudge/models/poisson.py
deleted file mode 100644
index 6636b94fe4d7ea90894bf92a86d30189fe93348e..0000000000000000000000000000000000000000
--- a/grudge/models/poisson.py
+++ /dev/null
@@ -1,261 +0,0 @@
-"""Operators for Poisson problems."""
-
-__copyright__ = "Copyright (C) 2007 Andreas Kloeckner"
-
-__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 grudge.models import Operator
-from grudge.second_order import LDGSecondDerivative
-import grudge.data
-import grudge.iterative
-
-
-class LaplacianOperatorBase:
- def sym_operator(self, apply_minv, u=None, dir_bc=None, neu_bc=None):
- """
- :param apply_minv: :class:`bool` specifying whether to compute a complete
- divergence operator. If False, the final application of the inverse
- mass operator is skipped. This is used in :meth:`op` in order to
- reduce the scheme :math:`M^{-1} S u = f` to :math:`S u = M f`, so
- that the mass operator only needs to be applied once, when preparing
- the right hand side in :meth:`prepare_rhs`.
-
- :class:`grudge.models.diffusion.DiffusionOperator` needs this.
- """
-
- from grudge.symbolic import Field, make_sym_vector
- from grudge.second_order import SecondDerivativeTarget
-
- if u is None:
- u = Field("u")
- if dir_bc is None:
- dir_bc = Field("dir_bc")
- if neu_bc is None:
- neu_bc = Field("neu_bc")
-
- # strong_form here allows IPDG to reuse the value of grad u.
- grad_tgt = SecondDerivativeTarget(
- self.dimensions, strong_form=True,
- operand=u)
-
- def grad_bc_getter(tag, expr):
- assert tag == self.dirichlet_tag
- return dir_bc
- self.scheme.grad(grad_tgt,
- bc_getter=grad_bc_getter,
- dirichlet_tags=[self.dirichlet_tag],
- neumann_tags=[self.neumann_tag])
-
- def apply_diff_tensor(v):
- if isinstance(self.diffusion_tensor, np.ndarray):
- sym_diff_tensor = self.diffusion_tensor
- else:
- sym_diff_tensor = (make_sym_vector(
- "diffusion", self.dimensions**2)
- .reshape(self.dimensions, self.dimensions))
-
- return np.dot(sym_diff_tensor, v)
-
- div_tgt = SecondDerivativeTarget(
- self.dimensions, strong_form=False,
- operand=apply_diff_tensor(grad_tgt.minv_all))
-
- def div_bc_getter(tag, expr):
- if tag == self.dirichlet_tag:
- return dir_bc
- elif tag == self.neumann_tag:
- return neu_bc
- else:
- assert False, "divergence bc getter " \
- "asked for '%s' BC for '%s'" % (tag, expr)
-
- self.scheme.div(div_tgt,
- div_bc_getter,
- dirichlet_tags=[self.dirichlet_tag],
- neumann_tags=[self.neumann_tag])
-
- if apply_minv:
- return div_tgt.minv_all
- else:
- return div_tgt.all
-
-
-class PoissonOperator(Operator, LaplacianOperatorBase):
- """Implements the Local Discontinuous Galerkin (LDG) Method for elliptic
- operators.
-
- See P. Castillo et al.,
- Local discontinuous Galerkin methods for elliptic problems",
- Communications in Numerical Methods in Engineering 18, no. 1 (2002): 69-75.
- """
-
- def __init__(self, dimensions, diffusion_tensor=None,
- dirichlet_bc=grudge.data.ConstantGivenFunction(),
- dirichlet_tag="dirichlet",
- neumann_bc=grudge.data.ConstantGivenFunction(),
- neumann_tag="neumann",
- scheme=LDGSecondDerivative()):
- self.dimensions = dimensions
-
- self.scheme = scheme
-
- self.dirichlet_bc = dirichlet_bc
- self.dirichlet_tag = dirichlet_tag
- self.neumann_bc = neumann_bc
- self.neumann_tag = neumann_tag
-
- if diffusion_tensor is None:
- diffusion_tensor = np.eye(dimensions)
- self.diffusion_tensor = diffusion_tensor
-
- # bound operator ----------------------------------------------------------
- def bind(self, discr):
- """Return a :class:`BoundPoissonOperator`."""
-
- assert self.dimensions == discr.dimensions
-
- from grudge.mesh import check_bc_coverage
- check_bc_coverage(discr.mesh, [self.dirichlet_tag, self.neumann_tag])
-
- return BoundPoissonOperator(self, discr)
-
-
-class BoundPoissonOperator(grudge.iterative.OperatorBase):
- """Returned by :meth:`PoissonOperator.bind`."""
-
- def __init__(self, poisson_op, discr):
- grudge.iterative.OperatorBase.__init__(self)
- self.discr = discr
-
- pop = self.poisson_op = poisson_op
-
- op = pop.sym_operator(
- apply_minv=False, dir_bc=0, neu_bc=0)
- bc_op = pop.sym_operator(apply_minv=False)
-
- self.compiled_op = discr.compile(op)
- self.compiled_bc_op = discr.compile(bc_op)
-
- if not isinstance(pop.diffusion_tensor, np.ndarray):
- self.diffusion = pop.diffusion_tensor.volume_interpolant(discr)
-
- # Check whether use of Poincaré mean-value method is required.
- # (for pure Neumann or pure periodic)
-
- from grudge.mesh import BTAG_ALL
- self.poincare_mean_value_hack = (
- len(self.discr.get_boundary(BTAG_ALL).nodes)
- == len(self.discr.get_boundary(poisson_op.neumann_tag).nodes))
-
- @property
- def dtype(self):
- return self.discr.default_scalar_type
-
- @property
- def shape(self):
- nodes = len(self.discr)
- return nodes, nodes
-
- def op(self, u):
- context = {"u": u}
- if not isinstance(self.poisson_op.diffusion_tensor, np.ndarray):
- context["diffusion"] = self.diffusion
-
- result = self.compiled_op(**context)
-
- if self.poincare_mean_value_hack:
- state_int = self.discr.integral(u)
- mean_state = state_int / self.discr.mesh_volume()
- return result - mean_state * self.discr._mass_ones()
- else:
- return result
-
- __call__ = op
-
- def prepare_rhs(self, rhs):
- """Prepare the right-hand side for the linear system op(u)=rhs(f).
-
- In matrix form, LDG looks like this:
-
- .. math::
- Mv = Cu + g
- Mf = Av + Bu + h
-
- where v is the auxiliary vector, u is the argument of the operator, f
- is the result of the grad operator, g and h are inhom boundary data, and
- A,B,C are some operator+lifting matrices.
-
- .. math::
-
- M f = A M^{-1}(Cu + g) + Bu + h
-
- so the linear system looks like
-
- .. math::
-
- M f = A M^{-1} Cu + A M^{-1} g + Bu + h
- M f - A M^{-1} g - h = (A M^{-1} C + B)u (*)
-
- So the right hand side we're putting together here is really
-
- .. math::
-
- M f - A M^{-1} g - h
-
- .. note::
-
- Resist the temptation to left-multiply by the inverse
- mass matrix, as this will result in a non-symmetric
- matrix, which (e.g.) a conjugate gradient Krylov
- solver will not take well.
- """
- pop = self.poisson_op
-
- from grudge.symbolic import MassOperator
- return (MassOperator().apply(self.discr, rhs)
- - self.compiled_bc_op(
- u=self.discr.volume_zeros(),
- dir_bc=pop.dirichlet_bc.boundary_interpolant(
- self.discr, pop.dirichlet_tag),
- neu_bc=pop.neumann_bc.boundary_interpolant(
- self.discr, pop.neumann_tag)))
-
-
-class HelmholtzOperator(PoissonOperator):
- def __init__(self, k, *args, **kwargs):
- PoissonOperator.__init__(self, *args, **kwargs)
- self.k = k
-
- def sym_operator(self, apply_minv, u=None, dir_bc=None, neu_bc=None):
- from grudge.symbolic import Field
- if u is None:
- u = Field("u")
-
- result = PoissonOperator.sym_operator(self,
- apply_minv, u, dir_bc, neu_bc)
-
- if apply_minv:
- return result + self.k**2 * u
- else:
- from grudge.symbolic import MassOperator
- return result + self.k**2 * MassOperator()(u)
diff --git a/grudge/models/second_order.py b/grudge/models/second_order.py
deleted file mode 100644
index 8f7251d329ea1cbbb960e368fe3c272f481ca36a..0000000000000000000000000000000000000000
--- a/grudge/models/second_order.py
+++ /dev/null
@@ -1,505 +0,0 @@
-"""Schemes for second-order derivatives."""
-
-__copyright__ = "Copyright (C) 2009 Andreas Kloeckner"
-
-__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
-import grudge.symbolic.mappers as mappers
-from grudge import sym
-
-
-# {{{ stabilization term generator
-
-class StabilizationTermGenerator(mappers.IdentityMapper):
- def __init__(self, flux_args):
- super().__init__()
- self.flux_args = flux_args
- self.flux_arg_lookup = {
- flux_arg: i for i, flux_arg in enumerate(flux_args)}
-
- def get_flux_arg_idx(self, expr, quad_above):
- from grudge.symbolic.mappers import QuadratureDetector
-
- quad_below = QuadratureDetector()(expr)
- if quad_above:
- if quad_below is not None:
- # Both the part of the expression above and below the
- # differentiation operator had quadrature upsamplers in it.
- # Since we're removing the differentiation operator, there are
- # now two layers of quadrature operators. We need to change the
- # inner layer to be the only layer.
-
- from grudge.symbolic.mappers import QuadratureUpsamplerChanger
- expr = QuadratureUpsamplerChanger(quad_above[0])(expr)
- else:
- # Only the part of the expression above the differentiation
- # operator had quadrature. Insert quadrature here, be done.
- expr = quad_above[0](expr)
- else:
- if quad_below is not None:
- # Only the part of the expression below the differentiation
- # operator had quadrature--the stuff above doesn't want it.
- # Get rid of it.
- from grudge.symbolic.mappers import QuadratureUpsamplerRemover
- expr = QuadratureUpsamplerRemover({}, do_warn=False)(expr)
- else:
- # No quadrature, no headaches.
- pass
-
- try:
- return self.flux_arg_lookup[expr]
- except KeyError:
- flux_arg_idx = len(self.flux_args)
- self.flux_arg_lookup[expr] = flux_arg_idx
- self.flux_args.append(expr)
- return flux_arg_idx
-
- def map_operator_binding(self, expr, quad_above=[]):
- from grudge.symbolic.operators import (
- DiffOperatorBase, FluxOperatorBase,
- InverseMassOperator,
- QuadratureInteriorFacesGridUpsampler)
-
- if isinstance(expr.op, DiffOperatorBase):
- flux_arg_idx = self.get_flux_arg_idx(expr.field, quad_above=quad_above)
-
- from grudge.symbolic import \
- WeakFormDiffOperatorBase, \
- StrongFormDiffOperatorBase
- if isinstance(expr.op, WeakFormDiffOperatorBase):
- factor = -1
- elif isinstance(expr.op, StrongFormDiffOperatorBase):
- factor = 1
- else:
- raise RuntimeError("unknown type of differentiation "
- "operator encountered by stab term generator")
-
- from grudge.flux import Normal, FluxScalarPlaceholder
- sph = FluxScalarPlaceholder(flux_arg_idx)
- return (factor
- * Normal(expr.op.xyz_axis)
- * (sph.int - sph.ext))
-
- elif isinstance(expr.op, FluxOperatorBase):
- return 0
- elif isinstance(expr.op, InverseMassOperator):
- return self.rec(expr.field, quad_above)
- elif isinstance(expr.op, QuadratureInteriorFacesGridUpsampler):
- if quad_above:
- raise RuntimeError("double quadrature upsampler found "
- "when generating stabilization term")
- return self.rec(expr.field, [expr.op])
- else:
- from grudge.symbolic.tools import pretty
- raise ValueError("stabilization term generator doesn't know "
- "what to do with '%s'" % pretty(expr))
-
- def map_variable(self, expr, quad_above=[]):
- from grudge.flux import FieldComponent
- return FieldComponent(
- self.get_flux_arg_idx(expr, quad_above),
- is_interior=True)
-
- def map_subscript(self, expr, quad_above=[]):
- from pymbolic.primitives import Variable
- assert isinstance(expr.aggregate, Variable)
- from grudge.flux import FieldComponent
- return FieldComponent(
- self.get_flux_arg_idx(expr, quad_above),
- is_interior=True)
-
-# }}}
-
-
-# {{{ neumann bc generator
-
-class NeumannBCGenerator(mappers.IdentityMapper):
- def __init__(self, tag, bc):
- super().__init__()
- self.tag = tag
- self.bc = bc
-
- def map_operator_binding(self, expr):
- if isinstance(expr.op, sym.DiffOperatorBase):
- from grudge.symbolic import \
- WeakFormDiffOperatorBase, \
- StrongFormDiffOperatorBase
- if isinstance(expr.op, WeakFormDiffOperatorBase):
- factor = -1
- elif isinstance(expr.op, StrongFormDiffOperatorBase):
- factor = 1
- else:
- raise RuntimeError("unknown type of differentiation "
- "operator encountered by stab term generator")
-
- from grudge.symbolic import BoundaryNormalComponent
- return (self.bc * factor
- * BoundaryNormalComponent(self.tag, expr.op.xyz_axis))
-
- elif isinstance(expr.op, sym.FluxOperatorBase):
- return 0
- elif isinstance(expr.op, sym.InverseMassOperator):
- return self.rec(expr.field)
- else:
- raise ValueError("neumann normal direction generator doesn't know "
- "what to do with '%s'" % expr)
-
-# }}}
-
-
-class IPDGDerivativeGenerator(mappers.IdentityMapper):
- def map_operator_binding(self, expr):
- if isinstance(expr.op, sym.DiffOperatorBase):
- from grudge.symbolic import (
- WeakFormDiffOperatorBase,
- StrongFormDiffOperatorBase)
-
- if isinstance(expr.op, WeakFormDiffOperatorBase):
- factor = -1
- elif isinstance(expr.op, StrongFormDiffOperatorBase):
- factor = 1
- else:
- raise RuntimeError("unknown type of differentiation "
- "operator encountered by stab term generator")
-
- from grudge.symbolic import DifferentiationOperator
- return factor*DifferentiationOperator(expr.op.xyz_axis)(expr.field)
-
- elif isinstance(expr.op, sym.FluxOperatorBase):
- return 0
- elif isinstance(expr.op, sym.InverseMassOperator):
- return self.rec(expr.field)
- elif isinstance(expr.op,
- sym.QuadratureInteriorFacesGridUpsampler):
- return super().map_operator_binding(expr)
- else:
- from grudge.symbolic.tools import pretty
- raise ValueError("IPDG derivative generator doesn't know "
- "what to do with '%s'" % pretty(expr))
-
-
-# {{{ second derivative target
-
-class SecondDerivativeTarget:
- def __init__(self, dimensions, strong_form, operand,
- int_flux_operand=None,
- bdry_flux_int_operand=None):
- """
- :param int_flux_operand: if not None, is used as the interior
- argument to the interior fluxes. This is useful e.g. if the boundary
- values are on a quadrature grid--in this case, *bdry_flux_int_operand*
- can be passed to also be on a boundary grid. If it is None, it defaults
- to *operand*.
-
- :param bdry_flux_int_operand: if not None, is used as the interior
- argument to the boundary fluxes. This is useful e.g. if the boundary
- values are on a quadrature grid--in this case, *bdry_flux_int_operand*
- can be passed to also be on a boundary grid. If it is None, it defaults
- to *int_flux_operand*.
- """
- self.dimensions = dimensions
- self.operand = operand
-
- if int_flux_operand is None:
- int_flux_operand = operand
- if bdry_flux_int_operand is None:
- bdry_flux_int_operand = int_flux_operand
-
- self.int_flux_operand = int_flux_operand
- self.bdry_flux_int_operand = bdry_flux_int_operand
-
- self.strong_form = strong_form
- if strong_form:
- self.strong_neg = -1
- else:
- self.strong_neg = 1
-
- self.local_derivatives = 0
- self.inner_fluxes = 0
- self.boundary_fluxes = 0
-
- def add_local_derivatives(self, expr):
- self.local_derivatives = self.local_derivatives \
- + (-self.strong_neg)*expr
-
- def vec_times(self, vec, operand):
- from pytools.obj_array import is_obj_array
-
- if is_obj_array(operand):
- if len(operand) != self.dimensions:
- raise ValueError("operand of vec_times must have %d dimensions"
- % self.dimensions)
-
- return np.dot(vec, operand)
- else:
- return vec*operand
-
- def normal_times_flux(self, flux):
- from grudge.flux import make_normal
- return self.vec_times(make_normal(self.dimensions), flux)
-
- def apply_diff(self, nabla, operand):
- from pytools.obj_array import make_obj_array, is_obj_array
- if is_obj_array(operand):
- if len(operand) != self.dimensions:
- raise ValueError("operand of apply_diff must have %d dimensions"
- % self.dimensions)
-
- return sum(nabla[i](operand[i]) for i in range(self.dimensions))
- else:
- return make_obj_array(
- [nabla[i](operand) for i in range(self.dimensions)])
-
- def _local_nabla(self):
- if self.strong_form:
- from grudge.symbolic import make_stiffness
- return make_stiffness(self.dimensions)
- else:
- from grudge.symbolic import make_stiffness_t
- return make_stiffness_t(self.dimensions)
-
- def add_derivative(self, operand=None):
- if operand is None:
- operand = self.operand
-
- self.add_local_derivatives(
- self.apply_diff(self._local_nabla(), operand))
-
- def add_inner_fluxes(self, flux, expr):
- from grudge.symbolic import get_flux_operator
- self.inner_fluxes = self.inner_fluxes \
- + get_flux_operator(self.strong_neg*flux)(expr)
-
- def add_boundary_flux(self, flux, volume_expr, bdry_expr, tag):
- from grudge.symbolic import BoundaryPair, get_flux_operator
- self.boundary_fluxes = self.boundary_fluxes + \
- get_flux_operator(self.strong_neg*flux)(BoundaryPair(
- volume_expr, bdry_expr, tag))
-
- @property
- def fluxes(self):
- return self.inner_fluxes + self.boundary_fluxes
-
- @property
- def all(self):
- return self.local_derivatives + self.fluxes
-
- @property
- def minv_all(self):
- from grudge.symbolic.primitives import make_common_subexpression as cse
- from grudge.symbolic.operators import InverseMassOperator
- return (cse(InverseMassOperator()(self.local_derivatives), "grad_loc")
- + cse(InverseMassOperator()(self.fluxes), "grad_flux"))
-
-# }}}
-
-
-# {{{ second derivative schemes
-
-class SecondDerivativeBase:
- def grad(self, tgt, bc_getter, dirichlet_tags, neumann_tags):
- """
- :param bc_getter: a function (tag, volume_expr) -> boundary expr.
- *volume_expr* will be None to query the Neumann condition.
- """
- n_times = tgt.normal_times_flux
-
- if tgt.strong_form:
- def adjust_flux(f):
- return n_times(u.int) - f
- else:
- def adjust_flux(f):
- return f
-
- from grudge.flux import FluxScalarPlaceholder
- u = FluxScalarPlaceholder()
-
- tgt.add_derivative()
- tgt.add_inner_fluxes(
- adjust_flux(self.grad_interior_flux(tgt, u)),
- tgt.int_flux_operand)
-
- for tag in dirichlet_tags:
- tgt.add_boundary_flux(
- adjust_flux(n_times(u.ext)),
- tgt.bdry_flux_int_operand, bc_getter(tag, tgt.operand), tag)
-
- for tag in neumann_tags:
- tgt.add_boundary_flux(
- adjust_flux(n_times(u.int)),
- tgt.bdry_flux_int_operand, 0, tag)
-
- def add_div_bcs(self, tgt, bc_getter, dirichlet_tags, neumann_tags,
- stab_term, adjust_flux, flux_v, flux_arg_int,
- grad_flux_arg_count):
- from pytools.obj_array import make_obj_array, join_fields
- n_times = tgt.normal_times_flux
-
- def unwrap_cse(expr):
- from pymbolic.primitives import CommonSubexpression
- if isinstance(expr, CommonSubexpression):
- return expr.child
- else:
- return expr
-
- for tag in dirichlet_tags:
- dir_bc_w = join_fields(
- [0]*grad_flux_arg_count,
- [bc_getter(tag, unwrap_cse(vol_expr)) for vol_expr in
- flux_arg_int[grad_flux_arg_count:]])
- tgt.add_boundary_flux(
- adjust_flux(n_times(flux_v.int-stab_term)),
- flux_arg_int, dir_bc_w, tag)
-
- loc_bc_vec = make_obj_array([0]*len(flux_arg_int))
-
- for tag in neumann_tags:
- neu_bc_w = join_fields(
- NeumannBCGenerator(tag, bc_getter(tag, None))(tgt.operand),
- [0]*len(flux_arg_int))
-
- tgt.add_boundary_flux(
- adjust_flux(n_times(flux_v.ext)),
- loc_bc_vec, neu_bc_w, tag)
-
-
-class LDGSecondDerivative(SecondDerivativeBase):
- def __init__(self, beta_value=0.5, stab_coefficient=1):
- self.beta_value = beta_value
- self.stab_coefficient = stab_coefficient
-
- def beta(self, tgt):
- return np.array([self.beta_value]*tgt.dimensions, dtype=np.float64)
-
- def grad_interior_flux(self, tgt, u):
- from grudge.symbolic.primitives import make_common_subexpression as cse
- n_times = tgt.normal_times_flux
- v_times = tgt.vec_times
-
- return n_times(
- cse(u.avg, "u_avg")
- - v_times(self.beta(tgt), n_times(u.int-u.ext)))
-
- def div(self, tgt, bc_getter, dirichlet_tags, neumann_tags):
- """
- :param bc_getter: a function (tag, volume_expr) -> boundary expr.
- *volume_expr* will be None to query the Neumann condition.
- """
-
- from grudge.symbolic.primitives import make_common_subexpression as cse
- from grudge.flux import FluxVectorPlaceholder, PenaltyTerm
-
- n_times = tgt.normal_times_flux
- v_times = tgt.vec_times
-
- if tgt.strong_form:
- def adjust_flux(f):
- return n_times(flux_v.int) - f
- else:
- def adjust_flux(f):
- return f
-
- flux_v = FluxVectorPlaceholder(tgt.dimensions)
-
- stab_term_generator = StabilizationTermGenerator(
- list(tgt.int_flux_operand))
- stab_term = (self.stab_coefficient * PenaltyTerm()
- * stab_term_generator(tgt.int_flux_operand))
-
- flux = n_times(flux_v.avg
- + v_times(self.beta(tgt),
- cse(n_times(flux_v.int - flux_v.ext), "jump_v"))
- - stab_term)
-
- from pytools.obj_array import make_obj_array
- flux_arg_int = cse(make_obj_array(stab_term_generator.flux_args))
-
- tgt.add_derivative(cse(tgt.operand))
- tgt.add_inner_fluxes(adjust_flux(flux), flux_arg_int)
-
- self.add_div_bcs(tgt, bc_getter, dirichlet_tags, neumann_tags,
- stab_term, adjust_flux, flux_v, flux_arg_int, tgt.dimensions)
-
-
-class StabilizedCentralSecondDerivative(LDGSecondDerivative):
- def __init__(self, stab_coefficient=1):
- LDGSecondDerivative.__init__(self, 0, stab_coefficient=stab_coefficient)
-
-
-class CentralSecondDerivative(LDGSecondDerivative):
- def __init__(self):
- LDGSecondDerivative.__init__(self, 0, 0)
-
-
-class IPDGSecondDerivative(SecondDerivativeBase):
- def __init__(self, stab_coefficient=1):
- self.stab_coefficient = stab_coefficient
-
- def grad_interior_flux(self, tgt, u):
- from grudge.symbolic.primitives import make_common_subexpression as cse
- n_times = tgt.normal_times_flux
- return n_times(cse(u.avg, "u_avg"))
-
- def div(self, tgt, bc_getter, dirichlet_tags, neumann_tags):
- """
- :param bc_getter: a function (tag, volume_expr) -> boundary expr.
- *volume_expr* will be None to query the Neumann condition.
- """
-
- from grudge.symbolic.primitives import make_common_subexpression as cse
- from grudge.flux import FluxVectorPlaceholder, PenaltyTerm
-
- n_times = tgt.normal_times_flux
-
- if tgt.strong_form:
- def adjust_flux(f):
- return n_times(flux_v.int) - f
- else:
- def adjust_flux(f):
- return f
-
- dim = tgt.dimensions
-
- flux_w = FluxVectorPlaceholder(2*tgt.dimensions)
- flux_v = flux_w[:dim]
- pure_diff_v = flux_w[dim:]
- flux_args = (
- list(tgt.int_flux_operand)
- + list(IPDGDerivativeGenerator()(tgt.int_flux_operand)))
-
- stab_term_generator = StabilizationTermGenerator(flux_args)
- stab_term = (self.stab_coefficient * PenaltyTerm()
- * stab_term_generator(tgt.int_flux_operand))
- flux = n_times(pure_diff_v.avg - stab_term)
-
- from pytools.obj_array import make_obj_array
- flux_arg_int = cse(make_obj_array(stab_term_generator.flux_args))
-
- tgt.add_derivative(cse(tgt.operand))
- tgt.add_inner_fluxes(adjust_flux(flux), flux_arg_int)
-
- self.add_div_bcs(tgt, bc_getter, dirichlet_tags, neumann_tags,
- stab_term, adjust_flux, flux_v, flux_arg_int, 2*tgt.dimensions)
-
-# }}}
-
-# vim: fdm=marker
diff --git a/grudge/op.py b/grudge/op.py
index c64c51b1d7a1dfb2576bdd06e07bab36ac5138ba..b3c402ca9c3f3c739f4e80378d5a284ad769f0d9 100644
--- a/grudge/op.py
+++ b/grudge/op.py
@@ -218,7 +218,7 @@ def weak_local_grad(dcoll, *args):
"""
if len(args) == 1:
vec, = args
- dd = dof_desc.DOFDesc("vol", dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc("vol", dof_desc.DISCR_TAG_BASE)
elif len(args) == 2:
dd, vec = args
else:
@@ -248,7 +248,7 @@ def weak_local_d_dx(dcoll, *args):
"""
if len(args) == 2:
xyz_axis, vec = args
- dd = dof_desc.DOFDesc("vol", dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc("vol", dof_desc.DISCR_TAG_BASE)
elif len(args) == 3:
dd, xyz_axis, vec = args
else:
@@ -273,7 +273,7 @@ def weak_local_div(dcoll, *args):
"""
if len(args) == 1:
vecs, = args
- dd = dof_desc.DOFDesc("vol", dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc("vol", dof_desc.DISCR_TAG_BASE)
elif len(args) == 2:
dd, vecs = args
else:
@@ -297,7 +297,7 @@ def _bound_mass(dcoll, dd):
def mass(dcoll, *args):
if len(args) == 1:
vec, = args
- dd = dof_desc.DOFDesc("vol", dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc("vol", dof_desc.DISCR_TAG_BASE)
elif len(args) == 2:
dd, vec = args
else:
@@ -333,7 +333,7 @@ def _bound_face_mass(dcoll, dd):
def face_mass(dcoll, *args):
if len(args) == 1:
vec, = args
- dd = dof_desc.DOFDesc("all_faces", dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc("all_faces", dof_desc.DISCR_TAG_BASE)
elif len(args) == 2:
dd, vec = args
else:
diff --git a/grudge/symbolic/compiler.py b/grudge/symbolic/compiler.py
index 742483aed9e02264c5ac65cfb40ee8ef0dbb8ee5..f5f7c66499b2adeadc470b9be416a06d5074048a 100644
--- a/grudge/symbolic/compiler.py
+++ b/grudge/symbolic/compiler.py
@@ -796,7 +796,7 @@ def aggregate_assignments(inf_mapper, instructions, result,
schedulable.sort()
if schedulable:
- for key, name, expr in schedulable:
+ for _key, name, expr in schedulable:
ordered_names_exprs.append((name, expr))
available_names.add(name)
else:
diff --git a/grudge/symbolic/dofdesc_inference.py b/grudge/symbolic/dofdesc_inference.py
index 9b5466190d9e7a303ea4ed967932e9a68ace1bdb..894a9bfc0e779e780b5da78e7a0aabb74cc925b2 100644
--- a/grudge/symbolic/dofdesc_inference.py
+++ b/grudge/symbolic/dofdesc_inference.py
@@ -48,8 +48,8 @@ def unify_dofdescs(dd_a, dd_b, expr=None):
raise ValueError("mismatched domain tags " + loc_str)
# domain tags match
- if dd_a.quadrature_tag != dd_b.quadrature_tag:
- raise ValueError("mismatched quadrature tags " + loc_str)
+ if dd_a.discretization_tag != dd_b.discretization_tag:
+ raise ValueError("mismatched discretization tags " + loc_str)
return dd_a
diff --git a/grudge/symbolic/mappers/__init__.py b/grudge/symbolic/mappers/__init__.py
index 99b4d0897f22c112e3c3af67a3b65f08bb4de2dc..84c11feb883db70b63bcac293c6b672c6e1af605 100644
--- a/grudge/symbolic/mappers/__init__.py
+++ b/grudge/symbolic/mappers/__init__.py
@@ -537,21 +537,21 @@ class OperatorSpecializer(CSECachingMapperMixin, IdentityMapper):
if isinstance(expr.op, op.MassOperator) and has_quad_operand:
return op.QuadratureMassOperator(
- field_repr_tag.quadrature_tag)(self.rec(expr.field))
+ field_repr_tag.discretization_tag)(self.rec(expr.field))
elif isinstance(expr.op, op.RefMassOperator) and has_quad_operand:
return op.RefQuadratureMassOperator(
- field_repr_tag.quadrature_tag)(self.rec(expr.field))
+ field_repr_tag.discretization_tag)(self.rec(expr.field))
elif (isinstance(expr.op, op.StiffnessTOperator) and has_quad_operand):
return op.QuadratureStiffnessTOperator(
- expr.op.xyz_axis, field_repr_tag.quadrature_tag)(
+ expr.op.xyz_axis, field_repr_tag.discretization_tag)(
self.rec(expr.field))
elif (isinstance(expr.op, op.RefStiffnessTOperator)
and has_quad_operand):
return op.RefQuadratureStiffnessTOperator(
- expr.op.xyz_axis, field_repr_tag.quadrature_tag)(
+ expr.op.xyz_axis, field_repr_tag.discretization_tag)(
self.rec(expr.field))
elif (isinstance(expr.op, op.QuadratureGridUpsampler)
@@ -560,11 +560,11 @@ class OperatorSpecializer(CSECachingMapperMixin, IdentityMapper):
# if (isinstance(expr.field, OperatorBinding)
# and isinstance(expr.field.op, RestrictToBoundary)):
# return QuadratureRestrictToBoundary(
- # expr.field.op.tag, expr.op.quadrature_tag)(
+ # expr.field.op.tag, expr.op.discretization_tag)(
# self.rec(expr.field.field))
return op.QuadratureBoundaryGridUpsampler(
- expr.op.quadrature_tag, field_type.boundary_tag)(expr.field)
+ expr.op.discretization_tag, field_type.boundary_tag)(expr.field)
# }}}
elif isinstance(expr.op, op.RestrictToBoundary) and has_quad_operand:
@@ -612,7 +612,7 @@ class GlobalToReferenceMapper(CSECachingMapperMixin, IdentityMapper):
jac_in = sym.area_element(self.ambient_dim, dim, dd=dd_in)
jac_noquad = sym.area_element(self.ambient_dim, dim,
- dd=dd_in.with_qtag(dof_desc.QTAG_NONE))
+ dd=dd_in.with_discr_tag(dof_desc.DISCR_TAG_BASE))
def rewrite_derivative(ref_class, field, dd_in, with_jacobian=True):
def imd(rst):
@@ -713,12 +713,12 @@ class StringifyMapper(pymbolic.mapper.stringifier.StringifyMapper):
else:
result = fmt(dd.domain_tag)
- if dd.quadrature_tag is None:
+ if dd.discretization_tag is None:
pass
- elif dd.quadrature_tag is dof_desc.QTAG_NONE:
+ elif dd.discretization_tag is dof_desc.DISCR_TAG_BASE:
result += "q"
else:
- result += "Q"+fmt(dd.quadrature_tag)
+ result += "Q"+fmt(dd.discretization_tag)
return result
@@ -868,23 +868,24 @@ class QuadratureCheckerAndRemover(CSECachingMapperMixin, IdentityMapper):
"""Checks whether all quadratu
"""
- def __init__(self, quad_tag_to_group_factory):
+ def __init__(self, discr_tag_to_group_factory):
IdentityMapper.__init__(self)
CSECachingMapperMixin.__init__(self)
- self.quad_tag_to_group_factory = quad_tag_to_group_factory
+ self.discr_tag_to_group_factory = discr_tag_to_group_factory
map_common_subexpression_uncached = \
IdentityMapper.map_common_subexpression
def _process_dd(self, dd, location_descr):
- if dd.quadrature_tag is not dof_desc.QTAG_NONE:
- if dd.quadrature_tag not in self.quad_tag_to_group_factory:
+
+ if dd.discretization_tag is not dof_desc.DISCR_TAG_BASE:
+ if dd.discretization_tag not in self.discr_tag_to_group_factory:
raise ValueError("found unknown quadrature tag '%s' in '%s'"
- % (dd.quadrature_tag, location_descr))
+ % (dd.discretization_tag, location_descr))
- grp_factory = self.quad_tag_to_group_factory[dd.quadrature_tag]
+ grp_factory = self.discr_tag_to_group_factory[dd.discretization_tag]
if grp_factory is None:
- dd = dof_desc.DOFDesc(dd.domain_tag, dof_desc.QTAG_NONE)
+ dd = dof_desc.DOFDesc(dd.domain_tag, dof_desc.DISCR_TAG_BASE)
return dd
diff --git a/grudge/symbolic/operators.py b/grudge/symbolic/operators.py
index 701661ca9b09a99da74bd2319b7139c43fcfc782..031ec60e0fab4ea8847ea4a4c6929269d98658d5 100644
--- a/grudge/symbolic/operators.py
+++ b/grudge/symbolic/operators.py
@@ -264,7 +264,7 @@ class DiffOperatorBase(Operator):
dd_in = dof_desc.DD_VOLUME
if dd_out is None:
- dd_out = dd_in.with_qtag(dof_desc.QTAG_NONE)
+ dd_out = dd_in.with_discr_tag(dof_desc.DISCR_TAG_BASE)
else:
dd_out = dof_desc.as_dofdesc(dd_out)
@@ -322,7 +322,7 @@ class RefDiffOperatorBase(ElementwiseLinearOperator):
dd_in = dof_desc.DD_VOLUME
if dd_out is None:
- dd_out = dd_in.with_qtag(dof_desc.QTAG_NONE)
+ dd_out = dd_in.with_discr_tag(dof_desc.DISCR_TAG_BASE)
if dd_out.uses_quadrature():
raise ValueError("differentiation outputs are not on "
@@ -633,7 +633,7 @@ class FaceMassOperatorBase(ElementwiseLinearOperator):
dd_in = dof_desc.DOFDesc(FACE_RESTR_ALL, None)
if dd_out is None or dd_out == "vol":
- dd_out = dof_desc.DOFDesc("vol", dof_desc.QTAG_NONE)
+ dd_out = dof_desc.DOFDesc("vol", dof_desc.DISCR_TAG_BASE)
if dd_out.uses_quadrature():
raise ValueError("face mass operator outputs are not on "
diff --git a/grudge/symbolic/primitives.py b/grudge/symbolic/primitives.py
index 4fc28a733daa9d5578395f14f9c622b86516173a..c880d9d0925d0b8c48a01b788cf81249160865e9 100644
--- a/grudge/symbolic/primitives.py
+++ b/grudge/symbolic/primitives.py
@@ -35,12 +35,6 @@ from pymbolic.primitives import (
make_common_subexpression as cse)
from pymbolic.geometric_algebra import MultiVector
-# FIXME: Need to import these for backwards compatibility with
-# Mirgecom. Need to remove importing the following from
-# `grudge.symbolic.primitives` in https://github.com/illinois-ceesd/mirgecom
-from grudge.dof_desc import \
- DOFDesc, as_dofdesc, DTAG_BOUNDARY, QTAG_NONE # noqa: F401
-
class ExpressionBase(prim.Expression):
def make_stringifier(self, originating_stringifier=None):
@@ -368,7 +362,7 @@ def forward_metric_nth_derivative(xyz_axis, ref_axes, dd=None):
if dd is None:
dd = dof_desc.DD_VOLUME
- inner_dd = dd.with_qtag(dof_desc.QTAG_NONE)
+ inner_dd = dd.with_discr_tag(dof_desc.DISCR_TAG_BASE)
from pytools import flatten
flat_ref_axes = flatten([rst_axis] * n for rst_axis, n in ref_axes)
diff --git a/setup.cfg b/setup.cfg
index e8791cc48f86a67e726046a2f7510f32a39a98a5..38dcd9de2c01e21191ca1612388f36ecbc251b70 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -9,11 +9,12 @@ exclude=
grudge/models/nd_calculus.py,
grudge/dt_finding.py
-
inline-quotes = "
docstring-quotes = """
multiline-quotes = """
+# enable-flake8-bugbear
+
[mypy]
ignore_missing_imports = True
diff --git a/test/test_grudge.py b/test/test_grudge.py
index 1b375bdfe7241422142b4d8511235738b214253c..f246a895c2bf2b1229986c52c152d405aebf68b4 100644
--- a/test/test_grudge.py
+++ b/test/test_grudge.py
@@ -93,8 +93,9 @@ def test_inverse_metric(actx_factory, dim):
# {{{ mass operator trig integration
@pytest.mark.parametrize("ambient_dim", [1, 2, 3])
-@pytest.mark.parametrize("quad_tag", [dof_desc.QTAG_NONE, "OVSMP"])
-def test_mass_mat_trig(actx_factory, ambient_dim, quad_tag):
+@pytest.mark.parametrize("discr_tag", [dof_desc.DISCR_TAG_BASE,
+ dof_desc.DISCR_TAG_QUAD])
+def test_mass_mat_trig(actx_factory, ambient_dim, discr_tag):
"""Check the integral of some trig functions on an interval using the mass
matrix.
"""
@@ -108,19 +109,21 @@ def test_mass_mat_trig(actx_factory, ambient_dim, quad_tag):
true_integral = 13*np.pi/2 * (b - a)**(ambient_dim - 1)
from meshmode.discretization.poly_element import QuadratureSimplexGroupFactory
- dd_quad = dof_desc.DOFDesc(dof_desc.DTAG_VOLUME_ALL, quad_tag)
- if quad_tag is dof_desc.QTAG_NONE:
- quad_tag_to_group_factory = {}
+ dd_quad = dof_desc.DOFDesc(dof_desc.DTAG_VOLUME_ALL, discr_tag)
+ if discr_tag is dof_desc.DISCR_TAG_BASE:
+ discr_tag_to_group_factory = {}
else:
- quad_tag_to_group_factory = {
- quad_tag: QuadratureSimplexGroupFactory(order=2*order)
- }
+ discr_tag_to_group_factory = {
+ discr_tag: QuadratureSimplexGroupFactory(order=2*order)
+ }
mesh = mgen.generate_regular_rect_mesh(
a=(a,)*ambient_dim, b=(b,)*ambient_dim,
nelements_per_axis=(nel_1d,)*ambient_dim, order=1)
- discr = DiscretizationCollection(actx, mesh, order=order,
- quad_tag_to_group_factory=quad_tag_to_group_factory)
+ discr = DiscretizationCollection(
+ actx, mesh, order=order,
+ discr_tag_to_group_factory=discr_tag_to_group_factory
+ )
def _get_variables_on(dd):
sym_f = sym.var("f", dd=dd)
@@ -147,7 +150,7 @@ def test_mass_mat_trig(actx_factory, ambient_dim, quad_tag):
err_2 = abs(num_integral_2 - true_integral)
assert err_2 < 1.0e-9, err_2
- if quad_tag is dof_desc.QTAG_NONE:
+ if discr_tag is dof_desc.DISCR_TAG_BASE:
# NOTE: `integral` always makes a square mass matrix and
# `QuadratureSimplexGroupFactory` does not have a `mass_matrix` method.
num_integral_3 = bind(discr,
@@ -568,17 +571,19 @@ def test_surface_divergence_theorem(actx_factory, mesh_name, visualize=False):
from meshmode.discretization.poly_element import \
QuadratureSimplexGroupFactory
- discr = DiscretizationCollection(actx, mesh, order=builder.order,
- quad_tag_to_group_factory={
+ discr = DiscretizationCollection(
+ actx, mesh, order=builder.order,
+ discr_tag_to_group_factory={
"product": QuadratureSimplexGroupFactory(2 * builder.order)
- })
+ }
+ )
volume = discr.discr_from_dd(dof_desc.DD_VOLUME)
logger.info("ndofs: %d", volume.ndofs)
logger.info("nelements: %d", volume.mesh.nelements)
dd = dof_desc.DD_VOLUME
- dq = dd.with_qtag("product")
+ dq = dd.with_discr_tag("product")
df = dof_desc.as_dofdesc(FACE_RESTR_ALL)
ambient_dim = discr.ambient_dim
dim = discr.dim
@@ -910,14 +915,16 @@ def test_improvement_quadrature(actx_factory, order):
order=order)
if use_quad:
- quad_tag_to_group_factory = {
+ discr_tag_to_group_factory = {
"product": QuadratureSimplexGroupFactory(order=4*order)
- }
+ }
else:
- quad_tag_to_group_factory = {"product": None}
+ discr_tag_to_group_factory = {"product": None}
- discr = DiscretizationCollection(actx, mesh, order=order,
- quad_tag_to_group_factory=quad_tag_to_group_factory)
+ discr = DiscretizationCollection(
+ actx, mesh, order=order,
+ discr_tag_to_group_factory=discr_tag_to_group_factory
+ )
bound_op = bind(discr, op.sym_operator())
fields = bind(discr, gaussian_mode())(actx, t=0)
diff --git a/test/test_modal_connections.py b/test/test_modal_connections.py
index 539a78069a9363323d07d1f26c1f3dcebf5c88eb..eeb9a9782f71ffc99cc41093ed84296a3155f599 100644
--- a/test/test_modal_connections.py
+++ b/test/test_modal_connections.py
@@ -66,8 +66,8 @@ def test_inverse_modal_connections(actx_factory, nodal_group_factory):
dcoll = DiscretizationCollection(
actx, mesh,
- quad_tag_to_group_factory={
- dof_desc.QTAG_NONE: nodal_group_factory(order)
+ discr_tag_to_group_factory={
+ dof_desc.DISCR_TAG_BASE: nodal_group_factory(order)
}
)
@@ -106,15 +106,16 @@ def test_inverse_modal_connections_quadgrid(actx_factory):
dcoll = DiscretizationCollection(
actx, mesh,
- quad_tag_to_group_factory={
- dof_desc.QTAG_NONE: PolynomialWarpAndBlendGroupFactory(order),
- "quad": QuadratureSimplexGroupFactory(2*order)
+ discr_tag_to_group_factory={
+ dof_desc.DISCR_TAG_BASE: PolynomialWarpAndBlendGroupFactory(order),
+ dof_desc.DISCR_TAG_QUAD: QuadratureSimplexGroupFactory(2*order)
}
)
# Use dof descriptors on the quadrature grid
dd_modal = dof_desc.DD_VOLUME_MODAL
- dd_quad = dof_desc.DOFDesc(dof_desc.DTAG_VOLUME_ALL, "quad")
+ dd_quad = dof_desc.DOFDesc(dof_desc.DTAG_VOLUME_ALL,
+ dof_desc.DISCR_TAG_QUAD)
x_quad = thaw(actx, dcoll.discr_from_dd(dd_quad).nodes()[0])
quad_f = f(x_quad)