Commit 6e0972a3 authored by Célestin Marot's avatar Célestin Marot
Browse files

chirotope -> orient3d (easier to understand)

parent 4ec7fad3
......@@ -137,12 +137,12 @@ static inline void get_compressed_index(unsigned v0,
// *index = v0 + v1*(v1-1)/2 + v2*(v2-1)*(v2-2)/6 + v3*(v3-1)*(v3-2)*(v3-3)/24;
}
#ifdef SPR_SAVE_CHIROTOPES
static SPRNOINLINE uint32_t add_chirotope(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
#ifdef SPR_SAVE_ORIENT3D
static SPRNOINLINE uint32_t add_orient3d(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
{
const unsigned config[12][4] = {
{0, 1, 2, 3},{0, 2, 3, 1},{0, 3, 1, 2},
......@@ -153,12 +153,12 @@ static SPRNOINLINE uint32_t add_chirotope(SPRCavity* SPR,
const unsigned v[4] = {v0, v1, v2, v3};
// compute the chirotope
double orientation = orient3d(&SPR->values,
SPR->points.array[v0].coord,
SPR->points.array[v1].coord,
SPR->points.array[v2].coord,
SPR->points.array[v3].coord);
// compute the orientation
double orientation = orient3d(SPR->values,
SPR->points.array[v0].coord,
SPR->points.array[v1].coord,
SPR->points.array[v2].coord,
SPR->points.array[v3].coord);
uint32_t val = 2 + (orientation>0.0) - (orientation<0.0);
uint32_t oppval = 2 + (orientation<0.0) - (orientation>0.0);
......@@ -176,24 +176,24 @@ static SPRNOINLINE uint32_t add_chirotope(SPRCavity* SPR,
ID = v0*n3+v1*n2+v2*n1+v3;
index = ID/16;
bit = ID%16 * 2;
SPR->map.chirotopes[index] &= ~(3U<<bit); // clear bits
SPR->map.chirotopes[index] |= val<<bit; // set bits
SPR->map.orient3d[index] &= ~(3U<<bit); // clear bits
SPR->map.orient3d[index] |= val<<bit; // set bits
ID = v0*n3+v1*n2+v3*n1+v2;
index = ID/16;
bit = ID%16 * 2;
SPR->map.chirotopes[index] &= ~(3U<<bit); // clear bits
SPR->map.chirotopes[index] |= oppval<<bit; // set bits
SPR->map.orient3d[index] &= ~(3U<<bit); // clear bits
SPR->map.orient3d[index] |= oppval<<bit; // set bits
}
return val;
}
static int get_chirotope(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
static int get_orient3d(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
{
if(v0==v1 || v0==v2 || v0==v3 || v1==v2 || v1==v3 || v2==v3)
return 0;
......@@ -201,9 +201,9 @@ static int get_chirotope(SPRCavity* SPR,
unsigned ID = v0*(n*n*n)+v1*(n*n)+v2*n+v3;
unsigned index = ID/16;
unsigned bit = ID%16 * 2;
uint32_t val = SPR->map.chirotopes[index]>>bit & 3;
uint32_t val = SPR->map.orient3d[index]>>bit & 3;
if(val==0) {
return (int)add_chirotope(SPR, v0, v1, v2 ,v3) - 2;
return (int)add_orient3d(SPR, v0, v1, v2 ,v3) - 2;
}
else {
return (int)val - 2;
......@@ -211,17 +211,17 @@ static int get_chirotope(SPRCavity* SPR,
}
#else
static int get_chirotope(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
static int get_orient3d(SPRCavity* SPR,
unsigned v0,
unsigned v1,
unsigned v2,
unsigned v3)
{
if(v0==v1 || v0==v2 || v0==v3 || v1==v2 || v1==v3 || v2==v3)
return 0;
return orient3d_inline(&SPR->values,
return orient3d_inline(SPR->values,
SPR->points.array[v0].coord,
SPR->points.array[v1].coord,
SPR->points.array[v2].coord,
......@@ -230,7 +230,7 @@ static int get_chirotope(SPRCavity* SPR,
#endif
// triangle-edge intersection following chirotopes
// triangle-edge intersection following orient3ds
// !! the case where v0 or v1 lie on the triangle
// must be checked independently !!
static inline int tri_edge_intersection(int t012v0,
......@@ -255,7 +255,7 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
uint8_t* tri = SPR->triangles.array[triangleID].node;
#ifndef NDEBUG
int v = get_chirotope(SPR, tet[0], tet[1], tet[2], tet[3]);
int v = get_orient3d(SPR, tet[0], tet[1], tet[2], tet[3]);
if(v>=0) {
fprintf(stderr, "Fatal Error: degenerated tetrahedron created %d\n",v);
exit(EXIT_FAILURE);
......@@ -325,9 +325,9 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
if(tri0_to_tet!=-1 && tri1_to_tet!=-1 && tri2_to_tet!=-1) // the triangle is there
return 0;
int f0v0 = get_chirotope(SPR, tri[0], tet[1], tet[2], tet[3]);
int f0v1 = get_chirotope(SPR, tri[1], tet[1], tet[2], tet[3]);
int f0v2 = get_chirotope(SPR, tri[2], tet[1], tet[2], tet[3]);
int f0v0 = get_orient3d(SPR, tri[0], tet[1], tet[2], tet[3]);
int f0v1 = get_orient3d(SPR, tri[1], tet[1], tet[2], tet[3]);
int f0v2 = get_orient3d(SPR, tri[2], tet[1], tet[2], tet[3]);
// if the triangle is entirely on one side of the tet: no intersection
if((f0v0 > 0 || (f0v0==0 && tri0_to_tet!=-1)) &&
......@@ -335,9 +335,9 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
(f0v2 > 0 || (f0v2==0 && tri2_to_tet!=-1)) )
return 0;
int f1v0 = get_chirotope(SPR, tet[0], tri[0], tet[2], tet[3]);
int f1v1 = get_chirotope(SPR, tet[0], tri[1], tet[2], tet[3]);
int f1v2 = get_chirotope(SPR, tet[0], tri[2], tet[2], tet[3]);
int f1v0 = get_orient3d(SPR, tet[0], tri[0], tet[2], tet[3]);
int f1v1 = get_orient3d(SPR, tet[0], tri[1], tet[2], tet[3]);
int f1v2 = get_orient3d(SPR, tet[0], tri[2], tet[2], tet[3]);
// we repeat that with the 3 other faces
if((f1v0 > 0 || (f1v0==0 && tri0_to_tet!=-1)) &&
......@@ -345,18 +345,18 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
(f1v2 > 0 || (f1v2==0 && tri2_to_tet!=-1)) )
return 0;
int f2v0 = get_chirotope(SPR, tet[0], tet[1], tri[0], tet[3]);
int f2v1 = get_chirotope(SPR, tet[0], tet[1], tri[1], tet[3]);
int f2v2 = get_chirotope(SPR, tet[0], tet[1], tri[2], tet[3]);
int f2v0 = get_orient3d(SPR, tet[0], tet[1], tri[0], tet[3]);
int f2v1 = get_orient3d(SPR, tet[0], tet[1], tri[1], tet[3]);
int f2v2 = get_orient3d(SPR, tet[0], tet[1], tri[2], tet[3]);
if((f2v0 > 0 || (f2v0==0 && tri0_to_tet!=-1)) &&
(f2v1 > 0 || (f2v1==0 && tri1_to_tet!=-1)) &&
(f2v2 > 0 || (f2v2==0 && tri2_to_tet!=-1)) )
return 0;
int f3v0 = get_chirotope(SPR, tet[0], tet[1], tet[2], tri[0]);
int f3v1 = get_chirotope(SPR, tet[0], tet[1], tet[2], tri[1]);
int f3v2 = get_chirotope(SPR, tet[0], tet[1], tet[2], tri[2]);
int f3v0 = get_orient3d(SPR, tet[0], tet[1], tet[2], tri[0]);
int f3v1 = get_orient3d(SPR, tet[0], tet[1], tet[2], tri[1]);
int f3v2 = get_orient3d(SPR, tet[0], tet[1], tet[2], tri[2]);
if((f3v0 > 0 || (f3v0==0 && tri0_to_tet!=-1)) &&
(f3v1 > 0 || (f3v1==0 && tri1_to_tet!=-1)) &&
......@@ -371,23 +371,23 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
int side = 0;
int outside = 1;
int v0f = get_chirotope(SPR, tet[0], tri[0], tri[1], tri[2]);
int v0f = get_orient3d(SPR, tet[0], tri[0], tri[1], tri[2]);
side = v0f;
int v1f = get_chirotope(SPR, tet[1], tri[0], tri[1], tri[2]);
int v1f = get_orient3d(SPR, tet[1], tri[0], tri[1], tri[2]);
if(v1f!=0) {
if(side==0)
side = v1f;
else if(v1f*side<0) // not on the same side
outside = 0;
}
int v2f = get_chirotope(SPR, tet[2], tri[0], tri[1], tri[2]);
int v2f = get_orient3d(SPR, tet[2], tri[0], tri[1], tri[2]);
if(v2f!=0) {
if(side==0)
side = v2f;
else if(v2f*side<0) // not on the same side
outside = 0;
}
int v3f = get_chirotope(SPR, tet[3], tri[0], tri[1], tri[2]);
int v3f = get_orient3d(SPR, tet[3], tri[0], tri[1], tri[2]);
if(v3f!=0) {
if(side==0)
side = v3f;
......@@ -404,12 +404,12 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
int e0e1, e1e1, e2e1, e3e1, e4e1, e5e1;
int e0e2, e1e2, e2e2, e3e2, e4e2, e5e2;
e0e0 = get_chirotope(SPR, tet[0], tet[1], tri[0], tri[1]);
e1e0 = get_chirotope(SPR, tet[0], tet[2], tri[0], tri[1]);
e2e0 = get_chirotope(SPR, tet[0], tet[3], tri[0], tri[1]);
e3e0 = get_chirotope(SPR, tet[1], tet[2], tri[0], tri[1]);
e4e0 = get_chirotope(SPR, tet[1], tet[3], tri[0], tri[1]);
e5e0 = get_chirotope(SPR, tet[2], tet[3], tri[0], tri[1]);
e0e0 = get_orient3d(SPR, tet[0], tet[1], tri[0], tri[1]);
e1e0 = get_orient3d(SPR, tet[0], tet[2], tri[0], tri[1]);
e2e0 = get_orient3d(SPR, tet[0], tet[3], tri[0], tri[1]);
e3e0 = get_orient3d(SPR, tet[1], tet[2], tri[0], tri[1]);
e4e0 = get_orient3d(SPR, tet[1], tet[3], tri[0], tri[1]);
e5e0 = get_orient3d(SPR, tet[2], tet[3], tri[0], tri[1]);
// if a point of the edge is coplanar with the facet of
// the tet, we don't really care trying to check the intersecton
......@@ -425,12 +425,12 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
return 1;
e0e1 = get_chirotope(SPR, tet[0], tet[1], tri[1], tri[2]);
e1e1 = get_chirotope(SPR, tet[0], tet[2], tri[1], tri[2]);
e2e1 = get_chirotope(SPR, tet[0], tet[3], tri[1], tri[2]);
e3e1 = get_chirotope(SPR, tet[1], tet[2], tri[1], tri[2]);
e4e1 = get_chirotope(SPR, tet[1], tet[3], tri[1], tri[2]);
e5e1 = get_chirotope(SPR, tet[2], tet[3], tri[1], tri[2]);
e0e1 = get_orient3d(SPR, tet[0], tet[1], tri[1], tri[2]);
e1e1 = get_orient3d(SPR, tet[0], tet[2], tri[1], tri[2]);
e2e1 = get_orient3d(SPR, tet[0], tet[3], tri[1], tri[2]);
e3e1 = get_orient3d(SPR, tet[1], tet[2], tri[1], tri[2]);
e4e1 = get_orient3d(SPR, tet[1], tet[3], tri[1], tri[2]);
e5e1 = get_orient3d(SPR, tet[2], tet[3], tri[1], tri[2]);
if( tri_edge_intersection( f0v1, f0v2, -e3e1, -e5e1, e4e1)
|| tri_edge_intersection( f1v1, f1v2, -e2e1, e5e1, e1e1)
......@@ -438,12 +438,12 @@ static inline int tet_tri_intersection(SPRCavity* SPR,
|| tri_edge_intersection( f3v1, f3v2, -e1e1, e3e1, e0e1))
return 1;
e0e2 = get_chirotope(SPR, tet[0], tet[1], tri[2], tri[0]);
e1e2 = get_chirotope(SPR, tet[0], tet[2], tri[2], tri[0]);
e2e2 = get_chirotope(SPR, tet[0], tet[3], tri[2], tri[0]);
e3e2 = get_chirotope(SPR, tet[1], tet[2], tri[2], tri[0]);
e4e2 = get_chirotope(SPR, tet[1], tet[3], tri[2], tri[0]);
e5e2 = get_chirotope(SPR, tet[2], tet[3], tri[2], tri[0]);
e0e2 = get_orient3d(SPR, tet[0], tet[1], tri[2], tri[0]);
e1e2 = get_orient3d(SPR, tet[0], tet[2], tri[2], tri[0]);
e2e2 = get_orient3d(SPR, tet[0], tet[3], tri[2], tri[0]);
e3e2 = get_orient3d(SPR, tet[1], tet[2], tri[2], tri[0]);
e4e2 = get_orient3d(SPR, tet[1], tet[3], tri[2], tri[0]);
e5e2 = get_orient3d(SPR, tet[2], tet[3], tri[2], tri[0]);
if( tri_edge_intersection( f0v2, f0v0, -e3e2, -e5e2, e4e2)
|| tri_edge_intersection( f1v2, f1v0, -e2e2, e5e2, e1e2)
......@@ -515,30 +515,30 @@ static SPRNOINLINE int tet_edge_intersection(SPRCavity* SPR,
if(l0_to_tet!=-1 && l1_to_tet!=-1) // the line is there: no interseection
return 0;
int f0v0 = get_chirotope(SPR, edge[0], tet[1], tet[2], tet[3]);
int f0v1 = get_chirotope(SPR, edge[1], tet[1], tet[2], tet[3]);
int f0v0 = get_orient3d(SPR, edge[0], tet[1], tet[2], tet[3]);
int f0v1 = get_orient3d(SPR, edge[1], tet[1], tet[2], tet[3]);
// if the triangle is entirely on one side of the tet: no intersection
if((f0v0 > 0 || (f0v0==0 && l0_to_tet!=-1)) &&
(f0v1 > 0 || (f0v1==0 && l1_to_tet!=-1)) )
return 0;
int f1v0 = get_chirotope(SPR, tet[0], edge[0], tet[2], tet[3]);
int f1v1 = get_chirotope(SPR, tet[0], edge[1], tet[2], tet[3]);
int f1v0 = get_orient3d(SPR, tet[0], edge[0], tet[2], tet[3]);
int f1v1 = get_orient3d(SPR, tet[0], edge[1], tet[2], tet[3]);
if((f1v0 > 0 || (f1v0==0 && l0_to_tet!=-1)) &&
(f1v1 > 0 || (f1v1==0 && l1_to_tet!=-1)) )
return 0;
int f2v0 = get_chirotope(SPR, tet[0], tet[1], edge[0], tet[3]);
int f2v1 = get_chirotope(SPR, tet[0], tet[1], edge[1], tet[3]);
int f2v0 = get_orient3d(SPR, tet[0], tet[1], edge[0], tet[3]);
int f2v1 = get_orient3d(SPR, tet[0], tet[1], edge[1], tet[3]);
if((f2v0 > 0 || (f2v0==0 && l0_to_tet!=-1)) &&
(f2v1 > 0 || (f2v1==0 && l1_to_tet!=-1)) )
return 0;
int f3v0 = get_chirotope(SPR, tet[0], tet[1], tet[2], edge[0]);
int f3v1 = get_chirotope(SPR, tet[0], tet[1], tet[2], edge[1]);
int f3v0 = get_orient3d(SPR, tet[0], tet[1], tet[2], edge[0]);
int f3v1 = get_orient3d(SPR, tet[0], tet[1], tet[2], edge[1]);
if((f3v0 > 0 || (f3v0==0 && l0_to_tet!=-1)) &&
(f3v1 > 0 || (f3v1==0 && l1_to_tet!=-1)) )
......@@ -558,12 +558,12 @@ static SPRNOINLINE int tet_edge_intersection(SPRCavity* SPR,
return 0;
// the line with each edge of the tetrahedra
int c0 = get_chirotope(SPR, tet[0], tet[1], edge[0], edge[1]);
int c1 = get_chirotope(SPR, tet[0], tet[2], edge[0], edge[1]);
int c2 = get_chirotope(SPR, tet[0], tet[3], edge[0], edge[1]);
int c3 = get_chirotope(SPR, tet[1], tet[2], edge[0], edge[1]);
int c4 = get_chirotope(SPR, tet[1], tet[3], edge[0], edge[1]);
int c5 = get_chirotope(SPR, tet[2], tet[3], edge[0], edge[1]);
int c0 = get_orient3d(SPR, tet[0], tet[1], edge[0], edge[1]);
int c1 = get_orient3d(SPR, tet[0], tet[2], edge[0], edge[1]);
int c2 = get_orient3d(SPR, tet[0], tet[3], edge[0], edge[1]);
int c3 = get_orient3d(SPR, tet[1], tet[2], edge[0], edge[1]);
int c4 = get_orient3d(SPR, tet[1], tet[3], edge[0], edge[1]);
int c5 = get_orient3d(SPR, tet[2], tet[3], edge[0], edge[1]);
// if one point of the edge is coplanar with the facet of
// the tet, we don't really care trying to check the intersecton
......@@ -594,10 +594,10 @@ static int tet_contains_pt(SPRCavity* SPR,
pt_coord[0] > tetBbox->max[0] ||
pt_coord[1] > tetBbox->max[1] ||
pt_coord[2] > tetBbox->max[2] ||
get_chirotope(SPR, tet[1], tet[3], tet[2], pt)>0 ||
get_chirotope(SPR, tet[2], tet[3], tet[0], pt)>0 ||
get_chirotope(SPR, tet[3], tet[1], tet[0], pt)>0 ||
get_chirotope(SPR, tet[0], tet[1], tet[2], pt)>0)
get_orient3d(SPR, tet[1], tet[3], tet[2], pt)>0 ||
get_orient3d(SPR, tet[2], tet[3], tet[0], pt)>0 ||
get_orient3d(SPR, tet[3], tet[1], tet[0], pt)>0 ||
get_orient3d(SPR, tet[0], tet[1], tet[2], pt)>0)
return 0;
return 1;
......@@ -665,7 +665,7 @@ static double get_quality_map(SPRCavity* SPR,
unsigned v2,
unsigned v3)
{
int orientation = get_chirotope(SPR, v0, v1, v2, v3);
int orientation = get_orient3d(SPR, v0, v1, v2, v3);
if(orientation > 0)
return -1.0;
else if(orientation==0)
......@@ -697,7 +697,7 @@ static double get_quality_map(SPRCavity* SPR,
unsigned v2,
unsigned v3)
{
int orientation = get_chirotope(SPR, v0, v1, v2, v3);
int orientation = get_orient3d(SPR, v0, v1, v2, v3);
if(orientation > 0)
return -1.0;
else if(orientation==0)
......@@ -1040,11 +1040,11 @@ static inline int SPR_init(SPRCavity* SPR)
}
#endif
#ifdef SPR_SAVE_CHIROTOPES
{ // allocating & setting to zero map.chirotopes
#ifdef SPR_SAVE_ORIENT3D
{ // allocating & setting to zero map.orient3d
size_t chir_len = n*n*n*n/16+1;
SPR->map.chirotopes = calloc(chir_len, sizeof(uint32_t));
if(SPR->map.chirotopes==NULL)
SPR->map.orient3d = calloc(chir_len, sizeof(uint32_t));
if(SPR->map.orient3d==NULL)
goto allocation_error;
}
#endif
......@@ -1058,7 +1058,7 @@ static inline int SPR_init(SPRCavity* SPR)
SPRBboxAddOne(&bbox, p->coord);
}
exactinit(&SPR->values,
exactinit(SPR->values,
bbox.max[0] - bbox.min[0],
bbox.max[1] - bbox.min[1],
bbox.max[2] - bbox.min[2]);
......@@ -1097,8 +1097,8 @@ static inline void SPR_terminate(SPRCavity* SPR)
#ifdef SPR_SAVE_QUALITIES
free(SPR->map.qualities);
#endif
#ifdef SPR_SAVE_CHIROTOPES
free(SPR->map.chirotopess);
#ifdef SPR_SAVE_ORIENT3D
free(SPR->map.orient3d);
#endif
}
......@@ -1127,6 +1127,9 @@ static inline void SPR_terminate(SPRCavity* SPR)
int hxtSPR(SPRCavity* SPR)
{
SPRExact valuesOnStack;
SPR->values = &valuesOnStack;
if(SPR_init(SPR))
return 2;
......@@ -1149,6 +1152,8 @@ int hxtSPR(SPRCavity* SPR)
rewind = 1;
}
else if(step->numCandidates<0){
// choosing the starting triangle
double qualMax;
const uint16_t triangleID = best_face_heuristic(SPR, &qualMax);
......@@ -1159,7 +1164,7 @@ int hxtSPR(SPRCavity* SPR)
else {
step->tet = SPR->triangles.array[triangleID];
// create the candidate array
// create (and order) the candidate array
compute_candidates(SPR, step, qualMax, triangleID);
// remove the face from the map
......@@ -1176,6 +1181,7 @@ int hxtSPR(SPRCavity* SPR)
step->tet.node[2],
step->tet.node[3]);
// best quality could have changed if there was a rewind
if(qual <= SPR->tetrahedra.quality)
continue;
......@@ -1186,7 +1192,7 @@ int hxtSPR(SPRCavity* SPR)
uint16_t index = SPR_get_faceMap(SPR, p0, p1, p2);
if(index==UINT16_MAX) {
if(index == UINT16_MAX) {
// add the face to the map
add_face(SPR, p0, p2, p1);
}
......@@ -1223,7 +1229,7 @@ int hxtSPR(SPRCavity* SPR)
uint16_t index = SPR_get_faceMap(SPR, p0, p2, p1);
if(index==UINT16_MAX) {
if(index == UINT16_MAX) {
// add the face to the map
add_face(SPR, p0, p1, p2);
}
......@@ -1233,8 +1239,7 @@ int hxtSPR(SPRCavity* SPR)
}
}
if(step->quality <= SPR->tetrahedra.quality ||
SPR->num_search_nodes >= SPR->max_search_nodes) {
if(step->quality <= SPR->tetrahedra.quality) {
step->numCandidates = 0;
}
}
......
......@@ -28,13 +28,12 @@ extern "C" {
#include <stdlib.h>
#include <stdio.h>
#include <stdint.h>
#include "predicates.h"
// defines to set/unset the main optimizations
#define SPR_SORT_CANDIDATES
#define SPR_SAVE_QUALITIES
#define SPR_USE_FACE_HEURISTIC
// #define SPR_SAVE_CHIROTOPES
// #define SPR_SAVE_ORIENT3D
// #define SPR_INPUT_VERIFICATION // not implemented yet
#define SPR_MAX_PTS 48 // you can change only this value
......@@ -92,6 +91,8 @@ typedef struct {
double quality; // quality before inserting this tet
} SPRStep;
typedef struct SPRExact_s SPRExact; // an opaque structure...
/* This is the main structure defining your cavity.
* You might want to allocate this heavy structure on the heap ! */
......@@ -169,13 +170,13 @@ typedef struct {
#ifdef SPR_SAVE_QUALITIES
double* qualities;
#endif
#ifdef SPR_SAVE_CHIROTOPES
uint32_t* chirotopes;
#ifdef SPR_SAVE_ORIENT3D
uint32_t* orient3d;
#endif
uint16_t* faces;
} map;
SPRExact values;
SPRExact* values;
} SPRCavity;
......
CFLAGS=-Wall -O3 -g -DDEBUG -Wextra -fno-omit-frame-pointer -Wno-missing-field-initializers -Wpointer-arith -Wcast-align -Wcast-qual -Wstrict-aliasing -Wpointer-arith -Winit-self -Wredundant-decls -fsanitize=address -fsanitize=leak -fsanitize=undefined
LDFLAGS=-lm -fsanitize=address -fsanitize=leak -fsanitize=undefined
CFLAGS=-Wall -Wextra -O3 #add '-fPIC' for shared library
LDFLAGS=-lm #add '-shared' for shared library
TARGET=./test #replace by 'libHXTSPR.so' for shared library
EXEC=./test
SRC= $(wildcard *.c)
OBJ= $(SRC:.c=.o)
SRC=$(wildcard *.c) #replace by '$(filter-out src/test.c,$(wildcard *.c))' for shared library
OBJ=$(SRC:.c=.o)
all: $(EXEC)
$(EXEC): $(OBJ)
$(TARGET): $(OBJ)
$(CC) -o $@ $^ $(LDFLAGS)
predicates.o: predicates.c predicates.h
......@@ -16,4 +14,9 @@ test.o: test.c HXTSPR.h predicates.h
$(OBJ):
# @echo $@ depends on $^
$(CC) -o $@ -c $< $(CFLAGS) -lm
$(CC) $(CFLAGS) -o $@ -c $< $(LDFLAGS)
help:
less README.md
.PHONY: all help
......@@ -222,7 +222,7 @@
/* */
/*****************************************************************************/
void exactinit(SPRExact* values, REAL maxx, REAL maxy, REAL maxz)
void exactinit(struct SPRExact_s* values, REAL maxx, REAL maxy, REAL maxz)
{
#ifdef CPU86
#ifdef SINGLE
......@@ -435,7 +435,7 @@ int fast_expansion_sum_zeroelim(const int elen,
/* */
/*****************************************************************************/
static int scale_expansion_zeroelim(const SPRExact* const values,
static int scale_expansion_zeroelim(const struct SPRExact_s* const values,
const int elen,
const REAL* const __restrict__ e,
const REAL b,
......@@ -530,7 +530,7 @@ static REAL estimate(const int elen, const REAL* const e)
/*****************************************************************************/
static REAL orient3dadapt(const SPRExact* const values,
static REAL orient3dadapt(const struct SPRExact_s* const values,
const REAL* const __restrict__ pa,
const REAL* const __restrict__ pb,
const REAL* const __restrict__ pc,
......@@ -947,7 +947,7 @@ static REAL orient3dadapt(const SPRExact* const values,
}
REAL orient3d(const SPRExact* const values,
REAL orient3d(const struct SPRExact_s* const values,
const REAL* const __restrict__ pa,
const REAL* const __restrict__ pb,
const REAL* const __restrict__ pc,
......
......@@ -30,7 +30,7 @@ extern "C" {
/* we used a structure instead of the global variables used originally */
typedef struct {
struct SPRExact_s{
/* splitter = 2^ceiling(p / 2) + 1. Used to split floats in half. */
REAL splitter;
/* A set of coefficients used to calculate maximum roundoff errors. */
......@@ -40,15 +40,15 @@ typedef struct {
// Static filters for orient3d() and insphere().
// Added by H. Si, 2012-08-23.
REAL o3dstaticfilter, o3derrboundA;
} SPRExact;
};
void exactinit(SPRExact* values, REAL maxx, REAL maxy, REAL maxz);
void exactinit(struct SPRExact_s* values, REAL maxx, REAL maxy, REAL maxz);
REAL orient3d(const SPRExact* const values,
REAL orient3d(const struct SPRExact_s* const values,
const REAL* const __restrict__ pa,
const REAL* const __restrict__ pb,
const REAL* const __restrict__ pc,
......@@ -56,7 +56,7 @@ REAL orient3d(const SPRExact* const values,
// remember to call exactinit with the right ranges beforehand !!!
static inline int orient3d_inline(const SPRExact* const values,
static inline int orient3d_inline(const struct SPRExact_s* const values,
const REAL* const __restrict__ pa,
const REAL* const __restrict__ pb,
const REAL* const __restrict__ pc,
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment