Actual source code: absolute.c
1: #include <petsc/private/vecimpl.h>
2: #include "../src/vec/vec/utils/tagger/impls/simple.h"
4: static PetscErrorCode VecTaggerComputeBoxes_Absolute(VecTagger tagger, Vec vec, PetscInt *numBoxes, VecTaggerBox **boxes, PetscBool *listed)
5: {
6: VecTagger_Simple *smpl = (VecTagger_Simple *)tagger->data;
7: PetscInt bs, i;
8: VecTaggerBox *bxs;
10: PetscFunctionBegin;
11: PetscCall(VecTaggerGetBlockSize(tagger, &bs));
12: *numBoxes = 1;
13: PetscCall(PetscMalloc1(bs, &bxs));
14: for (i = 0; i < bs; i++) {
15: bxs[i].min = smpl->box[i].min;
16: bxs[i].max = smpl->box[i].max;
17: }
18: *boxes = bxs;
19: if (listed) *listed = PETSC_TRUE;
20: PetscFunctionReturn(PETSC_SUCCESS);
21: }
23: /*@C
24: VecTaggerAbsoluteSetBox - Set the box defining the values to be tagged by the tagger.
26: Logically Collective
28: Input Parameters:
29: + tagger - the `VecTagger` context
30: - box - the box: a blocksize array of `VecTaggerBox` boxes
32: Level: advanced
34: .seealso: `VecTagger`, `VecTaggerBox`, `VecTaggerAbsoluteGetBox()`
35: @*/
36: PetscErrorCode VecTaggerAbsoluteSetBox(VecTagger tagger, VecTaggerBox box[])
37: {
38: PetscFunctionBegin;
39: PetscCall(VecTaggerSetBox_Simple(tagger, box));
40: PetscFunctionReturn(PETSC_SUCCESS);
41: }
43: /*@C
44: VecTaggerAbsoluteGetBox - Get the box defining the values to be tagged by the tagger.
46: Logically Collective
48: Input Parameter:
49: . tagger - the `VecTagger` context
51: Output Parameter:
52: . box - the box: a blocksize array of `VecTaggerBox` boxes
54: Level: advanced
56: .seealso: `VecTagger`, `VecTaggerBox`, `VecTaggerAbsoluteSetBox()`
57: @*/
58: PetscErrorCode VecTaggerAbsoluteGetBox(VecTagger tagger, const VecTaggerBox *box[])
59: {
60: PetscFunctionBegin;
61: PetscCall(VecTaggerGetBox_Simple(tagger, box));
62: PetscFunctionReturn(PETSC_SUCCESS);
63: }
65: PETSC_INTERN PetscErrorCode VecTaggerCreate_Absolute(VecTagger tagger)
66: {
67: PetscFunctionBegin;
68: PetscCall(VecTaggerCreate_Simple(tagger));
69: tagger->ops->computeboxes = VecTaggerComputeBoxes_Absolute;
70: PetscFunctionReturn(PETSC_SUCCESS);
71: }