Actual source code: ex21.c
1: static const char help[] = "Test VecScatterCopy() on an SF with duplicated leaves \n\n";
3: #include <petscvec.h>
4: #include <petscsf.h>
6: /*
7: Contributed-by: "Hammond, Glenn E" <glenn.hammond@pnnl.gov>
8: */
9: int main(int argc,char* argv[])
10: {
11: PetscMPIInt size;
12: PetscInt n;
13: PetscInt *indices;
14: Vec vec;
15: Vec vec2;
16: IS is;
17: IS is2;
18: VecScatter scatter;
19: VecScatter scatter2;
21: PetscInitialize(&argc,&argv,NULL,help);
22: MPI_Comm_size(PETSC_COMM_WORLD, &size);
25: n = 4;
26: PetscMalloc1(n,&indices);
27: indices[0] = 0;
28: indices[1] = 1;
29: indices[2] = 2;
30: indices[3] = 3;
31: ISCreateGeneral(PETSC_COMM_WORLD,n,indices,PETSC_COPY_VALUES,&is);
32: PetscFree(indices);
33: VecCreateMPI(PETSC_COMM_WORLD,n,n,&vec);
35: n = 4;
36: PetscMalloc1(n,&indices);
37: indices[0] = 0;
38: indices[1] = 0;
39: indices[2] = 1;
40: indices[3] = 1;
41: ISCreateGeneral(PETSC_COMM_WORLD,n,indices,PETSC_COPY_VALUES,&is2);
42: PetscFree(indices);
43: VecCreateMPI(PETSC_COMM_WORLD,n/2,n/2,&vec2);
45: VecScatterCreate(vec,is,vec2,is2,&scatter);
46: ISDestroy(&is);
47: ISDestroy(&is2);
49: VecScatterCopy(scatter,&scatter2);
51: VecDestroy(&vec);
52: VecDestroy(&vec2);
53: VecScatterDestroy(&scatter);
54: VecScatterDestroy(&scatter2);
55: PetscFinalize();
56: }
58: /*TEST
59: test:
60: nsize: 1
61: output_file: output/empty.out
62: TEST*/