1 static char help[] = "Test PetscSFCompose() against some corner cases \n\n";
2
3 #include <petscsf.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 PetscMPIInt size;
8 PetscSF sfA0, sfA1, sfA2, sfB;
9 PetscInt nroots, nleaves;
10 PetscInt *ilocalA0, *ilocalA1, *ilocalA2, *ilocalB;
11 PetscSFNode *iremoteA0, *iremoteA1, *iremoteA2, *iremoteB;
12
13 PetscFunctionBeginUser;
14 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
16 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Only coded for one MPI process");
17 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA0));
18 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA1));
19 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA2));
20 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfB));
21 /* sfA0 */
22 nroots = 1;
23 nleaves = 0;
24 PetscCall(PetscMalloc1(nleaves, &ilocalA0));
25 PetscCall(PetscMalloc1(nleaves, &iremoteA0));
26 PetscCall(PetscSFSetGraph(sfA0, nroots, nleaves, ilocalA0, PETSC_OWN_POINTER, iremoteA0, PETSC_OWN_POINTER));
27 PetscCall(PetscSFSetUp(sfA0));
28 PetscCall(PetscObjectSetName((PetscObject)sfA0, "sfA0"));
29 PetscCall(PetscSFView(sfA0, NULL));
30 /* sfA1 */
31 nroots = 1;
32 nleaves = 1;
33 PetscCall(PetscMalloc1(nleaves, &ilocalA1));
34 PetscCall(PetscMalloc1(nleaves, &iremoteA1));
35 ilocalA1[0] = 1;
36 iremoteA1[0].rank = 0;
37 iremoteA1[0].index = 0;
38 PetscCall(PetscSFSetGraph(sfA1, nroots, nleaves, ilocalA1, PETSC_OWN_POINTER, iremoteA1, PETSC_OWN_POINTER));
39 PetscCall(PetscSFSetUp(sfA1));
40 PetscCall(PetscObjectSetName((PetscObject)sfA1, "sfA1"));
41 PetscCall(PetscSFView(sfA1, NULL));
42 /* sfA2 */
43 nroots = 1;
44 nleaves = 1;
45 PetscCall(PetscMalloc1(nleaves, &ilocalA2));
46 PetscCall(PetscMalloc1(nleaves, &iremoteA2));
47 ilocalA2[0] = 0;
48 iremoteA2[0].rank = 0;
49 iremoteA2[0].index = 0;
50 PetscCall(PetscSFSetGraph(sfA2, nroots, nleaves, ilocalA2, PETSC_OWN_POINTER, iremoteA2, PETSC_OWN_POINTER));
51 PetscCall(PetscSFSetUp(sfA2));
52 PetscCall(PetscObjectSetName((PetscObject)sfA2, "sfA2"));
53 PetscCall(PetscSFView(sfA2, NULL));
54 /* sfB */
55 nroots = 2;
56 nleaves = 2;
57 PetscCall(PetscMalloc1(nleaves, &ilocalB));
58 PetscCall(PetscMalloc1(nleaves, &iremoteB));
59 ilocalB[0] = 100;
60 iremoteB[0].rank = 0;
61 iremoteB[0].index = 0;
62 ilocalB[1] = 101;
63 iremoteB[1].rank = 0;
64 iremoteB[1].index = 1;
65 PetscCall(PetscSFSetGraph(sfB, nroots, nleaves, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER));
66 PetscCall(PetscSFSetUp(sfB));
67 PetscCall(PetscObjectSetName((PetscObject)sfB, "sfB"));
68 PetscCall(PetscSFView(sfB, NULL));
69 /* Test 0 */
70 {
71 PetscSF sfC;
72
73 PetscCall(PetscSFCompose(sfA0, sfB, &sfC));
74 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA0, sfB)"));
75 PetscCall(PetscSFView(sfC, NULL));
76 PetscCall(PetscSFDestroy(&sfC));
77 }
78 /* Test 1 */
79 {
80 PetscSF sfC;
81
82 PetscCall(PetscSFCompose(sfA1, sfB, &sfC));
83 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA1, sfB)"));
84 PetscCall(PetscSFView(sfC, NULL));
85 PetscCall(PetscSFDestroy(&sfC));
86 }
87 /* Test 2 */
88 {
89 PetscSF sfC;
90
91 PetscCall(PetscSFCompose(sfA2, sfB, &sfC));
92 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA2, sfB)"));
93 PetscCall(PetscSFView(sfC, NULL));
94 PetscCall(PetscSFDestroy(&sfC));
95 }
96 PetscCall(PetscSFDestroy(&sfA0));
97 PetscCall(PetscSFDestroy(&sfA1));
98 PetscCall(PetscSFDestroy(&sfA2));
99 PetscCall(PetscSFDestroy(&sfB));
100 PetscCall(PetscFinalize());
101 return 0;
102 }
103
104 /*TEST
105
106 test:
107 suffix: 0
108
109 TEST*/
110