pymor.grids package

Submodules

_unstructured module

boundaryinfos module


class pymor.grids.boundaryinfos.AllDirichletBoundaryInfo(grid)[source]

Bases: pymor.grids.interfaces.BoundaryInfoInterface

BoundaryInfo where BoundaryType('dirichlet') is attached to each boundary entity.

mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.


class pymor.grids.boundaryinfos.BoundaryInfoFromIndicators(grid, indicators, assert_unique_type=None, assert_some_type=None)[source]

Bases: pymor.grids.interfaces.BoundaryInfoInterface

BoundaryInfo where the BoundaryTypes are determined by indicator functions.

Parameters

grid
The Grid to which the BoundaryInfo is associated.
indicators
Dict where each key is a BoundaryType and the corresponding value is a boolean valued function defined on the analytical domain which indicates if a point belongs to a boundary of the given BoundaryType (the indicator functions must be vectorized).
mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.


class pymor.grids.boundaryinfos.EmptyBoundaryInfo(grid)[source]

Bases: pymor.grids.interfaces.BoundaryInfoInterface

BoundaryInfo with no BoundaryTypes attached to any boundary.

mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.


class pymor.grids.boundaryinfos.SubGridBoundaryInfo(subgrid, grid, grid_boundary_info, new_boundary_type=None)[source]

Bases: pymor.grids.interfaces.BoundaryInfoInterface

Derives a BoundaryInfo for a SubGrid.

Parameters

subrid
The SubGrid for which a BoundaryInfo is created.
grid
The parent Grid.
grid_boundary_info
The BoundaryInfo of the parent Grid from which to derive the BoundaryInfo
new_boundary_type
The BoundaryType which is assigned to the new boundaries of subgrid. If None, no BoundaryType is assigned.
mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.

constructions module


pymor.grids.constructions.flatten_grid(grid)[source]

This method is used by our visualizers to render n-dimensional grids which cannot be embedded into R^n by duplicating vertices which would have to be mapped to multiple points at once (think of grids on rectangular domains with identified edges).

Parameters

grid
The Grid to flatten.

Returns

subentities
The subentities(0, grid.dim) relation for the flattened grid.
coordinates
The coordinates of the codim-grid.dim entities.
entity_map
Maps the indices of the codim-grid.dim entities of the flattened grid to the indices of the corresponding entities in the original grid.

defaultimpl module


class pymor.grids.defaultimpl.AffineGridDefaultImplementations[source]

Bases: object

Provides default implementations for AffineGrids.


class pymor.grids.defaultimpl.ConformalTopologicalGridDefaultImplementations[source]

Bases: object

Provides default informations for ConformalTopologicalGrids.


class pymor.grids.defaultimpl.ReferenceElementDefaultImplementations[source]

Bases: object

Provides default implementations for ReferenceElements.

gmsh module


class pymor.grids.gmsh.GmshBoundaryInfo(grid, sections)[source]

Bases: pymor.grids.interfaces.BoundaryInfoInterface

BoundaryInfo for a GmshGrid.

Parameters

grid
The corresponding GmshGrid.
sections
Parsed sections of the MSH-file as returned by load_gmsh.
mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.


class pymor.grids.gmsh.GmshGrid(sections)[source]

Bases: pymor.grids.unstructured.UnstructuredTriangleGrid

An UnstructuredTriangleGrid built from an existing Gmsh MSH-file.

Parameters

sections
Parsed sections of the MSH-file as returned by load_gmsh.
__str__() <==> str(x)[source]

pymor.grids.gmsh.load_gmsh(gmsh_file)[source]

Parse a Gmsh file and create a corresponding GmshGrid and GmshBoundaryInfo.

Parameters

gmsh_file
File handle of the Gmsh MSH-file.

Returns

grid
The generated GmshGrid.
boundary_info
The generated GmshBoundaryInfo.

interfaces module


class pymor.grids.interfaces.AffineGridInterface[source]

Bases: pymor.grids.defaultimpl.AffineGridDefaultImplementations, pymor.grids.interfaces.ConformalTopologicalGridInterface

Topological grid with geometry where each codim-0 entity is affinely mapped to the same ReferenceElement.

The grid is completely determined via the subentity relation given by subentities and the embeddings given by embeddings. In addition, only size and reference_element have to be implemented. Cached default implementations for all other methods are provided by AffineGridDefaultImplementations.

dim_outer

The dimension of the space into which the grid is embedded.

bounding_box()[source]

returns a (2, dim_outer)-shaped array containing lower/upper bounding box coordinates.

centers(codim)[source]

retval[e] is the barycenter of the codim-codim entity with global index e.

diameters(codim)[source]

retval[e] is the diameter of the codim-codim entity with global index e.

embeddings(codim)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

integration_elements(codim)[source]

retval[e] is given as sqrt(det(A^T*A)), where A = embeddings(codim)[0][e].

jacobian_inverse_transposed(codim)[source]

retval[e] is the transposed (pseudo-)inverse of the Jacobian of embeddings(codim)[e].

quadrature_points(codim, order=None, npoints=None, quadrature_type='default')[source]

retval[e] is an array of quadrature points in global coordinates for the codim-codim entity with global index e.

The quadrature is of order order or has npoints integration points. To integrate a function f over e one has to form

np.dot(f(quadrature_points(codim, order)[e]), reference_element(codim).quadrature(order)[1]) *
integration_elements(codim)[e].  # NOQA
reference_element(codim)[source]

The ReferenceElement of the codim-codim entities.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

unit_outer_normals()[source]

retval[e,i] is the unit outer normal to the i-th codim-1 subentity of the codim-0 entitiy with global index e.

volumes(codim)[source]

retval[e] is the (dim-codim)-dimensional volume of the codim-codim entity with global index e.

volumes_inverse(codim)[source]

retval[e] = 1 / volumes(codim)[e].


class pymor.grids.interfaces.AffineGridWithOrthogonalCentersInterface[source]

Bases: pymor.grids.interfaces.AffineGridInterface

AffineGrid with an additional orthogonal_centers method.

orthogonal_centers()[source]

retval[e] is a point inside the codim-0 entity with global index e such that the line segment from retval[e] to retval[e2] is always orthogonal to the codim-1 entity shared by the codim-0 entities with global index e and e2.

(This is mainly useful for gradient approximation in finite volume schemes.)


class pymor.grids.interfaces.BoundaryInfoInterface[source]

Bases: pymor.core.cache.CacheableInterface

Provides BoundaryTypes for the boundaries of a given ConformalTopologicalGrid.

For every BoundaryType and codimension a mask is provided, marking grid entities of the respective type and codimension by their global index.

boundary_types

set of all BoundaryTypes the grid has.

mask(boundary_type, codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to the BoundaryType boundary_type.

no_boundary_type_mask(codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to no BoundaryType.

unique_boundary_type_mask(codim)[source]

retval[i] is True if the codim-codim entity of global index i is associated to one and only one BoundaryType.


class pymor.grids.interfaces.ConformalTopologicalGridInterface[source]

Bases: pymor.grids.defaultimpl.ConformalTopologicalGridDefaultImplementations, pymor.core.cache.CacheableInterface

A topological grid without hanging nodes.

The grid is completely determined via the subentity relation given by subentities. In addition, only size has to be implemented, cached default implementations for all other methods are provided by ConformalTopologicalGridDefaultImplementations.

dim

The dimension of the grid.

boundaries(codim)[source]

Returns the global indices of all codim-codim boundary entities.

By definition, a codim-1 entity is a boundary entity if it has only one codim-0 superentity. For codim != 1, a codim-codim entity is a boundary entity if it has a codim-1 sub/super-entity.

boundary_mask(codim)[source]

retval[e] is true iff the codim-codim entity with global index e is a boundary entity.

By definition, a codim-1 entity is a boundary entity if it has only one codim-0 superentity. For codim != 1, a codim-codim entity is a boundary entity if it has a codim-1 sub/super-entity.

neighbours(codim, neighbour_codim, intersection_codim=None)[source]

retval[e,n] is the global index of the n-th codim-neighbour_codim entitiy of the codim-codim entity e that shares with e a subentity of codimension intersection_codim.

If intersection_codim == None, it is set to codim + 1 if codim == neighbour_codim and to min(codim, neighbour_codim) otherwise.

The default implementation is to compute the result from subentities(codim, intersection_codim) and superentities(intersection_codim, neihbour_codim).

size(codim)[source]

The number of entities of codimension codim.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

Only subentities(codim, codim+1) has to be implemented; a default implementation is provided which evaluates subentities(codim, subentity_codim) by computing the transitive closure of subentities(codim, codim+1).

superentities(codim, superentity_codim)[source]

retval[e,s] is the global index of the s-th codim-superentity_codim superentity of the codim-codim entity with global index e.

retval[e] is sorted by global index.

The default implementation is to compute the result from subentities(superentity_codim, codim).

superentity_indices(codim, superentity_codim)[source]

retval[e,s] is the local index of the codim-codim entity e in the codim-superentity_codim superentity superentities(codim, superentity_codim)[e,s].


class pymor.grids.interfaces.ReferenceElementInterface[source]

Bases: pymor.grids.defaultimpl.ReferenceElementDefaultImplementations, pymor.core.cache.CacheableInterface

Defines a reference element.

All reference elements have the property that all subentities of a given codimension are of the same type. I.e. a three-dimensional reference element cannot have triangles and rectangles as faces at the same time.

dim

The dimension of the reference element

volume

The volume of the reference element

__call__(codim)[source]

Returns the reference element of the codim-codim subentities.

center()[source]

Coordinates of the barycenter.

mapped_diameter(A)[source]

The diameter of the reference element after transforming it with the matrix A (vectorized).

quadrature(order=None, npoints=None, quadrature_type='default')[source]

Returns tuple (P, W) where P is an array of quadrature points with corresponding weights W.

The quadrature is of order order or has npoints integration points.

quadrature_info()[source]

Returns a tuple of dicts (O, N) where O[quadrature_type] is a list of orders which are implemented for quadrature_type and N[quadrature_type] is a list of the corresponding numbers of integration points.

size(codim)[source]

Number of subentities of codimension codim.

sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

subentities(codim, subentity_codim)[source]

subentities(c,sc)[i,j] is, with respect to the indexing inside the reference element, the index of the j-th codim-subentity_codim subentity of the i-th codim-codim subentity of the reference element.

subentity_embedding(subentity_codim)[source]

Returns a tuple (A, B) which defines the embedding of the codim-subentity_codim subentities into the reference element.

For subentity_codim > 1', the embedding is by default given recursively via `subentity_embedding(subentity_codim - 1) and sub_reference_element(subentity_codim - 1).subentity_embedding(1) choosing always the superentity with smallest index.

unit_outer_normals()[source]

retval[e] is the unit outer-normal vector to the codim-1 subentity with index e.

oned module


class pymor.grids.oned.OnedGrid(domain=(0, 1), num_intervals=4, identify_left_right=False)[source]

Bases: pymor.grids.interfaces.AffineGridWithOrthogonalCentersInterface

One-dimensional Grid on an interval.

Parameters

domain
Tuple (left, right) containing the left and right boundary of the domain.
num_intervals
The number of codim-0 entities.
__reduce__()[source]

helper for pickle

__str__() <==> str(x)[source]
bounding_box()[source]

returns a (2, dim_outer)-shaped array containing lower/upper bounding box coordinates.

embeddings(codim)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

orthogonal_centers()[source]

retval[e] is a point inside the codim-0 entity with global index e such that the line segment from retval[e] to retval[e2] is always orthogonal to the codim-1 entity shared by the codim-0 entities with global index e and e2.

(This is mainly useful for gradient approximation in finite volume schemes.)

size(codim=0)[source]

The number of entities of codimension codim.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

visualize(U, codim=2, **kwargs)[source]

Visualize scalar data associated to the grid as a patch plot.

Parameters

U
NumPy array of the data to visualize. If U.dim == 2 and len(U) > 1, the data is visualized as a time series of plots. Alternatively, a tuple of NumPy arrays can be provided, in which case a subplot is created for each entry of the tuple. The lengths of all arrays have to agree.
codim
The codimension of the entities the data in U is attached to (either 0 or 2).
kwargs
See visualize_patch

rect module


class pymor.grids.rect.RectGrid(num_intervals=(2, 2), domain=([0, 0], [1, 1]), identify_left_right=False, identify_bottom_top=False)[source]

Bases: pymor.grids.interfaces.AffineGridWithOrthogonalCentersInterface

Basic implementation of a rectangular Grid on a rectangular domain.

The global face, edge and vertex indices are given as follows

x1
^
|

6--10---7--11---8
|       |       |
3   2   4   3   5
|       |       |
3---8---4---9---5
|       |       |
0   0   1   1   2
|       |       |
0---6---1---7---2  --> x0

Parameters

num_intervals
Tuple (n0, n1) determining a grid with n0 x n1 codim-0 entities.
domain
Tuple (ll, ur) where ll defines the lower left and ur the upper right corner of the domain.
__reduce__()[source]

helper for pickle

__str__() <==> str(x)[source]
bounding_box()[source]

returns a (2, dim_outer)-shaped array containing lower/upper bounding box coordinates.

embeddings(codim=0)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

global_to_structured(codim)[source]

Returns an array which maps global codim-codim indices to structured indices.

I.e. if GTS = global_to_structured(codim) and STG = structured_to_global(codim), then STG[GTS[:, 0], GTS[:, 1]] == numpy.arange(size(codim)).

orthogonal_centers()[source]

retval[e] is a point inside the codim-0 entity with global index e such that the line segment from retval[e] to retval[e2] is always orthogonal to the codim-1 entity shared by the codim-0 entities with global index e and e2.

(This is mainly useful for gradient approximation in finite volume schemes.)

size(codim=0)[source]

The number of entities of codimension codim.

structured_to_global(codim)[source]

Returns an NumPy array which maps structured indices to global codim-codim indices.

In other words structured_to_global(codim)[i, j] is the global index of the i-th in x0-direction and j-th in x1-direction codim-codim entity of the grid.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

vertex_coordinates(dim)[source]

Returns an array of the x_dim coordinates of the grid vertices.

I.e.

centers(2)[structured_to_global(2)[i, j]] == np.array([vertex_coordinates(0)[i], vertex_coordinates(1)[j]])
visualize(U, codim=2, **kwargs)[source]

Visualize scalar data associated to the grid as a patch plot.

Parameters

U
NumPy array of the data to visualize. If U.dim == 2 and len(U) > 1, the data is visualized as a time series of plots. Alternatively, a tuple of NumPy arrays can be provided, in which case a subplot is created for each entry of the tuple. The lengths of all arrays have to agree.
codim
The codimension of the entities the data in U is attached to (either 0 or 2).
kwargs
See visualize_patch

referenceelements module


class pymor.grids.referenceelements.Line[source]

Bases: pymor.grids.interfaces.ReferenceElementInterface

center()[source]

Coordinates of the barycenter.

mapped_diameter(A)[source]

The diameter of the reference element after transforming it with the matrix A (vectorized).

quadrature(order=None, npoints=None, quadrature_type='default')[source]

Returns tuple (P, W) where P is an array of quadrature points with corresponding weights W.

The quadrature is of order order or has npoints integration points.

quadrature_info()[source]

Returns a tuple of dicts (O, N) where O[quadrature_type] is a list of orders which are implemented for quadrature_type and N[quadrature_type] is a list of the corresponding numbers of integration points.

size(codim)[source]

Number of subentities of codimension codim.

sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

subentities(codim, subentity_codim)[source]

subentities(c,sc)[i,j] is, with respect to the indexing inside the reference element, the index of the j-th codim-subentity_codim subentity of the i-th codim-codim subentity of the reference element.

subentity_embedding(subentity_codim)[source]

Returns a tuple (A, B) which defines the embedding of the codim-subentity_codim subentities into the reference element.

For subentity_codim > 1', the embedding is by default given recursively via `subentity_embedding(subentity_codim - 1) and sub_reference_element(subentity_codim - 1).subentity_embedding(1) choosing always the superentity with smallest index.

unit_outer_normals()[source]

retval[e] is the unit outer-normal vector to the codim-1 subentity with index e.


class pymor.grids.referenceelements.Point[source]

Bases: pymor.grids.interfaces.ReferenceElementInterface

center()[source]

Coordinates of the barycenter.

mapped_diameter(A)[source]

The diameter of the reference element after transforming it with the matrix A (vectorized).

quadrature(order=None, npoints=None, quadrature_type='default')[source]

Returns tuple (P, W) where P is an array of quadrature points with corresponding weights W.

The quadrature is of order order or has npoints integration points.

quadrature_info()[source]

Returns a tuple of dicts (O, N) where O[quadrature_type] is a list of orders which are implemented for quadrature_type and N[quadrature_type] is a list of the corresponding numbers of integration points.

size(codim)[source]

Number of subentities of codimension codim.

sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

subentities(codim, subentity_codim)[source]

subentities(c,sc)[i,j] is, with respect to the indexing inside the reference element, the index of the j-th codim-subentity_codim subentity of the i-th codim-codim subentity of the reference element.

subentity_embedding(subentity_codim)[source]

Returns a tuple (A, B) which defines the embedding of the codim-subentity_codim subentities into the reference element.

For subentity_codim > 1', the embedding is by default given recursively via `subentity_embedding(subentity_codim - 1) and sub_reference_element(subentity_codim - 1).subentity_embedding(1) choosing always the superentity with smallest index.

unit_outer_normals()[source]

retval[e] is the unit outer-normal vector to the codim-1 subentity with index e.


class pymor.grids.referenceelements.Square[source]

Bases: pymor.grids.interfaces.ReferenceElementInterface

center()[source]

Coordinates of the barycenter.

mapped_diameter(A)[source]

The diameter of the reference element after transforming it with the matrix A (vectorized).

quadrature(order=None, npoints=None, quadrature_type='default')[source]

Returns tuple (P, W) where P is an array of quadrature points with corresponding weights W.

The quadrature is of order order or has npoints integration points.

quadrature_info()[source]

Returns a tuple of dicts (O, N) where O[quadrature_type] is a list of orders which are implemented for quadrature_type and N[quadrature_type] is a list of the corresponding numbers of integration points.

size(codim)[source]

Number of subentities of codimension codim.

sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

subentities(codim, subentity_codim)[source]

subentities(c,sc)[i,j] is, with respect to the indexing inside the reference element, the index of the j-th codim-subentity_codim subentity of the i-th codim-codim subentity of the reference element.

subentity_embedding(subentity_codim)[source]

Returns a tuple (A, B) which defines the embedding of the codim-subentity_codim subentities into the reference element.

For subentity_codim > 1', the embedding is by default given recursively via `subentity_embedding(subentity_codim - 1) and sub_reference_element(subentity_codim - 1).subentity_embedding(1) choosing always the superentity with smallest index.

unit_outer_normals()[source]

retval[e] is the unit outer-normal vector to the codim-1 subentity with index e.


class pymor.grids.referenceelements.Triangle[source]

Bases: pymor.grids.interfaces.ReferenceElementInterface

center()[source]

Coordinates of the barycenter.

mapped_diameter(A)[source]

The diameter of the reference element after transforming it with the matrix A (vectorized).

quadrature(order=None, npoints=None, quadrature_type='default')[source]

Returns tuple (P, W) where P is an array of quadrature points with corresponding weights W.

The quadrature is of order order or has npoints integration points.

quadrature_info()[source]

Returns a tuple of dicts (O, N) where O[quadrature_type] is a list of orders which are implemented for quadrature_type and N[quadrature_type] is a list of the corresponding numbers of integration points.

size(codim)[source]

Number of subentities of codimension codim.

sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

subentities(codim, subentity_codim)[source]

subentities(c,sc)[i,j] is, with respect to the indexing inside the reference element, the index of the j-th codim-subentity_codim subentity of the i-th codim-codim subentity of the reference element.

subentity_embedding(subentity_codim)[source]

Returns a tuple (A, B) which defines the embedding of the codim-subentity_codim subentities into the reference element.

For subentity_codim > 1', the embedding is by default given recursively via `subentity_embedding(subentity_codim - 1) and sub_reference_element(subentity_codim - 1).subentity_embedding(1) choosing always the superentity with smallest index.

unit_outer_normals()[source]

retval[e] is the unit outer-normal vector to the codim-1 subentity with index e.

subgrid module


class pymor.grids.subgrid.SubGrid(grid, entities)[source]

Bases: pymor.grids.interfaces.AffineGridInterface

A subgrid of a Grid.

Given a Grid and a list of codim-0 entities we construct the minimal subgrid of the grid, containing all the given entities.

Parameters

grid
Grid of which a subgrid is to be created.
entities
NumPy array of global indices of the codim-0 entities which are to be contained in the subgrid.
parent_grid

The Grid from which the subgrid was constructed. Subgrid only stores a weakref to the grid, so accessing this property might return None if the original grid has been destroyed.

embeddings(codim)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

indices_from_parent_indices(ind, codim)[source]

Maps a NumPy array of indicies of codim-codim entites of the parent grid to indicies of the subgrid.

Raises

ValueError
Not all provided indices correspond to entities contained in the subgrid.
parent_indices(codim)[source]

retval[e] is the index of the e-th codim-codim entity in the parent grid.

size(codim)[source]

The number of entities of codimension codim.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

tria module


class pymor.grids.tria.TriaGrid(num_intervals=(2, 2), domain=([0, 0], [1, 1]), identify_left_right=False, identify_bottom_top=False)[source]

Bases: pymor.grids.interfaces.AffineGridWithOrthogonalCentersInterface

Basic implementation of a triangular grid on a rectangular domain.

The global face, edge and vertex indices are given as follows

6---------10----------7---------11----------8
| \                 / | \                 / |
|    22   10    18    |    23   11    19    |
|       \     /       |       \     /       |
3    14   11     6    4    15   12     7    5
|       /     \       |       /     \       |
|    14    2    26    |    15    3    27    |
| /                 \ | /                 \ |
3----------8----------4----------9----------5
| \                 / | \                 / |
|    20    8    16    |    21    9    17    |
|       \     /       |       \     /       |
0    12    9     4    1    13   10     5    2
|       /     \       |       /     \       |
|    12    0    24    |    13    1    25    |
| /                 \ | /                 \ |
0----------6----------1----------7----------2

Parameters

num_intervals
Tuple (n0, n1) determining a grid with n0 x n1 codim-0 entities.
domain
Tuple (ll, ur) where ll defines the lower left and ur the upper right corner of the domain.
__reduce__()[source]

helper for pickle

__str__() <==> str(x)[source]
bounding_box()[source]

returns a (2, dim_outer)-shaped array containing lower/upper bounding box coordinates.

embeddings(codim=0)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

size(codim=0)[source]

The number of entities of codimension codim.

subentities(codim, subentity_codim)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

visualize(U, codim=2, **kwargs)[source]

Visualize scalar data associated to the grid as a patch plot.

Parameters

U
NumPy array of the data to visualize. If U.dim == 2 and len(U) > 1, the data is visualized as a time series of plots. Alternatively, a tuple of NumPy arrays can be provided, in which case a subplot is created for each entry of the tuple. The lengths of all arrays have to agree.
codim
The codimension of the entities the data in U is attached to (either 0 or 2).
kwargs
See visualize_patch

unstructured module


class pymor.grids.unstructured.UnstructuredTriangleGrid(vertices, faces)[source]

Bases: pymor.grids.interfaces.AffineGridInterface

A generic unstructured, triangular grid.

Parameters

vertices
A (num_vertices, 2)-shaped NumPy array containing the coordinates of all vertices in the grid. The row numbers in the array will be the global indices of the given vertices (codim 2 entities).
faces
A (num_faces, 3)-shaped NumPy array containing the global indices of the vertices which define a given triangle in the grid. The row numbers in the array will be the global indices of the given triangles (codim 0 entities).
__str__() <==> str(x)[source]
embeddings(codim=0)[source]

Returns tuple (A, B) where A[e] and B[e] are the linear part and the translation part of the map from the reference element of e to e.

For codim > 0, we provide a default implementation by taking the embedding of the codim-1 parent entity e_0 of e with lowest global index and composing it with the subentity_embedding of e into e_0 determined by the reference element.

size(codim=0)[source]

The number of entities of codimension codim.

subentities(codim=0, subentity_codim=None)[source]

retval[e,s] is the global index of the s-th codim-subentity_codim subentity of the codim-codim entity with global index e.

The ordering of subentities(0, subentity_codim)[e] has to correspond, w.r.t. the embedding of e, to the local ordering inside the reference element.

For codim > 0, we provide a default implementation by calculating the subentities of e as follows:

  1. Find the codim-1 parent entity e_0 of e with minimal global index
  2. Lookup the local indices of the subentities of e inside e_0 using the reference element.
  3. Map these local indices to global indices using subentities(codim - 1, subentity_codim).

This procedures assures that subentities(codim, subentity_codim)[e] has the right ordering w.r.t. the embedding determined by e_0, which agrees with what is returned by embeddings(codim)

visualize(U, codim=2, **kwargs)[source]

Visualize scalar data associated to the grid as a patch plot.

Parameters

U
NumPy array of the data to visualize. If U.dim == 2 and len(U) > 1, the data is visualized as a time series of plots. Alternatively, a tuple of NumPy arrays can be provided, in which case a subplot is created for each entry of the tuple. The lengths of all arrays have to agree.
codim
The codimension of the entities the data in U is attached to (either 0 or 2).
kwargs
See visualize_patch