Actual source code: relative.c
2: #include <petsc/private/vecimpl.h>
3: #include "../src/vec/vec/utils/tagger/impls/simple.h"
5: static PetscErrorCode VecTaggerComputeBoxes_Relative(VecTagger tagger, Vec vec, PetscInt *numBoxes, VecTaggerBox **boxes, PetscBool *listed)
6: {
7: VecTagger_Simple *smpl = (VecTagger_Simple *)tagger->data;
8: PetscInt bs, i, j, k, n;
9: VecTaggerBox *bxs;
10: const PetscScalar *vArray;
12: PetscFunctionBegin;
13: PetscCall(VecTaggerGetBlockSize(tagger, &bs));
14: *numBoxes = 1;
15: PetscCall(PetscMalloc1(bs, &bxs));
16: PetscCall(VecGetLocalSize(vec, &n));
17: n /= bs;
18: for (i = 0; i < bs; i++) {
19: #if !defined(PETSC_USE_COMPLEX)
20: bxs[i].min = PETSC_MAX_REAL;
21: bxs[i].max = PETSC_MIN_REAL;
22: #else
23: bxs[i].min = PetscCMPLX(PETSC_MAX_REAL, PETSC_MAX_REAL);
24: bxs[i].max = PetscCMPLX(PETSC_MIN_REAL, PETSC_MIN_REAL);
25: #endif
26: }
27: PetscCall(VecGetArrayRead(vec, &vArray));
28: for (i = 0, k = 0; i < n; i++) {
29: for (j = 0; j < bs; j++, k++) {
30: #if !defined(PETSC_USE_COMPLEX)
31: bxs[j].min = PetscMin(bxs[j].min, vArray[k]);
32: bxs[j].max = PetscMax(bxs[j].max, vArray[k]);
33: #else
34: bxs[j].min = PetscCMPLX(PetscMin(PetscRealPart(bxs[j].min), PetscRealPart(vArray[k])), PetscMin(PetscImaginaryPart(bxs[j].min), PetscImaginaryPart(vArray[k])));
35: bxs[j].max = PetscCMPLX(PetscMax(PetscRealPart(bxs[j].max), PetscRealPart(vArray[k])), PetscMax(PetscImaginaryPart(bxs[j].max), PetscImaginaryPart(vArray[k])));
36: #endif
37: }
38: }
39: for (i = 0; i < bs; i++) bxs[i].max = -bxs[i].max;
40: PetscCall(VecRestoreArrayRead(vec, &vArray));
41: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, (PetscReal *)bxs, 2 * (sizeof(PetscScalar) / sizeof(PetscReal)) * bs, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)tagger)));
42: for (i = 0; i < bs; i++) {
43: PetscScalar mins = bxs[i].min;
44: PetscScalar difs = -bxs[i].max - mins;
45: #if !defined(PETSC_USE_COMPLEX)
46: bxs[i].min = mins + smpl->box[i].min * difs;
47: bxs[i].max = mins + smpl->box[i].max * difs;
48: #else
49: bxs[i].min = mins + PetscCMPLX(PetscRealPart(smpl->box[i].min) * PetscRealPart(difs), PetscImaginaryPart(smpl->box[i].min) * PetscImaginaryPart(difs));
50: bxs[i].max = mins + PetscCMPLX(PetscRealPart(smpl->box[i].max) * PetscRealPart(difs), PetscImaginaryPart(smpl->box[i].max) * PetscImaginaryPart(difs));
51: #endif
52: }
53: *boxes = bxs;
54: if (listed) *listed = PETSC_TRUE;
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }
58: /*@C
59: VecTaggerRelativeSetBox - Set the relative box defining the values to be tagged by the tagger, where relative boxes are subsets of [0,1] (or [0,1]+[0,1]i for complex scalars), where 0 indicates the smallest value present in the vector and 1 indicates the largest.
61: Logically Collective
63: Input Parameters:
64: + tagger - the VecTagger context
65: - box - a blocksize list of VecTaggerBox boxes
67: Level: advanced
69: .seealso: `VecTaggerRelativeGetBox()`
70: @*/
71: PetscErrorCode VecTaggerRelativeSetBox(VecTagger tagger, VecTaggerBox *box)
72: {
73: PetscFunctionBegin;
74: PetscCall(VecTaggerSetBox_Simple(tagger, box));
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: /*@C
79: VecTaggerRelativeGetBox - Get the relative box defining the values to be tagged by the tagger, where relative boxess are subsets of [0,1] (or [0,1]+[0,1]i for complex scalars), where 0 indicates the smallest value present in the vector and 1 indicates the largest.
81: Logically Collective
83: Input Parameter:
84: . tagger - the VecTagger context
86: Output Parameter:
87: . box - a blocksize list of VecTaggerBox boxes
89: Level: advanced
91: .seealso: `VecTaggerRelativeSetBox()`
92: @*/
93: PetscErrorCode VecTaggerRelativeGetBox(VecTagger tagger, const VecTaggerBox **box)
94: {
95: PetscFunctionBegin;
96: PetscCall(VecTaggerGetBox_Simple(tagger, box));
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: PETSC_INTERN PetscErrorCode VecTaggerCreate_Relative(VecTagger tagger)
101: {
102: PetscFunctionBegin;
103: PetscCall(VecTaggerCreate_Simple(tagger));
104: tagger->ops->computeboxes = VecTaggerComputeBoxes_Relative;
105: PetscFunctionReturn(PETSC_SUCCESS);
106: }