Actual source code: ex52.c
petsc-3.14.6 2021-03-30
1: static char help[] = "A benchmark for testing PetscSortInt(), PetscSortIntSemiOrdered(), PetscSortIntWithArrayPair(), PetscIntSortSemiOrderedWithArray(), and PetscSortIntWithArray()\n\
2: The array is filled with random numbers, but one can control average duplicates for each unique integer with the -d option.\n\
3: Usage:\n\
4: mpirun -n 1 ./ex52 -n <length of the array to sort>, default=100 \n\
5: -r <repeat times for each sort>, default=10 \n\
6: -d <average duplicates for each unique integer>, default=1, i.e., no duplicates \n\n";
8: #include <petscsys.h>
9: #include <petsctime.h>
10: #include <petscviewer.h>
11: #include <petscvec.h>
12: int main(int argc,char **argv)
13: {
15: PetscInt i,l,n=100,r=10,d=1,vsize=1;
16: PetscInt *X,*X1,*XR,*XSO,*W,*Y,*Z,*XP,*X1P;
17: PetscReal val,norm1,nreal;
18: PetscRandom rdm,rdm2;
19: PetscLogDouble time, time1, time2;
20: PetscMPIInt size;
21: PetscViewer vwr;
22: Vec x;
23: unsigned long seedr, seedo;
24: PetscBool order=PETSC_FALSE;
26: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
27: MPI_Comm_size(PETSC_COMM_WORLD,&size);
28: if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_ARG_OUTOFRANGE,"This is a uniprocessor example only!");
30: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
31: PetscOptionsGetInt(NULL,NULL,"-r",&r,NULL);
32: PetscOptionsGetInt(NULL,NULL,"-d",&d,NULL);
33: PetscOptionsGetInt(NULL,NULL,"-vsize",&vsize,NULL);
34: PetscOptionsGetBool(NULL,NULL,"-order",NULL,&order);
35: PetscOptionsGetViewer(PETSC_COMM_WORLD,NULL,NULL,"-array_view",&vwr,NULL,NULL);
36: if (n<1 || r<1 || d<1 || d>n) SETERRQ3(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Wrong input n=%D,r=%D,d=%D. They must be >=1 and n>=d\n",n,r,d);
38: PetscCalloc6(n,&X,n,&X1,n,&XR,n,&XSO,n,&Y,n,&Z);
39: PetscRandomCreate(PETSC_COMM_SELF,&rdm);
40: PetscRandomSetFromOptions(rdm);
41: PetscRandomGetSeed(rdm, &seedr);
43: for (i=0; i<n; ++i) {
44: PetscRandomGetValueReal(rdm,&val);
45: XR[i] = val*PETSC_MAX_INT;
46: if (d > 1) XR[i] = XR[i] % (n/d);
47: XSO[i] = i;
48: if (d > 1) XSO[i] = XSO[i] % (n/d);
49: }
51: nreal = (PetscReal) n;
52: PetscRandomCreate(PETSC_COMM_SELF,&rdm2);
53: PetscRandomGetSeed(rdm, &seedo);
54: PetscRandomSetInterval(rdm2,0,nreal);
55: for (i = 0; i < n/10; ++i) {
56: PetscInt swapi, t;
57: PetscRandomGetValueReal(rdm2,&val);
58: swapi = (PetscInt) val;
59: t = XSO[swapi-1];
60: XSO[swapi-1] = XSO[swapi];
61: XSO[swapi] = t;
62: }
63: PetscRandomDestroy(&rdm2);
65: if (vwr) {PetscIntView(n, order ? XSO : XR, vwr);}
66: PetscViewerDestroy(&vwr);
67: VecCreate(PETSC_COMM_WORLD,&x);
68: VecSetSizes(x,PETSC_DECIDE,vsize);
69: VecSetFromOptions(x);
70: VecSetRandom(x,rdm);
71: time = 0.0;
72: time1 = 0.0;
73: for (l=0; l<r; l++) { /* r loops */
74: PetscArraycpy(X,order ? XSO : XR,n);
75: PetscArraycpy(X1,order ? XSO : XR,n);
77: VecNorm(x,NORM_1,&norm1);
78: PetscTimeSubtract(&time1);
79: PetscIntSortSemiOrdered(n,X1);
80: PetscTimeAdd(&time1);
82: VecNorm(x,NORM_1,&norm1);
83: PetscTimeSubtract(&time);
84: PetscSortInt(n,X);
85: PetscTimeAdd(&time);
87: for (i=0; i<n-1; i++) {if (X[i] > X[i+1]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscSortInt() produced wrong results!");}
88: for (i=0; i<n; i++) {if (X[i] != X1[i]) SETERRQ7(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscIntSortSemiOrdered() rep %D X1[%D]:%D does not match PetscSortInt() X[%D]:%D! randomSeed %lu, orderedSeed %lu",l,i,X1[i],i,X[i],seedr,seedo);}
89: for (i=0; i<n-1; i++) {if (X1[i] > X1[i+1]) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscIntSortSemiOrdered() produced wrong results! randomSeed %lu orderedSeed %lu",seedr,seedo);}
90: PetscArrayzero(X,n);
91: PetscArrayzero(X1,n);
92: }
93: PetscPrintf(PETSC_COMM_SELF,"PetscSortInt() with %D integers, %D duplicate(s) per unique value took %g seconds\n",n,d,time/r);
94: PetscPrintf(PETSC_COMM_SELF,"PetscIntSortSemiOrdered() with %D integers, %D duplicate(s) per unique value took %g seconds\n",n,d,time1/r);
95: PetscPrintf(PETSC_COMM_SELF,"Speedup of PetscIntSortSemiOrdered() was %g (0:1 = slower, >1 means faster)\n",time/time1);
98: for (i=0; i<n; i++) { /* Init X[] */
99: PetscRandomGetValueReal(rdm,&val);
100: X[i] = val*PETSC_MAX_INT;
101: if (d > 1) X[i] = X[i] % (n/d);
102: }
103: PetscCalloc3(n,&XP,n,&X1P,n,&W);
105: time = 0.0;
106: time1 = 0.0;
107: time2 = 0.0;
108: for (l=0; l<r; l++) { /* r loops */
109: PetscArraycpy(X1, order ? XSO : XR,n);
110: PetscArraycpy(X1P,order ? XSO : XR,n);
111: PetscArraycpy(X, order ? XSO : XR,n);
112: PetscArraycpy(XP, order ? XSO : XR,n);
113: PetscArraycpy(W, order ? XSO : XR,n);
115: VecNorm(x,NORM_1,&norm1);
116: PetscTimeSubtract(&time1);
117: PetscIntSortSemiOrderedWithArray(n,X1,X1P);
118: PetscTimeAdd(&time1);
120: VecNorm(x,NORM_1,&norm1);
121: PetscTimeSubtract(&time2);
122: PetscSortIntWithArray(n,X,XP);
123: PetscTimeAdd(&time2);
125: PetscTimeSubtract(&time);
126: PetscSortIntWithArrayPair(n,W,Y,Z);
127: PetscTimeAdd(&time);
129: for (i=0; i<n-1; i++) {if (Y[i] > Y[i+1]) {PetscIntView(n,Y,0);SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscSortIntWithArray() produced wrong results!");}}
130: for (i=0; i<n-1; i++) {if (W[i] > W[i+1]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscSortIntWithArrayPair() produced wrong results!");}
131: for (i=0; i<n; i++) {if (X1P[i] != X[i]) SETERRQ7(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscIntSortSemiOrdered() rep %D X1[%D]:%D does not match PetscSortIntWithArray() X[%D]:%D! randomSeed %lu, orderedSeed %lu",l,i,X1[i],i,X[i],seedr,seedo);}
132: for (i=0; i<n-1; i++) {if (X1[i] > X1[i+1]) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscIntSortSemiOrdered() produced wrong results! randomSeed %lu orderedSeed %lu",seedr,seedo);}
133: PetscArrayzero(X1,n);
134: PetscArrayzero(X1P,n);
135: PetscArrayzero(X,n);
136: PetscArrayzero(XP,n);
137: PetscArrayzero(W,n);
138: }
139: PetscPrintf(PETSC_COMM_SELF,"PetscSortIntWithArrayPair() with %D integers, %D duplicate(s) per unique value took %g seconds\n",n,d,time/r);
140: PetscPrintf(PETSC_COMM_SELF,"PetscSortIntWithArray() with %D integers, %D duplicate(s) per unique value took %g seconds\n",n,d,time2/r);
141: PetscPrintf(PETSC_COMM_SELF,"PetscIntSortSemiOrderedWithArray() with %D integers, %D duplicate(s) per unique value took %g seconds\n",n,d,time1/r);
142: PetscPrintf(PETSC_COMM_SELF,"Speedup of PetscIntSortSemiOrderedWithArray() was %g (0:1 = slower, >1 means faster)\n",time2/time1);
143: PetscPrintf(PETSC_COMM_SELF,"SUCCEEDED\n");
145: VecDestroy(&x);
146: PetscRandomDestroy(&rdm);
147: PetscFree3(XP,X1P,W);
148: PetscFree6(X,X1,XR,XSO,Y,Z);
149: PetscFinalize();
150: return ierr;
151: }
153: /*TEST
155: test:
156: args: -n 1000 -r 10 -d 1
157: # Do not need to output timing results for test
158: filter: grep -vE "per unique value took|Speedup of "
159: TEST*/