From f8c1fdd9a99451df74cddb64b662fa209f13696e Mon Sep 17 00:00:00 2001
From: Matt Wala <wala1@illinois.edu>
Date: Thu, 28 Jul 2016 14:47:52 -0500
Subject: [PATCH] Split out make_refinement_connection() into its own
submodule, and add it to the top level documentation.
---
.../discretization/connection/__init__.py | 159 +-------------
.../discretization/connection/refinement.py | 195 ++++++++++++++++++
2 files changed, 202 insertions(+), 152 deletions(-)
create mode 100644 meshmode/discretization/connection/refinement.py
diff --git a/meshmode/discretization/connection/__init__.py b/meshmode/discretization/connection/__init__.py
index 80a9cfe3..ba158358 100644
--- a/meshmode/discretization/connection/__init__.py
+++ b/meshmode/discretization/connection/__init__.py
@@ -1,4 +1,3 @@
-# -*- coding: utf-8 -*-
from __future__ import division, print_function, absolute_import
__copyright__ = "Copyright (C) 2014 Andreas Kloeckner"
@@ -37,6 +36,9 @@ from meshmode.discretization.connection.face import (
make_face_restriction, make_face_to_all_faces_embedding)
from meshmode.discretization.connection.opposite_face import \
make_opposite_face_connection
+from meshmode.discretization.connection.refinement import \
+ make_refinement_connection
+
import logging
logger = logging.getLogger(__name__)
@@ -48,7 +50,8 @@ __all__ = [
"FRESTR_INTERIOR_FACES", "FRESTR_ALL_FACES",
"make_face_restriction",
"make_face_to_all_faces_embedding",
- "make_opposite_face_connection"
+ "make_opposite_face_connection",
+ "make_refinement_connection"
]
__doc__ = """
@@ -63,6 +66,8 @@ __doc__ = """
.. autofunction:: make_opposite_face_connection
+.. autofunction:: make_refinement_connection
+
Implementation details
^^^^^^^^^^^^^^^^^^^^^^
@@ -435,154 +440,4 @@ def check_connection(connection):
# }}}
-# {{{ refinement connection
-
-
-def _map_unit_nodes_to_children(unit_nodes, tesselation):
- """
- Given a collection of unit nodes, return the coordinates of the
- unit nodes mapped onto each of the children of the reference
- element.
-
- The tesselation should follow the format of
- :func:`meshmode.mesh.tesselate.tesselatetri()` or
- :func:`meshmode.mesh.tesselate.tesselatetet()`.
-
- `unit_nodes` should be relative to the unit simplex coordinates in
- :module:`modepy`.
-
- :arg unit_nodes: shaped `(dim, nunit_nodes)`
- :arg tesselation: With attributes `ref_vertices`, `children`
- """
- ref_vertices = np.array(tesselation.ref_vertices, dtype=np.float)
-
- for child_element in tesselation.children:
- center = np.vstack(ref_vertices[child_element[0]])
- # Scale by 1/2 since sides in the tesselation have length 2.
- aff_mat = (ref_vertices.T[:, child_element[1:]] - center) / 2
- # (-1, -1, ...) in unit_nodes = (0, 0, ...) in ref_vertices.
- # Hence the translation by +/- 1.
- yield aff_mat.dot(unit_nodes + 1) + center - 1
-
-
-def _build_interpolation_batches_for_group(
- queue, group_idx, coarse_discr_group, fine_discr_group, record):
- """
- To map between discretizations, we sort each of the fine mesh
- elements into an interpolation batch. Which batch they go
- into is determined by where the refined unit nodes live
- relative to the coarse reference element.
-
- For instance, consider the following refinement:
-
- ______ ______
- |\ | |\ e|
- | \ | |d\ |
- | \ | |__\ |
- | \ | => |\c|\ |
- | \ | |a\|b\ |
- | \| | | \|
- ‾‾‾‾‾‾ ‾‾‾‾‾‾
-
- Here, the discretization unit nodes for elements a,b,c,d,e
- will each have different positions relative to the reference
- element, so each element gets its own batch. On the other
- hand, for
-
- ______ ______
- |\ | |\ f|\e|
- | \ | |d\ |g\|
- | \ | |__\|__|
- | \ | => |\c|\ |
- | \ | |a\|b\h|
- | \| | | \|
- ‾‾‾‾‾‾ ‾‾‾‾‾‾
-
- the pairs {a,e}, {b,f}, {c,g}, {d,h} can share interpolation
- batches because their unit nodes are mapped from the same part
- of the reference element.
- """
- num_children = len(record.tesselation.children) \
- if record.tesselation else 0
- from_bins = [[] for i in range(1 + num_children)]
- to_bins = [[] for i in range(1 + num_children)]
- for elt_idx, refinement_result in enumerate(record.element_mapping):
- if len(refinement_result) == 1:
- # Not refined -> interpolates to self
- from_bins[0].append(elt_idx)
- to_bins[0].append(refinement_result[0])
- else:
- assert len(refinement_result) == num_children
- # Refined -> interpolates to children
- for from_bin, to_bin, child_idx in zip(
- from_bins[1:], to_bins[1:], refinement_result):
- from_bin.append(elt_idx)
- to_bin.append(child_idx)
-
- fine_unit_nodes = fine_discr_group.unit_nodes
- mapped_unit_nodes = _map_unit_nodes_to_children(
- fine_unit_nodes, record.tesselation)
-
- from itertools import chain
- for from_bin, to_bin, unit_nodes in zip(
- from_bins, to_bins,
- chain([fine_unit_nodes], mapped_unit_nodes)):
- if not from_bin:
- continue
- yield InterpolationBatch(
- from_group_index=group_idx,
- from_element_indices=cl.array.to_device(queue, np.asarray(from_bin)),
- to_element_indices=cl.array.to_device(queue, np.asarray(to_bin)),
- result_unit_nodes=unit_nodes,
- to_element_face=None)
-
-
-def make_refinement_connection(refiner, coarse_discr, group_factory):
- """
- :arg refiner: An instance of :class:`meshmode.mesh.refinement.Refiner`
-
- :arg coarse_discr: An instance of
- :class:`meshmode.mesh.discretization.Discretization`
-
- :arg group_factory: An instance of
- :class:`meshmode.mesh.discretization.ElementGroupFactory`. Used for
- discretizing the fine mesh.
-
- :return: A :class:`DiscretizationConnection` mapping `coarse_discr` to a
- discretization on the fine mesh
- """
- coarse_mesh = refiner.get_previous_mesh()
- fine_mesh = refiner.last_mesh
- assert coarse_discr.mesh is coarse_mesh
-
- from meshmode.discretization import Discretization
- fine_discr = Discretization(
- coarse_discr.cl_context,
- fine_mesh,
- group_factory,
- real_dtype=coarse_discr.real_dtype)
-
- logger.info("building refinement connection: start")
-
- groups = []
- with cl.CommandQueue(fine_discr.cl_context) as queue:
- for group_idx, (coarse_discr_group, fine_discr_group, record) in \
- enumerate(zip(coarse_discr.groups, fine_discr.groups,
- refiner.group_refinement_records)):
- groups.append(
- DiscretizationConnectionElementGroup(
- list(_build_interpolation_batches_for_group(
- queue, group_idx, coarse_discr_group,
- fine_discr_group, record))))
-
- logger.info("building refinement connection: done")
-
- return DiscretizationConnection(
- from_discr=coarse_discr,
- to_discr=fine_discr,
- groups=groups,
- is_surjective=True)
-
-# }}}
-
# vim: foldmethod=marker
diff --git a/meshmode/discretization/connection/refinement.py b/meshmode/discretization/connection/refinement.py
new file mode 100644
index 00000000..75604ec0
--- /dev/null
+++ b/meshmode/discretization/connection/refinement.py
@@ -0,0 +1,195 @@
+# -*- coding: utf-8 -*-
+from __future__ import division, print_function, absolute_import
+
+__copyright__ = "Copyright (C) 2016 Matt Wala"
+
+__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 pyopencl as cl
+import pyopencl.array # noqa
+
+import logging
+logger = logging.getLogger(__name__)
+
+
+# {{{ Map unit nodes to children
+
+def _map_unit_nodes_to_children(unit_nodes, tesselation):
+ """
+ Given a collection of unit nodes, return the coordinates of the
+ unit nodes mapped onto each of the children of the reference
+ element.
+
+ The tesselation should follow the format of
+ :func:`meshmode.mesh.tesselate.tesselatetri()` or
+ :func:`meshmode.mesh.tesselate.tesselatetet()`.
+
+ `unit_nodes` should be relative to the unit simplex coordinates in
+ :module:`modepy`.
+
+ :arg unit_nodes: shaped `(dim, nunit_nodes)`
+ :arg tesselation: With attributes `ref_vertices`, `children`
+ """
+ ref_vertices = np.array(tesselation.ref_vertices, dtype=np.float)
+
+ for child_element in tesselation.children:
+ center = np.vstack(ref_vertices[child_element[0]])
+ # Scale by 1/2 since sides in the tesselation have length 2.
+ aff_mat = (ref_vertices.T[:, child_element[1:]] - center) / 2
+ # (-1, -1, ...) in unit_nodes = (0, 0, ...) in ref_vertices.
+ # Hence the translation by +/- 1.
+ yield aff_mat.dot(unit_nodes + 1) + center - 1
+
+# }}}
+
+
+# {{{ Build interpolation batches for group
+
+def _build_interpolation_batches_for_group(
+ queue, group_idx, coarse_discr_group, fine_discr_group, record):
+ """
+ To map between discretizations, we sort each of the fine mesh
+ elements into an interpolation batch. Which batch they go
+ into is determined by where the refined unit nodes live
+ relative to the coarse reference element.
+
+ For instance, consider the following refinement:
+
+ ______ ______
+ |\ | |\ e|
+ | \ | |d\ |
+ | \ | |__\ |
+ | \ | => |\c|\ |
+ | \ | |a\|b\ |
+ | \| | | \|
+ ‾‾‾‾‾‾ ‾‾‾‾‾‾
+
+ Here, the discretization unit nodes for elements a,b,c,d,e
+ will each have different positions relative to the reference
+ element, so each element gets its own batch. On the other
+ hand, for
+
+ ______ ______
+ |\ | |\ f|\e|
+ | \ | |d\ |g\|
+ | \ | |__\|__|
+ | \ | => |\c|\ |
+ | \ | |a\|b\h|
+ | \| | | \|
+ ‾‾‾‾‾‾ ‾‾‾‾‾‾
+
+ the pairs {a,e}, {b,f}, {c,g}, {d,h} can share interpolation
+ batches because their unit nodes are mapped from the same part
+ of the reference element.
+ """
+ from meshmode.discretization.connection import InterpolationBatch
+
+ num_children = len(record.tesselation.children) \
+ if record.tesselation else 0
+ from_bins = [[] for i in range(1 + num_children)]
+ to_bins = [[] for i in range(1 + num_children)]
+ for elt_idx, refinement_result in enumerate(record.element_mapping):
+ if len(refinement_result) == 1:
+ # Not refined -> interpolates to self
+ from_bins[0].append(elt_idx)
+ to_bins[0].append(refinement_result[0])
+ else:
+ assert len(refinement_result) == num_children
+ # Refined -> interpolates to children
+ for from_bin, to_bin, child_idx in zip(
+ from_bins[1:], to_bins[1:], refinement_result):
+ from_bin.append(elt_idx)
+ to_bin.append(child_idx)
+
+ fine_unit_nodes = fine_discr_group.unit_nodes
+ mapped_unit_nodes = _map_unit_nodes_to_children(
+ fine_unit_nodes, record.tesselation)
+
+ from itertools import chain
+ for from_bin, to_bin, unit_nodes in zip(
+ from_bins, to_bins,
+ chain([fine_unit_nodes], mapped_unit_nodes)):
+ if not from_bin:
+ continue
+ yield InterpolationBatch(
+ from_group_index=group_idx,
+ from_element_indices=cl.array.to_device(queue, np.asarray(from_bin)),
+ to_element_indices=cl.array.to_device(queue, np.asarray(to_bin)),
+ result_unit_nodes=unit_nodes,
+ to_element_face=None)
+
+# }}}
+
+
+def make_refinement_connection(refiner, coarse_discr, group_factory):
+ """Return a
+ :class:`meshmode.discretization.connection.DiscretizationConnection`
+ connecting `coarse_discr` to a discretization on the fine mesh.
+
+ :arg refiner: An instance of
+ :class:`meshmode.mesh.refinement.Refiner`
+
+ :arg coarse_discr: An instance of
+ :class:`meshmode.mesh.discretization.Discretization` associated
+ with the mesh given to the refiner
+
+ :arg group_factory: An instance of
+ :class:`meshmode.mesh.discretization.ElementGroupFactory`. Used for
+ discretizing the fine mesh.
+ """
+ from meshmode.discretization.connection import (
+ DiscretizationConnectionElementGroup, DiscretizationConnection)
+
+ coarse_mesh = refiner.get_previous_mesh()
+ fine_mesh = refiner.last_mesh
+ assert coarse_discr.mesh is coarse_mesh
+
+ from meshmode.discretization import Discretization
+ fine_discr = Discretization(
+ coarse_discr.cl_context,
+ fine_mesh,
+ group_factory,
+ real_dtype=coarse_discr.real_dtype)
+
+ logger.info("building refinement connection: start")
+
+ groups = []
+ with cl.CommandQueue(fine_discr.cl_context) as queue:
+ for group_idx, (coarse_discr_group, fine_discr_group, record) in \
+ enumerate(zip(coarse_discr.groups, fine_discr.groups,
+ refiner.group_refinement_records)):
+ groups.append(
+ DiscretizationConnectionElementGroup(
+ list(_build_interpolation_batches_for_group(
+ queue, group_idx, coarse_discr_group,
+ fine_discr_group, record))))
+
+ logger.info("building refinement connection: done")
+
+ return DiscretizationConnection(
+ from_discr=coarse_discr,
+ to_discr=fine_discr,
+ groups=groups,
+ is_surjective=True)
+
+
+# vim: foldmethod=marker
--
GitLab