Actual source code: cdf.c

  1: #include <petsc/private/vecimpl.h>
  2: #include "../src/vec/vec/utils/tagger/impls/simple.h"

  4: const char *const VecTaggerCDFMethods[VECTAGGER_CDF_NUM_METHODS] = {"gather", "iterative"};

  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

 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;

 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;

 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: }

 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

 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;

 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));

 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;

 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));

 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: }

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;

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:   PetscCallMPI(MPI_Bcast((PetscScalar *)boxes, (PetscMPIInt)(2 * bs), MPIU_SCALAR, 0, PetscObjectComm((PetscObject)vec)));
118:   PetscCall(VecDestroy(&gVec));
119:   PetscFunctionReturn(PETSC_SUCCESS);
120: }

122: typedef struct _n_CDFStats {
123:   PetscReal min;
124:   PetscReal max;
125:   PetscReal moment[3];
126: } CDFStats;

128: static void MPIAPI VecTaggerCDFStatsReduce(void *a, void *b, int *len, MPI_Datatype *datatype)
129: {
130:   PetscInt  i, j, N = *len;
131:   CDFStats *A = (CDFStats *)a;
132:   CDFStats *B = (CDFStats *)b;

134:   for (i = 0; i < N; i++) {
135:     B[i].min = PetscMin(A[i].min, B[i].min);
136:     B[i].max = PetscMax(A[i].max, B[i].max);
137:     for (j = 0; j < 3; j++) B[i].moment[j] += A[i].moment[j];
138:   }
139: }

141: static PetscErrorCode CDFUtilInverseEstimate(const CDFStats *stats, PetscReal cdfTarget, PetscReal *absEst)
142: {
143:   PetscReal min, max;

145:   PetscFunctionBegin;
146:   min     = stats->min;
147:   max     = stats->max;
148:   *absEst = min + cdfTarget * (max - min);
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }

152: static PetscErrorCode VecTaggerComputeBox_CDF_SortedArray_Iterative(VecTagger tagger, MPI_Datatype statType, MPI_Op statReduce, const PetscReal *cArray, PetscInt m, const VecTaggerBoxReal *cdfBox, VecTaggerBoxReal *absBox)
153: {
154:   MPI_Comm       comm;
155:   VecTagger_CDF *cdf;
156:   PetscInt       maxit, i, j, k, l, M;
157:   PetscInt       bounds[2][2];
158:   PetscInt       offsets[2];
159:   PetscReal      intervalLen = cdfBox->max - cdfBox->min;
160:   PetscReal      rtol, atol;

162:   PetscFunctionBegin;
163:   comm  = PetscObjectComm((PetscObject)tagger);
164:   cdf   = (VecTagger_CDF *)tagger->data;
165:   maxit = cdf->maxit;
166:   rtol  = cdf->rtol;
167:   atol  = cdf->atol;
168:   /* local range of sorted values that can contain the sought radix */
169:   offsets[0]   = 0;
170:   offsets[1]   = 0;
171:   bounds[0][0] = 0;
172:   bounds[0][1] = m;
173:   bounds[1][0] = 0;
174:   bounds[1][1] = m;
175:   PetscCall(VecTaggerComputeBox_CDF_SortedArray(cArray, m, cdfBox, absBox)); /* compute a local estimate of the interval */
176:   {
177:     CDFStats stats[3];

179:     for (i = 0; i < 2; i++) { /* compute statistics of those local estimates */
180:       PetscReal val = i ? absBox->max : absBox->min;

182:       stats[i].min       = m ? val : PETSC_MAX_REAL;
183:       stats[i].max       = m ? val : PETSC_MIN_REAL;
184:       stats[i].moment[0] = m;
185:       stats[i].moment[1] = m * val;
186:       stats[i].moment[2] = m * val * val;
187:     }
188:     stats[2].min = PETSC_MAX_REAL;
189:     stats[2].max = PETSC_MAX_REAL;
190:     for (i = 0; i < 3; i++) stats[2].moment[i] = 0.;
191:     for (i = 0; i < m; i++) {
192:       PetscReal val = cArray[i];

194:       stats[2].min = PetscMin(stats[2].min, val);
195:       stats[2].max = PetscMax(stats[2].max, val);
196:       stats[2].moment[0]++;
197:       stats[2].moment[1] += val;
198:       stats[2].moment[2] += val * val;
199:     }
200:     /* reduce those statistics */
201:     PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, stats, 3, statType, statReduce, comm));
202:     M = (PetscInt)stats[2].moment[0];
203:     /* use those initial statistics to get the initial (globally agreed-upon) choices for the absolute box bounds */
204:     for (i = 0; i < 2; i++) PetscCall(CDFUtilInverseEstimate(&stats[i], i ? cdfBox->max : cdfBox->min, i ? &absBox->max : &absBox->min));
205:   }
206:   /* refine the estimates by computing how close they come to the desired box and refining */
207:   for (k = 0; k < maxit; k++) {
208:     PetscReal maxDiff = 0.;

210:     CDFStats stats[2][2];
211:     PetscInt newBounds[2][2][2];
212:     for (i = 0; i < 2; i++) {
213:       for (j = 0; j < 2; j++) {
214:         stats[i][j].min = PETSC_MAX_REAL;
215:         stats[i][j].max = PETSC_MIN_REAL;
216:         for (l = 0; l < 3; l++) stats[i][j].moment[l] = 0.;
217:         newBounds[i][j][0] = PetscMax(bounds[i][0], bounds[i][1]);
218:         newBounds[i][j][1] = PetscMin(bounds[i][0], bounds[i][1]);
219:       }
220:     }
221:     for (i = 0; i < 2; i++) {
222:       for (j = 0; j < bounds[i][1] - bounds[i][0]; j++) {
223:         PetscInt  thisInd = bounds[i][0] + j;
224:         PetscReal val     = cArray[thisInd];
225:         PetscInt  section;
226:         if (!i) {
227:           section = (val < absBox->min) ? 0 : 1;
228:         } else {
229:           section = (val <= absBox->max) ? 0 : 1;
230:         }
231:         stats[i][section].min = PetscMin(stats[i][section].min, val);
232:         stats[i][section].max = PetscMax(stats[i][section].max, val);
233:         stats[i][section].moment[0]++;
234:         stats[i][section].moment[1] += val;
235:         stats[i][section].moment[2] += val * val;
236:         newBounds[i][section][0] = PetscMin(newBounds[i][section][0], thisInd);
237:         newBounds[i][section][1] = PetscMax(newBounds[i][section][0], thisInd + 1);
238:       }
239:     }
240:     PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, stats, 4, statType, statReduce, comm));
241:     for (i = 0; i < 2; i++) {
242:       PetscInt  totalLessThan = offsets[i] + stats[i][0].moment[0];
243:       PetscReal cdfOfAbs      = (PetscReal)totalLessThan / (PetscReal)M;
244:       PetscReal diff;
245:       PetscInt  section;

247:       if (cdfOfAbs == (i ? cdfBox->max : cdfBox->min)) {
248:         offsets[i]   = totalLessThan;
249:         bounds[i][0] = bounds[i][1] = 0;
250:         continue;
251:       }
252:       if (cdfOfAbs > (i ? cdfBox->max : cdfBox->min)) { /* the correct absolute value lies in the lower section */
253:         section = 0;
254:       } else {
255:         section    = 1;
256:         offsets[i] = totalLessThan;
257:       }
258:       for (j = 0; j < 2; j++) bounds[i][j] = newBounds[i][section][j];
259:       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));
260:       diff    = PetscAbs(cdfOfAbs - (i ? cdfBox->max : cdfBox->min));
261:       maxDiff = PetscMax(maxDiff, diff);
262:     }
263:     if (!maxDiff) PetscFunctionReturn(PETSC_SUCCESS);
264:     if ((atol || rtol) && ((!atol) || (maxDiff <= atol)) && ((!rtol) || (maxDiff <= rtol * intervalLen))) break;
265:   }
266:   PetscFunctionReturn(PETSC_SUCCESS);
267: }

269: static PetscErrorCode VecTaggerComputeBoxes_CDF_Iterative(VecTagger tagger, Vec vec, PetscInt bs, VecTaggerBox *boxes)
270: {
271:   VecTagger_CDF    *cdf  = (VecTagger_CDF *)tagger->data;
272:   VecTagger_Simple *smpl = &cdf->smpl;
273:   Vec               vComp;
274:   PetscInt          i, N, M, n, m, rstart;
275: #if defined(PETSC_USE_COMPLEX)
276:   PetscReal *cReal, *cImag;
277: #endif
278:   MPI_Comm     comm;
279:   MPI_Datatype statType;
280:   MPI_Op       statReduce;

282:   PetscFunctionBegin;
283:   comm = PetscObjectComm((PetscObject)vec);
284:   PetscCall(VecGetSize(vec, &N));
285:   PetscCall(VecGetLocalSize(vec, &n));
286:   M = N / bs;
287:   m = n / bs;
288:   PetscCall(VecCreateMPI(comm, m, M, &vComp));
289:   PetscCall(VecSetUp(vComp));
290:   PetscCall(VecGetOwnershipRange(vComp, &rstart, NULL));
291: #if defined(PETSC_USE_COMPLEX)
292:   PetscCall(PetscMalloc2(m, &cReal, m, &cImag));
293: #endif
294:   PetscCallMPI(MPI_Type_contiguous(5, MPIU_REAL, &statType));
295:   PetscCallMPI(MPI_Type_commit(&statType));
296:   PetscCallMPI(MPI_Op_create(VecTaggerCDFStatsReduce, 1, &statReduce));
297:   for (i = 0; i < bs; i++) {
298:     IS           isStride;
299:     VecScatter   vScat;
300:     PetscScalar *cArray;

302:     PetscCall(ISCreateStride(comm, m, bs * rstart + i, bs, &isStride));
303:     PetscCall(VecScatterCreate(vec, isStride, vComp, NULL, &vScat));
304:     PetscCall(VecScatterBegin(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
305:     PetscCall(VecScatterEnd(vScat, vec, vComp, INSERT_VALUES, SCATTER_FORWARD));
306:     PetscCall(VecScatterDestroy(&vScat));
307:     PetscCall(ISDestroy(&isStride));

309:     PetscCall(VecGetArray(vComp, &cArray));
310: #if !defined(PETSC_USE_COMPLEX)
311:     PetscCall(PetscSortReal(m, cArray));
312:     PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cArray, m, &smpl->box[i], &boxes[i]));
313: #else
314:     {
315:       PetscInt         j;
316:       VecTaggerBoxReal realBxs, imagBxs;
317:       VecTaggerBoxReal realBoxes, imagBoxes;

319:       for (j = 0; j < m; j++) {
320:         cReal[j] = PetscRealPart(cArray[j]);
321:         cImag[j] = PetscImaginaryPart(cArray[j]);
322:       }
323:       PetscCall(PetscSortReal(m, cReal));
324:       PetscCall(PetscSortReal(m, cImag));

326:       realBxs.min = PetscRealPart(smpl->box[i].min);
327:       realBxs.max = PetscRealPart(smpl->box[i].max);
328:       imagBxs.min = PetscImaginaryPart(smpl->box[i].min);
329:       imagBxs.max = PetscImaginaryPart(smpl->box[i].max);
330:       PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cReal, m, &realBxs, &realBoxes));
331:       PetscCall(VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger, statType, statReduce, cImag, m, &imagBxs, &imagBoxes));
332:       boxes[i].min = PetscCMPLX(realBoxes.min, imagBoxes.min);
333:       boxes[i].max = PetscCMPLX(realBoxes.max, imagBoxes.max);
334:     }
335: #endif
336:     PetscCall(VecRestoreArray(vComp, &cArray));
337:   }
338:   PetscCallMPI(MPI_Op_free(&statReduce));
339:   PetscCallMPI(MPI_Type_free(&statType));
340: #if defined(PETSC_USE_COMPLEX)
341:   PetscCall(PetscFree2(cReal, cImag));
342: #endif
343:   PetscCall(VecDestroy(&vComp));
344:   PetscFunctionReturn(PETSC_SUCCESS);
345: }

347: static PetscErrorCode VecTaggerComputeBoxes_CDF(VecTagger tagger, Vec vec, PetscInt *numBoxes, VecTaggerBox **boxes, PetscBool *listed)
348: {
349:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
350:   PetscMPIInt    size;
351:   PetscInt       bs;
352:   VecTaggerBox  *bxs;

354:   PetscFunctionBegin;
355:   PetscCall(VecTaggerGetBlockSize(tagger, &bs));
356:   *numBoxes = 1;
357:   PetscCall(PetscMalloc1(bs, &bxs));
358:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)tagger), &size));
359:   if (size == 1) {
360:     PetscCall(VecTaggerComputeBoxes_CDF_Serial(tagger, vec, bs, bxs));
361:     *boxes = bxs;
362:     PetscFunctionReturn(PETSC_SUCCESS);
363:   }
364:   switch (cuml->method) {
365:   case VECTAGGER_CDF_GATHER:
366:     PetscCall(VecTaggerComputeBoxes_CDF_Gather(tagger, vec, bs, bxs));
367:     break;
368:   case VECTAGGER_CDF_ITERATIVE:
369:     PetscCall(VecTaggerComputeBoxes_CDF_Iterative(tagger, vec, bs, bxs));
370:     break;
371:   default:
372:     SETERRQ(PetscObjectComm((PetscObject)tagger), PETSC_ERR_SUP, "Unknown CDF calculation/estimation method.");
373:   }
374:   *boxes = bxs;
375:   if (listed) *listed = PETSC_TRUE;
376:   PetscFunctionReturn(PETSC_SUCCESS);
377: }

379: static PetscErrorCode VecTaggerView_CDF(VecTagger tagger, PetscViewer viewer)
380: {
381:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
382:   PetscBool      iascii;
383:   PetscMPIInt    size;

385:   PetscFunctionBegin;
386:   PetscCall(VecTaggerView_Simple(tagger, viewer));
387:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)tagger), &size));
388:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
389:   if (size > 1 && iascii) {
390:     PetscCall(PetscViewerASCIIPrintf(viewer, "CDF computation method: %s\n", VecTaggerCDFMethods[cuml->method]));
391:     if (cuml->method == VECTAGGER_CDF_ITERATIVE) {
392:       PetscCall(PetscViewerASCIIPushTab(viewer));
393:       PetscCall(PetscViewerASCIIPrintf(viewer, "max its: %" PetscInt_FMT ", abs tol: %g, rel tol %g\n", cuml->maxit, (double)cuml->atol, (double)cuml->rtol));
394:       PetscCall(PetscViewerASCIIPopTab(viewer));
395:     }
396:   }
397:   PetscFunctionReturn(PETSC_SUCCESS);
398: }

400: static PetscErrorCode VecTaggerSetFromOptions_CDF(VecTagger tagger, PetscOptionItems *PetscOptionsObject)
401: {
402:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
403:   PetscInt       method;
404:   PetscBool      set;

406:   PetscFunctionBegin;
407:   PetscCall(VecTaggerSetFromOptions_Simple(tagger, PetscOptionsObject));
408:   PetscOptionsHeadBegin(PetscOptionsObject, "VecTagger options for CDF boxes");
409:   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));
410:   if (set) cuml->method = (VecTaggerCDFMethod)method;
411:   PetscCall(PetscOptionsInt("-vec_tagger_cdf_max_it", "Maximum iterations for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->maxit, &cuml->maxit, NULL));
412:   PetscCall(PetscOptionsReal("-vec_tagger_cdf_rtol", "Maximum relative tolerance for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->rtol, &cuml->rtol, NULL));
413:   PetscCall(PetscOptionsReal("-vec_tagger_cdf_atol", "Maximum absolute tolerance for iterative computation of absolute boxes from CDF boxes", "VecTaggerCDFIterativeSetTolerances()", cuml->atol, &cuml->atol, NULL));
414:   PetscOptionsHeadEnd();
415:   PetscFunctionReturn(PETSC_SUCCESS);
416: }

418: /*@
419:   VecTaggerCDFSetMethod - Set the method used to compute absolute boxes from CDF boxes

421:   Logically Collective

423:   Input Parameters:
424: + tagger - the `VecTagger` context
425: - method - the method

427:   Level: advanced

429: .seealso: `Vec`, `VecTagger`, `VecTaggerCDFMethod`
430: @*/
431: PetscErrorCode VecTaggerCDFSetMethod(VecTagger tagger, VecTaggerCDFMethod method)
432: {
433:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;

435:   PetscFunctionBegin;
438:   cuml->method = method;
439:   PetscFunctionReturn(PETSC_SUCCESS);
440: }

442: /*@
443:   VecTaggerCDFGetMethod - Get the method used to compute absolute boxes from CDF boxes

445:   Logically Collective

447:   Input Parameter:
448: . tagger - the `VecTagger` context

450:   Output Parameter:
451: . method - the method

453:   Level: advanced

455: .seealso: `Vec`, `VecTagger`, `VecTaggerCDFMethod`
456: @*/
457: PetscErrorCode VecTaggerCDFGetMethod(VecTagger tagger, VecTaggerCDFMethod *method)
458: {
459:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;

461:   PetscFunctionBegin;
463:   PetscAssertPointer(method, 2);
464:   *method = cuml->method;
465:   PetscFunctionReturn(PETSC_SUCCESS);
466: }

468: /*@C
469:   VecTaggerCDFIterativeSetTolerances - Set the tolerances for iterative computation of absolute boxes from CDF boxes.

471:   Logically Collective

473:   Input Parameters:
474: + tagger - the `VecTagger` context
475: . 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,
476:           and standard deviation of the box endpoints.
477: . rtol   - the acceptable relative tolerance in the absolute values from the initial guess
478: - atol   - the acceptable absolute tolerance in the absolute values from the initial guess

480:   Level: advanced

482: .seealso: `VecTagger`, `VecTaggerCDFSetMethod()`
483: @*/
484: PetscErrorCode VecTaggerCDFIterativeSetTolerances(VecTagger tagger, PetscInt maxit, PetscReal rtol, PetscReal atol)
485: {
486:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;

488:   PetscFunctionBegin;
493:   cuml->maxit = maxit;
494:   cuml->rtol  = rtol;
495:   cuml->atol  = atol;
496:   PetscFunctionReturn(PETSC_SUCCESS);
497: }

499: /*@
500:   VecTaggerCDFIterativeGetTolerances - Get the tolerances for iterative computation of absolute boxes from CDF boxes.

502:   Logically Collective

504:   Input Parameter:
505: . tagger - the `VecTagger` context

507:   Output Parameters:
508: + maxit - the maximum number of iterations: 0 indicates the absolute values will be estimated from an initial guess based only on the minimum, maximum,
509:           mean, and standard deviation of the box endpoints.
510: . rtol  - the acceptable relative tolerance in the absolute values from the initial guess
511: - atol  - the acceptable absolute tolerance in the absolute values from the initial guess

513:   Level: advanced

515: .seealso: `VecTagger`, `VecTaggerCDFSetMethod()`
516: @*/
517: PetscErrorCode VecTaggerCDFIterativeGetTolerances(VecTagger tagger, PetscInt *maxit, PetscReal *rtol, PetscReal *atol)
518: {
519:   VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;

521:   PetscFunctionBegin;
523:   if (maxit) *maxit = cuml->maxit;
524:   if (rtol) *rtol = cuml->rtol;
525:   if (atol) *atol = cuml->atol;
526:   PetscFunctionReturn(PETSC_SUCCESS);
527: }

529: /*@C
530:   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:   Logically Collective

534:   Input Parameters:
535: + tagger - the `VecTagger` context
536: - box    - a blocksize array of `VecTaggerBox` boxes

538:   Level: advanced

540: .seealso: `VecTagger`, `VecTaggerCDFGetBox()`, `VecTaggerBox`
541: @*/
542: PetscErrorCode VecTaggerCDFSetBox(VecTagger tagger, VecTaggerBox box[])
543: {
544:   PetscFunctionBegin;
545:   PetscCall(VecTaggerSetBox_Simple(tagger, box));
546:   PetscFunctionReturn(PETSC_SUCCESS);
547: }

549: /*@C
550:   VecTaggerCDFGetBox - Get the cumulative box (multi-dimensional box) defining the values to be tagged by the tagger, where cumulative boxes
551:   are subsets of [0,1], where 0 indicates the smallest value present in the vector and 1 indicates the largest.

553:   Logically Collective

555:   Input Parameter:
556: . tagger - the `VecTagger` context

558:   Output Parameter:
559: . box - a blocksize array of `VecTaggerBox` boxes

561:   Level: advanced

563: .seealso: `VecTagger`, `VecTaggerCDFSetBox()`, `VecTaggerBox`
564: @*/
565: PetscErrorCode VecTaggerCDFGetBox(VecTagger tagger, const VecTaggerBox *box[])
566: {
567:   PetscFunctionBegin;
568:   PetscCall(VecTaggerGetBox_Simple(tagger, box));
569:   PetscFunctionReturn(PETSC_SUCCESS);
570: }

572: PETSC_INTERN PetscErrorCode VecTaggerCreate_CDF(VecTagger tagger)
573: {
574:   VecTagger_CDF *cuml;

576:   PetscFunctionBegin;
577:   PetscCall(VecTaggerCreate_Simple(tagger));
578:   PetscCall(PetscNew(&cuml));
579:   PetscCall(PetscMemcpy(&cuml->smpl, tagger->data, sizeof(VecTagger_Simple)));
580:   PetscCall(PetscFree(tagger->data));
581:   tagger->data                = cuml;
582:   tagger->ops->view           = VecTaggerView_CDF;
583:   tagger->ops->setfromoptions = VecTaggerSetFromOptions_CDF;
584:   tagger->ops->computeboxes   = VecTaggerComputeBoxes_CDF;
585:   PetscFunctionReturn(PETSC_SUCCESS);
586: }