1 #include "../src/dm/impls/swarm/data_bucket.h"
2
3 /* string helpers */
DMSwarmDataFieldStringInList(const char name[],const PetscInt N,const DMSwarmDataField gfield[],PetscBool * val)4 PetscErrorCode DMSwarmDataFieldStringInList(const char name[], const PetscInt N, const DMSwarmDataField gfield[], PetscBool *val)
5 {
6 PetscInt i;
7
8 PetscFunctionBegin;
9 *val = PETSC_FALSE;
10 for (i = 0; i < N; ++i) {
11 PetscBool flg;
12 PetscCall(PetscStrcmp(name, gfield[i]->name, &flg));
13 if (flg) {
14 *val = PETSC_TRUE;
15 PetscFunctionReturn(PETSC_SUCCESS);
16 }
17 }
18 PetscFunctionReturn(PETSC_SUCCESS);
19 }
20
DMSwarmDataFieldStringFindInList(const char name[],const PetscInt N,const DMSwarmDataField gfield[],PetscInt * index)21 PetscErrorCode DMSwarmDataFieldStringFindInList(const char name[], const PetscInt N, const DMSwarmDataField gfield[], PetscInt *index)
22 {
23 PetscInt i;
24
25 PetscFunctionBegin;
26 *index = -1;
27 for (i = 0; i < N; ++i) {
28 PetscBool flg;
29 PetscCall(PetscStrcmp(name, gfield[i]->name, &flg));
30 if (flg) {
31 *index = i;
32 PetscFunctionReturn(PETSC_SUCCESS);
33 }
34 }
35 PetscFunctionReturn(PETSC_SUCCESS);
36 }
37
DMSwarmDataFieldCreate(const char registration_function[],const char name[],const size_t size,const PetscInt L,DMSwarmDataField * DF)38 PetscErrorCode DMSwarmDataFieldCreate(const char registration_function[], const char name[], const size_t size, const PetscInt L, DMSwarmDataField *DF)
39 {
40 DMSwarmDataField df;
41
42 PetscFunctionBegin;
43 PetscCall(PetscNew(&df));
44 PetscCall(PetscStrallocpy(registration_function, &df->registration_function));
45 PetscCall(PetscStrallocpy(name, &df->name));
46 df->atomic_size = size;
47 df->L = L;
48 df->bs = 1;
49 /* allocate something so we don't have to reallocate */
50 PetscCall(PetscMalloc(size * L, &df->data));
51 PetscCall(PetscMemzero(df->data, size * L));
52 *DF = df;
53 PetscFunctionReturn(PETSC_SUCCESS);
54 }
55
DMSwarmDataFieldDestroy(DMSwarmDataField * DF)56 PetscErrorCode DMSwarmDataFieldDestroy(DMSwarmDataField *DF)
57 {
58 DMSwarmDataField df = *DF;
59
60 PetscFunctionBegin;
61 PetscCall(PetscFree(df->registration_function));
62 PetscCall(PetscFree(df->name));
63 PetscCall(PetscFree(df->data));
64 PetscCall(PetscFree(df));
65 *DF = NULL;
66 PetscFunctionReturn(PETSC_SUCCESS);
67 }
68
69 /* data bucket */
DMSwarmDataBucketCreate(DMSwarmDataBucket * DB)70 PetscErrorCode DMSwarmDataBucketCreate(DMSwarmDataBucket *DB)
71 {
72 DMSwarmDataBucket db;
73
74 PetscFunctionBegin;
75 PetscCall(PetscNew(&db));
76
77 db->finalised = PETSC_FALSE;
78 /* create empty spaces for fields */
79 db->L = -1;
80 db->buffer = 1;
81 db->allocated = 1;
82 db->nfields = 0;
83 PetscCall(PetscMalloc1(1, &db->field));
84 *DB = db;
85 PetscFunctionReturn(PETSC_SUCCESS);
86 }
87
DMSwarmDataBucketDestroy(DMSwarmDataBucket * DB)88 PetscErrorCode DMSwarmDataBucketDestroy(DMSwarmDataBucket *DB)
89 {
90 DMSwarmDataBucket db = *DB;
91 PetscInt f;
92
93 PetscFunctionBegin;
94 /* release fields */
95 for (f = 0; f < db->nfields; ++f) PetscCall(DMSwarmDataFieldDestroy(&db->field[f]));
96 /* this will catch the initially allocated objects in the event that no fields are registered */
97 if (db->field != NULL) PetscCall(PetscFree(db->field));
98 PetscCall(PetscFree(db));
99 *DB = NULL;
100 PetscFunctionReturn(PETSC_SUCCESS);
101 }
102
DMSwarmDataBucketQueryForActiveFields(DMSwarmDataBucket db,PetscBool * any_active_fields)103 PetscErrorCode DMSwarmDataBucketQueryForActiveFields(DMSwarmDataBucket db, PetscBool *any_active_fields)
104 {
105 PetscInt f;
106
107 PetscFunctionBegin;
108 *any_active_fields = PETSC_FALSE;
109 for (f = 0; f < db->nfields; ++f) {
110 if (db->field[f]->active) {
111 *any_active_fields = PETSC_TRUE;
112 PetscFunctionReturn(PETSC_SUCCESS);
113 }
114 }
115 PetscFunctionReturn(PETSC_SUCCESS);
116 }
117
DMSwarmDataBucketRegisterField(DMSwarmDataBucket db,const char registration_function[],const char field_name[],size_t atomic_size,DMSwarmDataField * _gfield)118 PetscErrorCode DMSwarmDataBucketRegisterField(DMSwarmDataBucket db, const char registration_function[], const char field_name[], size_t atomic_size, DMSwarmDataField *_gfield)
119 {
120 PetscBool val;
121 DMSwarmDataField fp;
122
123 PetscFunctionBegin;
124 /* check we haven't finalised the registration of fields */
125 /*
126 if (db->finalised==PETSC_TRUE) {
127 printf("ERROR: DMSwarmDataBucketFinalize() has been called. Cannot register more fields\n");
128 ERROR();
129 }
130 */
131 /* check for repeated name */
132 PetscCall(DMSwarmDataFieldStringInList(field_name, db->nfields, (const DMSwarmDataField *)db->field, &val));
133 PetscCheck(val != PETSC_TRUE, PETSC_COMM_SELF, PETSC_ERR_USER, "Field %s already exists. Cannot add same field twice", field_name);
134 /* create new space for data */
135 PetscCall(PetscRealloc(sizeof(DMSwarmDataField) * (db->nfields + 1), &db->field));
136 /* add field */
137 PetscCall(DMSwarmDataFieldCreate(registration_function, field_name, atomic_size, db->allocated, &fp));
138 db->field[db->nfields] = fp;
139 db->nfields++;
140 if (_gfield != NULL) *_gfield = fp;
141 PetscFunctionReturn(PETSC_SUCCESS);
142 }
143
144 /*
145 #define DMSwarmDataBucketRegisterField(db,name,size,k) {\
146 char *location;\
147 asprintf(&location,"Registered by %s() at line %d within file %s", __FUNCTION__, __LINE__, __FILE__);\
148 _DMSwarmDataBucketRegisterField( (db), location, (name), (size), (k));\
149 ierr = PetscFree(location);\
150 }
151 */
152
DMSwarmDataBucketGetDMSwarmDataFieldIdByName(DMSwarmDataBucket db,const char name[],PetscInt * idx)153 PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFieldIdByName(DMSwarmDataBucket db, const char name[], PetscInt *idx)
154 {
155 PetscBool found;
156
157 PetscFunctionBegin;
158 *idx = -1;
159 PetscCall(DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, &found));
160 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot find DMSwarmDataField with name %s", name);
161 PetscCall(DMSwarmDataFieldStringFindInList(name, db->nfields, (const DMSwarmDataField *)db->field, idx));
162 PetscFunctionReturn(PETSC_SUCCESS);
163 }
164
DMSwarmDataBucketGetDMSwarmDataFieldByName(DMSwarmDataBucket db,const char name[],DMSwarmDataField * gfield)165 PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFieldByName(DMSwarmDataBucket db, const char name[], DMSwarmDataField *gfield)
166 {
167 PetscInt idx;
168 PetscBool found;
169
170 PetscFunctionBegin;
171 PetscCall(DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, &found));
172 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot find DMSwarmDataField with name %s", name);
173 PetscCall(DMSwarmDataFieldStringFindInList(name, db->nfields, (const DMSwarmDataField *)db->field, &idx));
174 *gfield = db->field[idx];
175 PetscFunctionReturn(PETSC_SUCCESS);
176 }
177
DMSwarmDataBucketQueryDMSwarmDataFieldByName(DMSwarmDataBucket db,const char name[],PetscBool * found)178 PetscErrorCode DMSwarmDataBucketQueryDMSwarmDataFieldByName(DMSwarmDataBucket db, const char name[], PetscBool *found)
179 {
180 PetscFunctionBegin;
181 *found = PETSC_FALSE;
182 PetscCall(DMSwarmDataFieldStringInList(name, db->nfields, (const DMSwarmDataField *)db->field, found));
183 PetscFunctionReturn(PETSC_SUCCESS);
184 }
185
DMSwarmDataBucketFinalize(DMSwarmDataBucket db)186 PetscErrorCode DMSwarmDataBucketFinalize(DMSwarmDataBucket db)
187 {
188 PetscFunctionBegin;
189 db->finalised = PETSC_TRUE;
190 PetscFunctionReturn(PETSC_SUCCESS);
191 }
192
DMSwarmDataFieldGetNumEntries(DMSwarmDataField df,PetscInt * sum)193 PetscErrorCode DMSwarmDataFieldGetNumEntries(DMSwarmDataField df, PetscInt *sum)
194 {
195 PetscFunctionBegin;
196 *sum = df->L;
197 PetscFunctionReturn(PETSC_SUCCESS);
198 }
199
DMSwarmDataFieldSetBlockSize(DMSwarmDataField df,PetscInt blocksize)200 PetscErrorCode DMSwarmDataFieldSetBlockSize(DMSwarmDataField df, PetscInt blocksize)
201 {
202 PetscFunctionBegin;
203 df->bs = blocksize;
204 PetscFunctionReturn(PETSC_SUCCESS);
205 }
206
DMSwarmDataFieldSetSize(DMSwarmDataField df,const PetscInt new_L)207 PetscErrorCode DMSwarmDataFieldSetSize(DMSwarmDataField df, const PetscInt new_L)
208 {
209 PetscFunctionBegin;
210 PetscCheck(new_L >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot set size of DMSwarmDataField to be < 0");
211 if (new_L == df->L) PetscFunctionReturn(PETSC_SUCCESS);
212 if (new_L > df->L) {
213 PetscCall(PetscRealloc(df->atomic_size * (new_L), &df->data));
214 /* init new contents */
215 PetscCall(PetscMemzero(((char *)df->data) + df->L * df->atomic_size, (new_L - df->L) * df->atomic_size));
216 } else {
217 /* reallocate pointer list, add +1 in case new_L = 0 */
218 PetscCall(PetscRealloc(df->atomic_size * (new_L + 1), &df->data));
219 }
220 df->L = new_L;
221 PetscFunctionReturn(PETSC_SUCCESS);
222 }
223
DMSwarmDataFieldZeroBlock(DMSwarmDataField df,const PetscInt start,const PetscInt end)224 PetscErrorCode DMSwarmDataFieldZeroBlock(DMSwarmDataField df, const PetscInt start, const PetscInt end)
225 {
226 PetscFunctionBegin;
227 PetscCheck(start <= end, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot zero a block of entries if start(%" PetscInt_FMT ") > end(%" PetscInt_FMT ")", start, end);
228 PetscCheck(start >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot zero a block of entries if start(%" PetscInt_FMT ") < 0", start);
229 PetscCheck(end <= df->L, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot zero a block of entries if end(%" PetscInt_FMT ") >= array size(%" PetscInt_FMT ")", end, df->L);
230 PetscCall(PetscMemzero(((char *)df->data) + start * df->atomic_size, (end - start) * df->atomic_size));
231 PetscFunctionReturn(PETSC_SUCCESS);
232 }
233
234 /*
235 A negative buffer value will simply be ignored and the old buffer value will be used.
236 */
DMSwarmDataBucketSetSizes(DMSwarmDataBucket db,const PetscInt L,const PetscInt buffer)237 PetscErrorCode DMSwarmDataBucketSetSizes(DMSwarmDataBucket db, const PetscInt L, const PetscInt buffer)
238 {
239 PetscInt current_allocated, current_used, new_used, new_unused, new_buffer, new_allocated, f, end;
240 PetscBool any_active_fields;
241
242 PetscFunctionBegin;
243 PetscCheck(db->finalised != PETSC_FALSE, PETSC_COMM_SELF, PETSC_ERR_USER, "You must call DMSwarmDataBucketFinalize() before DMSwarmDataBucketSetSizes()");
244 PetscCall(DMSwarmDataBucketQueryForActiveFields(db, &any_active_fields));
245 PetscCheck(!any_active_fields, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot safely re-size as at least one DMSwarmDataField is currently being accessed");
246
247 current_allocated = db->allocated;
248 current_used = PetscMax(db->L, 0);
249 new_used = L;
250 new_unused = current_allocated - new_used;
251 new_buffer = db->buffer;
252 if (buffer >= 0) { /* update the buffer value */
253 new_buffer = buffer;
254 }
255 new_allocated = new_used + new_buffer;
256 /* action */
257 if (new_allocated > current_allocated) {
258 /* increase size to new_used + new_buffer and zero new space */
259 for (f = 0; f < db->nfields; f++) {
260 PetscCall(DMSwarmDataFieldSetSize(db->field[f], new_allocated));
261 PetscCall(DMSwarmDataFieldZeroBlock(db->field[f], current_allocated, new_allocated));
262 }
263 db->L = new_used;
264 db->buffer = new_buffer;
265 db->allocated = new_used + new_buffer;
266 } else {
267 if (new_unused > 2 * new_buffer) {
268 /* shrink array to new_used + new_buffer */
269 for (f = 0; f < db->nfields; ++f) PetscCall(DMSwarmDataFieldSetSize(db->field[f], new_allocated));
270 db->L = new_used;
271 db->buffer = new_buffer;
272 db->allocated = new_used + new_buffer;
273 } else {
274 db->L = new_used;
275 db->buffer = new_buffer;
276 }
277 }
278 /* if we shrunk, zero old entries from new_used to current_used or end of array */
279 end = PetscMin(current_used, new_allocated);
280 if (end > new_used) {
281 for (f = 0; f < db->nfields; ++f) {
282 DMSwarmDataField field = db->field[f];
283 PetscCall(DMSwarmDataFieldZeroBlock(field, new_used, end));
284 }
285 }
286 PetscFunctionReturn(PETSC_SUCCESS);
287 }
288
DMSwarmDataBucketSetInitialSizes(DMSwarmDataBucket db,const PetscInt L,const PetscInt buffer)289 PetscErrorCode DMSwarmDataBucketSetInitialSizes(DMSwarmDataBucket db, const PetscInt L, const PetscInt buffer)
290 {
291 PetscInt f;
292
293 PetscFunctionBegin;
294 PetscCall(DMSwarmDataBucketSetSizes(db, L, buffer));
295 for (f = 0; f < db->nfields; ++f) {
296 DMSwarmDataField field = db->field[f];
297 PetscCall(DMSwarmDataFieldZeroBlock(field, 0, db->allocated));
298 }
299 PetscFunctionReturn(PETSC_SUCCESS);
300 }
301
DMSwarmDataBucketGetSizes(DMSwarmDataBucket db,PetscInt * L,PetscInt * buffer,PetscInt * allocated)302 PetscErrorCode DMSwarmDataBucketGetSizes(DMSwarmDataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated)
303 {
304 PetscFunctionBegin;
305 if (L) *L = db->L;
306 if (buffer) *buffer = db->buffer;
307 if (allocated) *allocated = db->allocated;
308 PetscFunctionReturn(PETSC_SUCCESS);
309 }
310
DMSwarmDataBucketGetGlobalSizes(MPI_Comm comm,DMSwarmDataBucket db,PetscInt * L,PetscInt * buffer,PetscInt * allocated)311 PetscErrorCode DMSwarmDataBucketGetGlobalSizes(MPI_Comm comm, DMSwarmDataBucket db, PetscInt *L, PetscInt *buffer, PetscInt *allocated)
312 {
313 PetscFunctionBegin;
314 if (L) PetscCallMPI(MPIU_Allreduce(&db->L, L, 1, MPIU_INT, MPI_SUM, comm));
315 if (buffer) PetscCallMPI(MPIU_Allreduce(&db->buffer, buffer, 1, MPIU_INT, MPI_SUM, comm));
316 if (allocated) PetscCallMPI(MPIU_Allreduce(&db->allocated, allocated, 1, MPIU_INT, MPI_SUM, comm));
317 PetscFunctionReturn(PETSC_SUCCESS);
318 }
319
DMSwarmDataBucketGetDMSwarmDataFields(DMSwarmDataBucket db,PetscInt * L,DMSwarmDataField * fields[])320 PetscErrorCode DMSwarmDataBucketGetDMSwarmDataFields(DMSwarmDataBucket db, PetscInt *L, DMSwarmDataField *fields[])
321 {
322 PetscFunctionBegin;
323 if (L) *L = db->nfields;
324 if (fields) *fields = db->field;
325 PetscFunctionReturn(PETSC_SUCCESS);
326 }
327
DMSwarmDataFieldGetAccess(const DMSwarmDataField gfield)328 PetscErrorCode DMSwarmDataFieldGetAccess(const DMSwarmDataField gfield)
329 {
330 PetscFunctionBegin;
331 PetscCheck(!gfield->active, PETSC_COMM_SELF, PETSC_ERR_USER, "Field \"%s\" is already active. You must call DMSwarmDataFieldRestoreAccess()", gfield->name);
332 gfield->active = PETSC_TRUE;
333 PetscFunctionReturn(PETSC_SUCCESS);
334 }
335
DMSwarmDataFieldAccessPoint(const DMSwarmDataField gfield,const PetscInt pid,void ** ctx_p)336 PetscErrorCode DMSwarmDataFieldAccessPoint(const DMSwarmDataField gfield, const PetscInt pid, void **ctx_p)
337 {
338 PetscFunctionBegin;
339 *ctx_p = NULL;
340 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
341 /* debug mode */
342 /* check point is valid */
343 PetscCheck(pid >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
344 PetscCheck(pid < gfield->L, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, gfield->L);
345 PetscCheck(gfield->active != PETSC_FALSE, PETSC_COMM_SELF, PETSC_ERR_USER, "Field \"%s\" is not active. You must call DMSwarmDataFieldGetAccess() before point data can be retrivied", gfield->name);
346 #endif
347 *ctx_p = DMSWARM_DATAFIELD_point_access(gfield->data, pid, gfield->atomic_size);
348 PetscFunctionReturn(PETSC_SUCCESS);
349 }
350
DMSwarmDataFieldAccessPointOffset(const DMSwarmDataField gfield,const size_t offset,const PetscInt pid,void ** ctx_p)351 PetscErrorCode DMSwarmDataFieldAccessPointOffset(const DMSwarmDataField gfield, const size_t offset, const PetscInt pid, void **ctx_p)
352 {
353 PetscFunctionBegin;
354 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
355 /* debug mode */
356 /* check point is valid */
357 /* PetscCheck(offset >= 0,PETSC_COMM_SELF,PETSC_ERR_USER,"offset must be >= 0");*/
358 /* Note compiler realizes this can never happen with an unsigned PetscInt */
359 PetscCheck(offset < gfield->atomic_size, PETSC_COMM_SELF, PETSC_ERR_USER, "offset must be < %zu", gfield->atomic_size);
360 /* check point is valid */
361 PetscCheck(pid >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
362 PetscCheck(pid < gfield->L, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, gfield->L);
363 PetscCheck(gfield->active != PETSC_FALSE, PETSC_COMM_SELF, PETSC_ERR_USER, "Field \"%s\" is not active. You must call DMSwarmDataFieldGetAccess() before point data can be retrivied", gfield->name);
364 #endif
365 *ctx_p = DMSWARM_DATAFIELD_point_access_offset(gfield->data, pid, gfield->atomic_size, offset);
366 PetscFunctionReturn(PETSC_SUCCESS);
367 }
368
DMSwarmDataFieldRestoreAccess(DMSwarmDataField gfield)369 PetscErrorCode DMSwarmDataFieldRestoreAccess(DMSwarmDataField gfield)
370 {
371 PetscFunctionBegin;
372 PetscCheck(gfield->active != PETSC_FALSE, PETSC_COMM_SELF, PETSC_ERR_USER, "Field \"%s\" is not active. You must call DMSwarmDataFieldGetAccess()", gfield->name);
373 gfield->active = PETSC_FALSE;
374 PetscFunctionReturn(PETSC_SUCCESS);
375 }
376
DMSwarmDataFieldVerifyAccess(const DMSwarmDataField gfield,const size_t size)377 PetscErrorCode DMSwarmDataFieldVerifyAccess(const DMSwarmDataField gfield, const size_t size)
378 {
379 PetscFunctionBegin;
380 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
381 PetscCheck(gfield->atomic_size == size, PETSC_COMM_SELF, PETSC_ERR_USER, "Field \"%s\" must be mapped to %zu bytes, your intended structure is %zu bytes in length.", gfield->name, gfield->atomic_size, size);
382 #endif
383 PetscFunctionReturn(PETSC_SUCCESS);
384 }
385
DMSwarmDataFieldGetAtomicSize(const DMSwarmDataField gfield,size_t * size)386 PetscErrorCode DMSwarmDataFieldGetAtomicSize(const DMSwarmDataField gfield, size_t *size)
387 {
388 PetscFunctionBegin;
389 if (size) *size = gfield->atomic_size;
390 PetscFunctionReturn(PETSC_SUCCESS);
391 }
392
DMSwarmDataFieldGetEntries(const DMSwarmDataField gfield,void ** data)393 PetscErrorCode DMSwarmDataFieldGetEntries(const DMSwarmDataField gfield, void **data)
394 {
395 PetscFunctionBegin;
396 if (data) *data = gfield->data;
397 PetscFunctionReturn(PETSC_SUCCESS);
398 }
399
DMSwarmDataFieldRestoreEntries(const DMSwarmDataField gfield,void ** data)400 PetscErrorCode DMSwarmDataFieldRestoreEntries(const DMSwarmDataField gfield, void **data)
401 {
402 PetscFunctionBegin;
403 if (data) *data = NULL;
404 PetscFunctionReturn(PETSC_SUCCESS);
405 }
406
407 /* y = x */
DMSwarmDataBucketCopyPoint(const DMSwarmDataBucket xb,const PetscInt pid_x,const DMSwarmDataBucket yb,const PetscInt pid_y)408 PetscErrorCode DMSwarmDataBucketCopyPoint(const DMSwarmDataBucket xb, const PetscInt pid_x, const DMSwarmDataBucket yb, const PetscInt pid_y)
409 {
410 PetscInt f;
411
412 PetscFunctionBegin;
413 for (f = 0; f < xb->nfields; ++f) {
414 void *dest;
415 void *src;
416
417 PetscCall(DMSwarmDataFieldGetAccess(xb->field[f]));
418 if (xb != yb) PetscCall(DMSwarmDataFieldGetAccess(yb->field[f]));
419 PetscCall(DMSwarmDataFieldAccessPoint(xb->field[f], pid_x, &src));
420 PetscCall(DMSwarmDataFieldAccessPoint(yb->field[f], pid_y, &dest));
421 PetscCall(PetscMemcpy(dest, src, xb->field[f]->atomic_size));
422 PetscCall(DMSwarmDataFieldRestoreAccess(xb->field[f]));
423 if (xb != yb) PetscCall(DMSwarmDataFieldRestoreAccess(yb->field[f]));
424 }
425 PetscFunctionReturn(PETSC_SUCCESS);
426 }
427
DMSwarmDataBucketCreateFromSubset(DMSwarmDataBucket DBIn,const PetscInt N,const PetscInt list[],DMSwarmDataBucket * DB)428 PetscErrorCode DMSwarmDataBucketCreateFromSubset(DMSwarmDataBucket DBIn, const PetscInt N, const PetscInt list[], DMSwarmDataBucket *DB)
429 {
430 PetscInt nfields;
431 DMSwarmDataField *fields;
432 PetscInt f, L, buffer, allocated, p;
433
434 PetscFunctionBegin;
435 PetscCall(DMSwarmDataBucketCreate(DB));
436 /* copy contents of DBIn */
437 PetscCall(DMSwarmDataBucketGetDMSwarmDataFields(DBIn, &nfields, &fields));
438 PetscCall(DMSwarmDataBucketGetSizes(DBIn, &L, &buffer, &allocated));
439 for (f = 0; f < nfields; ++f) PetscCall(DMSwarmDataBucketRegisterField(*DB, "DMSwarmDataBucketCreateFromSubset", fields[f]->name, fields[f]->atomic_size, NULL));
440 PetscCall(DMSwarmDataBucketFinalize(*DB));
441 PetscCall(DMSwarmDataBucketSetSizes(*DB, L, buffer));
442 for (f = 0; f < nfields; ++f) {
443 DMSwarmDataField gfield;
444
445 PetscCall(DMSwarmDataBucketGetDMSwarmDataFieldByName(*DB, fields[f]->name, &gfield));
446 PetscCall(DMSwarmDataFieldSetBlockSize(gfield, fields[f]->bs));
447 gfield->petsc_type = fields[f]->petsc_type;
448 }
449 /* now copy the desired guys from DBIn => DB */
450 for (p = 0; p < N; ++p) PetscCall(DMSwarmDataBucketCopyPoint(DBIn, list[p], *DB, list[p]));
451 PetscFunctionReturn(PETSC_SUCCESS);
452 }
453
454 /* insert into an existing location */
DMSwarmDataFieldInsertPoint(const DMSwarmDataField field,const PetscInt index,const void * data)455 PetscErrorCode DMSwarmDataFieldInsertPoint(const DMSwarmDataField field, const PetscInt index, const void *data)
456 {
457 PetscFunctionBegin;
458 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
459 /* check point is valid */
460 PetscCheck(index >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
461 PetscCheck(index < field->L, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, field->L);
462 #endif
463 PetscCall(PetscMemcpy(DMSWARM_DATAFIELD_point_access(field->data, index, field->atomic_size), data, field->atomic_size));
464 PetscFunctionReturn(PETSC_SUCCESS);
465 }
466
467 /* remove data at index - replace with last point */
DMSwarmDataBucketRemovePointAtIndex(const DMSwarmDataBucket db,const PetscInt index)468 PetscErrorCode DMSwarmDataBucketRemovePointAtIndex(const DMSwarmDataBucket db, const PetscInt index)
469 {
470 PetscInt f;
471 PetscBool any_active_fields;
472
473 PetscFunctionBegin;
474 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
475 /* check point is valid */
476 PetscCheck(index >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
477 PetscCheck(index < db->allocated, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, db->L + db->buffer);
478 #endif
479 PetscCall(DMSwarmDataBucketQueryForActiveFields(db, &any_active_fields));
480 PetscCheck(!any_active_fields, PETSC_COMM_SELF, PETSC_ERR_USER, "Cannot safely remove point as at least one DMSwarmDataField is currently being accessed");
481 if (index >= db->L) { /* this point is not in the list - no need to error, but I will anyway */
482 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);
483 }
484 if (index != db->L - 1) { /* not last point in list */
485 for (f = 0; f < db->nfields; ++f) {
486 DMSwarmDataField field = db->field[f];
487
488 /* copy then remove */
489 PetscCall(DMSwarmDataFieldCopyPoint(db->L - 1, field, index, field));
490 /* DMSwarmDataFieldZeroPoint(field,index); */
491 }
492 }
493 /* decrement size */
494 /* this will zero out an crap at the end of the list */
495 PetscCall(DMSwarmDataBucketRemovePoint(db));
496 PetscFunctionReturn(PETSC_SUCCESS);
497 }
498
499 /* copy x into y */
DMSwarmDataFieldCopyPoint(const PetscInt pid_x,const DMSwarmDataField field_x,const PetscInt pid_y,const DMSwarmDataField field_y)500 PetscErrorCode DMSwarmDataFieldCopyPoint(const PetscInt pid_x, const DMSwarmDataField field_x, const PetscInt pid_y, const DMSwarmDataField field_y)
501 {
502 PetscFunctionBegin;
503 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
504 /* check point is valid */
505 PetscCheck(pid_x >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "(IN) index must be >= 0");
506 PetscCheck(pid_x < field_x->L, PETSC_COMM_SELF, PETSC_ERR_USER, "(IN) index must be < %" PetscInt_FMT, field_x->L);
507 PetscCheck(pid_y >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "(OUT) index must be >= 0");
508 PetscCheck(pid_y < field_y->L, PETSC_COMM_SELF, PETSC_ERR_USER, "(OUT) index must be < %" PetscInt_FMT, field_y->L);
509 PetscCheck(field_y->atomic_size == field_x->atomic_size, PETSC_COMM_SELF, PETSC_ERR_USER, "atomic size must match");
510 #endif
511 PetscCall(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));
512 PetscFunctionReturn(PETSC_SUCCESS);
513 }
514
515 /* zero only the datafield at this point */
DMSwarmDataFieldZeroPoint(const DMSwarmDataField field,const PetscInt index)516 PetscErrorCode DMSwarmDataFieldZeroPoint(const DMSwarmDataField field, const PetscInt index)
517 {
518 PetscFunctionBegin;
519 #if defined(DMSWARM_DATAFIELD_POINT_ACCESS_GUARD)
520 /* check point is valid */
521 PetscCheck(index >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
522 PetscCheck(index < field->L, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, field->L);
523 #endif
524 PetscCall(PetscMemzero(DMSWARM_DATAFIELD_point_access(field->data, index, field->atomic_size), field->atomic_size));
525 PetscFunctionReturn(PETSC_SUCCESS);
526 }
527
528 /* zero ALL data for this point */
DMSwarmDataBucketZeroPoint(const DMSwarmDataBucket db,const PetscInt index)529 PetscErrorCode DMSwarmDataBucketZeroPoint(const DMSwarmDataBucket db, const PetscInt index)
530 {
531 PetscInt f;
532
533 PetscFunctionBegin;
534 /* check point is valid */
535 PetscCheck(index >= 0, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be >= 0");
536 PetscCheck(index < db->allocated, PETSC_COMM_SELF, PETSC_ERR_USER, "index must be < %" PetscInt_FMT, db->allocated);
537 for (f = 0; f < db->nfields; ++f) {
538 DMSwarmDataField field = db->field[f];
539 PetscCall(DMSwarmDataFieldZeroPoint(field, index));
540 }
541 PetscFunctionReturn(PETSC_SUCCESS);
542 }
543
544 /* increment */
DMSwarmDataBucketAddPoint(DMSwarmDataBucket db)545 PetscErrorCode DMSwarmDataBucketAddPoint(DMSwarmDataBucket db)
546 {
547 PetscFunctionBegin;
548 PetscCall(DMSwarmDataBucketSetSizes(db, PetscMax(db->L, 0) + 1, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT));
549 PetscFunctionReturn(PETSC_SUCCESS);
550 }
551
552 /* decrement */
DMSwarmDataBucketRemovePoint(DMSwarmDataBucket db)553 PetscErrorCode DMSwarmDataBucketRemovePoint(DMSwarmDataBucket db)
554 {
555 PetscFunctionBegin;
556 PetscCheck(db->L > 0, PetscObjectComm((PetscObject)db), PETSC_ERR_ARG_WRONG, "Swarm has no points to be removed");
557 PetscCall(DMSwarmDataBucketSetSizes(db, db->L - 1, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT));
558 PetscFunctionReturn(PETSC_SUCCESS);
559 }
560
561 /* Should be redone to user PetscViewer */
DMSwarmDataBucketView_stdout(MPI_Comm comm,DMSwarmDataBucket db)562 static PetscErrorCode DMSwarmDataBucketView_stdout(MPI_Comm comm, DMSwarmDataBucket db)
563 {
564 PetscInt f;
565 double memory_usage_total = 0.0;
566
567 PetscFunctionBegin;
568 PetscCall(PetscPrintf(comm, "DMSwarmDataBucketView: \n"));
569 PetscCall(PetscPrintf(comm, " L = %" PetscInt_FMT " \n", db->L));
570 PetscCall(PetscPrintf(comm, " buffer = %" PetscInt_FMT " \n", db->buffer));
571 PetscCall(PetscPrintf(comm, " allocated = %" PetscInt_FMT " \n", db->allocated));
572 PetscCall(PetscPrintf(comm, " nfields registered = %" PetscInt_FMT " \n", db->nfields));
573
574 for (f = 0; f < db->nfields; ++f) {
575 double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
576 memory_usage_total += memory_usage_f;
577 }
578 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &memory_usage_total, 1, MPI_DOUBLE, MPI_SUM, comm));
579
580 for (f = 0; f < db->nfields; ++f) {
581 double memory_usage_f = (double)(db->field[f]->atomic_size * db->allocated) * 1.0e-6;
582 PetscCall(PetscPrintf(comm, " [%3" PetscInt_FMT "] %15s : Mem. usage = %1.2e (MB) [rank0]\n", f, db->field[f]->name, memory_usage_f));
583 PetscCall(PetscPrintf(comm, " blocksize = %" PetscInt_FMT " \n", db->field[f]->bs));
584 if (db->field[f]->bs != 1) {
585 PetscCall(PetscPrintf(comm, " atomic size = %zu [full block, bs=%" PetscInt_FMT "]\n", db->field[f]->atomic_size, db->field[f]->bs));
586 PetscCall(PetscPrintf(comm, " atomic size/item = %zu \n", (size_t)(db->field[f]->atomic_size / db->field[f]->bs)));
587 } else {
588 PetscCall(PetscPrintf(comm, " atomic size = %zu \n", db->field[f]->atomic_size));
589 }
590 }
591 PetscCall(PetscPrintf(comm, " Total mem. usage = %1.2e (MB) (collective)\n", memory_usage_total));
592 PetscFunctionReturn(PETSC_SUCCESS);
593 }
594
DMSwarmDataBucketView_Seq(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)595 static PetscErrorCode DMSwarmDataBucketView_Seq(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
596 {
597 PetscFunctionBegin;
598 switch (type) {
599 case DATABUCKET_VIEW_STDOUT:
600 PetscCall(DMSwarmDataBucketView_stdout(PETSC_COMM_SELF, db));
601 break;
602 case DATABUCKET_VIEW_ASCII:
603 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for ascii output");
604 case DATABUCKET_VIEW_BINARY:
605 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for binary output");
606 case DATABUCKET_VIEW_HDF5:
607 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for HDF5 output");
608 default:
609 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown viewer method requested");
610 }
611 PetscFunctionReturn(PETSC_SUCCESS);
612 }
613
DMSwarmDataBucketView_MPI(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)614 static PetscErrorCode DMSwarmDataBucketView_MPI(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
615 {
616 PetscFunctionBegin;
617 switch (type) {
618 case DATABUCKET_VIEW_STDOUT:
619 PetscCall(DMSwarmDataBucketView_stdout(comm, db));
620 break;
621 case DATABUCKET_VIEW_ASCII:
622 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for ascii output");
623 case DATABUCKET_VIEW_BINARY:
624 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for binary output");
625 case DATABUCKET_VIEW_HDF5:
626 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for HDF5 output");
627 default:
628 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Unknown viewer method requested");
629 }
630 PetscFunctionReturn(PETSC_SUCCESS);
631 }
632
DMSwarmDataBucketView(MPI_Comm comm,DMSwarmDataBucket db,const char filename[],DMSwarmDataBucketViewType type)633 PetscErrorCode DMSwarmDataBucketView(MPI_Comm comm, DMSwarmDataBucket db, const char filename[], DMSwarmDataBucketViewType type)
634 {
635 PetscMPIInt size;
636
637 PetscFunctionBegin;
638 PetscCallMPI(MPI_Comm_size(comm, &size));
639 if (size == 1) {
640 PetscCall(DMSwarmDataBucketView_Seq(comm, db, filename, type));
641 } else {
642 PetscCall(DMSwarmDataBucketView_MPI(comm, db, filename, type));
643 }
644 PetscFunctionReturn(PETSC_SUCCESS);
645 }
646
DMSwarmDataBucketDuplicateFields(DMSwarmDataBucket dbA,DMSwarmDataBucket * dbB)647 PetscErrorCode DMSwarmDataBucketDuplicateFields(DMSwarmDataBucket dbA, DMSwarmDataBucket *dbB)
648 {
649 DMSwarmDataBucket db2;
650 PetscInt f;
651
652 PetscFunctionBegin;
653 PetscCall(DMSwarmDataBucketCreate(&db2));
654 /* copy contents from dbA into db2 */
655 for (f = 0; f < dbA->nfields; ++f) {
656 DMSwarmDataField field;
657 size_t atomic_size;
658 char *name;
659
660 field = dbA->field[f];
661 atomic_size = field->atomic_size;
662 name = field->name;
663 PetscCall(DMSwarmDataBucketRegisterField(db2, "DMSwarmDataBucketDuplicateFields", name, atomic_size, NULL));
664 }
665 PetscCall(DMSwarmDataBucketFinalize(db2));
666 PetscCall(DMSwarmDataBucketSetInitialSizes(db2, 0, 1000));
667 *dbB = db2;
668 PetscFunctionReturn(PETSC_SUCCESS);
669 }
670
671 /*
672 Insert points from db2 into db1
673 db1 <<== db2
674 */
DMSwarmDataBucketInsertValues(DMSwarmDataBucket db1,DMSwarmDataBucket db2)675 PetscErrorCode DMSwarmDataBucketInsertValues(DMSwarmDataBucket db1, DMSwarmDataBucket db2)
676 {
677 PetscInt n_mp_points1, n_mp_points2;
678 PetscInt n_mp_points1_new, p;
679
680 PetscFunctionBegin;
681 PetscCall(DMSwarmDataBucketGetSizes(db1, &n_mp_points1, NULL, NULL));
682 PetscCall(DMSwarmDataBucketGetSizes(db2, &n_mp_points2, NULL, NULL));
683 n_mp_points1_new = n_mp_points1 + n_mp_points2;
684 PetscCall(DMSwarmDataBucketSetSizes(db1, n_mp_points1_new, DMSWARM_DATA_BUCKET_BUFFER_DEFAULT));
685 for (p = 0; p < n_mp_points2; ++p) {
686 /* db1 <<== db2 */
687 PetscCall(DMSwarmDataBucketCopyPoint(db2, p, db1, n_mp_points1 + p));
688 }
689 PetscFunctionReturn(PETSC_SUCCESS);
690 }
691
692 /* helpers for parallel send/recv */
DMSwarmDataBucketCreatePackedArray(DMSwarmDataBucket db,size_t * bytes,void ** buf)693 PetscErrorCode DMSwarmDataBucketCreatePackedArray(DMSwarmDataBucket db, size_t *bytes, void **buf)
694 {
695 PetscInt f;
696 size_t sizeof_marker_contents;
697 void *buffer;
698
699 PetscFunctionBegin;
700 sizeof_marker_contents = 0;
701 for (f = 0; f < db->nfields; ++f) {
702 DMSwarmDataField df = db->field[f];
703 sizeof_marker_contents += df->atomic_size;
704 }
705 PetscCall(PetscMalloc(sizeof_marker_contents, &buffer));
706 PetscCall(PetscMemzero(buffer, sizeof_marker_contents));
707 if (bytes) *bytes = sizeof_marker_contents;
708 if (buf) *buf = buffer;
709 PetscFunctionReturn(PETSC_SUCCESS);
710 }
711
DMSwarmDataBucketDestroyPackedArray(DMSwarmDataBucket db,void ** buf)712 PetscErrorCode DMSwarmDataBucketDestroyPackedArray(DMSwarmDataBucket db, void **buf)
713 {
714 PetscFunctionBegin;
715 if (buf) {
716 PetscCall(PetscFree(*buf));
717 *buf = NULL;
718 }
719 PetscFunctionReturn(PETSC_SUCCESS);
720 }
721
DMSwarmDataBucketFillPackedArray(DMSwarmDataBucket db,const PetscInt index,void * buf)722 PetscErrorCode DMSwarmDataBucketFillPackedArray(DMSwarmDataBucket db, const PetscInt index, void *buf)
723 {
724 PetscInt f;
725 void *data, *data_p;
726 size_t asize, offset;
727
728 PetscFunctionBegin;
729 offset = 0;
730 for (f = 0; f < db->nfields; ++f) {
731 DMSwarmDataField df = db->field[f];
732
733 asize = df->atomic_size;
734 data = df->data;
735 data_p = (void *)((char *)data + index * asize);
736 PetscCall(PetscMemcpy((void *)((char *)buf + offset), data_p, asize));
737 offset = offset + asize;
738 }
739 PetscFunctionReturn(PETSC_SUCCESS);
740 }
741
DMSwarmDataBucketInsertPackedArray(DMSwarmDataBucket db,const PetscInt idx,void * data)742 PetscErrorCode DMSwarmDataBucketInsertPackedArray(DMSwarmDataBucket db, const PetscInt idx, void *data)
743 {
744 PetscInt f;
745 void *data_p;
746 size_t offset;
747
748 PetscFunctionBegin;
749 offset = 0;
750 for (f = 0; f < db->nfields; ++f) {
751 DMSwarmDataField df = db->field[f];
752
753 data_p = (void *)((char *)data + offset);
754 PetscCall(DMSwarmDataFieldInsertPoint(df, idx, data_p));
755 offset = offset + df->atomic_size;
756 }
757 PetscFunctionReturn(PETSC_SUCCESS);
758 }
759