xref: /petsc/src/dm/impls/moab/dmmbmg.cxx (revision 9c36898578a91ec77675797fbfe6e623f4be5691)
1 #include <petsc/private/dmmbimpl.h>
2 #include <petscdmmoab.h>
3 #include <MBTagConventions.hpp>
4 #include <moab/NestedRefine.hpp>
5 
6 #undef __FUNCT__
7 #define __FUNCT__ "DMMoabGenerateHierarchy"
8 /*@
9   DMMoabGenerateHierarchy - Generate a multi-level uniform refinement hierarchy
10   by succesively refining a coarse mesh, already defined in the DM object
11   provided by the user.
12 
13   Collective on MPI_Comm
14 
15   Input Parameter:
16 + dmb  - The DMMoab object
17 
18   Output Parameter:
19 + nlevels   - The number of levels of refinement needed to generate the hierarchy
20 . ldegrees  - The degree of refinement at each level in the hierarchy
21 
22   Level: beginner
23 
24 .keywords: DMMoab, create, refinement
25 @*/
26 PetscErrorCode DMMoabGenerateHierarchy(DM dm, PetscInt nlevels, PetscInt *ldegrees)
27 {
28   DM_Moab        *dmmoab;
29   PetscErrorCode  ierr;
30   moab::ErrorCode merr;
31   PetscInt *pdegrees, ilevel;
32   std::vector<moab::EntityHandle> hsets;
33 
34   PetscFunctionBegin;
35   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
36   dmmoab = (DM_Moab*)(dm)->data;
37 
38   if (!ldegrees) {
39     ierr = PetscMalloc1(nlevels, &pdegrees);CHKERRQ(ierr);
40     for (ilevel = 0; ilevel < nlevels; ilevel++) pdegrees[ilevel] = 2; /* default = Degree 2 refinement */
41   }
42   else pdegrees = ldegrees;
43 
44   /* initialize set level refinement data for hierarchy */
45   dmmoab->nhlevels = nlevels;
46 
47   /* Instantiate the nested refinement class */
48 #ifdef MOAB_HAVE_MPI
49   dmmoab->hierarchy = new moab::NestedRefine(dynamic_cast<moab::Core*>(dmmoab->mbiface), dmmoab->pcomm, dmmoab->fileset);
50 #else
51   dmmoab->hierarchy = new moab::NestedRefine(dynamic_cast<moab::Core*>(dmmoab->mbiface), NULL, dmmoab->fileset);
52 #endif
53 
54   ierr = PetscMalloc1(nlevels + 1, &dmmoab->hsets);CHKERRQ(ierr);
55 
56   /* generate the mesh hierarchy */
57   merr = dmmoab->hierarchy->generate_mesh_hierarchy(nlevels, pdegrees, hsets, false);MBERRNM(merr);
58 
59 #ifdef MOAB_HAVE_MPI
60   if (dmmoab->pcomm->size() > 1) {
61     merr = dmmoab->hierarchy->exchange_ghosts(hsets, dmmoab->nghostrings);MBERRNM(merr);
62   }
63 #endif
64 
65   /* copy the mesh sets for nested refinement hierarchy */
66   for (ilevel = 0; ilevel <= nlevels; ilevel++)
67   {
68     if (ilevel) { /* Update material and other geometric tags from parent to child sets */
69       merr = dmmoab->hierarchy->update_special_tags(ilevel, hsets[ilevel]);MBERRNM(merr);
70     }
71 
72     dmmoab->hsets[ilevel] = hsets[ilevel];
73 #ifdef MOAB_HAVE_MPI
74     if (ilevel) {
75       merr = dmmoab->pcomm->assign_global_ids(hsets[ilevel], dmmoab->dim, 1, false, true, false);MBERRNM(merr);
76     }
77 #endif
78   }
79 
80   hsets.clear();
81   if (!ldegrees) {
82     ierr = PetscFree(pdegrees);CHKERRQ(ierr);
83   }
84   PetscFunctionReturn(0);
85 }
86 
87 #undef __FUNCT__
88 #define __FUNCT__ "DMRefineHierarchy_Moab"
89 /*@
90   DMRefineHierarchy_Moab - Generate a multi-level DM hierarchy
91   by succesively refining a coarse mesh.
92 
93   Collective on MPI_Comm
94 
95   Input Parameter:
96 + dm  - The DMMoab object
97 
98   Output Parameter:
99 + nlevels   - The number of levels of refinement needed to generate the hierarchy
100 . dmf  - The DM objects after successive refinement of the hierarchy
101 
102   Level: beginner
103 
104 .keywords: DMMoab, generate, refinement
105 @*/
106 PETSC_EXTERN PetscErrorCode  DMRefineHierarchy_Moab(DM dm, PetscInt nlevels, DM dmf[])
107 {
108   PetscErrorCode  ierr;
109   PetscInt        i;
110 
111   PetscFunctionBegin;
112 
113   ierr = DMRefine(dm, PetscObjectComm((PetscObject)dm), &dmf[0]);CHKERRQ(ierr);
114   for (i = 1; i < nlevels; i++) {
115     ierr = DMRefine(dmf[i - 1], PetscObjectComm((PetscObject)dm), &dmf[i]);CHKERRQ(ierr);
116   }
117   PetscFunctionReturn(0);
118 }
119 
120 #undef __FUNCT__
121 #define __FUNCT__ "DMCoarsenHierarchy_Moab"
122 /*@
123   DMCoarsenHierarchy_Moab - Generate a multi-level DM hierarchy
124   by succesively coarsening a refined mesh.
125 
126   Collective on MPI_Comm
127 
128   Input Parameter:
129 + dm  - The DMMoab object
130 
131   Output Parameter:
132 + nlevels   - The number of levels of refinement needed to generate the hierarchy
133 . dmc  - The DM objects after successive coarsening of the hierarchy
134 
135   Level: beginner
136 
137 .keywords: DMMoab, generate, coarsening
138 @*/
139 PETSC_EXTERN PetscErrorCode DMCoarsenHierarchy_Moab(DM dm, PetscInt nlevels, DM dmc[])
140 {
141   PetscErrorCode  ierr;
142   PetscInt        i;
143 
144   PetscFunctionBegin;
145 
146   ierr = DMCoarsen(dm, PetscObjectComm((PetscObject)dm), &dmc[0]);CHKERRQ(ierr);
147   for (i = 1; i < nlevels; i++) {
148     ierr = DMCoarsen(dmc[i - 1], PetscObjectComm((PetscObject)dm), &dmc[i]);CHKERRQ(ierr);
149   }
150   PetscFunctionReturn(0);
151 }
152 
153 PETSC_EXTERN PetscErrorCode DMMoab_Compute_NNZ_From_Connectivity(DM, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscBool);
154 
155 #undef __FUNCT__
156 #define __FUNCT__ "DMCreateInterpolation_Moab"
157 /*@
158   DMCreateInterpolation_Moab - Generate the interpolation operators to transform
159   operators (matrices, vectors) from parent level to child level as defined by
160   the DM inputs provided by the user.
161 
162   Collective on MPI_Comm
163 
164   Input Parameter:
165 + dm1  - The DMMoab object
166 - dm2  - the second, finer DMMoab object
167 
168   Output Parameter:
169 + interpl  - The interpolation operator for transferring data between the levels
170 - vec      - The scaling vector (optional)
171 
172   Level: developer
173 
174 .keywords: DMMoab, create, refinement
175 @*/
176 PETSC_EXTERN PetscErrorCode DMCreateInterpolation_Moab(DM dmp, DM dmc, Mat* interpl, Vec* vec)
177 {
178   DM_Moab         *dmbp, *dmbc;
179   PetscErrorCode   ierr;
180   moab::ErrorCode  merr;
181   PetscInt         dim;
182   PetscReal        factor;
183   PetscInt         innz, *nnz, ionz, *onz;
184   PetscInt         nlsizp, nlsizc, nlghsizp, ngsizp, ngsizc;
185   const PetscBool  use_consistent_bases=PETSC_TRUE;
186 
187   PetscFunctionBegin;
188   PetscValidHeaderSpecific(dmp, DM_CLASSID, 1);
189   PetscValidHeaderSpecific(dmc, DM_CLASSID, 2);
190   dmbp = (DM_Moab*)(dmp)->data;
191   dmbc = (DM_Moab*)(dmc)->data;
192   nlsizp = dmbp->nloc;// *dmb1->numFields;
193   nlsizc = dmbc->nloc;// *dmb2->numFields;
194   ngsizp = dmbp->n; // *dmb1->numFields;
195   ngsizc = dmbc->n; // *dmb2->numFields;
196   nlghsizp = (dmbp->nloc + dmbp->nghost); // *dmb1->numFields;
197 
198   // Columns = Parent DoFs ;  Rows = Child DoFs
199   // Interpolation matrix: \sum_{i=1}^P Owned(Child) * (Owned(Parent) + Ghosted(Parent))
200   // Size: nlsizc * nlghsizp
201   PetscInfo4(NULL, "Creating interpolation matrix %D X %D to apply transformation between levels %D -> %D.\n", ngsizc, nlghsizp, dmbp->hlevel, dmbc->hlevel);
202 
203   ierr = DMGetDimension(dmp, &dim);CHKERRQ(ierr);
204 
205   /* allocate the nnz, onz arrays based on block size and local nodes */
206   ierr = PetscCalloc2(nlsizc, &nnz, nlsizc, &onz);CHKERRQ(ierr);
207 
208   /* Loop through the local elements and compute the relation between the current parent and the refined_level. */
209   for (moab::Range::iterator iter = dmbc->vowned->begin(); iter != dmbc->vowned->end(); iter++) {
210 
211     const moab::EntityHandle vhandle = *iter;
212     /* define local variables */
213     moab::EntityHandle parent;
214     std::vector<moab::EntityHandle> adjs;
215     moab::Range     found;
216 
217     /* store the vertex DoF number */
218     const int ldof = dmbc->lidmap[vhandle - dmbc->seqstart];
219 
220     /* Get adjacency information for current vertex - i.e., all elements of dimension (dim) that connects
221        to the current vertex. We can then decipher if a vertex is ghosted or not and compute the
222        non-zero pattern accordingly. */
223     merr = dmbc->hierarchy->get_adjacencies(vhandle, dmbc->dim, adjs);MBERRNM(merr);
224 
225     /* loop over vertices and update the number of connectivity */
226     for (unsigned jter = 0; jter < adjs.size(); jter++) {
227 
228       const moab::EntityHandle jhandle = adjs[jter];
229 
230       /* Get the relation between the current (coarse) parent and its corresponding (finer) children elements */
231       merr = dmbc->hierarchy->child_to_parent(jhandle, dmbc->hlevel, dmbp->hlevel, &parent);MBERRNM(merr);
232 
233       /* Get connectivity information in canonical ordering for the local element */
234       std::vector<moab::EntityHandle> connp;
235       merr = dmbp->hierarchy->get_connectivity(parent, dmbp->hlevel, connp);MBERRNM(merr);
236 
237       for (unsigned ic=0; ic < connp.size(); ++ic) {
238 
239         /* loop over each element connected to the adjacent vertex and update as needed */
240         /* find the truly user-expected layer of ghosted entities to decipher NNZ pattern */
241         if (found.find(connp[ic]) != found.end()) continue; /* make sure we don't double count shared vertices */
242         if (dmbp->vghost->find(connp[ic]) != dmbp->vghost->end()) onz[ldof]++; /* update out-of-proc onz */
243         else nnz[ldof]++; /* else local vertex */
244         found.insert(connp[ic]);
245       }
246     }
247   }
248 
249   for (int i = 0; i < nlsizc; i++)
250     nnz[i] += 1; /* self count the node */
251 
252   ionz = onz[0];
253   innz = nnz[0];
254   for (int tc = 0; tc < nlsizc; tc++) {
255     // check for maximum allowed sparsity = fully dense
256     nnz[tc] = std::min(nlsizp, nnz[tc]);
257     onz[tc] = std::min(ngsizp - nlsizp, onz[tc]);
258 
259     PetscInfo3(NULL, "  %d: NNZ = %d, ONZ = %d\n", tc, nnz[tc], onz[tc]);
260 
261     innz = (innz < nnz[tc] ? nnz[tc] : innz);
262     ionz = (ionz < onz[tc] ? onz[tc] : ionz);
263   }
264 
265   /* create interpolation matrix */
266   ierr = MatCreate(PetscObjectComm((PetscObject)dmc), interpl);CHKERRQ(ierr);
267   ierr = MatSetSizes(*interpl, nlsizc, nlsizp, ngsizc, ngsizp);CHKERRQ(ierr);
268   ierr = MatSetType(*interpl, MATAIJ);CHKERRQ(ierr);
269   ierr = MatSetFromOptions(*interpl);CHKERRQ(ierr);
270 
271   ierr = MatSeqAIJSetPreallocation(*interpl, innz, nnz);CHKERRQ(ierr);
272   ierr = MatMPIAIJSetPreallocation(*interpl, innz, nnz, ionz, onz);CHKERRQ(ierr);
273 
274   /* clean up temporary memory */
275   ierr = PetscFree2(nnz, onz);CHKERRQ(ierr);
276 
277   /* set up internal matrix data-structures */
278   ierr = MatSetUp(*interpl);CHKERRQ(ierr);
279 
280   /* Define variables for assembly */
281   std::vector<moab::EntityHandle> children;
282   std::vector<moab::EntityHandle> connp, connc;
283   std::vector<PetscReal> pcoords, ccoords, values_phi;
284 
285   if (use_consistent_bases) {
286     const moab::EntityHandle ehandle = dmbp->elocal->front();
287 
288     merr = dmbp->hierarchy->parent_to_child(ehandle, dmbp->hlevel, dmbc->hlevel, children);MBERRNM(merr);
289 
290     /* Get connectivity and coordinates of the parent vertices */
291     merr = dmbp->hierarchy->get_connectivity(ehandle, dmbp->hlevel, connp);MBERRNM(merr);
292     merr = dmbc->mbiface->get_connectivity(&children[0], children.size(), connc);MBERRNM(merr);
293 
294     std::vector<PetscReal> natparam(3*connc.size(), 0.0);
295     pcoords.resize(connp.size() * 3);
296     ccoords.resize(connc.size() * 3);
297     values_phi.resize(connp.size()*connc.size());
298     /* Get coordinates for connectivity entities in canonical order for both coarse and finer levels */
299     merr = dmbp->hierarchy->get_coordinates(&connp[0], connp.size(), dmbp->hlevel, &pcoords[0]);MBERRNM(merr);
300     merr = dmbc->hierarchy->get_coordinates(&connc[0], connc.size(), dmbc->hlevel, &ccoords[0]);MBERRNM(merr);
301 
302     /* Set values: For each DOF in coarse grid cell, set the contribution or PHI evaluated at each fine grid DOF point */
303     for (unsigned tc = 0; tc < connc.size(); tc++) {
304       const PetscInt offset = tc * 3;
305 
306       /* Scale ccoords relative to pcoords */
307       ierr = DMMoabPToRMapping(dim, connp.size(), &pcoords[0], &ccoords[offset], &natparam[offset], &values_phi[connp.size()*tc]);CHKERRQ(ierr);
308     }
309   }
310   else {
311     factor = std::pow(2.0 /*degree_P_for_refinement*/, (dmbc->hlevel - dmbp->hlevel) * dmbp->dim * 1.0);
312   }
313 
314   /* TODO: Decipher the correct non-zero pattern. There is still some issue with onz allocation */
315   ierr = MatSetOption(*interpl, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE);CHKERRQ(ierr);
316 
317   /* Loop through the remaining vertices. These vertices appear only on the current refined_level. */
318   for (moab::Range::iterator iter = dmbp->elocal->begin(); iter != dmbp->elocal->end(); iter++) {
319 
320     const moab::EntityHandle ehandle = *iter;
321 
322     /* Get the relation between the current (coarse) parent and its corresponding (finer) children elements */
323     children.clear();
324     connc.clear();
325     merr = dmbp->hierarchy->parent_to_child(ehandle, dmbp->hlevel, dmbc->hlevel, children);MBERRNM(merr);
326 
327     /* Get connectivity and coordinates of the parent vertices */
328     merr = dmbp->hierarchy->get_connectivity(ehandle, dmbp->hlevel, connp);MBERRNM(merr);
329     merr = dmbc->mbiface->get_connectivity(&children[0], children.size(), connc);MBERRNM(merr);
330 
331     pcoords.resize(connp.size() * 3);
332     ccoords.resize(connc.size() * 3);
333     /* Get coordinates for connectivity entities in canonical order for both coarse and finer levels */
334     merr = dmbp->hierarchy->get_coordinates(&connp[0], connp.size(), dmbp->hlevel, &pcoords[0]);MBERRNM(merr);
335     merr = dmbc->hierarchy->get_coordinates(&connc[0], connc.size(), dmbc->hlevel, &ccoords[0]);MBERRNM(merr);
336 
337     std::vector<int> dofsp(connp.size()), dofsc(connc.size());
338     /* TODO: specific to scalar system - use GetDofs */
339     ierr = DMMoabGetDofsBlocked(dmp, connp.size(), &connp[0], &dofsp[0]);CHKERRQ(ierr);
340     ierr = DMMoabGetDofsBlocked(dmc, connc.size(), &connc[0], &dofsc[0]);CHKERRQ(ierr);
341 
342     /* Compute the actual interpolation weights when projecting solution/residual between levels */
343     if (use_consistent_bases) {
344 
345       /* Use the cached values of natural parameteric coordinates and basis pre-evaluated.
346          We are making an assumption here that UMR used in GMG to generate the hierarchy uses
347          the same template for all elements; This will fail for mixed element meshes (TRI/QUAD).
348 
349          TODO: Fix the above assumption by caching data for families (especially for Tets and mixed meshes)
350       */
351 
352       /* Set values: For each DOF in coarse grid cell, set the contribution or PHI evaluated at each fine grid DOF point */
353       for (unsigned tc = 0; tc < connc.size(); tc++) {
354         /* TODO: Check if we should be using INSERT_VALUES instead */
355         ierr = MatSetValues(*interpl, 1, &dofsc[tc], connp.size(), &dofsp[0], &values_phi[connp.size()*tc], ADD_VALUES);CHKERRQ(ierr);
356       }
357     }
358     else {
359       /* Compute the interpolation weights by determining distance of 1-ring
360          neighbor vertices from current vertex
361 
362          This should be used only when FEM basis is not used for the discretization.
363          Else, the consistent interface to compute the basis function for interpolation
364          between the levels should be evaluated correctly to preserve convergence of GMG.
365          Shephard's basis will be terrible for any unsmooth problems.
366       */
367       values_phi.resize(connp.size());
368       for (unsigned tc = 0; tc < connc.size(); tc++) {
369 
370         PetscReal normsum = 0.0;
371         for (unsigned tp = 0; tp < connp.size(); tp++) {
372           values_phi[tp] = 0.0;
373           for (unsigned k = 0; k < 3; k++)
374             values_phi[tp] += std::pow(pcoords[tp * 3 + k] - ccoords[k + tc * 3], dim);
375           if (values_phi[tp] < 1e-12) {
376             values_phi[tp] = 1e12;
377           }
378           else {
379             //values_phi[tp] = std::pow(values_phi[tp], -1.0/dim);
380             values_phi[tp] = std::pow(values_phi[tp], -1.0);
381             normsum += values_phi[tp];
382           }
383         }
384         for (unsigned tp = 0; tp < connp.size(); tp++) {
385           if (values_phi[tp] > 1e11)
386             values_phi[tp] = factor * 0.5 / connp.size();
387           else
388             values_phi[tp] = factor * values_phi[tp] * 0.5 / (connp.size() * normsum);
389         }
390         ierr = MatSetValues(*interpl, 1, &dofsc[tc], connp.size(), &dofsp[0], &values_phi[0], ADD_VALUES);CHKERRQ(ierr);
391       }
392     }
393   }
394   if (vec) *vec = NULL;
395   ierr = MatAssemblyBegin(*interpl, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
396   ierr = MatAssemblyEnd(*interpl, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
397   PetscFunctionReturn(0);
398 }
399 
400 #undef __FUNCT__
401 #define __FUNCT__ "DMCreateInjection_Moab"
402 /*@
403   DMCreateInjection_Moab - Generate a multi-level uniform refinement hierarchy
404   by succesively refining a coarse mesh, already defined in the DM object
405   provided by the user.
406 
407   Collective on MPI_Comm
408 
409   Input Parameter:
410 . dmb  - The DMMoab object
411 
412   Output Parameter:
413 . nlevels   - The number of levels of refinement needed to generate the hierarchy
414 + ldegrees  - The degree of refinement at each level in the hierarchy
415 
416   Level: beginner
417 
418 .keywords: DMMoab, create, refinement
419 @*/
420 PETSC_EXTERN PetscErrorCode DMCreateInjection_Moab(DM dm1, DM dm2, VecScatter* ctx)
421 {
422   //DM_Moab        *dmmoab;
423 
424   PetscFunctionBegin;
425   PetscValidHeaderSpecific(dm1, DM_CLASSID, 1);
426   PetscValidHeaderSpecific(dm2, DM_CLASSID, 2);
427   //dmmoab = (DM_Moab*)(dm1)->data;
428 
429   PetscPrintf(PETSC_COMM_WORLD, "[DMCreateInjection_Moab] :: Placeholder\n");
430   PetscFunctionReturn(0);
431 }
432 
433 
434 #undef __FUNCT__
435 #define __FUNCT__ "DM_UMR_Moab_Private"
436 PetscErrorCode  DM_UMR_Moab_Private(DM dm, MPI_Comm comm, PetscBool refine, DM *dmref)
437 {
438   PetscErrorCode  ierr;
439   PetscInt        i, dim;
440   DM              dm2;
441   moab::ErrorCode merr;
442   DM_Moab        *dmb = (DM_Moab*)dm->data, *dd2;
443 
444   PetscFunctionBegin;
445   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
446   PetscValidPointer(dmref, 4);
447 
448   if ( (dmb->hlevel == dmb->nhlevels && refine) || (dmb->hlevel == 0 && !refine) ) {
449     if (dmb->hlevel + 1 > dmb->nhlevels && refine) PetscInfo2(NULL, "Invalid multigrid refinement hierarchy level specified (%D). MOAB UMR max levels = %D. Creating a NULL object.\n", dmb->hlevel + 1, dmb->nhlevels);
450     if (dmb->hlevel - 1 < 0 && !refine) PetscInfo1(NULL, "Invalid multigrid coarsen hierarchy level specified (%D). Creating a NULL object.\n", dmb->hlevel - 1);
451     *dmref = PETSC_NULL;
452     PetscFunctionReturn(0);
453   }
454 
455   ierr = DMMoabCreate(PetscObjectComm((PetscObject)dm), &dm2);CHKERRQ(ierr);
456   dd2 = (DM_Moab*)dm2->data;
457 
458   dd2->mbiface = dmb->mbiface;
459 #ifdef MOAB_HAVE_MPI
460   dd2->pcomm = dmb->pcomm;
461 #endif
462   dd2->icreatedinstance = PETSC_FALSE;
463   dd2->nghostrings = dmb->nghostrings;
464 
465   /* set the new level based on refinement/coarsening */
466   if (refine) {
467     dd2->hlevel = dmb->hlevel + 1;
468   }
469   else {
470     dd2->hlevel = dmb->hlevel - 1;
471   }
472 
473   /* Copy the multilevel hierarchy pointers in MOAB */
474   dd2->hierarchy = dmb->hierarchy;
475   dd2->nhlevels = dmb->nhlevels;
476   ierr = PetscMalloc1(dd2->nhlevels + 1, &dd2->hsets);CHKERRQ(ierr);
477   for (i = 0; i <= dd2->nhlevels; i++) {
478     dd2->hsets[i] = dmb->hsets[i];
479   }
480   dd2->fileset = dd2->hsets[dd2->hlevel];
481 
482   /* do the remaining initializations for DMMoab */
483   dd2->bs = dmb->bs;
484   dd2->numFields = dmb->numFields;
485   dd2->rw_dbglevel = dmb->rw_dbglevel;
486   dd2->partition_by_rank = dmb->partition_by_rank;
487   ierr = PetscStrcpy(dd2->extra_read_options, dmb->extra_read_options);CHKERRQ(ierr);
488   ierr = PetscStrcpy(dd2->extra_write_options, dmb->extra_write_options);CHKERRQ(ierr);
489   dd2->read_mode = dmb->read_mode;
490   dd2->write_mode = dmb->write_mode;
491 
492   /* set global ID tag handle */
493   ierr = DMMoabSetLocalToGlobalTag(dm2, dmb->ltog_tag);CHKERRQ(ierr);
494 
495   merr = dd2->mbiface->tag_get_handle(MATERIAL_SET_TAG_NAME, dd2->material_tag);MBERRNM(merr);
496 
497   ierr = DMSetOptionsPrefix(dm2, ((PetscObject)dm)->prefix);CHKERRQ(ierr);
498   ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr);
499   ierr = DMSetDimension(dm2, dim);CHKERRQ(ierr);
500 
501   /* allow overloaded (user replaced) operations to be inherited by refinement clones */
502   dm2->ops->creatematrix = dm->ops->creatematrix;
503 
504   /* copy fill information if given */
505   ierr = DMMoabSetBlockFills(dm2, dmb->dfill, dmb->ofill);CHKERRQ(ierr);
506 
507   /* copy vector type information */
508   ierr = DMSetMatType(dm2, dm->mattype);CHKERRQ(ierr);
509   ierr = DMSetVecType(dm2, dm->vectype);CHKERRQ(ierr);
510   dd2->numFields = dmb->numFields;
511   if (dmb->numFields) {
512     ierr = DMMoabSetFieldNames(dm2, dmb->numFields, dmb->fieldNames);CHKERRQ(ierr);
513   }
514 
515   ierr = DMSetFromOptions(dm2);CHKERRQ(ierr);
516 
517   /* recreate Dof numbering for the refined DM and make sure the distribution is correctly populated */
518   ierr = DMSetUp(dm2);CHKERRQ(ierr);
519 
520   *dmref = dm2;
521   PetscFunctionReturn(0);
522 }
523 
524 
525 #undef __FUNCT__
526 #define __FUNCT__ "DMRefine_Moab"
527 /*@
528   DMRefine_Moab - Generate a multi-level uniform refinement hierarchy
529   by succesively refining a coarse mesh, already defined in the DM object
530   provided by the user.
531 
532   Collective on DM
533 
534   Input Parameter:
535 + dm  - The DMMoab object
536 - comm - the communicator to contain the new DM object (or MPI_COMM_NULL)
537 
538   Output Parameter:
539 . dmf - the refined DM, or NULL
540 
541   Note: If no refinement was done, the return value is NULL
542 
543   Level: developer
544 
545 .keywords: DMMoab, create, refinement
546 @*/
547 PETSC_EXTERN PetscErrorCode DMRefine_Moab(DM dm, MPI_Comm comm, DM* dmf)
548 {
549   PetscErrorCode  ierr;
550 
551   PetscFunctionBegin;
552   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
553 
554   ierr = DM_UMR_Moab_Private(dm, comm, PETSC_TRUE, dmf);CHKERRQ(ierr);
555   PetscFunctionReturn(0);
556 }
557 
558 #undef __FUNCT__
559 #define __FUNCT__ "DMCoarsen_Moab"
560 /*@
561   DMCoarsen_Moab - Generate a multi-level uniform refinement hierarchy
562   by succesively refining a coarse mesh, already defined in the DM object
563   provided by the user.
564 
565   Collective on DM
566 
567   Input Parameter:
568 . dm  - The DMMoab object
569 - comm - the communicator to contain the new DM object (or MPI_COMM_NULL)
570 
571   Output Parameter:
572 . dmf - the coarsened DM, or NULL
573 
574   Note: If no coarsening was done, the return value is NULL
575 
576   Level: developer
577 
578 .keywords: DMMoab, create, coarsening
579 @*/
580 PETSC_EXTERN PetscErrorCode DMCoarsen_Moab(DM dm, MPI_Comm comm, DM* dmc)
581 {
582   PetscErrorCode  ierr;
583 
584   PetscFunctionBegin;
585   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
586 
587   ierr = DM_UMR_Moab_Private(dm, comm, PETSC_FALSE, dmc);CHKERRQ(ierr);
588   PetscFunctionReturn(0);
589 }
590