diff --git a/Geo/Curvature.cpp b/Geo/Curvature.cpp index 4b33765d28e0ea7bb6cc80209e1c48064d0b897d..62a88f310fdf3fa21b9129c6e2706a3ccce19eb0 100644 --- a/Geo/Curvature.cpp +++ b/Geo/Curvature.cpp @@ -137,7 +137,7 @@ void Curvature::initializeMap() // face is a pointer to one surface of the group "FinalEntityList" GFace* face = _ptFinalEntityList[i]; - std::cout << "Face " << i << " has " << face->getNumMeshElements() << " elements" << std::endl; +// std::cout << "Face " << i << " has " << face->getNumMeshElements() << " elements" << std::endl; // Loop over the element all the element of the "myTag"-surface for (int iElem = 0; iElem < face->getNumMeshElements(); iElem++) @@ -488,7 +488,7 @@ void Curvature::computeRusinkiewiczNormals() void Curvature::computePointareas() { - std::cout << "Computing point areas... " << std::endl; +// std::cout << "Computing point areas... " << std::endl; // std::cout << "The mesh has " << _VertexToInt.size() << " nodes" << std::endl; SVector3 e[3]; @@ -565,12 +565,13 @@ void Curvature::computePointareas() } //End of loop over iElem - std::cout << "Done computing pointareas" << std::endl; // std::cout << "_pointareas.size = " << _pointareas.size() << std::endl; // std::cout << "_cornerareas.size = " << _cornerareas.size() << std::endl; } //End of loop over _ptFinalEntityList + //std::cout << "Done computing pointareas" << std::endl; + } //End of the method "computePointareas" @@ -918,6 +919,61 @@ void Curvature::computeCurvature_Rusinkiewicz() } + //======================================================================================================== + + void Curvature::elementNodalAbsoluteValues(MTriangle* triangle, double& c0, double& c1, double& c2) + { + MVertex* A = triangle->getVertex(0); + MVertex* B = triangle->getVertex(1); + MVertex* C = triangle->getVertex(2); + + int V0 = 0; + int V1 = 0; + int V2 = 0; + + std::map<int,int>::iterator vertexIterator; + + vertexIterator = _VertexToInt.find( A->getNum() ); + if ( vertexIterator != _VertexToInt.end() ) + { + V0 = (*vertexIterator).second; + } + else + { + std::cout << "Didn't find vertex with number " << A->getNum() << " in _VertextToInt !" << std::endl; + } + + vertexIterator = _VertexToInt.find( B->getNum() ); + if ( vertexIterator != _VertexToInt.end() ) + { + V1 = (*vertexIterator).second; + } + else + { + std::cout << "Didn't find vertex with number " << B->getNum() << " in _VertextToInt !" << std::endl; + } + + vertexIterator = _VertexToInt.find( C->getNum() ); + if ( vertexIterator != _VertexToInt.end() ) + { + V2 = (*vertexIterator).second; + } + else + { + std::cout << "Didn't find vertex with number " << C->getNum() << " in _VertextToInt !" << std::endl; + } + + + // const int V0 = _VertexToInt[A->getNum()]; + // const int V1 = _VertexToInt[B->getNum()]; + // const int V2 = _VertexToInt[C->getNum()]; + + c0 = std::abs(_VertexCurve[V0]); //Mean curvature in vertex 0 + c1 = std::abs(_VertexCurve[V1]); //Mean curvature in vertex 1 + c2 = std::abs(_VertexCurve[V2]); //Mean curvature in vertex 2 + + } + //======================================================================================================== void Curvature::writeToPosFile( const std::string & filename) diff --git a/Geo/Curvature.h b/Geo/Curvature.h index 9ed41086777e92b0bc6e83d97254a207e97d73ba..3134442043fffa6b6d6d6a4887e2b45f63e2858a 100644 --- a/Geo/Curvature.h +++ b/Geo/Curvature.h @@ -189,6 +189,8 @@ public: void elementNodalValues(MTriangle* triangle, double& c0, double& c1, double& c2); + void elementNodalAbsoluteValues(MTriangle* triangle, double& c0, double& c1, double& c2); + void writeToPosFile( const std::string & filename); void writeToVtkFile( const std::string & filename); diff --git a/Geo/GFaceCompound.cpp b/Geo/GFaceCompound.cpp index 349b01c9e5c560d447f5ec264532a5e4c3f37abc..b7c31552fa488724d38fe3b21c938cae42a55674 100644 --- a/Geo/GFaceCompound.cpp +++ b/Geo/GFaceCompound.cpp @@ -1373,7 +1373,7 @@ double GFaceCompound::curvatureMax(const SPoint2 ¶m) const if( !Curvature::valueAlreadyComputed() ) { - std::cout << "Need to compute curvature" << std::endl; + std::cout << "Need to compute discrete curvature" << std::endl; std::cout << "Getting instance of curvature" << std::endl; curvature.setGModel( model() ); @@ -1381,7 +1381,7 @@ double GFaceCompound::curvatureMax(const SPoint2 ¶m) const curvature.writeToPosFile("curvature.pos"); curvature.writeToVtkFile("curvature.vtk"); - std::cout << " ... finished" << std::endl; + std::cout << " ... computing curvature finished" << std::endl; } @@ -1414,15 +1414,29 @@ double GFaceCompound::curvatureMax(const SPoint2 ¶m) const } else if (lt->gf->geomType() == GEntity::DiscreteSurface) { - //std::cout << "I'm in DiscreteSurface" << std::endl; +// std::cout << "I'm in DiscreteSurface" << std::endl; double c0; double c1; double c2; - curvature.elementNodalValues(lt->tri,c0, c1, c2); + //curvature.elementNodalValues(lt->tri,c0, c1, c2); + curvature.elementNodalAbsoluteValues(lt->tri,c0, c1, c2); + double cv = (1-U-V)*c0 + U*c1 + V*c2; - //std::cin.get(); - //std::cout << "(" << c0 << "," << c1 << "," << c2 << ")" << std::endl; - //std::cout << "The curvature of the triangle " << lt->tri->getNum() << " is " << cv << std::endl; +// std::cout << "(" << c0 << "," << c1 << "," << c2 << ")" << std::endl; +// MVertex* V0 = lt->tri->getVertex(0); +// MVertex* V1 = lt->tri->getVertex(1); +// MVertex* V2 = lt->tri->getVertex(2); + +// std::cout << "=====================================================" << std::endl; +// std::cout << "Parametric coordinates: " << param.x() << "," << param.y() << std::endl; +// std::cout << "The coordinates of the triangle are:" << std::endl; +// std::cout << "\t[" << V0->x() << "," << V0->y() << "," << V0->z() << "]" << std::endl; +// std::cout << "\t[" << V1->x() << "," << V1->y() << "," << V1->z() << "]" << std::endl; +// std::cout << "\t[" << V2->x() << "," << V2->y() << "," << V2->z() << "]" << std::endl; +// std::cout << "The curvature of the triangle " << lt->tri->getNum() << " is " << cv << std::endl; +// std::cout << std::endl; +// std::cin.get(); +// return 1.0; return cv; // double curv= 0.; @@ -1430,7 +1444,6 @@ double GFaceCompound::curvatureMax(const SPoint2 ¶m) const // return curv; } - std::cin.get(); return 0.; }