Actual source code: ex6.c


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

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

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

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

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

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

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

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

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

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

 61: /*TEST

 63:    test:
 64:       suffix: 3
 65:       nsize: 3

 67: TEST*/