Actual source code: ex22.c

petsc-3.12.2 2019-11-22
Report Typos and Errors

  2: static char help[] = "Scatters from a parallel vector to a parallel vector.\n\n";

  4:  #include <petscvec.h>

  6: int main(int argc,char **argv)
  7: {
  9:   PetscInt       n = 5,N,i;
 10:   PetscMPIInt    size,rank;
 11:   PetscScalar    value,zero = 0.0;
 12:   Vec            x,y;
 13:   IS             is1,is2;
 14:   VecScatter     ctx = 0;

 16:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
 17:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 18:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);

 20:   /* create two vectors */
 21:   N    = size*n;
 22:   VecCreate(PETSC_COMM_WORLD,&y);
 23:   VecSetSizes(y,n,PETSC_DECIDE);
 24:   VecSetFromOptions(y);

 26:   VecCreate(PETSC_COMM_WORLD,&x);
 27:   VecSetSizes(x,n,PETSC_DECIDE);
 28:   VecSetFromOptions(x);

 30:   /* create two index sets */
 31:   ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&is1);
 32:   ISCreateStride(PETSC_COMM_WORLD,n,(n*(rank+1))%N,1,&is2);

 34:   /* fill local part of parallel vector x */
 35:   value = (PetscScalar)(rank+1);
 36:   for (i=n*rank; i<n*(rank+1); i++) {
 37:     VecSetValues(x,1,&i,&value,INSERT_VALUES);
 38:   }
 39:   VecAssemblyBegin(x);
 40:   VecAssemblyEnd(x);

 42:   VecSet(y,zero);

 44:   VecScatterCreate(x,is1,y,is2,&ctx);
 45:   VecScatterBegin(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
 46:   VecScatterEnd(ctx,x,y,ADD_VALUES,SCATTER_FORWARD);
 47:   VecScatterDestroy(&ctx);

 49:   VecView(y,PETSC_VIEWER_STDOUT_WORLD);

 51:   VecDestroy(&x);
 52:   VecDestroy(&y);
 53:   ISDestroy(&is1);
 54:   ISDestroy(&is2);

 56:   PetscFinalize();
 57:   return ierr;
 58: }



 62: /*TEST

 64:    testset:
 65:       nsize: 4
 66:       output_file: output/ex22_1.out
 67:       filter: grep -v "  type:"
 68:       test:
 69:         suffix: standard
 70:         args: -vec_type standard
 71:       test:
 72:         requires: cuda
 73:         suffix: cuda
 74:         args: -vec_type cuda
 75:       test:
 76:        requires: cuda
 77:         suffix: cuda_sf
 78:         args: -vec_type cuda -vecscatter_type sf
 79:       test:
 80:         requires: viennacl
 81:         suffix:  viennacl
 82:         args: -vec_type viennacl

 84: TEST*/