From bee34848c248b0a15b96665e8521eb067266f57a Mon Sep 17 00:00:00 2001
From: Andreas Kloeckner <inform@tiker.net>
Date: Sat, 4 Sep 2021 22:37:09 -0500
Subject: [PATCH] inverse_surface_metric_derivative_mat: add, use
times_area_element keyword arg
---
grudge/geometry/metrics.py | 16 +++++++++++++---
grudge/op.py | 32 ++++++++++++++------------------
2 files changed, 27 insertions(+), 21 deletions(-)
diff --git a/grudge/geometry/metrics.py b/grudge/geometry/metrics.py
index 1ce02ce2..24c5f8a2 100644
--- a/grudge/geometry/metrics.py
+++ b/grudge/geometry/metrics.py
@@ -414,7 +414,8 @@ def inverse_surface_metric_derivative(
def inverse_surface_metric_derivative_mat(
- actx: ArrayContext, dcoll: DiscretizationCollection, dd=None):
+ actx: ArrayContext, dcoll: DiscretizationCollection, dd=None,
+ *, times_area_element=False):
r"""Computes the matrix of inverse surface metric derivatives, indexed by
``(xyz_axis, rst_axis)``. It returns all values of
:func:`inverse_surface_metric_derivative_mat` in cached matrix form.
@@ -423,16 +424,25 @@ def inverse_surface_metric_derivative_mat(
:arg dd: a :class:`~grudge.dof_desc.DOFDesc`, or a value convertible to one.
Defaults to the base volume discretization.
+ :arg times_area_element: If *True*, each entry of the matrix is premultiplied
+ with the value of :func:`area_element`, reflecting the typical use
+ of the matrix in integrals evaluating weak derivatives.
:returns: a :class:`~meshmode.dof_array.DOFArray` containing the
inverse metric derivatives in per-group arrays of shape
``(xyz_dimension, rst_dimension, nelements, ndof)``.
"""
- @memoize_in(dcoll, (inverse_surface_metric_derivative_mat, dd))
+ @memoize_in(dcoll, (inverse_surface_metric_derivative_mat, dd,
+ times_area_element))
def _inv_surf_metric_deriv():
+ if times_area_element:
+ multiplier = area_element(actx, dcoll, dd=dd)
+ else:
+ multiplier = 1
+
mat = actx.np.stack([
actx.np.stack(
- [inverse_surface_metric_derivative(actx, dcoll,
+ [multiplier * inverse_surface_metric_derivative(actx, dcoll,
rst_axis, xyz_axis, dd=dd)
for rst_axis in range(dcoll.dim)])
for xyz_axis in range(dcoll.ambient_dim)])
diff --git a/grudge/op.py b/grudge/op.py
index 5176a710..4aac9881 100644
--- a/grudge/op.py
+++ b/grudge/op.py
@@ -306,32 +306,30 @@ def weak_local_grad(dcoll: DiscretizationCollection, *args, nested=False):
else:
return np.stack(grad, axis=0)
- from grudge.geometry import \
- inverse_surface_metric_derivative_mat, area_element
+ from grudge.geometry import inverse_surface_metric_derivative_mat
in_discr = dcoll.discr_from_dd(dd_in)
out_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME)
actx = vec.array_context
- area_elements = area_element(actx, dcoll, dd=dd_in)
- inverse_jac_mat = inverse_surface_metric_derivative_mat(actx, dcoll, dd=dd_in)
+ inverse_jac_mat = inverse_surface_metric_derivative_mat(actx, dcoll, dd=dd_in,
+ times_area_element=True)
per_group_grads = [
# r for rst axis
# x for xyz axis
- actx.einsum("rij,ej,ej,xrej->xei",
+ actx.einsum("rij,ej,xrej->xei",
reference_stiffness_transpose_matrix(
actx,
out_element_group=out_grp,
in_element_group=in_grp
),
- ae_i,
vec_i,
ijm_i,
- arg_names=("ref_stiffT_mat", "jac", "vec", "inv_jac_t"),
+ arg_names=("ref_stiffT_mat", "vec", "inv_jac_t"),
tagged=(FirstAxisIsElementsTag(),))
- for out_grp, in_grp, vec_i, ae_i, ijm_i in zip(
- out_discr.groups, in_discr.groups, vec, area_elements,
+ for out_grp, in_grp, vec_i, ijm_i in zip(
+ out_discr.groups, in_discr.groups, vec,
inverse_jac_mat)]
return make_obj_array([
@@ -374,33 +372,31 @@ def weak_local_d_dx(dcoll: DiscretizationCollection, *args):
else:
raise TypeError("invalid number of arguments")
- from grudge.geometry import \
- inverse_surface_metric_derivative_mat, area_element
+ from grudge.geometry import inverse_surface_metric_derivative_mat
in_discr = dcoll.discr_from_dd(dd_in)
out_discr = dcoll.discr_from_dd(dof_desc.DD_VOLUME)
actx = vec.array_context
- area_elements = area_element(actx, dcoll, dd=dd_in)
- inverse_jac_mat = inverse_surface_metric_derivative_mat(actx, dcoll, dd=dd_in)
+ inverse_jac_mat = inverse_surface_metric_derivative_mat(actx, dcoll, dd=dd_in,
+ times_area_element=True)
return DOFArray(
actx,
data=tuple(
- actx.einsum("dij,ej,ej,dej->ei",
+ actx.einsum("dij,ej,dej->ei",
reference_stiffness_transpose_matrix(
actx,
out_element_group=out_grp,
in_element_group=in_grp
),
- ae_i,
vec_i,
ijm_i[xyz_axis],
- arg_names=("ref_stiffT_mat", "jac", "vec", "inv_jac_t"),
+ arg_names=("ref_stiffT_mat", "vec", "inv_jac_t"),
tagged=(FirstAxisIsElementsTag(),))
- for out_grp, in_grp, vec_i, ae_i, ijm_i in zip(
- out_discr.groups, in_discr.groups, vec, area_elements,
+ for out_grp, in_grp, vec_i, ijm_i in zip(
+ out_discr.groups, in_discr.groups, vec,
inverse_jac_mat)))
--
GitLab