petsc4py.PETSc.DMPlex#

class petsc4py.PETSc.DMPlex#

Bases: DM

Encapsulate an unstructured mesh.

DMPlex encapsulates both topology and geometry. It is capable of parallel refinement and coarsening (using Pragmatic or ParMmg) and parallel redistribution for load balancing.

Methods Summary

computeCellGeometryFVM(cell)

Compute the volume for a given cell.

computeGradientClementInterpolant(locX, locC)

Return the L2 projection of the cellwise gradient of a function onto P1.

constructGhostCells([labelName])

Construct ghost cells which connect to every boundary face.

coordinatesLoad(viewer, sfxc)

Load coordinates into this DMPlex object.

coordinatesView(viewer)

Save DMPlex coordinates into a file.

create([comm])

Create a DMPlex object, which encapsulates an unstructured mesh.

createBoxMesh(faces[, lower, upper, ...])

Create a mesh on the tensor product of intervals.

createBoxSurfaceMesh(faces[, lower, upper, ...])

Create a mesh on the surface of a box mesh using tensor cells.

createCGNS(cgid[, interpolate, comm])

Create a DMPlex mesh from a CGNS file.

createCGNSFromFile(filename[, interpolate, comm])

"Create a DMPlex mesh from a CGNS file.

createClosureIndex(sec)

Calculate an index for sec for the closure operation.

createCoarsePointIS()

Create an IS covering the coarse DMPlex chart with the fine points as data.

createCohesiveSubmesh(hasLagrange, value)

Extract the hypersurface defined by one face of the cohesive cells.

createExodus(exoid[, interpolate, comm])

Create a DMPlex mesh from an ExodusII file ID.

createExodusFromFile(filename[, ...])

Create a DMPlex mesh from an ExodusII file.

createFromCellList(dim, cells, coords[, ...])

Create a DMPlex from a list of vertices for each cell on process 0.

createFromFile(filename[, plexname, ...])

Create DMPlex from a file.

createGmsh(viewer[, interpolate, comm])

Create a DMPlex mesh from a Gmsh file viewer.

createPointNumbering()

Create a global numbering for all points.

createSection(numComp, numDof[, bcField, ...])

Create a Section based upon the DOF layout specification provided.

distribute([overlap])

Distribute the mesh and any associated sections.

distributeField(sf, sec, vec[, newsec, newvec])

Distribute field data with a with a given SF.

distributeGetDefault()

Return a flag indicating whether the DM should be distributed by default.

distributeOverlap([overlap])

Add partition overlap to a distributed non-overlapping DMPlex.

distributeSetDefault(flag)

Set flag indicating whether the DMPlex should be distributed by default.

distributionGetName()

Retrieve the name of the specific parallel distribution.

distributionSetName(name)

Set the name of the specific parallel distribution.

generate(boundary[, name, interpolate])

Generate a mesh.

getAdjacency(p)

Return all points adjacent to the given point.

getAdjacencyUseAnchors()

Query whether adjacency in the mesh uses the point-to-point constraints.

getCellNumbering()

Return a global cell numbering for all cells on this process.

getCellType(p)

Return the polytope type of a given cell.

getCellTypeLabel()

Return the DMLabel recording the polytope type of each cell.

getChart()

Return the interval for all mesh points [pStart, pEnd).

getCone(p)

Return the points on the in-edges for this point in the DAG.

getConeOrientation(p)

Return the orientations on the in-edges for this point in the DAG.

getConeSize(p)

Return the number of in-edges for this point in the DAG.

getDepth()

Return the depth of the DAG representing this mesh.

getDepthStratum(svalue)

Return the bounds [start, end) for all points at a certain depth.

getFullJoin(points)

Return an array for the join of the set of points.

getHeightStratum(svalue)

Return the bounds [start, end) for all points at a certain height.

getJoin(points)

Return an array for the join of the set of points.

getMaxSizes()

Return the maximum number of in-edges and out-edges of the DAG.

getMeet(points)

Return an array for the meet of the set of points.

getMinRadius()

Return the minimum distance from any cell centroid to a face.

getOrdering(otype)

Calculate a reordering of the mesh.

getPartitioner()

Return the mesh partitioner.

getPointDepth(point)

Return the depth of a given point.

getPointGlobal(point)

Return location of point data in global Vec.

getPointGlobalField(point, field)

Return location of point field data in global Vec.

getPointHeight(point)

Return the height of a given point.

getPointLocal(point)

Return location of point data in local Vec.

getPointLocalField(point, field)

Return location of point field data in local Vec.

getRefinementLimit()

Retrieve the maximum cell volume for refinement.

getRefinementUniform()

Retrieve the flag for uniform refinement.

getSubpointIS()

Return an IS covering the entire subdm chart.

getSubpointMap()

Return a DMLabel with point dimension as values.

getSupport(p)

Return the points on the out-edges for this point in the DAG.

getSupportSize(p)

Return the number of out-edges for this point in the DAG.

getTransitiveClosure(p[, useCone])

Return the points and orientations on the transitive closure of this point.

getVecClosure(sec, vec, point)

Return an array of the values on the closure of a point.

getVertexNumbering()

Return a global vertex numbering for all vertices on this process.

globalVectorLoad(viewer, sectiondm, sf, vec)

Load on-disk vector data into a global vector.

globalVectorView(viewer, sectiondm, vec)

Save a global vector.

insertCone(p, conePos, conePoint)

DMPlexInsertCone - Insert a point into the in-edges for the point p in the DAG.

insertConeOrientation(p, conePos, ...)

Insert a point orientation for the in-edge for the point p in the DAG.

interpolate()

Convert to a mesh with all intermediate faces, edges, etc.

isDistributed()

Return the flag indicating if the mesh is distributed.

isSimplex()

Return the flag indicating if the first cell is a simplex.

labelCohesiveComplete(label, bdlabel, ...)

Add all other mesh pieces to complete the surface.

labelComplete(label)

Add the transitive closure to the surface.

labelsLoad(viewer, sfxc)

Load labels into this DMPlex object.

labelsView(viewer)

Save DMPlex labels into a file.

localVectorLoad(viewer, sectiondm, sf, vec)

Load on-disk vector data into a local vector.

localVectorView(viewer, sectiondm, vec)

Save a local vector.

markBoundaryFaces(label[, value])

Mark all faces on the boundary.

metricAverage2(metric1, metric2, metricAvg)

Compute and return the unweighted average of two metrics.

metricAverage3(metric1, metric2, metric3, ...)

Compute and return the unweighted average of three metrics.

metricCreate([field])

Create a Riemannian metric field.

metricCreateIsotropic(indicator[, field])

Construct an isotropic metric from an error indicator.

metricCreateUniform(alpha[, field])

Construct a uniform isotropic metric.

metricDeterminantCreate([field])

Create the determinant field for a Riemannian metric.

metricEnforceSPD(metric, ometric, determinant)

Enforce symmetric positive-definiteness of a metric.

metricGetGradationFactor()

Return the metric gradation factor.

metricGetHausdorffNumber()

Return the metric Hausdorff number.

metricGetMaximumAnisotropy()

Return the maximum tolerated metric anisotropy.

metricGetMaximumMagnitude()

Return the maximum tolerated metric magnitude.

metricGetMinimumMagnitude()

Return the minimum tolerated metric magnitude.

metricGetNormalizationOrder()

Return the order p for L-p normalization.

metricGetNumIterations()

Return the number of parallel adaptation iterations.

metricGetTargetComplexity()

Return the target metric complexity.

metricGetVerbosity()

Return the verbosity of the mesh adaptation package.

metricIntersection2(metric1, metric2, metricInt)

Compute and return the intersection of two metrics.

metricIntersection3(metric1, metric2, ...)

Compute the intersection of three metrics.

metricIsIsotropic()

Return the flag indicating whether the metric is isotropic or not.

metricIsUniform()

Return the flag indicating whether the metric is uniform or not.

metricNoInsertion()

Return the flag indicating whether node insertion and deletion are turned off.

metricNoMovement()

Return the flag indicating whether node movement is turned off.

metricNoSurf()

Return the flag indicating whether surface modification is turned off.

metricNoSwapping()

Return the flag indicating whether facet swapping is turned off.

metricNormalize(metric, ometric, determinant)

Apply L-p normalization to a metric.

metricRestrictAnisotropyFirst()

Return true if anisotropy is restricted before normalization.

metricSetFromOptions()

Configure the object from the options database.

metricSetGradationFactor(beta)

Set the metric gradation factor.

metricSetHausdorffNumber(hausd)

Set the metric Hausdorff number.

metricSetIsotropic(isotropic)

Record whether the metric is isotropic or not.

metricSetMaximumAnisotropy(a_max)

Set the maximum tolerated metric anisotropy.

metricSetMaximumMagnitude(h_max)

Set the maximum tolerated metric magnitude.

metricSetMinimumMagnitude(h_min)

Set the minimum tolerated metric magnitude.

metricSetNoInsertion(noInsert)

Set the flag indicating whether node insertion should be turned off.

metricSetNoMovement(noMove)

Set the flag indicating whether node movement should be turned off.

metricSetNoSurf(noSurf)

Set the flag indicating whether surface modification should be turned off.

metricSetNoSwapping(noSwap)

Set the flag indicating whether facet swapping should be turned off.

metricSetNormalizationOrder(p)

Set the order p for L-p normalization.

metricSetNumIterations(numIter)

Set the number of parallel adaptation iterations.

metricSetRestrictAnisotropyFirst(...)

Record whether anisotropy is be restricted before normalization or after.

metricSetTargetComplexity(targetComplexity)

Set the target metric complexity.

metricSetUniform(uniform)

Record whether the metric is uniform or not.

metricSetVerbosity(verbosity)

Set the verbosity of the mesh adaptation package.

orient()

Give a consistent orientation to the input mesh.

permute(perm)

Reorder the mesh according to the input permutation.

rebalanceSharedPoints([entityDepth, ...])

Redistribute shared points in order to achieve better balancing.

reorderGetDefault()

Return flag indicating whether the DMPlex should be reordered by default.

reorderSetDefault(flag)

Set flag indicating whether the DM should be reordered by default.

sectionLoad(viewer, sectiondm, sfxc)

Load section into a DM.

sectionView(viewer, sectiondm)

Save a section associated with a DMPlex.

setAdjacencyUseAnchors([useAnchors])

Define adjacency in the mesh using the point-to-point constraints.

setCellType(p, ctype)

Set the polytope type of a given cell.

setChart(pStart, pEnd)

Set the interval for all mesh points [pStart, pEnd).

setCone(p, cone[, orientation])

Set the points on the in-edges for this point in the DAG.

setConeOrientation(p, orientation)

Set the orientations on the in-edges for this point in the DAG.

setConeSize(p, size)

Set the number of in-edges for this point in the DAG.

setMatClosure(sec, gsec, mat, point, values)

Set an array of the values on the closure of point.

setPartitioner(part)

Set the mesh partitioner.

setRefinementLimit(refinementLimit)

Set the maximum cell volume for refinement.

setRefinementUniform([refinementUniform])

Set the flag for uniform refinement.

setSupport(p, supp)

Set the points on the out-edges for this point in the DAG.

setSupportSize(p, size)

Set the number of out-edges for this point in the DAG.

setTetGenOptions(opts)

Set the options used for the Tetgen mesh generator.

setTriangleOptions(opts)

Set the options used for the Triangle mesh generator.

setVecClosure(sec, vec, point, values[, addv])

Set an array of the values on the closure of point.

stratify()

Calculate the strata of DAG.

symmetrize()

Create support (out-edge) information from cone (in-edge) information.

topologyLoad(viewer)

Load a topology into this DMPlex object.

topologyView(viewer)

Save a DMPlex topology into a file.

uninterpolate()

Convert to a mesh with only cells and vertices.

vecGetClosure(sec, vec, p)

Return an array of values on the closure of p.

Methods Documentation

computeCellGeometryFVM(cell)#

Compute the volume for a given cell.

Not collective.

Parameters:

cell (int) – The cell.

Returns:

  • volume (float) – The cell volume.

  • centroid (ArrayReal) – The cell centroid.

  • normal (ArrayReal) – The cell normal, if appropriate.

Return type:

tuple[float, ArrayReal, ArrayReal]

Source code at petsc4py/PETSc/DMPlex.pyx:2209

computeGradientClementInterpolant(locX, locC)#

Return the L2 projection of the cellwise gradient of a function onto P1.

Collective.

Parameters:
  • locX (Vec) – The coefficient vector of the function.

  • locC (Vec) – The output Vec which holds the Clement interpolant of the gradient.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:3143

constructGhostCells(labelName=None)#

Construct ghost cells which connect to every boundary face.

Collective.

Parameters:

labelName (str | None) – The name of the label specifying the boundary faces. Defaults to "Face Sets".

Returns:

numGhostCells – The number of ghost cells added to the DMPlex.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:2241

coordinatesLoad(viewer, sfxc)#

Load coordinates into this DMPlex object.

Collective.

Parameters:
Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3316

coordinatesView(viewer)#

Save DMPlex coordinates into a file.

Collective.

Parameters:

viewer (Viewer) – The Viewer for saving.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3183

create(comm=None)#

Create a DMPlex object, which encapsulates an unstructured mesh.

Collective.

Parameters:

comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:14

createBoxMesh(faces, lower=(0, 0, 0), upper=(1, 1, 1), simplex=True, periodic=False, interpolate=True, comm=None)#

Create a mesh on the tensor product of intervals.

Collective.

Parameters:
Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:90

createBoxSurfaceMesh(faces, lower=(0, 0, 0), upper=(1, 1, 1), interpolate=True, comm=None)#

Create a mesh on the surface of a box mesh using tensor cells.

Collective.

Parameters:
Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:141

createCGNS(cgid, interpolate=True, comm=None)#

Create a DMPlex mesh from a CGNS file.

Collective.

Parameters:
  • cgid (int) – The CG id associated with a file and obtained using cg_open.

  • interpolate (bool | None) – Create faces and edges in the mesh.

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:218

createCGNSFromFile(filename, interpolate=True, comm=None)#

“Create a DMPlex mesh from a CGNS file.

Collective.

Parameters:
  • filename (str) – The name of the CGNS file.

  • interpolate (bool | None) – Create faces and edges in the mesh.

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:246

createClosureIndex(sec)#

Calculate an index for sec for the closure operation.

Not collective.

Parameters:

sec (Section) – The Section describing the layout in the local vector, or None to use the default section.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2023

createCoarsePointIS()#

Create an IS covering the coarse DMPlex chart with the fine points as data.

Collective.

Returns:

fpointIS – The IS of all the fine points which exist in the original coarse mesh.

Return type:

IS

Source code at petsc4py/PETSc/DMPlex.pyx:1809

createCohesiveSubmesh(hasLagrange, value)#

Extract the hypersurface defined by one face of the cohesive cells.

Collective.

Parameters:
  • hasLagrange (bool) – Flag indicating whether the mesh has Lagrange dofs in the cohesive cells.

  • value (int) – A label value.

Return type:

DMPlex

Source code at petsc4py/PETSc/DMPlex.pyx:369

createExodus(exoid, interpolate=True, comm=None)#

Create a DMPlex mesh from an ExodusII file ID.

Collective.

Parameters:
  • exoid (int) – The ExodusII id associated with a file obtained using ex_open.

  • interpolate (bool | None) – Create faces and edges in the mesh,

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:304

createExodusFromFile(filename, interpolate=True, comm=None)#

Create a DMPlex mesh from an ExodusII file.

Collective.

Parameters:
  • filename (str) – The name of the ExodusII file.

  • interpolate (bool | None) – Create faces and edges in the mesh.

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:275

createFromCellList(dim, cells, coords, interpolate=True, comm=None)#

Create a DMPlex from a list of vertices for each cell on process 0.

Collective.

Parameters:
  • dim (int) – The topological dimension of the mesh.

  • cells (Sequence[int]) – An array of number of cells times number of vertices on each cell.

  • coords (Sequence[float]) – An array of number of vertices times spatial dimension for coordinates.

  • interpolate (bool | None) – Flag to interpolate the mesh.

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:35

createFromFile(filename, plexname='unnamed', interpolate=True, comm=None)#

Create DMPlex from a file.

Collective.

Parameters:
  • filename (str) – A file name.

  • plexname (str | None) – The name of the resulting DMPlex, also used for intra-datafile lookup by some formats.

  • interpolate (bool | None) – Flag to create intermediate mesh pieces (edges, faces).

  • comm (Comm | None) – MPI communicator, defaults to Sys.getDefaultComm.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:184

createGmsh(viewer, interpolate=True, comm=None)#

Create a DMPlex mesh from a Gmsh file viewer.

Collective.

Parameters:
Return type:

Self

Notes

-dm_plex_gmsh_hybrid forces triangular prisms to use tensor order.

-dm_plex_gmsh_periodic allows for reading Gmsh periodic section.

-dm_plex_gmsh_highorder allows for generating high-order coordinates.

-dm_plex_gmsh_project projects high-order coordinates to a different space, use the prefix -dm_plex_gmsh_project_ to define the space.

-dm_plex_gmsh_use_regions generates labels with region names.

-dm_plex_gmsh_mark_vertices adds vertices to generated labels.

-dm_plex_gmsh_multiple_tags allows multiple tags for default labels.

-dm_plex_gmsh_spacedim <d> embedding space dimension.

Source code at petsc4py/PETSc/DMPlex.pyx:331

createPointNumbering()#

Create a global numbering for all points.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:907

Return type:

IS

createSection(numComp, numDof, bcField=None, bcComps=None, bcPoints=None, perm=None)#

Create a Section based upon the DOF layout specification provided.

Not collective.

Parameters:
  • numComp (Sequence[int]) – An array of size numFields holding the number of components per field.

  • numDof (Sequence[int]) – An array of size numFields*(dim+1) holding the number of DOFs per field on a mesh piece of dimension dim.

  • bcField (Sequence[int] | None) – An array of size numBC giving the field number for each boundary condition, where numBC is the number of boundary conditions.

  • bcComps (Sequence[IS] | None) – An array of size numBC giving an IS holding the field components to which each boundary condition applies.

  • bcPoints (Sequence[IS] | None) – An array of size numBC giving an IS holding the DMPlex points to which each boundary condition applies.

  • perm (IS | None) – Permutation of the chart.

Return type:

Section

Source code at petsc4py/PETSc/DMPlex.pyx:1829

distribute(overlap=0)#

Distribute the mesh and any associated sections.

Collective.

Parameters:

overlap (int | None) – The overlap of partitions.

Returns:

sf – The SF used for point distribution, or None if not distributed.

Return type:

SF or None

Source code at petsc4py/PETSc/DMPlex.pyx:1553

distributeField(sf, sec, vec, newsec=None, newvec=None)#

Distribute field data with a with a given SF.

Collective.

Parameters:
  • sf (SF) – The SF describing the communication pattern.

  • sec (Section) – The Section for existing data layout.

  • vec (Vec) – The existing data in a local vector.

  • newsec (Section | None) – The SF describing the new data layout.

  • newvec (Vec | None) – The new data in a local vector.

Returns:

  • newSection (Section) – The SF describing the new data layout.

  • newVec (Vec) – The new data in a local vector.

Return type:

tuple[Section, Vec]

Source code at petsc4py/PETSc/DMPlex.pyx:1750

distributeGetDefault()#

Return a flag indicating whether the DM should be distributed by default.

Not collective.

Returns:

dist – Flag indicating whether the DMPlex should be distributed by default.

Return type:

bool

Source code at petsc4py/PETSc/DMPlex.pyx:1639

distributeOverlap(overlap=0)#

Add partition overlap to a distributed non-overlapping DMPlex.

Collective.

Parameters:

overlap (int | None) – The overlap of partitions (the same on all ranks).

Returns:

sf – The SF used for point distribution.

Return type:

SF

Source code at petsc4py/PETSc/DMPlex.pyx:1581

distributeSetDefault(flag)#

Set flag indicating whether the DMPlex should be distributed by default.

Logically collective.

Parameters:

flag (bool) – Flag indicating whether the DMPlex should be distributed by default.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1659

distributionGetName()#

Retrieve the name of the specific parallel distribution.

Not collective.

Returns:

name – The name of the specific parallel distribution.

Return type:

str

Source code at petsc4py/PETSc/DMPlex.pyx:1700

distributionSetName(name)#

Set the name of the specific parallel distribution.

Logically collective.

Parameters:

name (str) – The name of the specific parallel distribution.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1679

generate(boundary, name=None, interpolate=True)#

Generate a mesh.

Not collective.

Parameters:
  • boundary (DMPlex) – The DMPlex boundary object.

  • name (str | None) – The mesh generation package name.

  • interpolate (bool | None) – Flag to create intermediate mesh elements.

Return type:

Self

Source code at petsc4py/PETSc/DMPlex.pyx:1276

getAdjacency(p)#

Return all points adjacent to the given point.

Not collective.

Parameters:

p (int) – The point.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:1459

getAdjacencyUseAnchors()#

Query whether adjacency in the mesh uses the point-to-point constraints.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1444

Return type:

bool

getCellNumbering()#

Return a global cell numbering for all cells on this process.

Collective the first time it is called.

Source code at petsc4py/PETSc/DMPlex.pyx:877

Return type:

IS

getCellType(p)#

Return the polytope type of a given cell.

Not collective.

Parameters:

p (int) – The cell.

Return type:

PolytopeType

Source code at petsc4py/PETSc/DMPlex.pyx:677

getCellTypeLabel()#

Return the DMLabel recording the polytope type of each cell.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:698

Return type:

DMLabel

getChart()#

Return the interval for all mesh points [pStart, pEnd).

Not collective.

Returns:

  • pStart (int) – The first mesh point.

  • pEnd (int) – The upper bound for mesh points.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:392

getCone(p)#

Return the points on the in-edges for this point in the DAG.

Not collective.

Parameters:

p (int) – The point, which must lie in the chart set with DMPlex.setChart.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:483

getConeOrientation(p)#

Return the orientations on the in-edges for this point in the DAG.

Not collective.

Parameters:

p (int) – The point, which must lie in the chart set with DMPlex.setChart.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:598

getConeSize(p)#

Return the number of in-edges for this point in the DAG.

Not collective.

Parameters:

p (int) – The point, which must lie in the chart set with DMPlex.setChart.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:434

getDepth()#

Return the depth of the DAG representing this mesh.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:921

Return type:

int

getDepthStratum(svalue)#

Return the bounds [start, end) for all points at a certain depth.

Not collective.

Parameters:

svalue (int) – The requested depth.

Returns:

  • pStart (int) – The first stratum point.

  • pEnd (int) – The upper bound for stratum points.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:936

getFullJoin(points)#

Return an array for the join of the set of points.

Not collective.

Parameters:

points (Sequence[int]) – The input points.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:1084

getHeightStratum(svalue)#

Return the bounds [start, end) for all points at a certain height.

Not collective.

Parameters:

svalue (int) – The requested height.

Returns:

  • pStart (int) – The first stratum point.

  • pEnd (int) – The upper bound for stratum points.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:963

getJoin(points)#

Return an array for the join of the set of points.

Not collective.

Parameters:

points (Sequence[int]) – The input points.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:1058

getMaxSizes()#

Return the maximum number of in-edges and out-edges of the DAG.

Not collective.

Returns:

  • maxConeSize (int) – The maximum number of in-edges.

  • maxSupportSize (int) – The maximum number of out-edges.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:818

getMeet(points)#

Return an array for the meet of the set of points.

Not collective.

Parameters:

points (Sequence[int]) – The input points.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:1032

getMinRadius()#

Return the minimum distance from any cell centroid to a face.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1795

Return type:

float

getOrdering(otype)#

Calculate a reordering of the mesh.

Collective.

Parameters:

otype (OrderingType) – Type of reordering, see Mat.OrderingType.

Returns:

perm – The point permutation.

Return type:

IS

Source code at petsc4py/PETSc/DMPlex.pyx:2122

getPartitioner()#

Return the mesh partitioner.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1502

Return type:

Partitioner

getPointDepth(point)#

Return the depth of a given point.

Not collective.

Parameters:

point (int) – The point.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:990

getPointGlobal(point)#

Return location of point data in global Vec.

Not collective.

Parameters:

point (int) – The topological point.

Returns:

  • start (int) – Start of point data; returns -(globalStart+1) if point is not owned.

  • end (int) – End of point data; returns -(globalEnd+1) if point is not owned.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:1964

getPointGlobalField(point, field)#

Return location of point field data in global Vec.

Not collective.

Parameters:
  • point (int) – The topological point.

  • field (int) – The field number.

Returns:

  • start (int) – Start of point data; returns -(globalStart+1) if point is not owned.

  • end (int) – End of point data; returns -(globalEnd+1) if point is not owned.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:1992

getPointHeight(point)#

Return the height of a given point.

Not collective.

Parameters:

point (int) – The point.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:1011

getPointLocal(point)#

Return location of point data in local Vec.

Not collective.

Parameters:

point (int) – The topological point.

Returns:

  • start (int) – Start of point data.

  • end (int) – End of point data.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:1905

getPointLocalField(point, field)#

Return location of point field data in local Vec.

Not collective.

Parameters:
  • point (int) – The topological point.

  • field (int) – The field number.

Returns:

  • start (int) – Start of point data.

  • end (int) – End of point data.

Return type:

tuple[int, int]

Source code at petsc4py/PETSc/DMPlex.pyx:1933

getRefinementLimit()#

Retrieve the maximum cell volume for refinement.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2106

Return type:

float

getRefinementUniform()#

Retrieve the flag for uniform refinement.

Not collective.

Returns:

refinementUniform – The flag for uniform refinement.

Return type:

bool

Source code at petsc4py/PETSc/DMPlex.pyx:2065

getSubpointIS()#

Return an IS covering the entire subdm chart.

Not collective.

Returns:

iset – The IS containing subdm’s parent’s points.

Return type:

IS

Source code at petsc4py/PETSc/DMPlex.pyx:2270

getSubpointMap()#

Return a DMLabel with point dimension as values.

Not collective.

Returns:

label – The DMLabel whose values are subdm’s point dimensions.

Return type:

DMLabel

Source code at petsc4py/PETSc/DMPlex.pyx:2290

getSupport(p)#

Return the points on the out-edges for this point in the DAG.

Not collective.

Parameters:

p (int) – The point, which must lie in the chart set with DMPlex.setChart.

Return type:

ArrayInt

Source code at petsc4py/PETSc/DMPlex.pyx:763

getSupportSize(p)#

Return the number of out-edges for this point in the DAG.

Not collective.

Parameters:

p (int) – The point, which must lie in the chart set with DMPlex.setChart.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:714

getTransitiveClosure(p, useCone=True)#

Return the points and orientations on the transitive closure of this point.

Not collective.

Parameters:
  • p (int) – The mesh point.

  • useCone (bool | None) – True for the closure, otherwise return the star.

Returns:

  • points (ArrayInt) – The points.

  • orientations (ArrayInt) – The orientations.

Return type:

tuple[ArrayInt, ArrayInt]

Source code at petsc4py/PETSc/DMPlex.pyx:1110

getVecClosure(sec, vec, point)#

Return an array of the values on the closure of a point.

Not collective.

Parameters:
  • sec (Section) – The Section describing the layout in vec or None to use the default section.

  • vec (Vec) – The local vector.

  • point (int) – The point in the DMPlex.

Return type:

ArrayScalar

Source code at petsc4py/PETSc/DMPlex.pyx:1177

getVertexNumbering()#

Return a global vertex numbering for all vertices on this process.

Collective the first time it is called.

Source code at petsc4py/PETSc/DMPlex.pyx:892

Return type:

IS

globalVectorLoad(viewer, sectiondm, sf, vec)#

Load on-disk vector data into a global vector.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer that represents the on-disk vector data.

  • sectiondm (DM) – The DM that contains the global section on which vec is defined.

  • sf (SF) – The SF that migrates the on-disk vector data into vec.

  • vec (Vec) – The global vector to set values of.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3393

globalVectorView(viewer, sectiondm, vec)#

Save a global vector.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer to save data with.

  • sectiondm (DM) – The DM containing the global section on which vec is defined; may be the same as this DMPlex object.

  • vec (Vec) – The global vector to be saved.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3240

insertCone(p, conePos, conePoint)#

DMPlexInsertCone - Insert a point into the in-edges for the point p in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart.

  • conePos (int) – The local index in the cone where the point should be put.

  • conePoint (int) – The mesh point to insert.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:548

insertConeOrientation(p, conePos, coneOrientation)#

Insert a point orientation for the in-edge for the point p in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart

  • conePos (int) – The local index in the cone where the point should be put.

  • coneOrientation (int) – The point orientation to insert.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:573

interpolate()#

Convert to a mesh with all intermediate faces, edges, etc.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1720

Return type:

None

isDistributed()#

Return the flag indicating if the mesh is distributed.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1610

Return type:

bool

isSimplex()#

Return the flag indicating if the first cell is a simplex.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1624

Return type:

bool

labelCohesiveComplete(label, bdlabel, bdvalue, flip, split, subdm)#

Add all other mesh pieces to complete the surface.

Not collective.

Parameters:
  • label (DMLabel) – A DMLabel marking the surface.

  • bdlabel (DMLabel) – A DMLabel marking the vertices on the boundary which will not be duplicated.

  • bdvalue (int) – Value of DMLabel marking the vertices on the boundary.

  • flip (bool) – Flag to flip the submesh normal and replace points on the other side.

  • split (bool) – Flag to split faces incident on the surface boundary, rather than clamping those faces to the boundary

  • subdm (DMPlex) – The DMPlex associated with the label.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1388

labelComplete(label)#

Add the transitive closure to the surface.

Not collective.

Parameters:

label (DMLabel) – A DMLabel marking the surface points.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1371

labelsLoad(viewer, sfxc)#

Load labels into this DMPlex object.

Collective.

Parameters:
Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3336

labelsView(viewer)#

Save DMPlex labels into a file.

Collective.

Parameters:

viewer (Viewer) – The Viewer for saving.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3201

localVectorLoad(viewer, sectiondm, sf, vec)#

Load on-disk vector data into a local vector.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer that represents the on-disk vector data.

  • sectiondm (DM) – The DM that contains the local section on which vec is defined.

  • sf (SF) – The SF that migrates the on-disk vector data into vec.

  • vec (Vec) – The local vector to set values of.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3418

localVectorView(viewer, sectiondm, vec)#

Save a local vector.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer to save data with.

  • sectiondm (DM) – The DM that contains the local section on which vec is defined; may be the same as this DMPlex object.

  • vec (Vec) – The local vector to be saved.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3264

markBoundaryFaces(label, value=None)#

Mark all faces on the boundary.

Not collective.

Parameters:
  • value (int | None) – The marker value, or DETERMINE or None to use some value in the closure (or 1 if none are found).

  • label (str)

Return type:

DMLabel

Source code at petsc4py/PETSc/DMPlex.pyx:1344

metricAverage2(metric1, metric2, metricAvg)#

Compute and return the unweighted average of two metrics.

Collective.

Parameters:
  • metric1 (Vec) – The first metric to be averaged.

  • metric2 (Vec) – The second metric to be averaged.

  • metricAvg (Vec) – The output averaged metric.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:3051

metricAverage3(metric1, metric2, metric3, metricAvg)#

Compute and return the unweighted average of three metrics.

Collective.

Parameters:
  • metric1 (Vec) – The first metric to be averaged.

  • metric2 (Vec) – The second metric to be averaged.

  • metric3 (Vec) – The third metric to be averaged.

  • metricAvg (Vec) – The output averaged metric.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:3073

metricCreate(field=0)#

Create a Riemannian metric field.

Collective.

Parameters:

field (int | None) – The field number to use.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:2882

metricCreateIsotropic(indicator, field=0)#

Construct an isotropic metric from an error indicator.

Collective.

Parameters:
  • indicator (Vec) – The error indicator.

  • field (int | None) – The field number to use.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:2927

metricCreateUniform(alpha, field=0)#

Construct a uniform isotropic metric.

Collective.

Parameters:
  • alpha (float) – Scaling parameter for the diagonal.

  • field (int | None) – The field number to use.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:2903

metricDeterminantCreate(field=0)#

Create the determinant field for a Riemannian metric.

Collective.

Parameters:

field (int | None) – The field number to use.

Returns:

  • determinant (Vec) – The determinant field.

  • dmDet (DM) – The corresponding DM

Return type:

tuple[Vec, DM]

Source code at petsc4py/PETSc/DMPlex.pyx:2950

metricEnforceSPD(metric, ometric, determinant, restrictSizes=False, restrictAnisotropy=False)#

Enforce symmetric positive-definiteness of a metric.

Collective.

Parameters:
  • metric (Vec) – The metric.

  • ometric (Vec) – The output metric.

  • determinant (Vec) – The output determinant.

  • restrictSizes (bool | None) – Flag indicating whether maximum/minimum magnitudes should be enforced.

  • restrictAnisotropy (bool | None) – Flag indicating whether maximum anisotropy should be enforced.

Returns:

  • ometric (Vec) – The output metric.

  • determinant (Vec) – The output determinant.

Return type:

tuple[Vec, Vec]

Source code at petsc4py/PETSc/DMPlex.pyx:2979

metricGetGradationFactor()#

Return the metric gradation factor.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2833

Return type:

float

metricGetHausdorffNumber()#

Return the metric Hausdorff number.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2867

Return type:

float

metricGetMaximumAnisotropy()#

Return the maximum tolerated metric anisotropy.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2731

Return type:

float

metricGetMaximumMagnitude()#

Return the maximum tolerated metric magnitude.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2697

Return type:

float

metricGetMinimumMagnitude()#

Return the minimum tolerated metric magnitude.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2663

Return type:

float

metricGetNormalizationOrder()#

Return the order p for L-p normalization.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2799

Return type:

float

metricGetNumIterations()#

Return the number of parallel adaptation iterations.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2629

Return type:

int

metricGetTargetComplexity()#

Return the target metric complexity.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2765

Return type:

float

metricGetVerbosity()#

Return the verbosity of the mesh adaptation package.

Not collective.

Returns:

verbosity – The verbosity, where -1 is silent and 10 is maximum.

Return type:

int

Source code at petsc4py/PETSc/DMPlex.pyx:2590

metricIntersection2(metric1, metric2, metricInt)#

Compute and return the intersection of two metrics.

Collective.

Parameters:
  • metric1 (Vec) – The first metric to be intersected.

  • metric2 (Vec) – The second metric to be intersected.

  • metricInt (Vec) – The output intersected metric.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:3097

metricIntersection3(metric1, metric2, metric3, metricInt)#

Compute the intersection of three metrics.

Collective.

Parameters:
  • metric1 (Vec) – The first metric to be intersected.

  • metric2 (Vec) – The second metric to be intersected.

  • metric3 (Vec) – The third metric to be intersected.

  • metricInt (Vec) – The output intersected metric.

Return type:

Vec

Source code at petsc4py/PETSc/DMPlex.pyx:3119

metricIsIsotropic()#

Return the flag indicating whether the metric is isotropic or not.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2378

Return type:

bool

metricIsUniform()#

Return the flag indicating whether the metric is uniform or not.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2344

Return type:

bool

metricNoInsertion()#

Return the flag indicating whether node insertion and deletion are turned off.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2447

Return type:

bool

metricNoMovement()#

Return the flag indicating whether node movement is turned off.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2519

Return type:

bool

metricNoSurf()#

Return the flag indicating whether surface modification is turned off.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2555

Return type:

bool

metricNoSwapping()#

Return the flag indicating whether facet swapping is turned off.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2483

Return type:

bool

metricNormalize(metric, ometric, determinant, restrictSizes=True, restrictAnisotropy=True)#

Apply L-p normalization to a metric.

Collective.

Parameters:
  • metric (Vec) – The metric.

  • ometric (Vec) – The output metric.

  • determinant (Vec) – The output determinant.

  • restrictSizes (bool | None) – Flag indicating whether maximum/minimum magnitudes should be enforced.

  • restrictAnisotropy (bool | None) – Flag indicating whether maximum anisotropy should be enforced.

Returns:

  • ometric (Vec) – The output normalized metric.

  • determinant (Vec) – The output determinant.

Return type:

tuple[Vec, Vec]

Source code at petsc4py/PETSc/DMPlex.pyx:3015

metricRestrictAnisotropyFirst()#

Return true if anisotropy is restricted before normalization.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2412

Return type:

bool

metricSetFromOptions()#

Configure the object from the options database.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2312

Return type:

None

metricSetGradationFactor(beta)#

Set the metric gradation factor.

Logically collective.

Parameters:

beta (float) – The metric gradation factor.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2814

metricSetHausdorffNumber(hausd)#

Set the metric Hausdorff number.

Logically collective.

Parameters:

hausd (float) – The metric Hausdorff number.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2848

metricSetIsotropic(isotropic)#

Record whether the metric is isotropic or not.

Logically collective.

Parameters:

isotropic (bool) – Flag indicating whether the metric is isotropic or not.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2359

metricSetMaximumAnisotropy(a_max)#

Set the maximum tolerated metric anisotropy.

Logically collective.

Parameters:

a_max (float) – The maximum tolerated metric anisotropy.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2712

metricSetMaximumMagnitude(h_max)#

Set the maximum tolerated metric magnitude.

Logically collective.

Parameters:

h_max (float) – The maximum tolerated metric magnitude.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2678

metricSetMinimumMagnitude(h_min)#

Set the minimum tolerated metric magnitude.

Logically collective.

Parameters:

h_min (float) – The minimum tolerated metric magnitude.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2644

metricSetNoInsertion(noInsert)#

Set the flag indicating whether node insertion should be turned off.

Logically collective.

Parameters:

noInsert (bool) – Flag indicating whether node insertion and deletion should be turned off.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2427

metricSetNoMovement(noMove)#

Set the flag indicating whether node movement should be turned off.

Logically collective.

Parameters:

noMove (bool) – Flag indicating whether node movement should be turned off.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2499

metricSetNoSurf(noSurf)#

Set the flag indicating whether surface modification should be turned off.

Logically collective.

Parameters:

noSurf (bool) – Flag indicating whether surface modification should be turned off.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2535

metricSetNoSwapping(noSwap)#

Set the flag indicating whether facet swapping should be turned off.

Logically collective.

Parameters:

noSwap (bool) – Flag indicating whether facet swapping should be turned off.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2463

metricSetNormalizationOrder(p)#

Set the order p for L-p normalization.

Logically collective.

Parameters:

p (float) – The normalization order.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2780

metricSetNumIterations(numIter)#

Set the number of parallel adaptation iterations.

Logically collective.

Parameters:

numIter (int) – The number of parallel adaptation iterations.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2610

metricSetRestrictAnisotropyFirst(restrictAnisotropyFirst)#

Record whether anisotropy is be restricted before normalization or after.

Logically collective.

Parameters:

restrictAnisotropyFirst (bool) – Flag indicating if anisotropy is restricted before normalization or after.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2393

metricSetTargetComplexity(targetComplexity)#

Set the target metric complexity.

Logically collective.

Parameters:

targetComplexity (float) – The target metric complexity.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2746

metricSetUniform(uniform)#

Record whether the metric is uniform or not.

Logically collective.

Parameters:

uniform (bool) – Flag indicating whether the metric is uniform or not.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2325

metricSetVerbosity(verbosity)#

Set the verbosity of the mesh adaptation package.

Logically collective.

Parameters:

verbosity (int) – The verbosity, where -1 is silent and 10 is maximum.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2571

orient()#

Give a consistent orientation to the input mesh.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:865

Return type:

None

permute(perm)#

Reorder the mesh according to the input permutation.

Collective.

Parameters:

perm (IS) – The point permutation, perm[old point number] = new point number.

Returns:

pdm – The permuted DMPlex.

Return type:

DMPlex

Source code at petsc4py/PETSc/DMPlex.pyx:2150

rebalanceSharedPoints(entityDepth=0, useInitialGuess=True, parallel=True)#

Redistribute shared points in order to achieve better balancing.

Collective.

Parameters:
  • entityDepth (int | None) – Depth of the entity to balance (e.g., 0 -> balance vertices).

  • useInitialGuess (bool | None) – Whether to use the current distribution as initial guess.

  • parallel (bool | None) – Whether to use ParMETIS and do the partition in parallel or gather the graph onto a single process.

Returns:

success – Whether the graph partitioning was successful or not. Unsuccessful simply means no change to the partitioning.

Return type:

bool

Source code at petsc4py/PETSc/DMPlex.pyx:1519

reorderGetDefault()#

Return flag indicating whether the DMPlex should be reordered by default.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:2174

Return type:

ReorderDefaultFlag

reorderSetDefault(flag)#

Set flag indicating whether the DM should be reordered by default.

Logically collective.

Parameters:
Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2188

sectionLoad(viewer, sectiondm, sfxc)#

Load section into a DM.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer that represents the on-disk section (sectionA).

  • sectiondm (DM) – The DM into which the on-disk section (sectionA) is migrated.

  • sfxc (SF) – The SF returned by topologyLoad.

Returns:

  • gsf (SF) – The SF that migrates any on-disk Vec data associated with sectionA into a global Vec associated with the sectiondm’s global section (None if not needed).

  • lsf (SF) – The SF that migrates any on-disk Vec data associated with sectionA into a local Vec associated with the sectiondm’s local section (None if not needed).

Return type:

tuple[SF, SF]

Source code at petsc4py/PETSc/DMPlex.pyx:3356

sectionView(viewer, sectiondm)#

Save a section associated with a DMPlex.

Collective.

Parameters:
  • viewer (Viewer) – The Viewer for saving.

  • sectiondm (DM) – The DM that contains the section to be saved.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3219

setAdjacencyUseAnchors(useAnchors=True)#

Define adjacency in the mesh using the point-to-point constraints.

Logically collective.

Parameters:

useAnchors (bool) – Flag to use the constraints. If True, then constrained points are omitted from DMPlex.getAdjacency, and their anchor points appear in their place.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1423

setCellType(p, ctype)#

Set the polytope type of a given cell.

Not collective.

Parameters:
  • p (int) – The cell.

  • ctype (PolytopeType) – The polytope type of the cell.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:655

setChart(pStart, pEnd)#

Set the interval for all mesh points [pStart, pEnd).

Not collective.

Parameters:
  • pStart (int) – The first mesh point.

  • pEnd (int) – The upper bound for mesh points.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:413

setCone(p, cone, orientation=None)#

Set the points on the in-edges for this point in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart.

  • cone (Sequence[int]) – An array of points which are on the in-edges for point p.

  • orientation (Sequence[int] | None) – An array of orientations, defaults to None.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:509

setConeOrientation(p, orientation)#

Set the orientations on the in-edges for this point in the DAG.

Not collective.

Parameters:
Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:624

setConeSize(p, size)#

Set the number of in-edges for this point in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart.

  • size (int) – The cone size for point p.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:458

setMatClosure(sec, gsec, mat, point, values, addv=None)#

Set an array of the values on the closure of point.

Not collective.

Parameters:
  • sec (Section) – The section describing the layout in mat, or None to use the default section.

  • gsec (Section) – The section describing the layout in mat, or None to use the default global section.

  • mat (Mat) – The matrix.

  • point (int) – The point in the DMPlex.

  • values (Sequence[Scalar]) – The array of values.

  • mode – The insertion mode.

  • addv (InsertModeSpec | None)

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1239

setPartitioner(part)#

Set the mesh partitioner.

Logically collective.

Parameters:

part (Partitioner) – The partitioner.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1484

setRefinementLimit(refinementLimit)#

Set the maximum cell volume for refinement.

Logically collective.

Parameters:

refinementLimit (float) – The maximum cell volume in the refined mesh.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2086

setRefinementUniform(refinementUniform=True)#

Set the flag for uniform refinement.

Logically collective.

Parameters:

refinementUniform (bool | None) – The flag for uniform refinement.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:2045

setSupport(p, supp)#

Set the points on the out-edges for this point in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart.

  • supp (Sequence[int]) – An array of points which are on the out-edges for point p.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:789

setSupportSize(p, size)#

Set the number of out-edges for this point in the DAG.

Not collective.

Parameters:
  • p (int) – The point, which must lie in the chart set with DMPlex.setChart.

  • size (int) – The support size for point p.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:738

setTetGenOptions(opts)#

Set the options used for the Tetgen mesh generator.

Not collective.

Parameters:

opts (str) – The command line options.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1324

setTriangleOptions(opts)#

Set the options used for the Triangle mesh generator.

Not collective.

Parameters:

opts (str) – The command line options.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1304

setVecClosure(sec, vec, point, values, addv=None)#

Set an array of the values on the closure of point.

Not collective.

Parameters:
  • sec (Section) – The section describing the layout in vec, or None to use the default section.

  • vec (Vec) – The local vector.

  • point (int) – The point in the DMPlex.

  • values (Sequence[Scalar]) – The array of values.

  • mode – The insertion mode.

  • addv (InsertModeSpec | None)

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:1207

stratify()#

Calculate the strata of DAG.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:853

Return type:

None

symmetrize()#

Create support (out-edge) information from cone (in-edge) information.

Not collective.

Source code at petsc4py/PETSc/DMPlex.pyx:840

Return type:

None

topologyLoad(viewer)#

Load a topology into this DMPlex object.

Collective.

Parameters:

viewer (Viewer) – The Viewer for the saved topology

Returns:

sfxc – The SF that pushes points in [0, N) to the associated points in the loaded DMPlex, where N is the global number of points.

Return type:

SF

Source code at petsc4py/PETSc/DMPlex.pyx:3290

topologyView(viewer)#

Save a DMPlex topology into a file.

Collective.

Parameters:

viewer (Viewer) – The Viewer for saving.

Return type:

None

Source code at petsc4py/PETSc/DMPlex.pyx:3165

uninterpolate()#

Convert to a mesh with only cells and vertices.

Collective.

Source code at petsc4py/PETSc/DMPlex.pyx:1735

Return type:

None

vecGetClosure(sec, vec, p)#

Return an array of values on the closure of p.

Not collective.

Parameters:
  • sec (Section) – The section describing the layout in vec.

  • vec (Vec) – The local vector.

  • p (int) – The point in the DMPlex.

Return type:

ArrayScalar

Source code at petsc4py/PETSc/DMPlex.pyx:1149