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: {
40: PetscErrorCode ierr;
43: if (!tag && (!range || range->empty())) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Both tag and range cannot be null.");
45: DMCreateVector_Moab_Private(dm, tag, range, is_global_vec, destroy_tag, vec);
46: return(0);
47: }
49: /*@C
50: DMMoabGetVecTag - Get the MOAB tag associated with this Vec
52: Input Parameter:
53: . vec - Vec being queried
55: Output Parameter:
56: . tag - Tag associated with this Vec. NULL if vec is a native PETSc Vec.
58: Level: beginner
60: .seealso: DMMoabCreateVector(), DMMoabGetVecRange()
61: @*/
62: PetscErrorCode DMMoabGetVecTag(Vec vec, moab::Tag *tag)
63: {
64: PetscContainer moabdata;
65: Vec_MOAB *vmoab;
66: PetscErrorCode ierr;
71: /* Get the MOAB private data */
72: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
73: PetscContainerGetPointer(moabdata, (void**) &vmoab);
75: *tag = vmoab->tag;
76: return(0);
77: }
79: /*@C
80: DMMoabGetVecRange - Get the MOAB entities associated with this Vec
82: Input Parameter:
83: . vec - Vec being queried
85: Output Parameter:
86: . range - Entities associated with this Vec. NULL if vec is a native PETSc Vec.
88: Level: beginner
90: .seealso: DMMoabCreateVector(), DMMoabGetVecTag()
91: @*/
92: PetscErrorCode DMMoabGetVecRange(Vec vec, moab::Range *range)
93: {
94: PetscContainer moabdata;
95: Vec_MOAB *vmoab;
96: PetscErrorCode ierr;
101: /* Get the MOAB private data handle */
102: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
103: PetscContainerGetPointer(moabdata, (void**) &vmoab);
105: *range = *vmoab->tag_range;
106: return(0);
107: }
109: /*@C
110: 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
112: Collective
114: Input Parameters:
115: + dm - The DMMoab object being set
116: - vec - The Vector whose underlying data is requested
118: Output Parameter:
119: . array - The local data array
121: Level: intermediate
123: .seealso: DMMoabVecRestoreArray(), DMMoabVecGetArrayRead(), DMMoabVecRestoreArrayRead()
124: @*/
125: PetscErrorCode DMMoabVecGetArray(DM dm, Vec vec, void* array)
126: {
127: DM_Moab *dmmoab;
128: moab::ErrorCode merr;
129: PetscErrorCode ierr;
130: PetscInt count, i, f;
131: moab::Tag vtag;
132: PetscScalar **varray;
133: PetscScalar *marray;
134: PetscContainer moabdata;
135: Vec_MOAB *vmoab, *xmoab;
141: dmmoab = (DM_Moab*)dm->data;
143: /* Get the Vec_MOAB struct for the original vector */
144: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
145: PetscContainerGetPointer(moabdata, (void**)&vmoab);
147: /* Get the real scalar array handle */
148: varray = reinterpret_cast<PetscScalar**>(array);
150: if (vmoab->is_native_vec) {
152: /* get the local representation of the arrays from Vectors */
153: DMGetLocalVector(dm, &vmoab->local);
154: DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
155: DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);
157: /* Get the Vec_MOAB struct for the original vector */
158: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
159: PetscContainerGetPointer(moabdata, (void**)&xmoab);
161: /* get the local representation of the arrays from Vectors */
162: VecGhostGetLocalForm(vmoab->local, &xmoab->local);
163: VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
164: VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
166: VecGetArray(xmoab->local, varray);
167: }
168: else {
170: /* Get the MOAB private data */
171: DMMoabGetVecTag(vec, &vtag);
173: #ifdef MOAB_HAVE_MPI
174: /* exchange the data into ghost cells first */
175: merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal); MBERRNM(merr);
176: #endif
178: PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);
180: /* Get the array data for local entities */
181: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
182: if (count != dmmoab->nloc + dmmoab->nghost) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mismatch between local vertices and tag partition for Vec. %D != %D.", count, dmmoab->nloc + dmmoab->nghost);
184: i = 0;
185: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
186: for (f = 0; f < dmmoab->numFields; f++, i++)
187: (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart]*dmmoab->numFields + f] = marray[i];
188: //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
189: }
190: }
191: return(0);
192: }
194: /*@C
195: DMMoabVecRestoreArray - Restores the writable direct access array obtained via DMMoabVecGetArray
197: Collective
199: Input Parameters:
200: + dm - The DMMoab object being set
201: . vec - The Vector whose underlying data is requested
202: - array - The local data array
204: Level: intermediate
206: .seealso: DMMoabVecGetArray(), DMMoabVecGetArrayRead(), DMMoabVecRestoreArrayRead()
207: @*/
208: PetscErrorCode DMMoabVecRestoreArray(DM dm, Vec vec, void* array)
209: {
210: DM_Moab *dmmoab;
211: moab::ErrorCode merr;
212: PetscErrorCode ierr;
213: moab::Tag vtag;
214: PetscInt count, i, f;
215: PetscScalar **varray;
216: PetscScalar *marray;
217: PetscContainer moabdata;
218: Vec_MOAB *vmoab, *xmoab;
224: dmmoab = (DM_Moab*)dm->data;
226: /* Get the Vec_MOAB struct for the original vector */
227: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
228: PetscContainerGetPointer(moabdata, (void**)&vmoab);
230: /* Get the real scalar array handle */
231: varray = reinterpret_cast<PetscScalar**>(array);
233: if (vmoab->is_native_vec) {
235: /* Get the Vec_MOAB struct for the original vector */
236: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
237: PetscContainerGetPointer(moabdata, (void**)&xmoab);
239: /* get the local representation of the arrays from Vectors */
240: VecRestoreArray(xmoab->local, varray);
241: VecGhostUpdateBegin(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
242: VecGhostUpdateEnd(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
243: VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);
245: /* restore local pieces */
246: DMLocalToGlobalBegin(dm, vmoab->local, INSERT_VALUES, vec);
247: DMLocalToGlobalEnd(dm, vmoab->local, INSERT_VALUES, vec);
248: DMRestoreLocalVector(dm, &vmoab->local);
249: }
250: else {
252: /* Get the MOAB private data */
253: DMMoabGetVecTag(vec, &vtag);
255: /* Get the array data for local entities */
256: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
257: if (count != dmmoab->nloc + dmmoab->nghost) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mismatch between local vertices and tag partition for Vec. %D != %D.", count, dmmoab->nloc + dmmoab->nghost);
259: i = 0;
260: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
261: for (f = 0; f < dmmoab->numFields; f++, i++)
262: marray[i] = (*varray)[dmmoab->lidmap[(PetscInt) * iter - dmmoab->seqstart] * dmmoab->numFields + f];
263: //marray[i] = (*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]];
264: }
266: #ifdef MOAB_HAVE_MPI
267: /* reduce the tags correctly -> should probably let the user choose how to reduce in the future
268: For all FEM residual based assembly calculations, MPI_SUM should serve well */
269: merr = dmmoab->pcomm->reduce_tags(vtag, MPI_SUM, *dmmoab->vlocal); MBERRV(dmmoab->mbiface, merr);
270: #endif
271: PetscFree(*varray);
272: }
273: return(0);
274: }
276: /*@C
277: 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
279: Collective
281: Input Parameters:
282: + dm - The DMMoab object being set
283: - vec - The Vector whose underlying data is requested
285: Output Parameter:
286: . array - The local data array
288: Level: intermediate
290: .seealso: DMMoabVecRestoreArrayRead(), DMMoabVecGetArray(), DMMoabVecRestoreArray()
291: @*/
292: PetscErrorCode DMMoabVecGetArrayRead(DM dm, Vec vec, void* array)
293: {
294: DM_Moab *dmmoab;
295: moab::ErrorCode merr;
296: PetscErrorCode ierr;
297: PetscInt count, i, f;
298: moab::Tag vtag;
299: PetscScalar **varray;
300: PetscScalar *marray;
301: PetscContainer moabdata;
302: Vec_MOAB *vmoab, *xmoab;
308: dmmoab = (DM_Moab*)dm->data;
310: /* Get the Vec_MOAB struct for the original vector */
311: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
312: PetscContainerGetPointer(moabdata, (void**)&vmoab);
314: /* Get the real scalar array handle */
315: varray = reinterpret_cast<PetscScalar**>(array);
317: if (vmoab->is_native_vec) {
318: /* get the local representation of the arrays from Vectors */
319: DMGetLocalVector(dm, &vmoab->local);
320: DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
321: DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);
323: /* Get the Vec_MOAB struct for the original vector */
324: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
325: PetscContainerGetPointer(moabdata, (void**)&xmoab);
327: /* get the local representation of the arrays from Vectors */
328: VecGhostGetLocalForm(vmoab->local, &xmoab->local);
329: VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
330: VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
331: VecGetArray(xmoab->local, varray);
332: }
333: else {
334: /* Get the MOAB private data */
335: DMMoabGetVecTag(vec, &vtag);
337: #ifdef MOAB_HAVE_MPI
338: /* exchange the data into ghost cells first */
339: merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal); MBERRNM(merr);
340: #endif
341: PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);
343: /* Get the array data for local entities */
344: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
345: if (count != dmmoab->nloc + dmmoab->nghost) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mismatch between local vertices and tag partition for Vec. %D != %D.", count, dmmoab->nloc + dmmoab->nghost);
347: i = 0;
348: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
349: for (f = 0; f < dmmoab->numFields; f++, i++)
350: (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart]*dmmoab->numFields + f] = marray[i];
351: //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
352: }
353: }
354: return(0);
355: }
357: /*@C
358: DMMoabVecRestoreArrayRead - Restores the read-only direct access array obtained via DMMoabVecGetArray
360: Collective
362: Input Parameters:
363: + dm - The DMMoab object being set
364: . vec - The Vector whose underlying data is requested
365: - array - The local data array
367: Level: intermediate
369: .seealso: DMMoabVecGetArrayRead(), DMMoabVecGetArray(), DMMoabVecRestoreArray()
370: @*/
371: PetscErrorCode DMMoabVecRestoreArrayRead(DM dm, Vec vec, void* array)
372: {
373: PetscErrorCode ierr;
374: PetscScalar **varray;
375: PetscContainer moabdata;
376: Vec_MOAB *vmoab, *xmoab;
383: /* Get the Vec_MOAB struct for the original vector */
384: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
385: PetscContainerGetPointer(moabdata, (void**)&vmoab);
387: /* Get the real scalar array handle */
388: varray = reinterpret_cast<PetscScalar**>(array);
390: if (vmoab->is_native_vec) {
391: /* Get the Vec_MOAB struct for the original vector */
392: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
393: PetscContainerGetPointer(moabdata, (void**)&xmoab);
395: /* restore the local representation of the arrays from Vectors */
396: VecRestoreArray(xmoab->local, varray);
397: VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);
399: /* restore local pieces */
400: DMRestoreLocalVector(dm, &vmoab->local);
401: }
402: else {
403: /* Nothing to do but just free the memory allocated before */
404: PetscFree(*varray);
406: }
407: return(0);
408: }
410: PetscErrorCode DMCreateVector_Moab_Private(DM dm, moab::Tag tag, const moab::Range* userrange, PetscBool is_global_vec, PetscBool destroy_tag, Vec *vec)
411: {
412: PetscErrorCode ierr;
413: moab::ErrorCode merr;
414: PetscBool is_newtag;
415: const moab::Range *range;
416: PetscInt count, lnative_vec, gnative_vec;
417: std::string ttname;
418: PetscScalar *data_ptr, *defaultvals;
420: Vec_MOAB *vmoab;
421: DM_Moab *dmmoab = (DM_Moab*)dm->data;
422: #ifdef MOAB_HAVE_MPI
423: moab::ParallelComm *pcomm = dmmoab->pcomm;
424: #endif
425: moab::Interface *mbiface = dmmoab->mbiface;
428: if (sizeof(PetscReal) != sizeof(PetscScalar)) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "MOAB tags only support Real types (Complex-type unsupported)");
429: if (!userrange) range = dmmoab->vowned;
430: else range = userrange;
431: if (!range) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Input range cannot be empty or call DMSetUp first.");
433: #ifndef USE_NATIVE_PETSCVEC
434: /* If the tag data is in a single sequence, use PETSc native vector since tag_iterate isn't useful anymore */
435: lnative_vec = (range->psize() - 1);
436: #else
437: lnative_vec = 1; /* NOTE: Testing PETSc vector: will force to create native vector all the time */
438: // lnative_vec=0; /* NOTE: Testing MOAB vector: will force to create MOAB tag_iterate based vector all the time */
439: #endif
441: #ifdef MOAB_HAVE_MPI
442: MPIU_Allreduce(&lnative_vec, &gnative_vec, 1, MPI_INT, MPI_MAX, (((PetscObject)dm)->comm));
443: #else
444: gnative_vec = lnative_vec;
445: #endif
447: /* Create the MOAB internal data object */
448: PetscNew(&vmoab);
449: vmoab->is_native_vec = (gnative_vec > 0 ? PETSC_TRUE : PETSC_FALSE);
451: if (!vmoab->is_native_vec) {
452: merr = moab::MB_SUCCESS;
453: if (tag != 0) merr = mbiface->tag_get_name(tag, ttname);
454: if (!ttname.length() || merr != moab::MB_SUCCESS) {
455: /* get the new name for the anonymous MOABVec -> the tag_name will be destroyed along with Tag */
456: char *tag_name = NULL;
457: #ifdef MOAB_HAVE_MPI
458: DMVecCreateTagName_Moab_Private(mbiface,pcomm,&tag_name);
459: #else
460: DMVecCreateTagName_Moab_Private(mbiface,&tag_name);
461: #endif
462: is_newtag = PETSC_TRUE;
464: /* Create the default value for the tag (all zeros) */
465: PetscCalloc1(dmmoab->numFields, &defaultvals);
467: /* Create the tag */
468: merr = mbiface->tag_get_handle(tag_name, dmmoab->numFields, moab::MB_TYPE_DOUBLE, tag,
469: moab::MB_TAG_DENSE | moab::MB_TAG_CREAT, defaultvals); MBERRNM(merr);
470: PetscFree(tag_name);
471: PetscFree(defaultvals);
472: }
473: else {
474: /* Make sure the tag data is of type "double" */
475: moab::DataType tag_type;
476: merr = mbiface->tag_get_data_type(tag, tag_type); MBERRNM(merr);
477: if (tag_type != moab::MB_TYPE_DOUBLE) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Tag data type must be MB_TYPE_DOUBLE");
478: is_newtag = destroy_tag;
479: }
481: vmoab->tag = tag;
482: vmoab->new_tag = is_newtag;
483: }
484: vmoab->mbiface = mbiface;
485: #ifdef MOAB_HAVE_MPI
486: vmoab->pcomm = pcomm;
487: #endif
488: vmoab->is_global_vec = is_global_vec;
489: vmoab->tag_size = dmmoab->bs;
491: if (vmoab->is_native_vec) {
493: /* Create the PETSc Vector directly and attach our functions accordingly */
494: if (!is_global_vec) {
495: /* This is an MPI Vector with ghosted padding */
496: VecCreateGhostBlock((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc,
497: dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], vec);
498: }
499: else {
500: /* This is an MPI/SEQ Vector */
501: VecCreate((((PetscObject)dm)->comm), vec);
502: VecSetSizes(*vec, dmmoab->numFields * dmmoab->nloc, PETSC_DECIDE);
503: VecSetBlockSize(*vec, dmmoab->bs);
504: VecSetType(*vec, VECMPI);
505: }
506: }
507: else {
508: /* Call tag_iterate. This will cause MOAB to allocate memory for the
509: tag data if it hasn't already happened */
510: merr = mbiface->tag_iterate(tag, range->begin(), range->end(), count, (void*&)data_ptr); MBERRNM(merr);
512: /* set the reference for vector range */
513: vmoab->tag_range = new moab::Range(*range);
514: merr = mbiface->tag_get_length(tag, dmmoab->numFields); MBERRNM(merr);
516: /* Create the PETSc Vector
517: Query MOAB mesh to check if there are any ghosted entities
518: -> if we do, create a ghosted vector to map correctly to the same layout
519: -> else, create a non-ghosted parallel vector */
520: if (!is_global_vec) {
521: /* This is an MPI Vector with ghosted padding */
522: VecCreateGhostBlockWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc,
523: dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], data_ptr, vec);
524: }
525: else {
526: /* This is an MPI Vector without ghosted padding */
527: VecCreateMPIWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * range->size(),
528: PETSC_DECIDE, data_ptr, vec);
529: }
530: }
531: VecSetFromOptions(*vec);
533: /* create a container and store the internal MOAB data for faster access based on Entities etc */
534: PetscContainer moabdata;
535: PetscContainerCreate(PETSC_COMM_WORLD, &moabdata);
536: PetscContainerSetPointer(moabdata, vmoab);
537: PetscContainerSetUserDestroy(moabdata, DMVecUserDestroy_Moab);
538: PetscObjectCompose((PetscObject) * vec, "MOABData", (PetscObject)moabdata);
539: (*vec)->ops->duplicate = DMVecDuplicate_Moab;
540: PetscContainerDestroy(&moabdata);
542: /* Vector created, manually set local to global mapping */
543: if (dmmoab->ltog_map) {
544: VecSetLocalToGlobalMapping(*vec, dmmoab->ltog_map);
545: }
547: /* set the DM reference to the vector */
548: VecSetDM(*vec, dm);
549: return(0);
550: }
552: /* DMVecCreateTagName_Moab_Private
553: *
554: * Creates a unique tag name that will be shared across processes. If
555: * pcomm is NULL, then this is a serial vector. A unique tag name
556: * will be returned in tag_name in either case.
557: *
558: * The tag names have the format _PETSC_VEC_N where N is some integer.
559: *
560: * NOTE: The tag_name is allocated in this routine; The user needs to free
561: * the character array.
562: */
563: #ifdef MOAB_HAVE_MPI
564: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, moab::ParallelComm *pcomm, char** tag_name)
565: #else
566: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, char** tag_name)
567: #endif
568: {
569: moab::ErrorCode mberr;
570: PetscErrorCode ierr;
571: PetscInt n, global_n;
572: moab::Tag indexTag;
575: const char* PVEC_PREFIX = "__PETSC_VEC_";
576: PetscMalloc1(PETSC_MAX_PATH_LEN, tag_name);
578: moab::EntityHandle rootset = mbiface->get_root_set();
580: /* Check to see if there are any PETSc vectors defined */
581: /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */
582: mberr = mbiface->tag_get_handle("__PETSC_VECS__", 1, moab::MB_TYPE_INTEGER, indexTag,
583: moab::MB_TAG_SPARSE | moab::MB_TAG_CREAT, 0); MBERRNM(mberr);
584: mberr = mbiface->tag_get_data(indexTag, &rootset, 1, &n);
585: if (mberr == moab::MB_TAG_NOT_FOUND) n = 0; /* this is the first temporary vector */
586: else MBERRNM(mberr);
588: /* increment the new value of n */
589: ++n;
591: #ifdef MOAB_HAVE_MPI
592: /* Make sure that n is consistent across all processes */
593: MPIU_Allreduce(&n, &global_n, 1, MPI_INT, MPI_MAX, pcomm->comm());
594: #else
595: global_n = n;
596: #endif
598: /* Set the new name accordingly and return */
599: PetscSNPrintf(*tag_name, PETSC_MAX_PATH_LEN - 1, "%s_%D", PVEC_PREFIX, global_n);
600: mberr = mbiface->tag_set_data(indexTag, &rootset, 1, (const void*)&global_n); MBERRNM(mberr);
601: return(0);
602: }
604: PETSC_EXTERN PetscErrorCode DMCreateGlobalVector_Moab(DM dm, Vec *gvec)
605: {
606: PetscErrorCode ierr;
607: DM_Moab *dmmoab = (DM_Moab*)dm->data;
612: DMCreateVector_Moab_Private(dm,NULL,dmmoab->vowned,PETSC_TRUE,PETSC_TRUE,gvec);
613: return(0);
614: }
616: PETSC_EXTERN PetscErrorCode DMCreateLocalVector_Moab(DM dm, Vec *lvec)
617: {
618: PetscErrorCode ierr;
619: DM_Moab *dmmoab = (DM_Moab*)dm->data;
624: DMCreateVector_Moab_Private(dm,NULL,dmmoab->vlocal,PETSC_FALSE,PETSC_TRUE,lvec);
625: return(0);
626: }
628: PetscErrorCode DMVecDuplicate_Moab(Vec x, Vec *y)
629: {
631: DM dm;
632: PetscContainer moabdata;
633: Vec_MOAB *vmoab;
639: /* Get the Vec_MOAB struct for the original vector */
640: PetscObjectQuery((PetscObject)x, "MOABData", (PetscObject*) &moabdata);
641: PetscContainerGetPointer(moabdata, (void**)&vmoab);
643: VecGetDM(x, &dm);
646: DMCreateVector_Moab_Private(dm,NULL,vmoab->tag_range,vmoab->is_global_vec,PETSC_TRUE,y);
647: VecSetDM(*y, dm);
648: return(0);
649: }
651: PetscErrorCode DMVecUserDestroy_Moab(void *user)
652: {
653: Vec_MOAB *vmoab = (Vec_MOAB*)user;
654: PetscErrorCode ierr;
655: moab::ErrorCode merr;
658: if (vmoab->new_tag && vmoab->tag) {
659: /* Tag was created via a call to VecDuplicate, delete the underlying tag in MOAB */
660: merr = vmoab->mbiface->tag_delete(vmoab->tag); MBERRNM(merr);
661: }
662: delete vmoab->tag_range;
663: vmoab->tag = NULL;
664: vmoab->mbiface = NULL;
665: #ifdef MOAB_HAVE_MPI
666: vmoab->pcomm = NULL;
667: #endif
668: PetscFree(vmoab);
669: return(0);
670: }
672: PETSC_EXTERN PetscErrorCode DMGlobalToLocalBegin_Moab(DM dm, Vec g, InsertMode mode, Vec l)
673: {
674: PetscErrorCode ierr;
675: DM_Moab *dmmoab = (DM_Moab*)dm->data;
678: VecScatterBegin(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
679: return(0);
680: }
682: PETSC_EXTERN PetscErrorCode DMGlobalToLocalEnd_Moab(DM dm, Vec g, InsertMode mode, Vec l)
683: {
684: PetscErrorCode ierr;
685: DM_Moab *dmmoab = (DM_Moab*)dm->data;
688: VecScatterEnd(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
689: return(0);
690: }
692: PETSC_EXTERN PetscErrorCode DMLocalToGlobalBegin_Moab(DM dm, Vec l, InsertMode mode, Vec g)
693: {
694: PetscErrorCode ierr;
695: DM_Moab *dmmoab = (DM_Moab*)dm->data;
698: VecScatterBegin(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
699: return(0);
700: }
702: PETSC_EXTERN PetscErrorCode DMLocalToGlobalEnd_Moab(DM dm, Vec l, InsertMode mode, Vec g)
703: {
704: PetscErrorCode ierr;
705: DM_Moab *dmmoab = (DM_Moab*)dm->data;
708: VecScatterEnd(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
709: return(0);
710: }