Actual source code: dapreallocate.c

  1: #include <petsc/private/dmdaimpl.h>
  2: #include <petsc/private/isimpl.h>
  3: #include <petscsf.h>

  5: /*@
  6:   DMDASetPreallocationCenterDimension - Determine the topology used to determine adjacency

  8:   Input Parameters:
  9: + dm - The DM object
 10: - preallocCenterDim - The dimension of points which connect adjacent entries

 12:   Level: developer

 14:   Notes:
 15: $     FEM:   Two points p and q are adjacent if q \in closure(star(p)), preallocCenterDim = dim
 16: $     FVM:   Two points p and q are adjacent if q \in star(cone(p)),    preallocCenterDim = dim-1
 17: $     FVM++: Two points p and q are adjacent if q \in star(closure(p)), preallocCenterDim = 0

 19: .seealso: `DMCreateMatrix()`, `DMDAPreallocateOperator()`
 20: @*/
 21: PetscErrorCode DMDASetPreallocationCenterDimension(DM dm, PetscInt preallocCenterDim)
 22: {
 23:   DM_DA *mesh = (DM_DA *)dm->data;

 26:   mesh->preallocCenterDim = preallocCenterDim;
 27:   return 0;
 28: }

 30: /*@
 31:   DMDAGetPreallocationCenterDimension - Return the topology used to determine adjacency

 33:   Input Parameter:
 34: . dm - The DM object

 36:   Output Parameter:
 37: . preallocCenterDim - The dimension of points which connect adjacent entries

 39:   Level: developer

 41:   Notes:
 42: $     FEM:   Two points p and q are adjacent if q \in closure(star(p)), preallocCenterDim = dim
 43: $     FVM:   Two points p and q are adjacent if q \in star(cone(p)),    preallocCenterDim = dim-1
 44: $     FVM++: Two points p and q are adjacent if q \in star(closure(p)), preallocCenterDim = 0

 46: .seealso: `DMCreateMatrix()`, `DMDAPreallocateOperator()`, `DMDASetPreallocationCenterDimension()`
 47: @*/
 48: PetscErrorCode DMDAGetPreallocationCenterDimension(DM dm, PetscInt *preallocCenterDim)
 49: {
 50:   DM_DA *mesh = (DM_DA *)dm->data;

 54:   *preallocCenterDim = mesh->preallocCenterDim;
 55:   return 0;
 56: }