From 345da99b0447457424210839198003dd1bd591bd Mon Sep 17 00:00:00 2001
From: Thomas Gibson <gibsonthomas1120@hotmail.com>
Date: Mon, 6 Dec 2021 22:38:35 -0600
Subject: [PATCH] Patch: numpify 1-d metric computations and check for 0-dim
discretizations
---
grudge/geometry/metrics.py | 26 +++++++++++++++++++-------
1 file changed, 19 insertions(+), 7 deletions(-)
diff --git a/grudge/geometry/metrics.py b/grudge/geometry/metrics.py
index af04e09f..34e84e51 100644
--- a/grudge/geometry/metrics.py
+++ b/grudge/geometry/metrics.py
@@ -488,14 +488,19 @@ def _signed_face_ones(
signed_face_ones = dcoll.discr_from_dd(dd).zeros(
actx, dtype=dcoll.real_dtype
) + 1
+
+ from arraycontext import to_numpy, from_numpy, thaw
+
+ _signed_face_ones_numpy = to_numpy(signed_face_ones, actx)
+
for igrp, grp in enumerate(all_faces_conn.groups):
for batch in grp.batches:
- i = actx.thaw(batch.to_element_indices)
- grp_field = signed_face_ones[igrp].reshape(-1)
+ i = to_numpy(thaw(batch.to_element_indices, actx), actx)
+ grp_field = _signed_face_ones_numpy[igrp].reshape(-1)
grp_field[i] = \
(2.0 * (batch.to_element_face % 2) - 1.0) * grp_field[i]
- return signed_face_ones
+ return from_numpy(_signed_face_ones_numpy, actx)
def parametrization_derivative(
@@ -573,8 +578,12 @@ def area_element(
result = actx.np.sqrt(
pseudoscalar(actx, dcoll, dd=dd).norm_squared())
- if _use_geoderiv_connection:
- result = dcoll._base_to_geoderiv_connection(dd)(result)
+ # NOTE: Don't interpolate 0-dim geometric factors to the
+ # "geometry discretization"
+ # (the metrics are just single scalars for facets in 1-D)
+ if dcoll.discr_from_dd(dd).dim != 0:
+ if _use_geoderiv_connection:
+ result = dcoll._base_to_geoderiv_connection(dd)(result)
return freeze(result, actx)
@@ -672,8 +681,11 @@ def mv_normal(
result = mv / actx.np.sqrt(mv.norm_squared())
- if _use_geoderiv_connection:
- result = dcoll._base_to_geoderiv_connection(dd)(result)
+ # NOTE: Don't interpolate normals to the "geometry discretization"
+ # in 1-D (just a single value for every facet: +1 or -1)
+ if dim != 0:
+ if _use_geoderiv_connection:
+ result = dcoll._base_to_geoderiv_connection(dd)(result)
return freeze(result, actx)
--
GitLab