diff --git a/meshmode/interop/FInAT/__init__.py b/meshmode/interop/FInAT/__init__.py deleted file mode 100644 index a79ad2ae3c12bc9b82afb3d17d15b26c81e906d7..0000000000000000000000000000000000000000 --- a/meshmode/interop/FInAT/__init__.py +++ /dev/null @@ -1,25 +0,0 @@ -__copyright__ = "Copyright (C) 2020 Benjamin Sepanski" - -__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 meshmode.interop.FInAT.lagrange_element import FinatLagrangeElementImporter - -__all__ = ['FinatLagrangeElementImporter'] diff --git a/meshmode/interop/FInAT/lagrange_element.py b/meshmode/interop/FInAT/lagrange_element.py deleted file mode 100644 index 4b3d977c04ba8a0db2517f78290e7b7147d275b0..0000000000000000000000000000000000000000 --- a/meshmode/interop/FInAT/lagrange_element.py +++ /dev/null @@ -1,218 +0,0 @@ -__copyright__ = "Copyright (C) 2020 Benjamin Sepanski" - -__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. -""" - -import numpy as np -import numpy.linalg as la - -from meshmode.interop import ExternalImportHandler -from meshmode.interop.fiat import FIATSimplexCellImporter - - -__doc__ = """ -.. autoclass:: FinatLagrangeElementImporter - :members: -""" - - -class FinatLagrangeElementImporter(ExternalImportHandler): - """ - An importer for a FInAT element, usually instantiated from - ``some_instantiated_firedrake_function_space.finat_element`` - """ - def __init__(self, finat_element): - """ - :param finat_element: A FInAT element of type - :class:`finat.fiat_elements.Lagrange` or - :class:`finat.fiat_elements.DiscontinuousLagrange` - (or :class:`finat.spectral.GaussLegendre` or - :class:`finat.spectral.GaussLobattoLegendre` in 1D) - which uses affine mapping of the basis functions - (i.e. ``finat_element.mapping`` must be - ``"affine"``) - - :raises TypeError: If :param:`finat_element` is not of type - :class:`finat.fiat_elements.Lagrange` or - :class:`finat.fiat_elements.DiscontinuousLagrange` - (or the 1D classes listed above) - - :raises ValueError: If :param:`finat_element` does not - use affine mappings of the basis functions - """ - # {{{ Check and store input - - # Check types - from finat.fiat_elements import DiscontinuousLagrange, Lagrange - from finat.spectral import GaussLegendre, GaussLobattoLegendre - valid_types = (Lagrange, DiscontinuousLagrange, - GaussLegendre, GaussLobattoLegendre) - if not isinstance(finat_element, valid_types): - raise TypeError(":param:`finat_element` must be of type" - " `finat.fiat_elements.Lagrange`," - " `finat.fiat_elements.DiscontinuousLagrange`," - " or `finat.spectral.GaussLegendre` or" - " `finat.spectral.GaussLobattoLegendre` in 1D," - " but is instead of type `%s`" % type(finat_element)) - - if finat_element.mapping != 'affine': - raise ValueError("FInAT element must use affine mappings" - " of the bases") - # }}} - - super(FinatLagrangeElementImporter, self).__init__(finat_element) - - self.cell_importer = FIATSimplexCellImporter(finat_element.cell) - - # computed and stored once :meth:`unit_nodes`, :meth:`unit_vertices`, - # and :meth:`flip_matrix` are called - self._unit_nodes = None - self._unit_vertex_indices = None - self._flip_matrix = None - - def _compute_unit_vertex_indices_and_nodes(self): - """ - Compute the unit nodes, as well as the unit vertex indices, - if they have not already been computed. - """ - if self._unit_nodes is None or self._unit_vertex_indices is None: - # FIXME : This should work, but uses some private info - # {{{ Compute unit nodes - - # each point evaluator is a function p(f) evaluates f at a node, - # so we need to evaluate each point evaluator at the identity to - # recover the nodes - point_evaluators = self.data._element.dual.nodes - unit_nodes = [p(lambda x: x) for p in point_evaluators] - unit_nodes = np.array(unit_nodes).T - self._unit_nodes = \ - self.cell_importer.affinely_map_firedrake_to_meshmode(unit_nodes) - - # Is this safe?, I think so bc on a reference element - close = 1e-8 - # Get vertices as (dim, nunit_vertices) - unit_vertices = np.array(self.data.cell.vertices).T - unit_vertices = \ - self.cell_importer.affinely_map_firedrake_to_meshmode(unit_vertices) - self._unit_vertex_indices = [] - for n_ndx in range(self._unit_nodes.shape[1]): - for v_ndx in range(unit_vertices.shape[1]): - diff = self._unit_nodes[:, n_ndx] - unit_vertices[:, v_ndx] - if np.max(np.abs(diff)) < close: - self._unit_vertex_indices.append(n_ndx) - break - - self._unit_vertex_indices = np.array(self._unit_vertex_indices) - - # }}} - - def dim(self): - """ - :return: The dimension of the FInAT element's cell - """ - return self.cell_importer.data.get_dimension() - - def unit_vertex_indices(self): - """ - :return: A numpy integer array of indices - so that *self.unit_nodes()[self.unit_vertex_indices()]* - are the nodes of the reference element which coincide - with its vertices (this is possibly empty). - """ - self._compute_unit_vertex_indices_and_nodes() - return self._unit_vertex_indices - - def unit_nodes(self): - """ - :return: The unit nodes used by the FInAT element mapped - onto the appropriate :mod:`modepy` `reference - element <https://documen.tician.de/modepy/nodes.html>`_ - as an array of shape *(dim, nunit_nodes)*. - """ - self._compute_unit_vertex_indices_and_nodes() - return self._unit_nodes - - def nunit_nodes(self): - """ - :return: The number of unit nodes. - """ - return self.unit_nodes().shape[1] - - def flip_matrix(self): - """ - :return: The matrix which should be applied to the - *(dim, nunitnodes)*-shaped array of nodes corresponding to - an element in order to change orientation - <-> +. - - The matrix will be *(dim, dim)* and orthogonal with - *np.float64* type entries. - """ - if self._flip_matrix is None: - # This is very similar to :mod:`meshmode` in processing.py - # the function :function:`from_simplex_element_group`, but - # we needed to use firedrake nodes - - from modepy.tools import barycentric_to_unit, unit_to_barycentric - - # Generate a resampling matrix that corresponds to the - # first two barycentric coordinates being swapped. - - bary_unit_nodes = unit_to_barycentric(self.unit_nodes()) - - flipped_bary_unit_nodes = bary_unit_nodes.copy() - flipped_bary_unit_nodes[0, :] = bary_unit_nodes[1, :] - flipped_bary_unit_nodes[1, :] = bary_unit_nodes[0, :] - flipped_unit_nodes = barycentric_to_unit(flipped_bary_unit_nodes) - - from modepy import resampling_matrix, simplex_best_available_basis - - flip_matrix = resampling_matrix( - simplex_best_available_basis(self.dim(), self.data.degree), - flipped_unit_nodes, self.unit_nodes()) - - flip_matrix[np.abs(flip_matrix) < 1e-15] = 0 - - # Flipping twice should be the identity - assert la.norm( - np.dot(flip_matrix, flip_matrix) - - np.eye(len(flip_matrix))) < 1e-13 - - self._flip_matrix = flip_matrix - - return self._flip_matrix - - def make_resampling_matrix(self, element_grp): - """ - :param element_grp: A - :class:`meshmode.discretization.InterpolatoryElementGroupBase` whose - basis functions span the same space as the FInAT element. - :return: A matrix which resamples a function sampled at - the firedrake unit nodes to a function sampled at - *element_grp.unit_nodes()* (by matrix multiplication) - """ - from meshmode.discretization import InterpolatoryElementGroupBase - assert isinstance(element_grp, InterpolatoryElementGroupBase), \ - "element group must be an interpolatory element group so that" \ - " can redistribute onto its nodes" - - from modepy import resampling_matrix - return resampling_matrix(element_grp.basis(), - new_nodes=element_grp.unit_nodes, - old_nodes=self.unit_nodes()) diff --git a/meshmode/interop/__init__.py b/meshmode/interop/__init__.py index 2c1bf0727a0becb3a13d9d6dbad12b751d60db39..ba465860a4545ee4aeeaec84c9df953efd10db76 100644 --- a/meshmode/interop/__init__.py +++ b/meshmode/interop/__init__.py @@ -19,54 +19,3 @@ 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. """ - -__doc__ = """ -Development Interface ---------------------- -.. autoclass:: ExternalDataHandler -.. autoclass:: ExternalExportHandler -.. autoclass:: ExternalImportHandler -""" - - -# {{{ Generic, most abstract class for transporting meshmode <-> external - -class ExternalDataHandler: - """ - A data handler takes data from meshmode and facilitates its use - in another package or the reverse: takes data from another package - and facilitates its use in meshmode. - - .. attribute:: data - - The object which needs to be interfaced either into meshmode or - out of meshmode. - Should not be modified after creation. - """ - def __init__(self, data): - self.data = data - -# }}} - - -# {{{ Define specific classes for meshmode -> external and meshmode <- external - -class ExternalExportHandler(ExternalDataHandler): - """ - Subclass of :class:`ExternalDataHandler` for meshmode -> external - data transfer - """ - pass - - -class ExternalImportHandler(ExternalDataHandler): - """ - Subclass of :class:`ExternalDataHandler` for external -> meshmode - data transfer - """ - pass - -# }}} - - -# vim: fdm=marker diff --git a/meshmode/interop/fiat/__init__.py b/meshmode/interop/fiat/__init__.py deleted file mode 100644 index 57a42bfe47fb1ed78cb8be603326b222761add2f..0000000000000000000000000000000000000000 --- a/meshmode/interop/fiat/__init__.py +++ /dev/null @@ -1,25 +0,0 @@ -__copyright__ = "Copyright (C) 2020 Benjamin Sepanski" - -__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 meshmode.interop.fiat.simplex_cell import FIATSimplexCellImporter - -__all__ = ['FIATSimplexCellImporter'] diff --git a/meshmode/interop/fiat/simplex_cell.py b/meshmode/interop/fiat/simplex_cell.py deleted file mode 100644 index c6c0ca2ca4e4241a0baa99860a083214645b4de5..0000000000000000000000000000000000000000 --- a/meshmode/interop/fiat/simplex_cell.py +++ /dev/null @@ -1,122 +0,0 @@ -__copyright__ = "Copyright (C) 2020 Benjamin Sepanski" - -__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. -""" - -import numpy as np -import numpy.linalg as la - -from meshmode.interop import ExternalImportHandler - - -__doc__ = """ -.. autoclass:: FIATSimplexCellImporter -""" - - -# {{{ Compute an affine mapping from given input/outputs - -def get_affine_mapping(reference_vects, vects): - r""" - Returns (mat, shift), - a matrix *mat* and vector *shift* which maps the - *i*th vector in *reference_vects* to the *i*th vector in *vects by - - ..math:: - - A ri + b -> vi, \qquad A = mat, b = shift - - :arg reference_vects: An np.array of *n* vectors of dimension *ref_dim* - :arg vects: An np.array of *n* vectors of dimension *dim*, with - *ref_dim* <= *dim*. - - NOTE : Should be shape *(ref_dim, nvectors)*, *(dim, nvectors)* respectively. - - *mat* will have shape *(dim, ref_dim)*, *shift* will have shape *(dim,)* - """ - # Make sure both have same number of vectors - ref_dim, num_vects = reference_vects.shape - assert num_vects == vects.shape[1] - - # Make sure d1 <= d2 (see docstring) - dim = vects.shape[0] - assert ref_dim <= dim - - # If there is only one vector, set M = I, b = vect - reference - if num_vects == 1: - mat = np.eye(dim, ref_dim) - shift = vects[:, 0] - np.matmul(mat, reference_vects[:, 0]) - else: - ref_span_vects = reference_vects[:, 1:] - reference_vects[:, 0, np.newaxis] - span_vects = vects[:, 1:] - vects[:, 0, np.newaxis] - mat = la.solve(ref_span_vects, span_vects) - shift = -np.matmul(mat, reference_vects[:, 0]) + vects[:, 0] - - return mat, shift - -# }}} - - -# {{{ Interoperator for FIAT's reference_element.Simplex - -class FIATSimplexCellImporter(ExternalImportHandler): - """ - Import handler for a :mod:`FIAT` simplex cell. - - .. attribute:: data - - An instance of :class:`fiat.FIAT.reference_element.Simplex`. - """ - def __init__(self, cell): - """ - :arg cell: a :class:`fiat.FIAT.reference_element.Simplex`. - """ - # Ensure this cell is actually a simplex - from FIAT.reference_element import Simplex - assert isinstance(cell, Simplex) - - super(FIATSimplexCellImporter, self).__init__(cell) - - # Stored as (dim, nunit_vertices) - from modepy.tools import unit_vertices - self._unit_vertices = unit_vertices(cell.get_dimension()).T - - # Maps FIAT reference vertices to :mod:`meshmode` - # unit vertices by x -> Ax + b, where A is :attr:`_mat` - # and b is :attr:`_shift` - reference_vertices = np.array(cell.vertices).T - self._mat, self._shift = get_affine_mapping(reference_vertices, - self._unit_vertices) - - def affinely_map_firedrake_to_meshmode(self, points): - """ - Map points on the firedrake reference simplex to - :mod:`modepy` - `unit coordinates <https://documen.tician.de/modepy/nodes.html>`_. - - :arg points: *n* points on the reference simplex - as a numpy array of shape *(dim, n)* - :return: A numpy array of shape *(dim, n)* wherein the - firedrake refernece simplex has been affinely mapped - to the modepy reference simplex - """ - return np.matmul(self._mat, points) + self._shift[:, np.newaxis] - -# }}}