Merge pull request #283 from coprigent/feature/pmmg-edge-tag-consistency

Fix storage of tags in `MMG5_swap23`
MmgTools · Sep 24, 2024 · 8a84488 · 8a84488
2 parents f6b3408 + d97c353
commit 8a84488
Show file tree

Hide file tree

Showing 2 changed files with 41 additions and 16 deletions.
diff --git a/src/mmg3d/chkmsh_3d.c b/src/mmg3d/chkmsh_3d.c
@@ -241,8 +241,10 @@ void MMG3D_chkmeshedgestags(MMG5_pMesh mesh) {
 
         if ( !ier ) {
            /* First time we meet the edge: store the its tag from the current
-            * tetra in the hash table */
-          int ier2 = MMG5_hEdge ( mesh,&hash,ip1,ip2,0,pxt->tag[i]);
+            * tetra in the hash table. Ignore OLDPARBDY tag because it is not
+            * consistent through meshes inside ParMmg and forbid the use of the
+            * current function to check tag consistency if not ignored. */
+          int ier2 = MMG5_hEdge ( mesh,&hash,ip1,ip2,0,(pxt->tag[i] & ~MG_OLDPARBDY));
           if ( !ier2 ) {
             /* Realloc error */
             fprintf(stderr,"Error: %s: %d: Unable to add to hash table the edge "
@@ -252,8 +254,12 @@ void MMG3D_chkmeshedgestags(MMG5_pMesh mesh) {
           }
         }
         else {
-          /* Edge tag has been stored from another tet: check consistency */
-          if ( tag != pxt->tag[i] ) {
+          /* Edge tag has been stored from another tet: check consistency.
+             Ignore OLDPARBDY tag because it is not
+            * consistent through meshes inside ParMmg and forbid the use of the
+            * current function to check tag consistency if not ignored.
+            */
+          if ( tag != (pxt->tag[i] & ~MG_OLDPARBDY) ) {
             fprintf(stderr,"Error: %s: %d: Non consistency at tet %" MMG5_PRId
                     " (%" MMG5_PRId "), edge %d:%" MMG5_PRId "--%" MMG5_PRId "\n ",
                   __func__,__LINE__,k,MMG3D_indElt(mesh,k),i,ip1,ip2);

diff --git a/src/mmg3d/swap_3d.c b/src/mmg3d/swap_3d.c
@@ -739,9 +739,20 @@ int MMG3D_swap23(MMG5_pMesh mesh,MMG5_pSol met,MMG5_int k,int8_t metRidTyp,
   }
 
   /** Swap */
+  /* Store useful information from pt1 before overwrite by memcpy*/
   xt1 = pt1->xt;
 
   np    = pt1->v[tau1[0]];
+
+  MMG5_int ref[3] = {0};
+  uint16_t  tag[3] = {0};
+  for (i=0;i<3;i++) {
+    if ( !MMG3D_get_shellEdgeTag(mesh,k1,taued1[i],&tag[i],&ref[i]) ) {
+      fprintf(stderr,"\n  ## Error: %s: %d. unable to get edge info.\n",__func__,i);
+      return 0;
+    }
+  }
+
   memcpy(pt1,pt0,sizeof(MMG5_Tetra));
 
   iel = MMG3D_newElt(mesh);
@@ -879,14 +890,20 @@ int MMG3D_swap23(MMG5_pMesh mesh,MMG5_pSol met,MMG5_int k,int8_t metRidTyp,
     pxt1 = &mesh->xtetra[xt1];
 
     /* Assignation of the xt fields to the appropriate tets */
+    /* Warning: after collapses, some boundary edges not connected to boundary
+     * faces may have a 0 tag inside a xtetra (see \ref MMG5_colver when a
+     * xtetra is assigned to one of the neighbours of the tetra of the edge
+     * shell). In consequence, we cannot simply use the stored tags. */
+
     /* xt[0] */
     xt[0].tag[taued0[0]] = 0;
-    xt[0].tag[taued0[3]] = pxt1->tag[taued1[2]];
-    xt[0].tag[taued0[4]] = pxt1->tag[taued1[1]];
+
+    xt[0].tag[taued0[3]] = tag[2];
+    xt[0].tag[taued0[4]] = tag[1];
 
     xt[0].edg[taued0[0]] = 0;
-    xt[0].edg[taued0[3]] =  pxt1->edg[taued1[2]];
-    xt[0].edg[taued0[4]] =  pxt1->edg[taued1[1]];
+    xt[0].edg[taued0[3]] = ref[2];
+    xt[0].edg[taued0[4]] = ref[1];
 
     xt[0].ref[ tau0[0]] = pxt1->ref[tau1[1]];
     xt[0].ref[ tau0[2]] = 0;
@@ -901,12 +918,13 @@ int MMG3D_swap23(MMG5_pMesh mesh,MMG5_pSol met,MMG5_int k,int8_t metRidTyp,
 
     /* xt[1] */
     xt[1].tag[taued0[1]] = 0;
-    xt[1].tag[taued0[3]] = pxt1->tag[taued1[0]];
-    xt[1].tag[taued0[5]] = pxt1->tag[taued1[2]];
+
+    xt[1].tag[taued0[3]] = tag[0];
+    xt[1].tag[taued0[5]] = tag[1];
 
     xt[1].edg[taued0[1]] = 0;
-    xt[1].edg[taued0[3]] =  pxt1->edg[taued1[0]];
-    xt[1].edg[taued0[5]] =  pxt1->edg[taued1[2]];
+    xt[1].edg[taued0[3]] = ref[0];
+    xt[1].edg[taued0[5]] = ref[1];
 
     xt[1].ref[ tau0[0]] = pxt1->ref[tau1[3]];
     xt[1].ref[ tau0[1]] = 0;
@@ -921,12 +939,13 @@ int MMG3D_swap23(MMG5_pMesh mesh,MMG5_pSol met,MMG5_int k,int8_t metRidTyp,
 
     /* xt[2] */
     xt[2].tag[taued0[2]] = 0;
-    xt[2].tag[taued0[4]] = pxt1->tag[taued1[0]];
-    xt[2].tag[taued0[5]] = pxt1->tag[taued1[1]];
+
+    xt[2].tag[taued0[4]] = tag[0];
+    xt[2].tag[taued0[5]] = tag[2];
 
     xt[2].edg[taued0[2]] = 0;
-    xt[2].edg[taued0[4]] = pxt1->edg[taued1[0]];
-    xt[2].edg[taued0[5]] = pxt1->edg[taued1[1]];
+    xt[2].edg[taued0[4]] = ref[0];
+    xt[2].edg[taued0[5]] = ref[2];
 
     xt[2].ref[ tau0[0]] = pxt1->ref[tau1[2]];
     xt[2].ref[ tau0[1]] = 0;