Actual source code: ex45.c

  1: static char help[] = "Demonstrates VecStrideSubSetScatter() and VecStrideSubSetGather().\n\n";

  3: /*
  4:    Allows one to easily pull out some components of a multi-component vector and put them in another vector.

  6:    Note that these are special cases of VecScatter
  7: */

  9: /*
 10:   Include "petscvec.h" so that we can use vectors.  Note that this file
 11:   automatically includes:
 12:      petscsys.h       - base PETSc routines   petscis.h     - index sets
 13:      petscviewer.h - viewers
 14: */

 16: #include <petscvec.h>

 18: int main(int argc, char **argv)
 19: {
 20:   Vec            v, s;
 21:   PetscInt       i, start, end, n = 8;
 22:   PetscScalar    value;
 23:   const PetscInt vidx[] = {1, 2}, sidx[] = {1, 0};
 24:   PetscInt       miidx[2];
 25:   PetscReal      mvidx[2];

 27:   PetscFunctionBeginUser;
 28:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
 29:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));

 31:   /*
 32:       Create multi-component vector with 4 components
 33:   */
 34:   PetscCall(VecCreate(PETSC_COMM_WORLD, &v));
 35:   PetscCall(VecSetSizes(v, PETSC_DECIDE, n));
 36:   PetscCall(VecSetBlockSize(v, 4));
 37:   PetscCall(VecSetFromOptions(v));

 39:   /*
 40:       Create double-component vectors
 41:   */
 42:   PetscCall(VecCreate(PETSC_COMM_WORLD, &s));
 43:   PetscCall(VecSetSizes(s, PETSC_DECIDE, n / 2));
 44:   PetscCall(VecSetBlockSize(s, 2));
 45:   PetscCall(VecSetFromOptions(s));

 47:   /*
 48:      Set the vector values
 49:   */
 50:   PetscCall(VecGetOwnershipRange(v, &start, &end));
 51:   for (i = start; i < end; i++) {
 52:     value = i;
 53:     PetscCall(VecSetValues(v, 1, &i, &value, INSERT_VALUES));
 54:   }
 55:   PetscCall(VecAssemblyBegin(v));
 56:   PetscCall(VecAssemblyEnd(v));

 58:   /*
 59:      Get the components from the large multi-component vector to the small multi-component vector,
 60:      scale the smaller vector and then move values back to the large vector
 61:   */
 62:   PetscCall(VecStrideSubSetGather(v, PETSC_DETERMINE, vidx, NULL, s, INSERT_VALUES));
 63:   PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD));
 64:   PetscCall(VecScale(s, 100.0));

 66:   PetscCall(VecStrideSubSetScatter(s, PETSC_DETERMINE, NULL, vidx, v, ADD_VALUES));
 67:   PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD));

 69:   /*
 70:      Get the components from the large multi-component vector to the small multi-component vector,
 71:      scale the smaller vector and then move values back to the large vector
 72:   */
 73:   PetscCall(VecStrideSubSetGather(v, 2, vidx, sidx, s, INSERT_VALUES));
 74:   PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD));
 75:   PetscCall(VecScale(s, 100.0));

 77:   PetscCall(VecStrideSubSetScatter(s, 2, sidx, vidx, v, ADD_VALUES));
 78:   PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD));

 80:   PetscCall(VecStrideMax(v, 1, &miidx[0], &mvidx[0]));
 81:   PetscCall(VecStrideMin(v, 1, &miidx[1], &mvidx[1]));
 82:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Min/Max: %" PetscInt_FMT " %g, %" PetscInt_FMT " %g\n", miidx[0], (double)mvidx[0], miidx[1], (double)mvidx[1]));
 83:   /*
 84:      Free work space.  All PETSc objects should be destroyed when they
 85:      are no longer needed.
 86:   */
 87:   PetscCall(VecDestroy(&v));
 88:   PetscCall(VecDestroy(&s));
 89:   PetscCall(PetscFinalize());
 90:   return 0;
 91: }

 93: /*TEST

 95:    test:
 96:       filter: grep -v type | grep -v " MPI process" | grep -v Process
 97:       diff_args: -j
 98:       nsize: 2

100:    test:
101:       filter: grep -v type | grep -v " MPI process" | grep -v Process
102:       output_file: output/ex45_1.out
103:       diff_args: -j
104:       suffix: 2
105:       nsize: 1
106:       args: -vec_type {{seq mpi}}

108: TEST*/