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++) { 64 sum += vars[j][i].x; 65 } 66 } 67 PetscCall(PetscPrintf(PETSC_COMM_SELF,"The sum of the x coordinates is %g\n",(double)PetscRealPart(sum))); 68 PetscCall(DMDAVecRestoreArray(sda,sxy,&vars)); 69 } 70 71 PetscCall(VecDestroy(&sxy)); 72 PetscCall(PetscSFDestroy(&sf)); 73 PetscCall(DMDestroy(&sda)); 74 PetscCall(DMDestroy(&da)); 75 PetscCall(PetscFinalize()); 76 return 0; 77 } 78 79 /*TEST 80 81 test: 82 83 test: 84 nsize: 4 85 suffix: 2 86 87 TEST*/ 88