pymor.discretizers.builtin.grids package

Submodules

_unstructured module

boundaryinfos module


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

Bases: pymor.discretizers.builtin.grids.interfaces.BoundaryInfo

BoundaryInfo where the boundary type ‘dirichlet’ is attached to each boundary entity.

Methods

AllDirichletBoundaryInfo

mask

BoundaryInfo

boundaries, check_boundary_types, dirichlet_boundaries, dirichlet_mask, neumann_boundaries, neumann_mask, no_boundary_type_mask, robin_boundaries, robin_mask, unique_boundary_type_mask

CacheableObject

cached_method_call, disable_caching, enable_caching

ImmutableObject

with_, __setattr__

BasicObject

disable_logging, enable_logging

mask(boundary_type, codim)[source]

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


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

Bases: pymor.discretizers.builtin.grids.interfaces.BoundaryInfo

BoundaryInfo with no boundary types attached to any boundary.

Methods

EmptyBoundaryInfo

mask

BoundaryInfo

boundaries, check_boundary_types, dirichlet_boundaries, dirichlet_mask, neumann_boundaries, neumann_mask, no_boundary_type_mask, robin_boundaries, robin_mask, unique_boundary_type_mask

CacheableObject

cached_method_call, disable_caching, enable_caching

ImmutableObject

with_, __setattr__

BasicObject

disable_logging, enable_logging

mask(boundary_type, codim)[source]

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


class pymor.discretizers.builtin.grids.boundaryinfos.GenericBoundaryInfo(grid, masks, assert_unique_type=(1), assert_some_type=())[source]

Bases: pymor.discretizers.builtin.grids.interfaces.BoundaryInfo

Generic BoundaryInfo storing entity masks per boundary type.

Methods

GenericBoundaryInfo

from_indicators, mask

BoundaryInfo

boundaries, check_boundary_types, dirichlet_boundaries, dirichlet_mask, neumann_boundaries, neumann_mask, no_boundary_type_mask, robin_boundaries, robin_mask, unique_boundary_type_mask

CacheableObject

cached_method_call, disable_caching, enable_caching

ImmutableObject

with_, __setattr__

BasicObject

disable_logging, enable_logging

classmethod from_indicators(grid, indicators, assert_unique_type=None, assert_some_type=None)[source]

Create BoundaryInfo from indicator functions.

Parameters

grid

The Grid to which the BoundaryInfo is associated.

indicators

Dict where each key is a boundary type 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 boundary type (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 boundary type boundary_type.

constructions module


pymor.discretizers.builtin.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.

gmsh module


pymor.discretizers.builtin.grids.gmsh.load_gmsh(filename)[source]

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

Parameters

filename

Path of the Gmsh MSH-file.

Returns

grid

The generated GmshGrid.

boundary_info

The generated GmshBoundaryInfo.

interfaces module


class pymor.discretizers.builtin.grids.interfaces.BoundaryInfo[source]

Bases: pymor.core.cache.CacheableObject

Provides boundary types for the boundaries of a given Grid.

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

Methods

BoundaryInfo

boundaries, check_boundary_types, dirichlet_boundaries, dirichlet_mask, mask, neumann_boundaries, neumann_mask, no_boundary_type_mask, robin_boundaries, robin_mask, unique_boundary_type_mask

CacheableObject

cached_method_call, disable_caching, enable_caching

ImmutableObject

with_, __setattr__

BasicObject

disable_logging, enable_logging

boundary_types

set of all boundary types 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 boundary type 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 boundary type.

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 boundary type.


class pymor.discretizers.builtin.grids.interfaces.Grid[source]

Bases: pymor.core.cache.CacheableObject

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.

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.

bounding_box()[source]

returns a (2, dim)-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.

abstract 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].

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).

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
abstract reference_element(codim)[source]

The ReferenceElement of the codim-codim entities.

abstract size(codim)[source]

The number of entities of codimension codim.

abstract 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)

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].

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.discretizers.builtin.grids.interfaces.GridWithOrthogonalCenters[source]

Bases: pymor.discretizers.builtin.grids.interfaces.Grid

Grid with an additional orthogonal_centers method.

abstract 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.discretizers.builtin.grids.interfaces.ReferenceElement[source]

Bases: pymor.core.cache.CacheableObject

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.

abstract center()[source]

Coordinates of the barycenter.

abstract mapped_diameter(A)[source]

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

abstract 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.

abstract 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.

abstract size(codim)[source]

Number of subentities of codimension codim.

abstract sub_reference_element(codim)[source]

Returns the reference element of the codim-codim subentities.

abstract 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.

abstract 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.

abstract unit_outer_normals()[source]

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

oned module


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

Bases: pymor.discretizers.builtin.grids.interfaces.GridWithOrthogonalCenters

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__()[source]

Return str(self).

bounding_box()[source]

returns a (2, dim)-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.discretizers.builtin.grids.rect.RectGrid(num_intervals=(2, 2), domain=([0, 0], [1, 1]), identify_left_right=False, identify_bottom_top=False)[source]

Bases: pymor.discretizers.builtin.grids.interfaces.GridWithOrthogonalCenters

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.

identify_left_right

If True, the left and right boundaries are identified, i.e. the left-most codim-0 entities become neighbors of the right-most codim-0 entities.

identify_bottom_top

If True, the bottom and top boundaries are identified, i.e. the bottom-most codim-0 entities become neighbors of the top-most codim-0 entities.

__reduce__()[source]

Helper for pickle.

__str__()[source]

Return str(self).

bounding_box()[source]

returns a (2, dim)-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.discretizers.builtin.grids.referenceelements.Line[source]

Bases: pymor.discretizers.builtin.grids.interfaces.ReferenceElement

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.discretizers.builtin.grids.referenceelements.Point[source]

Bases: pymor.discretizers.builtin.grids.interfaces.ReferenceElement

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.discretizers.builtin.grids.referenceelements.Square[source]

Bases: pymor.discretizers.builtin.grids.interfaces.ReferenceElement

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.discretizers.builtin.grids.referenceelements.Triangle[source]

Bases: pymor.discretizers.builtin.grids.interfaces.ReferenceElement

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.discretizers.builtin.grids.subgrid.SubGrid(parent_grid, parent_entities)[source]

Bases: pymor.discretizers.builtin.grids.interfaces.Grid

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

parent_grid

Grid of which a subgrid is to be created.

parent_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)


pymor.discretizers.builtin.grids.subgrid.make_sub_grid_boundary_info(sub_grid, parent_grid, parent_grid_boundary_info, new_boundary_type=None)[source]

Derives a BoundaryInfo for a SubGrid.

Parameters

sub_grid

The SubGrid for which a BoundaryInfo is created.

parent_grid

The parent Grid.

parent_grid_boundary_info

The BoundaryInfo of the parent Grid from which to derive the BoundaryInfo

new_boundary_type

The boundary type which is assigned to the new boundaries of subgrid. If None, no boundary type is assigned.

Returns

BoundaryInfo associated with sub_grid.

tria module


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

Bases: pymor.discretizers.builtin.grids.interfaces.GridWithOrthogonalCenters

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.

identify_left_right

If True, the left and right boundaries are identified, i.e. the left-most codim-0 entities become neighbors of the right-most codim-0 entities.

identify_bottom_top

If True, the bottom and top boundaries are identified, i.e. the bottom-most codim-0 entities become neighbors of the top-most codim-0 entities.

__reduce__()[source]

Helper for pickle.

__str__()[source]

Return str(self).

bounding_box()[source]

returns a (2, dim)-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.

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

unstructured module


class pymor.discretizers.builtin.grids.unstructured.UnstructuredTriangleGrid(sizes, subentity_data, embedding_data)[source]

Bases: pymor.discretizers.builtin.grids.interfaces.Grid

A generic unstructured, triangular grid.

Use from_vertices to instantiate the grid from vertex coordinates and connectivity data.

__str__()[source]

Return str(self).

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.

classmethod from_vertices(vertices, faces)[source]

Instantiate grid from vertex coordinates and connectivity data.

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).

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

vtkio module


pymor.discretizers.builtin.grids.vtkio.write_vtk(grid, data, filename_base, codim=2, binary_vtk=True, last_step=None)[source]

Output grid-associated data in (legacy) vtk format

Parameters

grid

A Grid with triangular or rectilinear reference element.

data

VectorArray with either cell (ie one datapoint per codim 0 entity) or vertex (ie one datapoint per codim 2 entity) data in each array element.

codim

the codimension associated with the data

filename_base

common component for output files in timeseries

binary_vtk

if false, output files contain human readable inline ascii data, else appended binary

last_step

if set must be <= len(data) to restrict output of timeseries