Actual source code: ex5.c
1: static char help[]= "Test PetscSFFCompose when the ilocal arrays are not identity nor dense\n\n";
3: #include <petsc.h>
4: #include <petscsf.h>
6: int main(int argc, char **argv)
7: {
8: PetscSF sfA, sfB, sfBA, sfAAm, sfBBm, sfAm, sfBm;
9: PetscInt nrootsA, nleavesA, nrootsB, nleavesB;
10: PetscInt *ilocalA, *ilocalB;
11: PetscSFNode *iremoteA, *iremoteB;
12: PetscMPIInt rank,size;
13: PetscInt i,m,n,k,nl = 2,mA,mB,nldataA,nldataB;
14: PetscInt *rdA,*rdB,*ldA,*ldB;
15: PetscBool inverse = PETSC_FALSE;
17: PetscInitialize(&argc,&argv,NULL,help);
18: PetscOptionsGetInt(NULL,NULL,"-nl",&nl,NULL);
19: PetscOptionsGetBool(NULL,NULL,"-explicit_inverse",&inverse,NULL);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
23: PetscSFCreate(PETSC_COMM_WORLD, &sfA);
24: PetscSFCreate(PETSC_COMM_WORLD, &sfB);
25: PetscSFSetFromOptions(sfA);
26: PetscSFSetFromOptions(sfB);
28: n = 4*nl*size;
29: m = 2*nl;
30: k = nl;
32: nldataA = rank == 0 ? n : 0;
33: nldataB = 3*nl;
35: nrootsA = m;
36: nleavesA = rank == 0 ? size*m : 0;
37: nrootsB = rank == 0 ? n : 0;
38: nleavesB = k;
40: PetscMalloc1(nleavesA, &ilocalA);
41: PetscMalloc1(nleavesA, &iremoteA);
42: PetscMalloc1(nleavesB, &ilocalB);
43: PetscMalloc1(nleavesB, &iremoteB);
45: /* sf A bcast is equivalent to a sparse gather on process 0
46: process 0 receives data in the middle [nl,3*nl] of the leaf data array for A */
47: for (i = 0; i < nleavesA; i++) {
48: iremoteA[i].rank = i/m;
49: iremoteA[i].index = i%m;
50: ilocalA[i] = nl + i/m * 4*nl + i%m;
51: }
53: /* sf B bcast is equivalent to a sparse scatter from process 0
54: process 0 sends data from [nl,2*nl] of the leaf data array for A
55: each process receives, in reverse order, in the middle [nl,2*nl] of the leaf data array for B */
56: for (i = 0; i < nleavesB; i++) {
57: iremoteB[i].rank = 0;
58: iremoteB[i].index = rank * 4*nl + nl + i%m;
59: ilocalB[i] = 2*nl - i - 1;
60: }
61: PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER);
62: PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);
63: PetscSFSetUp(sfA);
64: PetscSFSetUp(sfB);
65: PetscObjectSetName((PetscObject)sfA, "sfA");
66: PetscObjectSetName((PetscObject)sfB, "sfB");
67: PetscSFViewFromOptions(sfA, NULL, "-view");
68: PetscSFViewFromOptions(sfB, NULL, "-view");
70: PetscSFGetLeafRange(sfA, NULL, &mA);
71: PetscSFGetLeafRange(sfB, NULL, &mB);
72: PetscMalloc2(nrootsA, &rdA, nldataA, &ldA);
73: PetscMalloc2(nrootsB, &rdB, nldataB, &ldB);
74: for (i = 0; i < nrootsA; i++) rdA[i] = m*rank + i;
75: for (i = 0; i < nldataA; i++) ldA[i] = -1;
76: for (i = 0; i < nldataB; i++) ldB[i] = -1;
78: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastB(BcastA)\n");
79: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: root data\n");
80: PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD);
81: PetscSFBcastBegin(sfA, MPIU_INT, rdA, ldA,MPI_REPLACE);
82: PetscSFBcastEnd(sfA, MPIU_INT, rdA, ldA,MPI_REPLACE);
83: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: leaf data (all)\n");
84: PetscIntView(nldataA, ldA, PETSC_VIEWER_STDOUT_WORLD);
85: PetscSFBcastBegin(sfB, MPIU_INT, ldA, ldB,MPI_REPLACE);
86: PetscSFBcastEnd(sfB, MPIU_INT, ldA, ldB,MPI_REPLACE);
87: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "B: leaf data (all)\n");
88: PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD);
90: PetscSFCompose(sfA, sfB, &sfBA);
91: PetscSFSetFromOptions(sfBA);
92: PetscSFSetUp(sfBA);
93: PetscObjectSetName((PetscObject)sfBA, "sfBA");
94: PetscSFViewFromOptions(sfBA, NULL, "-view");
96: for (i = 0; i < nldataB; i++) ldB[i] = -1;
97: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastBA\n");
98: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: root data\n");
99: PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD);
100: PetscSFBcastBegin(sfBA, MPIU_INT, rdA, ldB,MPI_REPLACE);
101: PetscSFBcastEnd(sfBA, MPIU_INT, rdA, ldB,MPI_REPLACE);
102: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: leaf data (all)\n");
103: PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD);
105: PetscSFCreateInverseSF(sfA, &sfAm);
106: PetscSFSetFromOptions(sfAm);
107: PetscObjectSetName((PetscObject)sfAm, "sfAm");
108: PetscSFViewFromOptions(sfAm, NULL, "-view");
110: if (!inverse) {
111: PetscSFComposeInverse(sfA, sfA, &sfAAm);
112: } else {
113: PetscSFCompose(sfA, sfAm, &sfAAm);
114: }
115: PetscSFSetFromOptions(sfAAm);
116: PetscSFSetUp(sfAAm);
117: PetscObjectSetName((PetscObject)sfAAm, "sfAAm");
118: PetscSFViewFromOptions(sfAAm, NULL, "-view");
120: PetscSFCreateInverseSF(sfB, &sfBm);
121: PetscSFSetFromOptions(sfBm);
122: PetscObjectSetName((PetscObject)sfBm, "sfBm");
123: PetscSFViewFromOptions(sfBm, NULL, "-view");
125: if (!inverse) {
126: PetscSFComposeInverse(sfB, sfB, &sfBBm);
127: } else {
128: PetscSFCompose(sfB, sfBm, &sfBBm);
129: }
130: PetscSFSetFromOptions(sfBBm);
131: PetscSFSetUp(sfBBm);
132: PetscObjectSetName((PetscObject)sfBBm, "sfBBm");
133: PetscSFViewFromOptions(sfBBm, NULL, "-view");
135: PetscFree2(rdA, ldA);
136: PetscFree2(rdB, ldB);
138: PetscSFDestroy(&sfA);
139: PetscSFDestroy(&sfB);
140: PetscSFDestroy(&sfBA);
141: PetscSFDestroy(&sfAm);
142: PetscSFDestroy(&sfBm);
143: PetscSFDestroy(&sfAAm);
144: PetscSFDestroy(&sfBBm);
146: PetscFinalize();
147: return 0;
148: }
150: /*TEST
152: test:
153: suffix: 1
154: args: -view -explicit_inverse {{0 1}}
156: test:
157: nsize: 7
158: filter: grep -v "type" | grep -v "sort"
159: suffix: 2
160: args: -view -nl 5 -explicit_inverse {{0 1}}
162: # we cannot test for -sf_window_flavor dynamic because SFCompose with sparse leaves may change the root data pointer only locally, and this is not supported by the dynamic case
163: test:
164: nsize: 7
165: suffix: 2_window
166: filter: grep -v "type" | grep -v "sort"
167: output_file: output/ex5_2.out
168: args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor {{create allocate}}
169: requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW)
171: # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes
172: test:
173: nsize: 7
174: suffix: 2_window_shared
175: filter: grep -v "type" | grep -v "sort"
176: output_file: output/ex5_2.out
177: args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor shared
178: requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH_NUMVERSION) defined(PETSC_HAVE_MPI_ONE_SIDED)
180: TEST*/