Actual source code: ex3.c
petsc-3.11.4 2019-09-28
1: static char help[] = "Spot check DMStag Compatibility Checks";
3: #include <petscdm.h>
4: #include <petscdmstag.h>
6: #define NDMS 4
8: int main(int argc,char **argv)
9: {
11: DM dms[NDMS];
12: PetscInt i;
14: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16: /* Two 3d DMs, with all the same parameters */
17: DMStagCreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,4,3,2,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,2,3,4,5,DMSTAG_STENCIL_BOX,1,NULL,NULL,NULL,&dms[0]);
18: DMSetUp(dms[0]);
19: DMStagCreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,4,3,2,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,2,3,4,5,DMSTAG_STENCIL_BOX,1,NULL,NULL,NULL,&dms[1]);
20: DMSetUp(dms[1]);
22: /* A derived 3d DM, with a different section */
23: DMStagCreateCompatibleDMStag(dms[0],0,1,0,1,&dms[2]);
25: /* A DM expected to be incompatible (different stencil width) */
26: DMStagCreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,4,3,2,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,2,3,4,5,DMSTAG_STENCIL_BOX,2,NULL,NULL,NULL,&dms[3]);
28: /* Check expected self-compatibility */
29: for (i=0; i<NDMS; ++i) {
30: PetscBool compatible,set;
31: DMGetCompatibility(dms[i],dms[i],&compatible,&set);
32: if (!set || !compatible) SETERRQ1(PetscObjectComm((PetscObject)dms[i]),PETSC_ERR_PLIB,"DM %D not determined compatible with itself",i);
33: }
35: /* Check expected compatibility */
36: for (i=1; i<=2; ++i) {
37: PetscBool compatible,set;
38: DMGetCompatibility(dms[0],dms[i],&compatible,&set);
39: if (!set || !compatible) SETERRQ2(PetscObjectComm((PetscObject)dms[i]),PETSC_ERR_PLIB,"DM %D not determined compatible with DM %D",i,0);
40: }
42: /* Check expected incompatibility */
43: {
44: PetscBool compatible,set;
45: DMGetCompatibility(dms[0],dms[3],&compatible,&set);
46: if (!set || compatible) SETERRQ2(PetscObjectComm((PetscObject)dms[i]),PETSC_ERR_PLIB,"DM %D not determined incompatible with DM %D",i,0);
47: }
49: for (i=0; i<NDMS; ++i) {
50: DMDestroy(&dms[i]);
51: }
52: PetscFinalize();
53: return ierr;
54: }
56: /*TEST
58: test:
59: nsize: 1
60: suffix: 1
62: test:
63: nsize: 3
64: suffix: 2
66: TEST*/