Actual source code: ex6.c

  1: static char help[] = "Tests ISComplement().\n\n";

  3: #include <petscis.h>
  4: #include <petscviewer.h>

  6: int main(int argc, char **argv)
  7: {
  8:   PetscMPIInt rank, size;
  9:   PetscInt    i, j, n, cnt = 0, rstart, rend;
 10:   PetscBool   flg;
 11:   IS          is[2], isc;

 13:   PetscFunctionBeginUser;
 14:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
 15:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
 16:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));

 18:   n      = 3 * size;                    /* Number of local indices, same on each process. */
 19:   rstart = 3 * (size + 2) * rank;       /* start of local range */
 20:   rend   = 3 * (size + 2) * (rank + 1); /* end of local range */
 21:   for (i = 0; i < 2; i++) {
 22:     PetscCall(ISCreate(PETSC_COMM_WORLD, &is[i]));
 23:     PetscCall(ISSetType(is[i], ISGENERAL));
 24:   }
 25:   {
 26:     PetscBool *mask;

 28:     PetscCall(PetscCalloc1(rend - rstart, &mask));
 29:     for (i = 0; i < 3; i++) {
 30:       for (j = 0; j < size; j++) mask[i * (size + 2) + j] = PETSC_TRUE;
 31:     }
 32:     PetscCall(ISGeneralSetIndicesFromMask(is[0], rstart, rend, mask));
 33:     PetscCall(PetscFree(mask));
 34:   }
 35:   {
 36:     PetscInt *indices;

 38:     PetscCall(PetscMalloc1(n, &indices));
 39:     for (i = 0; i < 3; i++) {
 40:       for (j = 0; j < size; j++) indices[cnt++] = rstart + i * (size + 2) + j;
 41:     }
 42:     PetscCheck(cnt == n, PETSC_COMM_SELF, PETSC_ERR_PLIB, "inconsistent count");
 43:     PetscCall(ISGeneralSetIndices(is[1], n, indices, PETSC_COPY_VALUES));
 44:     PetscCall(PetscFree(indices));
 45:   }

 47:   PetscCall(ISEqual(is[0], is[1], &flg));
 48:   PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "is[0] should be equal to is[1]");

 50:   PetscCall(ISComplement(is[0], rstart, rend, &isc));
 51:   PetscCall(ISView(is[0], PETSC_VIEWER_STDOUT_WORLD));
 52:   PetscCall(ISView(isc, PETSC_VIEWER_STDOUT_WORLD));

 54:   for (i = 0; i < 2; i++) PetscCall(ISDestroy(&is[i]));
 55:   PetscCall(ISDestroy(&isc));
 56:   PetscCall(PetscFinalize());
 57:   return 0;
 58: }

 60: /*TEST

 62:    test:
 63:       suffix: 3
 64:       nsize: 3

 66: TEST*/