xref: /petsc/src/dm/impls/moab/dmmoab.cxx (revision aa768e4c55e7f3620f7ce62e5a108dc2c69d8893)
1 #include <petsc-private/dmimpl.h> /*I  "petscdm.h"   I*/
2 #include <petsc-private/vecimpl.h> /*I  "petscdm.h"   I*/
3 
4 #include <petscdmmoab.h>
5 #include "MBTagConventions.hpp"
6 
7 typedef struct {
8   PetscInt bs; /* Number of degrees of freedom on each entity, aka tag size in moab */
9   PetscBool icreatedinstance; /* true if DM created moab instance internally, will destroy instance in DMDestroy */
10   moab::ParallelComm *pcomm;
11   moab::Interface *mbiface;
12   moab::Tag ltog_tag; /* moab supports "global id" tags, which are usually local to global numbering */
13   moab::Range range;
14 } DM_Moab;
15 
16 typedef struct {
17   moab::Interface    *mbiface;
18   moab::ParallelComm *pcomm;
19   moab::Range         tag_range; /* entities to which this tag applies */
20   moab::Tag           tag;
21   moab::Tag           ltog_tag;
22   PetscInt            tag_size;
23   PetscBool           new_tag;
24   PetscBool           serial;
25 
26 } Vec_MOAB;
27 
28 #undef __FUNCT__
29 #define __FUNCT__ "DMCreateGlobalVector_Moab"
30 PetscErrorCode DMCreateGlobalVector_Moab(DM dm,Vec *gvec)
31 {
32   PetscErrorCode  ierr;
33   DM_Moab         *dmmoab = (DM_Moab*)dm->data;
34 
35   PetscFunctionBegin;
36   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
37   PetscValidPointer(gvec,2);
38   PetscInt block_size = ((DM_Moab*)dm->data)->bs;
39   moab::Tag tag = 0;
40   ierr = DMMoabCreateVector(dm,tag,block_size,dmmoab->range,PETSC_FALSE,PETSC_TRUE,gvec);CHKERRQ(ierr);
41   PetscFunctionReturn(0);
42 }
43 
44 
45 #undef __FUNCT__
46 #define __FUNCT__ "DMCreateLocalVector_Moab"
47 PetscErrorCode DMCreateLocalVector_Moab(DM dm,Vec *gvec)
48 {
49   PetscErrorCode  ierr;
50   DM_Moab         *dmmoab = (DM_Moab*)dm->data;
51 
52   PetscFunctionBegin;
53   PetscInt bs = 1;
54   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
55   PetscValidPointer(gvec,2);
56   moab::Tag tag = 0;
57   ierr = DMMoabCreateVector(dm,tag,bs,dmmoab->range,PETSC_TRUE,PETSC_TRUE,gvec);CHKERRQ(ierr);
58   PetscFunctionReturn(0);
59 }
60 
61 EXTERN_C_BEGIN
62 #undef __FUNCT__
63 #define __FUNCT__ "DMCreate_Moab"
64 PetscErrorCode DMCreate_Moab(DM dm)
65 {
66   DM_Moab        *moab;
67   PetscErrorCode ierr;
68 
69   PetscFunctionBegin;
70   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
71   ierr     = PetscNewLog(dm, DM_Moab, &moab);CHKERRQ(ierr);
72   dm->data = moab;
73 
74   PetscFunctionReturn(0);
75 }
76 EXTERN_C_END
77 
78 #undef __FUNCT__
79 #define __FUNCT__ "DMDestroy_Moab"
80 PetscErrorCode DMDestroy_Moab(DM dm)
81 {
82   PetscFunctionBegin;
83   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
84 
85   // Delete the DM_Moab:
86   if(dm->data) {
87     if (((DM_Moab*)dm->data)->icreatedinstance) {
88       delete ((DM_Moab*)dm->data)->mbiface;
89       ((DM_Moab*)dm->data)->mbiface = NULL;
90       ((DM_Moab*)dm->data)->pcomm = NULL;
91     }
92     delete (DM_Moab*)dm->data;
93     dm->data = NULL;
94   }
95   PetscFunctionReturn(0);
96 }
97 
98 
99 #undef __FUNCT__
100 #define __FUNCT__ "DMMoabCreate"
101 /*@
102   DMMoabCreate - Creates a DMMoab object, which encapsulates a moab instance
103 
104   Collective on MPI_Comm
105 
106   Input Parameter:
107 . comm - The communicator for the DMMoab object
108 
109   Output Parameter:
110 . moab  - The DMMoab object
111 
112   Level: beginner
113 
114 .keywords: DMMoab, create
115 @*/
116 PetscErrorCode DMMoabCreate(MPI_Comm comm, DM *moab)
117 {
118   PetscErrorCode ierr;
119 
120   PetscFunctionBegin;
121   PetscValidPointer(moab,2);
122   ierr = DMCreate(comm, moab);CHKERRQ(ierr);
123   ierr = DMSetType(*moab, DMMOAB);CHKERRQ(ierr);
124   PetscFunctionReturn(0);
125 }
126 
127 #undef __FUNCT__
128 #define __FUNCT__ "DMMoabCreateMoab"
129 /*@
130   DMMoabCreate - Creates a DMMoab object, optionally from an instance and other data
131 
132   Collective on MPI_Comm
133 
134   Input Parameter:
135 . comm - The communicator for the DMMoab object
136 . moab - (ptr to) the MOAB Instance; if passed in NULL, MOAB instance is created inside PETSc, and destroyed
137          along with the DMMoab
138 . pcomm - (ptr to) a ParallelComm; if NULL, creates one internally for the whole communicator
139 . ltog_tag - A tag to use to retrieve global id for an entity; if 0, will use GLOBAL_ID_TAG_NAME/tag
140 . range - If non-NULL, contains range of entities to which DOFs will be assigned
141 
142   Output Parameter:
143 . moab  - The DMMoab object
144 
145   Level: beginner
146 
147 .keywords: DMMoab, create
148 @*/
149 PetscErrorCode DMMoabCreateMoab(MPI_Comm comm, moab::Interface *mbiface, moab::ParallelComm *pcomm, moab::Tag ltog_tag, moab::Range *range, DM *moab)
150 {
151   PetscErrorCode ierr;
152 
153   PetscFunctionBegin;
154   PetscValidPointer(moab,2);
155   ierr = DMCreate(comm, moab);CHKERRQ(ierr);
156   ierr = DMSetType(*moab, DMMOAB);CHKERRQ(ierr);
157 
158   if (!mbiface) {
159     mbiface = new moab::Core();
160     ((DM_Moab*)(*moab)->data)->icreatedinstance = PETSC_TRUE;
161   }
162   if (!pcomm) {
163     PetscInt rank, nprocs;
164     MPI_Comm_rank(comm, &rank);
165     MPI_Comm_size(comm, &nprocs);
166     pcomm = new moab::ParallelComm(mbiface, comm);
167   }
168 
169     // do the initialization of the DM
170   DM_Moab *dmmoab = new DM_Moab;
171   (*moab)->data      = dmmoab;
172   dmmoab->bs       = 0;
173   dmmoab->pcomm    = pcomm;
174   dmmoab->mbiface    = mbiface;
175   dmmoab->ltog_tag = ltog_tag;
176 
177     // initialize various functions
178   (*moab)->ops->createglobalvector              = DMCreateGlobalVector_Moab;
179   (*moab)->ops->createlocalvector               = DMCreateLocalVector_Moab;
180   (*moab)->ops->destroy                         = DMDestroy_Moab;
181 
182   ierr = DMMoabSetInterface(*moab, mbiface);CHKERRQ(ierr);
183   if (!pcomm) pcomm = new moab::ParallelComm(mbiface, comm);
184   ierr = DMMoabSetParallelComm(*moab, pcomm);CHKERRQ(ierr);
185   if (!ltog_tag) {
186     moab::ErrorCode merr = mbiface->tag_get_handle(GLOBAL_ID_TAG_NAME, ltog_tag);MBERRNM(merr);
187   }
188   if (ltog_tag) {
189     ierr = DMMoabSetLocalToGlobalTag(*moab, ltog_tag);CHKERRQ(ierr);
190   }
191   if (range) {
192     ierr = DMMoabSetRange(*moab, *range);CHKERRQ(ierr);
193   }
194   PetscFunctionReturn(0);
195 }
196 
197 #undef __FUNCT__
198 #define __FUNCT__ "DMMoabSetParallelComm"
199 /*@
200   DMMoabSetParallelComm - Set the ParallelComm used with this DMMoab
201 
202   Collective on MPI_Comm
203 
204   Input Parameter:
205 . dm    - The DMMoab object being set
206 . pcomm - The ParallelComm being set on the DMMoab
207 
208   Level: beginner
209 
210 .keywords: DMMoab, create
211 @*/
212 PetscErrorCode DMMoabSetParallelComm(DM dm,moab::ParallelComm *pcomm)
213 {
214   PetscFunctionBegin;
215   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
216   ((DM_Moab*)dm->data)->pcomm = pcomm;
217   ((DM_Moab*)dm->data)->mbiface = pcomm->get_moab();
218   PetscFunctionReturn(0);
219 }
220 
221 
222 #undef __FUNCT__
223 #define __FUNCT__ "DMMoabGetParallelComm"
224 /*@
225   DMMoabGetParallelComm - Get the ParallelComm used with this DMMoab
226 
227   Collective on MPI_Comm
228 
229   Input Parameter:
230 . dm    - The DMMoab object being set
231 
232   Output Parameter:
233 . pcomm - The ParallelComm for the DMMoab
234 
235   Level: beginner
236 
237 .keywords: DMMoab, create
238 @*/
239 PetscErrorCode DMMoabGetParallelComm(DM dm,moab::ParallelComm **pcomm)
240 {
241   PetscFunctionBegin;
242   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
243   *pcomm = ((DM_Moab*)dm->data)->pcomm;
244   PetscFunctionReturn(0);
245 }
246 
247 
248 #undef __FUNCT__
249 #define __FUNCT__ "DMMoabSetInterface"
250 /*@
251   DMMoabSetInterface - Set the MOAB instance used with this DMMoab
252 
253   Collective on MPI_Comm
254 
255   Input Parameter:
256 . dm      - The DMMoab object being set
257 . mbiface - The MOAB instance being set on this DMMoab
258 
259   Level: beginner
260 
261 .keywords: DMMoab, create
262 @*/
263 PetscErrorCode DMMoabSetInterface(DM dm,moab::Interface *mbiface)
264 {
265   PetscFunctionBegin;
266   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
267   ((DM_Moab*)dm->data)->pcomm = NULL;
268   ((DM_Moab*)dm->data)->mbiface = mbiface;
269   PetscFunctionReturn(0);
270 }
271 
272 
273 #undef __FUNCT__
274 #define __FUNCT__ "DMMoabGetInterface"
275 /*@
276   DMMoabGetInterface - Get the MOAB instance used with this DMMoab
277 
278   Collective on MPI_Comm
279 
280   Input Parameter:
281 . dm      - The DMMoab object being set
282 
283   Output Parameter:
284 . mbiface - The MOAB instance set on this DMMoab
285 
286   Level: beginner
287 
288 .keywords: DMMoab, create
289 @*/
290 PetscErrorCode DMMoabGetInterface(DM dm,moab::Interface **mbiface)
291 {
292   PetscFunctionBegin;
293   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
294   *mbiface = ((DM_Moab*)dm->data)->mbiface;
295   PetscFunctionReturn(0);
296 }
297 
298 
299 #undef __FUNCT__
300 #define __FUNCT__ "DMMoabSetRange"
301 /*@
302   DMMoabSetRange - Set the entities having DOFs on this DMMoab
303 
304   Collective on MPI_Comm
305 
306   Input Parameter:
307 . dm    - The DMMoab object being set
308 . range - The entities treated by this DMMoab
309 
310   Level: beginner
311 
312 .keywords: DMMoab, create
313 @*/
314 PetscErrorCode DMMoabSetRange(DM dm,moab::Range range)
315 {
316   PetscFunctionBegin;
317   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
318   ((DM_Moab*)dm->data)->range = range;
319   PetscFunctionReturn(0);
320 }
321 
322 
323 #undef __FUNCT__
324 #define __FUNCT__ "DMMoabGetRange"
325 /*@
326   DMMoabGetRange - Get the entities having DOFs on this DMMoab
327 
328   Collective on MPI_Comm
329 
330   Input Parameter:
331 . dm    - The DMMoab object being set
332 
333   Output Parameter:
334 . range - The entities treated by this DMMoab
335 
336   Level: beginner
337 
338 .keywords: DMMoab, create
339 @*/
340 PetscErrorCode DMMoabGetRange(DM dm,moab::Range *range)
341 {
342   PetscFunctionBegin;
343   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
344   *range = ((DM_Moab*)dm->data)->range;
345   PetscFunctionReturn(0);
346 }
347 
348 #undef __FUNCT__
349 #define __FUNCT__ "DMMoabSetLocalToGlobalTag"
350 /*@
351   DMMoabSetLocalToGlobalTag - Set the tag used for local to global numbering
352 
353   Collective on MPI_Comm
354 
355   Input Parameter:
356 . dm      - The DMMoab object being set
357 . ltogtag - The MOAB tag used for local to global ids
358 
359   Level: beginner
360 
361 .keywords: DMMoab, create
362 @*/
363 PetscErrorCode DMMoabSetLocalToGlobalTag(DM dm,moab::Tag ltogtag)
364 {
365   PetscFunctionBegin;
366   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
367   ((DM_Moab*)dm->data)->ltog_tag = ltogtag;
368   PetscFunctionReturn(0);
369 }
370 
371 
372 #undef __FUNCT__
373 #define __FUNCT__ "DMMoabGetLocalToGlobalTag"
374 /*@
375   DMMoabGetLocalToGlobalTag - Get the tag used for local to global numbering
376 
377   Collective on MPI_Comm
378 
379   Input Parameter:
380 . dm      - The DMMoab object being set
381 
382   Output Parameter:
383 . ltogtag - The MOAB tag used for local to global ids
384 
385   Level: beginner
386 
387 .keywords: DMMoab, create
388 @*/
389 PetscErrorCode DMMoabGetLocalToGlobalTag(DM dm,moab::Tag *ltog_tag)
390 {
391   PetscFunctionBegin;
392   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
393   *ltog_tag = ((DM_Moab*)dm->data)->ltog_tag;
394   PetscFunctionReturn(0);
395 }
396 
397 
398 #undef __FUNCT__
399 #define __FUNCT__ "DMMoabSetBlockSize"
400 /*@
401   DMMoabSetBlockSize - Set the block size used with this DMMoab
402 
403   Collective on MPI_Comm
404 
405   Input Parameter:
406 . dm - The DMMoab object being set
407 . bs - The block size used with this DMMoab
408 
409   Level: beginner
410 
411 .keywords: DMMoab, create
412 @*/
413 PetscErrorCode DMMoabSetBlockSize(DM dm,PetscInt bs)
414 {
415   PetscFunctionBegin;
416   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
417   ((DM_Moab*)dm->data)->bs = bs;
418   PetscFunctionReturn(0);
419 }
420 
421 
422 #undef __FUNCT__
423 #define __FUNCT__ "DMMoabGetBlockSize"
424 /*@
425   DMMoabGetBlockSize - Get the block size used with this DMMoab
426 
427   Collective on MPI_Comm
428 
429   Input Parameter:
430 . dm - The DMMoab object being set
431 
432   Output Parameter:
433 . bs - The block size used with this DMMoab
434 
435   Level: beginner
436 
437 .keywords: DMMoab, create
438 @*/
439 PetscErrorCode DMMoabGetBlockSize(DM dm,PetscInt *bs)
440 {
441   PetscFunctionBegin;
442   PetscValidHeaderSpecific(dm,DM_CLASSID,1);
443   *bs = ((DM_Moab*)dm->data)->bs;
444   PetscFunctionReturn(0);
445 }
446 
447 
448 // declare for use later but before they're defined
449 PetscErrorCode DMMoab_VecUserDestroy(void *user);
450 PetscErrorCode DMMoab_VecDuplicate(Vec x,Vec *y);
451 PetscErrorCode DMMoab_CreateTagName(const moab::ParallelComm *pcomm,std::string& tag_name);
452 PetscErrorCode DMMoab_CreateVector(moab::Interface *iface,moab::ParallelComm *pcomm,moab::Tag tag,PetscInt tag_size,moab::Tag ltog_tag,moab::Range range,PetscBool serial, PetscBool destroy_tag,Vec *vec);
453 
454 #undef __FUNCT__
455 #define __FUNCT__ "DMMoabCreateVector"
456 /*@
457   DMMoabCreateVector - Create a Vec from either an existing tag, or a specified tag size, and a range of entities
458 
459   Collective on MPI_Comm
460 
461   Input Parameter:
462 . dm          - The DMMoab object being set
463 . tag         - If non-zero, block size will be taken from the tag size
464 . tag_size    - If tag was zero, this parameter specifies the block size; unique tag name will be generated automatically
465 . range       - If non-empty, Vec corresponds to these entities, otherwise to the entities set on the DMMoab
466 . serial      - If true, this is a serial Vec, otherwise a parallel one
467 . destroy_tag - If true, MOAB tag is destroyed with Vec, otherwise it is left on MOAB
468 
469   Output Parameter:
470 . vec         - The created vector
471 
472   Level: beginner
473 
474 .keywords: DMMoab, create
475 @*/
476 PetscErrorCode DMMoabCreateVector(DM dm,moab::Tag tag,PetscInt tag_size,moab::Range range,PetscBool serial, PetscBool destroy_tag,Vec *vec)
477 {
478   PetscErrorCode     ierr;
479 
480   PetscFunctionBegin;
481 
482   DM_Moab *dmmoab = (DM_Moab*)dm->data;
483   moab::ParallelComm *pcomm = dmmoab->pcomm;
484   moab::Interface *mbiface = dmmoab->mbiface;
485   moab::Tag ltog_tag = dmmoab->ltog_tag;
486 
487   if (!tag && !tag_size) {
488     PetscFunctionReturn(PETSC_ERR_ARG_WRONG);
489   }
490   else if (!tag && tag_size) {
491     ierr = DMMoab_CreateVector(mbiface,pcomm,tag,tag_size,ltog_tag,range,serial,destroy_tag,vec);CHKERRQ(ierr);
492   }
493   PetscFunctionReturn(0);
494 }
495 
496 
497 #undef __FUNCT__
498 #define __FUNCT__ "DMMoab_CreateVector"
499 PetscErrorCode DMMoab_CreateVector(moab::Interface *mbiface,moab::ParallelComm *pcomm,moab::Tag tag,PetscInt tag_size,moab::Tag ltog_tag,moab::Range range,PetscBool serial, PetscBool destroy_tag,Vec *vec)
500 {
501   PetscErrorCode     ierr;
502   moab::ErrorCode    merr;
503 
504   PetscFunctionBegin;
505 
506   if (!tag) {
507     std::string tag_name;
508     ierr = DMMoab_CreateTagName(pcomm,tag_name);CHKERRQ(ierr);
509 
510       // Create the default value for the tag (all zeros):
511     std::vector<PetscScalar> default_value(tag_size, 0.0);
512 
513       // Create the tag:
514     merr = mbiface->tag_get_handle(tag_name.c_str(),tag_size,moab::MB_TYPE_DOUBLE,tag,
515                                    moab::MB_TAG_DENSE | moab::MB_TAG_CREAT,default_value.data());MBERRNM(merr);
516   }
517   else {
518 
519       // Make sure the tag data is of type "double":
520     moab::DataType tag_type;
521     merr = mbiface->tag_get_data_type(tag, tag_type);MBERRNM(merr);
522     if(tag_type != moab::MB_TYPE_DOUBLE) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Tag data type must be MB_TYPE_DOUBLE");
523   }
524 
525     // Create the MOAB internal data object
526   Vec_MOAB *vmoab;
527   ierr = PetscMalloc(sizeof(Vec_MOAB),&vmoab);CHKERRQ(ierr);
528   new (vmoab) Vec_MOAB();
529   vmoab->tag = tag;
530   vmoab->ltog_tag = ltog_tag;
531   vmoab->mbiface = mbiface;
532   vmoab->pcomm = pcomm;
533   vmoab->tag_range = range;
534   vmoab->new_tag = destroy_tag;
535   vmoab->serial = serial;
536   merr = mbiface->tag_get_length(tag,vmoab->tag_size);MBERR("tag_get_size", merr);
537 
538     // Call tag_iterate. This will cause MOAB to allocate memory for the
539     // tag data if it hasn't already happened:
540   int  count;
541   void *void_ptr;
542   merr = mbiface->tag_iterate(tag,range.begin(),range.end(),count,void_ptr);MBERRNM(merr);
543 
544     // Check to make sure the tag data is in a single sequence:
545   if ((unsigned)count != range.size()) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Can only create MOAB Vector for single sequence");
546   PetscScalar *data_ptr = (PetscScalar*)void_ptr;
547 
548     // Create the PETSc Vector:
549   if(!serial) {
550       // This is an MPI Vector:
551     ierr = VecCreateMPIWithArray(vmoab->pcomm->comm(),vmoab->tag_size,vmoab->tag_size*range.size(),
552                                  PETSC_DECIDE,data_ptr,vec);CHKERRXX(ierr);
553 
554       // Vector created, manually set local to global mapping:
555     ISLocalToGlobalMapping ltog;
556     PetscInt               *gindices = new PetscInt[range.size()];
557     PetscInt               count = 0;
558     moab::Range::iterator  iter;
559     for(iter = range.begin(); iter != range.end(); iter++) {
560       int dof;
561       merr = mbiface->tag_get_data(ltog_tag,&(*iter),1,&dof);MBERRNM(merr);
562       gindices[count] = dof;
563       count++;
564     }
565 
566     ierr = ISLocalToGlobalMappingCreate(PETSC_COMM_SELF,range.size(),gindices,
567                                         PETSC_COPY_VALUES,&ltog);CHKERRQ(ierr);
568     ierr = VecSetLocalToGlobalMappingBlock(*vec,ltog);CHKERRQ(ierr);
569 
570       // Clean up:
571     ierr = ISLocalToGlobalMappingDestroy(&ltog);CHKERRQ(ierr);
572     delete [] gindices;
573   } else {
574       // This is a serial vector:
575     ierr = VecCreateSeqWithArray(PETSC_COMM_SELF,vmoab->tag_size,vmoab->tag_size*range.size(),data_ptr,vec);CHKERRXX(ierr);
576   }
577 
578 
579   PetscContainer moabdata;
580   ierr = PetscContainerCreate(PETSC_COMM_SELF,&moabdata);CHKERRQ(ierr);
581   ierr = PetscContainerSetPointer(moabdata,vmoab);CHKERRQ(ierr);
582   ierr = PetscContainerSetUserDestroy(moabdata,DMMoab_VecUserDestroy);CHKERRQ(ierr);
583   ierr = PetscObjectCompose((PetscObject)*vec,"MOABData",(PetscObject)moabdata);CHKERRQ(ierr);
584   (*vec)->ops->duplicate = DMMoab_VecDuplicate;
585 
586   ierr = PetscContainerDestroy(&moabdata);CHKERRQ(ierr);
587   PetscFunctionReturn(0);
588 }
589 
590 #undef __FUNCT__
591 #define __FUNCT__ "DMMoabGetVecTag"
592 /*@
593   DMMoabGetVecTag - Get the MOAB tag associated with this Vec
594 
595   Collective on MPI_Comm
596 
597   Input Parameter:
598 . vec - Vec being queried
599 
600   Output Parameter:
601 . tag - Tag associated with this Vec
602 
603   Level: beginner
604 
605 .keywords: DMMoab, create
606 @*/
607 PetscErrorCode DMMoabGetVecTag(Vec vec,moab::Tag *tag)
608 {
609   PetscContainer  moabdata;
610   Vec_MOAB        *vmoab;
611   PetscErrorCode  ierr;
612 
613   PetscFunctionBegin;
614 
615   // Get the MOAB private data:
616   ierr = PetscObjectQuery((PetscObject)vec,"MOABData", (PetscObject*) &moabdata);CHKERRQ(ierr);
617   ierr = PetscContainerGetPointer(moabdata, (void**) &vmoab);CHKERRQ(ierr);
618 
619   *tag = vmoab->tag;
620 
621   PetscFunctionReturn(0);
622 }
623 
624 
625 #undef __FUNCT__
626 #define __FUNCT__ "DMMoabGetVecRange"
627 /*@
628   DMMoabGetVecRange - Get the MOAB entities associated with this Vec
629 
630   Collective on MPI_Comm
631 
632   Input Parameter:
633 . vec   - Vec being queried
634 
635   Output Parameter:
636 . range - Entities associated with this Vec
637 
638   Level: beginner
639 
640 .keywords: DMMoab, create
641 @*/
642 PetscErrorCode DMMoabGetVecRange(Vec vec,moab::Range *range)
643 {
644   PetscContainer  moabdata;
645   Vec_MOAB        *vmoab;
646   PetscErrorCode  ierr;
647 
648   PetscFunctionBegin;
649 
650   // Get the MOAB private data:
651   ierr = PetscObjectQuery((PetscObject)vec,"MOABData", (PetscObject*) &moabdata);CHKERRQ(ierr);
652   ierr = PetscContainerGetPointer(moabdata, (void**) &vmoab);CHKERRQ(ierr);
653 
654   *range = vmoab->tag_range;
655 
656   PetscFunctionReturn(0);
657 }
658 
659 
660 #undef __FUNCT__
661 #define __FUNCT__ "DMMoab_VecDuplicate"
662 PetscErrorCode DMMoab_VecDuplicate(Vec x,Vec *y)
663 {
664   PetscErrorCode ierr;
665   PetscFunctionBegin;
666   PetscValidHeaderSpecific(x,VEC_CLASSID,1);
667   PetscValidPointer(y,2);
668 
669   // Get the Vec_MOAB struct for the original vector:
670   PetscContainer  moabdata;
671   Vec_MOAB        *vmoab;
672   ierr = PetscObjectQuery((PetscObject)x,"MOABData", (PetscObject*) &moabdata);CHKERRQ(ierr);
673   ierr = PetscContainerGetPointer(moabdata, (void**)&vmoab);CHKERRQ(ierr);
674 
675   ierr = DMMoab_CreateVector(vmoab->mbiface,vmoab->pcomm,vmoab->tag, vmoab->tag_size,vmoab->ltog_tag,vmoab->tag_range,vmoab->serial,PETSC_TRUE,y);CHKERRQ(ierr);
676   PetscFunctionReturn(0);
677 }
678 
679 
680 #undef __FUNCT__
681 #define __FUNCT__ "DMMoab_CreateTagName"
682 /*  DMMoab_CreateTagName
683  *
684  *  Creates a unique tag name that will be shared across processes. If
685  *  pcomm is NULL, then this is a serial vector. A unique tag name
686  *  will be returned in tag_name in either case.
687  *
688  *  The tag names have the format _PETSC_VEC_N where N is some integer.
689  */
690 PetscErrorCode DMMoab_CreateTagName(const moab::ParallelComm *pcomm,std::string& tag_name)
691 {
692   moab::ErrorCode mberr;
693   PetscErrorCode  ierr;
694 
695   PetscFunctionBegin;
696   const std::string PVEC_PREFIX      = "_PETSC_VEC_";
697   const PetscInt    PVEC_PREFIX_SIZE = PVEC_PREFIX.size();
698 
699   // Check to see if there are any PETSc vectors defined:
700   const moab::Interface  *mbiface = pcomm->get_moab();
701   std::vector<moab::Tag> tags;
702   PetscInt               n = 0;
703   mberr = mbiface->tag_get_tags(tags);MBERRNM(mberr);
704   for(unsigned i = 0; i < tags.size(); i++) {
705     std::string s;
706     mberr = mbiface->tag_get_name(tags[i],s);MBERRNM(mberr);
707     if(s.find(PVEC_PREFIX) != std::string::npos){
708       // This tag represents a PETSc vector. Find the vector number:
709       PetscInt m;
710       std::istringstream(s.substr(PVEC_PREFIX_SIZE)) >> m;
711       if(m >= n) n = m+1;
712     }
713   }
714 
715   // Make sure that n is consistent across all processes:
716   PetscInt global_n;
717   MPI_Comm comm = PETSC_COMM_SELF;
718   if(pcomm) comm = pcomm->comm();
719   ierr = MPI_Allreduce(&n,&global_n,1,MPI_INT,MPI_MAX,comm);CHKERRQ(ierr);
720 
721   // Set the answer and return:
722   std::stringstream ss;
723   ss << PVEC_PREFIX << global_n;
724   tag_name = ss.str();
725   PetscFunctionReturn(0);
726 }
727 
728 
729 #undef __FUNCT__
730 #define __FUNCT__ "DMMoab_VecUserDestroy"
731 PetscErrorCode DMMoab_VecUserDestroy(void *user)
732 {
733   Vec_MOAB        *vmoab;
734   PetscErrorCode  ierr;
735   moab::ErrorCode merr;
736 
737   PetscFunctionBegin;
738   vmoab = (Vec_MOAB*)user;
739   if(vmoab->new_tag == PETSC_TRUE) {
740     // Tag created via a call to VecDuplicate, delete the underlying tag in MOAB...
741     merr = vmoab->mbiface->tag_delete(vmoab->tag);MBERRNM(merr);
742   }
743 
744   ierr = PetscFree(vmoab);CHKERRQ(ierr);
745   PetscFunctionReturn(0);
746 }
747 
748