Actual source code: sfbasic.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  2: #include <petsc/private/sfimpl.h> /*I "petscsf.h" I*/

  4: typedef struct _n_PetscSFBasicPack *PetscSFBasicPack;
  5: struct _n_PetscSFBasicPack {
  6:   void (*Pack)(PetscInt,PetscInt,const PetscInt*,const void*,void*);
  7:   void (*UnpackInsert)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
  8:   void (*UnpackAdd)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
  9:   void (*UnpackMin)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
 10:   void (*UnpackMax)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
 11:   void (*UnpackMinloc)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
 12:   void (*UnpackMaxloc)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
 13:   void (*UnpackMult)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 14:   void (*UnpackLAND)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 15:   void (*UnpackBAND)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 16:   void (*UnpackLOR)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 17:   void (*UnpackBOR)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 18:   void (*UnpackLXOR)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 19:   void (*UnpackBXOR)(PetscInt,PetscInt,const PetscInt*,void*,const void *);
 20:   void (*FetchAndInsert)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 21:   void (*FetchAndAdd)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 22:   void (*FetchAndMin)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 23:   void (*FetchAndMax)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 24:   void (*FetchAndMinloc)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 25:   void (*FetchAndMaxloc)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 26:   void (*FetchAndMult)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 27:   void (*FetchAndLAND)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 28:   void (*FetchAndBAND)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 29:   void (*FetchAndLOR)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 30:   void (*FetchAndBOR)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 31:   void (*FetchAndLXOR)(PetscInt,PetscInt,const PetscInt*,void*,void*);
 32:   void (*FetchAndBXOR)(PetscInt,PetscInt,const PetscInt*,void*,void*);

 34:   MPI_Datatype     unit;
 35:   size_t           unitbytes;   /* Number of bytes in a unit */
 36:   PetscInt         bs;          /* Number of basic units in a unit */
 37:   const void       *key;        /* Array used as key for operation */
 38:   char             *root;       /* Packed root data, contiguous by leaf rank */
 39:   char             *leaf;       /* Packed leaf data, contiguous by root rank */
 40:   MPI_Request      *requests;   /* Array of root requests followed by leaf requests */
 41:   PetscSFBasicPack next;
 42: };

 44: typedef struct {
 45:   PetscMPIInt      tag;
 46:   PetscMPIInt      niranks;     /* Number of incoming ranks (ranks accessing my roots) */
 47:   PetscMPIInt      *iranks;     /* Array of ranks that reference my roots */
 48:   PetscInt         itotal;      /* Total number of graph edges referencing my roots */
 49:   PetscInt         *ioffset;    /* Array of length niranks+1 holding offset in irootloc[] for each rank */
 50:   PetscInt         *irootloc;   /* Incoming roots referenced by ranks starting at ioffset[rank] */
 51:   PetscSFBasicPack avail;       /* One or more entries per MPI Datatype, lazily constructed */
 52:   PetscSFBasicPack inuse;       /* Buffers being used for transactions that have not yet completed */
 53: } PetscSF_Basic;

 55: #if !defined(PETSC_HAVE_MPI_TYPE_DUP) /* Danger: type is not reference counted; subject to ABA problem */
 56: PETSC_STATIC_INLINE PetscErrorCode MPI_Type_dup(MPI_Datatype datatype,MPI_Datatype *newtype)
 57: {
 58:   *newtype = datatype;
 59:   return 0;
 60: }
 61: #endif

 63: /*
 64:  * MPI_Reduce_local is not really useful because it can't handle sparse data and it vectorizes "in the wrong direction",
 65:  * therefore we pack data types manually. This section defines packing routines for the standard data types.
 66:  */

 68: #define CPPJoin2_exp(a,b) a ## b
 69: #define CPPJoin2(a,b) CPPJoin2_exp(a,b)
 70: #define CPPJoin3_exp_(a,b,c) a ## b ## _ ## c
 71: #define CPPJoin3_(a,b,c) CPPJoin3_exp_(a,b,c)

 73: /* Basic types without addition */
 74: #define DEF_PackNoInit(type,BS)                                         \
 75:   static void CPPJoin3_(Pack_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,const void *unpacked,void *packed) { \
 76:     const type *u = (const type*)unpacked;                              \
 77:     type *p = (type*)packed;                                            \
 78:     PetscInt i,j,k;                                                     \
 79:     for (i=0; i<n; i++)                                                 \
 80:       for (j=0; j<bs; j+=BS)                                            \
 81:         for (k=j; k<j+BS; k++)                                          \
 82:           p[i*bs+k] = u[idx[i]*bs+k];                                   \
 83:   }                                                                     \
 84:   static void CPPJoin3_(UnpackInsert_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
 85:     type *u = (type*)unpacked;                                          \
 86:     const type *p = (const type*)packed;                                \
 87:     PetscInt i,j,k;                                                     \
 88:     for (i=0; i<n; i++)                                                 \
 89:       for (j=0; j<bs; j+=BS)                                            \
 90:         for (k=j; k<j+BS; k++)                                          \
 91:           u[idx[i]*bs+k] = p[i*bs+k];                                   \
 92:   }                                                                     \
 93:   static void CPPJoin3_(FetchAndInsert_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
 94:     type *u = (type*)unpacked;                                          \
 95:     type *p = (type*)packed;                                            \
 96:     PetscInt i,j,k;                                                     \
 97:     for (i=0; i<n; i++) {                                               \
 98:       PetscInt ii = idx[i];                                             \
 99:       for (j=0; j<bs; j+=BS)                                            \
100:         for (k=j; k<j+BS; k++) {                                        \
101:           type t = u[ii*bs+k];                                          \
102:           u[ii*bs+k] = p[i*bs+k];                                       \
103:           p[i*bs+k] = t;                                                \
104:         }                                                               \
105:     }                                                                   \
106:   }

108: /* Basic types defining addition */
109: #define DEF_PackAddNoInit(type,BS)                                      \
110:   DEF_PackNoInit(type,BS)                                               \
111:   static void CPPJoin3_(UnpackAdd_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
112:     type *u = (type*)unpacked;                                          \
113:     const type *p = (const type*)packed;                                \
114:     PetscInt i,j,k;                                                     \
115:     for (i=0; i<n; i++)                                                 \
116:       for (j=0; j<bs; j+=BS)                                            \
117:         for (k=j; k<j+BS; k++)                                          \
118:           u[idx[i]*bs+k] += p[i*bs+k];                                  \
119:   }                                                                     \
120:   static void CPPJoin3_(FetchAndAdd_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
121:     type *u = (type*)unpacked;                                          \
122:     type *p = (type*)packed;                                            \
123:     PetscInt i,j,k;                                                     \
124:     for (i=0; i<n; i++) {                                               \
125:       PetscInt ii = idx[i];                                             \
126:       for (j=0; j<bs; j+=BS)                                            \
127:         for (k=j; k<j+BS; k++) {                                        \
128:           type t = u[ii*bs+k];                                          \
129:           u[ii*bs+k] = t + p[i*bs+k];                                   \
130:           p[i*bs+k] = t;                                                \
131:         }                                                               \
132:     }                                                                   \
133:   }                                                                     \
134:   static void CPPJoin3_(UnpackMult_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
135:     type *u = (type*)unpacked;                                          \
136:     const type *p = (const type*)packed;                                \
137:     PetscInt i,j,k;                                                     \
138:     for (i=0; i<n; i++)                                                 \
139:       for (j=0; j<bs; j+=BS)                                            \
140:         for (k=j; k<j+BS; k++)                                          \
141:           u[idx[i]*bs+k] *= p[i*bs+k];                                  \
142:   }                                                                     \
143:   static void CPPJoin3_(FetchAndMult_,type,BS)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
144:     type *u = (type*)unpacked;                                          \
145:     type *p = (type*)packed;                                            \
146:     PetscInt i,j,k;                                                     \
147:     for (i=0; i<n; i++) {                                               \
148:       PetscInt ii = idx[i];                                             \
149:       for (j=0; j<bs; j+=BS)                                            \
150:         for (k=j; k<j+BS; k++) {                                        \
151:           type t = u[ii*bs+k];                                          \
152:           u[ii*bs+k] = t * p[i*bs+k];                                   \
153:           p[i*bs+k] = t;                                                \
154:         }                                                               \
155:     }                                                                   \
156:   }
157: #define DEF_Pack(type,BS)                                               \
158:   DEF_PackAddNoInit(type,BS)                                            \
159:   static void CPPJoin3_(PackInit_,type,BS)(PetscSFBasicPack link) {     \
160:     link->Pack = CPPJoin3_(Pack_,type,BS);                              \
161:     link->UnpackInsert = CPPJoin3_(UnpackInsert_,type,BS);              \
162:     link->UnpackAdd = CPPJoin3_(UnpackAdd_,type,BS);                    \
163:     link->UnpackMult = CPPJoin3_(UnpackMult_,type,BS);                  \
164:     link->FetchAndInsert = CPPJoin3_(FetchAndInsert_,type,BS);          \
165:     link->FetchAndAdd = CPPJoin3_(FetchAndAdd_,type,BS);                \
166:     link->FetchAndMult = CPPJoin3_(FetchAndMult_,type,BS);              \
167:     link->unitbytes = sizeof(type);                                     \
168:   }
169: /* Comparable types */
170: #define DEF_PackCmp(type)                                               \
171:   DEF_PackAddNoInit(type,1)                                             \
172:   static void CPPJoin2(UnpackMax_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
173:     type *u = (type*)unpacked;                                          \
174:     const type *p = (const type*)packed;                                \
175:     PetscInt i;                                                         \
176:     for (i=0; i<n; i++) {                                               \
177:       type v = u[idx[i]];                                               \
178:       u[idx[i]] = PetscMax(v,p[i]);                                     \
179:     }                                                                   \
180:   }                                                                     \
181:   static void CPPJoin2(UnpackMin_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
182:     type *u = (type*)unpacked;                                          \
183:     const type *p = (const type*)packed;                                \
184:     PetscInt i;                                                         \
185:     for (i=0; i<n; i++) {                                               \
186:       type v = u[idx[i]];                                               \
187:       u[idx[i]] = PetscMin(v,p[i]);                                     \
188:     }                                                                   \
189:   }                                                                     \
190:   static void CPPJoin2(FetchAndMax_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
191:     type *u = (type*)unpacked;                                          \
192:     type *p = (type*)packed;                                            \
193:     PetscInt i;                                                         \
194:     for (i=0; i<n; i++) {                                               \
195:       PetscInt j = idx[i];                                              \
196:       type v = u[j];                                                    \
197:       u[j] = PetscMax(v,p[i]);                                          \
198:       p[i] = v;                                                         \
199:     }                                                                   \
200:   }                                                                     \
201:   static void CPPJoin2(FetchAndMin_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
202:     type *u = (type*)unpacked;                                          \
203:     type *p = (type*)packed;                                            \
204:     PetscInt i;                                                         \
205:     for (i=0; i<n; i++) {                                               \
206:       PetscInt j = idx[i];                                              \
207:       type v = u[j];                                                    \
208:       u[j] = PetscMin(v,p[i]);                                          \
209:       p[i] = v;                                                         \
210:     }                                                                   \
211:   }                                                                     \
212:   static void CPPJoin2(PackInit_,type)(PetscSFBasicPack link) {         \
213:     link->Pack = CPPJoin3_(Pack_,type,1);                               \
214:     link->UnpackInsert = CPPJoin3_(UnpackInsert_,type,1);               \
215:     link->UnpackAdd  = CPPJoin3_(UnpackAdd_,type,1);                    \
216:     link->UnpackMax  = CPPJoin2(UnpackMax_,type);                       \
217:     link->UnpackMin  = CPPJoin2(UnpackMin_,type);                       \
218:     link->UnpackMult = CPPJoin3_(UnpackMult_,type,1);                   \
219:     link->FetchAndInsert = CPPJoin3_(FetchAndInsert_,type,1);           \
220:     link->FetchAndAdd = CPPJoin3_(FetchAndAdd_ ,type,1);                \
221:     link->FetchAndMax = CPPJoin2(FetchAndMax_ ,type);                   \
222:     link->FetchAndMin = CPPJoin2(FetchAndMin_ ,type);                   \
223:     link->FetchAndMult = CPPJoin3_(FetchAndMult_,type,1);               \
224:     link->unitbytes = sizeof(type);                                     \
225:   }

227: /* Logical Types */
228: #define DEF_PackLog(type)                                               \
229:   static void CPPJoin2(UnpackLAND_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
230:     type *u = (type*)unpacked;                                          \
231:     const type *p = (const type*)packed;                                \
232:     PetscInt i;                                                         \
233:     for (i=0; i<n; i++) {                                               \
234:       type v = u[idx[i]];                                               \
235:       u[idx[i]] = v && p[i];                                            \
236:     }                                                                   \
237:   }                                                                     \
238:   static void CPPJoin2(UnpackLOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
239:     type *u = (type*)unpacked;                                          \
240:     const type *p = (const type*)packed;                                \
241:     PetscInt i;                                                         \
242:     for (i=0; i<n; i++) {                                               \
243:       type v = u[idx[i]];                                               \
244:       u[idx[i]] = v || p[i];                                            \
245:     }                                                                   \
246:   }                                                                     \
247:   static void CPPJoin2(UnpackLXOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
248:     type *u = (type*)unpacked;                                          \
249:     const type *p = (const type*)packed;                                \
250:     PetscInt i;                                                         \
251:     for (i=0; i<n; i++) {                                               \
252:       type v = u[idx[i]];                                               \
253:       u[idx[i]] = (!v)!=(!p[i]);                                        \
254:     }                                                                   \
255:   }                                                                     \
256:   static void CPPJoin2(FetchAndLAND_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
257:     type *u = (type*)unpacked;                                          \
258:     type *p = (type*)packed;                                            \
259:     PetscInt i;                                                         \
260:     for (i=0; i<n; i++) {                                               \
261:       PetscInt j = idx[i];                                              \
262:       type v = u[j];                                                    \
263:       u[j] = v && p[i];                                                 \
264:       p[i] = v;                                                         \
265:     }                                                                   \
266:   }                                                                     \
267:   static void CPPJoin2(FetchAndLOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
268:     type *u = (type*)unpacked;                                          \
269:     type *p = (type*)packed;                                            \
270:     PetscInt i;                                                         \
271:     for (i=0; i<n; i++) {                                               \
272:       PetscInt j = idx[i];                                              \
273:       type v = u[j];                                                    \
274:       u[j] = v || p[i];                                                 \
275:       p[i] = v;                                                         \
276:     }                                                                   \
277:   }                                                                     \
278:   static void CPPJoin2(FetchAndLXOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
279:     type *u = (type*)unpacked;                                          \
280:     type *p = (type*)packed;                                            \
281:     PetscInt i;                                                         \
282:     for (i=0; i<n; i++) {                                               \
283:       PetscInt j = idx[i];                                              \
284:       type v = u[j];                                                    \
285:       u[j] = (!v)!=(!p[i]);                                             \
286:       p[i] = v;                                                         \
287:     }                                                                   \
288:   }                                                                     \
289:   static void CPPJoin2(PackInit_Logical_,type)(PetscSFBasicPack link) { \
290:     link->UnpackLAND = CPPJoin2(UnpackLAND_,type);                      \
291:     link->UnpackLOR  = CPPJoin2(UnpackLOR_,type);                       \
292:     link->UnpackLXOR = CPPJoin2(UnpackLXOR_,type);                      \
293:     link->FetchAndLAND = CPPJoin2(FetchAndLAND_,type);                  \
294:     link->FetchAndLOR  = CPPJoin2(FetchAndLOR_,type);                   \
295:     link->FetchAndLXOR = CPPJoin2(FetchAndLXOR_,type);                  \
296:   }


299: /* Bitwise Types */
300: #define DEF_PackBit(type)                                               \
301:   static void CPPJoin2(UnpackBAND_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
302:     type *u = (type*)unpacked;                                          \
303:     const type *p = (const type*)packed;                                \
304:     PetscInt i;                                                         \
305:     for (i=0; i<n; i++) {                                               \
306:       type v = u[idx[i]];                                               \
307:       u[idx[i]] = v & p[i];                                             \
308:     }                                                                   \
309:   }                                                                     \
310:   static void CPPJoin2(UnpackBOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
311:     type *u = (type*)unpacked;                                          \
312:     const type *p = (const type*)packed;                                \
313:     PetscInt i;                                                         \
314:     for (i=0; i<n; i++) {                                               \
315:       type v = u[idx[i]];                                               \
316:       u[idx[i]] = v | p[i];                                             \
317:     }                                                                   \
318:   }                                                                     \
319:   static void CPPJoin2(UnpackBXOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
320:     type *u = (type*)unpacked;                                          \
321:     const type *p = (const type*)packed;                                \
322:     PetscInt i;                                                         \
323:     for (i=0; i<n; i++) {                                               \
324:       type v = u[idx[i]];                                               \
325:       u[idx[i]] = v^p[i];                                               \
326:     }                                                                   \
327:   }                                                                     \
328:   static void CPPJoin2(FetchAndBAND_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
329:     type *u = (type*)unpacked;                                          \
330:     type *p = (type*)packed;                                            \
331:     PetscInt i;                                                         \
332:     for (i=0; i<n; i++) {                                               \
333:       PetscInt j = idx[i];                                              \
334:       type v = u[j];                                                    \
335:       u[j] = v & p[i];                                                  \
336:       p[i] = v;                                                         \
337:     }                                                                   \
338:   }                                                                     \
339:   static void CPPJoin2(FetchAndBOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
340:     type *u = (type*)unpacked;                                          \
341:     type *p = (type*)packed;                                            \
342:     PetscInt i;                                                         \
343:     for (i=0; i<n; i++) {                                               \
344:       PetscInt j = idx[i];                                              \
345:       type v = u[j];                                                    \
346:       u[j] = v | p[i];                                                  \
347:       p[i] = v;                                                         \
348:     }                                                                   \
349:   }                                                                     \
350:   static void CPPJoin2(FetchAndBXOR_,type)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
351:     type *u = (type*)unpacked;                                          \
352:     type *p = (type*)packed;                                            \
353:     PetscInt i;                                                         \
354:     for (i=0; i<n; i++) {                                               \
355:       PetscInt j = idx[i];                                              \
356:       type v = u[j];                                                    \
357:       u[j] = v^p[i];                                                    \
358:       p[i] = v;                                                         \
359:     }                                                                   \
360:   }                                                                     \
361:   static void CPPJoin2(PackInit_Bitwise_,type)(PetscSFBasicPack link) { \
362:     link->UnpackBAND = CPPJoin2(UnpackBAND_,type);                      \
363:     link->UnpackBOR  = CPPJoin2(UnpackBOR_,type);                       \
364:     link->UnpackBXOR = CPPJoin2(UnpackBXOR_,type);                      \
365:     link->FetchAndBAND = CPPJoin2(FetchAndBAND_,type);                  \
366:     link->FetchAndBOR  = CPPJoin2(FetchAndBOR_,type);                   \
367:     link->FetchAndBXOR = CPPJoin2(FetchAndBXOR_,type);                  \
368:   }

370: /* Pair types */
371: #define CPPJoinloc_exp(base,op,t1,t2) base ## op ## loc_ ## t1 ## _ ## t2
372: #define CPPJoinloc(base,op,t1,t2) CPPJoinloc_exp(base,op,t1,t2)
373: #define PairType(type1,type2) CPPJoin3_(_pairtype_,type1,type2)
374: #define DEF_UnpackXloc(type1,type2,locname,op)                              \
375:   static void CPPJoinloc(Unpack,locname,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
376:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;        \
377:     const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
378:     PetscInt i;                                                         \
379:     for (i=0; i<n; i++) {                                               \
380:       PetscInt j = idx[i];                                              \
381:       if (p[i].a op u[j].a) {                                           \
382:         u[j].a = p[i].a;                                                \
383:         u[j].b = p[i].b;                                                \
384:       } else if (u[j].a == p[i].a) {                                    \
385:         u[j].b = PetscMin(u[j].b,p[i].b);                               \
386:       }                                                                 \
387:     }                                                                   \
388:   }                                                                     \
389:   static void CPPJoinloc(FetchAnd,locname,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
390:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;        \
391:     PairType(type1,type2) *p = (PairType(type1,type2)*)packed;          \
392:     PetscInt i;                                                         \
393:     for (i=0; i<n; i++) {                                               \
394:       PetscInt j = idx[i];                                              \
395:       PairType(type1,type2) v;                                          \
396:       v.a = u[j].a;                                                     \
397:       v.b = u[j].b;                                                     \
398:       if (p[i].a op u[j].a) {                                           \
399:         u[j].a = p[i].a;                                                \
400:         u[j].b = p[i].b;                                                \
401:       } else if (u[j].a == p[i].a) {                                    \
402:         u[j].b = PetscMin(u[j].b,p[i].b);                               \
403:       }                                                                 \
404:       p[i].a = v.a;                                                     \
405:       p[i].b = v.b;                                                     \
406:     }                                                                   \
407:   }
408: #define DEF_PackPair(type1,type2)                                       \
409:   typedef struct {type1 a; type2 b;} PairType(type1,type2);             \
410:   static void CPPJoin3_(Pack_,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,const void *unpacked,void *packed) { \
411:     const PairType(type1,type2) *u = (const PairType(type1,type2)*)unpacked; \
412:     PairType(type1,type2) *p = (PairType(type1,type2)*)packed;          \
413:     PetscInt i;                                                         \
414:     for (i=0; i<n; i++) {                                               \
415:       p[i].a = u[idx[i]].a;                                             \
416:       p[i].b = u[idx[i]].b;                                             \
417:     }                                                                   \
418:   }                                                                     \
419:   static void CPPJoin3_(UnpackInsert_,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
420:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;       \
421:     const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
422:     PetscInt i;                                                         \
423:     for (i=0; i<n; i++) {                                               \
424:       u[idx[i]].a = p[i].a;                                             \
425:       u[idx[i]].b = p[i].b;                                             \
426:     }                                                                   \
427:   }                                                                     \
428:   static void CPPJoin3_(UnpackAdd_,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,const void *packed) { \
429:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;       \
430:     const PairType(type1,type2) *p = (const PairType(type1,type2)*)packed; \
431:     PetscInt i;                                                         \
432:     for (i=0; i<n; i++) {                                               \
433:       u[idx[i]].a += p[i].a;                                            \
434:       u[idx[i]].b += p[i].b;                                            \
435:     }                                                                   \
436:   }                                                                     \
437:   static void CPPJoin3_(FetchAndInsert_,type1,type2)(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
438:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;        \
439:     PairType(type1,type2) *p = (PairType(type1,type2)*)packed;          \
440:     PetscInt i;                                                         \
441:     for (i=0; i<n; i++) {                                               \
442:       PetscInt j = idx[i];                                              \
443:       PairType(type1,type2) v;                                          \
444:       v.a = u[j].a;                                                     \
445:       v.b = u[j].b;                                                     \
446:       u[j].a = p[i].a;                                                  \
447:       u[j].b = p[i].b;                                                  \
448:       p[i].a = v.a;                                                     \
449:       p[i].b = v.b;                                                     \
450:     }                                                                   \
451:   }                                                                     \
452:   static void FetchAndAdd_ ## type1 ## _ ## type2(PetscInt n,PetscInt bs,const PetscInt *idx,void *unpacked,void *packed) { \
453:     PairType(type1,type2) *u = (PairType(type1,type2)*)unpacked;       \
454:     PairType(type1,type2) *p = (PairType(type1,type2)*)packed;         \
455:     PetscInt i;                                                         \
456:     for (i=0; i<n; i++) {                                               \
457:       PetscInt j = idx[i];                                              \
458:       PairType(type1,type2) v;                                          \
459:       v.a = u[j].a;                                                     \
460:       v.b = u[j].b;                                                     \
461:       u[j].a = v.a + p[i].a;                                            \
462:       u[j].b = v.b + p[i].b;                                            \
463:       p[i].a = v.a;                                                     \
464:       p[i].b = v.b;                                                     \
465:     }                                                                   \
466:   }                                                                     \
467:   DEF_UnpackXloc(type1,type2,Max,>)                                     \
468:   DEF_UnpackXloc(type1,type2,Min,<)                                     \
469:   static void CPPJoin3_(PackInit_,type1,type2)(PetscSFBasicPack link) { \
470:     link->Pack = CPPJoin3_(Pack_,type1,type2);                          \
471:     link->UnpackInsert = CPPJoin3_(UnpackInsert_,type1,type2);          \
472:     link->UnpackAdd = CPPJoin3_(UnpackAdd_,type1,type2);                \
473:     link->UnpackMaxloc = CPPJoin3_(UnpackMaxloc_,type1,type2);          \
474:     link->UnpackMinloc = CPPJoin3_(UnpackMinloc_,type1,type2);          \
475:     link->FetchAndInsert = CPPJoin3_(FetchAndInsert_,type1,type2);      \
476:     link->FetchAndAdd = CPPJoin3_(FetchAndAdd_,type1,type2);            \
477:     link->FetchAndMaxloc = CPPJoin3_(FetchAndMaxloc_,type1,type2);      \
478:     link->FetchAndMinloc = CPPJoin3_(FetchAndMinloc_,type1,type2);      \
479:     link->unitbytes = sizeof(PairType(type1,type2));                    \
480:   }

482: /* Currently only dumb blocks of data */
483: #define BlockType(unit,count) CPPJoin3_(_blocktype_,unit,count)
484: #define DEF_Block(unit,count)                                           \
485:   typedef struct {unit v[count];} BlockType(unit,count);                \
486:   DEF_PackNoInit(BlockType(unit,count),1)                               \
487:   static void CPPJoin3_(PackInit_block_,unit,count)(PetscSFBasicPack link) { \
488:     link->Pack = CPPJoin3_(Pack_,BlockType(unit,count),1);               \
489:     link->UnpackInsert = CPPJoin3_(UnpackInsert_,BlockType(unit,count),1); \
490:     link->FetchAndInsert = CPPJoin3_(FetchAndInsert_,BlockType(unit,count),1); \
491:     link->unitbytes = sizeof(BlockType(unit,count));                    \
492:   }

494: DEF_PackCmp(int)
495: DEF_PackBit(int)
496: DEF_PackLog(int)
497: DEF_PackCmp(PetscInt)
498: DEF_PackBit(PetscInt)
499: DEF_PackLog(PetscInt)
500: DEF_Pack(PetscInt,2)
501: DEF_Pack(PetscInt,3)
502: DEF_Pack(PetscInt,4)
503: DEF_Pack(PetscInt,5)
504: DEF_Pack(PetscInt,7)
505: DEF_PackCmp(PetscReal)
506: DEF_PackLog(PetscReal)
507: DEF_Pack(PetscReal,2)
508: DEF_Pack(PetscReal,3)
509: DEF_Pack(PetscReal,4)
510: DEF_Pack(PetscReal,5)
511: DEF_Pack(PetscReal,7)
512: #if defined(PETSC_HAVE_COMPLEX)
513: DEF_Pack(PetscComplex,1)
514: DEF_Pack(PetscComplex,2)
515: DEF_Pack(PetscComplex,3)
516: DEF_Pack(PetscComplex,4)
517: DEF_Pack(PetscComplex,5)
518: DEF_Pack(PetscComplex,7)
519: #endif
520: DEF_PackPair(int,int)
521: DEF_PackPair(PetscInt,PetscInt)
522: DEF_Block(int,1)
523: DEF_Block(int,2)
524: DEF_Block(int,3)
525: DEF_Block(int,4)
526: DEF_Block(int,5)
527: DEF_Block(int,6)
528: DEF_Block(int,7)
529: DEF_Block(int,8)

533: static PetscErrorCode PetscSFSetUp_Basic(PetscSF sf)
534: {
535:   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
537:   PetscInt *rlengths,*ilengths,i;
538:   MPI_Comm comm;
539:   MPI_Request *rootreqs,*leafreqs;

542:   PetscObjectGetComm((PetscObject)sf,&comm);
543:   PetscObjectGetNewTag((PetscObject)sf,&bas->tag);
544:   /*
545:    * Inform roots about how many leaves and from which ranks
546:    */
547:   PetscMalloc1(sf->nranks,&rlengths);
548:   /* Determine number, sending ranks, and length of incoming  */
549:   for (i=0; i<sf->nranks; i++) {
550:     rlengths[i] = sf->roffset[i+1] - sf->roffset[i]; /* Number of roots referenced by my leaves; for rank sf->ranks[i] */
551:   }
552:   PetscCommBuildTwoSided(comm,1,MPIU_INT,sf->nranks,sf->ranks,rlengths,&bas->niranks,&bas->iranks,(void**)&ilengths);
553:   PetscFree(rlengths);

555:   /* Send leaf identities to roots */
556:   for (i=0,bas->itotal=0; i<bas->niranks; i++) bas->itotal += ilengths[i];
557:   PetscMalloc2(bas->niranks+1,&bas->ioffset,bas->itotal,&bas->irootloc);
558:   PetscMalloc2(bas->niranks,&rootreqs,sf->nranks,&leafreqs);
559:   bas->ioffset[0] = 0;
560:   for (i=0; i<bas->niranks; i++) {
561:     bas->ioffset[i+1] = bas->ioffset[i] + ilengths[i];
562:     MPI_Irecv(bas->irootloc+bas->ioffset[i],ilengths[i],MPIU_INT,bas->iranks[i],bas->tag,comm,&rootreqs[i]);
563:   }
564:   for (i=0; i<sf->nranks; i++) {
565:     PetscMPIInt npoints;
566:     PetscMPIIntCast(sf->roffset[i+1] - sf->roffset[i],&npoints);
567:     MPI_Isend(sf->rremote+sf->roffset[i],npoints,MPIU_INT,sf->ranks[i],bas->tag,comm,&leafreqs[i]);
568:   }
569:   MPI_Waitall(bas->niranks,rootreqs,MPI_STATUSES_IGNORE);
570:   MPI_Waitall(sf->nranks,leafreqs,MPI_STATUSES_IGNORE);
571:   PetscFree(ilengths);
572:   PetscFree2(rootreqs,leafreqs);
573:   return(0);
574: }

578: static PetscErrorCode PetscSFBasicPackTypeSetup(PetscSFBasicPack link,MPI_Datatype unit)
579: {
581:   PetscBool      isInt,isPetscInt,isPetscReal,is2Int,is2PetscInt;
582:   PetscInt       nPetscIntContig,nPetscRealContig;
583: #if defined(PETSC_HAVE_COMPLEX)
584:   PetscBool isPetscComplex;
585:   PetscInt nPetscComplexContig;
586: #endif

589:   MPIPetsc_Type_compare(unit,MPI_INT,&isInt);
590:   MPIPetsc_Type_compare(unit,MPIU_INT,&isPetscInt);
591:   MPIPetsc_Type_compare_contig(unit,MPIU_INT,&nPetscIntContig);
592:   MPIPetsc_Type_compare(unit,MPIU_REAL,&isPetscReal);
593:   MPIPetsc_Type_compare_contig(unit,MPIU_REAL,&nPetscRealContig);
594: #if defined(PETSC_HAVE_COMPLEX)
595:   MPIPetsc_Type_compare(unit,MPIU_COMPLEX,&isPetscComplex);
596:   MPIPetsc_Type_compare_contig(unit,MPIU_COMPLEX,&nPetscComplexContig);
597: #endif
598:   MPIPetsc_Type_compare(unit,MPI_2INT,&is2Int);
599:   MPIPetsc_Type_compare(unit,MPIU_2INT,&is2PetscInt);
600:   link->bs = 1;
601:   if (isInt) {PackInit_int(link); PackInit_Logical_int(link); PackInit_Bitwise_int(link);}
602:   else if (isPetscInt) {PackInit_PetscInt(link); PackInit_Logical_PetscInt(link); PackInit_Bitwise_PetscInt(link);}
603:   else if (isPetscReal) {PackInit_PetscReal(link); PackInit_Logical_PetscReal(link);}
604: #if defined(PETSC_HAVE_COMPLEX)
605:   else if (isPetscComplex) PackInit_PetscComplex_1(link);
606: #endif
607:   else if (is2Int) PackInit_int_int(link);
608:   else if (is2PetscInt) PackInit_PetscInt_PetscInt(link);
609:   else if (nPetscIntContig) {
610:     if (nPetscIntContig%7 == 0) PackInit_PetscInt_7(link);
611:     else if (nPetscIntContig%5 == 0) PackInit_PetscInt_5(link);
612:     else if (nPetscIntContig%4 == 0) PackInit_PetscInt_4(link);
613:     else if (nPetscIntContig%3 == 0) PackInit_PetscInt_3(link);
614:     else if (nPetscIntContig%2 == 0) PackInit_PetscInt_2(link);
615:     else PackInit_PetscInt(link);
616:     link->bs = nPetscIntContig;
617:     link->unitbytes *= nPetscIntContig;
618:   } else if (nPetscRealContig) {
619:     if (nPetscRealContig%7 == 0) PackInit_PetscReal_7(link);
620:     else if (nPetscRealContig%5 == 0) PackInit_PetscReal_5(link);
621:     else if (nPetscRealContig%4 == 0) PackInit_PetscReal_4(link);
622:     else if (nPetscRealContig%3 == 0) PackInit_PetscReal_3(link);
623:     else if (nPetscRealContig%2 == 0) PackInit_PetscReal_2(link);
624:     else PackInit_PetscReal(link);
625:     link->bs = nPetscRealContig;
626:     link->unitbytes *= nPetscRealContig;
627: #if defined(PETSC_HAVE_COMPLEX)
628:   } else if (nPetscComplexContig) {
629:     if (nPetscComplexContig%7 == 0) PackInit_PetscComplex_7(link);
630:     else if (nPetscComplexContig%5 == 0) PackInit_PetscComplex_5(link);
631:     else if (nPetscComplexContig%4 == 0) PackInit_PetscComplex_4(link);
632:     else if (nPetscComplexContig%3 == 0) PackInit_PetscComplex_3(link);
633:     else if (nPetscComplexContig%2 == 0) PackInit_PetscComplex_2(link);
634:     else PackInit_PetscComplex_1(link);
635:     link->bs = nPetscComplexContig;
636:     link->unitbytes *= nPetscComplexContig;
637: #endif
638:   } else {
639:     MPI_Aint lb,bytes;
640:     MPI_Type_get_extent(unit,&lb,&bytes);
641:     if (lb != 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Datatype with nonzero lower bound %ld\n",(long)lb);
642:     if (bytes % sizeof(int)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for type size not divisible by %D",sizeof(int));
643:     switch (bytes / sizeof(int)) {
644:     case 1: PackInit_block_int_1(link); break;
645:     case 2: PackInit_block_int_2(link); break;
646:     case 3: PackInit_block_int_3(link); break;
647:     case 4: PackInit_block_int_4(link); break;
648:     case 5: PackInit_block_int_5(link); break;
649:     case 6: PackInit_block_int_6(link); break;
650:     case 7: PackInit_block_int_7(link); break;
651:     case 8: PackInit_block_int_8(link); break;
652:     default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for arbitrary block sizes");
653:     }
654:   }
655:   MPI_Type_dup(unit,&link->unit);
656:   return(0);
657: }

661: static PetscErrorCode PetscSFBasicPackGetUnpackOp(PetscSF sf,PetscSFBasicPack link,MPI_Op op,void (**UnpackOp)(PetscInt,PetscInt,const PetscInt*,void*,const void*))
662: {
664:   *UnpackOp = NULL;
665:   if (op == MPIU_REPLACE) *UnpackOp = link->UnpackInsert;
666:   else if (op == MPI_SUM || op == MPIU_SUM) *UnpackOp = link->UnpackAdd;
667:   else if (op == MPI_PROD) *UnpackOp = link->UnpackMult;
668:   else if (op == MPI_MAX || op == MPIU_MAX) *UnpackOp = link->UnpackMax;
669:   else if (op == MPI_MIN || op == MPIU_MIN) *UnpackOp = link->UnpackMin;
670:   else if (op == MPI_LAND) *UnpackOp = link->UnpackLAND;
671:   else if (op == MPI_BAND) *UnpackOp = link->UnpackBAND;
672:   else if (op == MPI_LOR) *UnpackOp = link->UnpackLOR;
673:   else if (op == MPI_BOR) *UnpackOp = link->UnpackBOR;
674:   else if (op == MPI_LXOR) *UnpackOp = link->UnpackLXOR;
675:   else if (op == MPI_BXOR) *UnpackOp = link->UnpackBXOR;
676:   else if (op == MPI_MAXLOC) *UnpackOp = link->UnpackMaxloc;
677:   else if (op == MPI_MINLOC) *UnpackOp = link->UnpackMinloc;
678:   else *UnpackOp = NULL;
679:   return(0);
680: }
683: static PetscErrorCode PetscSFBasicPackGetFetchAndOp(PetscSF sf,PetscSFBasicPack link,MPI_Op op,void (**FetchAndOp)(PetscInt,PetscInt,const PetscInt*,void*,void*))
684: {
686:   *FetchAndOp = NULL;
687:   if (op == MPIU_REPLACE) *FetchAndOp = link->FetchAndInsert;
688:   else if (op == MPI_SUM || op == MPIU_SUM) *FetchAndOp = link->FetchAndAdd;
689:   else if (op == MPI_MAX || op == MPIU_MAX) *FetchAndOp = link->FetchAndMax;
690:   else if (op == MPI_MIN || op == MPIU_MIN) *FetchAndOp = link->FetchAndMin;
691:   else if (op == MPI_MAXLOC) *FetchAndOp = link->FetchAndMaxloc;
692:   else if (op == MPI_MINLOC) *FetchAndOp = link->FetchAndMinloc;
693:   else if (op == MPI_PROD)   *FetchAndOp = link->FetchAndMult;
694:   else if (op == MPI_LAND)   *FetchAndOp = link->FetchAndLAND;
695:   else if (op == MPI_BAND)   *FetchAndOp = link->FetchAndBAND;
696:   else if (op == MPI_LOR)    *FetchAndOp = link->FetchAndLOR;
697:   else if (op == MPI_BOR)    *FetchAndOp = link->FetchAndBOR;
698:   else if (op == MPI_LXOR)   *FetchAndOp = link->FetchAndLXOR;
699:   else if (op == MPI_BXOR)   *FetchAndOp = link->FetchAndBXOR;
700:   else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_SUP,"No support for MPI_Op");
701:   return(0);
702: }

706: static PetscErrorCode PetscSFBasicPackGetReqs(PetscSF sf,PetscSFBasicPack link,MPI_Request **rootreqs,MPI_Request **leafreqs)
707: {
708:   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;

711:   if (rootreqs) *rootreqs = link->requests;
712:   if (leafreqs) *leafreqs = link->requests + bas->niranks;
713:   return(0);
714: }

718: static PetscErrorCode PetscSFBasicPackWaitall(PetscSF sf,PetscSFBasicPack link)
719: {
720:   PetscSF_Basic  *bas = (PetscSF_Basic*)sf->data;

724:   MPI_Waitall(bas->niranks+sf->nranks,link->requests,MPI_STATUSES_IGNORE);
725:   return(0);
726: }

730: static PetscErrorCode PetscSFBasicGetRootInfo(PetscSF sf,PetscInt *nrootranks,const PetscMPIInt **rootranks,const PetscInt **rootoffset,const PetscInt **rootloc)
731: {
732:   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;

735:   if (nrootranks) *nrootranks = bas->niranks;
736:   if (rootranks)  *rootranks  = bas->iranks;
737:   if (rootoffset) *rootoffset = bas->ioffset;
738:   if (rootloc)    *rootloc    = bas->irootloc;
739:   return(0);
740: }

744: static PetscErrorCode PetscSFBasicGetLeafInfo(PetscSF sf,PetscInt *nleafranks,const PetscMPIInt **leafranks,const PetscInt **leafoffset,const PetscInt **leafloc)
745: {
747:   if (nleafranks) *nleafranks = sf->nranks;
748:   if (leafranks)  *leafranks  = sf->ranks;
749:   if (leafoffset) *leafoffset = sf->roffset;
750:   if (leafloc)    *leafloc    = sf->rmine;
751:   return(0);
752: }

756: static PetscErrorCode PetscSFBasicGetPack(PetscSF sf,MPI_Datatype unit,const void *key,PetscSFBasicPack *mylink)
757: {
758:   PetscSF_Basic    *bas = (PetscSF_Basic*)sf->data;
759:   PetscErrorCode   ierr;
760:   PetscSFBasicPack link,*p;
761:   PetscInt         nrootranks,nleafranks;
762:   const PetscInt   *rootoffset,*leafoffset;

765:   /* Look for types in cache */
766:   for (p=&bas->avail; (link=*p); p=&link->next) {
767:     PetscBool match;
768:     MPIPetsc_Type_compare(unit,link->unit,&match);
769:     if (match) {
770:       *p = link->next;          /* Remove from available list */
771:       goto found;
772:     }
773:   }

775:   /* Create new composite types for each send rank */
776:   PetscSFBasicGetRootInfo(sf,&nrootranks,NULL,&rootoffset,NULL);
777:   PetscSFBasicGetLeafInfo(sf,&nleafranks,NULL,&leafoffset,NULL);
778:   PetscNew(&link);
779:   PetscSFBasicPackTypeSetup(link,unit);
780:   PetscCalloc2(rootoffset[nrootranks]*link->unitbytes,&link->root,leafoffset[nleafranks]*link->unitbytes,&link->leaf);
781:   PetscCalloc1(nrootranks+nleafranks,&link->requests);

783: found:
784:   link->key  = key;
785:   link->next = bas->inuse;
786:   bas->inuse = link;

788:   *mylink = link;
789:   return(0);
790: }

794: static PetscErrorCode PetscSFBasicGetPackInUse(PetscSF sf,MPI_Datatype unit,const void *key,PetscCopyMode cmode,PetscSFBasicPack *mylink)
795: {
796:   PetscSF_Basic    *bas = (PetscSF_Basic*)sf->data;
797:   PetscErrorCode   ierr;
798:   PetscSFBasicPack link,*p;

801:   /* Look for types in cache */
802:   for (p=&bas->inuse; (link=*p); p=&link->next) {
803:     PetscBool match;
804:     MPIPetsc_Type_compare(unit,link->unit,&match);
805:     if (match && (key == link->key)) {
806:       switch (cmode) {
807:       case PETSC_OWN_POINTER: *p = link->next; break; /* Remove from inuse list */
808:       case PETSC_USE_POINTER: break;
809:       default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"invalid cmode");
810:       }
811:       *mylink = link;
812:       return(0);
813:     }
814:   }
815:   SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Could not find pack");
816:   return(0);
817: }

821: static PetscErrorCode PetscSFBasicReclaimPack(PetscSF sf,PetscSFBasicPack *link)
822: {
823:   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;

826:   (*link)->key  = NULL;
827:   (*link)->next = bas->avail;
828:   bas->avail    = *link;
829:   *link         = NULL;
830:   return(0);
831: }

835: static PetscErrorCode PetscSFSetFromOptions_Basic(PetscOptionItems *PetscOptionsObject,PetscSF sf)
836: {

840:   PetscOptionsHead(PetscOptionsObject,"PetscSF Basic options");
841:   PetscOptionsTail();
842:   return(0);
843: }

847: static PetscErrorCode PetscSFReset_Basic(PetscSF sf)
848: {
849:   PetscSF_Basic    *bas = (PetscSF_Basic*)sf->data;
850:   PetscErrorCode   ierr;
851:   PetscSFBasicPack link,next;

854:   PetscFree(bas->iranks);
855:   PetscFree2(bas->ioffset,bas->irootloc);
856:   if (bas->inuse) SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Outstanding operation has not been completed");
857:   for (link=bas->avail; link; link=next) {
858:     next = link->next;
859: #if defined(PETSC_HAVE_MPI_TYPE_DUP)
860:     MPI_Type_free(&link->unit);
861: #endif
862:     PetscFree2(link->root,link->leaf);
863:     PetscFree(link->requests);
864:     PetscFree(link);
865:   }
866:   bas->avail = NULL;
867:   return(0);
868: }

872: static PetscErrorCode PetscSFDestroy_Basic(PetscSF sf)
873: {

877:   PetscSFReset_Basic(sf);
878:   PetscFree(sf->data);
879:   return(0);
880: }

884: static PetscErrorCode PetscSFView_Basic(PetscSF sf,PetscViewer viewer)
885: {
886:   /* PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; */
888:   PetscBool      iascii;

891:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
892:   if (iascii) {
893:     PetscViewerASCIIPrintf(viewer,"  sort=%s\n",sf->rankorder ? "rank-order" : "unordered");
894:   }
895:   return(0);
896: }

900: /* Send from roots to leaves */
901: static PetscErrorCode PetscSFBcastBegin_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata)
902: {
903:   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
904:   PetscErrorCode    ierr;
905:   PetscSFBasicPack  link;
906:   PetscInt          i,nrootranks,nleafranks;
907:   const PetscInt    *rootoffset,*leafoffset,*rootloc,*leafloc;
908:   const PetscMPIInt *rootranks,*leafranks;
909:   MPI_Request       *rootreqs,*leafreqs;
910:   size_t            unitbytes;

913:   PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
914:   PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
915:   PetscSFBasicGetPack(sf,unit,rootdata,&link);

917:   unitbytes = link->unitbytes;

919:   PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
920:   /* Eagerly post leaf receives */
921:   for (i=0; i<nleafranks; i++) {
922:     PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
923:     MPI_Irecv(link->leaf+leafoffset[i]*unitbytes,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
924:   }
925:   /* Pack and send root data */
926:   for (i=0; i<nrootranks; i++) {
927:     PetscMPIInt n          = rootoffset[i+1] - rootoffset[i];
928:     void        *packstart = link->root+rootoffset[i]*unitbytes;
929:     (*link->Pack)(n,link->bs,rootloc+rootoffset[i],rootdata,packstart);
930:     MPI_Isend(packstart,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
931:   }
932:   return(0);
933: }

937: PetscErrorCode PetscSFBcastEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata)
938: {
939:   PetscErrorCode   ierr;
940:   PetscSFBasicPack link;
941:   PetscInt         i,nleafranks;
942:   const PetscInt   *leafoffset,*leafloc;

945:   PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
946:   PetscSFBasicPackWaitall(sf,link);
947:   PetscSFBasicGetLeafInfo(sf,&nleafranks,NULL,&leafoffset,&leafloc);
948:   for (i=0; i<nleafranks; i++) {
949:     PetscMPIInt n          = leafoffset[i+1] - leafoffset[i];
950:     const void  *packstart = link->leaf+leafoffset[i]*link->unitbytes;
951:     (*link->UnpackInsert)(n,link->bs,leafloc+leafoffset[i],leafdata,packstart);
952:   }
953:   PetscSFBasicReclaimPack(sf,&link);
954:   return(0);
955: }

959: /* leaf -> root with reduction */
960: PetscErrorCode PetscSFReduceBegin_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
961: {
962:   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
963:   PetscSFBasicPack  link;
964:   PetscErrorCode    ierr;
965:   PetscInt          i,nrootranks,nleafranks;
966:   const PetscInt    *rootoffset,*leafoffset,*rootloc,*leafloc;
967:   const PetscMPIInt *rootranks,*leafranks;
968:   MPI_Request       *rootreqs,*leafreqs;
969:   size_t            unitbytes;

972:   PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
973:   PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
974:   PetscSFBasicGetPack(sf,unit,rootdata,&link);

976:   unitbytes = link->unitbytes;

978:   PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
979:   /* Eagerly post root receives */
980:   for (i=0; i<nrootranks; i++) {
981:     PetscMPIInt n = rootoffset[i+1] - rootoffset[i];
982:     MPI_Irecv(link->root+rootoffset[i]*unitbytes,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
983:   }
984:   /* Pack and send leaf data */
985:   for (i=0; i<nleafranks; i++) {
986:     PetscMPIInt n          = leafoffset[i+1] - leafoffset[i];
987:     void        *packstart = link->leaf+leafoffset[i]*unitbytes;
988:     (*link->Pack)(n,link->bs,leafloc+leafoffset[i],leafdata,packstart);
989:     MPI_Isend(packstart,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
990:   }
991:   return(0);
992: }

996: static PetscErrorCode PetscSFReduceEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
997: {
998:   void             (*UnpackOp)(PetscInt,PetscInt,const PetscInt*,void*,const void*);
999:   PetscErrorCode   ierr;
1000:   PetscSFBasicPack link;
1001:   PetscInt         i,nrootranks;
1002:   PetscMPIInt      typesize = -1;
1003:   const PetscInt   *rootoffset,*rootloc;

1006:   PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
1007:   /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
1008:   PetscSFBasicPackWaitall(sf,link);
1009:   PetscSFBasicGetRootInfo(sf,&nrootranks,NULL,&rootoffset,&rootloc);
1010:   PetscSFBasicPackGetUnpackOp(sf,link,op,&UnpackOp);
1011:   if (UnpackOp) {
1012:     typesize = link->unitbytes;
1013:   }
1014:   else {
1015:     MPI_Type_size(unit,&typesize);
1016:   }
1017:   for (i=0; i<nrootranks; i++) {
1018:     PetscMPIInt n   = rootoffset[i+1] - rootoffset[i];
1019:     char *packstart = (char *) link->root+rootoffset[i]*typesize;

1021:     if (UnpackOp) {
1022:       (*UnpackOp)(n,link->bs,rootloc+rootoffset[i],rootdata,(const void *)packstart);
1023:     }
1024: #if PETSC_HAVE_MPI_REDUCE_LOCAL
1025:     else if (n) { /* the op should be defined to operate on the whole datatype, so we ignore link->bs */
1026:       PetscInt j;

1028:       for (j = 0; j < n; j++) {
1029:         MPI_Reduce_local(packstart+j*typesize,((char *) rootdata)+(rootloc[rootoffset[i]+j])*typesize,1,unit,op);
1030:       }
1031:     }
1032: #else
1033:     else {
1034:       SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No unpacking reduction operation for this MPI_Op");
1035:     }
1036: #endif
1037:   }
1038:   PetscSFBasicReclaimPack(sf,&link);
1039:   return(0);
1040: }

1044: static PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
1045: {

1049:   PetscSFReduceBegin_Basic(sf,unit,leafdata,rootdata,op);
1050:   return(0);
1051: }

1055: static PetscErrorCode PetscSFFetchAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
1056: {
1057:   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
1058:   void              (*FetchAndOp)(PetscInt,PetscInt,const PetscInt*,void*,void*);
1059:   PetscErrorCode    ierr;
1060:   PetscSFBasicPack  link;
1061:   PetscInt          i,nrootranks,nleafranks;
1062:   const PetscInt    *rootoffset,*leafoffset,*rootloc,*leafloc;
1063:   const PetscMPIInt *rootranks,*leafranks;
1064:   MPI_Request       *rootreqs,*leafreqs;
1065:   size_t            unitbytes;

1068:   PetscSFBasicGetPackInUse(sf,unit,rootdata,PETSC_OWN_POINTER,&link);
1069:   /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
1070:   PetscSFBasicPackWaitall(sf,link);
1071:   unitbytes = link->unitbytes;
1072:   PetscSFBasicGetRootInfo(sf,&nrootranks,&rootranks,&rootoffset,&rootloc);
1073:   PetscSFBasicGetLeafInfo(sf,&nleafranks,&leafranks,&leafoffset,&leafloc);
1074:   PetscSFBasicPackGetReqs(sf,link,&rootreqs,&leafreqs);
1075:   /* Post leaf receives */
1076:   for (i=0; i<nleafranks; i++) {
1077:     PetscMPIInt n = leafoffset[i+1] - leafoffset[i];
1078:     MPI_Irecv(link->leaf+leafoffset[i]*unitbytes,n,unit,leafranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&leafreqs[i]);
1079:   }
1080:   /* Process local fetch-and-op, post root sends */
1081:   PetscSFBasicPackGetFetchAndOp(sf,link,op,&FetchAndOp);
1082:   for (i=0; i<nrootranks; i++) {
1083:     PetscMPIInt n          = rootoffset[i+1] - rootoffset[i];
1084:     void        *packstart = link->root+rootoffset[i]*unitbytes;

1086:     (*FetchAndOp)(n,link->bs,rootloc+rootoffset[i],rootdata,packstart);
1087:     MPI_Isend(packstart,n,unit,rootranks[i],bas->tag,PetscObjectComm((PetscObject)sf),&rootreqs[i]);
1088:   }
1089:   PetscSFBasicPackWaitall(sf,link);
1090:   for (i=0; i<nleafranks; i++) {
1091:     PetscMPIInt n          = leafoffset[i+1] - leafoffset[i];
1092:     const void  *packstart = link->leaf+leafoffset[i]*unitbytes;
1093:     (*link->UnpackInsert)(n,link->bs,leafloc+leafoffset[i],leafupdate,packstart);
1094:   }
1095:   PetscSFBasicReclaimPack(sf,&link);
1096:   return(0);
1097: }

1101: PETSC_EXTERN PetscErrorCode PetscSFCreate_Basic(PetscSF sf)
1102: {
1103:   PetscSF_Basic  *bas = (PetscSF_Basic*)sf->data;

1107:   sf->ops->SetUp           = PetscSFSetUp_Basic;
1108:   sf->ops->SetFromOptions  = PetscSFSetFromOptions_Basic;
1109:   sf->ops->Reset           = PetscSFReset_Basic;
1110:   sf->ops->Destroy         = PetscSFDestroy_Basic;
1111:   sf->ops->View            = PetscSFView_Basic;
1112:   sf->ops->BcastBegin      = PetscSFBcastBegin_Basic;
1113:   sf->ops->BcastEnd        = PetscSFBcastEnd_Basic;
1114:   sf->ops->ReduceBegin     = PetscSFReduceBegin_Basic;
1115:   sf->ops->ReduceEnd       = PetscSFReduceEnd_Basic;
1116:   sf->ops->FetchAndOpBegin = PetscSFFetchAndOpBegin_Basic;
1117:   sf->ops->FetchAndOpEnd   = PetscSFFetchAndOpEnd_Basic;

1119:   PetscNewLog(sf,&bas);
1120:   sf->data = (void*)bas;
1121:   return(0);
1122: }