From b0b3afbccadb96e26708bdfdc07af8a437baa295 Mon Sep 17 00:00:00 2001
From: Andreas Kloeckner <inform@tiker.net>
Date: Wed, 26 Aug 2020 16:37:15 -0500
Subject: [PATCH] Fix var-velocity wave for 687a28b
---
examples/wave/var-propagation-speed.py | 2 +-
grudge/models/wave.py | 17 ++++++++---------
2 files changed, 9 insertions(+), 10 deletions(-)
diff --git a/examples/wave/var-propagation-speed.py b/examples/wave/var-propagation-speed.py
index f392b420..858e7254 100644
--- a/examples/wave/var-propagation-speed.py
+++ b/examples/wave/var-propagation-speed.py
@@ -56,7 +56,7 @@ def main(write_output=True, order=4):
sym_t = sym.ScalarVariable("t")
c = sym.If(sym.Comparison(
np.dot(sym_x, sym_x), "<", 0.15),
- np.float32(-0.1), np.float32(-0.2))
+ np.float32(0.1), np.float32(0.2))
from grudge.models.wave import VariableCoefficientWeakWaveOperator
from meshmode.mesh import BTAG_ALL, BTAG_NONE
diff --git a/grudge/models/wave.py b/grudge/models/wave.py
index 09a954e8..1be9abe2 100644
--- a/grudge/models/wave.py
+++ b/grudge/models/wave.py
@@ -370,18 +370,18 @@ class VariableCoefficientWeakWaveOperator(HyperbolicOperator):
v = w[2:]
normal = sym.normal(w.dd, self.ambient_dim)
+ flux_central_weak = -0.5 * flat_obj_array(
+ np.dot(v.int*c.int + v.ext*c.ext, normal),
+ (u.int * c.int + u.ext*c.ext) * normal)
+
if self.flux_type == "central":
- return -0.5 * flat_obj_array(
- np.dot(v.int*c.int + v.ext*c.ext, normal),
- (u.int * c.int + u.ext*c.ext) * normal)
+ return flux_central_weak
elif self.flux_type == "upwind":
- return -0.5 * flat_obj_array(
- np.dot(normal, c.ext * v.ext + c.int * v.int)
- + c.ext*u.ext - c.int * u.int,
+ return flux_central_weak - 0.5 * flat_obj_array(
+ c.ext*u.ext - c.int * u.int,
- normal * (np.dot(normal, c.ext * v.ext - c.int * v.int)
- + c.ext*u.ext + c.int * u.int))
+ normal * (np.dot(normal, c.ext * v.ext - c.int * v.int)))
else:
raise ValueError("invalid flux type '%s'" % self.flux_type)
@@ -443,7 +443,6 @@ class VariableCoefficientWeakWaveOperator(HyperbolicOperator):
-self.c*(sym.stiffness_t(self.ambient_dim)*u)
)
-
- sym.FaceMassOperator()(flux(sym.int_tpair(flux_w))
+ flux(sym.bv_tpair(self.dirichlet_tag, flux_w, dir_bc))
+ flux(sym.bv_tpair(self.neumann_tag, flux_w, neu_bc))
--
GitLab