Actual source code: ex32.c
petsc-3.3-p6 2013-02-11
1: #include <petscdmda.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: MPI_Allreduce(&nrm,&gnrm,1,MPIU_REAL,MPIU_MAX,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: DM 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=DMDACreate3d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,2,0,0,0,&Q2_da);
34: DMDASetUniformCoordinates(Q2_da,-1.0,1.0,-2.0,2.0,-3.0,3.0);
35: DMDAGetCoordinates(Q2_da,&coords);
36: DMDACreate3d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_STENCIL_BOX,mx,my,mz,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,3,1,0,0,0,&Q1_da);
37: DMDASetCoordinates(Q1_da,coords);
39: /* Get ghost coordinates one way */
40: DMDAGetGhostedCoordinates(Q1_da,&gcoords);
42: /* And another */
43: DMDAGetCoordinates(Q1_da,&coords);
44: DMDAGetCoordinateDA(Q1_da,&cda);
45: DMGetLocalVector(cda,&gcoords2);
46: DMGlobalToLocalBegin(cda,coords,INSERT_VALUES,gcoords2);
47: DMGlobalToLocalEnd(cda,coords,INSERT_VALUES,gcoords2);
49: CompareGhostedCoords(gcoords,gcoords2);
50: DMRestoreLocalVector(cda,&gcoords2);
52: VecScale(coords,10.0);
53: VecScale(gcoords,10.0);
54: DMDAGetGhostedCoordinates(Q1_da,&gcoords2);
55: CompareGhostedCoords(gcoords,gcoords2);
57: DMDestroy(&Q2_da);
58: DMDestroy(&Q1_da);
59: return(0);
60: }
64: int main(int argc,char **argv)
65: {
68: PetscInitialize(&argc,&argv,0,0);
69: TestQ2Q1DA();
70: PetscFinalize();
71: return 0;
72: }