Actual source code: ex16.c

  2: /* Program usage:  mpirun ex1 [-help] [all PETSc options] */

  4: static char help[] = "Demonstrates VecStrideScatter() and VecStrideGather() with subvectors that are also strided.\n\n";

  6: /*T
  7:    Concepts: vectors^sub-vectors;
  8:    Processors: n
  9: T*/

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

 18:  #include petscvec.h

 22: int main(int argc,char **argv)
 23: {
 24:   Vec           v,s,r,vecs[2];               /* vectors */
 25:   int           i,start,end,n = 20,ierr;
 26:   PetscScalar   one = 1.0,value;

 28:   PetscInitialize(&argc,&argv,(char*)0,help);
 29:   PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);

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

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

 48:   vecs[0] = s;
 49:   vecs[1] = r;
 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 multi-component vector to the other vectors
 61:   */
 62:   VecStrideGatherAll(v,vecs,INSERT_VALUES);

 64:   VecView(s,PETSC_VIEWER_STDOUT_WORLD);
 65:   VecView(r,PETSC_VIEWER_STDOUT_WORLD);

 67:   VecStrideScatterAll(vecs,v,ADD_VALUES);

 69:   VecView(v,PETSC_VIEWER_STDOUT_WORLD);

 71:   /* 
 72:      Free work space.  All PETSc objects should be destroyed when they
 73:      are no longer needed.
 74:   */
 75:   VecDestroy(v);
 76:   VecDestroy(s);
 77:   VecDestroy(r);
 78:   PetscFinalize();
 79:   return 0;
 80: }
 81: