1 2 static char help[] = "Takes a patch of a large DMDA vector to one process.\n\n"; 3 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 #include <petscdmpatch.h> 7 #include <petscsf.h> 8 9 typedef struct { 10 PetscScalar x, y; 11 } Field; 12 13 int main(int argc, char **argv) 14 { 15 Vec xy, sxy; 16 DM da, sda = NULL; 17 PetscSF sf; 18 PetscInt m = 10, n = 10, dof = 2; 19 MatStencil lower = {0, 3, 2, 0}, upper = {0, 7, 8, 0}; /* These are in the order of the z, y, x, logical coordinates, the fourth entry is ignored */ 20 MPI_Comm comm; 21 PetscMPIInt rank; 22 23 PetscFunctionBeginUser; 24 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 25 26 /* create the large DMDA and set coordinates (which we will copy down to the small DA). */ 27 PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, m, n, PETSC_DECIDE, PETSC_DECIDE, dof, 1, 0, 0, &da)); 28 PetscCall(DMSetFromOptions(da)); 29 PetscCall(DMSetUp(da)); 30 PetscCall(DMDASetUniformCoordinates(da, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0)); 31 /* Just as a simple example we use the coordinates as the variables in the vectors we wish to examine. */ 32 PetscCall(DMGetCoordinates(da, &xy)); 33 /* The vector entries are displayed in the "natural" ordering on the two dimensional grid; interlaced x and y with with the x variable increasing more rapidly than the y */ 34 PetscCall(VecView(xy, 0)); 35 36 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 37 if (rank == 0) comm = MPI_COMM_SELF; 38 else comm = MPI_COMM_NULL; 39 40 PetscCall(DMPatchZoom(da, lower, upper, comm, &sda, NULL, &sf)); 41 if (rank == 0) { 42 PetscCall(DMCreateGlobalVector(sda, &sxy)); 43 } else { 44 PetscCall(VecCreateSeq(PETSC_COMM_SELF, 0, &sxy)); 45 } 46 /* A PetscSF can also be used as a VecScatter context */ 47 PetscCall(VecScatterBegin(sf, xy, sxy, INSERT_VALUES, SCATTER_FORWARD)); 48 PetscCall(VecScatterEnd(sf, xy, sxy, INSERT_VALUES, SCATTER_FORWARD)); 49 /* Only rank == 0 has the entries of the patch, so run code only at that rank */ 50 if (rank == 0) { 51 Field **vars; 52 DMDALocalInfo info; 53 PetscInt i, j; 54 PetscScalar sum = 0; 55 56 /* The vector entries of the patch are displayed in the "natural" ordering on the two grid; interlaced x and y with with the x variable increasing more rapidly */ 57 PetscCall(VecView(sxy, PETSC_VIEWER_STDOUT_SELF)); 58 /* Compute some trivial statistic of the coordinates */ 59 PetscCall(DMDAGetLocalInfo(sda, &info)); 60 PetscCall(DMDAVecGetArray(sda, sxy, &vars)); 61 /* Loop over the patch of the entire domain */ 62 for (j = info.ys; j < info.ys + info.ym; j++) { 63 for (i = info.xs; i < info.xs + info.xm; i++) sum += vars[j][i].x; 64 } 65 PetscCall(PetscPrintf(PETSC_COMM_SELF, "The sum of the x coordinates is %g\n", (double)PetscRealPart(sum))); 66 PetscCall(DMDAVecRestoreArray(sda, sxy, &vars)); 67 } 68 69 PetscCall(VecDestroy(&sxy)); 70 PetscCall(PetscSFDestroy(&sf)); 71 PetscCall(DMDestroy(&sda)); 72 PetscCall(DMDestroy(&da)); 73 PetscCall(PetscFinalize()); 74 return 0; 75 } 76 77 /*TEST 78 79 test: 80 81 test: 82 nsize: 4 83 suffix: 2 84 85 TEST*/ 86