Actual source code: ex5.c

  1: static char help[] = "Demonstrates using ISLocalToGlobalMappings with block size.\n\n";

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

  6: int main(int argc, char **argv)
  7: {
  8:   PetscInt               i, n = 4, indices[] = {0, 3, 9, 12}, m = 2, input[] = {0, 2};
  9:   PetscInt               output[2], inglobals[13], outlocals[13];
 10:   ISLocalToGlobalMapping mapping;

 12:   PetscFunctionBeginUser;
 13:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));

 15:   /*
 16:       Create a local to global mapping. Each processor independently
 17:      creates a mapping
 18:   */
 19:   PetscCall(PetscIntView(n, indices, PETSC_VIEWER_STDOUT_WORLD));
 20:   PetscCall(ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 2, n, indices, PETSC_COPY_VALUES, &mapping));

 22:   /*
 23:      Map a set of local indices to their global values
 24:   */
 25:   PetscCall(PetscIntView(m, input, PETSC_VIEWER_STDOUT_WORLD));
 26:   PetscCall(ISLocalToGlobalMappingApply(mapping, m, input, output));
 27:   PetscCall(PetscIntView(m, output, PETSC_VIEWER_STDOUT_WORLD));

 29:   /*
 30:      Map some global indices to local, retaining the ones without a local index by -1
 31:   */
 32:   for (i = 0; i < 13; i++) inglobals[i] = i;
 33:   PetscCall(PetscIntView(13, inglobals, PETSC_VIEWER_STDOUT_WORLD));
 34:   PetscCall(ISGlobalToLocalMappingApply(mapping, IS_GTOLM_MASK, 13, inglobals, NULL, outlocals));
 35:   PetscCall(PetscIntView(13, outlocals, PETSC_VIEWER_STDOUT_WORLD));

 37:   /*
 38:      Map some block global indices to local, dropping the ones without a local index.
 39:   */
 40:   PetscCall(PetscIntView(13, inglobals, PETSC_VIEWER_STDOUT_WORLD));
 41:   PetscCall(ISGlobalToLocalMappingApplyBlock(mapping, IS_GTOLM_DROP, 13, inglobals, &m, outlocals));
 42:   PetscCall(PetscIntView(m, outlocals, PETSC_VIEWER_STDOUT_WORLD));

 44:   /*
 45:      Free the space used by the local to global mapping
 46:   */
 47:   PetscCall(ISLocalToGlobalMappingDestroy(&mapping));

 49:   PetscCall(PetscFinalize());
 50:   return 0;
 51: }

 53: /*TEST

 55:    test:

 57: TEST*/