Vijay Mahadevan avatar Vijay Mahadevan committed 741c5f9

Add boundary element storage for DM. Can instead set a bit tag in future to avoid storage. Use skinner and filter to get the necessary data for boundary entities since the previous method of find_skin alone fails in parallel due to internal faces also being returned.

Comments (0)

Files changed (3)


   moab::Range             *vlocal, *vowned, *vghost;
   moab::Range             *elocal, *eghost;
   moab::EntityHandle      fileset;
-  moab::Range             *bndyvtx,*bndyfaces;
+  moab::Range             *bndyvtx,*bndyfaces,*bndyelems;
   PetscInt                nfields;
   const char**            fields;


 PETSC_EXTERN PetscErrorCode DMMoabGetLocalToGlobalTag(DM dm,moab::Tag *ltog_tag);
 PETSC_EXTERN PetscErrorCode DMMoabSetBlockSize(DM dm,PetscInt bs);
 PETSC_EXTERN PetscErrorCode DMMoabGetBlockSize(DM dm,PetscInt *bs);
-PETSC_EXTERN PetscErrorCode DMMoabGetBoundaryEntities(DM dm,moab::Range *bdvtx,moab::Range* bdfaces);
+PETSC_EXTERN PetscErrorCode DMMoabGetBoundaryEntities(DM dm,moab::Range *bdvtx,moab::Range* bdfaces,moab::Range* bdelems);
 PETSC_EXTERN PetscErrorCode DMMoabVecGetArrayRead(DM,Vec,void*);
 PETSC_EXTERN PetscErrorCode DMMoabVecRestoreArrayRead(DM,Vec,void*);


   delete dmmoab->eghost;
   delete dmmoab->bndyvtx;
   delete dmmoab->bndyfaces;
+  delete dmmoab->bndyelems;
   ierr = PetscFree(dmmoab->gsindices);CHKERRQ(ierr);
   ierr = VecScatterDestroy(&dmmoab->ltog_sendrecv);CHKERRQ(ierr);
   PetscInt                totsize,local_min,local_max,global_min;
   PetscSection            section;
+  moab::Range adj;
   /* Get the local and shared vertices and cache it */
     /* find out the local and global minima of GLOBAL_ID */
-    for (i=1; i<totsize; ++i) {
+    for (i=0; i<totsize; ++i) {
 //      if (dmmoab->pcomm->rank())
 //        PetscPrintf(PETSC_COMM_SELF, "[%D] gsindices[%D] = %D\n", dmmoab->pcomm->rank(), i, dmmoab->gsindices[i]);
       if(local_min>dmmoab->gsindices[i]) local_min=dmmoab->gsindices[i];
           ierr = PetscSectionSetFieldOffset(section, locgid, i, (locgid-global_min)*dmmoab->nfields);
         else {
-          ierr = PetscSectionSetFieldDof(section, locgid, i, totsize*i+locgid-global_min);CHKERRQ(ierr);
-          ierr = PetscSectionSetFieldOffset(section, locgid, i, i*totsize);
-//          PetscPrintf(PETSC_COMM_SELF, "[%D] Local_GID = %D, FDOF = %D, OFF = %D.\n", dmmoab->pcomm->rank(), locgid, totsize*i+locgid-global_min, totsize );
+          ierr = PetscSectionSetFieldDof(section, locgid, i, dmmoab->n*i+locgid-global_min);CHKERRQ(ierr);
+          ierr = PetscSectionSetFieldOffset(section, locgid, i, i*dmmoab->n);
+          PetscPrintf(PETSC_COMM_SELF, "[%D] Local_GID = %D, FDOF = %D, OFF = %D.\n", dmmoab->pcomm->rank(), locgid, dmmoab->n*i+locgid-global_min, i*dmmoab->n );
       ierr = PetscSectionSetDof(section, locgid, dmmoab->nfields);CHKERRQ(ierr);
     moab::Skinner skinner(dmmoab->mbiface);
     dmmoab->bndyvtx = new moab::Range();
     dmmoab->bndyfaces = new moab::Range();
-    merr = skinner.find_skin(dmmoab->fileset, *dmmoab->elocal, true, *dmmoab->bndyvtx);MBERRNM(merr); // 'true' param indicates we want vertices back, not faces
-    merr = skinner.find_skin(dmmoab->fileset, *dmmoab->elocal, false, *dmmoab->bndyfaces);MBERRNM(merr); // 'false' param indicates we want faces back, not vertices
-    PetscInfo2(dm, "Found %D boundary vertices and %D faces.\n", dmmoab->bndyvtx->size(), dmmoab->bndyvtx->size());
+    dmmoab->bndyelems = new moab::Range();
+    merr = skinner.find_skin(dmmoab->fileset, *dmmoab->elocal, false, *dmmoab->bndyfaces, NULL, false, true, false, false);MBERRNM(merr); // 'false' param indicates we want faces back, not vertices
+//    merr = skinner.find_skin(dmmoab->fileset, *dmmoab->elocal, true, *dmmoab->bndyvtx, NULL, true, true, false, false);MBERRNM(merr); // 'true' param indicates we want vertices back, not faces
+    if (dmmoab->dim == 3) {
+      // get the edges from faces and then do the same if needed
+    }
+    merr = dmmoab->pcomm->filter_pstatus(*dmmoab->bndyfaces,PSTATUS_NOT_OWNED,PSTATUS_NOT);MBERRNM(merr);
+    merr = dmmoab->pcomm->filter_pstatus(*dmmoab->bndyfaces,PSTATUS_SHARED,PSTATUS_NOT);MBERRNM(merr);
+//    merr = dmmoab->pcomm->filter_pstatus(*dmmoab->bndyfaces,PSTATUS_MULTISHARED,PSTATUS_NOT);MBERRNM(merr);
+    merr = dmmoab->mbiface->get_connectivity(*dmmoab->bndyfaces, *dmmoab->bndyvtx, false);MBERRNM(ierr);
+    merr = dmmoab->mbiface->get_adjacencies(*dmmoab->bndyfaces, dmmoab->dim, false, *dmmoab->bndyelems, moab::Interface::UNION);MBERRNM(ierr);
+    dmmoab->bndyfaces->print(0);
+    dmmoab->bndyvtx->print(0);
+    dmmoab->bndyelems->print(0);
+//    merr = skinner.find_geometric_skin(dmmoab->fileset, *dmmoab->bndyelems, dmmoab->dim);MBERRV(dmmoab->mbiface,merr);
+//    merr = skinner.find_skin(dmmoab->fileset, *dmmoab->elocal, dmmoab->bndyvtx, dmmoab->bndyfaces, dmmoab->bndyelems, true, true);MBERRNM(merr);
+    PetscInfo3(dm, "Found %D boundary vertices, %D boundary faces and %D boundary elements.\n", dmmoab->bndyvtx->size(), dmmoab->bndyvtx->size(), dmmoab->bndyelems->size());
 #undef __FUNCT__
 #define __FUNCT__ "DMMoabGetBoundaryEntities"
-PetscErrorCode DMMoabGetBoundaryEntities(DM dm,moab::Range *bdvtx,moab::Range* bdfaces)
+PetscErrorCode DMMoabGetBoundaryEntities(DM dm,moab::Range *bdvtx,moab::Range* bdfaces,moab::Range* bdelems)
   DM_Moab        *dmmoab;
   if (bdvtx)  *bdvtx = *dmmoab->bndyvtx;
   if (bdfaces)  *bdfaces = *dmmoab->bndyfaces;
+  if (bdelems)  *bdfaces = *dmmoab->bndyelems;
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.