Actual source code: ex22.c


  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: {
  8:   PetscInt       n = 5,N,i;
  9:   PetscMPIInt    size,rank;
 10:   PetscScalar    value,zero = 0.0;
 11:   Vec            x,y;
 12:   IS             is1,is2;
 13:   VecScatter     ctx = 0;

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

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

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

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

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

 41:   VecSet(y,zero);

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

 48:   VecView(y,PETSC_VIEWER_STDOUT_WORLD);

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

 55:   PetscFinalize();
 56:   return 0;
 57: }

 59: /*TEST

 61:    testset:
 62:       nsize: 4
 63:       output_file: output/ex22_1.out
 64:       filter: grep -v "  type:"
 65:       diff_args: -j
 66:       test:
 67:         suffix: standard
 68:         args: -vec_type standard
 69:       test:
 70:         requires: cuda
 71:         suffix: cuda
 72:         args: -vec_type cuda
 73:       test:
 74:         requires: viennacl
 75:         suffix:  viennacl
 76:         args: -vec_type viennacl
 77:       test:
 78:         requires: !sycl kokkos_kernels
 79:         suffix: kokkos
 80:         args: -vec_type kokkos
 81:       test:
 82:         requires: hip
 83:         suffix: hip
 84:         args: -vec_type hip

 86: TEST*/