Actual source code: data_bucket.c

  1: #include "../src/dm/impls/swarm/data_bucket.h"

  3: /* string helpers */
  4: PetscErrorCode DMSwarmDataFieldStringInList(const char name[], const PetscInt N, const DMSwarmDataField gfield[], PetscBool *val)
  5: {
  6:   PetscInt i;

  8:   *val = PETSC_FALSE;
  9:   for (i = 0; i < N; ++i) {
 10:     PetscBool flg;
 11:     PetscStrcmp(name, gfield[i]->name, &flg);
 12:     if (flg) {
 13:       *val = PETSC_TRUE;
 14:       return 0;
 15:     }
 16:   }
 17:   return 0;
 18: }

 20: PetscErrorCode DMSwarmDataFieldStringFindInList(const char name[], const PetscInt N, const DMSwarmDataField gfield[], PetscInt *index)
 21: {
 22:   PetscInt i;

 24:   *index = -1;
 25:   for (i = 0; i < N; ++i) {
 26:     PetscBool flg;
 27:     PetscStrcmp(name, gfield[i]->name, &flg);
 28:     if (flg) {
 29:       *index = i;
 30:       return 0;
 31:     }
 32:   }
 33:   return 0;
 34: }

 36: PetscErrorCode DMSwarmDataFieldCreate(const char registration_function[], const char name[], const size_t size, const PetscInt L, DMSwarmDataField *DF)
 37: {
 38:   DMSwarmDataField df;

 40:   PetscNew(&df);
 41:   PetscStrallocpy(registration_function, &df->registration_function);
 42:   PetscStrallocpy(name, &df->name);
 43:   df->atomic_size = size;
 44:   df->L           = L;
 45:   df->bs          = 1;
 46:   /* allocate something so we don't have to reallocate */
 47:   PetscMalloc(size * L, &df->data);
 48:   PetscMemzero(df->data, size * L);
 49:   *DF = df;
 50:   return 0;
 51: }

 53: PetscErrorCode DMSwarmDataFieldDestroy(DMSwarmDataField *DF)
 54: {
 55:   DMSwarmDataField df = *DF;

 57:   PetscFree(df->registration_function);
 58:   PetscFree(df->name);
 59:   PetscFree(df->data);
 60:   PetscFree(df);
 61:   *DF = NULL;
 62:   return 0;
 63: }

 65: /* data bucket */
 66: PetscErrorCode DMSwarmDataBucketCreate(DMSwarmDataBucket *DB)
 67: {
 68:   DMSwarmDataBucket db;

 70:   PetscNew(&db);

 72:   db->finalised = PETSC_FALSE;
 73:   /* create empty spaces for fields */
 74:   db->L         = -1;
 75:   db->buffer    = 1;
 76:   db->allocated = 1;
 77:   db->nfields   = 0;
 78:   PetscMalloc1(1, &db->field);
 79:   *DB = db;
 80:   return 0;
 81: }

 83: PetscErrorCode DMSwarmDataBucketDestroy(DMSwarmDataBucket *DB)
 84: {
 85:   DMSwarmDataBucket db = *DB;
 86:   PetscInt          f;

 88:   /* release fields */
 89:   for (f = 0; f < db->nfields; ++f) DMSwarmDataFieldDestroy(&db->field[f]);
 90:   /* this will catch the initially allocated objects in the event that no fields are registered */
 91:   if (db->field != NULL) PetscFree(db->field);
 92:   PetscFree(db);
 93:   *DB = NULL;
 94:   return 0;
 95: }

 97: PetscErrorCode DMSwarmDataBucketQueryForActiveFields(DMSwarmDataBucket db, PetscBool *any_active_fields)
 98: {
 99:   PetscInt f;

101:   *any_active_fields = PETSC_FALSE;
102:   for (f = 0; f < db->nfields; ++f) {
103:     if (db->field[f]->active) {
104:       *any_active_fields = PETSC_TRUE;
105:       return 0;
106:     }
107:   }
108:   return 0;
109: }

111: PetscErrorCode DMSwarmDataBucketRegisterField(DMSwarmDataBucket db, const char registration_function[], const char field_name[], size_t atomic_size, DMSwarmDataField *_gfield)
112: {
113:   PetscBool        val;
114:   DMSwarmDataField fp;

116:   /* check we haven't finalised the registration of fields */
117:   /*
118:    if (db->finalised==PETSC_TRUE) {
119:    printf("ERROR: DMSwarmDataBucketFinalize() has been called. Cannot register more fields\n");
120:    ERROR();
121:    }
122:   */
123:   /* check for repeated name */
124:   DMSwarmDataFieldStringInList(field_name, db->nfields, (const DMSwarmDataField *)db->field, &val);
126:   /* create new space for data */
127:   PetscRealloc(sizeof(DMSwarmDataField) * (db->nfields + 1), &db->field);
128:   /* add field */
129:   DMSwarmDataFieldCreate(registration_function, field_name, atomic_size, db->allocated, &fp);
130:   db->field[db->nfields] = fp;
131:   db->nfields++;
132:   if (_gfield != NULL) *_gfield = fp;
133:   return 0;
134: }

136: /*
137:  #define DMSwarmDataBucketRegisterField(db,name,size,k) {\
138:  char *location;\
139:  asprintf(&location,"Registered by %s() at line %d within file %s", __FUNCTION__, __LINE__, __FILE__);\
140:  _DMSwarmDataBucketRegisterField( (db), location, (name), (size), (k));\
141:  PetscFree(location);\
142:  }
143:  */

145: PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFieldIdByName(DMSwarmDataBucket db, const char name[], PetscInt *idx)
146: {
147:   PetscBool found;

149:   *idx = -1;
150:   DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, &found);
152:   DMSwarmDataFieldStringFindInList(name, db->nfields, (const DMSwarmDataField *)db->field, idx);
153:   return 0;
154: }

156: PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFieldByName(DMSwarmDataBucket db, const char name[], DMSwarmDataField *gfield)
157: {
158:   PetscInt  idx;
159:   PetscBool found;

161:   DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, &found);
163:   DMSwarmDataFieldStringFindInList(name, db->nfields, (const DMSwarmDataField *)db->field, &idx);
164:   *gfield = db->field[idx];
165:   return 0;
166: }

168: PetscErrorCode DMSwarmDataBucketQueryDMSwarmDataFieldByName(DMSwarmDataBucket db, const char name[], PetscBool *found)
169: {
170:   *found = PETSC_FALSE;
171:   DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, found);
172:   return 0;
173: }

175: PetscErrorCode DMSwarmDataBucketFinalize(DMSwarmDataBucket db)
176: {
177:   db->finalised = PETSC_TRUE;
178:   return 0;
179: }

181: PetscErrorCode DMSwarmDataFieldGetNumEntries(DMSwarmDataField df, PetscInt *sum)
182: {
183:   *sum = df->L;
184:   return 0;
185: }

187: PetscErrorCode DMSwarmDataFieldSetBlockSize(DMSwarmDataField df, PetscInt blocksize)
188: {
189:   df->bs = blocksize;
190:   return 0;
191: }

193: PetscErrorCode DMSwarmDataFieldSetSize(DMSwarmDataField df, const PetscInt new_L)
194: {
196:   if (new_L == df->L) return 0;
197:   if (new_L > df->L) {
198:     PetscRealloc(df->atomic_size * (new_L), &df->data);
199:     /* init new contents */
200:     PetscMemzero((((char *)df->data) + df->L * df->atomic_size), (new_L - df->L) * df->atomic_size);
201:   } else {
202:     /* reallocate pointer list, add +1 in case new_L = 0 */
203:     PetscRealloc(df->atomic_size * (new_L + 1), &df->data);
204:   }
205:   df->L = new_L;
206:   return 0;
207: }

209: PetscErrorCode DMSwarmDataFieldZeroBlock(DMSwarmDataField df, const PetscInt start, const PetscInt end)
210: {
214:   PetscMemzero((((char *)df->data) + start * df->atomic_size), (end - start) * df->atomic_size);
215:   return 0;
216: }

218: /*
219:  A negative buffer value will simply be ignored and the old buffer value will be used.
220:  */
221: PetscErrorCode DMSwarmDataBucketSetSizes(DMSwarmDataBucket db, const PetscInt L, const PetscInt buffer)
222: {
223:   PetscInt  current_allocated, new_used, new_unused, new_buffer, new_allocated, f;
224:   PetscBool any_active_fields;

227:   DMSwarmDataBucketQueryForActiveFields(db, &any_active_fields);

230:   current_allocated = db->allocated;
231:   new_used          = L;
232:   new_unused        = current_allocated - new_used;
233:   new_buffer        = db->buffer;
234:   if (buffer >= 0) { /* update the buffer value */
235:     new_buffer = buffer;
236:   }
237:   new_allocated = new_used + new_buffer;
238:   /* action */
239:   if (new_allocated > current_allocated) {
240:     /* increase size to new_used + new_buffer */
241:     for (f = 0; f < db->nfields; f++) DMSwarmDataFieldSetSize(db->field[f], new_allocated);
242:     db->L         = new_used;
243:     db->buffer    = new_buffer;
244:     db->allocated = new_used + new_buffer;
245:   } else {
246:     if (new_unused > 2 * new_buffer) {
247:       /* shrink array to new_used + new_buffer */
248:       for (f = 0; f < db->nfields; ++f) DMSwarmDataFieldSetSize(db->field[f], new_allocated);
249:       db->L         = new_used;
250:       db->buffer    = new_buffer;
251:       db->allocated = new_used + new_buffer;
252:     } else {
253:       db->L      = new_used;
254:       db->buffer = new_buffer;
255:     }
256:   }
257:   /* zero all entries from db->L to db->allocated */
258:   for (f = 0; f < db->nfields; ++f) {
259:     DMSwarmDataField field = db->field[f];
260:     DMSwarmDataFieldZeroBlock(field, db->L, db->allocated);
261:   }
262:   return 0;
263: }

265: PetscErrorCode DMSwarmDataBucketSetInitialSizes(DMSwarmDataBucket db, const PetscInt L, const PetscInt buffer)
266: {
267:   PetscInt f;

269:   DMSwarmDataBucketSetSizes(db, L, buffer);
270:   for (f = 0; f < db->nfields; ++f) {
271:     DMSwarmDataField field = db->field[f];
272:     DMSwarmDataFieldZeroBlock(field, 0, db->allocated);
273:   }
274:   return 0;
275: }

277: PetscErrorCode DMSwarmDataBucketGetSizes(DMSwarmDataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated)
278: {
279:   if (L) *L = db->L;
280:   if (buffer) *buffer = db->buffer;
281:   if (allocated) *allocated = db->allocated;
282:   return 0;
283: }

285: PetscErrorCode DMSwarmDataBucketGetGlobalSizes(MPI_Comm comm, DMSwarmDataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated)
286: {
287:   if (L) MPI_Allreduce(&db->L, L, 1, MPIU_INT, MPI_SUM, comm);
288:   if (buffer) MPI_Allreduce(&db->buffer, buffer, 1, MPIU_INT, MPI_SUM, comm);
289:   if (allocated) MPI_Allreduce(&db->allocated, allocated, 1, MPIU_INT, MPI_SUM, comm);
290:   return 0;
291: }

293: PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFields(DMSwarmDataBucket db, PetscInt *L, DMSwarmDataField *fields[])
294: {
295:   if (L) *L = db->nfields;
296:   if (fields) *fields = db->field;
297:   return 0;
298: }

300: PetscErrorCode DMSwarmDataFieldGetAccess(const DMSwarmDataField gfield)
301: {
303:   gfield->active = PETSC_TRUE;
304:   return 0;
305: }

307: PetscErrorCode DMSwarmDataFieldAccessPoint(const DMSwarmDataField gfield, const PetscInt pid, void **ctx_p)
308: {
309:   *ctx_p = NULL;
310: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
311:   /* debug mode */
312:   /* check point is valid */
316: #endif
317:   *ctx_p = DMSWARM_DATAFIELD_point_access(gfield->data, pid, gfield->atomic_size);
318:   return 0;
319: }

321: PetscErrorCode DMSwarmDataFieldAccessPointOffset(const DMSwarmDataField gfield, const size_t offset, const PetscInt pid, void **ctx_p)
322: {
323: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
324:   /* debug mode */
325:   /* check point is valid */
327:   /* Note compiler realizes this can never happen with an unsigned PetscInt */
329:   /* check point is valid */
333: #endif
334:   *ctx_p = DMSWARM_DATAFIELD_point_access_offset(gfield->data, pid, gfield->atomic_size, offset);
335:   return 0;
336: }

338: PetscErrorCode DMSwarmDataFieldRestoreAccess(DMSwarmDataField gfield)
339: {
341:   gfield->active = PETSC_FALSE;
342:   return 0;
343: }

345: PetscErrorCode DMSwarmDataFieldVerifyAccess(const DMSwarmDataField gfield, const size_t size)
346: {
347: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
349: #endif
350:   return 0;
351: }

353: PetscErrorCode DMSwarmDataFieldGetAtomicSize(const DMSwarmDataField gfield, size_t *size)
354: {
355:   if (size) *size = gfield->atomic_size;
356:   return 0;
357: }

359: PetscErrorCode DMSwarmDataFieldGetEntries(const DMSwarmDataField gfield, void **data)
360: {
361:   if (data) *data = gfield->data;
362:   return 0;
363: }

365: PetscErrorCode DMSwarmDataFieldRestoreEntries(const DMSwarmDataField gfield, void **data)
366: {
367:   if (data) *data = NULL;
368:   return 0;
369: }

371: /* y = x */
372: PetscErrorCode DMSwarmDataBucketCopyPoint(const DMSwarmDataBucket xb, const PetscInt pid_x, const DMSwarmDataBucket yb, const PetscInt pid_y)
373: {
374:   PetscInt f;

376:   for (f = 0; f < xb->nfields; ++f) {
377:     void *dest;
378:     void *src;

380:     DMSwarmDataFieldGetAccess(xb->field[f]);
381:     if (xb != yb) DMSwarmDataFieldGetAccess(yb->field[f]);
382:     DMSwarmDataFieldAccessPoint(xb->field[f], pid_x, &src);
383:     DMSwarmDataFieldAccessPoint(yb->field[f], pid_y, &dest);
384:     PetscMemcpy(dest, src, xb->field[f]->atomic_size);
385:     DMSwarmDataFieldRestoreAccess(xb->field[f]);
386:     if (xb != yb) DMSwarmDataFieldRestoreAccess(yb->field[f]);
387:   }
388:   return 0;
389: }

391: PetscErrorCode DMSwarmDataBucketCreateFromSubset(DMSwarmDataBucket DBIn, const PetscInt N, const PetscInt list[], DMSwarmDataBucket *DB)
392: {
393:   PetscInt          nfields;
394:   DMSwarmDataField *fields;
395:   PetscInt          f, L, buffer, allocated, p;

397:   DMSwarmDataBucketCreate(DB);
398:   /* copy contents of DBIn */
399:   DMSwarmDataBucketGetDMSwarmDataFields(DBIn, &nfields, &fields);
400:   DMSwarmDataBucketGetSizes(DBIn, &L, &buffer, &allocated);
401:   for (f = 0; f < nfields; ++f) DMSwarmDataBucketRegisterField(*DB, "DMSwarmDataBucketCreateFromSubset", fields[f]->name, fields[f]->atomic_size, NULL);
402:   DMSwarmDataBucketFinalize(*DB);
403:   DMSwarmDataBucketSetSizes(*DB, L, buffer);
404:   /* now copy the desired guys from DBIn => DB */
405:   for (p = 0; p < N; ++p) DMSwarmDataBucketCopyPoint(DBIn, list[p], *DB, list[p]);
406:   return 0;
407: }

409: /* insert into an existing location */
410: PetscErrorCode DMSwarmDataFieldInsertPoint(const DMSwarmDataField field, const PetscInt index, const void *ctx)
411: {
412: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
413:   /* check point is valid */
416: #endif
417:   PetscMemcpy(DMSWARM_DATAFIELD_point_access(field->data, index, field->atomic_size), ctx, field->atomic_size);
418:   return 0;
419: }

421: /* remove data at index - replace with last point */
422: PetscErrorCode DMSwarmDataBucketRemovePointAtIndex(const DMSwarmDataBucket db, const PetscInt index)
423: {
424:   PetscInt  f;
425:   PetscBool any_active_fields;

427: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
428:   /* check point is valid */
431: #endif
432:   DMSwarmDataBucketQueryForActiveFields(db, &any_active_fields);
434:   if (index >= db->L) { /* this point is not in the list - no need to error, but I will anyway */
435:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_USER, "You should not be trying to remove point at index=%" PetscInt_FMT " since it's < db->L = %" PetscInt_FMT, index, db->L);
436:   }
437:   if (index != db->L - 1) { /* not last point in list */
438:     for (f = 0; f < db->nfields; ++f) {
439:       DMSwarmDataField field = db->field[f];

441:       /* copy then remove */
442:       DMSwarmDataFieldCopyPoint(db->L - 1, field, index, field);
443:       /* DMSwarmDataFieldZeroPoint(field,index); */
444:     }
445:   }
446:   /* decrement size */
447:   /* this will zero out an crap at the end of the list */
448:   DMSwarmDataBucketRemovePoint(db);
449:   return 0;
450: }

452: /* copy x into y */
453: PetscErrorCode DMSwarmDataFieldCopyPoint(const PetscInt pid_x, const DMSwarmDataField field_x, const PetscInt pid_y, const DMSwarmDataField field_y)
454: {
455: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
456:   /* check point is valid */
462: #endif
463:   PetscMemcpy(DMSWARM_DATAFIELD_point_access(field_y->data, pid_y, field_y->atomic_size), DMSWARM_DATAFIELD_point_access(field_x->data, pid_x, field_x->atomic_size), field_y->atomic_size);
464:   return 0;
465: }

467: /* zero only the datafield at this point */
468: PetscErrorCode DMSwarmDataFieldZeroPoint(const DMSwarmDataField field, const PetscInt index)
469: {
470: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
471:   /* check point is valid */
474: #endif
475:   PetscMemzero(DMSWARM_DATAFIELD_point_access(field->data, index, field->atomic_size), field->atomic_size);
476:   return 0;
477: }

479: /* zero ALL data for this point */
480: PetscErrorCode DMSwarmDataBucketZeroPoint(const DMSwarmDataBucket db, const PetscInt index)
481: {
482:   PetscInt f;

484:   /* check point is valid */
487:   for (f = 0; f < db->nfields; ++f) {
488:     DMSwarmDataField field = db->field[f];
489:     DMSwarmDataFieldZeroPoint(field, index);
490:   }
491:   return 0;
492: }

494: /* increment */
495: PetscErrorCode DMSwarmDataBucketAddPoint(DMSwarmDataBucket db)
496: {
497:   DMSwarmDataBucketSetSizes(db, db->L + 1, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
498:   return 0;
499: }

501: /* decrement */
502: PetscErrorCode DMSwarmDataBucketRemovePoint(DMSwarmDataBucket db)
503: {
504:   DMSwarmDataBucketSetSizes(db, db->L - 1, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
505:   return 0;
506: }

508: /*  Should be redone to user PetscViewer */
509: PetscErrorCode DMSwarmDataBucketView_stdout(MPI_Comm comm, DMSwarmDataBucket db)
510: {
511:   PetscInt f;
512:   double   memory_usage_total, memory_usage_total_local = 0.0;

514:   PetscPrintf(comm, "DMSwarmDataBucketView: \n");
515:   PetscPrintf(comm, "  L                  = %" PetscInt_FMT " \n", db->L);
516:   PetscPrintf(comm, "  buffer             = %" PetscInt_FMT " \n", db->buffer);
517:   PetscPrintf(comm, "  allocated          = %" PetscInt_FMT " \n", db->allocated);
518:   PetscPrintf(comm, "  nfields registered = %" PetscInt_FMT " \n", db->nfields);

520:   for (f = 0; f < db->nfields; ++f) {
521:     double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
522:     memory_usage_total_local += memory_usage_f;
523:   }
524:   MPI_Allreduce(&memory_usage_total_local, &memory_usage_total, 1, MPI_DOUBLE, MPI_SUM, comm);

526:   for (f = 0; f < db->nfields; ++f) {
527:     double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
528:     PetscPrintf(comm, "    [%3" PetscInt_FMT "] %15s : Mem. usage       = %1.2e (MB) [rank0]\n", f, db->field[f]->name, memory_usage_f);
529:     PetscPrintf(comm, "                            blocksize        = %" PetscInt_FMT " \n", db->field[f]->bs);
530:     if (db->field[f]->bs != 1) {
531:       PetscPrintf(comm, "                            atomic size      = %zu [full block, bs=%" PetscInt_FMT "]\n", db->field[f]->atomic_size, db->field[f]->bs);
532:       PetscPrintf(comm, "                            atomic size/item = %zu \n", (size_t)(db->field[f]->atomic_size / db->field[f]->bs));
533:     } else {
534:       PetscPrintf(comm, "                            atomic size      = %zu \n", db->field[f]->atomic_size);
535:     }
536:   }
537:   PetscPrintf(comm, "  Total mem. usage                           = %1.2e (MB) (collective)\n", memory_usage_total);
538:   return 0;
539: }

541: PetscErrorCode DMSwarmDataBucketView_Seq(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
542: {
543:   switch (type) {
544:   case DATABUCKET_VIEW_STDOUT:
545:     DMSwarmDataBucketView_stdout(PETSC_COMM_SELF, db);
546:     break;
547:   case DATABUCKET_VIEW_ASCII:
548:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for ascii output");
549:   case DATABUCKET_VIEW_BINARY:
550:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for binary output");
551:   case DATABUCKET_VIEW_HDF5:
552:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for HDF5 output");
553:   default:
554:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown viewer method requested");
555:   }
556:   return 0;
557: }

559: PetscErrorCode DMSwarmDataBucketView_MPI(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
560: {
561:   switch (type) {
562:   case DATABUCKET_VIEW_STDOUT:
563:     DMSwarmDataBucketView_stdout(comm, db);
564:     break;
565:   case DATABUCKET_VIEW_ASCII:
566:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for ascii output");
567:   case DATABUCKET_VIEW_BINARY:
568:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for binary output");
569:   case DATABUCKET_VIEW_HDF5:
570:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for HDF5 output");
571:   default:
572:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown viewer method requested");
573:   }
574:   return 0;
575: }

577: PetscErrorCode DMSwarmDataBucketView(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
578: {
579:   PetscMPIInt size;

581:   MPI_Comm_size(comm, &size);
582:   if (size == 1) {
583:     DMSwarmDataBucketView_Seq(comm, db, filename, type);
584:   } else {
585:     DMSwarmDataBucketView_MPI(comm, db, filename, type);
586:   }
587:   return 0;
588: }

590: PetscErrorCode DMSwarmDataBucketDuplicateFields(DMSwarmDataBucket dbA, DMSwarmDataBucket *dbB)
591: {
592:   DMSwarmDataBucket db2;
593:   PetscInt          f;

595:   DMSwarmDataBucketCreate(&db2);
596:   /* copy contents from dbA into db2 */
597:   for (f = 0; f < dbA->nfields; ++f) {
598:     DMSwarmDataField field;
599:     size_t           atomic_size;
600:     char            *name;

602:     field       = dbA->field[f];
603:     atomic_size = field->atomic_size;
604:     name        = field->name;
605:     DMSwarmDataBucketRegisterField(db2, "DMSwarmDataBucketDuplicateFields", name, atomic_size, NULL);
606:   }
607:   DMSwarmDataBucketFinalize(db2);
608:   DMSwarmDataBucketSetInitialSizes(db2, 0, 1000);
609:   *dbB = db2;
610:   return 0;
611: }

613: /*
614:  Insert points from db2 into db1
615:  db1 <<== db2
616:  */
617: PetscErrorCode DMSwarmDataBucketInsertValues(DMSwarmDataBucket db1, DMSwarmDataBucket db2)
618: {
619:   PetscInt n_mp_points1, n_mp_points2;
620:   PetscInt n_mp_points1_new, p;

622:   DMSwarmDataBucketGetSizes(db1, &n_mp_points1, NULL, NULL);
623:   DMSwarmDataBucketGetSizes(db2, &n_mp_points2, NULL, NULL);
624:   n_mp_points1_new = n_mp_points1 + n_mp_points2;
625:   DMSwarmDataBucketSetSizes(db1, n_mp_points1_new, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
626:   for (p = 0; p < n_mp_points2; ++p) {
627:     /* db1 <<== db2 */
628:     DMSwarmDataBucketCopyPoint(db2, p, db1, (n_mp_points1 + p));
629:   }
630:   return 0;
631: }

633: /* helpers for parallel send/recv */
634: PetscErrorCode DMSwarmDataBucketCreatePackedArray(DMSwarmDataBucket db, size_t *bytes, void **buf)
635: {
636:   PetscInt f;
637:   size_t   sizeof_marker_contents;
638:   void    *buffer;

640:   sizeof_marker_contents = 0;
641:   for (f = 0; f < db->nfields; ++f) {
642:     DMSwarmDataField df = db->field[f];
643:     sizeof_marker_contents += df->atomic_size;
644:   }
645:   PetscMalloc(sizeof_marker_contents, &buffer);
646:   PetscMemzero(buffer, sizeof_marker_contents);
647:   if (bytes) *bytes = sizeof_marker_contents;
648:   if (buf) *buf = buffer;
649:   return 0;
650: }

652: PetscErrorCode DMSwarmDataBucketDestroyPackedArray(DMSwarmDataBucket db, void **buf)
653: {
654:   if (buf) {
655:     PetscFree(*buf);
656:     *buf = NULL;
657:   }
658:   return 0;
659: }

661: PetscErrorCode DMSwarmDataBucketFillPackedArray(DMSwarmDataBucket db, const PetscInt index, void *buf)
662: {
663:   PetscInt f;
664:   void    *data, *data_p;
665:   size_t   asize, offset;

667:   offset = 0;
668:   for (f = 0; f < db->nfields; ++f) {
669:     DMSwarmDataField df = db->field[f];

671:     asize  = df->atomic_size;
672:     data   = (void *)(df->data);
673:     data_p = (void *)((char *)data + index * asize);
674:     PetscMemcpy((void *)((char *)buf + offset), data_p, asize);
675:     offset = offset + asize;
676:   }
677:   return 0;
678: }

680: PetscErrorCode DMSwarmDataBucketInsertPackedArray(DMSwarmDataBucket db, const PetscInt idx, void *data)
681: {
682:   PetscInt f;
683:   void    *data_p;
684:   size_t   offset;

686:   offset = 0;
687:   for (f = 0; f < db->nfields; ++f) {
688:     DMSwarmDataField df = db->field[f];

690:     data_p = (void *)((char *)data + offset);
691:     DMSwarmDataFieldInsertPoint(df, idx, (void *)data_p);
692:     offset = offset + df->atomic_size;
693:   }
694:   return 0;
695: }