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