Actual source code: relative.c
petsc-3.13.6 2020-09-29
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)
6: {
7: VecTagger_Simple *smpl = (VecTagger_Simple *)tagger->data;
8: PetscInt bs, i, j, k, n;
9: VecTaggerBox *bxs;
10: const PetscScalar *vArray;
11: PetscErrorCode ierr;
14: VecTaggerGetBlockSize(tagger,&bs);
15: *numBoxes = 1;
16: PetscMalloc1(bs,&bxs);
17: VecGetLocalSize(vec,&n);
18: n /= bs;
19: for (i = 0; i < bs; i++) {
20: #if !defined(PETSC_USE_COMPLEX)
21: bxs[i].min = PETSC_MAX_REAL;
22: bxs[i].max = PETSC_MIN_REAL;
23: #else
24: bxs[i].min = PetscCMPLX(PETSC_MAX_REAL,PETSC_MAX_REAL);
25: bxs[i].max = PetscCMPLX(PETSC_MIN_REAL,PETSC_MIN_REAL);
26: #endif
27: }
28: VecGetArrayRead(vec, &vArray);
29: for (i = 0, k = 0; i < n; i++) {
30: for (j = 0; j < bs; j++, k++) {
31: #if !defined(PETSC_USE_COMPLEX)
32: bxs[j].min = PetscMin(bxs[j].min,vArray[k]);
33: bxs[j].max = PetscMax(bxs[j].max,vArray[k]);
34: #else
35: bxs[j].min = PetscCMPLX(PetscMin(PetscRealPart(bxs[j].min),PetscRealPart(vArray[k])),PetscMin(PetscImaginaryPart(bxs[j].min),PetscImaginaryPart(vArray[k])));
36: bxs[j].max = PetscCMPLX(PetscMax(PetscRealPart(bxs[j].max),PetscRealPart(vArray[k])),PetscMax(PetscImaginaryPart(bxs[j].max),PetscImaginaryPart(vArray[k])));
37: #endif
38: }
39: }
40: for (i = 0; i < bs; i++) bxs[i].max = -bxs[i].max;
41: VecRestoreArrayRead(vec, &vArray);
42: MPI_Allreduce(MPI_IN_PLACE,(PetscReal *) bxs,2*(sizeof(PetscScalar)/sizeof(PetscReal))*bs,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)tagger));
43: for (i = 0; i < bs; i++) {
44: PetscScalar mins = bxs[i].min;
45: PetscScalar difs = -bxs[i].max - mins;
46: #if !defined(PETSC_USE_COMPLEX)
47: bxs[i].min = mins + smpl->box[i].min * difs;
48: bxs[i].max = mins + smpl->box[i].max * difs;
49: #else
50: bxs[i].min = mins + PetscCMPLX(PetscRealPart(smpl->box[i].min)*PetscRealPart(difs),PetscImaginaryPart(smpl->box[i].min)*PetscImaginaryPart(difs));
51: bxs[i].max = mins + PetscCMPLX(PetscRealPart(smpl->box[i].max)*PetscRealPart(difs),PetscImaginaryPart(smpl->box[i].max)*PetscImaginaryPart(difs));
52: #endif
53: }
54: *boxes = bxs;
55: return(0);
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 Arguments:
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: {
76: VecTaggerSetBox_Simple(tagger,box);
77: return(0);
78: }
80: /*@C
81: 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.
83: Logically Collective
85: Input Arguments:
86: . tagger - the VecTagger context
88: Output Arguments:
89: . box - a blocksize list of VecTaggerBox boxes
91: Level: advanced
93: .seealso: VecTaggerRelativeSetBox()
94: @*/
95: PetscErrorCode VecTaggerRelativeGetBox(VecTagger tagger,const VecTaggerBox **box)
96: {
100: VecTaggerGetBox_Simple(tagger,box);
101: return(0);
102: }
104: PETSC_INTERN PetscErrorCode VecTaggerCreate_Relative(VecTagger tagger)
105: {
106: PetscErrorCode ierr;
109: VecTaggerCreate_Simple(tagger);
110: tagger->ops->computeboxes = VecTaggerComputeBoxes_Relative;
111: return(0);
112: }