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 PetscCall(VecAssemblyBegin(v)); 57 PetscCall(VecAssemblyEnd(v)); 58 59 /* 60 Get the components from the large multi-component vector to the small multi-component vector, 61 scale the smaller vector and then move values back to the large vector 62 */ 63 PetscCall(VecStrideSubSetGather(v, PETSC_DETERMINE, vidx, NULL, s, INSERT_VALUES)); 64 PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD)); 65 PetscCall(VecScale(s, 100.0)); 66 67 PetscCall(VecStrideSubSetScatter(s, PETSC_DETERMINE, NULL, vidx, v, ADD_VALUES)); 68 PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD)); 69 70 /* 71 Get the components from the large multi-component vector to the small multi-component vector, 72 scale the smaller vector and then move values back to the large vector 73 */ 74 PetscCall(VecStrideSubSetGather(v, 2, vidx, sidx, s, INSERT_VALUES)); 75 PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD)); 76 PetscCall(VecScale(s, 100.0)); 77 78 PetscCall(VecStrideSubSetScatter(s, 2, sidx, vidx, v, ADD_VALUES)); 79 PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD)); 80 81 PetscCall(VecStrideMax(v, 1, &miidx[0], &mvidx[0])); 82 PetscCall(VecStrideMin(v, 1, &miidx[1], &mvidx[1])); 83 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Min/Max: %" PetscInt_FMT " %g, %" PetscInt_FMT " %g\n", miidx[0], (double)mvidx[0], miidx[1], (double)mvidx[1])); 84 /* 85 Free work space. All PETSc objects should be destroyed when they 86 are no longer needed. 87 */ 88 PetscCall(VecDestroy(&v)); 89 PetscCall(VecDestroy(&s)); 90 PetscCall(PetscFinalize()); 91 return 0; 92 } 93 94 /*TEST 95 96 test: 97 filter: grep -v type | grep -v " MPI process" | grep -v Process 98 diff_args: -j 99 nsize: 2 100 101 test: 102 filter: grep -v type | grep -v " MPI process" | grep -v Process 103 output_file: output/ex45_1.out 104 diff_args: -j 105 suffix: 2 106 nsize: 1 107 args: -vec_type {{seq mpi}} 108 109 TEST*/ 110