firedrake.mg package

Submodules

firedrake.mg.embedded module

class firedrake.mg.embedded.TransferManager(*, native_transfers=None, use_averaging=True)[source]

Bases: object

An object for managing transfers between levels in a multigrid hierarchy (possibly via embedding in DG spaces).

Parameters:
  • native_transfers – dict mapping UFL element to “natively supported” transfer operators. This should be a three-tuple of (prolong, restrict, inject).

  • use_averaging – Use averaging to approximate the projection out of the embedded DG space? If False, a global L2 projection will be performed.

class Cache(element)[source]

Bases: object

A caching object for work vectors and matrices.

Parameters:

element – The element to use for the caching.

DG_inv_mass(DG)[source]

Inverse DG mass matrix :arg DG: the DG space :returns: A PETSc Mat.

DG_work(V)[source]

A DG work Function matching V :arg V: a function space. :returns: A Function in the embedding DG space.

V_DG_mass(V, DG)[source]

Mass matrix from between V and DG spaces. :arg V: a function space :arg DG: the DG space :returns: A PETSc Mat mapping from V -> DG

V_approx_inv_mass(V, DG)[source]

Approximate inverse mass. Computes (cellwise) (V, V)^{-1} (V, DG). :arg V: a function space :arg DG: the DG space :returns: A PETSc Mat mapping from V -> DG.

V_dof_weights(V)[source]

Dof weights for averaging projection.

Parameters:

V – function space to compute weights for.

Returns:

A PETSc Vec.

V_inv_mass_ksp(V)[source]

A KSP inverting a mass matrix :arg V: a function space. :returns: A PETSc KSP for inverting (V, V).

cache(element)[source]
cache_dat_versions(element, transfer_op, source, target)[source]

Record the returned dat_versions of the source and target.

inject(uf, uc)[source]

Inject a function (primal restriction)

Parameters:
  • uf – The source (fine grid) function.

  • uc – The target (coarse grid) function.

is_native(element)[source]
op(source, target, transfer_op)[source]

Primal transfer (either prolongation or injection).

Parameters:
  • source – The source Function.

  • target – The target Function.

  • transfer_op – The transfer operation for the DG space.

prolong(uc, uf)[source]

Prolong a function.

Parameters:
  • uc – The source (coarse grid) function.

  • uf – The target (fine grid) function.

requires_transfer(element, transfer_op, source, target)[source]

Determine whether either the source or target have been modified since the last time a grid transfer was executed with them.

restrict(source, target)[source]

Restrict a dual function.

Parameters:
work_vec(V)[source]

A work Vec for V :arg V: a function space. :returns: A PETSc Vec for V.

firedrake.mg.interface module

firedrake.mg.interface.inject(fine, coarse)[source]
firedrake.mg.interface.prolong(coarse, fine)[source]
firedrake.mg.interface.restrict(fine_dual, coarse_dual)[source]

firedrake.mg.kernels module

class firedrake.mg.kernels.MacroKernelBuilder(scalar_type, num_entities)[source]

Bases: KernelBuilderBase

Kernel builder for integration on a macro-cell.

Parameters:

num_entities – the number of micro-entities to integrate over.

oriented = False
set_coefficients(coefficients)[source]
set_coordinates(domain)[source]

Prepare the coordinate field.

Parameters:

domainufl.AbstractDomain

firedrake.mg.kernels.compile_element(expression, dual_space=None, parameters=None, name='evaluate')[source]

Generate code for point evaluations.

Parameters:
  • expression – A UFL expression (may contain up to one coefficient, or one argument)

  • dual_space – if the expression has an argument, should we also distribute residual data?

Returns:

The generated code (loopy.TranslationUnit)

firedrake.mg.kernels.dg_injection_kernel(Vf, Vc, ncell)[source]
firedrake.mg.kernels.inject_kernel(Vf, Vc)[source]
firedrake.mg.kernels.prolong_kernel(expression)[source]
firedrake.mg.kernels.restrict_kernel(Vf, Vc)[source]
firedrake.mg.kernels.to_reference_coordinates(ufl_coordinate_element, parameters=None)[source]

firedrake.mg.mesh module

firedrake.mg.mesh.ExtrudedMeshHierarchy(base_hierarchy, height, base_layer=-1, refinement_ratio=2, layers=None, kernel=None, extrusion_type='uniform', gdim=None, mesh_builder=<cyfunction ExtrudedMesh>)[source]

Build a hierarchy of extruded meshes by extruding a hierarchy of meshes.

Parameters:
  • base_hierarchy – the unextruded base mesh hierarchy to extrude.

  • height – the height of the domain to extrude to. This is in contrast to the extrusion routines, which take in layer_height, the height of an individual layer. This is because when refining in the extruded dimension, the height of an individual layer will vary.

  • base_layer – the number of layers to use the extrusion of the coarsest grid.

  • refinement_ratio – the ratio by which base_layer should be increased on every refinement. refinement_ratio = 2 means standard uniform refinement. refinement_ratio = 1 means to not refine in the extruded dimension, i.e. the multigrid hierarchy will use semicoarsening.

  • layers – as an alternative to specifying base_layer and refinement_ratio, one may specify directly the number of layers to be used by each level in the extruded hierarchy. This option cannot be combined with base_layer and refinement_ratio. Note that the ratio of successive entries in this iterable must be an integer for the multigrid transfer operators to work.

  • mesh_builder – function used to turn a Mesh into an extruded mesh. Used by pyadjoint.

See ExtrudedMesh() for the meaning of the remaining parameters.

class firedrake.mg.mesh.HierarchyBase(meshes, coarse_to_fine_cells, fine_to_coarse_cells, refinements_per_level=1, nested=False)[source]

Bases: object

Create an encapsulation of an hierarchy of meshes.

Parameters:
  • meshes – list of meshes (coarse to fine)

  • coarse_to_fine_cells – list of numpy arrays for each level pair, mapping each coarse cell into fine cells it intersects.

  • fine_to_coarse_cells – list of numpy arrays for each level pair, mapping each fine cell into coarse cells it intersects.

  • refinements_per_level – number of mesh refinements each multigrid level should “see”.

  • nested – Is this mesh hierarchy nested?

Note

Most of the time, you do not need to create this object yourself, instead using MeshHierarchy(), ExtrudedMeshHierarchy(), or NonNestedHierarchy().

comm[source]
firedrake.mg.mesh.MeshHierarchy(mesh, refinement_levels, refinements_per_level=1, netgen_flags=False, reorder=None, distribution_parameters=None, callbacks=None, mesh_builder=<cyfunction Mesh>)[source]

Build a hierarchy of meshes by uniformly refining a coarse mesh.

Parameters:
  • mesh (MeshGeometry) – the coarse mesh to refine

  • refinement_levels (int) – the number of levels of refinement

  • refinements_per_level (int) – the number of refinements for each level in the hierarchy.

  • netgen_flags (bool, dict) – either a bool or a dictionary containing options for Netgen. If not False the hierachy is constructed using ngsPETSc, if None hierarchy constructed in a standard manner.

  • distribution_parameters (dict) – options controlling mesh distribution, see Mesh() for details. If None, use the same distribution parameters as were used to distribute the coarse mesh, otherwise, these options override the default.

  • reorder (bool) – optional flag indicating whether to reorder the refined meshes.

  • callbacks (tuple) – A 2-tuple of callbacks to call before and after refinement of the DM. The before callback receives the DM to be refined (and the current level), the after callback receives the refined DM (and the current level).

  • mesh_builder – Function to turn a DM into a Mesh. Used by pyadjoint.

Returns:

firedrake.mg.mesh.NonNestedHierarchy(*meshes)[source]
firedrake.mg.mesh.SemiCoarsenedExtrudedHierarchy(base_mesh, height, nref=1, base_layer=-1, refinement_ratio=2, layers=None, kernel=None, extrusion_type='uniform', gdim=None, mesh_builder=<cyfunction ExtrudedMesh>)[source]

Build a hierarchy of extruded meshes with refinement only in the extruded dimension.

Parameters:
  • base_mesh – the unextruded base mesh to extrude.

  • nref – Number of refinements.

  • height – the height of the domain to extrude to. This is in contrast to the extrusion routines, which take in layer_height, the height of an individual layer. This is because when refining in the extruded dimension, the height of an individual layer will vary.

  • base_layer – the number of layers to use the extrusion of the coarsest grid.

  • refinement_ratio – the ratio by which base_layer should be increased on every refinement. refinement_ratio = 2 means standard uniform refinement. refinement_ratio = 1 means to not refine in the extruded dimension, i.e. the multigrid hierarchy will use semicoarsening.

  • layers – as an alternative to specifying base_layer and refinement_ratio, one may specify directly the number of layers to be used by each level in the extruded hierarchy. This option cannot be combined with base_layer and refinement_ratio. Note that the ratio of successive entries in this iterable must be an integer for the multigrid transfer operators to work.

  • mesh_builder – function used to turn a Mesh into an extruded mesh. Used by pyadjoint.

See ExtrudedMesh() for the meaning of the remaining parameters.

See also ExtrudedMeshHierarchy() if you want to extruded a hierarchy of unstructured meshes.

firedrake.mg.opencascade_mh module

firedrake.mg.opencascade_mh.OpenCascadeMeshHierarchy(stepfile, element_size, levels, comm=<mpi4py.MPI.Intracomm object>, distribution_parameters=None, callbacks=None, order=1, mh_constructor=<function MeshHierarchy>, cache=True, verbose=True, gmsh='gmsh', project_refinements_to_cad=True, reorder=None)[source]

firedrake.mg.ufl_utils module

firedrake.mg.ufl_utils.coarsen(expr, self, coefficient_mapping=None)[source]
firedrake.mg.ufl_utils.coarsen(mesh: Mesh, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(expr: Expr, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(expr: BaseForm, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(form: Form, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(form: FormSum, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(bc: DirichletBC, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(V: WithGeometryBase, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(expr: Function, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(expr: Cofunction, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(problem: NonlinearVariationalProblem, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(basis: VectorSpaceBasis, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(mspbasis: MixedVectorSpaceBasis, self, coefficient_mapping=None)
firedrake.mg.ufl_utils.coarsen(context: _SNESContext, self, coefficient_mapping=None)

firedrake.mg.utils module

firedrake.mg.utils.coarse_cell_to_fine_node_map(Vc, Vf)[source]
firedrake.mg.utils.coarse_node_to_fine_node_map(Vc, Vf)[source]
firedrake.mg.utils.fine_node_to_coarse_node_map(Vf, Vc)[source]
firedrake.mg.utils.get_level(obj)[source]

Try and obtain hierarchy and level info from an object.

If no level info is available, return None, None.

firedrake.mg.utils.has_level(obj)[source]

Does the provided object have level info?

firedrake.mg.utils.physical_node_locations(V)[source]
firedrake.mg.utils.set_level(obj, hierarchy, level)[source]

Attach hierarchy and level info to an object.

Module contents