#include "petscdmplex.h" PetscErrorCode DMPlexVecGetClosure(DM dm, PetscSection section, Vec v, PetscInt point, PetscInt *csize, PetscScalar *values[])Not collective
dm | - The DM | |
section | - The section describing the layout in v, or NULL to use the default section | |
v | - The local vector | |
point | - The point in the DM | |
csize | - The size of the input values array, or NULL | |
values | - An array to use for the values, or NULL to have it allocated automatically |
csize | - The number of values in the closure | |
values | - The array of values. If the user provided NULL, it is a borrowed array and should not be freed |
Note that DMPlexVecGetClosure/DMPlexVecRestoreClosure only allocates the values array if it set to NULL in the
calling function. This is because DMPlexVecGetClosure() is typically called in the inner loop of a Vec or Mat
assembly function, and a user may already have allocated storage for this operation.
A typical use could be
values = NULL;
ierr = DMPlexVecGetClosure(dm, NULL, v, p, &clSize, &values);CHKERRQ(ierr);
for (cl = 0; cl < clSize; ++cl) {
<Compute on closure>
}
ierr = DMPlexVecRestoreClosure(dm, NULL, v, p, &clSize, &values);CHKERRQ(ierr);
or
PetscMalloc1(clMaxSize, &values);
for (p = pStart; p < pEnd; ++p) {
clSize = clMaxSize;
ierr = DMPlexVecGetClosure(dm, NULL, v, p, &clSize, &values);CHKERRQ(ierr);
for (cl = 0; cl < clSize; ++cl) {
<Compute on closure>
}
}
PetscFree(values);
The csize argument is not present in the Fortran 90 binding since it is internal to the array.