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 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 24 25 /* create the large DMDA and set coordinates (which we will copy down to the small DA). */ 26 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)); 27 PetscCall(DMSetFromOptions(da)); 28 PetscCall(DMSetUp(da)); 29 PetscCall(DMDASetUniformCoordinates(da,0.0,1.0,0.0,1.0,0.0,1.0)); 30 /* Just as a simple example we use the coordinates as the variables in the vectors we wish to examine. */ 31 PetscCall(DMGetCoordinates(da,&xy)); 32 /* 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 */ 33 PetscCall(VecView(xy,0)); 34 35 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 36 if (rank == 0) comm = MPI_COMM_SELF; 37 else comm = MPI_COMM_NULL; 38 39 PetscCall(DMPatchZoom(da,lower,upper,comm,&sda, NULL,&sf)); 40 if (rank == 0) { 41 PetscCall(DMCreateGlobalVector(sda,&sxy)); 42 } else { 43 PetscCall(VecCreateSeq(PETSC_COMM_SELF,0,&sxy)); 44 } 45 /* A PetscSF can also be used as a VecScatter context */ 46 PetscCall(VecScatterBegin(sf,xy,sxy,INSERT_VALUES,SCATTER_FORWARD)); 47 PetscCall(VecScatterEnd(sf,xy,sxy,INSERT_VALUES,SCATTER_FORWARD)); 48 /* Only rank == 0 has the entries of the patch, so run code only at that rank */ 49 if (rank == 0) { 50 Field **vars; 51 DMDALocalInfo info; 52 PetscInt i,j; 53 PetscScalar sum = 0; 54 55 /* 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 */ 56 PetscCall(VecView(sxy,PETSC_VIEWER_STDOUT_SELF)); 57 /* Compute some trivial statistic of the coordinates */ 58 PetscCall(DMDAGetLocalInfo(sda,&info)); 59 PetscCall(DMDAVecGetArray(sda,sxy,&vars)); 60 /* Loop over the patch of the entire domain */ 61 for (j=info.ys; j<info.ys+info.ym; j++) { 62 for (i=info.xs; i<info.xs+info.xm; i++) { 63 sum += vars[j][i].x; 64 } 65 } 66 PetscCall(PetscPrintf(PETSC_COMM_SELF,"The sum of the x coordinates is %g\n",(double)PetscRealPart(sum))); 67 PetscCall(DMDAVecRestoreArray(sda,sxy,&vars)); 68 } 69 70 PetscCall(VecDestroy(&sxy)); 71 PetscCall(PetscSFDestroy(&sf)); 72 PetscCall(DMDestroy(&sda)); 73 PetscCall(DMDestroy(&da)); 74 PetscCall(PetscFinalize()); 75 return 0; 76 } 77 78 /*TEST 79 80 test: 81 82 test: 83 nsize: 4 84 suffix: 2 85 86 TEST*/ 87