firedrake.cython package

Submodules

firedrake.cython.dmcommon module

firedrake.cython.dmcommon.cell_facet_labeling()

Computes a labeling for the facet numbers on a particular cell (interior and exterior facet labels with subdomain markers). The i-th local facet is represented as:

cell_facets[c, i]

If cell_facets[c, i, 0] is 0, then the local facet i is an exterior facet, otherwise if the result is 1 it is interior. cell_facets[c, i, 1] returns the subdomain marker for the local facet.

Parameters:
  • plex – The DMPlex object representing the mesh topology.

  • cell_numbering – PETSc.Section describing the global cell numbering

  • cell_closures – 2D array of ordered cell closures.

firedrake.cython.dmcommon.clear_adjacency_callback()

Clear the callback for DMPlexGetAdjacency.

Parameters:

dm – The DMPlex object

firedrake.cython.dmcommon.closure_ordering()

Apply Fenics local numbering to a cell closure.

Parameters:
  • dm – The DM object encapsulating the mesh topology

  • vertex_numbering – Section describing the universal vertex numbering

  • cell_numbering – Section describing the global cell numbering

  • entity_per_cell – List of the number of entity points in each dimension

Vertices := Ordered according to global/universal

vertex numbering

Edges/faces := Ordered according to lexicographical

ordering of non-incident vertices

firedrake.cython.dmcommon.complete_facet_labels()

Transfer label values from the facet labels to everything in the closure of the facets.

firedrake.cython.dmcommon.compute_point_cone_global_sizes()

Compute the total number of DMPlex points and the global sum of cone sizes for dm.

Parameters:

dm – The dm.

Returns:

A list of global number of points and global sum of cone sizes.

firedrake.cython.dmcommon.count_labelled_points()

Return the number of points in the chart [start, end) that are marked by the given label.

Note

This destroys any index on the label.

Parameters:
  • dm – The DM containing the label

  • name – The label name.

  • start – The smallest point to consider.

  • end – One past the largest point to consider.

firedrake.cython.dmcommon.create_cell_closure()

Create a map from FIAT local entity numbers to DMPlex point numbers for each cell.

Parameters:
  • dm – The DM object encapsulating the mesh topology

  • cell_numbering – Section describing the global cell numbering

  • _closureSize – Number of entities in the cell

firedrake.cython.dmcommon.create_section()

Create the section describing a global numbering.

Parameters:
  • mesh – The mesh.

  • nodes_per_entity – Number of nodes on each type of topological entity of the mesh. Or, if the mesh is extruded, the number of nodes on, and on top of, each topological entity in the base mesh.

  • on_base – If True, assume extruded space is actually Foo x Real.

  • block_size – The integer by which nodes_per_entity is uniformly multiplied to get the true data layout.

Returns:

A PETSc Section providing the number of dofs, and offset of each dof, on each mesh point.

firedrake.cython.dmcommon.entity_orientations()

Compute entity orientations.

Parameters:
  • mesh – The MeshTopology object encapsulating the mesh topology

  • cell_closure – The two-dimensional array, each row of which contains the closure of the associated cell

Returns:

A 2D array of the same shape as cell_closure, each row of which contains orientations of the entities in the closure of the associated cell

See entity_orientations() for details on the returned array.

See get_cell_nodes() for the usage of the returned array.

firedrake.cython.dmcommon.exchange_cell_orientations()

Halo exchange of cell orientations.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • section – Section describing the cell numbering

  • orientations – Cell orientations to exchange, values in the halo will be overwritten.

firedrake.cython.dmcommon.facet_closure_nodes()

Extract nodes in the closure of facets with a given marker.

This works fine for interior as well as exterior facets.

Parameters:
  • V – the function space

  • sub_domain – a tuple of mesh markers selecting the facets, or the magic string “on_boundary” indicating the entire boundary.

Returns:

a numpy array of unique nodes in the closure of facets with the given marker.

firedrake.cython.dmcommon.facet_numbering()

Compute the parent cell(s) and the local facet number within each parent cell for each given facet.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • kind – String indicating the facet kind (interior or exterior)

  • facets – Array of input facets

  • cell_numbering – Section describing the global cell numbering

  • cell_closures – 2D array of ordered cell closures

firedrake.cython.dmcommon.get_cell_markers()

Get the cells marked by a given subdomain_id.

Parameters:
  • dm – The DM for the mesh topology

  • cell_numbering – Section mapping dm cell points to firedrake cell indices.

  • subdomain_id – The subdomain_id to look for.

Raises:

ValueError – if the subdomain_id is not valid.

Returns:

A numpy array (possibly empty) of the cell ids.

firedrake.cython.dmcommon.get_cell_nodes()

Builds the DoF mapping.

Parameters:
  • mesh – The mesh

  • global_numbering – Section describing the global DoF numbering

  • entity_dofs – FInAT element entity dofs for the cell

  • entity_permutations – FInAT element entity permutations for the cell

  • offset – offsets for each entity dof walking up a column.

Preconditions: This function assumes that cell_closures contains mesh entities ordered by dimension, i.e. vertices first, then edges, faces, and finally the cell. For quadrilateral meshes, edges corresponding to dimension (0, 1) in the FInAT element must precede edges corresponding to dimension (1, 0) in the FInAT element.

firedrake.cython.dmcommon.get_cell_remote_ranks()

Returns an array assigning the rank of the owner to each locally visible cell. Locally owned cells have -1 assigned to them.

Parameters:

plex – The DMPlex object encapsulating the mesh topology

firedrake.cython.dmcommon.get_entity_classes()

Builds PyOP2 entity class offsets for all entity levels.

Parameters:

dm – The DM object encapsulating the mesh topology

firedrake.cython.dmcommon.get_facet_nodes()

Build to DoF mapping from facets.

Parameters:
  • mesh – The mesh.

  • cell_nodes – numpy array mapping from cells to function space nodes.

  • label – which set of facets to ask for (interior_facets or exterior_facets).

  • offset – optional offset (extruded only).

Returns:

numpy array mapping from facets to nodes in the closure of the support of that facet.

firedrake.cython.dmcommon.get_facet_ordering()

Builds a list of all facets ordered according to the given numbering.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • facet_numbering – A Section describing the global facet numbering

firedrake.cython.dmcommon.get_facets_by_class()

Builds a list of all facets ordered according to PyOP2 entity classes and computes the respective class offsets.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • ordering – An array giving the global traversal order of facets

  • label – Label string that marks the facets to order

firedrake.cython.dmcommon.get_topological_dimension()

Get the topological dimension of a DMPlex or DMSwarm. Assumes that a DMSwarm represents a mesh of vertices (tdim 0).

Parameters:

dm – A DMPlex or DMSwarm.

Returns:

For a DMPlex dm.getDimension(), for a DMSwarm 0.

firedrake.cython.dmcommon.label_facets()

Add labels to facets in the the plex

Facets on the boundary are marked with “exterior_facets” while all others are marked with “interior_facets”.

firedrake.cython.dmcommon.make_global_numbering()

Build an array of global numbers for local dofs

Parameters:
  • lsec – Section describing local dof layout and numbers.

  • gsec – Section describing global dof layout and numbers.

firedrake.cython.dmcommon.mark_entity_classes()

Mark all points in a given DM according to the PyOP2 entity classes:

core : owned and not in send halo owned : owned and in send halo ghost : in halo

by inspecting the pointSF graph.

Parameters:

dm – The DM object encapsulating the mesh topology

firedrake.cython.dmcommon.mark_entity_classes_using_cell_dm()

Mark all points in a given Particle in Cell (PIC) DMSwarm according to the PyOP2 entity classes (core, owned, ghost) using the markers of the parent DMPlex or DMSwarm cells (i.e. the cells within which the swarm points are located).

firedrake.cython.dmcommon.mark_points_with_function_array()

Marks points in a DMLabel using an indicator function array.

Parameters:
  • plex – DMPlex representing the mesh topology

  • section – Section describing the function space DoF layout and order

  • height – Height of marked points (0 for cells, 1 for facets)

  • array – Array representing the indicator function whose layout is defined by plex, section, and height

  • dmlabel – DMLabel that records marked entities

  • label_value – Value used in dmlabel

firedrake.cython.dmcommon.orientations_facet2cell()

Converts local quadrilateral facet orientations into global quadrilateral cell orientations.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • vertex_numbering – Section describing the universal vertex numbering

  • facet_orientations – Facet orientations (edge directions) relative to the local DMPlex ordering.

  • cell_numbering – Section describing the cell numbering

firedrake.cython.dmcommon.plex_renumbering()

Build a global node renumbering as a permutation of Plex points.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • entity_classes – Array of entity class offsets for each dimension.

  • reordering – A reordering from reordered to original plex points used to provide the traversal order of the cells (i.e. the inverse of the ordering obtained from DMPlexGetOrdering). Optional, if not provided (or None), no reordering is applied and the plex is traversed in original order.

The node permutation is derived from a depth-first traversal of the Plex graph over each entity class in turn. The returned IS is the Plex -> PyOP2 permutation.

firedrake.cython.dmcommon.prune_sf()

Prune an SF of roots referencing the local rank

Parameters:

sf – The PETSc SF to prune.

firedrake.cython.dmcommon.quadrilateral_closure_ordering()

Cellwise orders mesh entities according to the given cell orientations.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • vertex_numbering – Section describing the universal vertex numbering

  • cell_numbering – Section describing the cell numbering

  • cell_orientations – Specifies the starting vertex for each cell, and the order of traversal (CCW or CW).

firedrake.cython.dmcommon.quadrilateral_facet_orientations()

Returns globally synchronised facet orientations (edge directions) incident to locally owned quadrilateral cells.

Parameters:
  • plex – The DMPlex object encapsulating the mesh topology

  • vertex_numbering – Section describing the universal vertex numbering

  • cell_ranks – MPI rank of the owner of each (visible) non-owned cell, or -1 for (locally) owned cell.

firedrake.cython.dmcommon.reordered_coords()

Return coordinates for the dm, reordered according to the global numbering permutation for the coordinate function space.

Shape is a tuple of (mesh.num_vertices(), geometric_dim).

firedrake.cython.dmcommon.set_adjacency_callback()

Set the callback for DMPlexGetAdjacency.

Parameters:

dm – The DMPlex object.

This is used during DMPlexDistributeOverlap to determine where to grow the halos.

firedrake.cython.dmcommon.to_petsc_local_numbering()

Reorder a PETSc Vec corresponding to a Firedrake Function w.r.t. the PETSc natural numbering.

Parameters:
  • vec – the PETSc Vec to reorder; must be a global vector

  • V – the FunctionSpace of the Function which the Vec comes from

Ret out:

a copy of the Vec, ordered with the PETSc natural numbering

firedrake.cython.dmcommon.transform_vec_from_firedrake_to_petsc()

Transform Firedrake vec to PETSc vec.

Parameters:
  • dm (PETSc.DM) – PETSc.DM.

  • finat_element (finat.finiteelementbase.FiniteElementBase) – FInAT element.

  • firedrake_sec (PETSc.Section) – PETSc.Section representing the Firedrake DoF layout.

  • firedrake_vec (PETSc.Vec) – PETSc.Vec representing the Firedrake DoF vector.

  • petsc_sec (PETSc.Section) – PETSc.Section representing the PETSc DoF layout.

  • petsc_vec (PETSc.Vec) – PETSc.Vec representing the PETSc DoF vector.

  • use_transitive_closure_dg (bool) – Whether to use PETSc DG element whose DoFs are ordered according to the transitive closure of the cell.

firedrake.cython.dmcommon.validate_mesh()

Perform some validation of the input mesh.

Parameters:

dm – The DM object encapsulating the mesh topology.

firedrake.cython.extrusion_numbering module

Computation dof numberings for extruded meshes

On meshes with a constant number of cell layers (i.e. each column contains the same number of cells), it is possible to compute all the correct numberings by just lying to DMPlex about how many degrees of freedom there are on the base topological entities.

This ceases to be true as soon as we permit variable numbers of cells in each column, since now, although the number of degrees of freedom on a cell does not change from column to column, the number that are stacked up on each topological entity does change.

This module implements the necessary chicanery to deal with it.

Computation of topological layer extents

First, a picture.

Consider a one-dimensional mesh:

x---0---x---1---x---2---x

Extruded to form the following two-dimensional mesh:

                      x--------x
                      |        |
                      |        |
2                     |        |
                      |        |
    x--------x--------x--------x
    |        |        |
    |        |        |
1   |        |        |
    |        |        |
    x--------x--------x
    |        |
    |        |
0   |        |
    |        |
    x--------x

This is constructed by providing the number of cells in each column as well as the starting cell layer:

[[0, 2],
 [1, 1],
 [2, 1]]

We need to promote this cell layering to layering for all topological entities. Our solution to “interior” facets that only have one side is to require that they are geometrically zero sized, and then guarantee that we never iterate over them. We therefore need to keep track of two bits of information, the layer extent for allocation purposes and also the layer extent for iteration purposes.

We compute both by iterating over the cells and transferring cell layers to points in the closure of each cell. Allocation bounds use min-max on the cell bounds, iteration bounds use max-min.

To simplify some things, we require that the resulting mesh is not topologically disconnected anywhere. Offset cells must, at least, share a vertex with some other cell.

Computation of function space allocation size

With the layer extents computed, we need to compute the dof allocation. For this, we need the number of degrees of freedom on the base topological entity, and above it in each cell:

x-------x
|   o   |
o   o   o
o   o   o
|   o   |
o---o---o

This element has one degree of freedom on each base vertex and cell, two degrees of freedom “above” each vertex, and four above each cell. To compute the number of degrees of freedom on the column of topological entities we sum the number on the entity, multiplied by the number of layers with the number above, multiplied by the number of layers minus one (due to the fencepost error difference). This number of layers naturally changes from entity to entity, and so we can’t compute this up front, but must do it point by point, constructing the PETSc Section as we go.

Computation of function space maps

Now we need the maps from topological entities (cells and facets) to the function space nodes they can see. The allocation offsets that the numbering section gives us are wrong, because when we have a step in the column height, the offset will be wrong if we’re looking from the higher cell. Consider a vertex discretisation on the previous mesh, with a numbering:

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

The cell node map we get by just looking at allocation offsets is:

[[0, 1, 3, 4],
 [3, 4, 6, 7],
 [6, 7, 9, 10]]

note how the second and third cells have the wrong value for their “left” vertices. Instead, we need to shift the numbering we obtain from the allocation offset by the number of entities we’re skipping over, to result in:

[[0, 1, 3, 4],
 [4, 5, 6, 7],
 [7, 8, 9, 10]]

Now, when we iterate over cells, we ensure that we access the correct dofs. The same trick needs to be applied to facet maps too.

Computation of boundary nodes

For the top and bottom boundary nodes, we walk over the cells at, respectively, the top and bottom of the column and pull out those nodes whose entity height matches the appropriate cell height. As an example:

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

The bottom boundary nodes are:

[0, 3, 4, 6, 7, 9]

whereas the top are:

[2, 5, 7, 8, 10]

For these strange “interior” facets, we first walk over the cells, picking up the dofs in the closure of the base (ceiling) of the cell, then we walk over facets, picking up all the dofs in the closure of facets that are exposed (there may be more than one of these in the cell column). We don’t have to worry about any lower-dimensional entities, because if a co-dim 2 or greater entity is exposed in a column, then the co-dim 1 entity in its star is also exposed.

For the side boundary nodes, we can make a simplification: we know that the facet heights are always the same as the cell column heights (because there is only one cell in the support of the facet). Hence we just walk over the boundary facets of the base mesh, extract out the nodes on that facet on the bottom cell and walk up the column. This is guaranteed to pick up all the nodes in the closure of the facet column.

firedrake.cython.extrusion_numbering.entity_layers()

Compute the layers for a given entity type.

Parameters:
  • mesh – the extruded mesh to compute layers for.

  • height – the height of the entity to consider (in the DMPlex sense). e.g. 0 -> cells, 1 -> facets, etc…

  • label – optional label to select some subset of the points of the given height (may be None meaning select all points).

Returns:

a numpy array of shape (num_entities, 2) providing the layer extents for iteration on the requested entities.

firedrake.cython.extrusion_numbering.facet_closure_nodes()

Extract nodes in the closure of facets with a given marker.

This works fine for interior as well as exterior facets.

Note

Don’t call this function directly, but rather call facet_closure_nodes(), which will dispatch here if appropriate.

Parameters:
  • V – the function space

  • sub_domain – a mesh marker selecting the part of the boundary

Returns:

a numpy array of unique nodes on the boundary of the requested subdomain.

firedrake.cython.extrusion_numbering.layer_extents()

Compute the extents (start and stop layers) for an extruded mesh.

Parameters:
  • dm – The DMPlex.

  • cell_numbering – The cell numbering (plex points to Firedrake points).

  • cell_extents – The cell layers.

Returns:

a numpy array of shape (npoints, 4) where npoints is the number of mesh points in the base mesh. npoints[p, 0:2] gives the start and stop layers for allocation for mesh point p (in plex ordering), while npoints[p, 2:4] gives the start and stop layers for iteration over mesh point p (in plex ordering).

Warning

The indexing of this array uses DMPlex point ordering, not Firedrake ordering. So you always need to iterate over plex points and translate to Firedrake numbers if necessary.

firedrake.cython.extrusion_numbering.node_classes()

Compute the node classes for a given extruded mesh.

Parameters:
  • mesh – the extruded mesh.

  • nodes_per_entity – Number of nodes on, and on top of, each type of topological entity on the base mesh for a single cell layer. Multiplying up by the number of layers happens in this function.

Returns:

A numpy array of shape (3, ) giving the set entity sizes for the given nodes per entity.

firedrake.cython.extrusion_numbering.top_bottom_boundary_nodes()

Extract top or bottom boundary nodes from an extruded function space.

Parameters:
  • mesh – The extruded mesh.

  • cell_node_list – The map from cells to nodes.

  • masks – masks for dofs in the closure of the facets of the cell. First the vertical facets, then the horizontal facets (bottom then top).

  • offsets – Offsets to apply walking up the column.

  • kind – Whether we should select the bottom, or the top, nodes.

Returns:

a numpy array of unique indices of nodes on the bottom or top of the mesh.

firedrake.cython.hdf5interface module

firedrake.cython.hdf5interface.get_h5py_file()

Attempt to convert PETSc viewer file handle to h5py File.

Parameters:

vwr – The PETSc Viewer (must have type HDF5).

Warning

For this to work, h5py and PETSc must both have been compiled against the same HDF5 library (otherwise the file handles are not interchangeable). This is the likeliest reason for failure when attempting the conversion.

firedrake.cython.mgimpl module

firedrake.cython.mgimpl.coarse_to_fine_cells()

Return a map from (renumbered) cells in a coarse mesh to those in a refined fine mesh.

Parameters:
  • mc – the coarse mesh to create the map from.

  • mf – the fine mesh to map to.

  • clgmaps – coarse lgmaps (non-overlapped and overlapped)

  • flgmaps – fine lgmaps (non-overlapped and overlapped)

Returns:

Two arrays, one mapping coarse to fine cells, the second fine to coarse cells.

firedrake.cython.mgimpl.coarse_to_fine_nodes()
firedrake.cython.mgimpl.create_lgmap()

Create a local to global map for all points in the given DM.

Parameters:

dm – The DM to create the map for.

Returns a petsc4py LGMap.

firedrake.cython.mgimpl.filter_labels()

Remove labels from points that are not in keep. :arg dm: DM object with labels. :arg keep: subsection of the DMs chart on which to retain label values. :arg label_names: names of labels (strings) to clear. When refining, every point “underneath” the refined entity receives its label. But we typically have labels applied only to entities of a given stratum height (and rely on that elsewhere), so clear the labels from everything else.

firedrake.cython.mgimpl.fine_to_coarse_nodes()
firedrake.cython.mgimpl.get_entity_renumbering()

Given a section numbering a type of topological entity, return the renumberings from original plex numbers to new firedrake numbers (and vice versa)

Parameters:
  • plex – The DMPlex object

  • section – The Section defining the renumbering

  • entity_type – The type of entity (either "cell" or "vertex")

firedrake.cython.patchimpl module

firedrake.cython.patchimpl.set_patch_jacobian()
firedrake.cython.patchimpl.set_patch_residual()

firedrake.cython.spatialindex module

class firedrake.cython.spatialindex.SpatialIndex

Bases: object

Python class for holding a native spatial index object.

ctypes

Returns a ctypes pointer to the native spatial index.

firedrake.cython.spatialindex.bounding_boxes()

Given a spatial index and a point, return the bounding boxes the point is in.

Parameters:
  • sidx – the SpatialIndex

  • x – the point

Returns:

a numpy array of candidate bounding boxes.

firedrake.cython.spatialindex.from_regions()

Builds a spatial index from a set of maximum bounding regions (MBRs).

regions_lo and regions_hi must have the same size. regions_lo[i] and regions_hi[i] contain the coordinates of the diagonally opposite lower and higher corners of the i-th MBR, respectively.

firedrake.cython.supermeshimpl module

firedrake.cython.supermeshimpl.assemble_mixed_mass_matrix()
firedrake.cython.supermeshimpl.intersection_finder()

Module contents