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