Actual source code: dmmbvec.cxx

  1: #include <petsc/private/dmmbimpl.h>
  2: #include <petsc/private/vecimpl.h>

  4: #include <petscdmmoab.h>
  5: #include <MBTagConventions.hpp>

  7: #define USE_NATIVE_PETSCVEC

  9: /* declare some private DMMoab specific overrides */
 10: static PetscErrorCode DMCreateVector_Moab_Private(DM dm, moab::Tag tag, const moab::Range *userrange, PetscBool is_global_vec, PetscBool destroy_tag, Vec *vec);
 11: static PetscErrorCode DMVecUserDestroy_Moab(void *user);
 12: static PetscErrorCode DMVecDuplicate_Moab(Vec x, Vec *y);
 13: #ifdef MOAB_HAVE_MPI
 14: static PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, moab::ParallelComm *pcomm, char **tag_name);
 15: #else
 16: static PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, char **tag_name);
 17: #endif

 19: /*@C
 20:   DMMoabCreateVector - Create a Vec from either an existing tag, or a specified tag size, and a range of entities

 22:   Collective

 24:   Input Parameters:
 25: + dm              - The DMMoab object being set
 26: . tag             - If non-zero, block size will be taken from the tag size
 27: . range           - If non-empty, Vec corresponds to these entities, otherwise to the entities set on the DMMoab
 28: . is_global_vec   - If true, this is a local representation of the Vec (including ghosts in parallel), otherwise a truly parallel one
 29: - destroy_tag     - If true, MOAB tag is destroyed with Vec, otherwise it is left on MOAB

 31:   Output Parameter:
 32: . vec             - The created vector

 34:   Level: beginner

 36: .seealso: `VecCreate()`
 37: @*/
 38: PetscErrorCode DMMoabCreateVector(DM dm, moab::Tag tag, const moab::Range *range, PetscBool is_global_vec, PetscBool destroy_tag, Vec *vec)
 39: {

 42:   DMCreateVector_Moab_Private(dm, tag, range, is_global_vec, destroy_tag, vec);
 43:   return 0;
 44: }

 46: /*@C
 47:   DMMoabGetVecTag - Get the MOAB tag associated with this Vec

 49:   Input Parameter:
 50: . vec - Vec being queried

 52:   Output Parameter:
 53: . tag - Tag associated with this Vec. NULL if vec is a native PETSc Vec.

 55:   Level: beginner

 57: .seealso: `DMMoabCreateVector()`, `DMMoabGetVecRange()`
 58: @*/
 59: PetscErrorCode DMMoabGetVecTag(Vec vec, moab::Tag *tag)
 60: {
 61:   PetscContainer moabdata;
 62:   Vec_MOAB      *vmoab;


 66:   /* Get the MOAB private data */
 67:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
 68:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

 70:   *tag = vmoab->tag;
 71:   return 0;
 72: }

 74: /*@C
 75:   DMMoabGetVecRange - Get the MOAB entities associated with this Vec

 77:   Input Parameter:
 78: . vec   - Vec being queried

 80:   Output Parameter:
 81: . range - Entities associated with this Vec. NULL if vec is a native PETSc Vec.

 83:   Level: beginner

 85: .seealso: `DMMoabCreateVector()`, `DMMoabGetVecTag()`
 86: @*/
 87: PetscErrorCode DMMoabGetVecRange(Vec vec, moab::Range *range)
 88: {
 89:   PetscContainer moabdata;
 90:   Vec_MOAB      *vmoab;


 94:   /* Get the MOAB private data handle */
 95:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
 96:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

 98:   *range = *vmoab->tag_range;
 99:   return 0;
100: }

102: /*@C
103:   DMMoabVecGetArray - Returns the writable direct access array to the local representation of MOAB tag data for the underlying vector using locally owned+ghosted range of entities

105:   Collective

107:   Input Parameters:
108: + dm              - The DMMoab object being set
109: - vec             - The Vector whose underlying data is requested

111:   Output Parameter:
112: . array           - The local data array

114:   Level: intermediate

116: .seealso: `DMMoabVecRestoreArray()`, `DMMoabVecGetArrayRead()`, `DMMoabVecRestoreArrayRead()`
117: @*/
118: PetscErrorCode DMMoabVecGetArray(DM dm, Vec vec, void *array)
119: {
120:   DM_Moab        *dmmoab;
121:   moab::ErrorCode merr;
122:   PetscInt        count, i, f;
123:   moab::Tag       vtag;
124:   PetscScalar   **varray;
125:   PetscScalar    *marray;
126:   PetscContainer  moabdata;
127:   Vec_MOAB       *vmoab, *xmoab;

132:   dmmoab = (DM_Moab *)dm->data;

134:   /* Get the Vec_MOAB struct for the original vector */
135:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
136:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

138:   /* Get the real scalar array handle */
139:   varray = reinterpret_cast<PetscScalar **>(array);

141:   if (vmoab->is_native_vec) {
142:     /* get the local representation of the arrays from Vectors */
143:     DMGetLocalVector(dm, &vmoab->local);
144:     DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
145:     DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);

147:     /* Get the Vec_MOAB struct for the original vector */
148:     PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject *)&moabdata);
149:     PetscContainerGetPointer(moabdata, (void **)&xmoab);

151:     /* get the local representation of the arrays from Vectors */
152:     VecGhostGetLocalForm(vmoab->local, &xmoab->local);
153:     VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
154:     VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);

156:     VecGetArray(xmoab->local, varray);
157:   } else {
158:     /* Get the MOAB private data */
159:     DMMoabGetVecTag(vec, &vtag);

161: #ifdef MOAB_HAVE_MPI
162:     /* exchange the data into ghost cells first */
163:     merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal);
164:     MBERRNM(merr);
165: #endif

167:     PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);

169:     /* Get the array data for local entities */
170:     merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void *&>(marray), false);
171:     MBERRNM(merr);

174:     i = 0;
175:     for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
176:       for (f = 0; f < dmmoab->numFields; f++, i++) (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart] * dmmoab->numFields + f] = marray[i];
177:       //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
178:     }
179:   }
180:   return 0;
181: }

183: /*@C
184:   DMMoabVecRestoreArray - Restores the writable direct access array obtained via DMMoabVecGetArray

186:   Collective

188:   Input Parameters:
189: + dm              - The DMMoab object being set
190: . vec             - The Vector whose underlying data is requested
191: - array           - The local data array

193:   Level: intermediate

195: .seealso: `DMMoabVecGetArray()`, `DMMoabVecGetArrayRead()`, `DMMoabVecRestoreArrayRead()`
196: @*/
197: PetscErrorCode DMMoabVecRestoreArray(DM dm, Vec vec, void *array)
198: {
199:   DM_Moab        *dmmoab;
200:   moab::ErrorCode merr;
201:   moab::Tag       vtag;
202:   PetscInt        count, i, f;
203:   PetscScalar   **varray;
204:   PetscScalar    *marray;
205:   PetscContainer  moabdata;
206:   Vec_MOAB       *vmoab, *xmoab;

211:   dmmoab = (DM_Moab *)dm->data;

213:   /* Get the Vec_MOAB struct for the original vector */
214:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
215:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

217:   /* Get the real scalar array handle */
218:   varray = reinterpret_cast<PetscScalar **>(array);

220:   if (vmoab->is_native_vec) {
221:     /* Get the Vec_MOAB struct for the original vector */
222:     PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject *)&moabdata);
223:     PetscContainerGetPointer(moabdata, (void **)&xmoab);

225:     /* get the local representation of the arrays from Vectors */
226:     VecRestoreArray(xmoab->local, varray);
227:     VecGhostUpdateBegin(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
228:     VecGhostUpdateEnd(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
229:     VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);

231:     /* restore local pieces */
232:     DMLocalToGlobalBegin(dm, vmoab->local, INSERT_VALUES, vec);
233:     DMLocalToGlobalEnd(dm, vmoab->local, INSERT_VALUES, vec);
234:     DMRestoreLocalVector(dm, &vmoab->local);
235:   } else {
236:     /* Get the MOAB private data */
237:     DMMoabGetVecTag(vec, &vtag);

239:     /* Get the array data for local entities */
240:     merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void *&>(marray), false);
241:     MBERRNM(merr);

244:     i = 0;
245:     for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
246:       for (f = 0; f < dmmoab->numFields; f++, i++) marray[i] = (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart] * dmmoab->numFields + f];
247:       //marray[i] = (*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]];
248:     }

250: #ifdef MOAB_HAVE_MPI
251:     /* reduce the tags correctly -> should probably let the user choose how to reduce in the future
252:       For all FEM residual based assembly calculations, MPI_SUM should serve well */
253:     merr = dmmoab->pcomm->reduce_tags(vtag, MPI_SUM, *dmmoab->vlocal);
254:     MBERRV(dmmoab->mbiface, merr);
255: #endif
256:     PetscFree(*varray);
257:   }
258:   return 0;
259: }

261: /*@C
262:   DMMoabVecGetArrayRead - Returns the read-only direct access array to the local representation of MOAB tag data for the underlying vector using locally owned+ghosted range of entities

264:   Collective

266:   Input Parameters:
267: + dm              - The DMMoab object being set
268: - vec             - The Vector whose underlying data is requested

270:   Output Parameter:
271: . array           - The local data array

273:   Level: intermediate

275: .seealso: `DMMoabVecRestoreArrayRead()`, `DMMoabVecGetArray()`, `DMMoabVecRestoreArray()`
276: @*/
277: PetscErrorCode DMMoabVecGetArrayRead(DM dm, Vec vec, void *array)
278: {
279:   DM_Moab        *dmmoab;
280:   moab::ErrorCode merr;
281:   PetscInt        count, i, f;
282:   moab::Tag       vtag;
283:   PetscScalar   **varray;
284:   PetscScalar    *marray;
285:   PetscContainer  moabdata;
286:   Vec_MOAB       *vmoab, *xmoab;

291:   dmmoab = (DM_Moab *)dm->data;

293:   /* Get the Vec_MOAB struct for the original vector */
294:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
295:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

297:   /* Get the real scalar array handle */
298:   varray = reinterpret_cast<PetscScalar **>(array);

300:   if (vmoab->is_native_vec) {
301:     /* get the local representation of the arrays from Vectors */
302:     DMGetLocalVector(dm, &vmoab->local);
303:     DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
304:     DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);

306:     /* Get the Vec_MOAB struct for the original vector */
307:     PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject *)&moabdata);
308:     PetscContainerGetPointer(moabdata, (void **)&xmoab);

310:     /* get the local representation of the arrays from Vectors */
311:     VecGhostGetLocalForm(vmoab->local, &xmoab->local);
312:     VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
313:     VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
314:     VecGetArray(xmoab->local, varray);
315:   } else {
316:     /* Get the MOAB private data */
317:     DMMoabGetVecTag(vec, &vtag);

319: #ifdef MOAB_HAVE_MPI
320:     /* exchange the data into ghost cells first */
321:     merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal);
322:     MBERRNM(merr);
323: #endif
324:     PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);

326:     /* Get the array data for local entities */
327:     merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void *&>(marray), false);
328:     MBERRNM(merr);

331:     i = 0;
332:     for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
333:       for (f = 0; f < dmmoab->numFields; f++, i++) (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart] * dmmoab->numFields + f] = marray[i];
334:       //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
335:     }
336:   }
337:   return 0;
338: }

340: /*@C
341:   DMMoabVecRestoreArrayRead - Restores the read-only direct access array obtained via DMMoabVecGetArray

343:   Collective

345:   Input Parameters:
346: + dm              - The DMMoab object being set
347: . vec             - The Vector whose underlying data is requested
348: - array           - The local data array

350:   Level: intermediate

352: .seealso: `DMMoabVecGetArrayRead()`, `DMMoabVecGetArray()`, `DMMoabVecRestoreArray()`
353: @*/
354: PetscErrorCode DMMoabVecRestoreArrayRead(DM dm, Vec vec, void *array)
355: {
356:   PetscScalar  **varray;
357:   PetscContainer moabdata;
358:   Vec_MOAB      *vmoab, *xmoab;


364:   /* Get the Vec_MOAB struct for the original vector */
365:   PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject *)&moabdata);
366:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

368:   /* Get the real scalar array handle */
369:   varray = reinterpret_cast<PetscScalar **>(array);

371:   if (vmoab->is_native_vec) {
372:     /* Get the Vec_MOAB struct for the original vector */
373:     PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject *)&moabdata);
374:     PetscContainerGetPointer(moabdata, (void **)&xmoab);

376:     /* restore the local representation of the arrays from Vectors */
377:     VecRestoreArray(xmoab->local, varray);
378:     VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);

380:     /* restore local pieces */
381:     DMRestoreLocalVector(dm, &vmoab->local);
382:   } else {
383:     /* Nothing to do but just free the memory allocated before */
384:     PetscFree(*varray);
385:   }
386:   return 0;
387: }

389: PetscErrorCode DMCreateVector_Moab_Private(DM dm, moab::Tag tag, const moab::Range *userrange, PetscBool is_global_vec, PetscBool destroy_tag, Vec *vec)
390: {
391:   moab::ErrorCode    merr;
392:   PetscBool          is_newtag;
393:   const moab::Range *range;
394:   PetscInt           count, lnative_vec, gnative_vec;
395:   std::string        ttname;
396:   PetscScalar       *data_ptr, *defaultvals;

398:   Vec_MOAB *vmoab;
399:   DM_Moab  *dmmoab = (DM_Moab *)dm->data;
400: #ifdef MOAB_HAVE_MPI
401:   moab::ParallelComm *pcomm = dmmoab->pcomm;
402: #endif
403:   moab::Interface *mbiface = dmmoab->mbiface;

406:   if (!userrange) range = dmmoab->vowned;
407:   else range = userrange;

410: #ifndef USE_NATIVE_PETSCVEC
411:   /* If the tag data is in a single sequence, use PETSc native vector since tag_iterate isn't useful anymore */
412:   lnative_vec = (range->psize() - 1);
413: #else
414:   lnative_vec = 1; /* NOTE: Testing PETSc vector: will force to create native vector all the time */
415:                    //  lnative_vec=0; /* NOTE: Testing MOAB vector: will force to create MOAB tag_iterate based vector all the time */
416: #endif

418: #ifdef MOAB_HAVE_MPI
419:   MPIU_Allreduce(&lnative_vec, &gnative_vec, 1, MPI_INT, MPI_MAX, (((PetscObject)dm)->comm));
420: #else
421:   gnative_vec = lnative_vec;
422: #endif

424:   /* Create the MOAB internal data object */
425:   PetscNew(&vmoab);
426:   vmoab->is_native_vec = (gnative_vec > 0 ? PETSC_TRUE : PETSC_FALSE);

428:   if (!vmoab->is_native_vec) {
429:     merr = moab::MB_SUCCESS;
430:     if (tag != 0) merr = mbiface->tag_get_name(tag, ttname);
431:     if (!ttname.length() || merr != moab::MB_SUCCESS) {
432:       /* get the new name for the anonymous MOABVec -> the tag_name will be destroyed along with Tag */
433:       char *tag_name = NULL;
434: #ifdef MOAB_HAVE_MPI
435:       DMVecCreateTagName_Moab_Private(mbiface, pcomm, &tag_name);
436: #else
437:       DMVecCreateTagName_Moab_Private(mbiface, &tag_name);
438: #endif
439:       is_newtag = PETSC_TRUE;

441:       /* Create the default value for the tag (all zeros) */
442:       PetscCalloc1(dmmoab->numFields, &defaultvals);

444:       /* Create the tag */
445:       merr = mbiface->tag_get_handle(tag_name, dmmoab->numFields, moab::MB_TYPE_DOUBLE, tag, moab::MB_TAG_DENSE | moab::MB_TAG_CREAT, defaultvals);
446:       MBERRNM(merr);
447:       PetscFree(tag_name);
448:       PetscFree(defaultvals);
449:     } else {
450:       /* Make sure the tag data is of type "double" */
451:       moab::DataType tag_type;
452:       merr = mbiface->tag_get_data_type(tag, tag_type);
453:       MBERRNM(merr);
455:       is_newtag = destroy_tag;
456:     }

458:     vmoab->tag     = tag;
459:     vmoab->new_tag = is_newtag;
460:   }
461:   vmoab->mbiface = mbiface;
462: #ifdef MOAB_HAVE_MPI
463:   vmoab->pcomm = pcomm;
464: #endif
465:   vmoab->is_global_vec = is_global_vec;
466:   vmoab->tag_size      = dmmoab->bs;

468:   if (vmoab->is_native_vec) {
469:     /* Create the PETSc Vector directly and attach our functions accordingly */
470:     if (!is_global_vec) {
471:       /* This is an MPI Vector with ghosted padding */
472:       VecCreateGhostBlock((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc, dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], vec);
473:     } else {
474:       /* This is an MPI/SEQ Vector */
475:       VecCreate((((PetscObject)dm)->comm), vec);
476:       VecSetSizes(*vec, dmmoab->numFields * dmmoab->nloc, PETSC_DECIDE);
477:       VecSetBlockSize(*vec, dmmoab->bs);
478:       VecSetType(*vec, VECMPI);
479:     }
480:   } else {
481:     /* Call tag_iterate. This will cause MOAB to allocate memory for the
482:        tag data if it hasn't already happened */
483:     merr = mbiface->tag_iterate(tag, range->begin(), range->end(), count, (void *&)data_ptr);
484:     MBERRNM(merr);

486:     /* set the reference for vector range */
487:     vmoab->tag_range = new moab::Range(*range);
488:     merr             = mbiface->tag_get_length(tag, dmmoab->numFields);
489:     MBERRNM(merr);

491:     /* Create the PETSc Vector
492:       Query MOAB mesh to check if there are any ghosted entities
493:         -> if we do, create a ghosted vector to map correctly to the same layout
494:         -> else, create a non-ghosted parallel vector */
495:     if (!is_global_vec) {
496:       /* This is an MPI Vector with ghosted padding */
497:       VecCreateGhostBlockWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc, dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], data_ptr, vec);
498:     } else {
499:       /* This is an MPI Vector without ghosted padding */
500:       VecCreateMPIWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * range->size(), PETSC_DECIDE, data_ptr, vec);
501:     }
502:   }
503:   VecSetFromOptions(*vec);

505:   /* create a container and store the internal MOAB data for faster access based on Entities etc */
506:   PetscContainer moabdata;
507:   PetscContainerCreate(PETSC_COMM_WORLD, &moabdata);
508:   PetscContainerSetPointer(moabdata, vmoab);
509:   PetscContainerSetUserDestroy(moabdata, DMVecUserDestroy_Moab);
510:   PetscObjectCompose((PetscObject)*vec, "MOABData", (PetscObject)moabdata);
511:   (*vec)->ops->duplicate = DMVecDuplicate_Moab;
512:   PetscContainerDestroy(&moabdata);

514:   /* Vector created, manually set local to global mapping */
515:   if (dmmoab->ltog_map) VecSetLocalToGlobalMapping(*vec, dmmoab->ltog_map);

517:   /* set the DM reference to the vector */
518:   VecSetDM(*vec, dm);
519:   return 0;
520: }

522: /*  DMVecCreateTagName_Moab_Private
523:  *
524:  *  Creates a unique tag name that will be shared across processes. If
525:  *  pcomm is NULL, then this is a serial vector. A unique tag name
526:  *  will be returned in tag_name in either case.
527:  *
528:  *  The tag names have the format _PETSC_VEC_N where N is some integer.
529:  *
530:  *  NOTE: The tag_name is allocated in this routine; The user needs to free
531:  *        the character array.
532:  */
533: #ifdef MOAB_HAVE_MPI
534: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, moab::ParallelComm *pcomm, char **tag_name)
535: #else
536: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, char **tag_name)
537: #endif
538: {
539:   moab::ErrorCode mberr;
540:   PetscInt        n, global_n;
541:   moab::Tag       indexTag;

543:   const char *PVEC_PREFIX = "__PETSC_VEC_";
544:   PetscMalloc1(PETSC_MAX_PATH_LEN, tag_name);

546:   moab::EntityHandle rootset = mbiface->get_root_set();

548:   /* Check to see if there are any PETSc vectors defined */
549:   /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */
550:   mberr = mbiface->tag_get_handle("__PETSC_VECS__", 1, moab::MB_TYPE_INTEGER, indexTag, moab::MB_TAG_SPARSE | moab::MB_TAG_CREAT, 0);
551:   MBERRNM(mberr);
552:   mberr = mbiface->tag_get_data(indexTag, &rootset, 1, &n);
553:   if (mberr == moab::MB_TAG_NOT_FOUND) n = 0; /* this is the first temporary vector */
554:   else MBERRNM(mberr);

556:   /* increment the new value of n */
557:   ++n;

559: #ifdef MOAB_HAVE_MPI
560:   /* Make sure that n is consistent across all processes */
561:   MPIU_Allreduce(&n, &global_n, 1, MPI_INT, MPI_MAX, pcomm->comm());
562: #else
563:   global_n = n;
564: #endif

566:   /* Set the new name accordingly and return */
567:   PetscSNPrintf(*tag_name, PETSC_MAX_PATH_LEN - 1, "%s_%" PetscInt_FMT, PVEC_PREFIX, global_n);
568:   mberr = mbiface->tag_set_data(indexTag, &rootset, 1, (const void *)&global_n);
569:   MBERRNM(mberr);
570:   return 0;
571: }

573: PETSC_EXTERN PetscErrorCode DMCreateGlobalVector_Moab(DM dm, Vec *gvec)
574: {
575:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

579:   DMCreateVector_Moab_Private(dm, NULL, dmmoab->vowned, PETSC_TRUE, PETSC_TRUE, gvec);
580:   return 0;
581: }

583: PETSC_EXTERN PetscErrorCode DMCreateLocalVector_Moab(DM dm, Vec *lvec)
584: {
585:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

589:   DMCreateVector_Moab_Private(dm, NULL, dmmoab->vlocal, PETSC_FALSE, PETSC_TRUE, lvec);
590:   return 0;
591: }

593: PetscErrorCode DMVecDuplicate_Moab(Vec x, Vec *y)
594: {
595:   DM             dm;
596:   PetscContainer moabdata;
597:   Vec_MOAB      *vmoab;


602:   /* Get the Vec_MOAB struct for the original vector */
603:   PetscObjectQuery((PetscObject)x, "MOABData", (PetscObject *)&moabdata);
604:   PetscContainerGetPointer(moabdata, (void **)&vmoab);

606:   VecGetDM(x, &dm);

609:   DMCreateVector_Moab_Private(dm, NULL, vmoab->tag_range, vmoab->is_global_vec, PETSC_TRUE, y);
610:   VecSetDM(*y, dm);
611:   return 0;
612: }

614: PetscErrorCode DMVecUserDestroy_Moab(void *user)
615: {
616:   Vec_MOAB       *vmoab = (Vec_MOAB *)user;
617:   moab::ErrorCode merr;

619:   if (vmoab->new_tag && vmoab->tag) {
620:     /* Tag was created via a call to VecDuplicate, delete the underlying tag in MOAB */
621:     merr = vmoab->mbiface->tag_delete(vmoab->tag);
622:     MBERRNM(merr);
623:   }
624:   delete vmoab->tag_range;
625:   vmoab->tag     = NULL;
626:   vmoab->mbiface = NULL;
627: #ifdef MOAB_HAVE_MPI
628:   vmoab->pcomm = NULL;
629: #endif
630:   PetscFree(vmoab);
631:   return 0;
632: }

634: PETSC_EXTERN PetscErrorCode DMGlobalToLocalBegin_Moab(DM dm, Vec g, InsertMode mode, Vec l)
635: {
636:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

638:   VecScatterBegin(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
639:   return 0;
640: }

642: PETSC_EXTERN PetscErrorCode DMGlobalToLocalEnd_Moab(DM dm, Vec g, InsertMode mode, Vec l)
643: {
644:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

646:   VecScatterEnd(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
647:   return 0;
648: }

650: PETSC_EXTERN PetscErrorCode DMLocalToGlobalBegin_Moab(DM dm, Vec l, InsertMode mode, Vec g)
651: {
652:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

654:   VecScatterBegin(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
655:   return 0;
656: }

658: PETSC_EXTERN PetscErrorCode DMLocalToGlobalEnd_Moab(DM dm, Vec l, InsertMode mode, Vec g)
659: {
660:   DM_Moab *dmmoab = (DM_Moab *)dm->data;

662:   VecScatterEnd(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
663:   return 0;
664: }