1 static const char help[] = "Test VecScatterCopy() on an SF with duplicated leaves \n\n"; 2 3 #include <petscvec.h> 4 #include <petscsf.h> 5 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; 20 21 PetscFunctionBeginUser; 22 PetscCall(PetscInitialize(&argc,&argv,NULL,help)); 23 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 24 PetscCheck(size == 1,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This example requires 1 process"); 25 26 n = 4; 27 PetscCall(PetscMalloc1(n,&indices)); 28 indices[0] = 0; 29 indices[1] = 1; 30 indices[2] = 2; 31 indices[3] = 3; 32 PetscCall(ISCreateGeneral(PETSC_COMM_WORLD,n,indices,PETSC_COPY_VALUES,&is)); 33 PetscCall(PetscFree(indices)); 34 PetscCall(VecCreateMPI(PETSC_COMM_WORLD,n,n,&vec)); 35 36 n = 4; 37 PetscCall(PetscMalloc1(n,&indices)); 38 indices[0] = 0; 39 indices[1] = 0; 40 indices[2] = 1; 41 indices[3] = 1; 42 PetscCall(ISCreateGeneral(PETSC_COMM_WORLD,n,indices,PETSC_COPY_VALUES,&is2)); 43 PetscCall(PetscFree(indices)); 44 PetscCall(VecCreateMPI(PETSC_COMM_WORLD,n/2,n/2,&vec2)); 45 46 PetscCall(VecScatterCreate(vec,is,vec2,is2,&scatter)); 47 PetscCall(ISDestroy(&is)); 48 PetscCall(ISDestroy(&is2)); 49 50 PetscCall(VecScatterCopy(scatter,&scatter2)); 51 52 PetscCall(VecDestroy(&vec)); 53 PetscCall(VecDestroy(&vec2)); 54 PetscCall(VecScatterDestroy(&scatter)); 55 PetscCall(VecScatterDestroy(&scatter2)); 56 PetscCall(PetscFinalize()); 57 } 58 59 /*TEST 60 test: 61 nsize: 1 62 output_file: output/empty.out 63 TEST*/ 64