Actual source code: ex32.c

  1: static char help[] = "Tests DMDA ghost coordinates\n\n";

  3: #include <petscdm.h>
  4: #include <petscdmda.h>

  6: static PetscErrorCode CompareGhostedCoords(Vec gc1,Vec gc2)
  7: {
  9:   PetscReal      nrm,gnrm;
 10:   Vec            tmp;

 13:   VecDuplicate(gc1,&tmp);
 14:   VecWAXPY(tmp,-1.0,gc1,gc2);
 15:   VecNorm(tmp,NORM_INFINITY,&nrm);
 16:   MPI_Allreduce(&nrm,&gnrm,1,MPIU_REAL,MPIU_MAX,PETSC_COMM_WORLD);
 17:   PetscPrintf(PETSC_COMM_WORLD,"norm of difference of ghosted coordinates %8.2e\n",gnrm);
 18:   VecDestroy(&tmp);
 19:   return(0);
 20: }

 22: static PetscErrorCode TestQ2Q1DA(void)
 23: {
 24:   DM             Q2_da,Q1_da,cda;
 25:   PetscInt       mx,my,mz;
 26:   Vec            coords,gcoords,gcoords2;

 30:   mx   = 7;
 31:   my   = 11;
 32:   mz   = 13;
 33:   DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,2,0,0,0,&Q2_da);
 34:   DMSetFromOptions(Q2_da);
 35:   DMSetUp(Q2_da);
 36:   DMDASetUniformCoordinates(Q2_da,-1.0,1.0,-2.0,2.0,-3.0,3.0);
 37:   DMGetCoordinates(Q2_da,&coords);
 38:   DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,1,0,0,0,&Q1_da);
 39:   DMSetFromOptions(Q1_da);
 40:   DMSetUp(Q1_da);
 41:   DMSetCoordinates(Q1_da,coords);

 43:   /* Get ghost coordinates one way */
 44:   DMGetCoordinatesLocal(Q1_da,&gcoords);

 46:   /* And another */
 47:   DMGetCoordinates(Q1_da,&coords);
 48:   DMGetCoordinateDM(Q1_da,&cda);
 49:   DMGetLocalVector(cda,&gcoords2);
 50:   DMGlobalToLocalBegin(cda,coords,INSERT_VALUES,gcoords2);
 51:   DMGlobalToLocalEnd(cda,coords,INSERT_VALUES,gcoords2);

 53:   CompareGhostedCoords(gcoords,gcoords2);
 54:   DMRestoreLocalVector(cda,&gcoords2);

 56:   VecScale(coords,10.0);
 57:   VecScale(gcoords,10.0);
 58:   DMGetCoordinatesLocal(Q1_da,&gcoords2);
 59:   CompareGhostedCoords(gcoords,gcoords2);

 61:   DMDestroy(&Q2_da);
 62:   DMDestroy(&Q1_da);
 63:   return(0);
 64: }

 66: int main(int argc,char **argv)
 67: {

 70:   PetscInitialize(&argc,&argv,0,help);if (ierr) return ierr;
 71:   TestQ2Q1DA();
 72:   PetscFinalize();
 73:   return ierr;
 74: }

 76: /*TEST

 78:    test:
 79:       nsize: 2

 81: TEST*/