From 8ea19846b5b0de3144c9fc004f7d4056d8329fe6 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 10:37:08 +0300 Subject: [PATCH 01/16] pyproject: move project metadata from setup.py --- pyproject.toml | 86 +++++++++++++++++++++++++++++++++++++++++++++++--- pytest.ini | 3 -- setup.py | 58 ---------------------------------- 3 files changed, 81 insertions(+), 66 deletions(-) delete mode 100644 pytest.ini delete mode 100644 setup.py diff --git a/pyproject.toml b/pyproject.toml index 293bacf3b..51c740632 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,3 +1,81 @@ +[build-system] +build-backend = "setuptools.build_meta" +requires = [ + "setuptools>=63", +] + +[project] +name = "grudge" +version = "2021.1" +description = "Discretize discontinuous Galerkin operators quickly on heterogeneous hardware" +readme = "README.rst" +license = { text = "MIT" } +authors = [ + { name = "Andreas Kloeckner", email = "inform@tiker.net" }, +] +requires-python = ">=3.8" +classifiers = [ + "Development Status :: 3 - Alpha", + "Intended Audience :: Developers", + "Intended Audience :: Other Audience", + "Intended Audience :: Science/Research", + "License :: OSI Approved :: MIT License", + "Natural Language :: English", + "Programming Language :: Python", + "Programming Language :: Python :: 3", + "Topic :: Scientific/Engineering", + "Topic :: Scientific/Engineering :: Information Analysis", + "Topic :: Scientific/Engineering :: Mathematics", + "Topic :: Scientific/Engineering :: Visualization", + "Topic :: Software Development :: Libraries", + "Topic :: Utilities", +] +dependencies = [ + "arraycontext>=2021.1" + "immutables>=0.19", + "loopy>=2024.1", + "meshmode>=2021.2", + "modepy>=2021.1", + "pymbolic>=2022.2", + "pyopencl>=2022.1", + "pytools>=2024.1.3", +] + +[project.optional-dependencies] +all = [ + "dagrt>=2021.1", + "leap>=2021.1", + "meshpy>=2022.1", + "mpi4py", + "pymetis>=2023.1", + "pytato>=2021.1", + "pyvisfile>=2022.1", +] +doc = [ + "furo", + "sphinx-copybutton", + "sphinx>=4", +] +test = [ + "mypy", + "pytest", + "ruff" +] + +[project.urls] +Documentation = "https://documen.tician.de/grudge" +Homepage = "https://github.com/inducer/grudge" + +[tool.setuptools.packages.find] +include = [ + "grudge*", +] + +[tool.pytest.ini_options] +markers = [ + "mpi: mark a test using MPI", +] + [tool.ruff] preview = true @@ -7,14 +85,14 @@ extend-select = [ "C", # flake8-comprehensions "E", # pycodestyle "F", # pyflakes + "G", # flake8-logging-format "I", # flake8-isort "N", # pep8-naming "NPY", # numpy "Q", # flake8-quotes + "RUF", # ruff + "UP", # pyupgrade "W", # pycodestyle - "RUF", - "UP", - "G", ] extend-ignore = [ "C90", # McCabe complexity @@ -54,5 +132,3 @@ lines-after-imports = 2 [tool.mypy] warn_unused_ignores = true - - diff --git a/pytest.ini b/pytest.ini deleted file mode 100644 index 01a22d375..000000000 --- a/pytest.ini +++ /dev/null @@ -1,3 +0,0 @@ -[pytest] -markers = - mpi: marks tests as using MPI diff --git a/setup.py b/setup.py deleted file mode 100644 index 3c9e6488e..000000000 --- a/setup.py +++ /dev/null @@ -1,58 +0,0 @@ -#!/usr/bin/env python3 - - -def main(): - from setuptools import find_packages, setup - - version_dict = {} - init_filename = "grudge/version.py" - exec(compile(open(init_filename).read(), init_filename, "exec"), - version_dict) - - setup( - name="grudge", - version=version_dict["VERSION_TEXT"], - description=( - "Discretize discontinuous Galerkin operators quickly, " - "on heterogeneous hardware" - ), - long_description=open("README.rst").read(), - author="Andreas Kloeckner", - author_email="inform@tiker.net", - license="MIT", - url="https://github.com/inducer/grudge", - classifiers=[ - "Development Status :: 3 - Alpha", - "Intended Audience :: Developers", - "Intended Audience :: Other Audience", - "Intended Audience :: Science/Research", - "License :: OSI Approved :: MIT License", - "Natural Language :: English", - "Programming Language :: Python", - "Programming Language :: Python :: 3", - "Topic :: Scientific/Engineering", - "Topic :: Scientific/Engineering :: Information Analysis", - "Topic :: Scientific/Engineering :: Mathematics", - "Topic :: Scientific/Engineering :: Visualization", - "Topic :: Software Development :: Libraries", - "Topic :: Utilities", - ], - packages=find_packages(), - python_requires="~=3.8", - install_requires=[ - "pytest>=2.3", - "pytools>=2024.1.3", - "modepy>=2013.3", - "arraycontext>=2021.1", - "meshmode>=2020.2", - "pyopencl>=2013.1", - "pymbolic>=2013.2", - "loopy>=2020.2", - "cgen>=2013.1.2", - "immutabledict", - ], - ) - - -if __name__ == "__main__": - main() From e870a5fef7b7bc676a98820b0fc210cb1920469d Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 10:41:36 +0300 Subject: [PATCH 02/16] version: use importlib --- doc/conf.py | 43 ++++++++++++++----------------------------- grudge/version.py | 17 +++++++++++++++-- 2 files changed, 29 insertions(+), 31 deletions(-) diff --git a/doc/conf.py b/doc/conf.py index d21bde2a9..89ba77226 100644 --- a/doc/conf.py +++ b/doc/conf.py @@ -1,3 +1,5 @@ +import os +from importlib import metadata from urllib.request import urlopen @@ -9,39 +11,22 @@ extensions = globals()["extensions"] + [ "matplotlib.sphinxext.plot_directive"] -copyright = "2015-21, grudge contributors" -author = "grudge contributors" - - -def get_version(): - conf = {} - src = "../grudge/version.py" - exec( - compile(open(src).read(), src, "exec"), - conf) - return conf["VERSION_TEXT"] - - -version = get_version() - -# The full version, including alpha/beta/rc tags. -release = version +copyright = "2015-2024, Grudge contributors" +author = "Grudge contributors" +release = metadata.version("grudge") +version = ".".join(release.split(".")[:2]) intersphinx_mapping = { - "python": ("https://docs.python.org/3/", None), - "numpy": ("https://numpy.org/doc/stable/", None), - "pyopencl": ("https://documen.tician.de/pyopencl/", None), - "modepy": ("https://documen.tician.de/modepy/", None), - "pymbolic": ("https://documen.tician.de/pymbolic/", None), "arraycontext": ("https://documen.tician.de/arraycontext/", None), - "meshmode": ("https://documen.tician.de/meshmode/", None), "loopy": ("https://documen.tician.de/loopy/", None), + "meshmode": ("https://documen.tician.de/meshmode/", None), + "modepy": ("https://documen.tician.de/modepy/", None), "mpi4py": ("https://mpi4py.readthedocs.io/en/stable", None), - } + "numpy": ("https://numpy.org/doc/stable/", None), + "pymbolic": ("https://documen.tician.de/pymbolic/", None), + "pyopencl": ("https://documen.tician.de/pyopencl/", None), + "python": ("https://docs.python.org/3/", None), +} # index-page demo uses pyopencl via plot_directive -import os - - -# switch to "port:cpu" once we're firmly migrated to pocl 4.0 -os.environ["PYOPENCL_TEST"] = "port:0" +os.environ["PYOPENCL_TEST"] = "port:cpu" diff --git a/grudge/version.py b/grudge/version.py index 30c962914..ef4d357c8 100644 --- a/grudge/version.py +++ b/grudge/version.py @@ -1,2 +1,15 @@ -VERSION = (2021, 1) -VERSION_TEXT = ".".join(str(i) for i in VERSION) +from importlib import metadata +from typing import Tuple + + +def _parse_version(version: str) -> Tuple[Tuple[int, ...], str]: + import re + + m = re.match("^([0-9.]+)([a-z0-9]*?)$", VERSION_TEXT) + assert m is not None + + return tuple(int(nr) for nr in m.group(1).split(".")), m.group(2) + + +VERSION_TEXT = metadata.version("grudge") +VERSION, VERSION_STATUS = _parse_version(VERSION_TEXT) From 5945997aac89923419501bef7f7ed8db6931ae99 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 10:50:40 +0300 Subject: [PATCH 03/16] pyproject: remove and fix ruff ignores --- examples/wave/var-propagation-speed.py | 5 +++-- examples/wave/wave-min-mpi.py | 5 ++--- examples/wave/wave-op-mpi.py | 15 +++++++-------- examples/wave/wave-op-var-velocity.py | 8 +++----- pyproject.toml | 8 +------- test/test_euler_model.py | 2 +- test/test_grudge.py | 10 +++++----- test/test_op.py | 8 ++++---- 8 files changed, 26 insertions(+), 35 deletions(-) diff --git a/examples/wave/var-propagation-speed.py b/examples/wave/var-propagation-speed.py index 0fa11616c..48bd52762 100644 --- a/examples/wave/var-propagation-speed.py +++ b/examples/wave/var-propagation-speed.py @@ -141,8 +141,9 @@ def norm(u): step += 1 if step % 10 == 0: - logger.info(f"step: {step} t: {time()-t_last_step} " - f"L2: {norm(u=event.state_component[0])}") + logger.info("step: %d t: %.8e L2: %.8e", + step, time() - t_last_step, + norm(event.state_component[0])) if visualize: vis.write_vtk_file( f"fld-var-propogation-speed-{step:04d}.vtu", diff --git a/examples/wave/wave-min-mpi.py b/examples/wave/wave-min-mpi.py index c4e5c66a2..acdc47e4d 100644 --- a/examples/wave/wave-min-mpi.py +++ b/examples/wave/wave-min-mpi.py @@ -167,9 +167,8 @@ def norm(u): if step % 10 == 0: if comm.rank == 0: - logger.info(f"step: {step} " - f"t: {time()-t_last_step} " - f"L2: {l2norm}") + logger.info("step: %d t: %.8e L2: %.8e", + step, time() - t_last_step, l2norm) if visualize: vis.write_parallel_vtk_file( comm, diff --git a/examples/wave/wave-op-mpi.py b/examples/wave/wave-op-mpi.py index 367d0fff4..f72985608 100644 --- a/examples/wave/wave-op-mpi.py +++ b/examples/wave/wave-op-mpi.py @@ -281,8 +281,8 @@ def rhs(t, w): stop = time.time() if no_diagnostics: if comm.rank == 0: - logger.info(f"step: {istep} t: {t} " - f"wall: {stop-start} ") + logger.info("step: %d t: %.8e wall: %.8es", + istep, t, stop - start) else: l2norm = actx.to_numpy(op.norm(dcoll, fields.u, 2)) @@ -294,12 +294,11 @@ def rhs(t, w): nodalmax = actx.to_numpy(op.nodal_max(dcoll, "vol", fields.u)) nodalmin = actx.to_numpy(op.nodal_min(dcoll, "vol", fields.u)) if comm.rank == 0: - logger.info(f"step: {istep} t: {t} " - f"L2: {l2norm} " - f"Linf: {linfnorm} " - f"sol max: {nodalmax} " - f"sol min: {nodalmin} " - f"wall: {stop-start} ") + logger.info("step: %d t: %.8e L2: %.8e Linf: %.8e " + "sol max: %.8e sol min: %.8e wall: %.8e", + istep, t, l2norm, linfnorm, nodalmax, nodalmin, + stop - start) + if visualize: vis.write_parallel_vtk_file( comm, diff --git a/examples/wave/wave-op-var-velocity.py b/examples/wave/wave-op-var-velocity.py index 06456bcd4..051d13716 100644 --- a/examples/wave/wave-op-var-velocity.py +++ b/examples/wave/wave-op-var-velocity.py @@ -206,11 +206,9 @@ def rhs(t, w): linfnorm = actx.to_numpy(op.norm(dcoll, fields[0], np.inf)) nodalmax = actx.to_numpy(op.nodal_max(dcoll, "vol", fields[0])) nodalmin = actx.to_numpy(op.nodal_min(dcoll, "vol", fields[0])) - logger.info(f"step: {istep} t: {t} " - f"L2: {l2norm} " - f"Linf: {linfnorm} " - f"sol max: {nodalmax} " - f"sol min: {nodalmin}") + logger.info("step: %d t: %.8e L2: %.8e Linf: %.8e " + "sol max: %.8e sol min: %.8e", + istep, t, l2norm, linfnorm, nodalmax, nodalmin) if visualize: vis.write_vtk_file( f"fld-wave-eager-var-velocity-{istep:04d}.vtu", diff --git a/pyproject.toml b/pyproject.toml index 51c740632..ed9c9aea8 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -31,7 +31,7 @@ classifiers = [ "Topic :: Utilities", ] dependencies = [ - "arraycontext>=2021.1" + "arraycontext>=2021.1", "immutables>=0.19", "loopy>=2024.1", "meshmode>=2021.2", @@ -103,12 +103,6 @@ extend-ignore = [ "E402", # module level import not at the top of file ] -[tool.ruff.lint.per-file-ignores] -"test/test_grudge.py" = ["B023"] -"test/test_op.py" = ["B023"] -"test/test_euler_model.py" = ["B023"] -"examples/wave/*.py" = ["G004"] - [tool.ruff.lint.flake8-quotes] docstring-quotes = "double" inline-quotes = "double" diff --git a/test/test_euler_model.py b/test/test_euler_model.py index 5ef1c54dd..83d19a4e5 100644 --- a/test/test_euler_model.py +++ b/test/test_euler_model.py @@ -93,7 +93,7 @@ def test_euler_vortex_convergence(actx_factory, order): quadrature_tag=quad_tag ) - def rhs(t, q): + def rhs(t, q, euler_operator=euler_operator): return euler_operator.operator(t, q) compiled_rhs = actx.compile(rhs) diff --git a/test/test_grudge.py b/test/test_grudge.py index 9e78c0f71..103e445b0 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -810,7 +810,7 @@ def test_convergence_advec(actx_factory, mesh_name, mesh_pars, op_type, flux_typ def f(x): return actx.np.sin(10*x) - def u_analytic(x, t=0): + def u_analytic(x, t=0, v=v, norm_v=norm_v): return f(-v.dot(x)/norm_v + t*norm_v) from meshmode.mesh import BTAG_ALL @@ -824,7 +824,7 @@ def u_analytic(x, t=0): op_class = {"strong": StrongAdvectionOperator, "weak": WeakAdvectionOperator}[op_type] adv_operator = op_class(dcoll, v, - inflow_u=lambda t: u_analytic( + inflow_u=lambda t, dcoll=dcoll: u_analytic( actx.thaw(dcoll.nodes(dd=BTAG_ALL)), t=t ), @@ -833,7 +833,7 @@ def u_analytic(x, t=0): nodes = actx.thaw(dcoll.nodes()) u = u_analytic(nodes, t=0) - def rhs(t, u): + def rhs(t, u, adv_operator=adv_operator): return adv_operator.operator(t, u) compiled_rhs = actx.compile(rhs) @@ -938,7 +938,7 @@ def analytic_sol(x, t=0): ) maxwell_operator.check_bc_coverage(mesh) - def rhs(t, w): + def rhs(t, w, maxwell_operator=maxwell_operator): return maxwell_operator.operator(t, w) dt = actx.to_numpy(maxwell_operator.estimate_rk4_timestep(actx, dcoll)) @@ -1026,7 +1026,7 @@ def conv_test(descr, use_quad): nodes = actx.thaw(dcoll.nodes()) - def zero_inflow(dtag, t=0): + def zero_inflow(dtag, t=0, dcoll=dcoll): dd = dof_desc.DOFDesc(dtag, qtag) return dcoll.discr_from_dd(dd).zeros(actx) diff --git a/test/test_op.py b/test/test_op.py index 18faedc6c..1ad9ddb96 100644 --- a/test/test_op.py +++ b/test/test_op.py @@ -118,7 +118,7 @@ def vectorize_if_requested(vec): else: return vec - def get_flux(u_tpair): + def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd dd_allfaces = dd.with_domain_tag("all_faces") normal = geo.normal(actx, dcoll, dd) @@ -236,14 +236,14 @@ def test_divergence(actx_factory, form, dim, order, vectorize, nested, dcoll = DiscretizationCollection(actx, mesh, order=order) - def f(x): + def f(x, dcoll=dcoll): result = make_obj_array([dcoll.zeros(actx) + (i+1) for i in range(dim)]) for i in range(dim-1): result = result * actx.np.sin(np.pi*x[i]) result = result * actx.np.cos(np.pi/2*x[dim-1]) return result - def div_f(x): + def div_f(x, dcoll=dcoll): result = dcoll.zeros(actx) for i in range(dim-1): deriv = dcoll.zeros(actx) + (i+1) @@ -270,7 +270,7 @@ def div_f(x): else: u = f(x) - def get_flux(u_tpair): + def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd dd_allfaces = dd.with_dtag("all_faces") normal = geo.normal(actx, dcoll, dd) From aa4af9b626c32c051089cd8baadf2cbf3b44e830 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 11:02:22 +0300 Subject: [PATCH 04/16] test: use make_discretization_collection --- examples/advection/surface.py | 4 ++-- examples/advection/var-velocity.py | 4 ++-- examples/advection/weak.py | 4 ++-- examples/euler/acoustic_pulse.py | 4 ++-- examples/euler/vortex.py | 4 ++-- examples/geometry.py | 11 +++++----- examples/maxwell/cavities.py | 6 +++--- examples/wave/var-propagation-speed.py | 6 +++--- examples/wave/wave-op-var-velocity.py | 4 ++-- test/test_dt_utils.py | 12 +++++------ test/test_euler_model.py | 4 ++-- test/test_grudge.py | 28 +++++++++++++------------- test/test_metrics.py | 4 ++-- test/test_mpi_communication.py | 10 ++++----- test/test_op.py | 8 ++++---- test/test_reductions.py | 10 ++++----- test/test_trace_pair.py | 4 ++-- 17 files changed, 64 insertions(+), 63 deletions(-) diff --git a/examples/advection/surface.py b/examples/advection/surface.py index 554df7c85..34d59cb4a 100644 --- a/examples/advection/surface.py +++ b/examples/advection/surface.py @@ -158,9 +158,9 @@ def main(ctx_factory, dim=2, order=4, use_quad=False, visualize=False): discr_tag_to_group_factory[qtag] = \ QuadratureSimplexGroupFactory(order=4*order) - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory=discr_tag_to_group_factory ) diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index 9c2301068..211df3535 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -139,9 +139,9 @@ def main(ctx_factory, dim=2, order=4, use_quad=False, visualize=False, else: discr_tag_to_group_factory = {} - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, order=order, discr_tag_to_group_factory=discr_tag_to_group_factory ) diff --git a/examples/advection/weak.py b/examples/advection/weak.py index 270574fc4..ddd918322 100644 --- a/examples/advection/weak.py +++ b/examples/advection/weak.py @@ -131,9 +131,9 @@ def main(ctx_factory, dim=2, order=4, visualize=False): [np.linspace(-d/2, d/2, npoints) for _ in range(dim)], order=order) - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) # }}} diff --git a/examples/euler/acoustic_pulse.py b/examples/euler/acoustic_pulse.py index e126abe76..117399ddf 100644 --- a/examples/euler/acoustic_pulse.py +++ b/examples/euler/acoustic_pulse.py @@ -126,7 +126,7 @@ def run_acoustic_pulse(actx, default_simplex_group_factory, ) - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD exp_name = f"fld-acoustic-pulse-N{order}-K{resolution}" @@ -136,7 +136,7 @@ def run_acoustic_pulse(actx, else: quad_tag = None - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory={ DISCR_TAG_BASE: default_simplex_group_factory( diff --git a/examples/euler/vortex.py b/examples/euler/vortex.py index 07c23f31b..443526d6a 100644 --- a/examples/euler/vortex.py +++ b/examples/euler/vortex.py @@ -74,7 +74,7 @@ def run_vortex(actx, order=3, resolution=8, final_time=5, default_simplex_group_factory, ) - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD exp_name = f"fld-vortex-N{order}-K{resolution}-{flux_type}" @@ -85,7 +85,7 @@ def run_vortex(actx, order=3, resolution=8, final_time=5, else: quad_tag = None - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory={ DISCR_TAG_BASE: default_simplex_group_factory( diff --git a/examples/geometry.py b/examples/geometry.py index 63e66cf9d..5d8439122 100644 --- a/examples/geometry.py +++ b/examples/geometry.py @@ -32,8 +32,9 @@ from grudge.array_context import PyOpenCLArrayContext -from grudge import DiscretizationCollection, shortcuts -import grudge.geometry as geo +from grudge import shortcuts +from grudge import geometry +from grudge.discretization import make_discretization_collection def main(write_output=True): @@ -47,13 +48,13 @@ def main(write_output=True): from meshmode.mesh import BTAG_ALL from meshmode.mesh.generation import generate_warped_rect_mesh - mesh = generate_warped_rect_mesh(dim=2, order=4, nelements_side=6) - dcoll = DiscretizationCollection(actx, mesh, order=4) + mesh = generate_warped_rect_mesh(dim=2, order=4, nelements_side=6) + dcoll = make_discretization_collection(actx, mesh, order=4) nodes = actx.thaw(dcoll.nodes()) bdry_nodes = actx.thaw(dcoll.nodes(dd=BTAG_ALL)) - bdry_normals = geo.normal(actx, dcoll, dd=BTAG_ALL) + bdry_normals = geometry.normal(actx, dcoll, dd=BTAG_ALL) if write_output: vis = shortcuts.make_visualizer(dcoll) diff --git a/examples/maxwell/cavities.py b/examples/maxwell/cavities.py index 0376a0a00..c05272b2a 100644 --- a/examples/maxwell/cavities.py +++ b/examples/maxwell/cavities.py @@ -31,9 +31,9 @@ import pyopencl as cl import pyopencl.tools as cl_tools -import grudge.op as op -from grudge import DiscretizationCollection +from grudge import op from grudge.array_context import PyOpenCLArrayContext +from grudge.discretization import make_discretization_collection from grudge.models.em import get_rectangular_cavity_mode from grudge.shortcuts import set_up_rk4 @@ -56,7 +56,7 @@ def main(ctx_factory, dim=3, order=4, visualize=False): b=(1.0,)*dim, nelements_per_axis=(4,)*dim) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) if 0: epsilon0 = 8.8541878176e-12 # C**2 / (N m**2) diff --git a/examples/wave/var-propagation-speed.py b/examples/wave/var-propagation-speed.py index 48bd52762..3bbfe19db 100644 --- a/examples/wave/var-propagation-speed.py +++ b/examples/wave/var-propagation-speed.py @@ -32,9 +32,9 @@ import pyopencl.tools as cl_tools from pytools.obj_array import flat_obj_array -import grudge.op as op -from grudge import DiscretizationCollection +from grudge import op from grudge.array_context import PyOpenCLArrayContext +from grudge.discretization import make_discretization_collection from grudge.shortcuts import set_up_rk4 @@ -56,7 +56,7 @@ def main(ctx_factory, dim=2, order=4, visualize=False): b=(0.5,)*dim, nelements_per_axis=(20,)*dim) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) def source_f(actx, dcoll, t=0): source_center = np.array([0.1, 0.22, 0.33])[:dcoll.dim] diff --git a/examples/wave/wave-op-var-velocity.py b/examples/wave/wave-op-var-velocity.py index 051d13716..8b1b28a28 100644 --- a/examples/wave/wave-op-var-velocity.py +++ b/examples/wave/wave-op-var-velocity.py @@ -37,7 +37,7 @@ import grudge.geometry as geo import grudge.op as op from grudge.array_context import PyOpenCLArrayContext -from grudge.discretization import DiscretizationCollection +from grudge.discretization import make_discretization_collection from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD, DOFDesc from grudge.shortcuts import make_visualizer, rk4_step @@ -170,7 +170,7 @@ def main(ctx_factory, dim=2, order=3, visualize=False): QuadratureSimplexGroupFactory, default_simplex_group_factory, ) - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory={ DISCR_TAG_BASE: default_simplex_group_factory(base_dim=dim, order=order), diff --git a/test/test_dt_utils.py b/test/test_dt_utils.py index 5d99e5adb..729980917 100644 --- a/test/test_dt_utils.py +++ b/test/test_dt_utils.py @@ -42,7 +42,7 @@ import pytest import grudge.op as op -from grudge import DiscretizationCollection +from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) @@ -73,7 +73,7 @@ def test_geometric_factors_regular_refinement(actx_factory, name): min_factors = [] for resolution in builder.resolutions: mesh = builder.get_mesh(resolution, order) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) min_factors.append( actx.to_numpy( op.nodal_min(dcoll, "vol", actx.thaw(dt_geometric_factors(dcoll)))) @@ -87,7 +87,7 @@ def test_geometric_factors_regular_refinement(actx_factory, name): # Make sure it works with empty meshes mesh = builder.get_mesh(0) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) factors = actx.thaw(dt_geometric_factors(dcoll)) # noqa: F841 @@ -114,7 +114,7 @@ def test_non_geometric_factors(actx_factory, name): degrees = list(range(1, 8)) for degree in degrees: mesh = builder.get_mesh(1, degree) - dcoll = DiscretizationCollection(actx, mesh, order=degree) + dcoll = make_discretization_collection(actx, mesh, order=degree) factors.append(min(dt_non_geometric_factors(dcoll))) # Crude estimate, factors should behave like 1/N**2 @@ -133,7 +133,7 @@ def test_build_jacobian(actx_factory): mesh = mgen.generate_regular_rect_mesh(a=[0], b=[1], nelements_per_axis=(3,)) assert mesh.dim == 1 - dcoll = DiscretizationCollection(actx, mesh, order=1) + dcoll = make_discretization_collection(actx, mesh, order=1) def rhs(x): return 3*x**2 + 2*x + 5 @@ -162,7 +162,7 @@ def test_wave_dt_estimate(actx_factory, dim, degree, visualize=False): nelements_per_axis=(3,)*dim) assert mesh.dim == dim - dcoll = DiscretizationCollection(actx, mesh, order=degree) + dcoll = make_discretization_collection(actx, mesh, order=degree) from grudge.models.wave import WeakWaveOperator wave_op = WeakWaveOperator(dcoll, c=1) diff --git a/test/test_euler_model.py b/test/test_euler_model.py index 83d19a4e5..1540ef142 100644 --- a/test/test_euler_model.py +++ b/test/test_euler_model.py @@ -52,7 +52,7 @@ def test_euler_vortex_convergence(actx_factory, order): from meshmode.mesh.generation import generate_regular_rect_mesh from pytools.convergence import EOCRecorder - from grudge import DiscretizationCollection + from grudge.discretization import make_discretization_collection from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD from grudge.dt_utils import h_max_from_volume from grudge.models.euler import EulerOperator, vortex_initial_condition @@ -77,7 +77,7 @@ def test_euler_vortex_convergence(actx_factory, order): DISCR_TAG_QUAD: QuadratureSimplexGroupFactory(2*order) } - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory=discr_tag_to_group_factory ) diff --git a/test/test_grudge.py b/test/test_grudge.py index 103e445b0..3202435cd 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -52,7 +52,7 @@ import grudge.dof_desc as dof_desc import grudge.geometry as geo import grudge.op as op -from grudge import DiscretizationCollection, make_discretization_collection +from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) @@ -88,7 +88,7 @@ def test_mass_mat_trig(actx_factory, ambient_dim, discr_tag): mesh = mgen.generate_regular_rect_mesh( a=(a,)*ambient_dim, b=(b,)*ambient_dim, nelements_per_axis=(nel_1d,)*ambient_dim, order=1) - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, order=order, discr_tag_to_group_factory=discr_tag_to_group_factory ) @@ -193,7 +193,7 @@ def test_mass_surface_area(actx_factory, name): for resolution in builder.resolutions: mesh = builder.get_mesh(resolution, order) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) logger.info("ndofs: %d", volume_discr.ndofs) @@ -352,7 +352,7 @@ def test_face_normal_surface(actx_factory, mesh_name): order = 4 mesh = builder.get_mesh(builder.resolutions[0], order) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) logger.info("ndofs: %d", volume_discr.ndofs) @@ -443,7 +443,7 @@ def df(x, axis): mesh = mgen.generate_regular_rect_mesh(a=(-0.5,)*dim, b=(0.5,)*dim, nelements_per_axis=(n,)*dim, order=4) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) x = actx.thaw(volume_discr.nodes()) @@ -657,7 +657,7 @@ def f(x): from meshmode.discretization.poly_element import QuadratureSimplexGroupFactory qtag = dof_desc.DISCR_TAG_QUAD - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, order=order, discr_tag_to_group_factory={ qtag: QuadratureSimplexGroupFactory(2 * order) @@ -820,7 +820,7 @@ def u_analytic(x, t=0, v=v, norm_v=norm_v): WeakAdvectionOperator, ) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) op_class = {"strong": StrongAdvectionOperator, "weak": WeakAdvectionOperator}[op_type] adv_operator = op_class(dcoll, v, @@ -914,7 +914,7 @@ def test_convergence_maxwell(actx_factory, order): b=(1.0,)*dims, nelements_per_axis=(n,)*dims) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) epsilon = 1 mu = 1 @@ -1019,7 +1019,7 @@ def conv_test(descr, use_quad): else: discr_tag_to_group_factory = {} - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, order=order, discr_tag_to_group_factory=discr_tag_to_group_factory ) @@ -1072,7 +1072,7 @@ def test_bessel(actx_factory): b=(1.0,)*dims, nelements_per_axis=(8,)*dims) - dcoll = DiscretizationCollection(actx, mesh, order=3) + dcoll = make_discretization_collection(actx, mesh, order=3) nodes = actx.thaw(dcoll.nodes()) r = actx.np.sqrt(nodes[0]**2 + nodes[1]**2) @@ -1106,7 +1106,7 @@ def test_norm_real(actx_factory, p): mesh = mgen.generate_regular_rect_mesh( a=(0,)*dim, b=(1,)*dim, nelements_per_axis=(8,)*dim, order=1) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) nodes = actx.thaw(dcoll.nodes()) norm = op.norm(dcoll, nodes[0], p) @@ -1129,7 +1129,7 @@ def test_norm_complex(actx_factory, p): mesh = mgen.generate_regular_rect_mesh( a=(0,)*dim, b=(1,)*dim, nelements_per_axis=(8,)*dim, order=1) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) nodes = actx.thaw(dcoll.nodes()) norm = op.norm(dcoll, (1 + 1j)*nodes[0], p) @@ -1152,7 +1152,7 @@ def test_norm_obj_array(actx_factory, p): mesh = mgen.generate_regular_rect_mesh( a=(0,)*dim, b=(1,)*dim, nelements_per_axis=(8,)*dim, order=1) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) nodes = actx.thaw(dcoll.nodes()) norm = op.norm(dcoll, nodes, p) @@ -1183,7 +1183,7 @@ def test_empty_boundary(actx_factory): mesh = mgen.generate_regular_rect_mesh( a=(-0.5,)*dim, b=(0.5,)*dim, nelements_per_axis=(8,)*dim, order=4) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) normal = geo.normal(actx, dcoll, BTAG_NONE) from meshmode.dof_array import DOFArray for component in normal: diff --git a/test/test_metrics.py b/test/test_metrics.py index f26778e29..ec766342d 100644 --- a/test/test_metrics.py +++ b/test/test_metrics.py @@ -44,7 +44,7 @@ import meshmode.mesh.generation as mgen from meshmode.dof_array import flat_norm -from grudge import DiscretizationCollection +from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) @@ -85,7 +85,7 @@ def m(x): from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD, as_dofdesc - dcoll = DiscretizationCollection( + dcoll = make_discretization_collection( actx, mesh, discr_tag_to_group_factory={ DISCR_TAG_BASE: default_simplex_group_factory(base_dim=dim, order=order), diff --git a/test/test_mpi_communication.py b/test/test_mpi_communication.py index 9998a1b0a..f1ce00739 100644 --- a/test/test_mpi_communication.py +++ b/test/test_mpi_communication.py @@ -31,8 +31,6 @@ import numpy as np import pytest -import pyopencl as cl - from grudge.array_context import MPIPyOpenCLArrayContext, MPIPytatoArrayContext @@ -45,7 +43,7 @@ import grudge.dof_desc as dof_desc import grudge.op as op -from grudge import DiscretizationCollection +from grudge.discretization import make_discretization_collection from grudge.shortcuts import rk4_step @@ -81,6 +79,8 @@ def run_test_with_mpi_inner(): from pickle import loads f, (actx_class, *args) = loads(b64decode(os.environ["INVOCATION_INFO"].encode())) + import pyopencl as cl + cl_context = cl.create_some_context() queue = cl.CommandQueue(cl_context) @@ -136,7 +136,7 @@ def _test_func_comparison_mpi_communication_entrypoint(actx): else: local_mesh = comm.scatter(None) - dcoll = DiscretizationCollection(actx, local_mesh, order=5) + dcoll = make_discretization_collection(actx, local_mesh, order=5) x = actx.thaw(dcoll.nodes()) myfunc = actx.np.sin(np.dot(x, [2, 3])) @@ -214,7 +214,7 @@ def _test_mpi_wave_op_entrypoint(actx, visualize=False): else: local_mesh = comm.scatter(None) - dcoll = DiscretizationCollection(actx, local_mesh, order=order) + dcoll = make_discretization_collection(actx, local_mesh, order=order) def source_f(actx, dcoll, t=0): source_center = np.array([0.1, 0.22, 0.33])[:dcoll.dim] diff --git a/test/test_op.py b/test/test_op.py index 1ad9ddb96..79890cedb 100644 --- a/test/test_op.py +++ b/test/test_op.py @@ -33,7 +33,7 @@ from meshmode.mesh import BTAG_ALL from pytools.obj_array import make_obj_array -from grudge import DiscretizationCollection, geometry as geo, op +from grudge import geometry, op from grudge.array_context import PytestPyOpenCLArrayContextFactory from grudge.discretization import make_discretization_collection from grudge.dof_desc import ( @@ -121,7 +121,7 @@ def vectorize_if_requested(vec): def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd dd_allfaces = dd.with_domain_tag("all_faces") - normal = geo.normal(actx, dcoll, dd) + normal = geometry.normal(actx, dcoll, dd) u_avg = u_tpair.avg if vectorize: if nested: @@ -234,7 +234,7 @@ def test_divergence(actx_factory, form, dim, order, vectorize, nested, a=(-1,)*dim, b=(1,)*dim, nelements_per_axis=(n,)*dim) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) def f(x, dcoll=dcoll): result = make_obj_array([dcoll.zeros(actx) + (i+1) for i in range(dim)]) @@ -273,7 +273,7 @@ def div_f(x, dcoll=dcoll): def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd dd_allfaces = dd.with_dtag("all_faces") - normal = geo.normal(actx, dcoll, dd) + normal = geometry.normal(actx, dcoll, dd) flux = u_tpair.avg @ normal return op.project(dcoll, dd, dd_allfaces, flux) diff --git a/test/test_reductions.py b/test/test_reductions.py index 33f98b417..ac3e65d3e 100644 --- a/test/test_reductions.py +++ b/test/test_reductions.py @@ -48,7 +48,7 @@ from pytools.obj_array import make_obj_array import grudge.op as op -from grudge import DiscretizationCollection +from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) @@ -66,7 +66,7 @@ def test_nodal_reductions(actx_factory, mesh_size, with_initial): builder = mesh_data.BoxMeshBuilder1D() mesh = builder.get_mesh(mesh_size) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) x = actx.thaw(dcoll.nodes()) def f(x): @@ -132,7 +132,7 @@ def test_elementwise_reductions(actx_factory): nelements = 4 mesh = builder.get_mesh(nelements) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) x = actx.thaw(dcoll.nodes()) def f(x): @@ -192,7 +192,7 @@ def test_nodal_reductions_with_container(actx_factory): builder = mesh_data.BoxMeshBuilder2D() mesh = builder.get_mesh(4) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) x = actx.thaw(dcoll.nodes()) def f(x): @@ -239,7 +239,7 @@ def test_elementwise_reductions_with_container(actx_factory): nelements = 4 mesh = builder.get_mesh(nelements) - dcoll = DiscretizationCollection(actx, mesh, order=4) + dcoll = make_discretization_collection(actx, mesh, order=4) x = actx.thaw(dcoll.nodes()) def f(x): diff --git a/test/test_trace_pair.py b/test/test_trace_pair.py index 5bfbb0eb9..43dc6abfb 100644 --- a/test/test_trace_pair.py +++ b/test/test_trace_pair.py @@ -27,8 +27,8 @@ from arraycontext import pytest_generate_tests_for_array_contexts from meshmode.dof_array import DOFArray -from grudge import DiscretizationCollection from grudge.array_context import PytestPyOpenCLArrayContextFactory +from grudge.discretization import make_discretization_collection from grudge.trace_pair import TracePair @@ -52,7 +52,7 @@ def test_trace_pair(actx_factory): a=(-1,)*dim, b=(1,)*dim, nelements_per_axis=(n,)*dim) - dcoll = DiscretizationCollection(actx, mesh, order=order) + dcoll = make_discretization_collection(actx, mesh, order=order) rng = np.random.default_rng(1234) From 1aceaff08d75a664f099a406ba33bb1614f60ad8 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 11:09:34 +0300 Subject: [PATCH 05/16] test: fix import order --- test/test_euler_model.py | 11 +++---- test/test_grudge.py | 59 ++++++++++++++++------------------ test/test_metrics.py | 19 ++++------- test/test_mpi_communication.py | 22 +++++-------- test/test_op.py | 8 ++--- test/test_reductions.py | 23 +++++-------- test/test_trace_pair.py | 8 ++--- 7 files changed, 61 insertions(+), 89 deletions(-) diff --git a/test/test_euler_model.py b/test/test_euler_model.py index 1540ef142..5da3a0e7c 100644 --- a/test/test_euler_model.py +++ b/test/test_euler_model.py @@ -22,25 +22,22 @@ THE SOFTWARE. """ +import logging + import pytest from arraycontext import ( pytest_generate_tests_for_array_contexts, ) +from grudge import op from grudge.array_context import PytestPyOpenCLArrayContextFactory +logger = logging.getLogger(__name__) pytest_generate_tests = pytest_generate_tests_for_array_contexts( [PytestPyOpenCLArrayContextFactory]) -import logging - -import grudge.op as op - - -logger = logging.getLogger(__name__) - @pytest.mark.parametrize("order", [1, 2, 3]) def test_euler_vortex_convergence(actx_factory, order): diff --git a/test/test_grudge.py b/test/test_grudge.py index 3202435cd..0d8b6def7 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -23,39 +23,32 @@ THE SOFTWARE. """ +import logging + import mesh_data import numpy as np import numpy.linalg as la +import pytest +import meshmode.mesh.generation as mgen from arraycontext import pytest_generate_tests_for_array_contexts +from meshmode import _acf # noqa: F401 from meshmode.discretization.poly_element import ( InterpolatoryEdgeClusteredGroupFactory, QuadratureGroupFactory, ) -from meshmode.mesh import TensorProductElementGroup - -from grudge.array_context import PytestPyOpenCLArrayContextFactory - - -pytest_generate_tests = pytest_generate_tests_for_array_contexts( - [PytestPyOpenCLArrayContextFactory]) - -import logging - -import pytest - -import meshmode.mesh.generation as mgen -from meshmode import _acf # noqa: F401 from meshmode.dof_array import flat_norm +from meshmode.mesh import TensorProductElementGroup from pytools.obj_array import flat_obj_array -import grudge.dof_desc as dof_desc -import grudge.geometry as geo -import grudge.op as op +from grudge import dof_desc, geometry, op +from grudge.array_context import PytestPyOpenCLArrayContextFactory from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) +pytest_generate_tests = pytest_generate_tests_for_array_contexts( + [PytestPyOpenCLArrayContextFactory]) # {{{ mass operator trig integration @@ -77,6 +70,7 @@ def test_mass_mat_trig(actx_factory, ambient_dim, discr_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, discr_tag) if discr_tag is dof_desc.DISCR_TAG_BASE: discr_tag_to_group_factory = {} @@ -141,6 +135,7 @@ def _ellipse_surface_area(radius, aspect_ratio): ellip_e = 1.2110560275684594 else: from scipy.special import ellipe # pylint: disable=no-name-in-module + ellip_e = ellipe(eccentricity) return 4.0 * radius * ellip_e @@ -189,6 +184,7 @@ def test_mass_surface_area(actx_factory, name): # {{{ convergence from pytools.convergence import EOCRecorder + eoc = EOCRecorder() for resolution in builder.resolutions: @@ -274,6 +270,7 @@ def test_mass_operator_inverse(actx_factory, name): # {{{ inv(m) @ m == id from pytools.convergence import EOCRecorder + eoc = EOCRecorder() for resolution in builder.resolutions: @@ -363,8 +360,6 @@ def test_face_normal_surface(actx_factory, mesh_name): # {{{ Compute surface and face normals from meshmode.discretization.connection import FACE_RESTR_INTERIOR - from grudge.geometry import normal - dv = dof_desc.DD_VOLUME df = dof_desc.as_dofdesc(FACE_RESTR_INTERIOR) @@ -372,21 +367,21 @@ def test_face_normal_surface(actx_factory, mesh_name): surf_normal = op.project( dcoll, dv, df, - normal(actx, dcoll, dd=dv) + geometry.normal(actx, dcoll, dd=dv) ) surf_normal = surf_normal / actx.np.sqrt(sum(surf_normal**2)) - face_normal_i = geo.normal(actx, dcoll, df) + face_normal_i = geometry.normal(actx, dcoll, df) face_normal_e = dcoll.opposite_face_connection( dof_desc.BoundaryDomainTag( dof_desc.FACE_RESTR_INTERIOR, dof_desc.VTAG_ALL) )(face_normal_i) if ambient_dim == 3: - from grudge.geometry import area_element, pseudoscalar # NOTE: there's only one face tangent in 3d face_tangent = ( - pseudoscalar(actx, dcoll, dd=df) / area_element(actx, dcoll, dd=df) + geometry.pseudoscalar(actx, dcoll, dd=df) + / geometry.area_element(actx, dcoll, dd=df) ).as_vector(dtype=object) # }}} @@ -546,7 +541,7 @@ def f(x): int_1 = op.integral(dcoll, "vol", div_f) prj_f = op.project(dcoll, "vol", BTAG_ALL, f_volm) - normal = geo.normal(actx, dcoll, BTAG_ALL) + normal = geometry.normal(actx, dcoll, BTAG_ALL) f_dot_n = prj_f.dot(normal) int_2 = op.integral(dcoll, BTAG_ALL, f_dot_n) else: @@ -562,12 +557,13 @@ def f(x): int_1 = op.integral(dcoll, dd_quad_vol, div_f) prj_f = op.project(dcoll, "vol", dd_quad_bd, f_volm) - normal = geo.normal(actx, dcoll, dd_quad_bd) + normal = geometry.normal(actx, dcoll, dd_quad_bd) f_dot_n = prj_f.dot(normal) int_2 = op.integral(dcoll, dd_quad_bd, f_dot_n) if visualize: from grudge.shortcuts import make_boundary_visualizer, make_visualizer + vis = make_visualizer(dcoll) bvis = make_boundary_visualizer(dcoll) @@ -677,11 +673,9 @@ def f(x): f_num = f(actx.thaw(dcoll.nodes(dd=dd))) f_quad_num = f(actx.thaw(dcoll.nodes(dd=dq))) - from grudge.geometry import normal, summed_curvature - - kappa = summed_curvature(actx, dcoll, dd=dq) - normal = normal(actx, dcoll, dd=dq) - face_normal = geo.normal(actx, dcoll, df) + kappa = geometry.summed_curvature(actx, dcoll, dd=dq) + normal = geometry.normal(actx, dcoll, dd=dq) + face_normal = geometry.normal(actx, dcoll, df) face_f = op.project(dcoll, dd, df, f_num) # operators @@ -710,6 +704,7 @@ def f(x): if visualize: from grudge.shortcuts import make_visualizer + vis = make_visualizer(dcoll) filename = f"surface_divergence_theorem_{mesh_name}_{i:04d}.vtu" @@ -852,6 +847,7 @@ def rhs(t, u, adv_operator=adv_operator): dt = final_time/nsteps + tol from grudge.shortcuts import compiled_lsrk45_step, make_visualizer + vis = make_visualizer(dcoll) step = 0 @@ -946,6 +942,7 @@ def rhs(t, w, maxwell_operator=maxwell_operator): nsteps = int(final_t/dt) from grudge.shortcuts import set_up_rk4 + dt_stepper = set_up_rk4("w", dt, fields, rhs) logger.info("dt %.5e nsteps %5d", dt, nsteps) @@ -1184,7 +1181,7 @@ def test_empty_boundary(actx_factory): a=(-0.5,)*dim, b=(0.5,)*dim, nelements_per_axis=(8,)*dim, order=4) dcoll = make_discretization_collection(actx, mesh, order=4) - normal = geo.normal(actx, dcoll, BTAG_NONE) + normal = geometry.normal(actx, dcoll, BTAG_NONE) from meshmode.dof_array import DOFArray for component in normal: assert isinstance(component, DOFArray) diff --git a/test/test_metrics.py b/test/test_metrics.py index ec766342d..e6863c7f6 100644 --- a/test/test_metrics.py +++ b/test/test_metrics.py @@ -23,32 +23,27 @@ THE SOFTWARE. """ +import logging + import numpy as np +import pytest +import meshmode.mesh.generation as mgen from arraycontext import pytest_generate_tests_for_array_contexts +from meshmode.dof_array import flat_norm from grudge.array_context import ( PytestPyOpenCLArrayContextFactory, PytestPytatoPyOpenCLArrayContextFactory, ) +from grudge.discretization import make_discretization_collection +logger = logging.getLogger(__name__) pytest_generate_tests = pytest_generate_tests_for_array_contexts( [PytestPyOpenCLArrayContextFactory, PytestPytatoPyOpenCLArrayContextFactory]) -import logging - -import pytest - -import meshmode.mesh.generation as mgen -from meshmode.dof_array import flat_norm - -from grudge.discretization import make_discretization_collection - - -logger = logging.getLogger(__name__) - # {{{ inverse metric diff --git a/test/test_mpi_communication.py b/test/test_mpi_communication.py index f1ce00739..342c52676 100644 --- a/test/test_mpi_communication.py +++ b/test/test_mpi_communication.py @@ -31,22 +31,18 @@ import numpy as np import pytest -from grudge.array_context import MPIPyOpenCLArrayContext, MPIPytatoArrayContext - - -logger = logging.getLogger(__name__) -logging.basicConfig() -logger.setLevel(logging.INFO) - from meshmode.dof_array import flat_norm from pytools.obj_array import flat_obj_array -import grudge.dof_desc as dof_desc -import grudge.op as op +from grudge import dof_desc, op +from grudge.array_context import MPIPyOpenCLArrayContext, MPIPytatoArrayContext from grudge.discretization import make_discretization_collection from grudge.shortcuts import rk4_step +logger = logging.getLogger(__name__) + + class SimpleTag: pass @@ -141,11 +137,9 @@ def _test_func_comparison_mpi_communication_entrypoint(actx): x = actx.thaw(dcoll.nodes()) myfunc = actx.np.sin(np.dot(x, [2, 3])) - from grudge.dof_desc import as_dofdesc - - dd_int = as_dofdesc("int_faces") - dd_vol = as_dofdesc("vol") - dd_af = as_dofdesc("all_faces") + dd_int = dof_desc.as_dofdesc("int_faces") + dd_vol = dof_desc.as_dofdesc("vol") + dd_af = dof_desc.as_dofdesc("all_faces") all_faces_func = op.project(dcoll, dd_vol, dd_af, myfunc) int_faces_func = op.project(dcoll, dd_vol, dd_int, myfunc) diff --git a/test/test_op.py b/test/test_op.py index 79890cedb..125c69570 100644 --- a/test/test_op.py +++ b/test/test_op.py @@ -21,6 +21,8 @@ """ +import logging + import numpy as np import pytest @@ -46,14 +48,10 @@ from grudge.trace_pair import bv_trace_pair +logger = logging.getLogger(__name__) pytest_generate_tests = pytest_generate_tests_for_array_contexts( [PytestPyOpenCLArrayContextFactory]) -import logging - - -logger = logging.getLogger(__name__) - # {{{ gradient diff --git a/test/test_reductions.py b/test/test_reductions.py index ac3e65d3e..f0ef12dc5 100644 --- a/test/test_reductions.py +++ b/test/test_reductions.py @@ -22,36 +22,29 @@ THE SOFTWARE. """ +import logging from dataclasses import dataclass +import mesh_data import numpy as np +import pytest from arraycontext import ( dataclass_array_container, pytest_generate_tests_for_array_contexts, with_container_arithmetic, ) -from meshmode.dof_array import DOFArray - -from grudge.array_context import PytestPyOpenCLArrayContextFactory - - -pytest_generate_tests = pytest_generate_tests_for_array_contexts( - [PytestPyOpenCLArrayContextFactory]) - -import logging - -import mesh_data -import pytest - -from meshmode.dof_array import flatten +from meshmode.dof_array import DOFArray, flatten from pytools.obj_array import make_obj_array -import grudge.op as op +from grudge import op +from grudge.array_context import PytestPyOpenCLArrayContextFactory from grudge.discretization import make_discretization_collection logger = logging.getLogger(__name__) +pytest_generate_tests = pytest_generate_tests_for_array_contexts( + [PytestPyOpenCLArrayContextFactory]) @pytest.mark.parametrize(("mesh_size", "with_initial"), [ diff --git a/test/test_trace_pair.py b/test/test_trace_pair.py index 43dc6abfb..d9c9fa789 100644 --- a/test/test_trace_pair.py +++ b/test/test_trace_pair.py @@ -21,6 +21,8 @@ """ +import logging + import numpy as np import meshmode.mesh.generation as mgen @@ -32,14 +34,10 @@ from grudge.trace_pair import TracePair +logger = logging.getLogger(__name__) pytest_generate_tests = pytest_generate_tests_for_array_contexts( [PytestPyOpenCLArrayContextFactory]) -import logging - - -logger = logging.getLogger(__name__) - def test_trace_pair(actx_factory): """Simple smoke test for :class:`grudge.trace_pair.TracePair`.""" From 03fbe19881f27a0a6c03119270fd8c7f4da64924 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 11:13:54 +0300 Subject: [PATCH 06/16] test: replace deprecated DD_VOLUME --- examples/advection/surface.py | 6 +++--- examples/advection/var-velocity.py | 2 +- examples/advection/weak.py | 2 +- test/test_grudge.py | 24 ++++++++++++------------ 4 files changed, 17 insertions(+), 17 deletions(-) diff --git a/examples/advection/surface.py b/examples/advection/surface.py index 34d59cb4a..c1134de9c 100644 --- a/examples/advection/surface.py +++ b/examples/advection/surface.py @@ -61,7 +61,7 @@ def __init__(self, actx, dcoll, order, visualize=True): import matplotlib.pyplot as pt self.fig = pt.figure(figsize=(8, 8), dpi=300) - x = actx.thaw(dcoll.discr_from_dd(dof_desc.DD_VOLUME).nodes()) + x = actx.thaw(dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL).nodes()) self.x = actx.to_numpy(flatten(actx.np.arctan2(x[1], x[0]))) elif self.ambient_dim == 3: from grudge.shortcuts import make_visualizer @@ -165,7 +165,7 @@ def main(ctx_factory, dim=2, order=4, use_quad=False, visualize=False): discr_tag_to_group_factory=discr_tag_to_group_factory ) - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) logger.info("ndofs: %d", volume_discr.ndofs) logger.info("nelements: %d", volume_discr.mesh.nelements) @@ -196,7 +196,7 @@ def rhs(t, u): # check velocity is tangential from grudge.geometry import normal - surf_normal = normal(actx, dcoll, dd=dof_desc.DD_VOLUME) + surf_normal = normal(actx, dcoll, dd=dof_desc.DD_VOLUME_ALL) error = op.norm(dcoll, c.dot(surf_normal), 2) logger.info("u_dot_n: %.5e", error) diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index 211df3535..217df74d5 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -58,7 +58,7 @@ def __init__(self, actx, dcoll, order, visualize=True, ylim=None): self.fig = pt.figure(figsize=(8, 8), dpi=300) self.ylim = ylim - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) self.x = actx.to_numpy(flatten(actx.thaw(volume_discr.nodes()[0]))) else: from grudge.shortcuts import make_visualizer diff --git a/examples/advection/weak.py b/examples/advection/weak.py index ddd918322..204ef39f4 100644 --- a/examples/advection/weak.py +++ b/examples/advection/weak.py @@ -58,7 +58,7 @@ def __init__(self, actx, dcoll, order, visualize=True, ylim=None): self.fig = pt.figure(figsize=(8, 8), dpi=300) self.ylim = ylim - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) self.x = actx.to_numpy(flatten(actx.thaw(volume_discr.nodes()[0]))) else: from grudge.shortcuts import make_visualizer diff --git a/test/test_grudge.py b/test/test_grudge.py index 0d8b6def7..38f4af4fd 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -90,7 +90,7 @@ def test_mass_mat_trig(actx_factory, ambient_dim, discr_tag): def f(x): return actx.np.sin(x[0])**2 - volm_disc = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volm_disc = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) x_volm = actx.thaw(volm_disc.nodes()) f_volm = f(x_volm) ones_volm = volm_disc.zeros(actx) + 1 @@ -102,14 +102,14 @@ def f(x): mop_1 = op.mass(dcoll, dd_quad, f_quad) num_integral_1 = op.nodal_sum( - dcoll, dof_desc.DD_VOLUME, ones_volm * mop_1 + dcoll, dof_desc.DD_VOLUME_ALL, ones_volm * mop_1 ) err_1 = abs(num_integral_1 - true_integral) assert err_1 < 2e-9, err_1 mop_2 = op.mass(dcoll, dd_quad, ones_quad) - num_integral_2 = op.nodal_sum(dcoll, dof_desc.DD_VOLUME, f_volm * mop_2) + num_integral_2 = op.nodal_sum(dcoll, dof_desc.DD_VOLUME_ALL, f_volm * mop_2) err_2 = abs(num_integral_2 - true_integral) assert err_2 < 2e-9, err_2 @@ -117,7 +117,7 @@ def f(x): 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 = op.nodal_sum(dcoll, dof_desc.DD_VOLUME, f_quad * mop_2) + num_integral_3 = op.nodal_sum(dcoll, dof_desc.DD_VOLUME_ALL, f_quad * mop_2) err_3 = abs(num_integral_3 - true_integral) assert err_3 < 5e-10, err_3 @@ -190,14 +190,14 @@ def test_mass_surface_area(actx_factory, name): for resolution in builder.resolutions: mesh = builder.get_mesh(resolution, order) dcoll = make_discretization_collection(actx, mesh, order=order) - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) logger.info("ndofs: %d", volume_discr.ndofs) logger.info("nelements: %d", volume_discr.mesh.nelements) # {{{ compute surface area - dd = dof_desc.DD_VOLUME + dd = dof_desc.DD_VOLUME_ALL ones_volm = volume_discr.zeros(actx) + 1 approx_surface_area = actx.to_numpy(op.integral(dcoll, dd, ones_volm)) @@ -295,7 +295,7 @@ def f(x): x_volm = actx.thaw(volume_discr.nodes()) f_volm = f(x_volm) if not overintegrate: - dd = dof_desc.DD_VOLUME + dd = dof_desc.DD_VOLUME_ALL f_inv = op.inverse_mass( dcoll, op.mass(dcoll, dd, f_volm) ) @@ -351,7 +351,7 @@ def test_face_normal_surface(actx_factory, mesh_name): mesh = builder.get_mesh(builder.resolutions[0], order) dcoll = make_discretization_collection(actx, mesh, order=order) - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) logger.info("ndofs: %d", volume_discr.ndofs) logger.info("nelements: %d", volume_discr.mesh.nelements) @@ -360,7 +360,7 @@ def test_face_normal_surface(actx_factory, mesh_name): # {{{ Compute surface and face normals from meshmode.discretization.connection import FACE_RESTR_INTERIOR - dv = dof_desc.DD_VOLUME + dv = dof_desc.DD_VOLUME_ALL df = dof_desc.as_dofdesc(FACE_RESTR_INTERIOR) ambient_dim = mesh.ambient_dim @@ -439,7 +439,7 @@ def df(x, axis): nelements_per_axis=(n,)*dim, order=4) dcoll = make_discretization_collection(actx, mesh, order=4) - volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) x = actx.thaw(volume_discr.nodes()) for axis in range(dim): @@ -660,11 +660,11 @@ def f(x): } ) - volume = dcoll.discr_from_dd(dof_desc.DD_VOLUME) + volume = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) logger.info("ndofs: %d", volume.ndofs) logger.info("nelements: %d", volume.mesh.nelements) - dd = dof_desc.DD_VOLUME + dd = dof_desc.DD_VOLUME_ALL dq = dd.with_discr_tag(qtag) df = dof_desc.as_dofdesc(FACE_RESTR_ALL) ambient_dim = dcoll.ambient_dim From 4e6478965c1f361f4564baca3542b1be7b7460ad Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 14:56:43 +0300 Subject: [PATCH 07/16] replace sloppy uses of DOFDesc --- examples/advection/surface.py | 2 +- examples/advection/var-velocity.py | 2 +- examples/wave/wave-op-var-velocity.py | 6 +++--- grudge/dof_desc.py | 5 +++-- grudge/models/advection.py | 21 ++++++++++++--------- grudge/models/em.py | 7 +++++-- grudge/models/euler.py | 16 ++++++++++------ test/test_grudge.py | 2 +- test/test_op.py | 16 +++++++++++----- 9 files changed, 47 insertions(+), 30 deletions(-) diff --git a/examples/advection/surface.py b/examples/advection/surface.py index c1134de9c..05deb364c 100644 --- a/examples/advection/surface.py +++ b/examples/advection/surface.py @@ -236,7 +236,7 @@ def rhs(t, u): overwrite=True ) - df = dof_desc.DOFDesc(FACE_RESTR_INTERIOR) + df = dof_desc.as_dofdesc(FACE_RESTR_INTERIOR) face_discr = dcoll.discr_from_dd(df) face_normal = geo.normal(actx, dcoll, dd=df) diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index 217df74d5..9c73a1af3 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -160,7 +160,7 @@ def f_halfcircle(x): * (0.5+0.5*actx.np.tanh(500*(dist[0])))) def zero_inflow_bc(dtag, t=0): - dd = dof_desc.DOFDesc(dtag, qtag) + dd = dof_desc.as_dofdesc(dtag, qtag) return dcoll.discr_from_dd(dd).zeros(actx) from grudge.models.advection import VariableCoefficientAdvectionOperator diff --git a/examples/wave/wave-op-var-velocity.py b/examples/wave/wave-op-var-velocity.py index 8b1b28a28..486608a39 100644 --- a/examples/wave/wave-op-var-velocity.py +++ b/examples/wave/wave-op-var-velocity.py @@ -38,7 +38,7 @@ import grudge.op as op from grudge.array_context import PyOpenCLArrayContext from grudge.discretization import make_discretization_collection -from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD, DOFDesc +from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD, as_dofdesc from grudge.shortcuts import make_visualizer, rk4_step @@ -84,13 +84,13 @@ def wave_operator(actx, 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", DISCR_TAG_QUAD) + dd_quad = as_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", DISCR_TAG_QUAD) + dd_allfaces_quad = as_dofdesc("all_faces", DISCR_TAG_QUAD) return ( op.inverse_mass( diff --git a/grudge/dof_desc.py b/grudge/dof_desc.py index da1f4a41b..27c031b37 100644 --- a/grudge/dof_desc.py +++ b/grudge/dof_desc.py @@ -230,8 +230,9 @@ class DOFDesc: domain_tag: DomainTag discretization_tag: DiscretizationTag - def __init__(self, domain_tag: Any, - discretization_tag: Optional[type[DiscretizationTag]] = None): + def __init__(self, + domain_tag: Any, + discretization_tag: Optional[Type[DiscretizationTag]] = None): if ( not (isinstance(domain_tag, diff --git a/grudge/models/advection.py b/grudge/models/advection.py index c5cba7798..8bf2260d3 100644 --- a/grudge/models/advection.py +++ b/grudge/models/advection.py @@ -106,8 +106,10 @@ def flux(tpair): return op.project(dcoll, tpair.dd, "all_faces", self.flux(tpair)) if self.inflow_u is not None: + from grudge.dof_desc import as_dofdesc + inflow_flux = flux(op.bv_trace_pair(dcoll, - BTAG_ALL, + as_dofdesc(BTAG_ALL), interior=u, exterior=self.inflow_u(t))) else: @@ -149,8 +151,9 @@ def flux(tpair): return op.project(dcoll, tpair.dd, "all_faces", self.flux(tpair)) if self.inflow_u is not None: + from grudge.dof_desc import as_dofdesc inflow_flux = flux(op.bv_trace_pair(dcoll, - BTAG_ALL, + as_dofdesc(BTAG_ALL), interior=u, exterior=self.inflow_u(t))) else: @@ -225,11 +228,11 @@ def operator(self, t, u): from meshmode.discretization.connection import FACE_RESTR_ALL from meshmode.mesh import BTAG_ALL - from grudge.dof_desc import DD_VOLUME_ALL, DTAG_VOLUME_ALL, DOFDesc + from grudge.dof_desc import DD_VOLUME_ALL, DTAG_VOLUME_ALL, as_dofdesc - face_dd = DOFDesc(FACE_RESTR_ALL, self.quad_tag) - boundary_dd = DOFDesc(BTAG_ALL, self.quad_tag) - quad_dd = DOFDesc(DTAG_VOLUME_ALL, self.quad_tag) + face_dd = as_dofdesc(FACE_RESTR_ALL, self.quad_tag) + boundary_dd = as_dofdesc(BTAG_ALL, self.quad_tag) + quad_dd = as_dofdesc(DTAG_VOLUME_ALL, self.quad_tag) dcoll = self.dcoll @@ -340,10 +343,10 @@ def flux(self, u_tpair): def operator(self, t, u): from meshmode.discretization.connection import FACE_RESTR_ALL - from grudge.dof_desc import DD_VOLUME_ALL, DTAG_VOLUME_ALL, DOFDesc + from grudge.dof_desc import DD_VOLUME_ALL, DTAG_VOLUME_ALL, as_dofdesc - face_dd = DOFDesc(FACE_RESTR_ALL, self.quad_tag) - quad_dd = DOFDesc(DTAG_VOLUME_ALL, self.quad_tag) + face_dd = as_dofdesc(FACE_RESTR_ALL, self.quad_tag) + quad_dd = as_dofdesc(DTAG_VOLUME_ALL, self.quad_tag) dcoll = self.dcoll diff --git a/grudge/models/em.py b/grudge/models/em.py index c7c5a7256..4e32c9ba8 100644 --- a/grudge/models/em.py +++ b/grudge/models/em.py @@ -410,6 +410,8 @@ def operator(self, t, w): def flux(pair): return op.project(dcoll, pair.dd, "all_faces", self.flux(pair)) + from grudge.dof_desc import as_dofdesc + return ( - self.local_derivatives(w) - op.inverse_mass( @@ -417,7 +419,7 @@ def flux(pair): op.face_mass( dcoll, sum(flux(tpair) for tpair in op.interior_trace_pairs(dcoll, w)) - + sum(flux(op.bv_trace_pair(dcoll, tag, w, bc)) + + sum(flux(op.bv_trace_pair(dcoll, as_dofdesc(tag), w, bc)) for tag, bc in tags_and_bcs) ) ) @@ -462,7 +464,8 @@ def check_bc_coverage(self, mesh): self.pec_tag, self.pmc_tag, self.absorb_tag, - self.incident_tag]) + self.incident_tag, + ]) # }}} diff --git a/grudge/models/euler.py b/grudge/models/euler.py index a288b486e..8d182986c 100644 --- a/grudge/models/euler.py +++ b/grudge/models/euler.py @@ -180,7 +180,7 @@ def boundary_tpair( dd_bc: DOFDesc, state: ConservedEulerField, t=0): actx = state.array_context - dd_base = as_dofdesc("vol").with_discr_tag(DISCR_TAG_BASE) + dd_base = as_dofdesc("vol", DISCR_TAG_BASE) return TracePair( dd_bc, @@ -197,7 +197,7 @@ def boundary_tpair( dd_bc: DOFDesc, state: ConservedEulerField, t=0): actx = state.array_context - dd_base = as_dofdesc("vol").with_discr_tag(DISCR_TAG_BASE) + dd_base = as_dofdesc("vol", DISCR_TAG_BASE) nhat = geo.normal(actx, dcoll, dd_bc) interior = op.project(dcoll, dd_base, dd_bc, state) @@ -253,8 +253,12 @@ def euler_numerical_flux( dissipation. :returns: A :class:`ConservedEulerField` containing the interface fluxes. """ + from grudge.dof_desc import FACE_RESTR_ALL, VTAG_ALL, BoundaryDomainTag + dd_intfaces = tpair.dd - dd_allfaces = dd_intfaces.with_dtag("all_faces") + dd_allfaces = dd_intfaces.with_domain_tag( + BoundaryDomainTag(FACE_RESTR_ALL, VTAG_ALL) + ) q_ll = tpair.int q_rr = tpair.ext actx = q_ll.array_context @@ -310,8 +314,8 @@ def operator(self, t, q): dcoll = self.dcoll gamma = self.gamma qtag = self.qtag - dq = DOFDesc("vol", qtag) - df = DOFDesc("all_faces", qtag) + dq = as_dofdesc("vol", qtag) + df = as_dofdesc("all_faces", qtag) def interp_to_quad(u): return op.project(dcoll, "vol", dq, u) @@ -341,7 +345,7 @@ def interp_to_quad(u): dcoll, self.bdry_conditions[btag].boundary_tpair( dcoll, - as_dofdesc(btag).with_discr_tag(qtag), + as_dofdesc(btag, qtag), q, t=t ), diff --git a/test/test_grudge.py b/test/test_grudge.py index 38f4af4fd..8a5a5a678 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -1024,7 +1024,7 @@ def conv_test(descr, use_quad): nodes = actx.thaw(dcoll.nodes()) def zero_inflow(dtag, t=0, dcoll=dcoll): - dd = dof_desc.DOFDesc(dtag, qtag) + dd = dof_desc.as_dofdesc(dtag, qtag) return dcoll.discr_from_dd(dd).zeros(actx) adv_op = VariableCoefficientAdvectionOperator( diff --git a/test/test_op.py b/test/test_op.py index 125c69570..17f49a074 100644 --- a/test/test_op.py +++ b/test/test_op.py @@ -42,7 +42,9 @@ DISCR_TAG_BASE, DISCR_TAG_QUAD, DTAG_VOLUME_ALL, - DOFDesc, + FACE_RESTR_ALL, + VTAG_ALL, + BoundaryDomainTag, as_dofdesc, ) from grudge.trace_pair import bv_trace_pair @@ -118,7 +120,9 @@ def vectorize_if_requested(vec): def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd - dd_allfaces = dd.with_domain_tag("all_faces") + dd_allfaces = dd.with_domain_tag( + BoundaryDomainTag(FACE_RESTR_ALL, VTAG_ALL) + ) normal = geometry.normal(actx, dcoll, dd) u_avg = u_tpair.avg if vectorize: @@ -149,7 +153,7 @@ def get_flux(u_tpair, dcoll=dcoll): else: quad_discr_tag = DISCR_TAG_BASE - allfaces_dd_base = as_dofdesc("all_faces", quad_discr_tag) + allfaces_dd_base = as_dofdesc(FACE_RESTR_ALL, quad_discr_tag) vol_dd_base = as_dofdesc(DTAG_VOLUME_ALL) vol_dd_quad = vol_dd_base.with_discr_tag(quad_discr_tag) bdry_dd_quad = bdry_dd_base.with_discr_tag(quad_discr_tag) @@ -270,12 +274,14 @@ def div_f(x, dcoll=dcoll): def get_flux(u_tpair, dcoll=dcoll): dd = u_tpair.dd - dd_allfaces = dd.with_dtag("all_faces") + dd_allfaces = dd.with_domain_tag( + BoundaryDomainTag(FACE_RESTR_ALL, VTAG_ALL) + ) normal = geometry.normal(actx, dcoll, dd) flux = u_tpair.avg @ normal return op.project(dcoll, dd, dd_allfaces, flux) - dd_allfaces = DOFDesc("all_faces") + dd_allfaces = as_dofdesc(FACE_RESTR_ALL) if form == "strong": div_u = ( From 31e6690c7f32f0bf4eb6b687a1817361829d73c7 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 15:07:13 +0300 Subject: [PATCH 08/16] euler: add _cls_has_array_context_attr=True --- grudge/models/euler.py | 11 +++++++++-- 1 file changed, 9 insertions(+), 2 deletions(-) diff --git a/grudge/models/euler.py b/grudge/models/euler.py index 8d182986c..894683d3d 100644 --- a/grudge/models/euler.py +++ b/grudge/models/euler.py @@ -66,7 +66,10 @@ @with_container_arithmetic(bcast_obj_array=False, bcast_container_types=(DOFArray, np.ndarray), matmul=True, - rel_comparison=True) + rel_comparison=True, + # NOTE: this can be set to True because DOFArray + # also sets it to True and we just get it from there + _cls_has_array_context_attr=True) @dataclass_array_container @dataclass(frozen=True) class ConservedEulerField: @@ -76,7 +79,11 @@ class ConservedEulerField: @property def array_context(self): - return self.mass.array_context + return ( + self.mass.array_context + if isinstance(self.mass, DOFArray) + # NOTE: ConservedEulerField also holds the fluxes sometimes... + else self.mass[0].array_context) @property def dim(self): From d2714508885092ea5246346ea718c8c69649f6f5 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 15:08:27 +0300 Subject: [PATCH 09/16] op: use multi_vandermonde for gradients --- grudge/op.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/grudge/op.py b/grudge/op.py index 5052d61a7..0d76c180a 100644 --- a/grudge/op.py +++ b/grudge/op.py @@ -448,10 +448,11 @@ def get_ref_stiffness_transpose_mat(out_grp, in_grp): np.asarray( [dmat.T @ mmat.T for dmat in diff_matrices(out_grp)])))) - from modepy import vandermonde + from modepy import multi_vandermonde, vandermonde + basis = out_grp.basis_obj() vand = vandermonde(basis.functions, out_grp.unit_nodes) - grad_vand = vandermonde(basis.gradients, in_grp.unit_nodes) + grad_vand = multi_vandermonde(basis.gradients, in_grp.unit_nodes) vand_inv_t = np.linalg.inv(vand).T if not isinstance(grad_vand, tuple): From 5c6d752592a1a403326022fde8493004aab7af5f Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 15:09:23 +0300 Subject: [PATCH 10/16] replace deprecated with_dtag --- examples/wave/wave-op-var-velocity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/wave/wave-op-var-velocity.py b/examples/wave/wave-op-var-velocity.py index 486608a39..00099398a 100644 --- a/examples/wave/wave-op-var-velocity.py +++ b/examples/wave/wave-op-var-velocity.py @@ -68,7 +68,7 @@ def wave_flux(actx, dcoll, c, w_tpair): ) # FIXME this flux is only correct for continuous c - dd_allfaces_quad = dd_quad.with_dtag("all_faces") + dd_allfaces_quad = dd_quad.with_domain_tag("all_faces") c_quad = op.project(dcoll, "vol", dd_quad, c) flux_quad = op.project(dcoll, dd, dd_quad, flux_weak) From 8e0794e2ecca18704f6a9172e512eedff96c4937 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 15:18:21 +0300 Subject: [PATCH 11/16] port from deprecated meshmode flatten --- examples/advection/surface.py | 6 +++--- examples/advection/var-velocity.py | 6 +++--- examples/advection/weak.py | 6 +++--- test/test_reductions.py | 18 ++++++++++-------- 4 files changed, 19 insertions(+), 17 deletions(-) diff --git a/examples/advection/surface.py b/examples/advection/surface.py index 05deb364c..026bf9947 100644 --- a/examples/advection/surface.py +++ b/examples/advection/surface.py @@ -32,8 +32,8 @@ import pyopencl as cl import pyopencl.tools as cl_tools +from arraycontext import flatten from meshmode.discretization.connection import FACE_RESTR_INTERIOR -from meshmode.dof_array import flatten from pytools.obj_array import make_obj_array import grudge.dof_desc as dof_desc @@ -62,7 +62,7 @@ def __init__(self, actx, dcoll, order, visualize=True): self.fig = pt.figure(figsize=(8, 8), dpi=300) x = actx.thaw(dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL).nodes()) - self.x = actx.to_numpy(flatten(actx.np.arctan2(x[1], x[0]))) + self.x = actx.to_numpy(flatten(actx.np.arctan2(x[1], x[0]), self.actx)) elif self.ambient_dim == 3: from grudge.shortcuts import make_visualizer self.vis = make_visualizer(dcoll) @@ -74,7 +74,7 @@ def __call__(self, evt, basename, overwrite=True): return if self.ambient_dim == 2: - u = self.actx.to_numpy(flatten(evt.state_component)) + u = self.actx.to_numpy(flatten(evt.state_component, self.actx)) filename = f"{basename}.png" if not overwrite and os.path.exists(filename): diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index 9c73a1af3..dd266280a 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -30,7 +30,7 @@ import pyopencl as cl import pyopencl.tools as cl_tools -from meshmode.dof_array import flatten +from arraycontext import flatten from meshmode.mesh import BTAG_ALL from pytools.obj_array import flat_obj_array @@ -59,7 +59,7 @@ def __init__(self, actx, dcoll, order, visualize=True, ylim=None): self.ylim = ylim volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) - self.x = actx.to_numpy(flatten(actx.thaw(volume_discr.nodes()[0]))) + self.x = actx.to_numpy(flatten(volume_discr.nodes()[0], self.actx)) else: from grudge.shortcuts import make_visualizer self.vis = make_visualizer(dcoll) @@ -69,7 +69,7 @@ def __call__(self, evt, basename, overwrite=True): return if self.dim == 1: - u = self.actx.to_numpy(flatten(evt.state_component)) + u = self.actx.to_numpy(flatten(evt.state_component, self.actx)) filename = f"{basename}.png" if not overwrite and os.path.exists(filename): diff --git a/examples/advection/weak.py b/examples/advection/weak.py index 204ef39f4..eee9999da 100644 --- a/examples/advection/weak.py +++ b/examples/advection/weak.py @@ -31,7 +31,7 @@ import pyopencl as cl import pyopencl.tools as cl_tools -from meshmode.dof_array import flatten +from arraycontext import flatten from meshmode.mesh import BTAG_ALL import grudge.dof_desc as dof_desc @@ -59,7 +59,7 @@ def __init__(self, actx, dcoll, order, visualize=True, ylim=None): self.ylim = ylim volume_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME_ALL) - self.x = actx.to_numpy(flatten(actx.thaw(volume_discr.nodes()[0]))) + self.x = actx.to_numpy(flatten(volume_discr.nodes()[0], self.actx)) else: from grudge.shortcuts import make_visualizer self.vis = make_visualizer(dcoll) @@ -69,7 +69,7 @@ def __call__(self, evt, basename, overwrite=True): return if self.dim == 1: - u = self.actx.to_numpy(flatten(evt.state_component)) + u = self.actx.to_numpy(flatten(evt.state_component, self.actx)) filename = f"{basename}.png" if not overwrite and os.path.exists(filename): diff --git a/test/test_reductions.py b/test/test_reductions.py index f0ef12dc5..44850a54e 100644 --- a/test/test_reductions.py +++ b/test/test_reductions.py @@ -31,10 +31,11 @@ from arraycontext import ( dataclass_array_container, + flatten, pytest_generate_tests_for_array_contexts, with_container_arithmetic, ) -from meshmode.dof_array import DOFArray, flatten +from meshmode.dof_array import DOFArray from pytools.obj_array import make_obj_array from grudge import op @@ -73,9 +74,9 @@ def h(x): fields = make_obj_array([f(x), g(x), h(x)]) - f_ref = actx.to_numpy(flatten(fields[0])) - g_ref = actx.to_numpy(flatten(fields[1])) - h_ref = actx.to_numpy(flatten(fields[2])) + f_ref = actx.to_numpy(flatten(fields[0], actx)) + g_ref = actx.to_numpy(flatten(fields[1], actx)) + h_ref = actx.to_numpy(flatten(fields[2], actx)) concat_fields = np.concatenate([f_ref, g_ref, h_ref]) for grudge_op, np_op in [(op.nodal_max, np.max), @@ -206,10 +207,11 @@ def h(x): momentum=momentum, enthalpy=enthalpy) - mass_ref = actx.to_numpy(flatten(mass)) - momentum_ref = np.concatenate([actx.to_numpy(mom_i) - for mom_i in flatten(momentum)]) - enthalpy_ref = actx.to_numpy(flatten(enthalpy)) + mass_ref = actx.to_numpy(flatten(mass, actx)) + momentum_ref = np.concatenate([ + actx.to_numpy(mom_i) + for mom_i in flatten(momentum, actx, leaf_class=DOFArray)]) + enthalpy_ref = actx.to_numpy(flatten(enthalpy, actx)) concat_fields = np.concatenate([mass_ref, momentum_ref, enthalpy_ref]) for grudge_op, np_op in [(op.nodal_sum, np.sum), From b617e19160f4388f4c05d6a903499a9cd0d9ceba Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 11 Jul 2024 15:20:23 +0300 Subject: [PATCH 12/16] tests: fix deprecated call to TracePair --- test/test_trace_pair.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/test/test_trace_pair.py b/test/test_trace_pair.py index d9c9fa789..5cb060c7e 100644 --- a/test/test_trace_pair.py +++ b/test/test_trace_pair.py @@ -61,9 +61,10 @@ def rand(): rng.uniform(size=(grp.nelements, grp.nunit_dofs))) for grp in dcoll.discr_from_dd("vol").groups)) + from grudge.dof_desc import DD_VOLUME_ALL interior = rand() exterior = rand() - tpair = TracePair("vol", interior=interior, exterior=exterior) + tpair = TracePair(DD_VOLUME_ALL, interior=interior, exterior=exterior) import grudge.op as op assert op.norm(dcoll, tpair.avg - 0.5*(exterior + interior), np.inf) == 0 From fb971ab9ef855d364d8c0417686bc1b21a425e5b Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Sat, 13 Jul 2024 10:36:40 +0300 Subject: [PATCH 13/16] port to get_connected_parts --- grudge/discretization.py | 4 ++-- grudge/trace_pair.py | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/grudge/discretization.py b/grudge/discretization.py index 7d51e5145..3e66f283c 100644 --- a/grudge/discretization.py +++ b/grudge/discretization.py @@ -256,8 +256,8 @@ def _set_up_distributed_communication( boundary_connections = {} - from meshmode.distributed import get_connected_partitions - connected_parts = get_connected_partitions(self._volume_discrs[vtag].mesh) + from meshmode.distributed import get_connected_parts + connected_parts = get_connected_parts(self._volume_discrs[vtag].mesh) if connected_parts: if mpi_communicator is None: diff --git a/grudge/trace_pair.py b/grudge/trace_pair.py index 6480cadb7..6ac3ad632 100644 --- a/grudge/trace_pair.py +++ b/grudge/trace_pair.py @@ -361,8 +361,8 @@ def connected_ranks( if volume_dd is None: volume_dd = DD_VOLUME_ALL - from meshmode.distributed import get_connected_partitions - return get_connected_partitions( + from meshmode.distributed import get_connected_parts + return get_connected_parts( dcoll._volume_discrs[volume_dd.domain_tag.tag].mesh) From eba596dac40643ef2fdbca125042483ef1b328d9 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Sat, 13 Jul 2024 11:09:38 +0300 Subject: [PATCH 14/16] test: use MemoryPool in testing actx --- grudge/array_context.py | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/grudge/array_context.py b/grudge/array_context.py index 6fc44bf66..6c31e6727 100644 --- a/grudge/array_context.py +++ b/grudge/array_context.py @@ -164,6 +164,7 @@ def __init__(self, queue, allocator=None, warn("No memory allocator specified, please pass one. " "(Preferably a pyopencl.tools.MemoryPool in order " "to reduce device allocations)", stacklevel=2) + super().__init__(queue, allocator, compile_trace_callback=compile_trace_callback) @@ -509,11 +510,30 @@ class PytestPyOpenCLArrayContextFactory( _PytestPyOpenCLArrayContextFactoryWithClass): actx_class = PyOpenCLArrayContext + def __call__(self): + from pyopencl.tools import ImmediateAllocator, MemoryPool + + _ctx, queue = self.get_command_queue() + alloc = MemoryPool(ImmediateAllocator(queue)) + + return self.actx_class( + queue, + allocator=alloc, + force_device_scalars=self.force_device_scalars) + class PytestPytatoPyOpenCLArrayContextFactory( _PytestPytatoPyOpenCLArrayContextFactory): actx_class = PytatoPyOpenCLArrayContext + def __call__(self): + _ctx, queue = self.get_command_queue() + + from pyopencl.tools import ImmediateAllocator, MemoryPool + alloc = MemoryPool(ImmediateAllocator(queue)) + + return self.actx_class(queue, allocator=alloc) + # deprecated class PytestPyOpenCLArrayContextFactoryWithHostScalars( From d0921572c7f7288850cf00d45db8ccaf713e1b8f Mon Sep 17 00:00:00 2001 From: Andreas Kloeckner Date: Tue, 16 Jul 2024 12:20:08 -0500 Subject: [PATCH 15/16] Switch ConservedEulerField to not knowing its actx --- examples/euler/acoustic_pulse.py | 12 +++++---- examples/euler/vortex.py | 2 +- grudge/models/euler.py | 45 ++++++++++++++------------------ test/test_euler_model.py | 2 +- 4 files changed, 29 insertions(+), 32 deletions(-) diff --git a/examples/euler/acoustic_pulse.py b/examples/euler/acoustic_pulse.py index 117399ddf..c44a6322e 100644 --- a/examples/euler/acoustic_pulse.py +++ b/examples/euler/acoustic_pulse.py @@ -29,6 +29,7 @@ import pyopencl as cl import pyopencl.tools as cl_tools +from arraycontext import ArrayContext from meshmode.mesh import BTAG_ALL from pytools.obj_array import make_obj_array @@ -42,6 +43,7 @@ def gaussian_profile( + actx: ArrayContext, x_vec, t=0, rho0=1.0, rhoamp=1.0, p0=1.0, gamma=1.4, center=None, velocity=None): @@ -57,12 +59,11 @@ def gaussian_profile( rel_center = make_obj_array( [x_vec[i] - lump_loc[i] for i in range(dim)] ) - actx = x_vec[0].array_context r = actx.np.sqrt(np.dot(rel_center, rel_center)) expterm = rhoamp * actx.np.exp(1 - r ** 2) mass = expterm + rho0 - mom = velocity * mass + mom = velocity.astype(object) * mass energy = (p0 / (gamma - 1.0)) + np.dot(mom, mom) / (2.0 * mass) return ConservedEulerField(mass=mass, energy=energy, momentum=mom) @@ -81,11 +82,12 @@ def make_pulse(amplitude, r0, w, r): return amplitude * actx.np.exp(-.5 * r2) -def acoustic_pulse_condition(x_vec, t=0): +def acoustic_pulse_condition(actx: ArrayContext, x_vec, t=0): dim = len(x_vec) vel = np.zeros(shape=(dim,)) orig = np.zeros(shape=(dim,)) uniform_gaussian = gaussian_profile( + actx, x_vec, t=t, center=orig, velocity=vel, rhoamp=0.0) amplitude = 1.0 @@ -158,7 +160,7 @@ def run_acoustic_pulse(actx, ) def rhs(t, q): - return euler_operator.operator(t, q) + return euler_operator.operator(actx, t, q) compiled_rhs = actx.compile(rhs) @@ -168,7 +170,7 @@ def rhs(t, q): cn = 0.5*(order + 1)**2 dt = cfl * actx.to_numpy(h_min_from_volume(dcoll)) / cn - fields = acoustic_pulse_condition(actx.thaw(dcoll.nodes())) + fields = acoustic_pulse_condition(actx, actx.thaw(dcoll.nodes())) logger.info("Timestep size: %g", dt) diff --git a/examples/euler/vortex.py b/examples/euler/vortex.py index 443526d6a..974981efb 100644 --- a/examples/euler/vortex.py +++ b/examples/euler/vortex.py @@ -106,7 +106,7 @@ def run_vortex(actx, order=3, resolution=8, final_time=5, ) def rhs(t, q): - return euler_operator.operator(t, q) + return euler_operator.operator(actx, t, q) compiled_rhs = actx.compile(rhs) diff --git a/grudge/models/euler.py b/grudge/models/euler.py index 894683d3d..ef55ad63e 100644 --- a/grudge/models/euler.py +++ b/grudge/models/euler.py @@ -49,7 +49,11 @@ import numpy as np -from arraycontext import dataclass_array_container, with_container_arithmetic +from arraycontext import ( + ArrayContext, + dataclass_array_container, + with_container_arithmetic, +) from meshmode.dof_array import DOFArray from pytools.obj_array import make_obj_array @@ -67,9 +71,7 @@ bcast_container_types=(DOFArray, np.ndarray), matmul=True, rel_comparison=True, - # NOTE: this can be set to True because DOFArray - # also sets it to True and we just get it from there - _cls_has_array_context_attr=True) + ) @dataclass_array_container @dataclass(frozen=True) class ConservedEulerField: @@ -77,14 +79,6 @@ class ConservedEulerField: energy: DOFArray momentum: np.ndarray - @property - def array_context(self): - return ( - self.mass.array_context - if isinstance(self.mass, DOFArray) - # NOTE: ConservedEulerField also holds the fluxes sometimes... - else self.mass[0].array_context) - @property def dim(self): return len(self.momentum) @@ -146,7 +140,7 @@ def conservative_to_primitive_vars(cv_state: ConservedEulerField, gamma=1.4): return rho, u, p -def compute_wavespeed(cv_state: ConservedEulerField, gamma=1.4): +def compute_wavespeed(actx: ArrayContext, cv_state: ConservedEulerField, gamma=1.4): """Computes the total translational wavespeed. :arg cv_state: A :class:`ConservedEulerField` containing the conserved @@ -155,7 +149,6 @@ def compute_wavespeed(cv_state: ConservedEulerField, gamma=1.4): (default set to 1.4). :returns: A :class:`~meshmode.dof_array.DOFArray` containing local wavespeeds. """ - actx = cv_state.array_context rho, u, p = conservative_to_primitive_vars(cv_state, gamma=gamma) return actx.np.sqrt(np.dot(u, u)) + actx.np.sqrt(gamma * (p / rho)) @@ -173,6 +166,7 @@ def __init__(self, *, prescribed_state=None) -> None: @abstractmethod def boundary_tpair( self, + actx: ArrayContext, dcoll: DiscretizationCollection, dd_bc: DOFDesc, state: ConservedEulerField, t=0): @@ -183,6 +177,7 @@ class PrescribedBC(InviscidBCObject): def boundary_tpair( self, + actx: ArrayContext, dcoll: DiscretizationCollection, dd_bc: DOFDesc, state: ConservedEulerField, t=0): @@ -200,10 +195,10 @@ class InviscidWallBC(InviscidBCObject): def boundary_tpair( self, + actx: ArrayContext, dcoll: DiscretizationCollection, dd_bc: DOFDesc, state: ConservedEulerField, t=0): - actx = state.array_context dd_base = as_dofdesc("vol", DISCR_TAG_BASE) nhat = geo.normal(actx, dcoll, dd_bc) interior = op.project(dcoll, dd_base, dd_bc, state) @@ -243,11 +238,12 @@ def euler_volume_flux( return ConservedEulerField( mass=cv_state.momentum, energy=u * (cv_state.energy + p), - momentum=rho * outer(u, u) + np.eye(dcoll.dim) * p + momentum=rho * outer(u, u) + np.eye(dcoll.dim, dtype=object) * p ) def euler_numerical_flux( + actx: ArrayContext, dcoll: DiscretizationCollection, tpair: TracePair, gamma=1.4, lf_stabilization=False): """Computes the interface numerical flux for the Euler operator. @@ -268,7 +264,6 @@ def euler_numerical_flux( ) q_ll = tpair.int q_rr = tpair.ext - actx = q_ll.array_context flux_tpair = TracePair( tpair.dd, @@ -282,8 +277,8 @@ def euler_numerical_flux( from arraycontext import outer # Compute jump penalization parameter - lam = actx.np.maximum(compute_wavespeed(q_ll, gamma=gamma), - compute_wavespeed(q_rr, gamma=gamma)) + lam = actx.np.maximum(compute_wavespeed(actx, q_ll, gamma=gamma), + compute_wavespeed(actx, q_rr, gamma=gamma)) num_flux -= lam*outer(tpair.diff, normal)/2 return op.project(dcoll, dd_intfaces, dd_allfaces, num_flux @ normal) @@ -313,11 +308,11 @@ def __init__(self, dcoll: DiscretizationCollection, self.lf_stabilization = flux_type == "lf" self.qtag = quadrature_tag - def max_characteristic_velocity(self, actx, **kwargs): + def max_characteristic_velocity(self, actx: ArrayContext, **kwargs): state = kwargs["state"] - return compute_wavespeed(state, gamma=self.gamma) + return compute_wavespeed(actx, state, gamma=self.gamma) - def operator(self, t, q): + def operator(self, actx: ArrayContext, t, q): dcoll = self.dcoll gamma = self.gamma qtag = self.qtag @@ -337,7 +332,7 @@ def interp_to_quad(u): interface_fluxes = ( sum( euler_numerical_flux( - dcoll, + actx, dcoll, op.tracepair_with_discr_tag(dcoll, qtag, tpair), gamma=gamma, lf_stabilization=self.lf_stabilization @@ -349,9 +344,9 @@ def interp_to_quad(u): if self.bdry_conditions is not None: bc_fluxes = sum( euler_numerical_flux( - dcoll, + actx, dcoll, self.bdry_conditions[btag].boundary_tpair( - dcoll, + actx, dcoll, as_dofdesc(btag, qtag), q, t=t diff --git a/test/test_euler_model.py b/test/test_euler_model.py index 5da3a0e7c..6e9e9c1d2 100644 --- a/test/test_euler_model.py +++ b/test/test_euler_model.py @@ -91,7 +91,7 @@ def test_euler_vortex_convergence(actx_factory, order): ) def rhs(t, q, euler_operator=euler_operator): - return euler_operator.operator(t, q) + return euler_operator.operator(actx, t, q) compiled_rhs = actx.compile(rhs) From f44e6d25f3a159ce328f6fab5933c3b6ac084549 Mon Sep 17 00:00:00 2001 From: Andreas Kloeckner Date: Tue, 16 Jul 2024 12:20:47 -0500 Subject: [PATCH 16/16] ref_stiffness_transpose: Use matrix makers from modepy --- grudge/op.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/grudge/op.py b/grudge/op.py index 0d76c180a..f9e762f14 100644 --- a/grudge/op.py +++ b/grudge/op.py @@ -439,14 +439,13 @@ def _reference_stiffness_transpose_matrices( in_grp.discretization_key())) def get_ref_stiffness_transpose_mat(out_grp, in_grp): if in_grp == out_grp: - from meshmode.discretization.poly_element import diff_matrices, mass_matrix - - mmat = mass_matrix(out_grp) + mmat = mp.mass_matrix(out_grp.basis_obj(), out_grp.unit_nodes) + diff_matrices = mp.diff_matrices(out_grp.basis_obj(), out_grp.unit_nodes) return actx.freeze( actx.tag_axis(1, DiscretizationDOFAxisTag(), actx.from_numpy( np.asarray( - [dmat.T @ mmat.T for dmat in diff_matrices(out_grp)])))) + [dmat.T @ mmat.T for dmat in diff_matrices])))) from modepy import multi_vandermonde, vandermonde