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) {
143: /* get the local representation of the arrays from Vectors */
144: DMGetLocalVector(dm, &vmoab->local);
145: DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
146: DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);
148: /* Get the Vec_MOAB struct for the original vector */
149: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
150: PetscContainerGetPointer(moabdata, (void**)&xmoab);
152: /* get the local representation of the arrays from Vectors */
153: VecGhostGetLocalForm(vmoab->local, &xmoab->local);
154: VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
155: VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
157: VecGetArray(xmoab->local, varray);
158: }
159: else {
161: /* Get the MOAB private data */
162: DMMoabGetVecTag(vec, &vtag);
164: #ifdef MOAB_HAVE_MPI
165: /* exchange the data into ghost cells first */
166: merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal); MBERRNM(merr);
167: #endif
169: PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);
171: /* Get the array data for local entities */
172: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
175: i = 0;
176: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
177: for (f = 0; f < dmmoab->numFields; f++, i++)
178: (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart]*dmmoab->numFields + f] = marray[i];
179: //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
180: }
181: }
182: return 0;
183: }
185: /*@C
186: DMMoabVecRestoreArray - Restores the writable direct access array obtained via DMMoabVecGetArray
188: Collective
190: Input Parameters:
191: + dm - The DMMoab object being set
192: . vec - The Vector whose underlying data is requested
193: - array - The local data array
195: Level: intermediate
197: .seealso: DMMoabVecGetArray(), DMMoabVecGetArrayRead(), DMMoabVecRestoreArrayRead()
198: @*/
199: PetscErrorCode DMMoabVecRestoreArray(DM dm, Vec vec, void* array)
200: {
201: DM_Moab *dmmoab;
202: moab::ErrorCode merr;
203: moab::Tag vtag;
204: PetscInt count, i, f;
205: PetscScalar **varray;
206: PetscScalar *marray;
207: PetscContainer moabdata;
208: Vec_MOAB *vmoab, *xmoab;
213: dmmoab = (DM_Moab*)dm->data;
215: /* Get the Vec_MOAB struct for the original vector */
216: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
217: PetscContainerGetPointer(moabdata, (void**)&vmoab);
219: /* Get the real scalar array handle */
220: varray = reinterpret_cast<PetscScalar**>(array);
222: if (vmoab->is_native_vec) {
224: /* Get the Vec_MOAB struct for the original vector */
225: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
226: PetscContainerGetPointer(moabdata, (void**)&xmoab);
228: /* get the local representation of the arrays from Vectors */
229: VecRestoreArray(xmoab->local, varray);
230: VecGhostUpdateBegin(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
231: VecGhostUpdateEnd(vmoab->local, ADD_VALUES, SCATTER_REVERSE);
232: VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);
234: /* restore local pieces */
235: DMLocalToGlobalBegin(dm, vmoab->local, INSERT_VALUES, vec);
236: DMLocalToGlobalEnd(dm, vmoab->local, INSERT_VALUES, vec);
237: DMRestoreLocalVector(dm, &vmoab->local);
238: }
239: else {
241: /* Get the MOAB private data */
242: DMMoabGetVecTag(vec, &vtag);
244: /* Get the array data for local entities */
245: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
248: i = 0;
249: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
250: for (f = 0; f < dmmoab->numFields; f++, i++)
251: marray[i] = (*varray)[dmmoab->lidmap[(PetscInt) * iter - dmmoab->seqstart] * dmmoab->numFields + f];
252: //marray[i] = (*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]];
253: }
255: #ifdef MOAB_HAVE_MPI
256: /* reduce the tags correctly -> should probably let the user choose how to reduce in the future
257: For all FEM residual based assembly calculations, MPI_SUM should serve well */
258: merr = dmmoab->pcomm->reduce_tags(vtag, MPI_SUM, *dmmoab->vlocal); MBERRV(dmmoab->mbiface, merr);
259: #endif
260: PetscFree(*varray);
261: }
262: return 0;
263: }
265: /*@C
266: 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
268: Collective
270: Input Parameters:
271: + dm - The DMMoab object being set
272: - vec - The Vector whose underlying data is requested
274: Output Parameter:
275: . array - The local data array
277: Level: intermediate
279: .seealso: DMMoabVecRestoreArrayRead(), DMMoabVecGetArray(), DMMoabVecRestoreArray()
280: @*/
281: PetscErrorCode DMMoabVecGetArrayRead(DM dm, Vec vec, void* array)
282: {
283: DM_Moab *dmmoab;
284: moab::ErrorCode merr;
285: PetscInt count, i, f;
286: moab::Tag vtag;
287: PetscScalar **varray;
288: PetscScalar *marray;
289: PetscContainer moabdata;
290: Vec_MOAB *vmoab, *xmoab;
295: dmmoab = (DM_Moab*)dm->data;
297: /* Get the Vec_MOAB struct for the original vector */
298: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
299: PetscContainerGetPointer(moabdata, (void**)&vmoab);
301: /* Get the real scalar array handle */
302: varray = reinterpret_cast<PetscScalar**>(array);
304: if (vmoab->is_native_vec) {
305: /* get the local representation of the arrays from Vectors */
306: DMGetLocalVector(dm, &vmoab->local);
307: DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vmoab->local);
308: DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vmoab->local);
310: /* Get the Vec_MOAB struct for the original vector */
311: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
312: PetscContainerGetPointer(moabdata, (void**)&xmoab);
314: /* get the local representation of the arrays from Vectors */
315: VecGhostGetLocalForm(vmoab->local, &xmoab->local);
316: VecGhostUpdateBegin(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
317: VecGhostUpdateEnd(vmoab->local, INSERT_VALUES, SCATTER_FORWARD);
318: VecGetArray(xmoab->local, varray);
319: }
320: else {
321: /* Get the MOAB private data */
322: DMMoabGetVecTag(vec, &vtag);
324: #ifdef MOAB_HAVE_MPI
325: /* exchange the data into ghost cells first */
326: merr = dmmoab->pcomm->exchange_tags(vtag, *dmmoab->vlocal); MBERRNM(merr);
327: #endif
328: PetscMalloc1((dmmoab->nloc + dmmoab->nghost) * dmmoab->numFields, varray);
330: /* Get the array data for local entities */
331: merr = dmmoab->mbiface->tag_iterate(vtag, dmmoab->vlocal->begin(), dmmoab->vlocal->end(), count, reinterpret_cast<void*&>(marray), false); MBERRNM(merr);
334: i = 0;
335: for (moab::Range::iterator iter = dmmoab->vlocal->begin(); iter != dmmoab->vlocal->end(); iter++) {
336: for (f = 0; f < dmmoab->numFields; f++, i++)
337: (*varray)[dmmoab->lidmap[(PetscInt)*iter - dmmoab->seqstart]*dmmoab->numFields + f] = marray[i];
338: //(*varray)[dmmoab->llmap[dmmoab->lidmap[((PetscInt)*iter-dmmoab->seqstart)]*dmmoab->numFields+f]]=marray[i];
339: }
340: }
341: return 0;
342: }
344: /*@C
345: DMMoabVecRestoreArrayRead - Restores the read-only direct access array obtained via DMMoabVecGetArray
347: Collective
349: Input Parameters:
350: + dm - The DMMoab object being set
351: . vec - The Vector whose underlying data is requested
352: - array - The local data array
354: Level: intermediate
356: .seealso: DMMoabVecGetArrayRead(), DMMoabVecGetArray(), DMMoabVecRestoreArray()
357: @*/
358: PetscErrorCode DMMoabVecRestoreArrayRead(DM dm, Vec vec, void* array)
359: {
360: PetscScalar **varray;
361: PetscContainer moabdata;
362: Vec_MOAB *vmoab, *xmoab;
368: /* Get the Vec_MOAB struct for the original vector */
369: PetscObjectQuery((PetscObject)vec, "MOABData", (PetscObject*) &moabdata);
370: PetscContainerGetPointer(moabdata, (void**)&vmoab);
372: /* Get the real scalar array handle */
373: varray = reinterpret_cast<PetscScalar**>(array);
375: if (vmoab->is_native_vec) {
376: /* Get the Vec_MOAB struct for the original vector */
377: PetscObjectQuery((PetscObject)vmoab->local, "MOABData", (PetscObject*) &moabdata);
378: PetscContainerGetPointer(moabdata, (void**)&xmoab);
380: /* restore the local representation of the arrays from Vectors */
381: VecRestoreArray(xmoab->local, varray);
382: VecGhostRestoreLocalForm(vmoab->local, &xmoab->local);
384: /* restore local pieces */
385: DMRestoreLocalVector(dm, &vmoab->local);
386: }
387: else {
388: /* Nothing to do but just free the memory allocated before */
389: PetscFree(*varray);
391: }
392: return 0;
393: }
395: PetscErrorCode DMCreateVector_Moab_Private(DM dm, moab::Tag tag, const moab::Range* userrange, PetscBool is_global_vec, PetscBool destroy_tag, Vec *vec)
396: {
397: moab::ErrorCode merr;
398: PetscBool is_newtag;
399: const moab::Range *range;
400: PetscInt count, lnative_vec, gnative_vec;
401: std::string ttname;
402: PetscScalar *data_ptr, *defaultvals;
404: Vec_MOAB *vmoab;
405: DM_Moab *dmmoab = (DM_Moab*)dm->data;
406: #ifdef MOAB_HAVE_MPI
407: moab::ParallelComm *pcomm = dmmoab->pcomm;
408: #endif
409: moab::Interface *mbiface = dmmoab->mbiface;
412: if (!userrange) range = dmmoab->vowned;
413: else range = userrange;
416: #ifndef USE_NATIVE_PETSCVEC
417: /* If the tag data is in a single sequence, use PETSc native vector since tag_iterate isn't useful anymore */
418: lnative_vec = (range->psize() - 1);
419: #else
420: lnative_vec = 1; /* NOTE: Testing PETSc vector: will force to create native vector all the time */
421: // lnative_vec=0; /* NOTE: Testing MOAB vector: will force to create MOAB tag_iterate based vector all the time */
422: #endif
424: #ifdef MOAB_HAVE_MPI
425: MPIU_Allreduce(&lnative_vec, &gnative_vec, 1, MPI_INT, MPI_MAX, (((PetscObject)dm)->comm));
426: #else
427: gnative_vec = lnative_vec;
428: #endif
430: /* Create the MOAB internal data object */
431: PetscNew(&vmoab);
432: vmoab->is_native_vec = (gnative_vec > 0 ? PETSC_TRUE : PETSC_FALSE);
434: if (!vmoab->is_native_vec) {
435: merr = moab::MB_SUCCESS;
436: if (tag != 0) merr = mbiface->tag_get_name(tag, ttname);
437: if (!ttname.length() || merr != moab::MB_SUCCESS) {
438: /* get the new name for the anonymous MOABVec -> the tag_name will be destroyed along with Tag */
439: char *tag_name = NULL;
440: #ifdef MOAB_HAVE_MPI
441: DMVecCreateTagName_Moab_Private(mbiface,pcomm,&tag_name);
442: #else
443: DMVecCreateTagName_Moab_Private(mbiface,&tag_name);
444: #endif
445: is_newtag = PETSC_TRUE;
447: /* Create the default value for the tag (all zeros) */
448: PetscCalloc1(dmmoab->numFields, &defaultvals);
450: /* Create the tag */
451: merr = mbiface->tag_get_handle(tag_name, dmmoab->numFields, moab::MB_TYPE_DOUBLE, tag,
452: moab::MB_TAG_DENSE | moab::MB_TAG_CREAT, defaultvals); MBERRNM(merr);
453: PetscFree(tag_name);
454: PetscFree(defaultvals);
455: }
456: else {
457: /* Make sure the tag data is of type "double" */
458: moab::DataType tag_type;
459: merr = mbiface->tag_get_data_type(tag, tag_type); MBERRNM(merr);
461: is_newtag = destroy_tag;
462: }
464: vmoab->tag = tag;
465: vmoab->new_tag = is_newtag;
466: }
467: vmoab->mbiface = mbiface;
468: #ifdef MOAB_HAVE_MPI
469: vmoab->pcomm = pcomm;
470: #endif
471: vmoab->is_global_vec = is_global_vec;
472: vmoab->tag_size = dmmoab->bs;
474: if (vmoab->is_native_vec) {
476: /* Create the PETSc Vector directly and attach our functions accordingly */
477: if (!is_global_vec) {
478: /* This is an MPI Vector with ghosted padding */
479: VecCreateGhostBlock((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc,dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], vec);
480: } else {
481: /* This is an MPI/SEQ Vector */
482: VecCreate((((PetscObject)dm)->comm), vec);
483: VecSetSizes(*vec, dmmoab->numFields * dmmoab->nloc, PETSC_DECIDE);
484: VecSetBlockSize(*vec, dmmoab->bs);
485: VecSetType(*vec, VECMPI);
486: }
487: } else {
488: /* Call tag_iterate. This will cause MOAB to allocate memory for the
489: tag data if it hasn't already happened */
490: merr = mbiface->tag_iterate(tag, range->begin(), range->end(), count, (void*&)data_ptr); MBERRNM(merr);
492: /* set the reference for vector range */
493: vmoab->tag_range = new moab::Range(*range);
494: merr = mbiface->tag_get_length(tag, dmmoab->numFields); MBERRNM(merr);
496: /* Create the PETSc Vector
497: Query MOAB mesh to check if there are any ghosted entities
498: -> if we do, create a ghosted vector to map correctly to the same layout
499: -> else, create a non-ghosted parallel vector */
500: if (!is_global_vec) {
501: /* This is an MPI Vector with ghosted padding */
502: VecCreateGhostBlockWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * dmmoab->nloc,dmmoab->numFields * dmmoab->n, dmmoab->nghost, &dmmoab->gsindices[dmmoab->nloc], data_ptr, vec);
503: } else {
504: /* This is an MPI Vector without ghosted padding */
505: VecCreateMPIWithArray((((PetscObject)dm)->comm), dmmoab->bs, dmmoab->numFields * range->size(),PETSC_DECIDE, data_ptr, vec);
506: }
507: }
508: VecSetFromOptions(*vec);
510: /* create a container and store the internal MOAB data for faster access based on Entities etc */
511: PetscContainer moabdata;
512: PetscContainerCreate(PETSC_COMM_WORLD, &moabdata);
513: PetscContainerSetPointer(moabdata, vmoab);
514: PetscContainerSetUserDestroy(moabdata, DMVecUserDestroy_Moab);
515: PetscObjectCompose((PetscObject) * vec, "MOABData", (PetscObject)moabdata);
516: (*vec)->ops->duplicate = DMVecDuplicate_Moab;
517: PetscContainerDestroy(&moabdata);
519: /* Vector created, manually set local to global mapping */
520: if (dmmoab->ltog_map) {
521: VecSetLocalToGlobalMapping(*vec, dmmoab->ltog_map);
522: }
524: /* set the DM reference to the vector */
525: VecSetDM(*vec, dm);
526: return 0;
527: }
529: /* DMVecCreateTagName_Moab_Private
530: *
531: * Creates a unique tag name that will be shared across processes. If
532: * pcomm is NULL, then this is a serial vector. A unique tag name
533: * will be returned in tag_name in either case.
534: *
535: * The tag names have the format _PETSC_VEC_N where N is some integer.
536: *
537: * NOTE: The tag_name is allocated in this routine; The user needs to free
538: * the character array.
539: */
540: #ifdef MOAB_HAVE_MPI
541: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, moab::ParallelComm *pcomm, char** tag_name)
542: #else
543: PetscErrorCode DMVecCreateTagName_Moab_Private(moab::Interface *mbiface, char** tag_name)
544: #endif
545: {
546: moab::ErrorCode mberr;
547: PetscInt n, global_n;
548: moab::Tag indexTag;
550: const char* PVEC_PREFIX = "__PETSC_VEC_";
551: PetscMalloc1(PETSC_MAX_PATH_LEN, tag_name);
553: moab::EntityHandle rootset = mbiface->get_root_set();
555: /* Check to see if there are any PETSc vectors defined */
556: /* Create a tag in MOAB mesh to index and keep track of number of Petsc vec tags */
557: mberr = mbiface->tag_get_handle("__PETSC_VECS__", 1, moab::MB_TYPE_INTEGER, indexTag,
558: moab::MB_TAG_SPARSE | moab::MB_TAG_CREAT, 0); MBERRNM(mberr);
559: mberr = mbiface->tag_get_data(indexTag, &rootset, 1, &n);
560: if (mberr == moab::MB_TAG_NOT_FOUND) n = 0; /* this is the first temporary vector */
561: else MBERRNM(mberr);
563: /* increment the new value of n */
564: ++n;
566: #ifdef MOAB_HAVE_MPI
567: /* Make sure that n is consistent across all processes */
568: MPIU_Allreduce(&n, &global_n, 1, MPI_INT, MPI_MAX, pcomm->comm());
569: #else
570: global_n = n;
571: #endif
573: /* Set the new name accordingly and return */
574: PetscSNPrintf(*tag_name, PETSC_MAX_PATH_LEN - 1, "%s_%D", PVEC_PREFIX, global_n);
575: mberr = mbiface->tag_set_data(indexTag, &rootset, 1, (const void*)&global_n); MBERRNM(mberr);
576: return 0;
577: }
579: PETSC_EXTERN PetscErrorCode DMCreateGlobalVector_Moab(DM dm, Vec *gvec)
580: {
581: DM_Moab *dmmoab = (DM_Moab*)dm->data;
585: DMCreateVector_Moab_Private(dm,NULL,dmmoab->vowned,PETSC_TRUE,PETSC_TRUE,gvec);
586: return 0;
587: }
589: PETSC_EXTERN PetscErrorCode DMCreateLocalVector_Moab(DM dm, Vec *lvec)
590: {
591: DM_Moab *dmmoab = (DM_Moab*)dm->data;
595: DMCreateVector_Moab_Private(dm,NULL,dmmoab->vlocal,PETSC_FALSE,PETSC_TRUE,lvec);
596: return 0;
597: }
599: PetscErrorCode DMVecDuplicate_Moab(Vec x, Vec *y)
600: {
601: DM dm;
602: PetscContainer moabdata;
603: Vec_MOAB *vmoab;
608: /* Get the Vec_MOAB struct for the original vector */
609: PetscObjectQuery((PetscObject)x, "MOABData", (PetscObject*) &moabdata);
610: PetscContainerGetPointer(moabdata, (void**)&vmoab);
612: VecGetDM(x, &dm);
615: DMCreateVector_Moab_Private(dm,NULL,vmoab->tag_range,vmoab->is_global_vec,PETSC_TRUE,y);
616: VecSetDM(*y, dm);
617: return 0;
618: }
620: PetscErrorCode DMVecUserDestroy_Moab(void *user)
621: {
622: Vec_MOAB *vmoab = (Vec_MOAB*)user;
623: moab::ErrorCode merr;
625: if (vmoab->new_tag && vmoab->tag) {
626: /* Tag was created via a call to VecDuplicate, delete the underlying tag in MOAB */
627: merr = vmoab->mbiface->tag_delete(vmoab->tag); MBERRNM(merr);
628: }
629: delete vmoab->tag_range;
630: vmoab->tag = NULL;
631: vmoab->mbiface = NULL;
632: #ifdef MOAB_HAVE_MPI
633: vmoab->pcomm = NULL;
634: #endif
635: PetscFree(vmoab);
636: return 0;
637: }
639: PETSC_EXTERN PetscErrorCode DMGlobalToLocalBegin_Moab(DM dm, Vec g, InsertMode mode, Vec l)
640: {
641: DM_Moab *dmmoab = (DM_Moab*)dm->data;
643: VecScatterBegin(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
644: return 0;
645: }
647: PETSC_EXTERN PetscErrorCode DMGlobalToLocalEnd_Moab(DM dm, Vec g, InsertMode mode, Vec l)
648: {
649: DM_Moab *dmmoab = (DM_Moab*)dm->data;
651: VecScatterEnd(dmmoab->ltog_sendrecv, g, l, mode, SCATTER_REVERSE);
652: return 0;
653: }
655: PETSC_EXTERN PetscErrorCode DMLocalToGlobalBegin_Moab(DM dm, Vec l, InsertMode mode, Vec g)
656: {
657: DM_Moab *dmmoab = (DM_Moab*)dm->data;
659: VecScatterBegin(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
660: return 0;
661: }
663: PETSC_EXTERN PetscErrorCode DMLocalToGlobalEnd_Moab(DM dm, Vec l, InsertMode mode, Vec g)
664: {
665: DM_Moab *dmmoab = (DM_Moab*)dm->data;
667: VecScatterEnd(dmmoab->ltog_sendrecv, l, g, mode, SCATTER_FORWARD);
668: return 0;
669: }