Actual source code: ex45.c


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

  4: /*T
  5:    Concepts: vectors^sub-vectors;
  6:    Processors: n

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

 10:    Note that these are special cases of VecScatter
 11: T*/

 13: /*
 14:   Include "petscvec.h" so that we can use vectors.  Note that this file
 15:   automatically includes:
 16:      petscsys.h       - base PETSc routines   petscis.h     - index sets
 17:      petscviewer.h - viewers
 18: */

 20: #include <petscvec.h>

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

 31:   PetscInitialize(&argc,&argv,(char*)0,help);
 32:   PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);

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

 42:   /*
 43:       Create double-component vectors
 44:   */
 45:   VecCreate(PETSC_COMM_WORLD,&s);
 46:   VecSetSizes(s,PETSC_DECIDE,n/2);
 47:   VecSetBlockSize(s,2);
 48:   VecSetFromOptions(s);

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

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

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

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

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

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

 94: /*TEST

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

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

109: TEST*/