From 6b725f369f5b3347a4baa58d41d1e645ee1e329d Mon Sep 17 00:00:00 2001 From: Andreas Kloeckner Date: Tue, 6 Oct 2015 17:45:23 -0500 Subject: [PATCH] Wave example tweaks --- examples/wave/wave-min.py | 2 +- grudge/models/wave.py | 40 +++++++++++++++------------------------ grudge/shortcuts.py | 8 +++++--- 3 files changed, 21 insertions(+), 29 deletions(-) diff --git a/examples/wave/wave-min.py b/examples/wave/wave-min.py index 3b5f5d32..42f3041e 100644 --- a/examples/wave/wave-min.py +++ b/examples/wave/wave-min.py @@ -72,7 +72,7 @@ def main(write_output=True): #dt = op.estimate_rk4_timestep(discr, fields=fields) dt = 0.001 - dt_stepper = set_up_rk4(dt, fields) + dt_stepper = set_up_rk4(dt, fields, rhs) final_t = 10 nsteps = int(final_t/dt) diff --git a/grudge/models/wave.py b/grudge/models/wave.py index 53a0a848..f6276e70 100644 --- a/grudge/models/wave.py +++ b/grudge/models/wave.py @@ -108,48 +108,38 @@ class StrongWaveOperator(HyperbolicOperator): return -self.c*flux_strong def op_template(self): - from grudge.symbolic import \ - make_sym_vector, \ - BoundaryPair, \ - get_flux_operator, \ - make_nabla, \ - InverseMassOperator, \ - BoundarizeOperator - d = self.dimensions - w = make_sym_vector("w", d+1) + w = sym.make_sym_vector("w", d+1) u = w[0] v = w[1:] # boundary conditions ------------------------------------------------- # dirichlet BCs ------------------------------------------------------- - from grudge.symbolic import normal, Field - - dir_u = BoundarizeOperator(self.dirichlet_tag) * u - dir_v = BoundarizeOperator(self.dirichlet_tag) * v + dir_u = sym.BoundarizeOperator(self.dirichlet_tag) * u + dir_v = sym.BoundarizeOperator(self.dirichlet_tag) * v if self.dirichlet_bc_f: # FIXME from warnings import warn warn("Inhomogeneous Dirichlet conditions on the wave equation " "are still having issues.") - dir_g = Field("dir_bc_u") + dir_g = sym.Field("dir_bc_u") dir_bc = join_fields(2*dir_g - dir_u, dir_v) else: dir_bc = join_fields(-dir_u, dir_v) # neumann BCs --------------------------------------------------------- - neu_u = BoundarizeOperator(self.neumann_tag) * u - neu_v = BoundarizeOperator(self.neumann_tag) * v + neu_u = sym.BoundarizeOperator(self.neumann_tag) * u + neu_v = sym.BoundarizeOperator(self.neumann_tag) * v neu_bc = join_fields(neu_u, -neu_v) # radiation BCs ------------------------------------------------------- - rad_normal = normal(self.radiation_tag, d) + rad_normal = sym.normal(self.radiation_tag, d) - rad_u = BoundarizeOperator(self.radiation_tag) * u - rad_v = BoundarizeOperator(self.radiation_tag) * v + rad_u = sym.BoundarizeOperator(self.radiation_tag) * u + rad_v = sym.BoundarizeOperator(self.radiation_tag) * v rad_bc = join_fields( 0.5*(rad_u - self.sign*np.dot(rad_normal, rad_v)), @@ -157,8 +147,8 @@ class StrongWaveOperator(HyperbolicOperator): ) # entire operator ----------------------------------------------------- - nabla = make_nabla(d) - flux_op = get_flux_operator(self.flux()) + nabla = sym.make_nabla(d) + flux_op = sym.get_flux_operator(self.flux()) result = ( - join_fields( @@ -166,11 +156,11 @@ class StrongWaveOperator(HyperbolicOperator): -self.c*(nabla*u) ) + - InverseMassOperator() * ( + sym.InverseMassOperator() * ( flux_op(w) - + flux_op(BoundaryPair(w, dir_bc, self.dirichlet_tag)) - + flux_op(BoundaryPair(w, neu_bc, self.neumann_tag)) - + flux_op(BoundaryPair(w, rad_bc, self.radiation_tag)) + + flux_op(sym.BoundaryPair(w, dir_bc, self.dirichlet_tag)) + + flux_op(sym.BoundaryPair(w, neu_bc, self.neumann_tag)) + + flux_op(sym.BoundaryPair(w, rad_bc, self.radiation_tag)) )) result[0] += self.source_f diff --git a/grudge/shortcuts.py b/grudge/shortcuts.py index 6ed54d51..6a24bb23 100644 --- a/grudge/shortcuts.py +++ b/grudge/shortcuts.py @@ -39,13 +39,15 @@ def make_discretization(mesh, order, **kwargs): PolynomialWarpAndBlendGroupFactory(order=order)) -def set_up_rk4(dt, fields, t_start=0): +def set_up_rk4(dt, fields, rhs, t_start=0): from leap.method.rk import LSRK4Method from leap.vm.codegen import PythonCodeGenerator dt_method = LSRK4Method(component_id="y") dt_stepper = PythonCodeGenerator().get_class(dt_method.generate()) - dt_stepper.set_up(t_start=t_start, dt_start=dt, context={"y": fields}) + dt_stepper.set_up( + t_start=t_start, dt_start=dt, + context={dt_method.component_id: fields}) - return dt_stepper + return dt_stepper({""+dt_method.component_id: rhs}) -- GitLab