Actual source code: ex1.c
petsc-3.3-p7 2013-05-11
1: static const char help[] = "Test star forest communication (PetscSF)\n\n";
3: /*T
4: Description: A star forest is a simple tree with one root and zero or more leaves.
5: Many common communication patterns can be expressed as updates of rootdata using leafdata and vice-versa.
6: This example creates a star forest, communicates values using the graph (see options for types of communication), views the graph, then destroys it.
7: T*/
9: /*
10: Include petscsf.h so we can use PetscSF objects. Note that this automatically
11: includes petscsys.h.
12: */
13: #include <petscsf.h>
17: int main(int argc,char **argv)
18: {
20: PetscInt i,nroots,nleaves;
21: PetscSFNode *remote;
22: PetscMPIInt rank,size;
23: PetscSF sf;
24: PetscBool test_bcast,test_reduce,test_degree,test_fetchandop,test_gather,test_scatter,test_embed,test_invert;
26: PetscInitialize(&argc,&argv,(char*)0,help);
27: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
28: MPI_Comm_size(PETSC_COMM_WORLD,&size);
30: PetscOptionsBegin(PETSC_COMM_WORLD,"","PetscSF Test Options","none");
31: test_bcast = PETSC_FALSE;
32: PetscOptionsBool("-test_bcast","Test broadcast","",test_bcast,&test_bcast,PETSC_NULL);
33: test_reduce = PETSC_FALSE;
34: PetscOptionsBool("-test_reduce","Test reduction","",test_reduce,&test_reduce,PETSC_NULL);
35: test_degree = PETSC_FALSE;
36: PetscOptionsBool("-test_degree","Test computation of vertex degree","",test_degree,&test_degree,PETSC_NULL);
37: test_fetchandop = PETSC_FALSE;
38: PetscOptionsBool("-test_fetchandop","Test atomic Fetch-And-Op","",test_fetchandop,&test_fetchandop,PETSC_NULL);
39: test_gather = PETSC_FALSE;
40: PetscOptionsBool("-test_gather","Test point gather","",test_gather,&test_gather,PETSC_NULL);
41: test_scatter = PETSC_FALSE;
42: PetscOptionsBool("-test_scatter","Test point scatter","",test_scatter,&test_scatter,PETSC_NULL);
43: test_embed = PETSC_FALSE;
44: PetscOptionsBool("-test_embed","Test point embed","",test_embed,&test_embed,PETSC_NULL);
45: test_invert = PETSC_FALSE;
46: PetscOptionsBool("-test_invert","Test point invert","",test_invert,&test_invert,PETSC_NULL);
47: PetscOptionsEnd();
49: nroots = 2 + (PetscInt)(rank == 0);
50: nleaves = 2 + (PetscInt)(rank > 0);
51: PetscMalloc(nleaves*sizeof(*remote),&remote);
52: /* Left periodic neighbor */
53: remote[0].rank = (rank+size-1)%size;
54: remote[0].index = 1;
55: /* Right periodic neighbor */
56: remote[1].rank = (rank+1)%size;
57: remote[1].index = 0;
58: if (rank > 0) { /* All processes reference rank 0, index 1 */
59: remote[2].rank = 0;
60: remote[2].index = 2;
61: }
63: /* Create a star forest for communication. In this example, the leaf space is dense, so we pass PETSC_NULL. */
64: PetscSFCreate(PETSC_COMM_WORLD,&sf);
65: PetscSFSetFromOptions(sf);
66: PetscSFSetGraph(sf,nroots,nleaves,PETSC_NULL,PETSC_COPY_VALUES,remote,PETSC_OWN_POINTER);
68: /* View graph, mostly useful for debugging purposes. */
69: PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD);
72: if (test_bcast) { /* broadcast rootdata into leafdata */
73: PetscInt *rootdata,*leafdata;
74: /* Allocate space for send and recieve buffers. This example communicates PetscInt, but other types, including
75: * user-defined structures, could also be used. */
76: PetscMalloc2(nroots,PetscInt,&rootdata,nleaves,PetscInt,&leafdata);
77: /* Set rootdata buffer to be broadcast */
78: for (i=0; i<nroots; i++) rootdata[i] = 100*(rank+1) + i;
79: /* Initialize local buffer, these values are never used. */
80: for (i=0; i<nleaves; i++) leafdata[i] = -1;
81: /* Broadcast entries from rootdata to leafdata. Computation or other communication can be performed between the begin and end calls. */
82: PetscSFBcastBegin(sf,MPIU_INT,rootdata,leafdata);
83: PetscSFBcastEnd(sf,MPIU_INT,rootdata,leafdata);
84: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Bcast Rootdata\n");
85: PetscIntView(nroots,rootdata,PETSC_VIEWER_STDOUT_WORLD);
86: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Bcast Leafdata\n");
87: PetscIntView(nleaves,leafdata,PETSC_VIEWER_STDOUT_WORLD);
88: PetscFree2(rootdata,leafdata);
89: }
91: if (test_reduce) { /* Reduce leafdata into rootdata */
92: PetscInt *rootdata,*leafdata;
93: PetscMalloc2(nroots,PetscInt,&rootdata,nleaves,PetscInt,&leafdata);
94: /* Initialize rootdata buffer in which the result of the reduction will appear. */
95: for (i=0; i<nroots; i++) rootdata[i] = 100*(rank+1) + i;
96: /* Set leaf values to reduce. */
97: for (i=0; i<nleaves; i++) leafdata[i] = 1000*(rank+1) + 10*i;
98: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Pre-Reduce Rootdata\n");
99: PetscIntView(nroots,rootdata,PETSC_VIEWER_STDOUT_WORLD);
100: /* Perform reduction. Computation or other communication can be performed between the begin and end calls.
101: * This example sums the values, but other MPI_Ops can be used (e.g MPI_MAX, MPI_PROD). */
102: PetscSFReduceBegin(sf,MPIU_INT,leafdata,rootdata,MPIU_SUM);
103: PetscSFReduceEnd(sf,MPIU_INT,leafdata,rootdata,MPIU_SUM);
104: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Reduce Leafdata\n");
105: PetscIntView(nleaves,leafdata,PETSC_VIEWER_STDOUT_WORLD);
106: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Reduce Rootdata\n");
107: PetscIntView(nroots,rootdata,PETSC_VIEWER_STDOUT_WORLD);
108: PetscFree2(rootdata,leafdata);
109: }
111: if (test_degree) {
112: const PetscInt *degree;
113: PetscSFComputeDegreeBegin(sf,°ree);
114: PetscSFComputeDegreeEnd(sf,°ree);
115: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Root degrees\n");
116: PetscIntView(nroots,degree,PETSC_VIEWER_STDOUT_WORLD);
117: }
119: if (test_fetchandop) {
120: /* Cannot use text compare here because token ordering is not deterministic */
121: PetscInt *leafdata,*leafupdate,*rootdata;
122: PetscMalloc3(nleaves,PetscInt,&leafdata,nleaves,PetscInt,&leafupdate,nroots,PetscInt,&rootdata);
123: for (i=0; i<nleaves; i++) leafdata[i] = 1;
124: for (i=0; i<nroots; i++) rootdata[i] = 0;
125: PetscSFFetchAndOpBegin(sf,MPIU_INT,rootdata,leafdata,leafupdate,MPIU_SUM);
126: PetscSFFetchAndOpEnd(sf,MPIU_INT,rootdata,leafdata,leafupdate,MPIU_SUM);
127: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Rootdata (sum of 1 from each leaf)\n");
128: PetscIntView(nroots,rootdata,PETSC_VIEWER_STDOUT_WORLD);
129: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Leafupdate (value at roots prior to my atomic update)\n");
130: PetscIntView(nleaves,leafupdate,PETSC_VIEWER_STDOUT_WORLD);
131: PetscFree3(leafdata,leafupdate,rootdata);
132: }
134: if (test_gather) {
135: const PetscInt *degree;
136: PetscInt inedges,*indata,*outdata;
137: PetscSFComputeDegreeBegin(sf,°ree);
138: PetscSFComputeDegreeEnd(sf,°ree);
139: for (i=0,inedges=0; i<nroots; i++) inedges += degree[i];
140: PetscMalloc2(inedges,PetscInt,&indata,nleaves,PetscInt,&outdata);
141: for (i=0; i<nleaves; i++) outdata[i] = 1000*(rank+1) + i;
142: PetscSFGatherBegin(sf,MPIU_INT,outdata,indata);
143: PetscSFGatherEnd(sf,MPIU_INT,outdata,indata);
144: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Gathered data at multi-roots from leaves\n");
145: PetscIntView(inedges,indata,PETSC_VIEWER_STDOUT_WORLD);
146: PetscFree2(indata,outdata);
147: }
149: if (test_scatter) {
150: const PetscInt *degree;
151: PetscInt j,count,inedges,*indata,*outdata;
152: PetscSFComputeDegreeBegin(sf,°ree);
153: PetscSFComputeDegreeEnd(sf,°ree);
154: for (i=0,inedges=0; i<nroots; i++) inedges += degree[i];
155: PetscMalloc2(inedges,PetscInt,&indata,nleaves,PetscInt,&outdata);
156: for (i=0,count=0; i<nroots; i++) {
157: for (j=0; j<degree[i]; j++) indata[count++] = 1000*(rank+1) + 100*i + j;
158: }
159: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Data at multi-roots, to scatter to leaves\n");
160: PetscIntView(inedges,indata,PETSC_VIEWER_STDOUT_WORLD);
162: PetscSFScatterBegin(sf,MPIU_INT,indata,outdata);
163: PetscSFScatterEnd(sf,MPIU_INT,indata,outdata);
164: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Scattered data at leaves\n");
165: PetscIntView(nleaves,outdata,PETSC_VIEWER_STDOUT_WORLD);
166: PetscFree2(indata,outdata);
167: }
169: if (test_embed) {
170: const PetscInt nroots = 1 + (PetscInt)!rank,selected[] = {1,2};
171: PetscSF esf;
172: PetscSFCreateEmbeddedSF(sf,nroots,selected,&esf);
173: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Embedded PetscSF\n");
174: PetscSFView(esf,PETSC_VIEWER_STDOUT_WORLD);
175: PetscSFDestroy(&esf);
176: }
178: if (test_invert) {
179: PetscSF msf,imsf;
180: PetscSFGetMultiSF(sf,&msf);
181: PetscSFCreateInverseSF(msf,&imsf);
182: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Multi-SF\n");
183: PetscSFView(msf,PETSC_VIEWER_STDOUT_WORLD);
184: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Inverse of Multi-SF\n");
185: PetscSFView(imsf,PETSC_VIEWER_STDOUT_WORLD);
186: PetscSFDestroy(&imsf);
187: }
189: /* Clean storage for star forest. */
190: PetscSFDestroy(&sf);
191: PetscFinalize();
192: return 0;
193: }