1 static char help[] = "Spot check DMStag Compatibility Checks";
2
3 #include <petscdm.h>
4 #include <petscdmstag.h>
5
6 #define NDMS 4
7
main(int argc,char ** argv)8 int main(int argc, char **argv)
9 {
10 DM dms[NDMS];
11 PetscInt i;
12
13 PetscFunctionBeginUser;
14 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15
16 /* Two 3d DMs, with all the same parameters */
17 PetscCall(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 PetscCall(DMSetUp(dms[0]));
19 PetscCall(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 PetscCall(DMSetUp(dms[1]));
21
22 /* A derived 3d DM, with a different section */
23 PetscCall(DMStagCreateCompatibleDMStag(dms[0], 0, 1, 0, 1, &dms[2]));
24
25 /* A DM expected to be incompatible (different stencil width) */
26 PetscCall(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]));
27
28 /* Check expected self-compatibility */
29 for (i = 0; i < NDMS; ++i) {
30 PetscBool compatible, set;
31 PetscCall(DMGetCompatibility(dms[i], dms[i], &compatible, &set));
32 PetscCheck(set && compatible, PetscObjectComm((PetscObject)dms[i]), PETSC_ERR_PLIB, "DM %" PetscInt_FMT " not determined compatible with itself", i);
33 }
34
35 /* Check expected compatibility */
36 for (i = 1; i <= 2; ++i) {
37 PetscBool compatible, set;
38 PetscCall(DMGetCompatibility(dms[0], dms[i], &compatible, &set));
39 PetscCheck(set && compatible, PetscObjectComm((PetscObject)dms[i]), PETSC_ERR_PLIB, "DM %" PetscInt_FMT " not determined compatible with DM %d", i, 0);
40 }
41
42 /* Check expected incompatibility */
43 {
44 PetscBool compatible, set;
45 PetscCall(DMGetCompatibility(dms[0], dms[3], &compatible, &set));
46 PetscCheck(set && !compatible, PetscObjectComm((PetscObject)dms[i]), PETSC_ERR_PLIB, "DM %" PetscInt_FMT " not determined incompatible with DM %d", i, 0);
47 }
48
49 for (i = 0; i < NDMS; ++i) PetscCall(DMDestroy(&dms[i]));
50 PetscCall(PetscFinalize());
51 return 0;
52 }
53
54 /*TEST
55
56 test:
57 nsize: {{1 3}}
58 suffix: 1
59 output_file: output/empty.out
60
61 TEST*/
62