Actual source code: ex52.c
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: {
14: PetscInt i,l,n=100,r=10,d=1,vsize=1;
15: PetscInt *X,*X1,*XR,*XSO,*W,*Y,*Z,*XP,*X1P;
16: PetscReal val,norm1,nreal;
17: PetscRandom rdm,rdm2;
18: PetscLogDouble time, time1, time2;
19: PetscMPIInt size;
20: PetscViewer vwr;
21: Vec x;
22: unsigned long seedr, seedo;
23: PetscBool order=PETSC_FALSE;
25: PetscInitialize(&argc,&argv,(char*)0,help);
26: MPI_Comm_size(PETSC_COMM_WORLD,&size);
29: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
30: PetscOptionsGetInt(NULL,NULL,"-r",&r,NULL);
31: PetscOptionsGetInt(NULL,NULL,"-d",&d,NULL);
32: PetscOptionsGetInt(NULL,NULL,"-vsize",&vsize,NULL);
33: PetscOptionsGetBool(NULL,NULL,"-order",NULL,&order);
34: PetscOptionsGetViewer(PETSC_COMM_WORLD,NULL,NULL,"-array_view",&vwr,NULL,NULL);
37: PetscCalloc6(n,&X,n,&X1,n,&XR,n,&XSO,n,&Y,n,&Z);
38: PetscRandomCreate(PETSC_COMM_SELF,&rdm);
39: PetscRandomSetFromOptions(rdm);
40: PetscRandomGetSeed(rdm, &seedr);
42: for (i=0; i<n; ++i) {
43: PetscRandomGetValueReal(rdm,&val);
44: XR[i] = val*PETSC_MAX_INT;
45: if (d > 1) XR[i] = XR[i] % (n/d);
46: XSO[i] = i;
47: if (d > 1) XSO[i] = XSO[i] % (n/d);
48: }
50: nreal = (PetscReal) n;
51: PetscRandomCreate(PETSC_COMM_SELF,&rdm2);
52: PetscRandomGetSeed(rdm, &seedo);
53: PetscRandomSetInterval(rdm2,0,nreal);
54: for (i = 0; i < n/10; ++i) {
55: PetscInt swapi, t;
56: PetscRandomGetValueReal(rdm2,&val);
57: swapi = (PetscInt) val;
58: t = XSO[swapi-1];
59: XSO[swapi-1] = XSO[swapi];
60: XSO[swapi] = t;
61: }
62: PetscRandomDestroy(&rdm2);
64: if (vwr) PetscIntView(n, order ? XSO : XR, vwr);
65: PetscViewerDestroy(&vwr);
66: VecCreate(PETSC_COMM_WORLD,&x);
67: VecSetSizes(x,PETSC_DECIDE,vsize);
68: VecSetFromOptions(x);
69: VecSetRandom(x,rdm);
70: time = 0.0;
71: time1 = 0.0;
72: for (l=0; l<r; l++) { /* r loops */
73: PetscArraycpy(X,order ? XSO : XR,n);
74: PetscArraycpy(X1,order ? XSO : XR,n);
76: VecNorm(x,NORM_1,&norm1);
77: PetscTimeSubtract(&time1);
78: PetscIntSortSemiOrdered(n,X1);
79: PetscTimeAdd(&time1);
81: VecNorm(x,NORM_1,&norm1);
82: PetscTimeSubtract(&time);
83: PetscSortInt(n,X);
84: PetscTimeAdd(&time);
89: PetscArrayzero(X,n);
90: PetscArrayzero(X1,n);
91: }
92: PetscPrintf(PETSC_COMM_SELF,"PetscSortInt() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n",n,d,time/r);
93: PetscPrintf(PETSC_COMM_SELF,"PetscIntSortSemiOrdered() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n",n,d,time1/r);
94: PetscPrintf(PETSC_COMM_SELF,"Speedup of PetscIntSortSemiOrdered() was %g (0:1 = slower, >1 means faster)\n",time/time1);
96: for (i=0; i<n; i++) { /* Init X[] */
97: PetscRandomGetValueReal(rdm,&val);
98: X[i] = val*PETSC_MAX_INT;
99: if (d > 1) X[i] = X[i] % (n/d);
100: }
101: PetscCalloc3(n,&XP,n,&X1P,n,&W);
103: time = 0.0;
104: time1 = 0.0;
105: time2 = 0.0;
106: for (l=0; l<r; l++) { /* r loops */
107: PetscArraycpy(X1, order ? XSO : XR,n);
108: PetscArraycpy(X1P,order ? XSO : XR,n);
109: PetscArraycpy(X, order ? XSO : XR,n);
110: PetscArraycpy(XP, order ? XSO : XR,n);
111: PetscArraycpy(W, order ? XSO : XR,n);
113: VecNorm(x,NORM_1,&norm1);
114: PetscTimeSubtract(&time1);
115: PetscIntSortSemiOrderedWithArray(n,X1,X1P);
116: PetscTimeAdd(&time1);
118: VecNorm(x,NORM_1,&norm1);
119: PetscTimeSubtract(&time2);
120: PetscSortIntWithArray(n,X,XP);
121: PetscTimeAdd(&time2);
123: PetscTimeSubtract(&time);
124: PetscSortIntWithArrayPair(n,W,Y,Z);
125: PetscTimeAdd(&time);
127: 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!");}}
131: PetscArrayzero(X1,n);
132: PetscArrayzero(X1P,n);
133: PetscArrayzero(X,n);
134: PetscArrayzero(XP,n);
135: PetscArrayzero(W,n);
136: }
137: PetscPrintf(PETSC_COMM_SELF,"PetscSortIntWithArrayPair() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n",n,d,time/r);
138: PetscPrintf(PETSC_COMM_SELF,"PetscSortIntWithArray() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n",n,d,time2/r);
139: PetscPrintf(PETSC_COMM_SELF,"PetscIntSortSemiOrderedWithArray() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n",n,d,time1/r);
140: PetscPrintf(PETSC_COMM_SELF,"Speedup of PetscIntSortSemiOrderedWithArray() was %g (0:1 = slower, >1 means faster)\n",time2/time1);
141: PetscPrintf(PETSC_COMM_SELF,"SUCCEEDED\n");
143: VecDestroy(&x);
144: PetscRandomDestroy(&rdm);
145: PetscFree3(XP,X1P,W);
146: PetscFree6(X,X1,XR,XSO,Y,Z);
147: PetscFinalize();
148: return 0;
149: }
151: /*TEST
153: testset:
154: filter: grep -vE "per unique value took|Speedup of "
156: test:
157: suffix: small
158: args: -n 9 -r 1
160: test:
161: suffix: large
162: args: -n 1000 -r 10 -d 1
163: # Do not need to output timing results for test
165: TEST*/