Actual source code: cdf.c
petsc-3.9.4 2018-09-11
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];
133: } CDFStats;
135: static void VecTaggerCDFStatsReduce(void *a, void *b, int * len, MPI_Datatype *datatype)
136: {
137: PetscInt i, j, N = *len;
138: CDFStats *A = (CDFStats *) a;
139: CDFStats *B = (CDFStats *) b;
141: for (i = 0; i < N; i++) {
142: B[i].min = PetscMin(A[i].min,B[i].min);
143: B[i].max = PetscMax(A[i].max,B[i].max);
144: for (j = 0; j < 3; j++) {
145: B[i].moment[j] += A[i].moment[j];
146: }
147: }
148: }
150: static PetscErrorCode CDFUtilInverseEstimate(const CDFStats *stats,PetscReal cdfTarget,PetscReal *absEst)
151: {
152: PetscReal min, max;
155: min = stats->min;
156: max = stats->max;
157: *absEst = min + cdfTarget * (max - min);
158: return(0);
159: }
161: static PetscErrorCode VecTaggerComputeBox_CDF_SortedArray_Iterative(VecTagger tagger, MPI_Datatype statType, MPI_Op statReduce, const PetscReal *cArray, PetscInt m, const VecTaggerBoxReal *cdfBox, VecTaggerBoxReal *absBox)
162: {
163: MPI_Comm comm;
164: VecTagger_CDF *cdf;
165: PetscInt maxit, i, j, k, l, M;
166: PetscInt bounds[2][2];
167: PetscInt offsets[2];
168: PetscReal intervalLen = cdfBox->max - cdfBox->min;
169: PetscReal rtol, atol;
173: comm = PetscObjectComm((PetscObject)tagger);
174: cdf = (VecTagger_CDF *) tagger->data;
175: maxit = cdf->maxit;
176: rtol = cdf->rtol;
177: atol = cdf->atol;
178: /* local range of sorted values that can contain the sought radix */
179: offsets[0] = 0;
180: offsets[1] = 0;
181: bounds[0][0] = 0; bounds[0][1] = m;
182: bounds[1][0] = 0; bounds[1][1] = m;
183: VecTaggerComputeBox_CDF_SortedArray(cArray,m,cdfBox,absBox); /* compute a local estimate of the interval */
184: {
185: CDFStats stats[3];
187: for (i = 0; i < 2; i++) { /* compute statistics of those local estimates */
188: PetscReal val = i ? absBox->max : absBox->min;
190: stats[i].min = m ? val : PETSC_MAX_REAL;
191: stats[i].max = m ? val : PETSC_MIN_REAL;
192: stats[i].moment[0] = m;
193: stats[i].moment[1] = m * val;
194: stats[i].moment[2] = m * val * val;
195: }
196: stats[2].min = PETSC_MAX_REAL;
197: stats[2].max = PETSC_MAX_REAL;
198: for (i = 0; i < 3; i++) {
199: stats[2].moment[i] = 0.;
200: }
201: for (i = 0; i < m; i++) {
202: PetscReal val = cArray[i];
204: stats[2].min = PetscMin(stats[2].min, val);
205: stats[2].max = PetscMax(stats[2].max, val);
206: stats[2].moment[0] ++;
207: stats[2].moment[1] += val;
208: stats[2].moment[2] += val * val;
209: }
210: /* reduce those statistics */
211: MPI_Allreduce(MPI_IN_PLACE,stats,3,statType,statReduce,comm);
212: M = (PetscInt) stats[2].moment[0];
213: /* use those initial statistics to get the initial (globally agreed-upon) choices for the absolute box bounds */
214: for (i = 0; i < 2; i++) {
215: CDFUtilInverseEstimate(&stats[i],i ? cdfBox->max : cdfBox->min,(i ? &absBox->max : &absBox->min));
216: }
217: }
218: /* refine the estimates by computing how close they come to the desired box and refining */
219: for (k = 0; k < maxit; k++) {
220: PetscReal maxDiff = 0.;
222: CDFStats stats[2][2];
223: PetscInt newBounds[2][2][2];
224: for (i = 0; i < 2; i++) {
225: for (j = 0; j < 2; j++) {
226: stats[i][j].min = PETSC_MAX_REAL;
227: stats[i][j].max = PETSC_MIN_REAL;
228: for (l = 0; l < 3; l++) {
229: stats[i][j].moment[l] = 0.;
230: }
231: newBounds[i][j][0] = PetscMax(bounds[i][0],bounds[i][1]);
232: newBounds[i][j][1] = PetscMin(bounds[i][0],bounds[i][1]);
233: }
234: }
235: for (i = 0; i < 2; i++) {
236: for (j = 0; j < bounds[i][1] - bounds[i][0]; j++) {
237: PetscInt thisInd = bounds[i][0] + j;
238: PetscReal val = cArray[thisInd];
239: PetscInt section;
240: if (!i) {
241: section = (val < absBox->min) ? 0 : 1;
242: } else {
243: section = (val <= absBox->max) ? 0 : 1;
244: }
245: stats[i][section].min = PetscMin(stats[i][section].min,val);
246: stats[i][section].max = PetscMax(stats[i][section].max,val);
247: stats[i][section].moment[0] ++;
248: stats[i][section].moment[1] += val;
249: stats[i][section].moment[2] += val * val;
250: newBounds[i][section][0] = PetscMin(newBounds[i][section][0],thisInd);
251: newBounds[i][section][1] = PetscMax(newBounds[i][section][0],thisInd + 1);
252: }
253: }
254: MPI_Allreduce(MPI_IN_PLACE, stats, 4, statType, statReduce, comm);
255: for (i = 0; i < 2; i++) {
256: PetscInt totalLessThan = offsets[i] + stats[i][0].moment[0];
257: PetscReal cdfOfAbs = (PetscReal) totalLessThan / (PetscReal) M;
258: PetscReal diff;
259: PetscInt section;
261: if (cdfOfAbs == (i ? cdfBox->max : cdfBox->min)) {
262: offsets[i] = totalLessThan;
263: bounds[i][0] = bounds[i][1] = 0;
264: continue;
265: }
266: if (cdfOfAbs > (i ? cdfBox->max : cdfBox->min)) { /* the correct absolute value lies in the lower section */
267: section = 0;
268: } else {
269: section = 1;
270: offsets[i] = totalLessThan;
271: }
272: for (j = 0; j < 2; j++) {
273: bounds[i][j] = newBounds[i][section][j];
274: }
275: CDFUtilInverseEstimate(&stats[i][section],((i ? cdfBox->max : cdfBox->min) - ((PetscReal) offsets[i] / (PetscReal) M))/stats[i][section].moment[0],(i ? &absBox->max : &absBox->min));
276: diff = PetscAbs(cdfOfAbs - (i ? cdfBox->max : cdfBox->min));
277: maxDiff = PetscMax(maxDiff,diff);
278: }
279: if (!maxDiff) return(0);
280: if ((atol || rtol) && ((!atol) || (maxDiff <= atol)) && ((!rtol) || (maxDiff <= rtol * intervalLen))) {
281: break;
282: }
283: }
284: return(0);
285: }
287: static PetscErrorCode VecTaggerComputeBoxes_CDF_Iterative(VecTagger tagger,Vec vec,PetscInt bs,VecTaggerBox *boxes)
288: {
289: VecTagger_CDF *cdf = (VecTagger_CDF *) tagger->data;
290: VecTagger_Simple *smpl = &(cdf->smpl);
291: Vec vComp;
292: PetscInt i, N, M, n, m, rstart;
293: #if defined (PETSC_USE_COMPLEX)
294: PetscReal *cReal, *cImag;
295: #endif
296: MPI_Comm comm;
297: MPI_Datatype statType;
298: MPI_Op statReduce;
299: PetscErrorCode ierr;
302: comm = PetscObjectComm((PetscObject)vec);
303: VecGetSize(vec,&N);
304: VecGetLocalSize(vec,&n);
305: M = N/bs;
306: m = n/bs;
307: VecCreateMPI(comm,m,M,&vComp);
308: VecSetUp(vComp);
309: VecGetOwnershipRange(vComp,&rstart,NULL);
310: #if defined (PETSC_USE_COMPLEX)
311: PetscMalloc2(m,&cReal,m,&cImag);
312: #endif
313: MPI_Type_contiguous(5,MPIU_REAL,&statType);
314: MPI_Type_commit(&statType);
315: MPI_Op_create(VecTaggerCDFStatsReduce,1,&statReduce);
316: for (i = 0; i < bs; i++) {
317: IS isStride;
318: VecScatter vScat;
319: PetscScalar *cArray;
321: ISCreateStride(comm,m,bs * rstart + i,bs,&isStride);
322: VecScatterCreate(vec,isStride,vComp,NULL,&vScat);
323: VecScatterBegin(vScat,vec,vComp,INSERT_VALUES,SCATTER_FORWARD);
324: VecScatterEnd(vScat,vec,vComp,INSERT_VALUES,SCATTER_FORWARD);
325: VecScatterDestroy(&vScat);
326: ISDestroy(&isStride);
328: VecGetArray(vComp,&cArray);
329: #if !defined(PETSC_USE_COMPLEX)
330: PetscSortReal(m,cArray);
331: VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger,statType,statReduce,cArray,m,&smpl->box[i],&boxes[i]);
332: #else
333: {
334: PetscInt j;
335: VecTaggerBoxReal realBxs, imagBxs;
336: VecTaggerBoxReal realBoxes, imagBoxes;
338: for (j = 0; j < m; j++) {
339: cReal[j] = PetscRealPart(cArray[j]);
340: cImag[j] = PetscImaginaryPart(cArray[j]);
341: }
342: PetscSortReal(m,cReal);
343: PetscSortReal(m,cImag);
345: realBxs.min = PetscRealPart(smpl->box[i].min);
346: realBxs.max = PetscRealPart(smpl->box[i].max);
347: imagBxs.min = PetscImaginaryPart(smpl->box[i].min);
348: imagBxs.max = PetscImaginaryPart(smpl->box[i].max);
349: VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger,statType,statReduce,cReal,m,&realBxs,&realBoxes);
350: VecTaggerComputeBox_CDF_SortedArray_Iterative(tagger,statType,statReduce,cImag,m,&imagBxs,&imagBoxes);
351: boxes[i].min = PetscCMPLX(realBoxes.min,imagBoxes.min);
352: boxes[i].max = PetscCMPLX(realBoxes.max,imagBoxes.max);
353: }
354: #endif
355: VecRestoreArray(vComp,&cArray);
356: }
357: MPI_Op_free(&statReduce);
358: MPI_Type_free(&statType);
359: #if defined(PETSC_USE_COMPLEX)
360: PetscFree2(cReal,cImag);
361: #endif
362: VecDestroy(&vComp);
363: return(0);
364: }
366: static PetscErrorCode VecTaggerComputeBoxes_CDF(VecTagger tagger,Vec vec,PetscInt *numBoxes,VecTaggerBox **boxes)
367: {
368: VecTagger_CDF *cuml = (VecTagger_CDF *)tagger->data;
369: PetscMPIInt size;
370: PetscInt bs;
371: VecTaggerBox *bxs;
372: PetscErrorCode ierr;
375: VecTaggerGetBlockSize(tagger,&bs);
376: *numBoxes = 1;
377: PetscMalloc1(bs,&bxs);
378: MPI_Comm_size(PetscObjectComm((PetscObject)tagger),&size);
379: if (size == 1) {
380: VecTaggerComputeBoxes_CDF_Serial(tagger,vec,bs,bxs);
381: *boxes = bxs;
382: return(0);
383: }
384: switch (cuml->method) {
385: case VECTAGGER_CDF_GATHER:
386: VecTaggerComputeBoxes_CDF_Gather(tagger,vec,bs,bxs);
387: break;
388: case VECTAGGER_CDF_ITERATIVE:
389: VecTaggerComputeBoxes_CDF_Iterative(tagger,vec,bs,bxs);
390: break;
391: default:
392: SETERRQ(PetscObjectComm((PetscObject)tagger),PETSC_ERR_SUP,"Unknown CDF calculation/estimation method.");
393: }
394: *boxes = bxs;
395: return(0);
396: }
398: static PetscErrorCode VecTaggerView_CDF(VecTagger tagger,PetscViewer viewer)
399: {
400: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
401: PetscBool iascii;
402: PetscMPIInt size;
403: PetscErrorCode ierr;
406: VecTaggerView_Simple(tagger,viewer);
407: MPI_Comm_size(PetscObjectComm((PetscObject)tagger),&size);
408: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
409: if (size > 1 && iascii) {
410: PetscViewerASCIIPrintf(viewer,"CDF computation method: %s\n",VecTaggerCDFMethods[cuml->method]);
411: if (cuml->method == VECTAGGER_CDF_ITERATIVE) {
412: PetscViewerASCIIPushTab(viewer);
413: PetscViewerASCIIPrintf(viewer,"max its: %D, abs tol: %g, rel tol %g\n",cuml->maxit,(double) cuml->atol,(double) cuml->rtol);
414: PetscViewerASCIIPopTab(viewer);
415: }
416: }
417: return(0);
418: }
420: static PetscErrorCode VecTaggerSetFromOptions_CDF(PetscOptionItems *PetscOptionsObject,VecTagger tagger)
421: {
422: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
423: PetscInt method;
424: PetscBool set;
428: VecTaggerSetFromOptions_Simple(PetscOptionsObject,tagger);
429: PetscOptionsHead(PetscOptionsObject,"VecTagger options for CDF boxes");
430: PetscOptionsEList("-vec_tagger_cdf_method","Method for computing absolute boxes from CDF boxes","VecTaggerCDFSetMethod()",VecTaggerCDFMethods,VECTAGGER_CDF_NUM_METHODS,VecTaggerCDFMethods[cuml->method],&method,&set);
431: if (set) cuml->method = (VecTaggerCDFMethod) method;
432: PetscOptionsInt("-vec_tagger_cdf_max_it","Maximum iterations for iterative computation of absolute boxes from CDF boxes","VecTaggerCDFIterativeSetTolerances()",cuml->maxit,&cuml->maxit,NULL);
433: PetscOptionsReal("-vec_tagger_cdf_rtol","Maximum relative tolerance for iterative computation of absolute boxes from CDF boxes","VecTaggerCDFIterativeSetTolerances()",cuml->rtol,&cuml->rtol,NULL);
434: PetscOptionsReal("-vec_tagger_cdf_atol","Maximum absolute tolerance for iterative computation of absolute boxes from CDF boxes","VecTaggerCDFIterativeSetTolerances()",cuml->atol,&cuml->atol,NULL);
435: PetscOptionsTail();
437: return(0);
438: }
440: /*@C
441: VecTaggerCDFSetMethod - Set the method used to compute absolute boxes from CDF boxes
443: Logically Collective on VecTagger
445: Level: advanced
447: Input Parameters:
448: + tagger - the VecTagger context
449: - method - the method
451: .seealso VecTaggerCDFMethod
452: @*/
453: PetscErrorCode VecTaggerCDFSetMethod(VecTagger tagger, VecTaggerCDFMethod method)
454: {
455: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
460: cuml->method = method;
461: return(0);
462: }
464: /*@C
465: VecTaggerCDFGetMethod - Get the method used to compute absolute boxes from CDF boxes
467: Logically Collective on VecTagger
469: Level: advanced
471: Input Parameters:
472: . tagger - the VecTagger context
474: Output Parameters:
475: . method - the method
477: .seealso VecTaggerCDFMethod
478: @*/
479: PetscErrorCode VecTaggerCDFGetMethod(VecTagger tagger, VecTaggerCDFMethod *method)
480: {
481: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
486: *method = cuml->method;
487: return(0);
488: }
490: /*@C
491: VecTaggerCDFIterativeSetTolerances - Set the tolerances for iterative computation of absolute boxes from CDF boxes.
493: Logically Collective on VecTagger
495: Input Parameters:
496: + tagger - the VecTagger context
497: . 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.
498: . rtol - the acceptable relative tolerance in the absolute values from the initial guess
499: - atol - the acceptable absolute tolerance in the absolute values from the initial guess
501: Level: advanced
503: .seealso: VecTaggerCDFSetMethod()
504: @*/
505: PetscErrorCode VecTaggerCDFIterativeSetTolerances(VecTagger tagger, PetscInt maxit, PetscReal rtol, PetscReal atol)
506: {
507: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
514: cuml->maxit = maxit;
515: cuml->rtol = rtol;
516: cuml->atol = atol;
517: return(0);
518: }
520: /*@C
521: VecTaggerCDFIterativeGetTolerances - Get the tolerances for iterative computation of absolute boxes from CDF boxes.
523: Logically Collective on VecTagger
525: Input Parameters:
526: . tagger - the VecTagger context
528: Output Parameters:
529: + 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.
530: . rtol - the acceptable relative tolerance in the absolute values from the initial guess
531: - atol - the acceptable absolute tolerance in the absolute values from the initial guess
533: Level: advanced
535: .seealso: VecTaggerCDFSetMethod()
536: @*/
537: PetscErrorCode VecTaggerCDFIterativeGetTolerances(VecTagger tagger, PetscInt *maxit, PetscReal *rtol, PetscReal *atol)
538: {
539: VecTagger_CDF *cuml = (VecTagger_CDF *) tagger->data;
543: if (maxit) *maxit = cuml->maxit;
544: if (rtol) *rtol = cuml->rtol;
545: if (atol) *atol = cuml->atol;
546: return(0);
547: }
549: /*@C
550: 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.
552: Logically Collective
554: Input Arguments:
555: + tagger - the VecTagger context
556: - boxes - a blocksize array of VecTaggerBox boxes
558: Level: advanced
560: .seealso: VecTaggerCDFGetBox()
561: @*/
562: PetscErrorCode VecTaggerCDFSetBox(VecTagger tagger,VecTaggerBox *box)
563: {
567: VecTaggerSetBox_Simple(tagger,box);
568: return(0);
569: }
571: /*@C
572: 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.
574: Logically Collective
576: Input Arguments:
577: . tagger - the VecTagger context
579: Output Arguments:
580: . boxes - a blocksize array of VecTaggerBox boxes
582: Level: advanced
584: .seealso: VecTaggerCDFSetBox()
585: @*/
586: PetscErrorCode VecTaggerCDFGetBox(VecTagger tagger,const VecTaggerBox **box)
587: {
591: VecTaggerGetBox_Simple(tagger,box);
592: return(0);
593: }
595: PETSC_INTERN PetscErrorCode VecTaggerCreate_CDF(VecTagger tagger)
596: {
597: VecTagger_CDF *cuml;
598: PetscErrorCode ierr;
601: VecTaggerCreate_Simple(tagger);
602: PetscNewLog(tagger,&cuml);
603: PetscMemcpy(&(cuml->smpl),tagger->data,sizeof(VecTagger_Simple));
604: PetscFree(tagger->data);
605: tagger->data = cuml;
606: tagger->ops->view = VecTaggerView_CDF;
607: tagger->ops->setfromoptions = VecTaggerSetFromOptions_CDF;
608: tagger->ops->computeboxes = VecTaggerComputeBoxes_CDF;
609: return(0);
610: }