diff --git a/CMakeLists.txt b/CMakeLists.txt index 37932898c241f193655824e2a623d38a8c721801..8e30e1beb0048eb9cd7648d7af9e12e2cb0975e4 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -33,6 +33,7 @@ find_package(Boost COMPONENTS regex thread system REQUIRED) find_package(ZLIB REQUIRED) find_package(LibXml2 REQUIRED) find_package(GLEW REQUIRED) +find_package(SuiteSparse REQUIRED) IF (DEFINED ASSERTON) add_definitions(-DCGOGN_ASSERT_BOOL=${ASSERTON}) @@ -77,6 +78,7 @@ SET (CGoGN_EXT_INCLUDES # define libs for external libs SET (CGoGN_EXT_LIBS + nl ${OPENGL_LIBRARY} ${GLEW_LIBRARIES} ${ZLIB_LIBRARIES} @@ -84,6 +86,7 @@ SET (CGoGN_EXT_LIBS ${Boost_SYSTEM_LIBRARY} ${Boost_REGEX_LIBRARY} ${Boost_THREAD_LIBRARY} + ${SUITESPARSE_LIBRARIES} ) #optionnal libs diff --git a/apps_cmake.txt b/apps_cmake.txt index 72c6d1c9db0e2f3bd3e41f6925ad5e33adfdffa3..cf8c4d34599478cad63405eb0b35087fa5a1cfe1 100644 --- a/apps_cmake.txt +++ b/apps_cmake.txt @@ -32,6 +32,7 @@ find_package(Boost COMPONENTS regex thread system REQUIRED) find_package(ZLIB REQUIRED) find_package(LibXml2 REQUIRED) find_package(GLEW REQUIRED) +find_package(SuiteSparse REQUIRED) IF (DEFINED ASSERTON) add_definitions(-DCGOGN_ASSERT_BOOL=${ASSERTON}) @@ -53,6 +54,7 @@ SET(COMMON_INCLUDES # define libs for external libs SET (COMMON_LIBS + nl ${OPENGL_LIBRARY} ${GLEW_LIBRARIES} ${ZLIB_LIBRARIES} @@ -60,6 +62,7 @@ SET (COMMON_LIBS ${Boost_SYSTEM_LIBRARY} ${Boost_REGEX_LIBRARY} ${Boost_THREAD_LIBRARY} + ${SUITESPARSE_LIBRARIES} ) #optionnal libs diff --git a/include/Algo/Modelisation/subdivision3.hpp b/include/Algo/Modelisation/subdivision3.hpp index 148f891d68330db29a1b14ceeab58a3696d0acdc..dc31f615d5a6c765311d51a82aa6ab510dd4dbb1 100644 --- a/include/Algo/Modelisation/subdivision3.hpp +++ b/include/Algo/Modelisation/subdivision3.hpp @@ -672,7 +672,7 @@ void sqrt3Vol(typename PFP::MAP& map, VertexAttribute& posit Traversor3WF tWF(map, dit); for(Dart ditWF = tWF.begin() ; ditWF != tWF.end() ; ditWF = tWF.next()) { - if(!map.isBoundaryFace(ditWF)) + if(!map.isBoundaryFace(ditWF) && !m.isMarked(ditWF)) m.markOrbit(ditWF); } @@ -711,7 +711,7 @@ void sqrt3Vol(typename PFP::MAP& map, VertexAttribute& posit Traversor3WE tWE(map, dit); for(Dart ditWE = tWE.begin() ; ditWE != tWE.end() ; ditWE = tWE.next()) { - if(map.isBoundaryEdge(ditWE)) + if(map.isBoundaryEdge(ditWE) && !m.isMarked(ditWE)) m.markOrbit(ditWE); } @@ -728,45 +728,73 @@ void sqrt3Vol(typename PFP::MAP& map, VertexAttribute& posit } } - // - // edge-removal on all old boundary edges - // - TraversorE tE(map,selected); - for(Dart dit = tE.begin() ; dit != tE.end() ; dit = tE.next()) - { - if(m.isMarked(dit)) - { - m.unmarkOrbit(dit); - Dart d = map.phi2(map.phi3(map.findBoundaryFaceOfEdge(dit))); - Volume::Modelisation::Tetrahedralization::swapGen3To2(map, d); - - } - } - TraversorV tVg(map,selected); for(Dart dit = tVg.begin() ; dit != tVg.end() ; dit = tVg.next()) { if(map.isBoundaryVertex(dit) && !newBoundaryV.isMarked(dit)) { - typename PFP::VEC3 P = position[dit] ; + Dart db = map.findBoundaryFaceOfVertex(dit); + + typename PFP::VEC3 P = position[db] ; typename PFP::VEC3 newP(0) ; unsigned int val = 0 ; - Dart vit = dit ; + Dart vit = db ; do { newP += position[map.phi_1(map.phi2(map.phi1(vit)))] ; ++val ; vit = map.phi2(map.phi_1(vit)) ; - } while(vit != dit) ; + } while(vit != db) ; typename PFP::REAL K = sqrt3_K(val) ; newP *= typename PFP::REAL(3) ; newP -= typename PFP::REAL(val) * P ; newP *= K / typename PFP::REAL(2 * val) ; newP += (typename PFP::REAL(1) - K) * P ; - position[dit] = newP ; + position[db] = newP ; + } + } + + // + // edge-removal on all old boundary edges + // + TraversorE tE(map,selected); + for(Dart dit = tE.begin() ; dit != tE.end() ; dit = tE.next()) + { + if(m.isMarked(dit)) + { + m.unmarkOrbit(dit); + Dart d = map.phi2(map.phi3(map.findBoundaryFaceOfEdge(dit))); + Volume::Modelisation::Tetrahedralization::swapGen3To2(map, d); + } } +// TraversorV tVg(map,selected); +// for(Dart dit = tVg.begin() ; dit != tVg.end() ; dit = tVg.next()) +// { +// if(map.isBoundaryVertex(dit) && !newBoundaryV.isMarked(dit)) +// { +// Dart db = map.findBoundaryFaceOfVertex(dit); +// +// typename PFP::VEC3 P = position[db] ; +// typename PFP::VEC3 newP(0) ; +// unsigned int val = 0 ; +// Dart vit = db ; +// do +// { +// newP += position[map.phi_1(map.phi2(map.phi1(vit)))] ; +// ++val ; +// vit = map.phi2(map.phi_1(vit)) ; +// } while(vit != db) ; +// typename PFP::REAL K = sqrt3_K(val) ; +// newP *= typename PFP::REAL(3) ; +// newP -= typename PFP::REAL(val) * P ; +// newP *= K / typename PFP::REAL(2 * val) ; +// newP += (typename PFP::REAL(1) - K) * P ; +// position[db] = newP ; +// } +// } + //AutoVertexAttribute laplacian qui est une copie de position // TraversorV tVg2(map,selected); @@ -781,108 +809,48 @@ void sqrt3Vol(typename PFP::MAP& map, VertexAttribute& posit //echange lapaclian et position - VertexAutoAttribute diffCoord(map, "diffCoord"); - Algo::Volume::Geometry::computeLaplacianTopoVertices(map, position, diffCoord) ; - - VertexAutoAttribute vIndex(map, "vIndex"); - - unsigned int nb_vertices = map.template computeIndexCells(vIndex); - - - CellMarker lockingMarker(map); - - TraversorV tv(map); - for(Dart dit = tv.begin() ; dit != tv.end() ; dit = tv.next()) - { - if(!lockingMarker.isMarked(dit) && map.isBoundaryVertex(dit)) - lockingMarker.mark(dit); - } - - - NLContext nlContext = nlNewContext(); - nlSolverParameteri(NL_NB_VARIABLES, nb_vertices); - nlSolverParameteri(NL_LEAST_SQUARES, NL_TRUE); - nlSolverParameteri(NL_SOLVER, NL_CHOLMOD_EXT); - - - nlMakeCurrent(nlContext); - if(nlGetCurrentState() == NL_STATE_INITIAL) - nlBegin(NL_SYSTEM) ; - - for(int coord = 0; coord < 3; ++coord) - { - LinearSolving::setupVariables(map, vIndex, lockingMarker, position, coord) ; - nlBegin(NL_MATRIX) ; - LinearSolving::addRowsRHS_Laplacian_Topo(map, vIndex, diffCoord, coord) ; -// LinearSolving::addRowsRHS_Laplacian_Cotan(*map, perMap->vIndex, perMap->edgeWeight, perMap->vertexArea, perMap->diffCoord, coord) ; - nlEnd(NL_MATRIX) ; - nlEnd(NL_SYSTEM) ; - nlSolve() ; - LinearSolving::getResult(map, vIndex, position, coord) ; - nlReset(NL_TRUE) ; - } - - nlDeleteContext(nlContext); - - +// VertexAutoAttribute diffCoord(map, "diffCoord"); +// Algo::Volume::Geometry::computeLaplacianTopoVertices(map, position, diffCoord) ; // -// float weight = 1.0; - -// LinearSolving::LinearSolver ls(nbV); -// ls.set_least_squares(true); - -// for(unsigned int coord = 0 ; coord < 3 ; ++coord) +// VertexAutoAttribute vIndex(map, "vIndex"); +// +// unsigned int nb_vertices = map.template computeIndexCells(vIndex); +// +// +// CellMarker lockingMarker(map); +// +// TraversorV tv(map); +// for(Dart dit = tv.begin() ; dit != tv.end() ; dit = tv.next()) // { -// std::cout << "coord " << coord << std::flush; - -// TraversorV tv(map); -// for(Dart dit = tv.begin() ; dit != tv.end() ; dit = tv.next()) -// { -// ls.variable(indexV[dit]).set_value(position[dit][coord]); - -// if(map.isBoundaryVertex(dit)) -// ls.variable(indexV[dit]).lock(); -// } -// std::cout << "... variables set... " << std::flush; - -// ls.begin_system(); - -// TraversorV tv2(map); -// for(Dart dit = tv2.begin() ; dit != tv2.end() ; dit = tv2.next()) -// { -// if(!map.isBoundaryVertex(dit)) -// { -// ls.begin_row(); - -// float sum = 0; -// Traversor3VVaE tvvae(map, dit); -// for(Dart ditvvae = tvvae.begin() ; ditvvae != tvvae.end() ; ditvvae = tvvae.next()) -// { -// ls.add_coefficient(indexV[ditvvae],weight); -// sum += weight; -// } - -// ls.add_coefficient(indexV[dit],-sum); -// ls.normalize_row(); - -// ls.end_row(); -// } -// } - -// ls.end_system(); -// std::cout << "... system built... " << std::flush; -// ls.solve(); -// std::cout << "... system solved... " << std::flush; - -// TraversorV tv3(map); -// for(Dart dit = tv3.begin() ; dit != tv3.end() ; dit = tv3.next()) -// { -// position[dit][coord] = ls.variable(indexV[dit]).value(); -// } - -// ls.reset(false); -// std::cout << "... done" << std::endl; +// if(map.isBoundaryVertex(dit)) +// lockingMarker.mark(dit); +// } +// +// +// NLContext nlContext = nlNewContext(); +// nlSolverParameteri(NL_NB_VARIABLES, nb_vertices); +// nlSolverParameteri(NL_LEAST_SQUARES, NL_TRUE); +// nlSolverParameteri(NL_SOLVER, NL_CHOLMOD_EXT); +// +// +// nlMakeCurrent(nlContext); +// if(nlGetCurrentState() == NL_STATE_INITIAL) +// nlBegin(NL_SYSTEM) ; +// +// for(int coord = 0; coord < 3; ++coord) +// { +// LinearSolving::setupVariables(map, vIndex, lockingMarker, position, coord) ; +// nlBegin(NL_MATRIX) ; +// LinearSolving::addRowsRHS_Laplacian_Topo(map, vIndex, diffCoord, coord) ; +//// LinearSolving::addRowsRHS_Laplacian_Cotan(*map, perMap->vIndex, perMap->edgeWeight, perMap->vertexArea, perMap->diffCoord, coord) ; +// nlEnd(NL_MATRIX) ; +// nlEnd(NL_SYSTEM) ; +// nlSolve() ; +// LinearSolving::getResult(map, vIndex, position, coord) ; +// nlReset(NL_TRUE) ; // } +// +// nlDeleteContext(nlContext); } template diff --git a/include/Algo/Multiresolution/Map2MR/Filters/sqrt3.h b/include/Algo/Multiresolution/Map2MR/Filters/sqrt3.h index bbcbdf8ca990fac4db1affeb9f90fa28f2821a78..b79d9254eb1e7b64b40c9402453bc6d5b3b41df3 100644 --- a/include/Algo/Multiresolution/Map2MR/Filters/sqrt3.h +++ b/include/Algo/Multiresolution/Map2MR/Filters/sqrt3.h @@ -114,6 +114,11 @@ public: {} void operator() () + { + + } + + void operator() (bool filtering) { TraversorF trav(m_map) ; for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) @@ -145,10 +150,6 @@ public: } } - void operator() (bool filtering) - { - - } } ; template @@ -163,6 +164,11 @@ public: {} void operator() () + { + + } + + void operator() (bool filtering) { TraversorV trav(m_map) ; for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) @@ -234,10 +240,6 @@ public: } } - void operator() (bool filtering) - { - - } } ; template @@ -252,6 +254,11 @@ public: {} void operator() () + { + + } + + void operator() (bool filtering) { TraversorV trav(m_map) ; for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) @@ -289,10 +296,6 @@ public: } } - void operator() (bool filtering) - { - - } } ; template @@ -307,6 +310,11 @@ public: {} void operator() () + { + + } + + void operator() (bool filtering) { TraversorF trav(m_map) ; for (Dart d = trav.begin(); d != trav.end(); d = trav.next()) @@ -331,10 +339,6 @@ public: } } - void operator() (bool filtering) - { - - } } ; /*********************************************************************************