diff --git a/contrib/hxt/README.md b/contrib/hxt/README.md
index bab9f3cec87af6327c56c3015b4a21e8ac45f2ab..2fa03da350a54c5e93cd5b0df86eb58a2c4ab301 100644
--- a/contrib/hxt/README.md
+++ b/contrib/hxt/README.md
@@ -140,20 +140,18 @@ set_target_properties(my_target hxt_core hxt_predicates hxt_tetBR hxt_tetMesh
### TODO/future work
-core + others:
-
- - extend HXT to 40+ bits nodes, 32 bits colors
-
tetMesh:
- - compute Moore indices on the fly, only when few point to insert or few tetrahedra to optimize
- - try a mesh improvement parallelization based on Geogram point-associated locks
+ - optimize `respectEdgeConstraint` in the same way we did for `reshapeCavityIfNeeded`
+ - optimize the first tet optimization to avoid being stuck if all points are on the same plane
+ - compute Moore indices on the fly for points already in the mesh (find heuristic for when the old tactic should be employed)
+ - try a GSC parallelization based on Geogram point-associated locks
- add edge-collapse and point insertion
- modify surface mesh up to a certain distance allowed
- add a pass of smoothing aimed solely at mesh size control
- track the order of allocated blocks in reproducible Delaunay
- - simplify the check between Hxt colors and Gmsh colors
- code the boundary recovery, to replace the tetBR module
+ - extend HXT to 40+ bits nodes
## References
diff --git a/contrib/hxt/tetMesh/src/hxt_tetDelaunay.c b/contrib/hxt/tetMesh/src/hxt_tetDelaunay.c
index f5a2c687b3e8f087687d60999d320ac0eb7d9479..ebd561dceb46c829b7cfc34912be1c3072e82946 100644
--- a/contrib/hxt/tetMesh/src/hxt_tetDelaunay.c
+++ b/contrib/hxt/tetMesh/src/hxt_tetDelaunay.c
@@ -517,7 +517,6 @@ static inline HXTStatus walking2Cavity(HXTMesh* mesh, HXTPartition* partition, u
// we should pass in orient3d mode here :p
unsigned index = 4;
- unsigned outside = 0;
unsigned randomU = hxtReproducibleLCG(&seed);
for (unsigned j=0; j<4; j++)
{
@@ -535,7 +534,6 @@ static inline HXTStatus walking2Cavity(HXTMesh* mesh, HXTPartition* partition, u
// vtaCoord[0], vtaCoord[1], vtaCoord[2]);
// }
- outside = 1;
uint32_t node = mesh->tetrahedra.node[curNeigh[i]];
if(node==HXT_GHOST_VERTEX) {
@@ -543,10 +541,11 @@ static inline HXTStatus walking2Cavity(HXTMesh* mesh, HXTPartition* partition, u
return HXT_STATUS_OK;
}
- if(!vertexOutOfPartition(vertices, node, rel, partition->startDist)){
- index=i;
- break;
- }
+ if(vertexOutOfPartition(vertices, node, rel, partition->startDist))
+ return HXT_STATUS_CONFLICT;
+
+ index=i;
+ break;
}
}
}
@@ -556,13 +555,10 @@ static inline HXTStatus walking2Cavity(HXTMesh* mesh, HXTPartition* partition, u
const double* __restrict__ b = vertices[curNode[1]].coord;
const double* __restrict__ c = vertices[curNode[2]].coord;
const double* __restrict__ d = vertices[curNode[3]].coord;
- if(outside) {
- return HXT_STATUS_CONFLICT;
- }
- else if( (orient3d(a,b,c,vtaCoord)>=0) +
- (orient3d(a,b,vtaCoord,d)>=0) +
- (orient3d(a,vtaCoord,c,d)>=0) +
- (orient3d(vtaCoord,b,c,d)>=0)>2 ){
+ if( (orient3d(a,b,c,vtaCoord)>=0) +
+ (orient3d(a,b,vtaCoord,d)>=0) +
+ (orient3d(a,vtaCoord,c,d)>=0) +
+ (orient3d(vtaCoord,b,c,d)>=0)>2 ){
*curTet = nextTet;
return HXT_STATUS_DOUBLE_PT;
}
@@ -570,6 +566,12 @@ static inline HXTStatus walking2Cavity(HXTMesh* mesh, HXTPartition* partition, u
return HXT_STATUS_OK;
}
+#ifdef DEBUG
+ nstep++;
+ if(nstep > mesh->tetrahedra.num)
+ return HXT_ERROR_MSG(HXT_STATUS_ERROR, "walk stuck in an infinite loop");
+#endif
+
enteringFace = curNeigh[index]&3;
nextTet = curNeigh[index]/4;
}
@@ -928,6 +930,18 @@ static inline HXTStatus fillingACavity(HXTMesh* mesh, TetLocal* local, unsigned
mesh->tetrahedra.flag[newTet] = local->ball.array[i].flag;
+ #ifdef DEBUG
+ if((getFacetConstraint(mesh, neigh / 4, neigh % 4)!=0) ^ (getFacetConstraint(mesh, newTet, 0)!=0))
+ return HXT_ERROR_MSG(HXT_STATUS_ERROR, "constraints are not the same on both side of the new tet");
+ if(getFacetConstraint(mesh, newTet, 1) || getFacetConstraint(mesh, newTet, 2) || getFacetConstraint(mesh, newTet, 3))
+ return HXT_ERROR_MSG(HXT_STATUS_ERROR, "new tetrahedron flag are erroneous (face constraint on new face)");
+ if(getEdgeConstraint(mesh, newTet, getEdgeFromNodes(0, 1)) ||
+ getEdgeConstraint(mesh, newTet, getEdgeFromNodes(0, 2)) ||
+ getEdgeConstraint(mesh, newTet, getEdgeFromNodes(0, 3))
+ )
+ return HXT_ERROR_MSG(HXT_STATUS_ERROR, "new tetrahedron flag are erroneous (edge constraint on new edge)");
+ #endif
+
// update neighbor's neighbor
mesh->tetrahedra.neigh[neigh] = 4*newTet;
@@ -1450,8 +1464,6 @@ static HXTStatus parallelDelaunay3D(HXTMesh* mesh,
******************************************************/
for (uint32_t i=0; i<localN; i++)
{
- // if( i%1000==0)
- // printf("%f %%\n", 100.0*i/localN);
uint32_t passIndex = (localStart+i)%passLength;
uint32_t vta = nodeInfo[passStart + passIndex].node;
if(nodeInfo[passStart + passIndex].status==HXT_STATUS_TRYAGAIN){
diff --git a/contrib/hxt/tetMesh/src/hxt_tetDelaunayReshape.c b/contrib/hxt/tetMesh/src/hxt_tetDelaunayReshape.c
index 6d53b8fe997ab297aef9f55ee5c68405a608760c..d2307a3eea5b0b3ea84b2608bf7752574180de35 100644
--- a/contrib/hxt/tetMesh/src/hxt_tetDelaunayReshape.c
+++ b/contrib/hxt/tetMesh/src/hxt_tetDelaunayReshape.c
@@ -387,7 +387,6 @@ HXTStatus reshapeCavityIfNeeded(TetLocal* local, HXTMesh* mesh, const uint32_t v
uint64_t out = faces[4 * (in / 4) + f];
faces[out] = local->ball.num;
-
uint64_t* curNeigh = mesh->tetrahedra.neigh + tetToUndelete*4;
uint32_t* curNode = mesh->tetrahedra.node + tetToUndelete*4;
if(f==0) {
@@ -442,7 +441,9 @@ HXTStatus reshapeCavityIfNeeded(TetLocal* local, HXTMesh* mesh, const uint32_t v
}
local->deleted.num -= shift;
- // HXT_ASSERT( isStarShaped(local, mesh, vta, &blindFace)==HXT_STATUS_OK );
+#ifdef DEBUG
+ HXT_ASSERT( isStarShaped(local, mesh, vta, &blindFace)==HXT_STATUS_OK );
+#endif
return HXT_STATUS_OK;
}