Actual source code: ex23.c
petsc-3.12.0 2019-09-29
2: static char help[] = "Scatters from a parallel vector to a sequential vector.\n\
3: Using a blocked send and a strided receive.\n\n";
5: /*
6: 0 1 2 3 | 4 5 6 7 || 8 9 10 11
8: Scatter first and third block to first processor and
9: second and third block to second processor
10: */
12: #include <petscvec.h>
14: int main(int argc,char **argv)
15: {
17: PetscInt i,blocks[2],nlocal;
18: PetscMPIInt size,rank;
19: PetscScalar value;
20: Vec x,y;
21: IS is1,is2;
22: VecScatter ctx = 0;
23: PetscViewer subviewer;
25: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
26: MPI_Comm_size(PETSC_COMM_WORLD,&size);
27: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
29: if (size != 2) SETERRQ(PETSC_COMM_SELF,1,"Must run with 2 processors");
31: /* create two vectors */
32: if (!rank) nlocal = 8;
33: else nlocal = 4;
34: VecCreate(PETSC_COMM_WORLD,&x);
35: VecSetSizes(x,nlocal,12);
36: VecSetFromOptions(x);
37: VecCreate(PETSC_COMM_SELF,&y);
38: VecSetSizes(y,8,PETSC_DECIDE);
39: VecSetFromOptions(y);
41: /* create two index sets */
42: if (!rank) {
43: blocks[0] = 0; blocks[1] = 2;
44: } else {
45: blocks[0] = 1; blocks[1] = 2;
46: }
47: ISCreateBlock(PETSC_COMM_SELF,4,2,blocks,PETSC_COPY_VALUES,&is1);
48: ISCreateStride(PETSC_COMM_SELF,8,0,1,&is2);
50: for (i=0; i<12; i++) {
51: value = i;
52: VecSetValues(x,1,&i,&value,INSERT_VALUES);
53: }
54: VecAssemblyBegin(x);
55: VecAssemblyEnd(x);
57: VecScatterCreate(x,is1,y,is2,&ctx);
58: VecScatterBegin(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
59: VecScatterEnd(ctx,x,y,INSERT_VALUES,SCATTER_FORWARD);
60: VecScatterDestroy(&ctx);
62: PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&subviewer);
63: VecView(y,subviewer);
64: PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&subviewer);
66: VecDestroy(&x);
67: VecDestroy(&y);
68: ISDestroy(&is1);
69: ISDestroy(&is2);
71: PetscFinalize();
72: return ierr;
73: }
77: /*TEST
79: testset:
80: nsize: 2
81: output_file: output/ex23_1.out
82: filter: grep -v " type:"
83: test:
84: suffix: standard
85: args: -vec_type standard
86: test:
87: requires: cuda
88: suffix: cuda
89: args: -vec_type cuda
90: test:
91: requires: viennacl
92: suffix: viennacl
93: args: -vec_type viennacl
94: TEST*/