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;

 26:   PetscInitialize(&argc, &argv, (char *)0, help);
 27:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

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

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

 88:     for (i = 0; i < n; i++) {
 90:     }
 92:     PetscArrayzero(X, n);
 93:     PetscArrayzero(X1, n);
 94:   }
 95:   PetscPrintf(PETSC_COMM_SELF, "PetscSortInt()              with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n", n, d, time / r);
 96:   PetscPrintf(PETSC_COMM_SELF, "PetscIntSortSemiOrdered()   with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n", n, d, time1 / r);
 97:   PetscPrintf(PETSC_COMM_SELF, "Speedup of PetscIntSortSemiOrdered() was %g (0:1 = slower, >1 means faster)\n", time / time1);

 99:   for (i = 0; i < n; i++) { /* Init X[] */
100:     PetscRandomGetValueReal(rdm, &val);
101:     X[i] = val * PETSC_MAX_INT;
102:     if (d > 1) X[i] = X[i] % (n / d);
103:   }
104:   PetscCalloc3(n, &XP, n, &X1P, n, &W);

106:   time  = 0.0;
107:   time1 = 0.0;
108:   time2 = 0.0;
109:   for (l = 0; l < r; l++) { /* r loops */
110:     PetscArraycpy(X1, order ? XSO : XR, n);
111:     PetscArraycpy(X1P, order ? XSO : XR, n);
112:     PetscArraycpy(X, order ? XSO : XR, n);
113:     PetscArraycpy(XP, order ? XSO : XR, n);
114:     PetscArraycpy(W, order ? XSO : XR, n);

116:     VecNorm(x, NORM_1, &norm1);
117:     PetscTimeSubtract(&time1);
118:     PetscIntSortSemiOrderedWithArray(n, X1, X1P);
119:     PetscTimeAdd(&time1);

121:     VecNorm(x, NORM_1, &norm1);
122:     PetscTimeSubtract(&time2);
123:     PetscSortIntWithArray(n, X, XP);
124:     PetscTimeAdd(&time2);

126:     PetscTimeSubtract(&time);
127:     PetscSortIntWithArrayPair(n, W, Y, Z);
128:     PetscTimeAdd(&time);

130:     for (i = 0; i < n - 1; i++) {
131:       if (Y[i] > Y[i + 1]) {
132:         PetscIntView(n, Y, 0);
133:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscSortIntWithArray() produced wrong results!");
134:       }
135:     }
137:     for (i = 0; i < n; i++) {
139:     }
141:     PetscArrayzero(X1, n);
142:     PetscArrayzero(X1P, n);
143:     PetscArrayzero(X, n);
144:     PetscArrayzero(XP, n);
145:     PetscArrayzero(W, n);
146:   }
147:   PetscPrintf(PETSC_COMM_SELF, "PetscSortIntWithArrayPair()        with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n", n, d, time / r);
148:   PetscPrintf(PETSC_COMM_SELF, "PetscSortIntWithArray()            with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n", n, d, time2 / r);
149:   PetscPrintf(PETSC_COMM_SELF, "PetscIntSortSemiOrderedWithArray() with %" PetscInt_FMT " integers, %" PetscInt_FMT " duplicate(s) per unique value took %g seconds\n", n, d, time1 / r);
150:   PetscPrintf(PETSC_COMM_SELF, "Speedup of PetscIntSortSemiOrderedWithArray() was %g (0:1 = slower, >1 means faster)\n", time2 / time1);
151:   PetscPrintf(PETSC_COMM_SELF, "SUCCEEDED\n");

153:   VecDestroy(&x);
154:   PetscRandomDestroy(&rdm);
155:   PetscFree3(XP, X1P, W);
156:   PetscFree6(X, X1, XR, XSO, Y, Z);
157:   PetscFinalize();
158:   return 0;
159: }

161: /*TEST

163:    testset:
164:      filter: grep -vE "per unique value took|Speedup of "

166:      test:
167:        suffix: small
168:        args: -n 9 -r 1

170:      test:
171:        suffix: large
172:        args: -n 1000 -r 10 -d 1
173:        # Do not need to output timing results for test

175: TEST*/