xref: /petsc/src/vec/vec/utils/tagger/impls/cdf.c (revision bcda9346efad4e5ba2d553af84eb238771ba1e25)
1 #include <petsc/private/vecimpl.h> /*I "petscvec.h" I*/
2 #include "../src/vec/vec/utils/tagger/impls/simple.h"
3 
4 const char *const VecTaggerCDFMethods[VECTAGGER_CDF_NUM_METHODS] = {"gather", "iterative"};
5 
6 #if !defined(PETSC_USE_COMPLEX)
7 typedef VecTaggerBox VecTaggerBoxReal;
8 #else
9 typedef struct {
10   PetscReal min;
11   PetscReal max;
12 } VecTaggerBoxReal;
13 #endif
14 
15 typedef struct {
16   VecTagger_Simple   smpl;
17   PetscReal          atol;
18   PetscReal          rtol;
19   PetscInt           maxit;
20   PetscInt           numMoments;
21   VecTaggerCDFMethod method;
22 } VecTagger_CDF;
23 
VecTaggerComputeBox_CDF_SortedArray(const PetscReal * cArray,PetscInt m,const VecTaggerBoxReal * bxs,VecTaggerBoxReal * boxes)24 static PetscErrorCode VecTaggerComputeBox_CDF_SortedArray(const PetscReal *cArray, PetscInt m, const VecTaggerBoxReal *bxs, VecTaggerBoxReal *boxes)
25 {
26   PetscInt  minInd, maxInd;
27   PetscReal minCDF, maxCDF;
28 
29   PetscFunctionBegin;
30   minCDF     = PetscMax(0., bxs->min);
31   maxCDF     = PetscMin(1., bxs->max);
32   minInd     = (PetscInt)(minCDF * m);
33   maxInd     = (PetscInt)(maxCDF * m);
34   boxes->min = cArray[PetscMin(minInd, m - 1)];
35   boxes->max = cArray[PetscMax(minInd, maxInd - 1)];
36   PetscFunctionReturn(PETSC_SUCCESS);
37 }
38 
VecTaggerComputeBoxes_CDF_Serial(VecTagger tagger,Vec vec,PetscInt bs,VecTaggerBox * boxes)39 static PetscErrorCode VecTaggerComputeBoxes_CDF_Serial(VecTagger tagger, Vec vec, PetscInt bs, VecTaggerBox *boxes)
40 {
41   VecTagger_Simple *smpl = (VecTagger_Simple *)tagger->data;
42   Vec               vComp;
43   PetscInt          n, m;
44   PetscInt          i;
45 #if defined(PETSC_USE_COMPLEX)
46   PetscReal *cReal, *cImag;
47 #endif
48 
49   PetscFunctionBegin;
50   PetscCall(VecGetLocalSize(vec, &n));
51   m = n / bs;
52   PetscCall(VecCreateSeq(PETSC_COMM_SELF, m, &vComp));
53 #if defined(PETSC_USE_COMPLEX)
54   PetscCall(PetscMalloc2(m, &cReal, m, &cImag));
55 #endif
56   for (i = 0; i < bs; i++) {
57     IS           isStride;
58     VecScatter   vScat;
59     PetscScalar *cArray;
60 
61     PetscCall(ISCreateStride(PETSC_COMM_SELF, m, i, bs, &isStride));
62     PetscCall(VecScatterCreate(vec, isStride, vComp, NULL, &vScat));
63     PetscCall(VecScatterBegin(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
64     PetscCall(VecScatterEnd(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
65     PetscCall(VecScatterDestroy(&vScat));
66     PetscCall(ISDestroy(&isStride));
67 
68     PetscCall(VecGetArray(vComp, &cArray));
69 #if !defined(PETSC_USE_COMPLEX)
70     PetscCall(PetscSortReal(m, cArray));
71     PetscCall(VecTaggerComputeBox_CDF_SortedArray(cArray, m, &smpl->box[i], &boxes[i]));
72 #else
73     {
74       PetscInt         j;
75       VecTaggerBoxReal realBxs, imagBxs;
76       VecTaggerBoxReal realBoxes, imagBoxes;
77 
78       for (j = 0; j < m; j++) {
79         cReal[j] = PetscRealPart(cArray[j]);
80         cImag[j] = PetscImaginaryPart(cArray[j]);
81       }
82       PetscCall(PetscSortReal(m, cReal));
83       PetscCall(PetscSortReal(m, cImag));
84 
85       realBxs.min = PetscRealPart(smpl->box[i].min);
86       realBxs.max = PetscRealPart(smpl->box[i].max);
87       imagBxs.min = PetscImaginaryPart(smpl->box[i].min);
88       imagBxs.max = PetscImaginaryPart(smpl->box[i].max);
89       PetscCall(VecTaggerComputeBox_CDF_SortedArray(cReal, m, &realBxs, &realBoxes));
90       PetscCall(VecTaggerComputeBox_CDF_SortedArray(cImag, m, &imagBxs, &imagBoxes));
91       boxes[i].min = PetscCMPLX(realBoxes.min, imagBoxes.min);
92       boxes[i].max = PetscCMPLX(realBoxes.max, imagBoxes.max);
93     }
94 #endif
95     PetscCall(VecRestoreArray(vComp, &cArray));
96   }
97 #if defined(PETSC_USE_COMPLEX)
98   PetscCall(PetscFree2(cReal, cImag));
99 #endif
100   PetscCall(VecDestroy(&vComp));
101   PetscFunctionReturn(PETSC_SUCCESS);
102 }
103 
VecTaggerComputeBoxes_CDF_Gather(VecTagger tagger,Vec vec,PetscInt bs,VecTaggerBox * boxes)104 static PetscErrorCode VecTaggerComputeBoxes_CDF_Gather(VecTagger tagger, Vec vec, PetscInt bs, VecTaggerBox *boxes)
105 {
106   Vec         gVec = NULL;
107   VecScatter  vScat;
108   PetscMPIInt rank, bs2;
109 
110   PetscFunctionBegin;
111   PetscCall(VecScatterCreateToZero(vec, &vScat, &gVec));
112   PetscCall(VecScatterBegin(vScat, vec, gVec, INSERT_VALUES, SCATTER_FORWARD));
113   PetscCall(VecScatterEnd(vScat, vec, gVec, INSERT_VALUES, SCATTER_FORWARD));
114   PetscCall(VecScatterDestroy(&vScat));
115   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)vec), &rank));
116   if (rank == 0) PetscCall(VecTaggerComputeBoxes_CDF_Serial(tagger, gVec, bs, boxes));
117   PetscCall(PetscMPIIntCast(2 * bs, &bs2));
118   PetscCallMPI(MPI_Bcast(boxes, bs2, MPIU_SCALAR, 0, PetscObjectComm((PetscObject)vec)));
119   PetscCall(VecDestroy(&gVec));
120   PetscFunctionReturn(PETSC_SUCCESS);
121 }
122 
123 typedef struct _n_CDFStats {
124   PetscReal min;
125   PetscReal max;
126   PetscReal moment[3];
127 } CDFStats;
128 
VecTaggerCDFStatsReduce(void * a,void * b,int * len,MPI_Datatype * datatype)129 static void MPIAPI VecTaggerCDFStatsReduce(void *a, void *b, int *len, MPI_Datatype *datatype)
130 {
131   PetscInt  i, j, N = *len;
132   CDFStats *A = (CDFStats *)a;
133   CDFStats *B = (CDFStats *)b;
134 
135   for (i = 0; i < N; i++) {
136     B[i].min = PetscMin(A[i].min, B[i].min);
137     B[i].max = PetscMax(A[i].max, B[i].max);
138     for (j = 0; j < 3; j++) B[i].moment[j] += A[i].moment[j];
139   }
140 }
141 
CDFUtilInverseEstimate(const CDFStats * stats,PetscReal cdfTarget,PetscReal * absEst)142 static PetscErrorCode CDFUtilInverseEstimate(const CDFStats *stats, PetscReal cdfTarget, PetscReal *absEst)
143 {
144   PetscReal min, max;
145 
146   PetscFunctionBegin;
147   min     = stats->min;
148   max     = stats->max;
149   *absEst = min + cdfTarget * (max - min);
150   PetscFunctionReturn(PETSC_SUCCESS);
151 }
152 
VecTaggerComputeBox_CDF_SortedArray_Iterative(VecTagger tagger,MPI_Datatype statType,MPI_Op statReduce,const PetscReal * cArray,PetscInt m,const VecTaggerBoxReal * cdfBox,VecTaggerBoxReal * absBox)153 static PetscErrorCode VecTaggerComputeBox_CDF_SortedArray_Iterative(VecTagger tagger, MPI_Datatype statType, MPI_Op statReduce, const PetscReal *cArray, PetscInt m, const VecTaggerBoxReal *cdfBox, VecTaggerBoxReal *absBox)
154 {
155   MPI_Comm       comm;
156   VecTagger_CDF *cdf;
157   PetscInt       maxit, i, j, k, l, M;
158   PetscInt       bounds[2][2];
159   PetscInt       offsets[2];
160   PetscReal      intervalLen = cdfBox->max - cdfBox->min;
161   PetscReal      rtol, atol;
162 
163   PetscFunctionBegin;
164   comm  = PetscObjectComm((PetscObject)tagger);
165   cdf   = (VecTagger_CDF *)tagger->data;
166   maxit = cdf->maxit;
167   rtol  = cdf->rtol;
168   atol  = cdf->atol;
169   /* local range of sorted values that can contain the sought radix */
170   offsets[0]   = 0;
171   offsets[1]   = 0;
172   bounds[0][0] = 0;
173   bounds[0][1] = m;
174   bounds[1][0] = 0;
175   bounds[1][1] = m;
176   PetscCall(VecTaggerComputeBox_CDF_SortedArray(cArray, m, cdfBox, absBox)); /* compute a local estimate of the interval */
177   {
178     CDFStats stats[3];
179 
180     for (i = 0; i < 2; i++) { /* compute statistics of those local estimates */
181       PetscReal val = i ? absBox->max : absBox->min;
182 
183       stats[i].min       = m ? val : PETSC_MAX_REAL;
184       stats[i].max       = m ? val : PETSC_MIN_REAL;
185       stats[i].moment[0] = m;
186       stats[i].moment[1] = m * val;
187       stats[i].moment[2] = m * val * val;
188     }
189     stats[2].min = PETSC_MAX_REAL;
190     stats[2].max = PETSC_MAX_REAL;
191     for (i = 0; i < 3; i++) stats[2].moment[i] = 0.;
192     for (i = 0; i < m; i++) {
193       PetscReal val = cArray[i];
194 
195       stats[2].min = PetscMin(stats[2].min, val);
196       stats[2].max = PetscMax(stats[2].max, val);
197       stats[2].moment[0]++;
198       stats[2].moment[1] += val;
199       stats[2].moment[2] += val * val;
200     }
201     /* reduce those statistics */
202     PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, stats, 3, statType, statReduce, comm));
203     M = (PetscInt)stats[2].moment[0];
204     /* use those initial statistics to get the initial (globally agreed-upon) choices for the absolute box bounds */
205     for (i = 0; i < 2; i++) PetscCall(CDFUtilInverseEstimate(&stats[i], i ? cdfBox->max : cdfBox->min, i ? &absBox->max : &absBox->min));
206   }
207   /* refine the estimates by computing how close they come to the desired box and refining */
208   for (k = 0; k < maxit; k++) {
209     PetscReal maxDiff = 0.;
210 
211     CDFStats stats[2][2];
212     PetscInt newBounds[2][2][2];
213     for (i = 0; i < 2; i++) {
214       for (j = 0; j < 2; j++) {
215         stats[i][j].min = PETSC_MAX_REAL;
216         stats[i][j].max = PETSC_MIN_REAL;
217         for (l = 0; l < 3; l++) stats[i][j].moment[l] = 0.;
218         newBounds[i][j][0] = PetscMax(bounds[i][0], bounds[i][1]);
219         newBounds[i][j][1] = PetscMin(bounds[i][0], bounds[i][1]);
220       }
221     }
222     for (i = 0; i < 2; i++) {
223       for (j = 0; j < bounds[i][1] - bounds[i][0]; j++) {
224         PetscInt  thisInd = bounds[i][0] + j;
225         PetscReal val     = cArray[thisInd];
226         PetscInt  section;
227         if (!i) {
228           section = (val < absBox->min) ? 0 : 1;
229         } else {
230           section = (val <= absBox->max) ? 0 : 1;
231         }
232         stats[i][section].min = PetscMin(stats[i][section].min, val);
233         stats[i][section].max = PetscMax(stats[i][section].max, val);
234         stats[i][section].moment[0]++;
235         stats[i][section].moment[1] += val;
236         stats[i][section].moment[2] += val * val;
237         newBounds[i][section][0] = PetscMin(newBounds[i][section][0], thisInd);
238         newBounds[i][section][1] = PetscMax(newBounds[i][section][0], thisInd + 1);
239       }
240     }
241     PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, stats, 4, statType, statReduce, comm));
242     for (i = 0; i < 2; i++) {
243       PetscInt  totalLessThan = offsets[i] + stats[i][0].moment[0];
244       PetscReal cdfOfAbs      = (PetscReal)totalLessThan / (PetscReal)M;
245       PetscReal diff;
246       PetscInt  section;
247 
248       if (cdfOfAbs == (i ? cdfBox->max : cdfBox->min)) {
249         offsets[i]   = totalLessThan;
250         bounds[i][0] = bounds[i][1] = 0;
251         continue;
252       }
253       if (cdfOfAbs > (i ? cdfBox->max : cdfBox->min)) { /* the correct absolute value lies in the lower section */
254         section = 0;
255       } else {
256         section    = 1;
257         offsets[i] = totalLessThan;
258       }
259       for (j = 0; j < 2; j++) bounds[i][j] = newBounds[i][section][j];
260       PetscCall(CDFUtilInverseEstimate(&stats[i][section], ((i ? cdfBox->max : cdfBox->min) - ((PetscReal)offsets[i] / (PetscReal)M)) / stats[i][section].moment[0], i ? &absBox->max : &absBox->min));
261       diff    = PetscAbs(cdfOfAbs - (i ? cdfBox->max : cdfBox->min));
262       maxDiff = PetscMax(maxDiff, diff);
263     }
264     if (!maxDiff) PetscFunctionReturn(PETSC_SUCCESS);
265     if ((atol || rtol) && ((!atol) || (maxDiff <= atol)) && ((!rtol) || (maxDiff <= rtol * intervalLen))) break;
266   }
267   PetscFunctionReturn(PETSC_SUCCESS);
268 }
269 
VecTaggerComputeBoxes_CDF_Iterative(VecTagger tagger,Vec vec,PetscInt bs,VecTaggerBox * boxes)270 static PetscErrorCode VecTaggerComputeBoxes_CDF_Iterative(VecTagger tagger, Vec vec, PetscInt bs, VecTaggerBox *boxes)
271 {
272   VecTagger_CDF    *cdf  = (VecTagger_CDF *)tagger->data;
273   VecTagger_Simple *smpl = &cdf->smpl;
274   Vec               vComp;
275   PetscInt          i, N, M, n, m, rstart;
276 #if defined(PETSC_USE_COMPLEX)
277   PetscReal *cReal, *cImag;
278 #endif
279   MPI_Comm     comm;
280   MPI_Datatype statType;
281   MPI_Op       statReduce;
282 
283   PetscFunctionBegin;
284   comm = PetscObjectComm((PetscObject)vec);
285   PetscCall(VecGetSize(vec, &N));
286   PetscCall(VecGetLocalSize(vec, &n));
287   M = N / bs;
288   m = n / bs;
289   PetscCall(VecCreateMPI(comm, m, M, &vComp));
290   PetscCall(VecSetUp(vComp));
291   PetscCall(VecGetOwnershipRange(vComp, &rstart, NULL));
292 #if defined(PETSC_USE_COMPLEX)
293   PetscCall(PetscMalloc2(m, &cReal, m, &cImag));
294 #endif
295   PetscCallMPI(MPI_Type_contiguous(5, MPIU_REAL, &statType));
296   PetscCallMPI(MPI_Type_commit(&statType));
297   PetscCallMPI(MPI_Op_create(VecTaggerCDFStatsReduce, 1, &statReduce));
298   for (i = 0; i < bs; i++) {
299     IS           isStride;
300     VecScatter   vScat;
301     PetscScalar *cArray;
302 
303     PetscCall(ISCreateStride(comm, m, bs * rstart + i, bs, &isStride));
304     PetscCall(VecScatterCreate(vec, isStride, vComp, NULL, &vScat));
305     PetscCall(VecScatterBegin(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
306     PetscCall(VecScatterEnd(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
307     PetscCall(VecScatterDestroy(&vScat));
308     PetscCall(ISDestroy(&isStride));
309 
310     PetscCall(VecGetArray(vComp, &cArray));
311 #if !defined(PETSC_USE_COMPLEX)
312     PetscCall(PetscSortReal(m, cArray));
313     PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cArray, m, &smpl->box[i], &boxes[i]));
314 #else
315     {
316       PetscInt         j;
317       VecTaggerBoxReal realBxs, imagBxs;
318       VecTaggerBoxReal realBoxes, imagBoxes;
319 
320       for (j = 0; j < m; j++) {
321         cReal[j] = PetscRealPart(cArray[j]);
322         cImag[j] = PetscImaginaryPart(cArray[j]);
323       }
324       PetscCall(PetscSortReal(m, cReal));
325       PetscCall(PetscSortReal(m, cImag));
326 
327       realBxs.min = PetscRealPart(smpl->box[i].min);
328       realBxs.max = PetscRealPart(smpl->box[i].max);
329       imagBxs.min = PetscImaginaryPart(smpl->box[i].min);
330       imagBxs.max = PetscImaginaryPart(smpl->box[i].max);
331       PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cReal, m, &realBxs, &realBoxes));
332       PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cImag, m, &imagBxs, &imagBoxes));
333       boxes[i].min = PetscCMPLX(realBoxes.min, imagBoxes.min);
334       boxes[i].max = PetscCMPLX(realBoxes.max, imagBoxes.max);
335     }
336 #endif
337     PetscCall(VecRestoreArray(vComp, &cArray));
338   }
339   PetscCallMPI(MPI_Op_free(&statReduce));
340   PetscCallMPI(MPI_Type_free(&statType));
341 #if defined(PETSC_USE_COMPLEX)
342   PetscCall(PetscFree2(cReal, cImag));
343 #endif
344   PetscCall(VecDestroy(&vComp));
345   PetscFunctionReturn(PETSC_SUCCESS);
346 }
347 
VecTaggerComputeBoxes_CDF(VecTagger tagger,Vec vec,PetscInt * numBoxes,VecTaggerBox ** boxes,PetscBool * listed)348 static PetscErrorCode VecTaggerComputeBoxes_CDF(VecTagger tagger, Vec vec, PetscInt *numBoxes, VecTaggerBox **boxes, PetscBool *listed)
349 {
350   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
351   PetscMPIInt    size;
352   PetscInt       bs;
353   VecTaggerBox  *bxs;
354 
355   PetscFunctionBegin;
356   PetscCall(VecTaggerGetBlockSize(tagger, &bs));
357   *numBoxes = 1;
358   PetscCall(PetscMalloc1(bs, &bxs));
359   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)tagger), &size));
360   if (size == 1) {
361     PetscCall(VecTaggerComputeBoxes_CDF_Serial(tagger, vec, bs, bxs));
362     *boxes = bxs;
363     PetscFunctionReturn(PETSC_SUCCESS);
364   }
365   switch (cuml->method) {
366   case VECTAGGER_CDF_GATHER:
367     PetscCall(VecTaggerComputeBoxes_CDF_Gather(tagger, vec, bs, bxs));
368     break;
369   case VECTAGGER_CDF_ITERATIVE:
370     PetscCall(VecTaggerComputeBoxes_CDF_Iterative(tagger, vec, bs, bxs));
371     break;
372   default:
373     SETERRQ(PetscObjectComm((PetscObject)tagger), PETSC_ERR_SUP, "Unknown CDF calculation/estimation method.");
374   }
375   *boxes = bxs;
376   if (listed) *listed = PETSC_TRUE;
377   PetscFunctionReturn(PETSC_SUCCESS);
378 }
379 
VecTaggerView_CDF(VecTagger tagger,PetscViewer viewer)380 static PetscErrorCode VecTaggerView_CDF(VecTagger tagger, PetscViewer viewer)
381 {
382   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
383   PetscBool      isascii;
384   PetscMPIInt    size;
385 
386   PetscFunctionBegin;
387   PetscCall(VecTaggerView_Simple(tagger, viewer));
388   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)tagger), &size));
389   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
390   if (size > 1 && isascii) {
391     PetscCall(PetscViewerASCIIPrintf(viewer, "CDF computation method: %s\n", VecTaggerCDFMethods[cuml->method]));
392     if (cuml->method == VECTAGGER_CDF_ITERATIVE) {
393       PetscCall(PetscViewerASCIIPushTab(viewer));
394       PetscCall(PetscViewerASCIIPrintf(viewer, "max its: %" PetscInt_FMT ", abs tol: %g, rel tol %g\n", cuml->maxit, (double)cuml->atol, (double)cuml->rtol));
395       PetscCall(PetscViewerASCIIPopTab(viewer));
396     }
397   }
398   PetscFunctionReturn(PETSC_SUCCESS);
399 }
400 
VecTaggerSetFromOptions_CDF(VecTagger tagger,PetscOptionItems PetscOptionsObject)401 static PetscErrorCode VecTaggerSetFromOptions_CDF(VecTagger tagger, PetscOptionItems PetscOptionsObject)
402 {
403   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
404   PetscInt       method;
405   PetscBool      set;
406 
407   PetscFunctionBegin;
408   PetscCall(VecTaggerSetFromOptions_Simple(tagger, PetscOptionsObject));
409   PetscOptionsHeadBegin(PetscOptionsObject, "VecTagger options for CDF boxes");
410   PetscCall(PetscOptionsEList("-vec_tagger_cdf_method", "Method for computing absolute boxes from CDF boxes", "VecTaggerCDFSetMethod()", VecTaggerCDFMethods, VECTAGGER_CDF_NUM_METHODS, VecTaggerCDFMethods[cuml->method], &method, &set));
411   if (set) cuml->method = (VecTaggerCDFMethod)method;
412   PetscCall(PetscOptionsInt("-vec_tagger_cdf_max_it", "Maximum iterations for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->maxit, &cuml->maxit, NULL));
413   PetscCall(PetscOptionsReal("-vec_tagger_cdf_rtol", "Maximum relative tolerance for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->rtol, &cuml->rtol, NULL));
414   PetscCall(PetscOptionsReal("-vec_tagger_cdf_atol", "Maximum absolute tolerance for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->atol, &cuml->atol, NULL));
415   PetscOptionsHeadEnd();
416   PetscFunctionReturn(PETSC_SUCCESS);
417 }
418 
419 /*@
420   VecTaggerCDFSetMethod - Set the method used to compute absolute boxes from CDF boxes
421 
422   Logically Collective
423 
424   Input Parameters:
425 + tagger - the `VecTagger` context
426 - method - the method
427 
428   Level: advanced
429 
430 .seealso: `Vec`, `VecTagger`, `VecTaggerCDFMethod`
431 @*/
VecTaggerCDFSetMethod(VecTagger tagger,VecTaggerCDFMethod method)432 PetscErrorCode VecTaggerCDFSetMethod(VecTagger tagger, VecTaggerCDFMethod method)
433 {
434   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
435 
436   PetscFunctionBegin;
437   PetscValidHeaderSpecific(tagger, VEC_TAGGER_CLASSID, 1);
438   PetscValidLogicalCollectiveInt(tagger, method, 2);
439   cuml->method = method;
440   PetscFunctionReturn(PETSC_SUCCESS);
441 }
442 
443 /*@
444   VecTaggerCDFGetMethod - Get the method used to compute absolute boxes from CDF boxes
445 
446   Logically Collective
447 
448   Input Parameter:
449 . tagger - the `VecTagger` context
450 
451   Output Parameter:
452 . method - the method
453 
454   Level: advanced
455 
456 .seealso: `Vec`, `VecTagger`, `VecTaggerCDFMethod`
457 @*/
VecTaggerCDFGetMethod(VecTagger tagger,VecTaggerCDFMethod * method)458 PetscErrorCode VecTaggerCDFGetMethod(VecTagger tagger, VecTaggerCDFMethod *method)
459 {
460   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
461 
462   PetscFunctionBegin;
463   PetscValidHeaderSpecific(tagger, VEC_TAGGER_CLASSID, 1);
464   PetscAssertPointer(method, 2);
465   *method = cuml->method;
466   PetscFunctionReturn(PETSC_SUCCESS);
467 }
468 
469 /*@C
470   VecTaggerCDFIterativeSetTolerances - Set the tolerances for iterative computation of absolute boxes from CDF boxes.
471 
472   Logically Collective
473 
474   Input Parameters:
475 + tagger - the `VecTagger` context
476 . maxit  - the maximum number of iterations: 0 indicates the absolute values will be estimated from an initial guess based only on the minimum, maximum, mean,
477           and standard deviation of the box endpoints.
478 . rtol   - the acceptable relative tolerance in the absolute values from the initial guess
479 - atol   - the acceptable absolute tolerance in the absolute values from the initial guess
480 
481   Level: advanced
482 
483 .seealso: `VecTagger`, `VecTaggerCDFSetMethod()`
484 @*/
VecTaggerCDFIterativeSetTolerances(VecTagger tagger,PetscInt maxit,PetscReal rtol,PetscReal atol)485 PetscErrorCode VecTaggerCDFIterativeSetTolerances(VecTagger tagger, PetscInt maxit, PetscReal rtol, PetscReal atol)
486 {
487   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
488 
489   PetscFunctionBegin;
490   PetscValidHeaderSpecific(tagger, VEC_TAGGER_CLASSID, 1);
491   PetscValidLogicalCollectiveInt(tagger, maxit, 2);
492   PetscValidLogicalCollectiveReal(tagger, rtol, 3);
493   PetscValidLogicalCollectiveReal(tagger, atol, 4);
494   cuml->maxit = maxit;
495   cuml->rtol  = rtol;
496   cuml->atol  = atol;
497   PetscFunctionReturn(PETSC_SUCCESS);
498 }
499 
500 /*@
501   VecTaggerCDFIterativeGetTolerances - Get the tolerances for iterative computation of absolute boxes from CDF boxes.
502 
503   Logically Collective
504 
505   Input Parameter:
506 . tagger - the `VecTagger` context
507 
508   Output Parameters:
509 + maxit - the maximum number of iterations: 0 indicates the absolute values will be estimated from an initial guess based only on the minimum, maximum,
510           mean, and standard deviation of the box endpoints.
511 . rtol  - the acceptable relative tolerance in the absolute values from the initial guess
512 - atol  - the acceptable absolute tolerance in the absolute values from the initial guess
513 
514   Level: advanced
515 
516 .seealso: `VecTagger`, `VecTaggerCDFSetMethod()`
517 @*/
VecTaggerCDFIterativeGetTolerances(VecTagger tagger,PetscInt * maxit,PetscReal * rtol,PetscReal * atol)518 PetscErrorCode VecTaggerCDFIterativeGetTolerances(VecTagger tagger, PetscInt *maxit, PetscReal *rtol, PetscReal *atol)
519 {
520   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
521 
522   PetscFunctionBegin;
523   PetscValidHeaderSpecific(tagger, VEC_TAGGER_CLASSID, 1);
524   if (maxit) *maxit = cuml->maxit;
525   if (rtol) *rtol = cuml->rtol;
526   if (atol) *atol = cuml->atol;
527   PetscFunctionReturn(PETSC_SUCCESS);
528 }
529 
530 /*@C
531   VecTaggerCDFSetBox - Set the cumulative box defining the values to be tagged by the tagger, where cumulative boxes are subsets of [0,1], where 0 indicates the smallest value present in the vector and 1 indicates the largest.
532 
533   Logically Collective
534 
535   Input Parameters:
536 + tagger - the `VecTagger` context
537 - box    - a blocksize array of `VecTaggerBox` boxes
538 
539   Level: advanced
540 
541 .seealso: `VecTagger`, `VecTaggerCDFGetBox()`, `VecTaggerBox`
542 @*/
VecTaggerCDFSetBox(VecTagger tagger,VecTaggerBox box[])543 PetscErrorCode VecTaggerCDFSetBox(VecTagger tagger, VecTaggerBox box[])
544 {
545   PetscFunctionBegin;
546   PetscCall(VecTaggerSetBox_Simple(tagger, box));
547   PetscFunctionReturn(PETSC_SUCCESS);
548 }
549 
550 /*@C
551   VecTaggerCDFGetBox - Get the cumulative box (multi-dimensional box) defining the values to be tagged by the tagger, where cumulative boxes
552   are subsets of [0,1], where 0 indicates the smallest value present in the vector and 1 indicates the largest.
553 
554   Logically Collective
555 
556   Input Parameter:
557 . tagger - the `VecTagger` context
558 
559   Output Parameter:
560 . box - a blocksize array of `VecTaggerBox` boxes
561 
562   Level: advanced
563 
564 .seealso: `VecTagger`, `VecTaggerCDFSetBox()`, `VecTaggerBox`
565 @*/
VecTaggerCDFGetBox(VecTagger tagger,const VecTaggerBox * box[])566 PetscErrorCode VecTaggerCDFGetBox(VecTagger tagger, const VecTaggerBox *box[])
567 {
568   PetscFunctionBegin;
569   PetscCall(VecTaggerGetBox_Simple(tagger, box));
570   PetscFunctionReturn(PETSC_SUCCESS);
571 }
572 
VecTaggerCreate_CDF(VecTagger tagger)573 PETSC_INTERN PetscErrorCode VecTaggerCreate_CDF(VecTagger tagger)
574 {
575   VecTagger_CDF *cuml;
576 
577   PetscFunctionBegin;
578   PetscCall(VecTaggerCreate_Simple(tagger));
579   PetscCall(PetscNew(&cuml));
580   PetscCall(PetscMemcpy(&cuml->smpl, tagger->data, sizeof(VecTagger_Simple)));
581   PetscCall(PetscFree(tagger->data));
582   tagger->data                = cuml;
583   tagger->ops->view           = VecTaggerView_CDF;
584   tagger->ops->setfromoptions = VecTaggerSetFromOptions_CDF;
585   tagger->ops->computeboxes   = VecTaggerComputeBoxes_CDF;
586   PetscFunctionReturn(PETSC_SUCCESS);
587 }
588