Commit 4ec7fad3 authored by Célestin Marot's avatar Célestin Marot
Browse files

both iterative and recursive version works

+ rename files to correspond to HXTSort naming scheme
parent 5e9d8385
......@@ -17,7 +17,7 @@
*
Author: Célestin Marot (celestin.marot@uclouvain.be) */
#include "hxt_SPR.h"
#include "HXTSPR.h"
#include "predicates.h"
#include <float.h>
#include <math.h>
......@@ -1123,136 +1123,130 @@ static inline void SPR_terminate(SPRCavity* SPR)
|*| | | \ \ / -_) _/ _ \ ' \| ' \/ -_) _| _| / _ \ ' \ |*|
|*| |_| \_\\___\__\___/_||_|_||_\___\__|\__|_\___/_||_| |*|
\*|=======================================================|*/
// int hxtSPR(SPRCavity* SPR) {
// if(SPR_init(SPR))
// return 2;
#ifndef SPR_RECURSIVE_VERSION
// SPR->steps.num = 0;
// SPRStep* step = &SPR->steps.array[0];
// step->numCandidates = -1; // candidates not computed yet
// step->quality = DBL_MAX; // we can hope for the biggest quality at the moment
int hxtSPR(SPRCavity* SPR)
{
if(SPR_init(SPR))
return 2;
SPR->steps.num = 0;
SPRStep* step = &SPR->steps.array[0];
step->numCandidates = -1; // candidates not computed yet
step->quality = DBL_MAX; // we can hope for the biggest quality at the moment
while(1) {
if(SPR->num_search_nodes >= SPR->max_search_nodes){
SPR_terminate(SPR);
return 1;
}
SPR->num_search_nodes++;
int rewind = 0;
if(SPR->triangles.num==0) { // we found a new (better) tetrahedralization
new_tetrahedralization(SPR, step);
rewind = 1;
}
else if(step->numCandidates<0){
double qualMax;
const uint16_t triangleID = best_face_heuristic(SPR, &qualMax);
// while(1){
// if(SPR->num_search_nodes >= SPR->max_search_nodes){
// SPR_terminate(SPR);
// return 1;
// }
// SPR->num_search_nodes++;
if(triangleID == UINT16_MAX) {
step->numCandidates = -1;
rewind = 1;
}
else {
step->tet = SPR->triangles.array[triangleID];
// create the candidate array
compute_candidates(SPR, step, qualMax, triangleID);
// assert(step==&SPR->steps.array[SPR->steps.num]);
// remove the face from the map
remove_face(SPR, triangleID);
}
}
if(step->numCandidates>0) {
step->tet.node[3] = step->candidate[--step->numCandidates];
double qual = get_quality_map(SPR,
step->tet.node[0],
step->tet.node[1],
step->tet.node[2],
step->tet.node[3]);
// printf("step %d, %d candidates, quality %f\n", SPR->steps.num, step->numCandidates, step->quality);
if(qual <= SPR->tetrahedra.quality)
continue;
for (int f=0; f<3; f++) {
uint8_t p0 = step->tet.node[nodeFromFacet[0][f]];
uint8_t p1 = step->tet.node[nodeFromFacet[1][f]];
uint8_t p2 = step->tet.node[nodeFromFacet[2][f]];
uint16_t index = SPR_get_faceMap(SPR, p0, p1, p2);
if(index==UINT16_MAX) {
// add the face to the map
add_face(SPR, p0, p2, p1);
}
else {
// remove the face from the map
remove_face(SPR, index);
}
}
double oldQual = step->quality;
step = &SPR->steps.array[++SPR->steps.num];
step->quality = qual < oldQual ? qual : oldQual;
step->numCandidates = -1;
}
else if(step->numCandidates==0){
add_face(SPR,
step->tet.node[0],
step->tet.node[1],
step->tet.node[2]);
rewind = 1;
}
if(rewind) {
if(SPR->steps.num==0) {
SPR_terminate(SPR);
return 0;
}
step = &SPR->steps.array[--SPR->steps.num];
for (int f=0; f<3; f++) {
uint8_t p0 = step->tet.node[nodeFromFacet[0][f]];
uint8_t p1 = step->tet.node[nodeFromFacet[1][f]];
uint8_t p2 = step->tet.node[nodeFromFacet[2][f]];
uint16_t index = SPR_get_faceMap(SPR, p0, p2, p1);
if(index==UINT16_MAX) {
// add the face to the map
add_face(SPR, p0, p1, p2);
}
else {
// remove the face from the map
remove_face(SPR, index);
}
}
// int rewind = 0;
// if(SPR->triangles.num==0) { // we found a new (better) tetrahedralization
// printf("tetrahedralization found ! quality %f\n", step->quality);
// new_tetrahedralization(SPR, step);
// rewind = 1;
// }
// else if(step->quality <= SPR->tetrahedra.quality) {
// rewind = 1;
// }
// else if(step->numCandidates<0) {
// assert(step->quality > SPR->tetrahedra.quality);
// // we still have to compute the face and candidates
// double qualMax;
// const uint16_t triangleID = best_face_heuristic(SPR, &qualMax);
// if(triangleID != UINT16_MAX) {
// step->tet = SPR->triangles.array[triangleID];
// // create the candidate array
// compute_candidates(SPR, step, qualMax, triangleID);
// // there still can be 0 candidates because intersection are only computed here above
// }
// if(step->numCandidates>0) {
// // remove the face from the map
// remove_face(SPR, triangleID);
// }
// else {
// rewind = 1;
// }
// printf("found %d candidates\n", step->numCandidates);
// }
// else if(step->numCandidates==0) {
// add_face(SPR,
// step->tet.node[0],
// step->tet.node[1],
// step->tet.node[2]);
// rewind = 1;
// }
// if(rewind) {
// printf("rewinding\n");
// if(SPR->steps.num==0) {
// SPR_terminate(SPR);
// return 0;
// }
// step = &SPR->steps.array[--SPR->steps.num];
// for (int f=0; f<3; f++) {
// uint8_t p0 = step->tet.node[nodeFromFacet[0][f]];
// uint8_t p1 = step->tet.node[nodeFromFacet[1][f]];
// uint8_t p2 = step->tet.node[nodeFromFacet[2][f]];
// uint16_t index = SPR_get_faceMap(SPR, p0, p2, p1);
// if(index==UINT16_MAX) {
// // add the face to the map
// add_face(SPR, p0, p1, p2);
// }
// else {
// // remove the face from the map
// remove_face(SPR, index);
// }
// }
// continue;
// }
// assert(step->numCandidates>0);
// step->tet.node[3] = step->candidate[--step->numCandidates];
// double qual = get_quality_map(SPR,
// step->tet.node[0],
// step->tet.node[1],
// step->tet.node[2],
// step->tet.node[3]);
// // best.quality could have changed...
// if(qual>SPR->tetrahedra.quality) {
// for (int f=0; f<3; f++) {
// uint8_t p0 = step->tet.node[nodeFromFacet[0][f]];
// uint8_t p1 = step->tet.node[nodeFromFacet[1][f]];
// uint8_t p2 = step->tet.node[nodeFromFacet[2][f]];
// uint16_t index = SPR_get_faceMap(SPR, p0, p1, p2);
// if(index==UINT16_MAX) {
// // add the face to the map
// add_face(SPR, p0, p2, p1);
// }
// else {
// // remove the face from the map
// remove_face(SPR, index);
// }
// }
// SPRStep* nextStep = &SPR->steps.array[++SPR->steps.num];
// nextStep->quality = qual < step->quality ? qual : step->quality;
// nextStep->numCandidates = -1;
// step = nextStep;
// }
// }
// }
if(step->quality <= SPR->tetrahedra.quality ||
SPR->num_search_nodes >= SPR->max_search_nodes) {
step->numCandidates = 0;
}
}
}
}
#else
/*****************************************************************
* Here below is a recursive version. It works and is more easy *
* to understand (IMO) than the iterative one *
*****************************************************************/
void hxtSPR_recur(SPRCavity* SPR) {
if(SPR->num_search_nodes >= SPR->max_search_nodes){
return;
......@@ -1366,4 +1360,6 @@ int hxtSPR(SPRCavity* SPR) {
if(SPR->num_search_nodes>=SPR->max_search_nodes)
return 1;
return 0;
}
\ No newline at end of file
}
#endif
\ No newline at end of file
......@@ -18,8 +18,8 @@
Author: Célestin Marot (celestin.marot@uclouvain.be) */
#ifndef _HXT_SPR_
#define _HXT_SPR_
#ifndef _HXTSPR_
#define _HXTSPR_
#ifdef __cplusplus
extern "C" {
......
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
EXEC=./exec
EXEC=./test
SRC= $(wildcard *.c)
OBJ= $(SRC:.c=.o)
......@@ -11,8 +11,8 @@ $(EXEC): $(OBJ)
$(CC) -o $@ $^ $(LDFLAGS)
predicates.o: predicates.c predicates.h
hxt_SPR.o: hxt_SPR.c hxt_SPR.h predicates.h
test.o: test.c hxt_SPR.h predicates.h
HXTSPR.o: HXTSPR.c HXTSPR.h predicates.h
test.o: test.c HXTSPR.h predicates.h
$(OBJ):
# @echo $@ depends on $^
......
......@@ -17,7 +17,7 @@
*
Author: Célestin Marot (celestin.marot@uclouvain.be) */
#include "hxt_SPR.h"
#include "HXTSPR.h"
#include <math.h>
......
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