From 9cbf637a98d2752f1fcc16418166f25ff51356f3 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 23 Apr 2020 22:28:48 -0500 Subject: [PATCH 1/4] try and remove cyclic deps in symbolic module --- examples/advection/var-velocity.py | 5 +- examples/dagrt-fusion.py | 21 ++-- examples/wave/var-propagation-speed.py | 7 +- grudge/__init__.py | 11 +- grudge/models/advection.py | 10 +- grudge/models/em.py | 4 +- grudge/models/wave.py | 6 +- grudge/sym.py | 28 ----- grudge/symbolic/__init__.py | 3 + grudge/symbolic/operators.py | 104 +++++++++---------- grudge/symbolic/primitives.py | 136 ++++++++++++------------- test/test_grudge.py | 16 +-- 12 files changed, 161 insertions(+), 190 deletions(-) delete mode 100644 grudge/sym.py diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index bc720731..adf2f2fe 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -75,8 +75,9 @@ def main(write_output=True, order=4): / source_width**2) def f_step(x): - return sym.If( - sym.Comparison( + from pymbolic.primitives import If, Comparison + return If( + Comparison( np.dot(sym_source_center_dist, sym_source_center_dist), "<", (4*source_width)**2), diff --git a/examples/dagrt-fusion.py b/examples/dagrt-fusion.py index 27e35abd..b9ef892e 100755 --- a/examples/dagrt-fusion.py +++ b/examples/dagrt-fusion.py @@ -50,6 +50,7 @@ from pymbolic.mapper import Mapper from pymbolic.mapper.evaluator import EvaluationMapper \ as PymbolicEvaluationMapper from pytools import memoize +from pytools.obj_array import join_fields from grudge import sym, bind, DGDiscretizationWithBoundaries from leap.rk import LSRK4MethodBuilder @@ -245,7 +246,6 @@ class RK4TimeStepperBase(object): self.component_getter = component_getter def get_initial_context(self, fields, t_start, dt): - from pytools.obj_array import join_fields # Flatten fields. flattened_fields = [] @@ -286,7 +286,6 @@ class RK4TimeStepperBase(object): assert len(output_states) == num_fields assert len(output_states) == len(output_residuals) - from pytools.obj_array import join_fields flattened_results = join_fields(output_t, output_dt, *output_states) self.bound_op = bind( @@ -354,8 +353,6 @@ class RK4TimeStepper(RK4TimeStepperBase): + tuple( Variable(field_var_name, dd=sym.DD_VOLUME)[i] for i in range(num_fields))))) - - from pytools.obj_array import join_fields sym_rhs = join_fields(*(call[i] for i in range(num_fields))) self.queue = queue @@ -377,7 +374,6 @@ class RK4TimeStepper(RK4TimeStepperBase): self.component_getter = component_getter def _bound_op(self, queue, t, *args, profile_data=None): - from pytools.obj_array import join_fields context = { "t": t, self.field_var_name: join_fields(*args)} @@ -452,16 +448,15 @@ def dg_flux(c, tpair): dims = len(v) normal = sym.normal(tpair.dd, dims) - - flux_weak = sym.join_fields( + flux_weak = join_fields( np.dot(v.avg, normal), u.avg * normal) - flux_weak -= (1 if c > 0 else -1)*sym.join_fields( + flux_weak -= (1 if c > 0 else -1)*join_fields( 0.5*(u.int-u.ext), 0.5*(normal * np.dot(normal, v.int-v.ext))) - flux_strong = sym.join_fields( + flux_strong = join_fields( np.dot(v.int, normal), u.int * normal) - flux_weak @@ -507,13 +502,13 @@ def get_strong_wave_op_with_discr_direct(cl_ctx, dims=2, order=4): rad_u = sym.cse(sym.interp("vol", BTAG_ALL)(u)) rad_v = sym.cse(sym.interp("vol", BTAG_ALL)(v)) - rad_bc = sym.cse(sym.join_fields( + rad_bc = sym.cse(join_fields( 0.5*(rad_u - sign*np.dot(rad_normal, rad_v)), 0.5*rad_normal*(np.dot(rad_normal, rad_v) - sign*rad_u) ), "rad_bc") sym_operator = ( - - sym.join_fields( + - join_fields( -c*np.dot(sym.nabla(dims), v) - source_f, -c*(sym.nabla(dims)*u) ) @@ -550,7 +545,6 @@ def test_stepper_equivalence(ctx_factory, order=4): elif dims == 3: dt = 0.02 - from pytools.obj_array import join_fields ic = join_fields(discr.zeros(queue), [discr.zeros(queue) for i in range(discr.dim)]) @@ -804,7 +798,6 @@ def test_stepper_mem_ops(ctx_factory, use_fusion): dt = 0.04 t_end = 0.2 - from pytools.obj_array import join_fields ic = join_fields(discr.zeros(queue), [discr.zeros(queue) for i in range(discr.dim)]) @@ -975,7 +968,6 @@ def test_stepper_timing(ctx_factory, use_fusion): dt = 0.04 t_end = 0.1 - from pytools.obj_array import join_fields ic = join_fields(discr.zeros(queue), [discr.zeros(queue) for i in range(discr.dim)]) @@ -1040,7 +1032,6 @@ def get_example_stepper(queue, dims=2, order=3, use_fusion=True, exec_mapper_factory=exec_mapper_factory) if return_ic: - from pytools.obj_array import join_fields ic = join_fields(discr.zeros(queue), [discr.zeros(queue) for i in range(discr.dim)]) return stepper, ic diff --git a/examples/wave/var-propagation-speed.py b/examples/wave/var-propagation-speed.py index a7a634f7..937b87a0 100644 --- a/examples/wave/var-propagation-speed.py +++ b/examples/wave/var-propagation-speed.py @@ -51,9 +51,10 @@ def main(write_output=True, order=4): sym_x = sym.nodes(mesh.dim) sym_source_center_dist = sym_x - source_center 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)) + + from pymbolic.primitives import If, Comparison + c = If(Comparison(np.dot(sym_x, sym_x), "<", 0.15), + 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/__init__.py b/grudge/__init__.py index b854007a..17358dd4 100644 --- a/grudge/__init__.py +++ b/grudge/__init__.py @@ -22,6 +22,11 @@ OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. """ -import grudge.sym # noqa -from grudge.execution import bind # noqa -from grudge.discretization import DGDiscretizationWithBoundaries # noqa +import grudge.symbolic as sym +from grudge.execution import bind + +from grudge.discretization import DGDiscretizationWithBoundaries + +__all__ = [ + "sym", "bind", "DGDiscretizationWithBoundaries" +] diff --git a/grudge/models/advection.py b/grudge/models/advection.py index 091b9870..5ef1e5ce 100644 --- a/grudge/models/advection.py +++ b/grudge/models/advection.py @@ -59,9 +59,10 @@ class AdvectionOperatorBase(HyperbolicOperator): elif self.flux_type == "lf": return u.avg*v_dot_normal + 0.5*norm_v*(u.int - u.ext) elif self.flux_type == "upwind": + from pymbolic.primitives import If, Comparison return ( - v_dot_normal * sym.If( - sym.Comparison(v_dot_normal, ">", 0), + v_dot_normal * If( + Comparison(v_dot_normal, ">", 0), u.int, # outflow u.ext, # inflow )) @@ -155,9 +156,10 @@ class VariableCoefficientAdvectionOperator(HyperbolicOperator): elif self.flux_type == "lf": return u.avg*v_dot_normal + 0.5*norm_v*(u.int - u.ext) elif self.flux_type == "upwind": + from pymbolic.primitives import If, Comparison return ( - v_dot_normal * sym.If( - sym.Comparison(v_dot_normal, ">", 0), + v_dot_normal * If( + Comparison(v_dot_normal, ">", 0), u.int, # outflow u.ext, # inflow )) diff --git a/grudge/models/em.py b/grudge/models/em.py index 3422c3be..225b3b4a 100644 --- a/grudge/models/em.py +++ b/grudge/models/em.py @@ -437,7 +437,7 @@ def get_rectangular_cavity_mode(E_0, mode_indices): # noqa: N803 if dims == 2: tfac = sym.ScalarVariable("t") * omega - result = sym.join_fields( + result = join_fields( 0, 0, sym.sin(kx * x) * sym.sin(ky * y) * sym.cos(tfac), # ez @@ -449,7 +449,7 @@ def get_rectangular_cavity_mode(E_0, mode_indices): # noqa: N803 tdep = sym.exp(-1j * omega * sym.ScalarVariable("t")) gamma_squared = ky**2 + kx**2 - result = sym.join_fields( + result = join_fields( -kx * kz * E_0*cx*sy*sz*tdep / gamma_squared, # ex -ky * kz * E_0*sx*cy*sz*tdep / gamma_squared, # ey E_0 * sx*sy*cz*tdep, # ez diff --git a/grudge/models/wave.py b/grudge/models/wave.py index 0ac60a33..9b16c602 100644 --- a/grudge/models/wave.py +++ b/grudge/models/wave.py @@ -350,9 +350,9 @@ class VariableCoefficientWeakWaveOperator(HyperbolicOperator): self.ambient_dim = ambient_dim self.source_f = source_f - self.sign = sym.If(sym.Comparison( - self.c, ">", 0), - np.int32(1), np.int32(-1)) + from pymbolic.primitives import If, Comparison + self.sign = If(Comparison(self.c, ">", 0), + np.int32(1), np.int32(-1)) self.dirichlet_tag = dirichlet_tag self.neumann_tag = neumann_tag diff --git a/grudge/sym.py b/grudge/sym.py deleted file mode 100644 index bf9c7de6..00000000 --- a/grudge/sym.py +++ /dev/null @@ -1,28 +0,0 @@ -from __future__ import division, absolute_import - -__copyright__ = "Copyright (C) 2015 Andreas Kloeckner" - -__license__ = """ -Permission is hereby granted, free of charge, to any person obtaining a copy -of this software and associated documentation files (the "Software"), to deal -in the Software without restriction, including without limitation the rights -to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -copies of the Software, and to permit persons to whom the Software is -furnished to do so, subject to the following conditions: - -The above copyright notice and this permission notice shall be included in -all copies or substantial portions of the Software. - -THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN -THE SOFTWARE. -""" - -from grudge.symbolic.primitives import * # noqa -from grudge.symbolic.operators import * # noqa - -from grudge.symbolic.tools import pretty # noqa diff --git a/grudge/symbolic/__init__.py b/grudge/symbolic/__init__.py index acd34e35..b5182aa8 100644 --- a/grudge/symbolic/__init__.py +++ b/grudge/symbolic/__init__.py @@ -23,3 +23,6 @@ LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. """ + +from grudge.symbolic.primitives import * # noqa: F401,F403 +from grudge.symbolic.operators import * # noqa: F401,F403 diff --git a/grudge/symbolic/operators.py b/grudge/symbolic/operators.py index a1899ed2..7daea219 100644 --- a/grudge/symbolic/operators.py +++ b/grudge/symbolic/operators.py @@ -27,36 +27,28 @@ THE SOFTWARE. from six.moves import intern import numpy as np -import numpy.linalg as la # noqa import pymbolic.primitives -def _sym(): - # A hack to make referring to grudge.sym possible without - # circular imports and tons of typing. - - from grudge import sym - return sym - - # {{{ base classes class Operator(pymbolic.primitives.Expression): """ .. attribute:: dd_in - an instance of :class:`grudge.sym.DOFDesc` describing the + an instance of :class:`~grudge.sym.DOFDesc` describing the input discretization. .. attribute:: dd_out - an instance of :class:`grudge.sym.DOFDesc` describing the + an instance of :class:`~grudge.sym.DOFDesc` describing the output discretization. """ def __init__(self, dd_in, dd_out): - self.dd_in = _sym().as_dofdesc(dd_in) - self.dd_out = _sym().as_dofdesc(dd_out) + import grudge.symbolic.primitives as prim + self.dd_in = prim.as_dofdesc(dd_in) + self.dd_out = prim.as_dofdesc(dd_out) def stringifier(self): from grudge.symbolic.mappers import StringifyMapper @@ -100,8 +92,10 @@ class ElementwiseLinearOperator(Operator): class InterpolationOperator(Operator): def __init__(self, dd_in, dd_out): - official_dd_in = _sym().as_dofdesc(dd_in) - official_dd_out = _sym().as_dofdesc(dd_out) + import grudge.symbolic.primitives as prim + official_dd_in = prim.as_dofdesc(dd_in) + official_dd_out = prim.as_dofdesc(dd_out) + if official_dd_in == official_dd_out: raise ValueError("Interpolating from {} to {}" " does not do anything.".format(official_dd_in, official_dd_out)) @@ -123,7 +117,8 @@ class ElementwiseMaxOperator(Operator): class NodalReductionOperator(Operator): def __init__(self, dd_in, dd_out=None): if dd_out is None: - dd_out = _sym().DD_SCALAR + import grudge.symbolic.primitives as prim + dd_out = prim.DD_SCALAR assert dd_out.is_scalar() @@ -150,13 +145,14 @@ class NodalMin(NodalReductionOperator): class DiffOperatorBase(Operator): def __init__(self, xyz_axis, dd_in=None, dd_out=None): + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = _sym().DD_VOLUME + dd_in = prim.DD_VOLUME if dd_out is None: - dd_out = dd_in.with_qtag(_sym().QTAG_NONE) + dd_out = dd_in.with_qtag(prim.QTAG_NONE) else: - dd_out = _sym().as_dofdesc(dd_out) + dd_out = prim.as_dofdesc(dd_out) if dd_out.uses_quadrature(): raise ValueError("differentiation outputs are not on " @@ -207,10 +203,13 @@ class MInvSTOperator(WeakFormDiffOperatorBase): class RefDiffOperatorBase(ElementwiseLinearOperator): def __init__(self, rst_axis, dd_in=None, dd_out=None): + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = _sym().DD_VOLUME + dd_in = prim.DD_VOLUME + if dd_out is None: - dd_out = dd_in.with_qtag(_sym().QTAG_NONE) + dd_out = dd_in.with_qtag(prim.QTAG_NONE) + if dd_out.uses_quadrature(): raise ValueError("differentiation outputs are not on " "quadrature grids") @@ -275,8 +274,10 @@ class FilterOperator(ElementwiseLinearOperator): (For example an instance of :class:`ExponentialFilterResponseFunction`. """ + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = _sym().DD_VOLUME + dd_in = prim.DD_VOLUME + if dd_out is None: dd_out = dd_in @@ -314,8 +315,10 @@ class FilterOperator(ElementwiseLinearOperator): class AveragingOperator(ElementwiseLinearOperator): def __init__(self, dd_in=None, dd_out=None): + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = _sym().DD_VOLUME + dd_in = prim.DD_VOLUME + if dd_out is None: dd_out = dd_in @@ -362,8 +365,9 @@ class MassOperatorBase(Operator): """ def __init__(self, dd_in=None, dd_out=None): + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = _sym().DD_VOLUME + dd_in = prim.DD_VOLUME if dd_out is None: dd_out = dd_in @@ -429,15 +433,14 @@ class OppositeInteriorFaceSwap(Operator): """ def __init__(self, dd_in=None, dd_out=None, unique_id=None): - sym = _sym() - + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = sym.DOFDesc(sym.FACE_RESTR_INTERIOR, None) + dd_in = prim.DOFDesc(prim.FACE_RESTR_INTERIOR, None) if dd_out is None: dd_out = dd_in super(OppositeInteriorFaceSwap, self).__init__(dd_in, dd_out) - if self.dd_in.domain_tag is not sym.FACE_RESTR_INTERIOR: + if self.dd_in.domain_tag is not prim.FACE_RESTR_INTERIOR: raise ValueError("dd_in must be an interior faces domain") if self.dd_out != self.dd_in: raise ValueError("dd_out and dd_in must be identical") @@ -462,7 +465,7 @@ class OppositePartitionFaceSwap(Operator): MPI tag offset to keep different subexpressions apart in MPI traffic. """ def __init__(self, dd_in=None, dd_out=None, unique_id=None): - sym = _sym() + import grudge.symbolic.primitives as prim if dd_in is None and dd_out is None: raise ValueError("dd_in or dd_out must be specified") @@ -472,7 +475,7 @@ class OppositePartitionFaceSwap(Operator): dd_out = dd_in super(OppositePartitionFaceSwap, self).__init__(dd_in, dd_out) - if not isinstance(self.dd_in.domain_tag, sym.BTAG_PARTITION): + if not isinstance(self.dd_in.domain_tag, prim.BTAG_PARTITION): raise ValueError("dd_in must be a partition boundary faces domain") if self.dd_out != self.dd_in: raise ValueError("dd_out and dd_in must be identical") @@ -492,20 +495,20 @@ class OppositePartitionFaceSwap(Operator): class FaceMassOperatorBase(ElementwiseLinearOperator): def __init__(self, dd_in=None, dd_out=None): - sym = _sym() - + import grudge.symbolic.primitives as prim if dd_in is None: - dd_in = sym.DOFDesc(sym.FACE_RESTR_ALL, None) + dd_in = prim.DOFDesc(prim.FACE_RESTR_ALL, None) if dd_out is None or dd_out == "vol": - dd_out = sym.DOFDesc("vol", sym.QTAG_NONE) + dd_out = prim.DOFDesc("vol", prim.QTAG_NONE) + if dd_out.uses_quadrature(): raise ValueError("face mass operator outputs are not on " "quadrature grids") if not dd_out.is_volume(): raise ValueError("dd_out must be a volume domain") - if dd_in.domain_tag is not sym.FACE_RESTR_ALL: + if dd_in.domain_tag is not prim.FACE_RESTR_ALL: raise ValueError("dd_in must be an interior faces domain") super(FaceMassOperatorBase, self).__init__(dd_in, dd_out) @@ -575,43 +578,40 @@ def stiffness_t(dim, dd_in=None, dd_out=None): def integral(arg, dd=None): - sym = _sym() + import grudge.symbolic.primitives as prim if dd is None: - dd = sym.DD_VOLUME - - dd = sym.as_dofdesc(dd) + dd = prim.DD_VOLUME + dd = prim.as_dofdesc(dd) - return sym.NodalSum(dd)( - arg * sym.cse( - sym.MassOperator(dd_in=dd)(sym.Ones(dd)), + return NodalSum(dd)( + arg * prim.cse( + MassOperator(dd_in=dd)(prim.Ones(dd)), "mass_quad_weights", - sym.cse_scope.DISCRETIZATION)) + prim.cse_scope.DISCRETIZATION)) def norm(p, arg, dd=None): """ :arg arg: is assumed to be a vector, i.e. have shape ``(n,)``. """ - sym = _sym() + import grudge.symbolic.primitives as prim if dd is None: - dd = sym.DD_VOLUME - - dd = sym.as_dofdesc(dd) + dd = prim.DD_VOLUME + dd = prim.as_dofdesc(dd) if p == 2: - norm_squared = sym.NodalSum(dd_in=dd)( - sym.FunctionSymbol("fabs")( - arg * sym.MassOperator()(arg))) + norm_squared = NodalSum(dd_in=dd)( + prim.fabs(arg * MassOperator()(arg))) if isinstance(norm_squared, np.ndarray): norm_squared = norm_squared.sum() - return sym.FunctionSymbol("sqrt")(norm_squared) + return prim.sqrt(norm_squared) elif p == np.Inf: - result = sym.NodalMax(dd_in=dd)(sym.FunctionSymbol("fabs")(arg)) + result = NodalMax(dd_in=dd)(prim.fabs(arg)) from pymbolic.primitives import Max if isinstance(result, np.ndarray): diff --git a/grudge/symbolic/primitives.py b/grudge/symbolic/primitives.py index 0b9fac20..c43ca30d 100644 --- a/grudge/symbolic/primitives.py +++ b/grudge/symbolic/primitives.py @@ -27,37 +27,35 @@ THE SOFTWARE. from six.moves import range, intern import numpy as np -import pymbolic.primitives -from meshmode.mesh import BTAG_ALL, BTAG_REALLY_ALL, BTAG_NONE, BTAG_PARTITION # noqa -from meshmode.discretization.connection import ( # noqa - FACE_RESTR_ALL, FACE_RESTR_INTERIOR) - -from pymbolic.primitives import ( # noqa +from pytools.obj_array import make_obj_array + +from meshmode.mesh import ( + BTAG_ALL, + BTAG_REALLY_ALL, + BTAG_NONE, + BTAG_PARTITION) +from meshmode.discretization.connection import ( + FACE_RESTR_ALL, + FACE_RESTR_INTERIOR) + +import pymbolic.primitives as prim +from pymbolic.primitives import ( + Variable as VariableBase, + CommonSubexpression as CommonSubexpressionBase, cse_scope as cse_scope_base, - make_common_subexpression as cse, If, Comparison) + make_common_subexpression as cse) from pymbolic.geometric_algebra import MultiVector -from pytools.obj_array import join_fields, make_obj_array # noqa -class ExpressionBase(pymbolic.primitives.Expression): +class ExpressionBase(prim.Expression): def make_stringifier(self, originating_stringifier=None): from grudge.symbolic.mappers import StringifyMapper return StringifyMapper() -def _sym(): - # A hack to make referring to grudge.sym possible without - # circular imports and tons of typing. - - from grudge import sym - return sym - - __doc__ = """ -.. currentmodule:: grudge.sym - -.. autoclass:: If +.. currentmodule:: grudge.symbolic.primitives DOF description ^^^^^^^^^^^^^^^ @@ -66,6 +64,7 @@ DOF description .. autoclass:: DTAG_VOLUME_ALL .. autoclass:: QTAG_NONE .. autoclass:: DOFDesc + .. data:: DD_SCALAR .. data:: DD_VOLUME @@ -78,6 +77,7 @@ Symbols .. autoclass:: make_sym_mv .. autoclass:: FunctionSymbol +.. function :: fabs(arg) .. function :: sqrt(arg) .. function :: exp(arg) .. function :: sin(arg) @@ -111,20 +111,20 @@ Geometry data # {{{ DOF description -class DTAG_SCALAR: # noqa +class DTAG_SCALAR: # noqa: N801 pass -class DTAG_VOLUME_ALL: # noqa +class DTAG_VOLUME_ALL: # noqa: N801 pass -class DTAG_BOUNDARY: # noqa +class DTAG_BOUNDARY: # noqa: N801 def __init__(self, tag): self.tag = tag -class QTAG_NONE: # noqa +class QTAG_NONE: # noqa: N801 pass @@ -133,14 +133,18 @@ class DOFDesc(object): .. attribute:: domain_tag .. attribute:: quadrature_tag + .. automethod:: is_scalar .. automethod:: is_discretized .. automethod:: is_volume .. automethod:: is_boundary .. automethod:: is_trace + .. automethod:: uses_quadrature + .. automethod:: with_qtag .. automethod:: with_dtag + .. automethod:: __eq__ .. automethod:: __ne__ .. automethod:: __hash__ @@ -149,46 +153,32 @@ class DOFDesc(object): def __init__(self, domain_tag, quadrature_tag=None): """ :arg domain_tag: One of the following: - :class:`grudge.sym.DTAG_SCALAR` (or the string ``"scalar"``), - :class:`grudge.sym.DTAG_VOLUME_ALL` (or the string ``"vol"``) + :class:`DTAG_SCALAR` (or the string ``"scalar"``), + :class:`DTAG_VOLUME_ALL` (or the string ``"vol"``) for the default volume discretization, - :class:`meshmode.discretization.connection.FACE_RESTR_ALL` - (or the string ``"all_faces"``), - or - :class:`meshmode.discretization.connection.FACE_RESTR_INTERIOR` - (or the string ``"int_faces"``), - or one of - :class:`meshmode.discretization.BTAG_ALL`, - :class:`meshmode.discretization.BTAG_NONE`, - :class:`meshmode.discretization.BTAG_REALLY_ALL`, - :class:`meshmode.discretization.PARTITION`, - or :class + :class:`FACE_RESTR_ALL` (or the string ``"all_faces"``), + or :class:`FACE_RESTR_INTERIOR` (or the string ``"int_faces"``), + or one of :class:`BTAG_ALL`, :class:`BTAG_NONE`, + :class:`BTAG_REALLY_ALL`, :class:`BTAG_PARTITION`, or *None* to indicate that the geometry is not yet known. :arg quadrature_tag: - *None* to indicate that the quadrature grid is not known,or + *None* to indicate that the quadrature grid is not known, or :class:`QTAG_NONE` to indicate the use of the basic discretization grid, or a string to indicate the use of the thus-tagged quadratue grid. """ - if domain_tag == "scalar": - domain_tag = DTAG_SCALAR - elif domain_tag is DTAG_SCALAR: + + if domain_tag is None: + pass + elif domain_tag in [DTAG_SCALAR, "scalar"]: domain_tag = DTAG_SCALAR - elif domain_tag == "vol": + elif domain_tag in [DTAG_VOLUME_ALL, "vol"]: domain_tag = DTAG_VOLUME_ALL - elif domain_tag is DTAG_VOLUME_ALL: - pass - elif domain_tag == "all_faces": + elif domain_tag in [FACE_RESTR_ALL, "all_faces"]: domain_tag = FACE_RESTR_ALL - elif domain_tag is FACE_RESTR_ALL: - pass - elif domain_tag == "int_faces": + elif domain_tag in [FACE_RESTR_INTERIOR, "int_faces"]: domain_tag = FACE_RESTR_INTERIOR - elif domain_tag is FACE_RESTR_INTERIOR: - pass - elif domain_tag is None: - pass elif isinstance(domain_tag, BTAG_PARTITION): pass elif domain_tag in [BTAG_ALL, BTAG_REALLY_ALL, BTAG_NONE]: @@ -273,8 +263,7 @@ DD_VOLUME = DOFDesc(DTAG_VOLUME_ALL, None) def as_dofdesc(dd): if isinstance(dd, DOFDesc): return dd - else: - return DOFDesc(dd, None) + return DOFDesc(dd, quadrature_tag=None) # }}} @@ -314,11 +303,11 @@ class HasDOFDesc(object): # {{{ variables -class cse_scope(cse_scope_base): # noqa +class cse_scope(cse_scope_base): # noqa: N801 DISCRETIZATION = "grudge_discretization" -class Variable(HasDOFDesc, ExpressionBase, pymbolic.primitives.Variable): +class Variable(HasDOFDesc, ExpressionBase, VariableBase): """A user-supplied input variable with a known :class:`DOFDesc`. """ init_arg_names = ("name", "dd") @@ -347,7 +336,7 @@ def make_sym_array(name, shape, dd=None): def var_factory(name): return Variable(name, dd) - return pymbolic.primitives.make_sym_array(name, shape, var_factory) + return prim.make_sym_array(name, shape, var_factory) def make_sym_mv(name, dim, var_factory=None): @@ -355,15 +344,15 @@ def make_sym_mv(name, dim, var_factory=None): make_sym_array(name, dim, var_factory)) -class FunctionSymbol(ExpressionBase, pymbolic.primitives.Variable): +class FunctionSymbol(ExpressionBase, VariableBase): """A symbol to be used as the function argument of - :class:`pymbolic.primitives.Call`. - + :class:`~pymbolic.primitives.Call`. """ mapper_method = "map_function_symbol" +fabs = FunctionSymbol("fabs") sqrt = FunctionSymbol("sqrt") exp = FunctionSymbol("exp") sin = FunctionSymbol("sin") @@ -399,7 +388,7 @@ class OperatorBinding(ExpressionBase): return hash((self.__class__, self.op, obj_array_to_hashable(self.field))) -class PrioritizedSubexpression(pymbolic.primitives.CommonSubexpression): +class PrioritizedSubexpression(CommonSubexpressionBase): """When the optemplate-to-code transformation is performed, prioritized subexpressions work like common subexpression in that they are assigned their own separate identifier/register @@ -409,7 +398,7 @@ class PrioritizedSubexpression(pymbolic.primitives.CommonSubexpression): """ def __init__(self, child, priority=0): - pymbolic.primitives.CommonSubexpression.__init__(self, child) + super(PrioritizedSubexpression, self).__init__(child) self.priority = priority def __getinitargs__(self): @@ -477,12 +466,14 @@ def forward_metric_derivative(xyz_axis, rst_axis, dd=None): inner_dd = dd.with_qtag(QTAG_NONE) - diff_op = _sym().RefDiffOperator(rst_axis, inner_dd) + from grudge.symbolic.operators import RefDiffOperator + diff_op = RefDiffOperator(rst_axis, inner_dd) result = diff_op(NodeCoordinateComponent(xyz_axis, inner_dd)) if dd.uses_quadrature(): - result = _sym().interp(inner_dd, dd)(result) + from grudge.symbolic.operators import interp + result = interp(inner_dd, dd)(result) return cse( result, @@ -676,13 +667,15 @@ class TracePair: def int_tpair(expression, qtag=None): - i = _sym().interp("vol", "int_faces")(expression) - e = cse(_sym().OppositeInteriorFaceSwap()(i)) + from grudge.symbolic.operators import interp, OppositeInteriorFaceSwap + + i = interp("vol", "int_faces")(expression) + e = cse(OppositeInteriorFaceSwap()(i)) - if qtag is not None and qtag != _sym().QTAG_NONE: - q_dd = _sym().DOFDesc("int_faces", qtag) - i = cse(_sym().interp("int_faces", q_dd)(i)) - e = cse(_sym().interp("int_faces", q_dd)(e)) + if qtag is not None and qtag != QTAG_NONE: + q_dd = DOFDesc("int_faces", qtag) + i = cse(interp("int_faces", q_dd)(i)) + e = cse(interp("int_faces", q_dd)(e)) else: q_dd = "int_faces" @@ -710,7 +703,8 @@ def bv_tpair(dd, interior, exterior): representing the exterior value to be used for the flux. """ - interior = _sym().cse(_sym().interp("vol", dd)(interior)) + from grudge.symbolic.operators import interp + interior = cse(interp("vol", dd)(interior)) return TracePair(dd, interior, exterior) # }}} diff --git a/test/test_grudge.py b/test/test_grudge.py index e373e133..1d2732a9 100644 --- a/test/test_grudge.py +++ b/test/test_grudge.py @@ -23,11 +23,14 @@ THE SOFTWARE. """ -import numpy as np # noqa -import numpy.linalg as la # noqa -import pyopencl as cl # noqa -import pyopencl.array # noqa -import pyopencl.clmath # noqa +import numpy as np +import numpy.linalg as la + +import pyopencl as cl +import pyopencl.array +import pyopencl.clmath + +from pytools.obj_array import join_fields import pytest # noqa @@ -185,7 +188,7 @@ def test_2d_gauss_theorem(ctx_factory): discr = DGDiscretizationWithBoundaries(cl_ctx, mesh, order=2) def f(x): - return sym.join_fields( + return join_fields( sym.sin(3*x[0])+sym.cos(3*x[1]), sym.sin(2*x[0])+sym.cos(x[1])) @@ -408,7 +411,6 @@ def test_improvement_quadrature(ctx_factory, order): from meshmode.mesh.generation import generate_regular_rect_mesh from grudge.models.advection import VariableCoefficientAdvectionOperator from pytools.convergence import EOCRecorder - from pytools.obj_array import join_fields from meshmode.discretization.poly_element import QuadratureSimplexGroupFactory cl_ctx = cl.create_some_context() -- GitLab From 06650221dc4a1b67c0f5709957c9683c446e00dc Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Thu, 23 Apr 2020 23:57:22 -0500 Subject: [PATCH 2/4] add docs --- doc/conf.py | 4 +-- doc/index.rst | 11 ++----- doc/symbolic.rst | 14 +++++++++ grudge/symbolic/__init__.py | 3 ++ grudge/symbolic/operators.py | 58 +++++++++++++++++++++++++++++++++-- grudge/symbolic/primitives.py | 26 ++++++++++------ 6 files changed, 94 insertions(+), 22 deletions(-) create mode 100644 doc/symbolic.rst diff --git a/doc/conf.py b/doc/conf.py index 3becabe6..560efec5 100644 --- a/doc/conf.py +++ b/doc/conf.py @@ -320,11 +320,11 @@ texinfo_documents = [ # Example configuration for intersphinx: refer to the Python standard library. intersphinx_mapping = { 'http://docs.python.org/': None, - 'http://documen.tician.de/boxtree/': None, 'http://docs.scipy.org/doc/numpy/': None, - 'http://documen.tician.de/modepy/': None, 'http://documen.tician.de/pyopencl/': None, + 'http://documen.tician.de/modepy/': None, 'http://documen.tician.de/pymbolic/': None, + 'http://documen.tician.de/meshmode/': None, #'http://documen.tician.de/loopy/': None, } autoclass_content = "both" diff --git a/doc/index.rst b/doc/index.rst index a20ec509..b996c470 100644 --- a/doc/index.rst +++ b/doc/index.rst @@ -1,9 +1,4 @@ -.. grudge documentation master file, created by - sphinx-quickstart on Sun Sep 27 13:08:30 2015. - You can adapt this file completely to your liking, but it should at least - contain the root `toctree` directive. - -Welcome to grudge's documentation! +Welcome to grudge's Documentation! ================================== Contents: @@ -12,10 +7,10 @@ Contents: :maxdepth: 2 misc + symbolic - -Indices and tables +Indices and Tables ================== * :ref:`genindex` diff --git a/doc/symbolic.rst b/doc/symbolic.rst new file mode 100644 index 00000000..632b1f86 --- /dev/null +++ b/doc/symbolic.rst @@ -0,0 +1,14 @@ +Symbolic Operator Representation +================================ + +Based on :mod:`pymbolic`. + +Basic Objects +------------- + +.. automodule:: grudge.symbolic.primitives + +Operators +--------- + +.. automodule:: grudge.symbolic.operators diff --git a/grudge/symbolic/__init__.py b/grudge/symbolic/__init__.py index b5182aa8..59fa5b26 100644 --- a/grudge/symbolic/__init__.py +++ b/grudge/symbolic/__init__.py @@ -26,3 +26,6 @@ THE SOFTWARE. from grudge.symbolic.primitives import * # noqa: F401,F403 from grudge.symbolic.operators import * # noqa: F401,F403 +from grudge.symbolic.tools import pretty # noqa: F401 + +from pymbolic.primitives import If, Comparison # noqa: F401 diff --git a/grudge/symbolic/operators.py b/grudge/symbolic/operators.py index 7daea219..20f28e59 100644 --- a/grudge/symbolic/operators.py +++ b/grudge/symbolic/operators.py @@ -1,5 +1,3 @@ -"""Building blocks and mappers for operator expression trees.""" - from __future__ import division, absolute_import __copyright__ = "Copyright (C) 2008-2017 Andreas Kloeckner, Bogdan Enache" @@ -29,6 +27,62 @@ from six.moves import intern import numpy as np import pymbolic.primitives +__doc__ = """ + +Building blocks and mappers for operator expression trees. + +Basic Operators +^^^^^^^^^^^^^^^ + +.. autoclass:: Operator +.. autoclass:: ElementwiseLinearOperator +.. autoclass:: InterpolationOperator + +.. data:: interp + +Reductions +^^^^^^^^^^ + +.. autoclass:: ElementwiseMaxOperator + +.. autoclass:: NodalReductionOperator +.. autoclass:: NodalSum +.. autoclass:: NodalMax +.. autoclass:: NodalMin + +Differentiation +^^^^^^^^^^^^^^^ + +.. autoclass:: StrongFormDiffOperatorBase +.. autoclass:: WeakFormDiffOperatorBase +.. autoclass:: StiffnessOperator +.. autoclass:: DiffOperator +.. autoclass:: StiffnessTOperator +.. autoclass:: MInvSTOperator + +.. autoclass:: RefDiffOperator +.. autoclass:: RefStiffnessTOperator + +.. autofunction:: nabla +.. autofunction:: minv_stiffness_t +.. autofunction:: stiffness +.. autofunction:: stiffness_t + +Mass Operators +^^^^^^^^^^^^^^ + +.. autoclass:: MassOperatorBase + +.. autoclass:: MassOperator +.. autoclass:: InverseMassOperator +.. autoclass:: FaceMassOperator + +.. autoclass:: RefMassOperator +.. autoclass:: RefInverseMassOperator +.. autoclass:: RefFaceMassOperator + +""" + # {{{ base classes diff --git a/grudge/symbolic/primitives.py b/grudge/symbolic/primitives.py index c43ca30d..96279055 100644 --- a/grudge/symbolic/primitives.py +++ b/grudge/symbolic/primitives.py @@ -55,15 +55,16 @@ class ExpressionBase(prim.Expression): __doc__ = """ -.. currentmodule:: grudge.symbolic.primitives - DOF description ^^^^^^^^^^^^^^^ .. autoclass:: DTAG_SCALAR .. autoclass:: DTAG_VOLUME_ALL +.. autoclass:: DTAG_BOUNDARY .. autoclass:: QTAG_NONE + .. autoclass:: DOFDesc +.. autofunction:: as_dofdesc .. data:: DD_SCALAR .. data:: DD_VOLUME @@ -73,10 +74,11 @@ Symbols .. autoclass:: Variable .. autoclass:: ScalarVariable -.. autoclass:: make_sym_array -.. autoclass:: make_sym_mv .. autoclass:: FunctionSymbol +.. autofunction:: make_sym_array +.. autofunction:: make_sym_mv + .. function :: fabs(arg) .. function :: sqrt(arg) .. function :: exp(arg) @@ -90,11 +92,11 @@ Helpers .. autoclass:: OperatorBinding .. autoclass:: PrioritizedSubexpression -.. autoclass:: BoundaryPair .. autoclass:: Ones Geometry data ^^^^^^^^^^^^^ + .. autoclass:: NodeCoordinateComponent .. autofunction:: nodes .. autofunction:: mv_nodes @@ -156,10 +158,14 @@ class DOFDesc(object): :class:`DTAG_SCALAR` (or the string ``"scalar"``), :class:`DTAG_VOLUME_ALL` (or the string ``"vol"``) for the default volume discretization, - :class:`FACE_RESTR_ALL` (or the string ``"all_faces"``), - or :class:`FACE_RESTR_INTERIOR` (or the string ``"int_faces"``), - or one of :class:`BTAG_ALL`, :class:`BTAG_NONE`, - :class:`BTAG_REALLY_ALL`, :class:`BTAG_PARTITION`, + :class:`~meshmode.discretization.connection.FACE_RESTR_ALL` + (or the string ``"all_faces"``), or + :class:`~meshmode.discretization.connection.FACE_RESTR_INTERIOR` + (or the string ``"int_faces"``), or one of + :class:`~meshmode.mesh.BTAG_ALL`, + :class:`~meshmode.mesh.BTAG_NONE`, + :class:`~meshmode.mesh.BTAG_REALLY_ALL`, + :class:`~meshmode.mesh.BTAG_PARTITION`, or *None* to indicate that the geometry is not yet known. :arg quadrature_tag: @@ -586,7 +592,7 @@ def area_element(ambient_dim, dim=None, dd=None): def mv_normal(dd, ambient_dim, dim=None): - """Exterior unit normal as a MultiVector.""" + """Exterior unit normal as a :class:`~pymbolic.geometric_algebra.MultiVector`.""" dd = as_dofdesc(dd) -- GitLab From d36d43859947414bf8285d0a5d73f22a21b9f508 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Fri, 24 Apr 2020 15:36:28 -0500 Subject: [PATCH 3/4] fix doc links --- grudge/symbolic/operators.py | 4 ++-- grudge/symbolic/primitives.py | 6 +++--- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/grudge/symbolic/operators.py b/grudge/symbolic/operators.py index 20f28e59..c1360b16 100644 --- a/grudge/symbolic/operators.py +++ b/grudge/symbolic/operators.py @@ -90,12 +90,12 @@ class Operator(pymbolic.primitives.Expression): """ .. attribute:: dd_in - an instance of :class:`~grudge.sym.DOFDesc` describing the + an instance of :class:`~grudge.symbolic.primitives.DOFDesc` describing the input discretization. .. attribute:: dd_out - an instance of :class:`~grudge.sym.DOFDesc` describing the + an instance of :class:`~grudge.symbolic.primitives.DOFDesc` describing the output discretization. """ diff --git a/grudge/symbolic/primitives.py b/grudge/symbolic/primitives.py index 96279055..ea9e89f6 100644 --- a/grudge/symbolic/primitives.py +++ b/grudge/symbolic/primitives.py @@ -158,9 +158,9 @@ class DOFDesc(object): :class:`DTAG_SCALAR` (or the string ``"scalar"``), :class:`DTAG_VOLUME_ALL` (or the string ``"vol"``) for the default volume discretization, - :class:`~meshmode.discretization.connection.FACE_RESTR_ALL` + :data:`~meshmode.discretization.connection.FACE_RESTR_ALL` (or the string ``"all_faces"``), or - :class:`~meshmode.discretization.connection.FACE_RESTR_INTERIOR` + :data:`~meshmode.discretization.connection.FACE_RESTR_INTERIOR` (or the string ``"int_faces"``), or one of :class:`~meshmode.mesh.BTAG_ALL`, :class:`~meshmode.mesh.BTAG_NONE`, @@ -465,7 +465,7 @@ def forward_metric_derivative(xyz_axis, rst_axis, dd=None): .. math:: - \frac{d x_{\mathtt{xyz\_axis}} }{d r_{\mathtt{rst\_axis}} } + \frac{\partial x_{\mathrm{xyz\_axis}} }{\partial r_{\mathrm{rst\_axis}} } """ if dd is None: dd = DD_VOLUME -- GitLab From 9151bb64a8fe7114d69808453b14a74a38741601 Mon Sep 17 00:00:00 2001 From: Alexandru Fikl Date: Fri, 24 Apr 2020 15:49:08 -0500 Subject: [PATCH 4/4] remove some useless changes --- examples/advection/var-velocity.py | 5 ++--- examples/wave/var-propagation-speed.py | 7 +++---- grudge/models/advection.py | 10 ++++------ grudge/models/wave.py | 6 +++--- 4 files changed, 12 insertions(+), 16 deletions(-) diff --git a/examples/advection/var-velocity.py b/examples/advection/var-velocity.py index adf2f2fe..bc720731 100644 --- a/examples/advection/var-velocity.py +++ b/examples/advection/var-velocity.py @@ -75,9 +75,8 @@ def main(write_output=True, order=4): / source_width**2) def f_step(x): - from pymbolic.primitives import If, Comparison - return If( - Comparison( + return sym.If( + sym.Comparison( np.dot(sym_source_center_dist, sym_source_center_dist), "<", (4*source_width)**2), diff --git a/examples/wave/var-propagation-speed.py b/examples/wave/var-propagation-speed.py index 937b87a0..a7a634f7 100644 --- a/examples/wave/var-propagation-speed.py +++ b/examples/wave/var-propagation-speed.py @@ -51,10 +51,9 @@ def main(write_output=True, order=4): sym_x = sym.nodes(mesh.dim) sym_source_center_dist = sym_x - source_center sym_t = sym.ScalarVariable("t") - - from pymbolic.primitives import If, Comparison - c = If(Comparison(np.dot(sym_x, sym_x), "<", 0.15), - np.float32(-0.1), np.float32(-0.2)) + c = sym.If(sym.Comparison( + np.dot(sym_x, sym_x), "<", 0.15), + 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/advection.py b/grudge/models/advection.py index 5ef1e5ce..091b9870 100644 --- a/grudge/models/advection.py +++ b/grudge/models/advection.py @@ -59,10 +59,9 @@ class AdvectionOperatorBase(HyperbolicOperator): elif self.flux_type == "lf": return u.avg*v_dot_normal + 0.5*norm_v*(u.int - u.ext) elif self.flux_type == "upwind": - from pymbolic.primitives import If, Comparison return ( - v_dot_normal * If( - Comparison(v_dot_normal, ">", 0), + v_dot_normal * sym.If( + sym.Comparison(v_dot_normal, ">", 0), u.int, # outflow u.ext, # inflow )) @@ -156,10 +155,9 @@ class VariableCoefficientAdvectionOperator(HyperbolicOperator): elif self.flux_type == "lf": return u.avg*v_dot_normal + 0.5*norm_v*(u.int - u.ext) elif self.flux_type == "upwind": - from pymbolic.primitives import If, Comparison return ( - v_dot_normal * If( - Comparison(v_dot_normal, ">", 0), + v_dot_normal * sym.If( + sym.Comparison(v_dot_normal, ">", 0), u.int, # outflow u.ext, # inflow )) diff --git a/grudge/models/wave.py b/grudge/models/wave.py index 9b16c602..0ac60a33 100644 --- a/grudge/models/wave.py +++ b/grudge/models/wave.py @@ -350,9 +350,9 @@ class VariableCoefficientWeakWaveOperator(HyperbolicOperator): self.ambient_dim = ambient_dim self.source_f = source_f - from pymbolic.primitives import If, Comparison - self.sign = If(Comparison(self.c, ">", 0), - np.int32(1), np.int32(-1)) + self.sign = sym.If(sym.Comparison( + self.c, ">", 0), + np.int32(1), np.int32(-1)) self.dirichlet_tag = dirichlet_tag self.neumann_tag = neumann_tag -- GitLab