From 2bdcd938453a0d930ac4b329028a9c4337d45da5 Mon Sep 17 00:00:00 2001 From: Matt Wala Date: Sun, 11 Jun 2017 15:30:23 -0500 Subject: [PATCH] Extract a layer potential source superclass. This creates an abstract class LayerPotentialSourceBase with a few methods and properties and makes QBXLayerPotentialSource a subclass. The end goal is to share code and functionality with a Nystrom layer potential source (#33) which is not yet implemented. --- doc/discretization.rst | 13 +++-- pytential/qbx/__init__.py | 34 ++----------- pytential/source.py | 86 +++++++++++++++++++++++++++++++- pytential/symbolic/execution.py | 12 ++--- pytential/symbolic/mappers.py | 4 +- pytential/symbolic/primitives.py | 2 +- 6 files changed, 106 insertions(+), 45 deletions(-) diff --git a/doc/discretization.rst b/doc/discretization.rst index 6eb7a624..3229e267 100644 --- a/doc/discretization.rst +++ b/doc/discretization.rst @@ -1,7 +1,10 @@ Discretizations =============== -To compute a layer potential as an an end user, create a +QBX discretization +------------------ + +To compute a layer potential as an an end user, create a :class:`meshmode.discretization.Discretization` with a :class:`InterpolatoryQuadratureSimplexGroupFactory` as a discretization for the density. @@ -10,11 +13,13 @@ Then create :class:`pytential.qbx.QBXLayerPotentialSource`, :func:`pytential.bind` a layer potential operator to it, and you can start computing. -QBX discretization ------------------- - .. automodule:: pytential.qbx +Sources +------- + +.. automodule:: pytential.source + Targets ------- diff --git a/pytential/qbx/__init__.py b/pytential/qbx/__init__.py index 5292f7bb..2212b005 100644 --- a/pytential/qbx/__init__.py +++ b/pytential/qbx/__init__.py @@ -29,7 +29,7 @@ import numpy as np from pytools import memoize_method from meshmode.discretization import Discretization from pytential.qbx.target_assoc import QBXTargetAssociationFailedException -from pytential.source import PotentialSource +from pytential.source import LayerPotentialSourceBase import pyopencl as cl @@ -44,21 +44,13 @@ __doc__ = """ """ -class LayerPotentialSource(PotentialSource): - pass - - # {{{ QBX layer potential source -class QBXLayerPotentialSource(LayerPotentialSource): +class QBXLayerPotentialSource(LayerPotentialSourceBase): """A source discretization for a QBX layer potential. - .. attribute :: density_discr .. attribute :: qbx_order .. attribute :: fmm_order - .. attribute :: cl_context - .. automethod :: weights_and_area_elements - .. automethod :: with_refinement See :ref:`qbxguts` for some information on the inner workings of this. """ @@ -249,26 +241,6 @@ class QBXLayerPotentialSource(LayerPotentialSource): panel_sizes = self._panel_sizes("npanels").with_queue(queue) return np.asscalar(cl.array.max(panel_sizes).get()) - @property - def ambient_dim(self): - return self.density_discr.ambient_dim - - @property - def dim(self): - return self.density_discr.dim - - @property - def cl_context(self): - return self.density_discr.cl_context - - @property - def real_dtype(self): - return self.density_discr.real_dtype - - @property - def complex_dtype(self): - return self.density_discr.complex_dtype - # {{{ internal API @memoize_method @@ -455,7 +427,7 @@ class QBXLayerPotentialSource(LayerPotentialSource): len(target_discrs_and_qbx_sides) target_discr = bound_expr.places[o.target_name] - if isinstance(target_discr, LayerPotentialSource): + if isinstance(target_discr, LayerPotentialSourceBase): target_discr = target_discr.density_discr qbx_forced_limit = o.qbx_forced_limit diff --git a/pytential/source.py b/pytential/source.py index 87565060..f4dccda5 100644 --- a/pytential/source.py +++ b/pytential/source.py @@ -29,6 +29,13 @@ import six from pytools import memoize_method +__doc__ = """ +.. autoclass:: PotentialSource +.. autoclass:: PointPotentialSource +.. autoclass:: LayerPotentialSourceBase +""" + + class PotentialSource(object): """ .. method:: preprocess_optemplate(name, expr) @@ -43,12 +50,13 @@ class PotentialSource(object): """ +# {{{ point potential source + class PointPotentialSource(PotentialSource): """ ... attributes:: points An :class:`pyopencl.array.Array` of shape ``[ambient_dim, npoints]``. - """ def __init__(self, cl_context, points): @@ -130,3 +138,79 @@ class PointPotentialSource(PotentialSource): result.fill(1) return result.with_queue(None) + +# }}} + + +# {{{ layer potential source + +class LayerPotentialSourceBase(PotentialSource): + """A discretization of a layer potential using panel-based geometry, with + support for refinement and upsampling. + + .. rubric:: Discretizations + + .. attribute:: density_discr + .. attribute:: fine_density_discr + .. attribute:: resampler + .. method:: with_refinement + + .. rubric:: Discretization data + + .. attribute:: cl_context + .. attribute:: ambient_dim + .. attribute:: dim + .. attribute:: real_dtype + .. attribute:: complex_dtype + .. attribute:: h_max + + .. rubric:: Execution + + .. automethod:: weights_and_area_elements + .. method:: exec_compute_potential_insn + """ + + @property + def ambient_dim(self): + return self.density_discr.ambient_dim + + @property + def dim(self): + return self.density_discr.dim + + @property + def cl_context(self): + return self.density_discr.cl_context + + @property + def real_dtype(self): + return self.density_discr.real_dtype + + @property + def complex_dtype(self): + return self.density_discr.complex_dtype + + # {{{ weights and area elements + + @memoize_method + def weights_and_area_elements(self): + import pytential.symbolic.primitives as p + from pytential.symbolic.execution import bind + with cl.CommandQueue(self.cl_context) as queue: + # fine_density_discr is not guaranteed to be usable for + # interpolation/differentiation. Use density_discr to find + # area element instead, then upsample that. + + area_element = self.resampler(queue, + bind( + self.density_discr, + p.area_element(self.ambient_dim, self.dim) + )(queue)) + + qweight = bind(self.fine_density_discr, p.QWeight())(queue) + + return (area_element.with_queue(queue)*qweight).with_queue(None) + + # }}} + +# }}} diff --git a/pytential/symbolic/execution.py b/pytential/symbolic/execution.py index 5f67e808..8f30600a 100644 --- a/pytential/symbolic/execution.py +++ b/pytential/symbolic/execution.py @@ -246,8 +246,8 @@ class BoundExpression: def get_discretization(self, where): discr = self.places[where] - from pytential.qbx import LayerPotentialSource - if isinstance(discr, LayerPotentialSource): + from pytential.source import LayerPotentialSourceBase + if isinstance(discr, LayerPotentialSourceBase): discr = discr.density_discr return discr @@ -303,9 +303,9 @@ class BoundExpression: def prepare_places(places): from pytential.symbolic.primitives import DEFAULT_SOURCE, DEFAULT_TARGET from meshmode.discretization import Discretization - from pytential.qbx import LayerPotentialSource + from pytential.source import LayerPotentialSourceBase - if isinstance(places, LayerPotentialSource): + if isinstance(places, LayerPotentialSourceBase): places = { DEFAULT_SOURCE: places, DEFAULT_TARGET: places.density_discr, @@ -345,7 +345,7 @@ def prepare_expr(places, expr, auto_where=None): """ from pytential.symbolic.primitives import DEFAULT_SOURCE, DEFAULT_TARGET - from pytential.qbx import LayerPotentialSource + from pytential.source import LayerPotentialSourceBase from pytential.symbolic.mappers import ( ToTargetTagger, @@ -364,7 +364,7 @@ def prepare_expr(places, expr, auto_where=None): expr = DerivativeBinder()(expr) for name, place in six.iteritems(places): - if isinstance(place, LayerPotentialSource): + if isinstance(place, LayerPotentialSourceBase): expr = place.preprocess_optemplate(name, places, expr) return expr diff --git a/pytential/symbolic/mappers.py b/pytential/symbolic/mappers.py index 16679ac2..3eef6c4e 100644 --- a/pytential/symbolic/mappers.py +++ b/pytential/symbolic/mappers.py @@ -351,8 +351,8 @@ class QBXPreprocessor(IdentityMapper): source = self.places[self.source_name] target_discr = self.places[expr.target] - from pytential.qbx import LayerPotentialSource - if isinstance(target_discr, LayerPotentialSource): + from pytential.source import LayerPotentialSourceBase + if isinstance(target_discr, LayerPotentialSourceBase): target_discr = target_discr.density_discr if expr.qbx_forced_limit == 0: diff --git a/pytential/symbolic/primitives.py b/pytential/symbolic/primitives.py index c4c63b99..b6cf78a2 100644 --- a/pytential/symbolic/primitives.py +++ b/pytential/symbolic/primitives.py @@ -581,7 +581,7 @@ class IntG(Expression): to expressions that determine them) :arg source: The symbolic name of the source discretization. This name - is bound to a concrete :class:`pytential.qbx.QBXLayerPotentialSource` + is bound to a concrete :class:`pytential.source.LayerPotentialSourceBase` by :func:`pytential.bind`. :arg target: The symbolic name of the set of targets. This name gets -- GitLab