1 static char help[] = "Demonstrates use of VecCreateGhost().\n\n"; 2 3 /* 4 Description: Ghost padding is one way to handle local calculations that 5 involve values from other processors. VecCreateGhost() provides 6 a way to create vectors with extra room at the end of the vector 7 array to contain the needed ghost values from other processors, 8 vector computations are otherwise unaffected. 9 */ 10 11 /* 12 Include "petscvec.h" so that we can use vectors. Note that this file 13 automatically includes: 14 petscsys.h - base PETSc routines petscis.h - index sets 15 petscviewer.h - viewers 16 */ 17 #include <petscvec.h> 18 19 int main(int argc, char **argv) 20 { 21 PetscMPIInt rank, size; 22 PetscInt nlocal = 6, nghost = 2, ifrom[2], i, rstart, rend; 23 PetscBool flg, flg2, flg3, flg4, flg5; 24 PetscScalar value, *array, *tarray = 0; 25 Vec lx, gx, gxs; 26 IS ghost; 27 ISLocalToGlobalMapping mapping; 28 29 PetscFunctionBeginUser; 30 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 31 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 32 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 33 PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run example with two processors"); 34 35 /* 36 Construct a two dimensional graph connecting nlocal degrees of 37 freedom per processor. From this we will generate the global 38 indices of needed ghost values 39 40 For simplicity we generate the entire graph on each processor: 41 in real application the graph would stored in parallel, but this 42 example is only to demonstrate the management of ghost padding 43 with VecCreateGhost(). 44 45 In this example we consider the vector as representing 46 degrees of freedom in a one dimensional grid with periodic 47 boundary conditions. 48 49 ----Processor 1--------- ----Processor 2 -------- 50 0 1 2 3 4 5 6 7 8 9 10 11 51 |----| 52 |-------------------------------------------------| 53 54 */ 55 56 if (rank == 0) { 57 ifrom[0] = 11; 58 ifrom[1] = 6; 59 } else { 60 ifrom[0] = 0; 61 ifrom[1] = 5; 62 } 63 64 /* 65 Create the vector with two slots for ghost points. Note that both 66 the local vector (lx) and the global vector (gx) share the same 67 array for storing vector values. 68 */ 69 PetscCall(PetscOptionsHasName(NULL, NULL, "-allocate", &flg)); 70 PetscCall(PetscOptionsHasName(NULL, NULL, "-vecmpisetghost", &flg2)); 71 PetscCall(PetscOptionsHasName(NULL, NULL, "-minvalues", &flg3)); 72 if (flg) { 73 PetscCall(PetscMalloc1(nlocal + nghost, &tarray)); 74 PetscCall(VecCreateGhostWithArray(PETSC_COMM_WORLD, nlocal, PETSC_DECIDE, nghost, ifrom, tarray, &gxs)); 75 } else if (flg2) { 76 PetscCall(VecCreate(PETSC_COMM_WORLD, &gxs)); 77 PetscCall(VecSetType(gxs, VECMPI)); 78 PetscCall(VecSetSizes(gxs, nlocal, PETSC_DECIDE)); 79 PetscCall(VecMPISetGhost(gxs, nghost, ifrom)); 80 } else { 81 PetscCall(VecCreateGhost(PETSC_COMM_WORLD, nlocal, PETSC_DECIDE, nghost, ifrom, &gxs)); 82 } 83 84 /* 85 Test VecDuplicate() 86 */ 87 PetscCall(VecDuplicate(gxs, &gx)); 88 PetscCall(VecDestroy(&gxs)); 89 90 /* 91 Access the local representation 92 */ 93 PetscCall(VecGhostGetLocalForm(gx, &lx)); 94 95 /* 96 Set the values from 0 to 12 into the "global" vector 97 */ 98 PetscCall(VecGetOwnershipRange(gx, &rstart, &rend)); 99 for (i = rstart; i < rend; i++) { 100 value = (PetscScalar)i; 101 PetscCall(VecSetValues(gx, 1, &i, &value, INSERT_VALUES)); 102 } 103 PetscCall(VecAssemblyBegin(gx)); 104 PetscCall(VecAssemblyEnd(gx)); 105 106 PetscCall(VecGhostUpdateBegin(gx, INSERT_VALUES, SCATTER_FORWARD)); 107 PetscCall(VecGhostUpdateEnd(gx, INSERT_VALUES, SCATTER_FORWARD)); 108 109 /* 110 Print out each vector, including the ghost padding region. 111 */ 112 PetscCall(VecGetArray(lx, &array)); 113 for (i = 0; i < nlocal + nghost; i++) PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%" PetscInt_FMT " %g\n", i, (double)PetscRealPart(array[i]))); 114 PetscCall(VecRestoreArray(lx, &array)); 115 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT)); 116 PetscCall(VecGhostRestoreLocalForm(gx, &lx)); 117 118 /* Another test that sets ghost values and then accumulates onto the owning processors using MIN_VALUES */ 119 if (flg3) { 120 if (rank == 0) PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "\nTesting VecGhostUpdate with MIN_VALUES\n")); 121 PetscCall(VecGhostGetLocalForm(gx, &lx)); 122 PetscCall(VecGetArray(lx, &array)); 123 for (i = 0; i < nghost; i++) array[nlocal + i] = rank ? (PetscScalar)4 : (PetscScalar)8; 124 PetscCall(VecRestoreArray(lx, &array)); 125 PetscCall(VecGhostRestoreLocalForm(gx, &lx)); 126 127 PetscCall(VecGhostUpdateBegin(gx, MIN_VALUES, SCATTER_REVERSE)); 128 PetscCall(VecGhostUpdateEnd(gx, MIN_VALUES, SCATTER_REVERSE)); 129 130 PetscCall(VecGhostGetLocalForm(gx, &lx)); 131 PetscCall(VecGetArray(lx, &array)); 132 133 for (i = 0; i < nlocal + nghost; i++) PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "%" PetscInt_FMT " %g\n", i, (double)PetscRealPart(array[i]))); 134 PetscCall(VecRestoreArray(lx, &array)); 135 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT)); 136 PetscCall(VecGhostRestoreLocalForm(gx, &lx)); 137 } 138 139 PetscCall(PetscOptionsHasName(NULL, NULL, "-vecghostgetghostis", &flg4)); 140 if (flg4) { 141 PetscCall(VecGhostGetGhostIS(gx, &ghost)); 142 PetscCall(ISView(ghost, PETSC_VIEWER_STDOUT_WORLD)); 143 } 144 PetscCall(PetscOptionsHasName(NULL, NULL, "-getgtlmapping", &flg5)); 145 if (flg5) { 146 PetscCall(VecGetLocalToGlobalMapping(gx, &mapping)); 147 PetscCall(ISLocalToGlobalMappingView(mapping, NULL)); 148 } 149 150 PetscCall(VecDestroy(&gx)); 151 152 if (flg) PetscCall(PetscFree(tarray)); 153 PetscCall(PetscFinalize()); 154 return 0; 155 } 156 157 /*TEST 158 159 test: 160 nsize: 2 161 162 test: 163 suffix: 2 164 nsize: 2 165 args: -allocate 166 output_file: output/ex9_1.out 167 168 test: 169 suffix: 3 170 nsize: 2 171 args: -vecmpisetghost 172 output_file: output/ex9_1.out 173 174 test: 175 suffix: 4 176 nsize: 2 177 args: -minvalues 178 output_file: output/ex9_2.out 179 requires: !complex 180 181 test: 182 suffix: 5 183 nsize: 2 184 args: -vecghostgetghostis 185 186 test: 187 suffix: 6 188 nsize: 2 189 args: -getgtlmapping 190 191 TEST*/ 192