Actual source code: sfbasic.c
petsc-3.4.2 2013-07-02
1: #define PETSC_DESIRE_COMPLEX
2: #include <petsc-private/sfimpl.h> /*I "petscsf.h" I*/
4: typedef struct _n_PetscSFBasicPack *PetscSFBasicPack;
5: struct _n_PetscSFBasicPack {
6: void (*Pack)(PetscInt,const PetscInt*,const void*,void*);
7: void (*UnpackInsert)(PetscInt,const PetscInt*,void*,const void*);
8: void (*UnpackAdd)(PetscInt,const PetscInt*,void*,const void*);
9: void (*UnpackMin)(PetscInt,const PetscInt*,void*,const void*);
10: void (*UnpackMax)(PetscInt,const PetscInt*,void*,const void*);
11: void (*UnpackMinloc)(PetscInt,const PetscInt*,void*,const void*);
12: void (*UnpackMaxloc)(PetscInt,const PetscInt*,void*,const void*);
13: void (*FetchAndInsert)(PetscInt,const PetscInt*,void*,void*);
14: void (*FetchAndAdd)(PetscInt,const PetscInt*,void*,void*);
15: void (*FetchAndMin)(PetscInt,const PetscInt*,void*,void*);
16: void (*FetchAndMax)(PetscInt,const PetscInt*,void*,void*);
17: void (*FetchAndMinloc)(PetscInt,const PetscInt*,void*,void*);
18: void (*FetchAndMaxloc)(PetscInt,const PetscInt*,void*,void*);
20: MPI_Datatype unit;
21: size_t unitbytes; /* Number of bytes in a unit */
22: const void *key; /* Array used as key for operation */
23: char *root; /* Packed root data, contiguous by leaf rank */
24: char *leaf; /* Packed leaf data, contiguous by root rank */
25: MPI_Request *requests; /* Array of root requests followed by leaf requests */
26: PetscSFBasicPack next;
27: };
29: typedef struct {
30: PetscMPIInt tag;
31: PetscInt niranks; /* Number of incoming ranks (ranks accessing my roots) */
32: PetscMPIInt *iranks; /* Array of ranks that reference my roots */
33: PetscInt itotal; /* Total number of graph edges referencing my roots */
34: PetscInt *ioffset; /* Array of length niranks+1 holding offset in irootloc[] for each rank */
35: PetscInt *irootloc; /* Incoming roots referenced by ranks starting at ioffset[rank] */
36: PetscSFBasicPack avail; /* One or more entries per MPI Datatype, lazily constructed */
37: PetscSFBasicPack inuse; /* Buffers being used for transactions that have not yet completed */
38: } PetscSF_Basic;
40: #if !defined(PETSC_HAVE_MPI_TYPE_DUP) /* Danger: type is not reference counted; subject to ABA problem */
41: PETSC_STATIC_INLINE PetscErrorCode MPI_Type_dup(MPI_Datatype datatype,MPI_Datatype *newtype)
42: {
43: *newtype = datatype;
44: return 0;
45: }
46: #endif
48: /*
49: * MPI_Reduce_local is not really useful because it can't handle sparse data and it vectorizes "in the wrong direction",
50: * therefore we pack data types manually. This section defines packing routines for the standard data types.
51: */
53: #define CPPJoin2_exp(a,b) a ## b
54: #define CPPJoin2(a,b) CPPJoin2_exp(a,b)
55: #define CPPJoin3_exp_(a,b,c) a ## b ## _ ## c
56: #define CPPJoin3_(a,b,c) CPPJoin3_exp_(a,b,c)
58: /* Basic types without addition */
59: #define DEF_PackNoInit(type) \
60: static void CPPJoin2(Pack_,type)(PetscInt n,const PetscInt *idx,const void *unpacked,void *packed) { \
61: const type *u = (const type*)unpacked; \
62: type *p = (type*)packed; \
63: PetscInt i; \
64: for (i=0; i<n; i++) p[i] = u[idx[i]]; \
65: } \
66: static void CPPJoin2(UnpackInsert_,type)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
67: type *u = (type*)unpacked; \
68: const type *p = (const type*)packed; \
69: PetscInt i; \
70: for (i=0; i<n; i++) u[idx[i]] = p[i]; \
71: } \
72: static void CPPJoin2(FetchAndInsert_,type)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
73: type *u = (type*)unpacked; \
74: type *p = (type*)packed; \
75: PetscInt i; \
76: for (i=0; i<n; i++) { \
77: PetscInt j = idx[i]; \
78: type t = u[j]; \
79: u[j] = p[i]; \
80: p[i] = t; \
81: } \
82: }
84: /* Basic types defining addition */
85: #define DEF_PackAddNoInit(type) \
86: DEF_PackNoInit(type) \
87: static void CPPJoin2(UnpackAdd_,type)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
88: type *u = (type*)unpacked; \
89: const type *p = (const type*)packed; \
90: PetscInt i; \
91: for (i=0; i<n; i++) u[idx[i]] += p[i]; \
92: } \
93: static void CPPJoin2(FetchAndAdd_,type)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
94: type *u = (type*)unpacked; \
95: type *p = (type*)packed; \
96: PetscInt i; \
97: for (i=0; i<n; i++) { \
98: PetscInt j = idx[i]; \
99: type t = u[j]; \
100: u[j] = t + p[i]; \
101: p[i] = t; \
102: } \
103: }
104: #define DEF_Pack(type) \
105: DEF_PackAddNoInit(type) \
106: static void CPPJoin2(PackInit_,type)(PetscSFBasicPack link) { \
107: link->Pack = CPPJoin2(Pack_,type); \
108: link->UnpackInsert = CPPJoin2(UnpackInsert_,type); \
109: link->UnpackAdd = CPPJoin2(UnpackAdd_,type); \
110: link->FetchAndInsert = CPPJoin2(FetchAndInsert_,type); \
111: link->FetchAndAdd = CPPJoin2(FetchAndAdd_,type); \
112: link->unitbytes = sizeof(type); \
113: }
114: /* Comparable types */
115: #define DEF_PackCmp(type) \
116: DEF_PackAddNoInit(type) \
117: static void CPPJoin2(UnpackMax_,type)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
118: type *u = (type*)unpacked; \
119: const type *p = (const type*)packed; \
120: PetscInt i; \
121: for (i=0; i<n; i++) { \
122: type v = u[idx[i]]; \
123: u[idx[i]] = PetscMax(v,p[i]); \
124: } \
125: } \
126: static void CPPJoin2(UnpackMin_,type)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
127: type *u = (type*)unpacked; \
128: const type *p = (const type*)packed; \
129: PetscInt i; \
130: for (i=0; i<n; i++) { \
131: type v = u[idx[i]]; \
132: u[idx[i]] = PetscMin(v,p[i]); \
133: } \
134: } \
135: static void CPPJoin2(FetchAndMax_,type)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
136: type *u = (type*)unpacked; \
137: type *p = (type*)packed; \
138: PetscInt i; \
139: for (i=0; i<n; i++) { \
140: PetscInt j = idx[i]; \
141: type v = u[j]; \
142: u[j] = PetscMax(v,p[i]); \
143: p[i] = v; \
144: } \
145: } \
146: static void CPPJoin2(FetchAndMin_,type)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
147: type *u = (type*)unpacked; \
148: type *p = (type*)packed; \
149: PetscInt i; \
150: for (i=0; i<n; i++) { \
151: PetscInt j = idx[i]; \
152: type v = u[j]; \
153: u[j] = PetscMin(v,p[i]); \
154: p[i] = v; \
155: } \
156: } \
157: static void CPPJoin2(PackInit_,type)(PetscSFBasicPack link) { \
158: link->Pack = CPPJoin2(Pack_,type); \
159: link->UnpackInsert = CPPJoin2(UnpackInsert_,type); \
160: link->UnpackAdd = CPPJoin2(UnpackAdd_,type); \
161: link->UnpackMax = CPPJoin2(UnpackMax_,type); \
162: link->UnpackMin = CPPJoin2(UnpackMin_,type); \
163: link->FetchAndInsert = CPPJoin2(FetchAndInsert_,type); \
164: link->FetchAndAdd = CPPJoin2(FetchAndAdd_ ,type); \
165: link->FetchAndMax = CPPJoin2(FetchAndMax_ ,type); \
166: link->FetchAndMin = CPPJoin2(FetchAndMin_ ,type); \
167: link->unitbytes = sizeof(type); \
168: }
170: /* Pair types */
171: #define CPPJoinloc_exp(base,op,t1,t2) base ## op ## loc_ ## t1 ## _ ## t2
172: #define CPPJoinloc(base,op,t1,t2) CPPJoinloc_exp(base,op,t1,t2)
173: #define PairType(type1,type2) CPPJoin3_(_pairtype_,type1,type2)
174: #define DEF_UnpackXloc(type1,type2,locname,op) \
175: static void CPPJoinloc(Unpack,locname,type1,type2)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
176: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
177: const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
178: PetscInt i; \
179: for (i=0; i<n; i++) { \
180: PetscInt j = idx[i]; \
181: if (p[i].a op u[j].a) { \
182: u[j].a = p[i].a; \
183: u[j].b = p[i].b; \
184: } else if (u[j].a == p[i].a) { \
185: u[j].b = PetscMin(u[j].b,p[i].b); \
186: } \
187: } \
188: } \
189: static void CPPJoinloc(FetchAnd,locname,type1,type2)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
190: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
191: PairType(type1,type2) *p = (PairType(type1,type2)*)packed; \
192: PetscInt i; \
193: for (i=0; i<n; i++) { \
194: PetscInt j = idx[i]; \
195: PairType(type1,type2) v; \
196: v.a = u[j].a; \
197: v.b = u[j].b; \
198: if (p[i].a op u[j].a) { \
199: u[j].a = p[i].a; \
200: u[j].b = p[i].b; \
201: } else if (u[j].a == p[i].a) { \
202: u[j].b = PetscMin(u[j].b,p[i].b); \
203: } \
204: p[i].a = v.a; \
205: p[i].b = v.b; \
206: } \
207: }
208: #define DEF_PackPair(type1,type2) \
209: typedef struct {type1 a; type2 b;} PairType(type1,type2); \
210: static void CPPJoin3_(Pack_,type1,type2)(PetscInt n,const PetscInt *idx,const void *unpacked,void *packed) { \
211: const PairType(type1,type2) *u = (const PairType(type1,type2)*)unpacked; \
212: PairType(type1,type2) *p = (PairType(type1,type2)*)packed; \
213: PetscInt i; \
214: for (i=0; i<n; i++) { \
215: p[i].a = u[idx[i]].a; \
216: p[i].b = u[idx[i]].b; \
217: } \
218: } \
219: static void CPPJoin3_(UnpackInsert_,type1,type2)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
220: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
221: const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
222: PetscInt i; \
223: for (i=0; i<n; i++) { \
224: u[idx[i]].a = p[i].a; \
225: u[idx[i]].b = p[i].b; \
226: } \
227: } \
228: static void CPPJoin3_(UnpackAdd_,type1,type2)(PetscInt n,const PetscInt *idx,void *unpacked,const void *packed) { \
229: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
230: const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
231: PetscInt i; \
232: for (i=0; i<n; i++) { \
233: u[idx[i]].a += p[i].a; \
234: u[idx[i]].b += p[i].b; \
235: } \
236: } \
237: static void CPPJoin3_(FetchAndInsert_,type1,type2)(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
238: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
239: PairType(type1,type2) *p = (PairType(type1,type2)*)packed; \
240: PetscInt i; \
241: for (i=0; i<n; i++) { \
242: PetscInt j = idx[i]; \
243: PairType(type1,type2) v; \
244: v.a = u[j].a; \
245: v.b = u[j].b; \
246: u[j].a = p[i].a; \
247: u[j].b = p[i].b; \
248: p[i].a = v.a; \
249: p[i].b = v.b; \
250: } \
251: } \
252: static void FetchAndAdd_ ## type1 ## _ ## type2(PetscInt n,const PetscInt *idx,void *unpacked,void *packed) { \
253: PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked; \
254: PairType(type1,type2) *p = (PairType(type1,type2)*)packed; \
255: PetscInt i; \
256: for (i=0; i<n; i++) { \
257: PetscInt j = idx[i]; \
258: PairType(type1,type2) v; \
259: v.a = u[j].a; \
260: v.b = u[j].b; \
261: u[j].a = v.a + p[i].a; \
262: u[j].b = v.b + p[i].b; \
263: p[i].a = v.a; \
264: p[i].b = v.b; \
265: } \
266: } \
267: DEF_UnpackXloc(type1,type2,Max,>) \
268: DEF_UnpackXloc(type1,type2,Min,<) \
269: static void CPPJoin3_(PackInit_,type1,type2)(PetscSFBasicPack link) { \
270: link->Pack = CPPJoin3_(Pack_,type1,type2); \
271: link->UnpackInsert = CPPJoin3_(UnpackInsert_,type1,type2); \
272: link->UnpackAdd = CPPJoin3_(UnpackAdd_,type1,type2); \
273: link->UnpackMaxloc = CPPJoin3_(UnpackMaxloc_,type1,type2); \
274: link->UnpackMinloc = CPPJoin3_(UnpackMinloc_,type1,type2); \
275: link->FetchAndInsert = CPPJoin3_(FetchAndInsert_,type1,type2); \
276: link->FetchAndAdd = CPPJoin3_(FetchAndAdd_,type1,type2); \
277: link->FetchAndMaxloc = CPPJoin3_(FetchAndMaxloc_,type1,type2); \
278: link->FetchAndMinloc = CPPJoin3_(FetchAndMinloc_,type1,type2); \
279: link->unitbytes = sizeof(PairType(type1,type2)); \
280: }
282: /* Currently only dumb blocks of data */
283: #define BlockType(unit,count) CPPJoin3_(_blocktype_,unit,count)
284: #define DEF_Block(unit,count) \
285: typedef struct {unit v[count];} BlockType(unit,count); \
286: DEF_PackNoInit(BlockType(unit,count)) \
287: static void CPPJoin3_(PackInit_block_,unit,count)(PetscSFBasicPack link) { \
288: link->Pack = CPPJoin2(Pack_,BlockType(unit,count)); \
289: link->UnpackInsert = CPPJoin2(UnpackInsert_,BlockType(unit,count)); \
290: link->FetchAndInsert = CPPJoin2(FetchAndInsert_,BlockType(unit,count)); \
291: link->unitbytes = sizeof(BlockType(unit,count)); \
292: }
294: DEF_PackCmp(int)
295: DEF_PackCmp(PetscInt)
296: DEF_PackCmp(PetscReal)
297: #if defined(PETSC_HAVE_COMPLEX)
298: DEF_Pack(PetscComplex)
299: #endif
300: DEF_PackPair(int,int)
301: DEF_PackPair(PetscInt,PetscInt)
302: DEF_Block(int,1)
303: DEF_Block(int,2)
304: DEF_Block(int,3)
305: DEF_Block(int,4)
306: DEF_Block(int,5)
307: DEF_Block(int,6)
308: DEF_Block(int,7)
309: DEF_Block(int,8)
313: static PetscErrorCode PetscSFSetUp_Basic(PetscSF sf)
314: {
315: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
317: PetscInt *rlengths,*ilengths,i;
318: MPI_Comm comm;
319: MPI_Request *rootreqs,*leafreqs;
322: PetscObjectGetComm((PetscObject)sf,&comm);
323: PetscObjectGetNewTag((PetscObject)sf,&bas->tag);
324: /*
325: * Inform roots about how many leaves and from which ranks
326: */
327: PetscMalloc(sf->nranks*sizeof(PetscInt),&rlengths);
328: /* Determine number, sending ranks, and length of incoming */
329: for (i=0; i<sf->nranks; i++) {
330: rlengths[i] = sf->roffset[i+1] - sf->roffset[i]; /* Number of roots referenced by my leaves; for rank sf->ranks[i] */
331: }
332: PetscCommBuildTwoSided(comm,1,MPIU_INT,sf->nranks,sf->ranks,rlengths,&bas->niranks,&bas->iranks,(void**)&ilengths);
333: PetscFree(rlengths);
335: /* Send leaf identities to roots */
336: for (i=0,bas->itotal=0; i<bas->niranks; i++) bas->itotal += ilengths[i];
337: PetscMalloc2(bas->niranks+1,PetscInt,&bas->ioffset,bas->itotal,PetscInt,&bas->irootloc);
338: PetscMalloc((bas->niranks+sf->nranks)*sizeof(MPI_Request),&rootreqs);
340: leafreqs = rootreqs + bas->niranks;
341: bas->ioffset[0] = 0;
342: for (i=0; i<bas->niranks; i++) {
343: bas->ioffset[i+1] = bas->ioffset[i] + ilengths[i];
344: MPI_Irecv(bas->irootloc+bas->ioffset[i],ilengths[i],MPIU_INT,bas->iranks[i],bas->tag,comm,&rootreqs[i]);
345: }
346: for (i=0; i<sf->nranks; i++) {
347: PetscMPIInt npoints;
348: PetscMPIIntCast(sf->roffset[i+1] - sf->roffset[i],&npoints);
349: MPI_Isend(sf->rremote+sf->roffset[i],npoints,MPIU_INT,sf->ranks[i],bas->tag,comm,&leafreqs[i]);
350: }
351: MPI_Waitall(sf->nranks+bas->niranks,rootreqs,MPI_STATUSES_IGNORE);
352: PetscFree(ilengths);
353: PetscFree(rootreqs);
354: return(0);
355: }
359: static PetscErrorCode PetscSFBasicPackTypeSetup(PetscSFBasicPack link,MPI_Datatype unit)
360: {
362: PetscBool isInt,isPetscInt,isPetscReal,is2Int,is2PetscInt;
363: #if defined(PETSC_HAVE_COMPLEX)
364: PetscBool isPetscComplex;
365: #endif
368: MPIPetsc_Type_compare(unit,MPI_INT,&isInt);
369: MPIPetsc_Type_compare(unit,MPIU_INT,&isPetscInt);
370: MPIPetsc_Type_compare(unit,MPIU_REAL,&isPetscReal);
371: #if defined(PETSC_HAVE_COMPLEX)
372: MPIPetsc_Type_compare(unit,MPIU_COMPLEX,&isPetscComplex);
373: #endif
374: MPIPetsc_Type_compare(unit,MPI_2INT,&is2Int);
375: MPIPetsc_Type_compare(unit,MPIU_2INT,&is2PetscInt);
376: if (isInt) PackInit_int(link);
377: else if (isPetscInt) PackInit_PetscInt(link);
378: else if (isPetscReal) PackInit_PetscReal(link);
379: #if defined(PETSC_HAVE_COMPLEX)
380: else if (isPetscComplex) PackInit_PetscComplex(link);
381: #endif
382: else if (is2Int) PackInit_int_int(link);
383: else if (is2PetscInt) PackInit_PetscInt_PetscInt(link);
384: else {
385: PetscMPIInt bytes;
386: MPI_Type_size(unit,&bytes);
387: if (bytes % sizeof(int)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for type size not divisible by %D",sizeof(int));
388: switch (bytes / sizeof(int)) {
389: case 1: PackInit_block_int_1(link); break;
390: case 2: PackInit_block_int_2(link); break;
391: case 3: PackInit_block_int_3(link); break;
392: case 4: PackInit_block_int_4(link); break;
393: case 5: PackInit_block_int_5(link); break;
394: case 6: PackInit_block_int_6(link); break;
395: case 7: PackInit_block_int_7(link); break;
396: case 8: PackInit_block_int_8(link); break;
397: default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for arbitrary block sizes");
398: }
399: }
400: MPI_Type_dup(unit,&link->unit);
401: return(0);
402: }
406: static PetscErrorCode PetscSFBasicPackGetUnpackOp(PetscSF sf,PetscSFBasicPack link,MPI_Op op,void (**UnpackOp)(PetscInt,const PetscInt*,void*,const void*))
407: {
409: *UnpackOp = NULL;
410: if (op == MPIU_REPLACE) *UnpackOp = link->UnpackInsert;
411: else if (op == MPI_SUM || op == MPIU_SUM) *UnpackOp = link->UnpackAdd;
412: else if (op == MPI_MAX || op == MPIU_MAX) *UnpackOp = link->UnpackMax;
413: else if (op == MPI_MIN || op == MPIU_MIN) *UnpackOp = link->UnpackMin;
414: else if (op == MPI_MAXLOC) *UnpackOp = link->UnpackMaxloc;
415: else if (op == MPI_MINLOC) *UnpackOp = link->UnpackMinloc;
416: else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_SUP,"No support for MPI_Op");
417: return(0);
418: }
421: static PetscErrorCode PetscSFBasicPackGetFetchAndOp(PetscSF sf,PetscSFBasicPack link,MPI_Op op,void (**FetchAndOp)(PetscInt,const PetscInt*,void*,void*))
422: {
424: *FetchAndOp = NULL;
425: if (op == MPIU_REPLACE) *FetchAndOp = link->FetchAndInsert;
426: else if (op == MPI_SUM || op == MPIU_SUM) *FetchAndOp = link->FetchAndAdd;
427: else if (op == MPI_MAX || op == MPIU_MAX) *FetchAndOp = link->FetchAndMax;
428: else if (op == MPI_MIN || op == MPIU_MIN) *FetchAndOp = link->FetchAndMin;
429: else if (op == MPI_MAXLOC) *FetchAndOp = link->FetchAndMaxloc;
430: else if (op == MPI_MINLOC) *FetchAndOp = link->FetchAndMinloc;
431: else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_SUP,"No support for MPI_Op");
432: return(0);
433: }
437: static PetscErrorCode PetscSFBasicPackGetReqs(PetscSF sf,PetscSFBasicPack link,MPI_Request **rootreqs,MPI_Request **leafreqs)
438: {
439: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
442: if (rootreqs) *rootreqs = link->requests;
443: if (leafreqs) *leafreqs = link->requests + bas->niranks;
444: return(0);
445: }
449: static PetscErrorCode PetscSFBasicPackWaitall(PetscSF sf,PetscSFBasicPack link)
450: {
451: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
455: MPI_Waitall(bas->niranks+sf->nranks,link->requests,MPI_STATUSES_IGNORE);
456: return(0);
457: }
461: static PetscErrorCode PetscSFBasicGetRootInfo(PetscSF sf,PetscInt *nrootranks,const PetscMPIInt **rootranks,const PetscInt **rootoffset,const PetscInt **rootloc)
462: {
463: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
466: if (nrootranks) *nrootranks = bas->niranks;
467: if (rootranks) *rootranks = bas->iranks;
468: if (rootoffset) *rootoffset = bas->ioffset;
469: if (rootloc) *rootloc = bas->irootloc;
470: return(0);
471: }
475: static PetscErrorCode PetscSFBasicGetLeafInfo(PetscSF sf,PetscInt *nleafranks,const PetscMPIInt **leafranks,const PetscInt **leafoffset,const PetscInt **leafloc)
476: {
478: if (nleafranks) *nleafranks = sf->nranks;
479: if (leafranks) *leafranks = sf->ranks;
480: if (leafoffset) *leafoffset = sf->roffset;
481: if (leafloc) *leafloc = sf->rmine;
482: return(0);
483: }
487: static PetscErrorCode PetscSFBasicGetPack(PetscSF sf,MPI_Datatype unit,const void *key,PetscSFBasicPack *mylink)
488: {
489: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
490: PetscErrorCode ierr;
491: PetscSFBasicPack link,*p;
492: PetscInt nrootranks,nleafranks;
493: const PetscInt *rootoffset,*leafoffset;
496: /* Look for types in cache */
497: for (p=&bas->avail; (link=*p); p=&link->next) {
498: PetscBool match;
499: MPIPetsc_Type_compare(unit,link->unit,&match);
500: if (match) {
501: *p = link->next; /* Remove from available list */
502: goto found;
503: }
504: }
506: /* Create new composite types for each send rank */
507: PetscSFBasicGetRootInfo(sf,&nrootranks,NULL,&rootoffset,NULL);
508: PetscSFBasicGetLeafInfo(sf,&nleafranks,NULL,&leafoffset,NULL);
509: PetscNew(struct _n_PetscSFBasicPack,&link);
510: PetscSFBasicPackTypeSetup(link,unit);
511: PetscMalloc2(rootoffset[nrootranks]*link->unitbytes,char,&link->root,leafoffset[nleafranks]*link->unitbytes,char,&link->leaf);
512: PetscMalloc((nrootranks+nleafranks)*sizeof(MPI_Request),&link->requests);
514: found:
515: link->key = key;
516: link->next = bas->inuse;
517: bas->inuse = link;
519: *mylink = link;
520: return(0);
521: }
525: static PetscErrorCode PetscSFBasicGetPackInUse(PetscSF sf,MPI_Datatype unit,const void *key,PetscCopyMode cmode,PetscSFBasicPack *mylink)
526: {
527: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
528: PetscErrorCode ierr;
529: PetscSFBasicPack link,*p;
532: /* Look for types in cache */
533: for (p=&bas->inuse; (link=*p); p=&link->next) {
534: PetscBool match;
535: MPIPetsc_Type_compare(unit,link->unit,&match);
536: if (match) {
537: switch (cmode) {
538: case PETSC_OWN_POINTER: *p = link->next; break; /* Remove from inuse list */
539: case PETSC_USE_POINTER: break;
540: default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"invalid cmode");
541: }
542: *mylink = link;
543: return(0);
544: }
545: }
546: SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Could not find pack");
547: return(0);
548: }
552: static PetscErrorCode PetscSFBasicReclaimPack(PetscSF sf,PetscSFBasicPack *link)
553: {
554: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
557: (*link)->key = NULL;
558: (*link)->next = bas->avail;
559: bas->avail = *link;
560: *link = NULL;
561: return(0);
562: }
566: static PetscErrorCode PetscSFSetFromOptions_Basic(PetscSF sf)
567: {
571: PetscOptionsHead("PetscSF Basic options");
572: PetscOptionsTail();
573: return(0);
574: }
578: static PetscErrorCode PetscSFReset_Basic(PetscSF sf)
579: {
580: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
581: PetscErrorCode ierr;
582: PetscSFBasicPack link,next;
585: PetscFree(bas->iranks);
586: PetscFree2(bas->ioffset,bas->irootloc);
587: if (bas->inuse) SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Outstanding operation has not been completed");
588: for (link=bas->avail; link; link=next) {
589: next = link->next;
590: #if defined(PETSC_HAVE_MPI_TYPE_DUP)
591: MPI_Type_free(&link->unit);
592: #endif
593: PetscFree2(link->root,link->leaf);
594: PetscFree(link->requests);
595: PetscFree(link);
596: }
597: bas->avail = NULL;
598: return(0);
599: }
603: static PetscErrorCode PetscSFDestroy_Basic(PetscSF sf)
604: {
608: PetscSFReset_Basic(sf);
609: PetscFree(sf->data);
610: return(0);
611: }
615: static PetscErrorCode PetscSFView_Basic(PetscSF sf,PetscViewer viewer)
616: {
617: /* PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; */
619: PetscBool iascii;
622: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
623: if (iascii) {
624: PetscViewerASCIIPrintf(viewer," sort=%s\n",sf->rankorder ? "rank-order" : "unordered");
625: }
626: return(0);
627: }
631: /* Send from roots to leaves */
632: static PetscErrorCode PetscSFBcastBegin_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata)
633: {
634: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
635: PetscErrorCode ierr;
636: PetscSFBasicPack link;
637: PetscInt i,nrootranks,nleafranks;
638: const PetscInt *rootoffset,*leafoffset,*rootloc,*leafloc;
639: const PetscMPIInt *rootranks,*leafranks;
640: MPI_Request *rootreqs,*leafreqs;
641: size_t unitbytes;
644: PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
645: PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
646: PetscSFBasicGetPack(sf,unit,rootdata,&link);
648: unitbytes = link->unitbytes;
650: PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
651: /* Eagerly post leaf receives */
652: for (i=0; i<nleafranks; i++) {
653: PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
654: MPI_Irecv(link->leaf+leafoffset[i]*unitbytes,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
655: }
656: /* Pack and send root data */
657: for (i=0; i<nrootranks; i++) {
658: PetscMPIInt n = rootoffset[i+1] - rootoffset[i];
659: void *packstart = link->root+rootoffset[i]*unitbytes;
660: (*link->Pack)(n,rootloc+rootoffset[i],rootdata,packstart);
661: MPI_Isend(packstart,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
662: }
663: return(0);
664: }
668: PetscErrorCode PetscSFBcastEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata)
669: {
670: PetscErrorCode ierr;
671: PetscSFBasicPack link;
672: PetscInt i,nleafranks;
673: const PetscInt *leafoffset,*leafloc;
676: PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
677: PetscSFBasicPackWaitall(sf,link);
678: PetscSFBasicGetLeafInfo(sf,&nleafranks,NULL,&leafoffset,&leafloc);
679: for (i=0; i<nleafranks; i++) {
680: PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
681: const void *packstart = link->leaf+leafoffset[i]*link->unitbytes;
682: (*link->UnpackInsert)(n,leafloc+leafoffset[i],leafdata,packstart);
683: }
684: PetscSFBasicReclaimPack(sf,&link);
685: return(0);
686: }
690: /* leaf -> root with reduction */
691: PetscErrorCode PetscSFReduceBegin_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
692: {
693: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
694: PetscSFBasicPack link;
695: PetscErrorCode ierr;
696: PetscInt i,nrootranks,nleafranks;
697: const PetscInt *rootoffset,*leafoffset,*rootloc,*leafloc;
698: const PetscMPIInt *rootranks,*leafranks;
699: MPI_Request *rootreqs,*leafreqs;
700: size_t unitbytes;
703: PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
704: PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
705: PetscSFBasicGetPack(sf,unit,rootdata,&link);
707: unitbytes = link->unitbytes;
709: PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
710: /* Eagerly post root receives */
711: for (i=0; i<nrootranks; i++) {
712: PetscMPIInt n = rootoffset[i+1] - rootoffset[i];
713: MPI_Irecv(link->root+rootoffset[i]*unitbytes,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
714: }
715: /* Pack and send leaf data */
716: for (i=0; i<nleafranks; i++) {
717: PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
718: void *packstart = link->leaf+leafoffset[i]*unitbytes;
719: (*link->Pack)(n,leafloc+leafoffset[i],leafdata,packstart);
720: MPI_Isend(packstart,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
721: }
722: return(0);
723: }
727: static PetscErrorCode PetscSFReduceEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
728: {
729: void (*UnpackOp)(PetscInt,const PetscInt*,void*,const void*);
730: PetscErrorCode ierr;
731: PetscSFBasicPack link;
732: PetscInt i,nrootranks;
733: const PetscInt *rootoffset,*rootloc;
736: PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
737: /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
738: PetscSFBasicPackWaitall(sf,link);
739: PetscSFBasicGetRootInfo(sf,&nrootranks,NULL,&rootoffset,&rootloc);
740: PetscSFBasicPackGetUnpackOp(sf,link,op,&UnpackOp);
741: for (i=0; i<nrootranks; i++) {
742: PetscMPIInt n = rootoffset[i+1] - rootoffset[i];
743: const void *packstart = link->root+rootoffset[i]*link->unitbytes;
745: (*UnpackOp)(n,rootloc+rootoffset[i],rootdata,packstart);
746: }
747: PetscSFBasicReclaimPack(sf,&link);
748: return(0);
749: }
753: static PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
754: {
758: PetscSFReduceBegin_Basic(sf,unit,leafdata,rootdata,op);
759: return(0);
760: }
764: static PetscErrorCode PetscSFFetchAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
765: {
766: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
767: void (*FetchAndOp)(PetscInt,const PetscInt*,void*,void*);
768: PetscErrorCode ierr;
769: PetscSFBasicPack link;
770: PetscInt i,nrootranks,nleafranks;
771: const PetscInt *rootoffset,*leafoffset,*rootloc,*leafloc;
772: const PetscMPIInt *rootranks,*leafranks;
773: MPI_Request *rootreqs,*leafreqs;
774: size_t unitbytes;
777: PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
778: /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
779: PetscSFBasicPackWaitall(sf,link);
780: unitbytes = link->unitbytes;
781: PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
782: PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
783: PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
784: /* Post leaf receives */
785: for (i=0; i<nleafranks; i++) {
786: PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
787: MPI_Irecv(link->leaf+leafoffset[i]*unitbytes,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
788: }
789: /* Process local fetch-and-op, post root sends */
790: PetscSFBasicPackGetFetchAndOp(sf,link,op,&FetchAndOp);
791: for (i=0; i<nrootranks; i++) {
792: PetscMPIInt n = rootoffset[i+1] - rootoffset[i];
793: void *packstart = link->root+rootoffset[i]*unitbytes;
795: (*FetchAndOp)(n,rootloc+rootoffset[i],rootdata,packstart);
796: MPI_Isend(packstart,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
797: }
798: PetscSFBasicPackWaitall(sf,link);
799: for (i=0; i<nleafranks; i++) {
800: PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
801: const void *packstart = link->leaf+leafoffset[i]*unitbytes;
802: (*link->UnpackInsert)(n,leafloc+leafoffset[i],leafupdate,packstart);
803: }
804: PetscSFBasicReclaimPack(sf,&link);
805: return(0);
806: }
810: PETSC_EXTERN PetscErrorCode PetscSFCreate_Basic(PetscSF sf)
811: {
812: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
816: sf->ops->SetUp = PetscSFSetUp_Basic;
817: sf->ops->SetFromOptions = PetscSFSetFromOptions_Basic;
818: sf->ops->Reset = PetscSFReset_Basic;
819: sf->ops->Destroy = PetscSFDestroy_Basic;
820: sf->ops->View = PetscSFView_Basic;
821: sf->ops->BcastBegin = PetscSFBcastBegin_Basic;
822: sf->ops->BcastEnd = PetscSFBcastEnd_Basic;
823: sf->ops->ReduceBegin = PetscSFReduceBegin_Basic;
824: sf->ops->ReduceEnd = PetscSFReduceEnd_Basic;
825: sf->ops->FetchAndOpBegin = PetscSFFetchAndOpBegin_Basic;
826: sf->ops->FetchAndOpEnd = PetscSFFetchAndOpEnd_Basic;
828: PetscNewLog(sf,PetscSF_Basic,&bas);
829: sf->data = (void*)bas;
830: return(0);
831: }