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) {
90: DMSwarmDataFieldDestroy(&db->field[f]);
91: }
92: /* this will catch the initially allocated objects in the event that no fields are registered */
93: if (db->field != NULL) {
94: PetscFree(db->field);
95: }
96: PetscFree(db);
97: *DB = NULL;
98: return 0;
99: }
101: PetscErrorCode DMSwarmDataBucketQueryForActiveFields(DMSwarmDataBucket db,PetscBool *any_active_fields)
102: {
103: PetscInt f;
105: *any_active_fields = PETSC_FALSE;
106: for (f = 0; f < db->nfields; ++f) {
107: if (db->field[f]->active) {
108: *any_active_fields = PETSC_TRUE;
109: return 0;
110: }
111: }
112: return 0;
113: }
115: PetscErrorCode DMSwarmDataBucketRegisterField(DMSwarmDataBucket db,const char registration_function[],const char field_name[],size_t atomic_size, DMSwarmDataField *_gfield)
116: {
117: PetscBool val;
118: DMSwarmDataField fp;
120: /* check we haven't finalised the registration of fields */
121: /*
122: if (db->finalised==PETSC_TRUE) {
123: printf("ERROR: DMSwarmDataBucketFinalize() has been called. Cannot register more fields\n");
124: ERROR();
125: }
126: */
127: /* check for repeated name */
128: DMSwarmDataFieldStringInList(field_name, db->nfields, (const DMSwarmDataField*) db->field, &val);
130: /* create new space for data */
131: PetscRealloc(sizeof(DMSwarmDataField)*(db->nfields+1), &db->field);
132: /* add field */
133: DMSwarmDataFieldCreate(registration_function, field_name, atomic_size, db->allocated, &fp);
134: db->field[db->nfields] = fp;
135: db->nfields++;
136: if (_gfield != NULL) {
137: *_gfield = fp;
138: }
139: return 0;
140: }
142: /*
143: #define DMSwarmDataBucketRegisterField(db,name,size,k) {\
144: char *location;\
145: asprintf(&location,"Registered by %s() at line %d within file %s", __FUNCTION__, __LINE__, __FILE__);\
146: _DMSwarmDataBucketRegisterField( (db), location, (name), (size), (k));\
147: PetscFree(location);\
148: }
149: */
151: PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFieldByName(DMSwarmDataBucket db,const char name[],DMSwarmDataField *gfield)
152: {
153: PetscInt idx;
154: PetscBool found;
156: DMSwarmDataFieldStringInList(name,db->nfields,(const DMSwarmDataField*)db->field,&found);
158: DMSwarmDataFieldStringFindInList(name,db->nfields,(const DMSwarmDataField*)db->field,&idx);
159: *gfield = db->field[idx];
160: return 0;
161: }
163: PetscErrorCode DMSwarmDataBucketQueryDMSwarmDataFieldByName(DMSwarmDataBucket db,const char name[],PetscBool *found)
164: {
165: *found = PETSC_FALSE;
166: DMSwarmDataFieldStringInList(name,db->nfields,(const DMSwarmDataField*)db->field,found);
167: return 0;
168: }
170: PetscErrorCode DMSwarmDataBucketFinalize(DMSwarmDataBucket db)
171: {
172: db->finalised = PETSC_TRUE;
173: return 0;
174: }
176: PetscErrorCode DMSwarmDataFieldGetNumEntries(DMSwarmDataField df,PetscInt *sum)
177: {
178: *sum = df->L;
179: return 0;
180: }
182: PetscErrorCode DMSwarmDataFieldSetBlockSize(DMSwarmDataField df,PetscInt blocksize)
183: {
184: df->bs = blocksize;
185: return 0;
186: }
188: PetscErrorCode DMSwarmDataFieldSetSize(DMSwarmDataField df,const PetscInt new_L)
189: {
191: if (new_L == df->L) return 0;
192: if (new_L > df->L) {
193: PetscRealloc(df->atomic_size * (new_L), &df->data);
194: /* init new contents */
195: PetscMemzero(( ((char*)df->data)+df->L*df->atomic_size), (new_L-df->L)*df->atomic_size);
196: } else {
197: /* reallocate pointer list, add +1 in case new_L = 0 */
198: PetscRealloc(df->atomic_size * (new_L+1), &df->data);
199: }
200: df->L = new_L;
201: return 0;
202: }
204: PetscErrorCode DMSwarmDataFieldZeroBlock(DMSwarmDataField df,const PetscInt start,const PetscInt end)
205: {
209: PetscMemzero((((char*)df->data)+start*df->atomic_size), (end-start)*df->atomic_size);
210: return 0;
211: }
213: /*
214: A negative buffer value will simply be ignored and the old buffer value will be used.
215: */
216: PetscErrorCode DMSwarmDataBucketSetSizes(DMSwarmDataBucket db,const PetscInt L,const PetscInt buffer)
217: {
218: PetscInt current_allocated,new_used,new_unused,new_buffer,new_allocated,f;
219: PetscBool any_active_fields;
222: DMSwarmDataBucketQueryForActiveFields(db,&any_active_fields);
225: current_allocated = db->allocated;
226: new_used = L;
227: new_unused = current_allocated - new_used;
228: new_buffer = db->buffer;
229: if (buffer >= 0) { /* update the buffer value */
230: new_buffer = buffer;
231: }
232: new_allocated = new_used + new_buffer;
233: /* action */
234: if (new_allocated > current_allocated) {
235: /* increase size to new_used + new_buffer */
236: for (f=0; f<db->nfields; f++) {
237: DMSwarmDataFieldSetSize(db->field[f], new_allocated);
238: }
239: db->L = new_used;
240: db->buffer = new_buffer;
241: db->allocated = new_used + new_buffer;
242: } else {
243: if (new_unused > 2 * new_buffer) {
244: /* shrink array to new_used + new_buffer */
245: for (f = 0; f < db->nfields; ++f) {
246: DMSwarmDataFieldSetSize(db->field[f], new_allocated);
247: }
248: db->L = new_used;
249: db->buffer = new_buffer;
250: db->allocated = new_used + new_buffer;
251: } else {
252: db->L = new_used;
253: db->buffer = new_buffer;
254: }
255: }
256: /* zero all entries from db->L to db->allocated */
257: for (f = 0; f < db->nfields; ++f) {
258: DMSwarmDataField field = db->field[f];
259: DMSwarmDataFieldZeroBlock(field, db->L,db->allocated);
260: }
261: return 0;
262: }
264: PetscErrorCode DMSwarmDataBucketSetInitialSizes(DMSwarmDataBucket db,const PetscInt L,const PetscInt buffer)
265: {
266: PetscInt f;
268: DMSwarmDataBucketSetSizes(db,L,buffer);
269: for (f = 0; f < db->nfields; ++f) {
270: DMSwarmDataField field = db->field[f];
271: DMSwarmDataFieldZeroBlock(field,0,db->allocated);
272: }
273: return 0;
274: }
276: PetscErrorCode DMSwarmDataBucketGetSizes(DMSwarmDataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated)
277: {
278: if (L) {*L = db->L;}
279: if (buffer) {*buffer = db->buffer;}
280: if (allocated) {*allocated = db->allocated;}
281: return 0;
282: }
284: PetscErrorCode DMSwarmDataBucketGetGlobalSizes(MPI_Comm comm,DMSwarmDataBucket db,PetscInt *L,PetscInt *buffer,PetscInt *allocated)
285: {
286: if (L) MPI_Allreduce(&db->L,L,1,MPIU_INT,MPI_SUM,comm);
287: if (buffer) MPI_Allreduce(&db->buffer,buffer,1,MPIU_INT,MPI_SUM,comm);
288: if (allocated) MPI_Allreduce(&db->allocated,allocated,1,MPIU_INT,MPI_SUM,comm);
289: return 0;
290: }
292: PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFields(DMSwarmDataBucket db,PetscInt *L,DMSwarmDataField *fields[])
293: {
294: if (L) {*L = db->nfields;}
295: if (fields) {*fields = db->field;}
296: return 0;
297: }
299: PetscErrorCode DMSwarmDataFieldGetAccess(const DMSwarmDataField gfield)
300: {
302: gfield->active = PETSC_TRUE;
303: return 0;
304: }
306: PetscErrorCode DMSwarmDataFieldAccessPoint(const DMSwarmDataField gfield,const PetscInt pid,void **ctx_p)
307: {
308: *ctx_p = NULL;
309: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
310: /* debug mode */
311: /* check point is valid */
315: #endif
316: *ctx_p = DMSWARM_DATAFIELD_point_access(gfield->data,pid,gfield->atomic_size);
317: return 0;
318: }
320: PetscErrorCode DMSwarmDataFieldAccessPointOffset(const DMSwarmDataField gfield,const size_t offset,const PetscInt pid,void **ctx_p)
321: {
322: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
323: /* debug mode */
324: /* check point is valid */
326: /* Note compiler realizes this can never happen with an unsigned PetscInt */
328: /* check point is valid */
332: #endif
333: *ctx_p = DMSWARM_DATAFIELD_point_access_offset(gfield->data,pid,gfield->atomic_size,offset);
334: return 0;
335: }
337: PetscErrorCode DMSwarmDataFieldRestoreAccess(DMSwarmDataField gfield)
338: {
340: gfield->active = PETSC_FALSE;
341: return 0;
342: }
344: PetscErrorCode DMSwarmDataFieldVerifyAccess(const DMSwarmDataField gfield,const size_t size)
345: {
346: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
348: #endif
349: return 0;
350: }
352: PetscErrorCode DMSwarmDataFieldGetAtomicSize(const DMSwarmDataField gfield,size_t *size)
353: {
354: if (size) {*size = gfield->atomic_size;}
355: return 0;
356: }
358: PetscErrorCode DMSwarmDataFieldGetEntries(const DMSwarmDataField gfield,void **data)
359: {
360: if (data) {*data = gfield->data;}
361: return 0;
362: }
364: PetscErrorCode DMSwarmDataFieldRestoreEntries(const DMSwarmDataField gfield,void **data)
365: {
366: if (data) {*data = NULL;}
367: return 0;
368: }
370: /* y = x */
371: PetscErrorCode DMSwarmDataBucketCopyPoint(const DMSwarmDataBucket xb,const PetscInt pid_x,const DMSwarmDataBucket yb,const PetscInt pid_y)
372: {
373: PetscInt f;
375: for (f = 0; f < xb->nfields; ++f) {
376: void *dest;
377: void *src;
379: DMSwarmDataFieldGetAccess(xb->field[f]);
380: if (xb != yb) DMSwarmDataFieldGetAccess( yb->field[f]);
381: DMSwarmDataFieldAccessPoint(xb->field[f],pid_x, &src);
382: DMSwarmDataFieldAccessPoint(yb->field[f],pid_y, &dest);
383: PetscMemcpy(dest, src, xb->field[f]->atomic_size);
384: DMSwarmDataFieldRestoreAccess(xb->field[f]);
385: if (xb != yb) DMSwarmDataFieldRestoreAccess(yb->field[f]);
386: }
387: return 0;
388: }
390: PetscErrorCode DMSwarmDataBucketCreateFromSubset(DMSwarmDataBucket DBIn,const PetscInt N,const PetscInt list[],DMSwarmDataBucket *DB)
391: {
392: PetscInt nfields;
393: DMSwarmDataField *fields;
394: PetscInt f,L,buffer,allocated,p;
396: DMSwarmDataBucketCreate(DB);
397: /* copy contents of DBIn */
398: DMSwarmDataBucketGetDMSwarmDataFields(DBIn,&nfields,&fields);
399: DMSwarmDataBucketGetSizes(DBIn,&L,&buffer,&allocated);
400: for (f = 0; f < nfields; ++f) {
401: DMSwarmDataBucketRegisterField(*DB,"DMSwarmDataBucketCreateFromSubset",fields[f]->name,fields[f]->atomic_size,NULL);
402: }
403: DMSwarmDataBucketFinalize(*DB);
404: DMSwarmDataBucketSetSizes(*DB,L,buffer);
405: /* now copy the desired guys from DBIn => DB */
406: for (p = 0; p < N; ++p) {
407: DMSwarmDataBucketCopyPoint(DBIn,list[p], *DB,list[p]);
408: }
409: return 0;
410: }
412: /* insert into an exisitng location */
413: PetscErrorCode DMSwarmDataFieldInsertPoint(const DMSwarmDataField field,const PetscInt index,const void *ctx)
414: {
415: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
416: /* check point is valid */
419: #endif
420: PetscMemcpy(DMSWARM_DATAFIELD_point_access(field->data,index,field->atomic_size), ctx, field->atomic_size);
421: return 0;
422: }
424: /* remove data at index - replace with last point */
425: PetscErrorCode DMSwarmDataBucketRemovePointAtIndex(const DMSwarmDataBucket db,const PetscInt index)
426: {
427: PetscInt f;
428: PetscBool any_active_fields;
430: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
431: /* check point is valid */
434: #endif
435: DMSwarmDataBucketQueryForActiveFields(db,&any_active_fields);
437: if (index >= db->L) { /* this point is not in the list - no need to error, but I will anyway */
438: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"You should not be trying to remove point at index=%D since it's < db->L = %D", index, db->L);
439: }
440: if (index != db->L-1) { /* not last point in list */
441: for (f = 0; f < db->nfields; ++f) {
442: DMSwarmDataField field = db->field[f];
444: /* copy then remove */
445: DMSwarmDataFieldCopyPoint(db->L-1, field, index, field);
446: /* DMSwarmDataFieldZeroPoint(field,index); */
447: }
448: }
449: /* decrement size */
450: /* this will zero out an crap at the end of the list */
451: DMSwarmDataBucketRemovePoint(db);
452: return 0;
453: }
455: /* copy x into y */
456: PetscErrorCode DMSwarmDataFieldCopyPoint(const PetscInt pid_x,const DMSwarmDataField field_x,const PetscInt pid_y,const DMSwarmDataField field_y)
457: {
458: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
459: /* check point is valid */
465: #endif
466: 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);
467: return 0;
468: }
470: /* zero only the datafield at this point */
471: PetscErrorCode DMSwarmDataFieldZeroPoint(const DMSwarmDataField field,const PetscInt index)
472: {
473: #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
474: /* check point is valid */
477: #endif
478: PetscMemzero(DMSWARM_DATAFIELD_point_access(field->data,index,field->atomic_size), field->atomic_size);
479: return 0;
480: }
482: /* zero ALL data for this point */
483: PetscErrorCode DMSwarmDataBucketZeroPoint(const DMSwarmDataBucket db,const PetscInt index)
484: {
485: PetscInt f;
487: /* check point is valid */
490: for (f = 0; f < db->nfields; ++f) {
491: DMSwarmDataField field = db->field[f];
492: DMSwarmDataFieldZeroPoint(field,index);
493: }
494: return 0;
495: }
497: /* increment */
498: PetscErrorCode DMSwarmDataBucketAddPoint(DMSwarmDataBucket db)
499: {
500: DMSwarmDataBucketSetSizes(db,db->L+1,DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
501: return 0;
502: }
504: /* decrement */
505: PetscErrorCode DMSwarmDataBucketRemovePoint(DMSwarmDataBucket db)
506: {
507: DMSwarmDataBucketSetSizes(db,db->L-1,DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
508: return 0;
509: }
511: /* Should be redone to user PetscViewer */
512: PetscErrorCode DMSwarmDataBucketView_stdout(MPI_Comm comm,DMSwarmDataBucket db)
513: {
514: PetscInt f;
515: double memory_usage_total,memory_usage_total_local = 0.0;
517: PetscPrintf(comm,"DMSwarmDataBucketView: \n");
518: PetscPrintf(comm," L = %D \n", db->L);
519: PetscPrintf(comm," buffer = %D \n", db->buffer);
520: PetscPrintf(comm," allocated = %D \n", db->allocated);
521: PetscPrintf(comm," nfields registered = %D \n", db->nfields);
523: for (f = 0; f < db->nfields; ++f) {
524: double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
525: memory_usage_total_local += memory_usage_f;
526: }
527: MPI_Allreduce(&memory_usage_total_local,&memory_usage_total,1,MPI_DOUBLE,MPI_SUM,comm);
529: for (f = 0; f < db->nfields; ++f) {
530: double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
531: PetscPrintf(comm," [%3D] %15s : Mem. usage = %1.2e (MB) [rank0]\n", f, db->field[f]->name, memory_usage_f);
532: PetscPrintf(comm," blocksize = %D \n", db->field[f]->bs);
533: if (db->field[f]->bs != 1) {
534: PetscPrintf(comm," atomic size = %zu [full block, bs=%D]\n", db->field[f]->atomic_size,db->field[f]->bs);
535: PetscPrintf(comm," atomic size/item = %zu \n", db->field[f]->atomic_size/db->field[f]->bs);
536: } else {
537: PetscPrintf(comm," atomic size = %zu \n", db->field[f]->atomic_size);
538: }
539: }
540: PetscPrintf(comm," Total mem. usage = %1.2e (MB) (collective)\n", memory_usage_total);
541: return 0;
542: }
544: PetscErrorCode DMSwarmDataBucketView_Seq(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)
545: {
546: switch (type) {
547: case DATABUCKET_VIEW_STDOUT:
548: DMSwarmDataBucketView_stdout(PETSC_COMM_SELF,db);
549: break;
550: case DATABUCKET_VIEW_ASCII:
551: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for ascii output");
552: case DATABUCKET_VIEW_BINARY:
553: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for binary output");
554: case DATABUCKET_VIEW_HDF5:
555: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for HDF5 output");
556: default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unknown viewer method requested");
557: }
558: return 0;
559: }
561: PetscErrorCode DMSwarmDataBucketView_MPI(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)
562: {
563: switch (type) {
564: case DATABUCKET_VIEW_STDOUT:
565: DMSwarmDataBucketView_stdout(comm,db);
566: break;
567: case DATABUCKET_VIEW_ASCII:
568: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for ascii output");
569: case DATABUCKET_VIEW_BINARY:
570: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for binary output");
571: case DATABUCKET_VIEW_HDF5:
572: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for HDF5 output");
573: default: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unknown viewer method requested");
574: }
575: return 0;
576: }
578: PetscErrorCode DMSwarmDataBucketView(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)
579: {
580: PetscMPIInt size;
582: MPI_Comm_size(comm,&size);
583: if (size == 1) {
584: DMSwarmDataBucketView_Seq(comm,db,filename,type);
585: } else {
586: DMSwarmDataBucketView_MPI(comm,db,filename,type);
587: }
588: return 0;
589: }
591: PetscErrorCode DMSwarmDataBucketDuplicateFields(DMSwarmDataBucket dbA,DMSwarmDataBucket *dbB)
592: {
593: DMSwarmDataBucket db2;
594: PetscInt f;
596: DMSwarmDataBucketCreate(&db2);
597: /* copy contents from dbA into db2 */
598: for (f = 0; f < dbA->nfields; ++f) {
599: DMSwarmDataField field;
600: size_t atomic_size;
601: char *name;
603: field = dbA->field[f];
604: atomic_size = field->atomic_size;
605: name = field->name;
606: DMSwarmDataBucketRegisterField(db2,"DMSwarmDataBucketDuplicateFields",name,atomic_size,NULL);
607: }
608: DMSwarmDataBucketFinalize(db2);
609: DMSwarmDataBucketSetInitialSizes(db2,0,1000);
610: *dbB = db2;
611: return 0;
612: }
614: /*
615: Insert points from db2 into db1
616: db1 <<== db2
617: */
618: PetscErrorCode DMSwarmDataBucketInsertValues(DMSwarmDataBucket db1,DMSwarmDataBucket db2)
619: {
620: PetscInt n_mp_points1,n_mp_points2;
621: PetscInt n_mp_points1_new,p;
623: DMSwarmDataBucketGetSizes(db1,&n_mp_points1,NULL,NULL);
624: DMSwarmDataBucketGetSizes(db2,&n_mp_points2,NULL,NULL);
625: n_mp_points1_new = n_mp_points1 + n_mp_points2;
626: DMSwarmDataBucketSetSizes(db1,n_mp_points1_new,DMSWARM_DATA_BUCKET_BUFFER_DEFAULT);
627: for (p = 0; p < n_mp_points2; ++p) {
628: /* db1 <<== db2 */
629: DMSwarmDataBucketCopyPoint(db2,p, db1,(n_mp_points1 + p));
630: }
631: return 0;
632: }
634: /* helpers for parallel send/recv */
635: PetscErrorCode DMSwarmDataBucketCreatePackedArray(DMSwarmDataBucket db,size_t *bytes,void **buf)
636: {
637: PetscInt f;
638: size_t sizeof_marker_contents;
639: void *buffer;
641: sizeof_marker_contents = 0;
642: for (f = 0; f < db->nfields; ++f) {
643: DMSwarmDataField df = db->field[f];
644: sizeof_marker_contents += df->atomic_size;
645: }
646: PetscMalloc(sizeof_marker_contents, &buffer);
647: PetscMemzero(buffer, sizeof_marker_contents);
648: if (bytes) {*bytes = sizeof_marker_contents;}
649: if (buf) {*buf = buffer;}
650: return 0;
651: }
653: PetscErrorCode DMSwarmDataBucketDestroyPackedArray(DMSwarmDataBucket db,void **buf)
654: {
655: if (buf) {
656: PetscFree(*buf);
657: *buf = NULL;
658: }
659: return 0;
660: }
662: PetscErrorCode DMSwarmDataBucketFillPackedArray(DMSwarmDataBucket db,const PetscInt index,void *buf)
663: {
664: PetscInt f;
665: void *data, *data_p;
666: size_t asize, offset;
668: offset = 0;
669: for (f = 0; f < db->nfields; ++f) {
670: DMSwarmDataField df = db->field[f];
672: asize = df->atomic_size;
673: data = (void*)( df->data);
674: data_p = (void*)( (char*)data + index*asize);
675: PetscMemcpy((void*)((char*)buf + offset), data_p, asize);
676: offset = offset + asize;
677: }
678: return 0;
679: }
681: PetscErrorCode DMSwarmDataBucketInsertPackedArray(DMSwarmDataBucket db,const PetscInt idx,void *data)
682: {
683: PetscInt f;
684: void *data_p;
685: size_t offset;
687: offset = 0;
688: for (f = 0; f < db->nfields; ++f) {
689: DMSwarmDataField df = db->field[f];
691: data_p = (void*)( (char*)data + offset);
692: DMSwarmDataFieldInsertPoint(df, idx, (void*)data_p);
693: offset = offset + df->atomic_size;
694: }
695: return 0;
696: }