Actual source code: ex32.c

  1: #include "petscda.h"

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

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

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

 30:   mx=7;
 31:   my=11;
 32:   mz=13;
 33:   ierr=DACreate3d(PETSC_COMM_WORLD,DA_NONPERIODIC,DA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,2,0,0,0,&Q2_da);
 34:   DASetUniformCoordinates(Q2_da,-1.0,1.0,-2.0,2.0,-3.0,3.0);
 35:   DAGetCoordinates(Q2_da,&coords);
 36:   DACreate3d(PETSC_COMM_WORLD,DA_NONPERIODIC,DA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,1,0,0,0,&Q1_da);
 37:   DASetCoordinates(Q1_da,coords);
 38:   VecDestroy(coords);

 40:   /* Get ghost coordinates one way */
 41:   DAGetGhostedCoordinates(Q1_da,&gcoords);

 43:   /* And another */
 44:   DAGetCoordinates(Q1_da,&coords);
 45:   DAGetCoordinateDA(Q1_da,&cda);
 46:   DAGetLocalVector(cda,&gcoords2);
 47:   DAGlobalToLocalBegin(cda,coords,INSERT_VALUES,gcoords2);
 48:   DAGlobalToLocalEnd(cda,coords,INSERT_VALUES,gcoords2);

 50:   CompareGhostedCoords(gcoords,gcoords2);
 51:   DARestoreLocalVector(cda,&gcoords2);
 52:   DADestroy(cda);

 54:   VecScale(coords,10.0);
 55:   VecScale(gcoords,10.0);
 56:   DAGetGhostedCoordinates(Q1_da,&gcoords2);
 57:   CompareGhostedCoords(gcoords,gcoords2);
 58:   VecDestroy(coords);
 59:   VecDestroy(gcoords2);

 61:   VecDestroy(gcoords);
 62:   DADestroy(Q2_da);
 63:   DADestroy(Q1_da);
 64:   return(0);
 65: }

 69: int main(int argc,char **argv)
 70: {

 73:   PetscInitialize(&argc,&argv,0,0);
 74:   TestQ2Q1DA();
 75:   PetscFinalize();
 76:   return 0;
 77: }