Actual source code: ex7.c
petsc-3.4.2 2013-07-02

      2: static char help[] = "Tests DMDALocalToLocalxxx().\n\n";

      4: #include <petscdmda.h>

      8: int main(int argc,char **argv)
      9: {
     10:   PetscMPIInt      rank;
     11:   PetscInt         M=8,dof=1,stencil_width=1,i,start,end,P=5,N = 6,m=PETSC_DECIDE,n=PETSC_DECIDE,p=PETSC_DECIDE,pt = 0,st = 0;
     12:   PetscErrorCode   ierr;
     13:   PetscBool        flg = PETSC_FALSE,flg2,flg3;
     14:   DMDABoundaryType periodic;
     15:   DMDAStencilType  stencil_type;
     16:   DM               da;
     17:   Vec              local,global,local_copy;
     18:   PetscScalar      value;
     19:   PetscReal        norm,work;
     20:   PetscViewer      viewer;
     21:   char             filename[64];
     22:   FILE             *file;

     24:   PetscInitialize(&argc,&argv,(char*)0,help);

     26:   PetscOptionsGetInt(NULL,"-M",&M,NULL);
     27:   PetscOptionsGetInt(NULL,"-N",&N,NULL);
     28:   PetscOptionsGetInt(NULL,"-dof",&dof,NULL);
     29:   PetscOptionsGetInt(NULL,"-stencil_width",&stencil_width,NULL);
     30:   PetscOptionsGetInt(NULL,"-periodic",&pt,NULL);

     32:   periodic = (DMDABoundaryType) pt;

     34:   PetscOptionsGetInt(NULL,"-stencil_type",&st,NULL);

     36:   stencil_type = (DMDAStencilType) st;

     38:   PetscOptionsHasName(NULL,"-2d",&flg2);
     39:   PetscOptionsHasName(NULL,"-3d",&flg3);
     40:   if (flg2) {
     41:     DMDACreate2d(PETSC_COMM_WORLD,periodic,periodic,stencil_type,M,N,m,n,dof,stencil_width,
     42:                         NULL,NULL,&da);
     43:   } else if (flg3) {
     44:     DMDACreate3d(PETSC_COMM_WORLD,periodic,periodic,periodic,stencil_type,M,N,P,m,n,p,dof,stencil_width,
     45:                         NULL,NULL,NULL,&da);
     46:   } else {
     47:     DMDACreate1d(PETSC_COMM_WORLD,periodic,M,dof,stencil_width,NULL,&da);
     48:   }

     50:   DMCreateGlobalVector(da,&global);
     51:   DMCreateLocalVector(da,&local);
     52:   VecDuplicate(local,&local_copy);

     55:   /* zero out vectors so that ghostpoints are zero */
     56:   value = 0;
     57:   VecSet(local,value);
     58:   VecSet(local_copy,value);

     60:   VecGetOwnershipRange(global,&start,&end);
     61:   for (i=start; i<end; i++) {
     62:     value = i + 1;
     63:     VecSetValues(global,1,&i,&value,INSERT_VALUES);
     64:   }
     65:   VecAssemblyBegin(global);
     66:   VecAssemblyEnd(global);

     68:   DMGlobalToLocalBegin(da,global,INSERT_VALUES,local);
     69:   DMGlobalToLocalEnd(da,global,INSERT_VALUES,local);

     72:   DMDALocalToLocalBegin(da,local,INSERT_VALUES,local_copy);
     73:   DMDALocalToLocalEnd(da,local,INSERT_VALUES,local_copy);

     76:   PetscOptionsGetBool(NULL,"-save",&flg,NULL);
     77:   if (flg) {
     78:     MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
     79:     sprintf(filename,"local.%d",rank);
     80:     PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&viewer);
     81:     PetscViewerASCIIGetPointer(viewer,&file);
     82:     VecView(local,viewer);
     83:     fprintf(file,"Vector with correct ghost points\n");
     84:     VecView(local_copy,viewer);
     85:     PetscViewerDestroy(&viewer);
     86:   }

     88:   VecAXPY(local_copy,-1.0,local);
     89:   VecNorm(local_copy,NORM_MAX,&work);
     90:   MPI_Allreduce(&work,&norm,1,MPIU_REAL,MPIU_MAX,PETSC_COMM_WORLD);
     91:   PetscPrintf(PETSC_COMM_WORLD,"Norm of difference %G should be zero\n",norm);

     93:   VecDestroy(&local_copy);
     94:   VecDestroy(&local);
     95:   VecDestroy(&global);
     96:   DMDestroy(&da);
     97:   PetscFinalize();
     98:   return 0;
     99: }