1 2 static char help[] = "Demonstrates VecStrideSubSetScatter() and VecStrideSubSetGather().\n\n"; 3 4 /* 5 Allows one to easily pull out some components of a multi-component vector and put them in another vector. 6 7 Note that these are special cases of VecScatter 8 */ 9 10 /* 11 Include "petscvec.h" so that we can use vectors. Note that this file 12 automatically includes: 13 petscsys.h - base PETSc routines petscis.h - index sets 14 petscviewer.h - viewers 15 */ 16 17 #include <petscvec.h> 18 19 int main(int argc, char **argv) 20 { 21 Vec v, s; 22 PetscInt i, start, end, n = 8; 23 PetscScalar value; 24 const PetscInt vidx[] = {1, 2}, sidx[] = {1, 0}; 25 PetscInt miidx[2]; 26 PetscReal mvidx[2]; 27 28 PetscFunctionBeginUser; 29 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 30 PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL)); 31 32 /* 33 Create multi-component vector with 4 components 34 */ 35 PetscCall(VecCreate(PETSC_COMM_WORLD, &v)); 36 PetscCall(VecSetSizes(v, PETSC_DECIDE, n)); 37 PetscCall(VecSetBlockSize(v, 4)); 38 PetscCall(VecSetFromOptions(v)); 39 40 /* 41 Create double-component vectors 42 */ 43 PetscCall(VecCreate(PETSC_COMM_WORLD, &s)); 44 PetscCall(VecSetSizes(s, PETSC_DECIDE, n / 2)); 45 PetscCall(VecSetBlockSize(s, 2)); 46 PetscCall(VecSetFromOptions(s)); 47 48 /* 49 Set the vector values 50 */ 51 PetscCall(VecGetOwnershipRange(v, &start, &end)); 52 for (i = start; i < end; i++) { 53 value = i; 54 PetscCall(VecSetValues(v, 1, &i, &value, INSERT_VALUES)); 55 } 56 57 /* 58 Get the components from the large multi-component vector to the small multi-component vector, 59 scale the smaller vector and then move values back to the large vector 60 */ 61 PetscCall(VecStrideSubSetGather(v, PETSC_DETERMINE, vidx, NULL, s, INSERT_VALUES)); 62 PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD)); 63 PetscCall(VecScale(s, 100.0)); 64 65 PetscCall(VecStrideSubSetScatter(s, PETSC_DETERMINE, NULL, vidx, v, ADD_VALUES)); 66 PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD)); 67 68 /* 69 Get the components from the large multi-component vector to the small multi-component vector, 70 scale the smaller vector and then move values back to the large vector 71 */ 72 PetscCall(VecStrideSubSetGather(v, 2, vidx, sidx, s, INSERT_VALUES)); 73 PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD)); 74 PetscCall(VecScale(s, 100.0)); 75 76 PetscCall(VecStrideSubSetScatter(s, 2, sidx, vidx, v, ADD_VALUES)); 77 PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD)); 78 79 PetscCall(VecStrideMax(v, 1, &miidx[0], &mvidx[0])); 80 PetscCall(VecStrideMin(v, 1, &miidx[1], &mvidx[1])); 81 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Min/Max: %" PetscInt_FMT " %g, %" PetscInt_FMT " %g\n", miidx[0], (double)mvidx[0], miidx[1], (double)mvidx[1])); 82 /* 83 Free work space. All PETSc objects should be destroyed when they 84 are no longer needed. 85 */ 86 PetscCall(VecDestroy(&v)); 87 PetscCall(VecDestroy(&s)); 88 PetscCall(PetscFinalize()); 89 return 0; 90 } 91 92 /*TEST 93 94 test: 95 filter: grep -v type | grep -v " MPI process" | grep -v Process 96 diff_args: -j 97 nsize: 2 98 99 test: 100 filter: grep -v type | grep -v " MPI process" | grep -v Process 101 output_file: output/ex45_1.out 102 diff_args: -j 103 suffix: 2 104 nsize: 1 105 args: -vec_type {{seq mpi}} 106 107 TEST*/ 108